% define desired region
region_name='global';
nx=2160;

% get and save daily-averaged regional fields
pn='~dmenemen/llc_2160/MITgcm/';
pout=['~dmenemen/llc_2160/regions/' region_name '/'];
for fnm={'U','V'}
    eval(['mkdir ' pout fnm{1} '3d_daily'])
   
    eval(['cd ' pout fnm{1} '_daily'])
    
   
    for  ts=92160+90*80*24:(80*24):92160+90*80*24%1586400;
         
         if (ts<1198080)
            pin=[pn 'run_day49_624/'];
         else
            pin=[pn 'run/'];
         end
        
        dy=ts2dte(ts,45,2011,1,17,30);
             
        fout=[fnm{1} '_100mAVG_' int2str((nx)) 'x' int2str((nx*13)) '.' dy];

        if mod(ts,144*24*10) 
            disp(fout);
        end
        
        fld=0; n=0; fld2 = 0; pwr = 0; fld_mn = 0; fld_mn2 = 0; n2 = 0;
        
        for h=0:23 %, mydisp(h)
            fin=[pin fnm{1} '.' myint2str(ts+h*80,10) '.data'];
            
            if exist(fin)
                for k = 1 : 23 % up to 100m, 23rd layer
                    
                    skip=((k)-1)*nx*nx*13;
                    tmp = readbin(fin,[nx, nx*13],1,'real*4',skip);
                    if sum(sum(tmp)) == 0
                        warning(['ZEROs FILE ' fin])
                        break
                    else
                        fld=fld+tmp;
                        %fld_mn = fld_mn + tmp.^2;
                        n=n+1
                    end
                end
            else
                warning(['MISSING FILE ' fin])
                break
            end
            
         

        end
        if n==24
            
            writebin([pout fnm{1} '_100mAVG_daily/' fout],fld./n);
%             fout_msq = [pout fnm{1} '_mnsq_daily/' fnm{1} '_mnsq_' int2str((nx)) 'x' int2str((nx*13)) '.' dy];
%             writebin(fout_msq,fld_mn./n);
            
        end
      
        
    end
end