%make a movie of surface Fe in the southern ocean

%define input
 idir='../dat/mod/daily/';
 fils = dir([idir '*.nc']);

%plotting details
 xlims=[0 360];
 ylims=[-78.15 -20];
 in=[0:0.025:0.5 0.6:0.1:1];

%get grid once
 x=ncread([idir fils(1).name],'xt');
 y=ncread([idir fils(1).name],'yt');
 [X,Y]=meshgrid(x,y); 

%some ini
 n=0; m=0; l=1; tmp=0;
 mbin = linspace(0,365,13);
 SFE=zeros(length(x),length(y),12);

for i=1:length(fils)
 ssf=ncread([idir fils(i).name],'ssf'); 
 t=ncread([idir fils(i).name],'time');
 I=find(ssf<-1e10); ssf(I)=NaN; 
 ssf=double(ssf).*1e9;
 for k=1:length(t)
    n=n+1; m=m+1;
    if n==366; break; end
    tmp = tmp+squee(ssf(:,:,k));
    if n>=mbin(l+1)
     SFE(:,:,l)=tmp./m;
     m=0; l=l+1; tmp=0;
    end
 end
 if n==366; break; end
end

unit='mumolFe/m^3';
save ../dat/mod/processed/feMonth_1918.mat x y SFE unit

