% timatch_30
% Xcruise 2009 :: 2009-09-07 :: LBB
% Aggregates/averages data in 30-minute windows.
disp(['timatch_30: ', cruise]);

fclose('all');

clear ztsg zt zsnk zo y x w vv vrsp uu usr un5 un1 un ue5 ue1 ue ts_tsg ts_snk ts tn tm_s tm_ot tm_n tm_in tcol taub ta_im ta_etl ta t sogy sogx
clear sogm5 sog sody5 sody sodx5 sodx sm skip_it sig_sp shp_spd shp_hed sal_tsg s3 s1 rh_im rh_etl rf relsp reldir rdcon qs qcol qa_im qa_etl qa psp_im psp1 psp prt_it prt_f pir1 pir org odec5 ode npz
clear np n minn m lrgsm lrgs5 lrgs lrgm5 lrgcm lrgc5 lrgc lrg lon lm lat k jl jk jj jg jf jdx jd jazz jazy5 jazy1 jazy jazx5 jazx1 jazx jaz5 jaz1 jaz jax jaq jadson jadscs jadmns jad5 
clear jad0 j imuy imux imum5 imu imdm5 imd ii i hsb hr hnet hm hlb g frac flist3 flist2 flist f e dt dir1 d cogm5 cog c by bx ay ax ans V U Sm Rns
clear Rnl Lm Dt

% inherit saveit from calling flux_eval program
%saveit=true;    %set true to one to print data to file
clearit=false;   %set true to clear variables at end

% jazx=zeros(6,1440); % never used
% jaz=zeros(16,1440);
% jazy=zeros(12,1440);

jd=sprintf('%03i',ddd);
f=[way_proc_data_flux cruise '_proc_son_1min_' sprintf('%04i_%02i_%02i_%03i',Vdate(1),Vdate(2),Vdate(3),ddd) '.txt'];
jazx=load(f)';

%if scs file is present or made, then use it. Otherwise will set matrix to
%NaN and use COG, SOG and heading from PSD GPS system.
if exist([way_proc_data_flux cruise '_proc_scs_1min_' sprintf('%04i_%02i_%02i_%03i',Vdate(1),Vdate(2),Vdate(3),ddd)  '.txt'],'file')
    e=[way_proc_data_flux cruise '_proc_scs_1min_' sprintf('%04i_%02i_%02i_%03i',Vdate(1),Vdate(2),Vdate(3),ddd)  '.txt'];
    jaz=load(e)';
    latm=jaz(2,:);
    lonm=jaz(3,:);
% % % % jaz(1,1:npz)=jd_scs;	%pcode gps time,
% % % % jaz(2,1:npz)=latm(1:npz)';	%decimal latitude, deg
% % % % jaz(3,1:npz)=lonm(1:npz)';	%longitude, deg
% % % % jaz(4,1:npz)=cogm(1:npz)';	% pcode cog, deg
% % % % jaz(5,1:npz)=sogm(1:npz)';	% pcode sog, m/s
% % % % jaz(6,1:npz)=tsgm(1:npz)';	%thermosalinograph T, C
% % % % jaz(7,1:npz)=tssm(1:npz)';	%salinity, psu
% % % % jaz(8,1:npz)=rhm(1:npz)';	%rel hum,%
% % % % jaz(9,1:npz)=tam(1:npz)';	%air t, C
% % % % jaz(10,1:npz)=solarm(1:npz)';	%solar flux, w/m^2
% % % % jaz(11,1:npz)=lrgm(1:npz);	%laser gyro, deg
% % % % jaz(12,1:npz)=imum(1:npz)';	%imet true wind speed, m/s
% % % % jaz(13,1:npz)=imdm(1:npz)';	%imet true wind dir, deg
% % % % jaz(14,1:npz)=rainm(1:npz)'; %rain, mm 
% % % % jaz(15,1:npz)=odecm(1:npz)'; %odec 
% % % % jaz(16,1:npz)=irm(1:npz)'; %imet ir
else   %if no SCS data file, use COG, SOG and heading from PSD GPS system
    jaz=ones(16,1440)*NaN;
    e=load([way_proc_data_flux cruise '_proc_gpsnav_1min_' sprintf('%04i_%02i_%02i_%03i',Vdate(1),Vdate(2),Vdate(3),ddd)  '.txt'])';
    jaz(1,:)=e(1,:);
    jaz(2,:)=e(2,:);
    jaz(3,:)=e(3,:);
    jaz(4,:)=e(5,:);
    jaz(5,:)=e(4,:);
    jaz(11,:)=e(6,:);
    
    latm=jaz(2,:);
    lonm=jaz(3,:);
end;


g=[way_proc_data_flux cruise '_proc_met_1min_' sprintf('%04i_%02i_%02i_%03i',Vdate(1),Vdate(2),Vdate(3),ddd) '.txt'];
jazy=load(g)';

zsnk=ones(1,48);
ztsg=ones(1,48);

[ax,bx]=size(jazx);			%jazx
[c,d]=size(jaz);			%jaz
[ay,by]=size(jazy);
tn=60*24;														%number of minutes in a day
n=1:tn;
jad0=str2double(jd)+(n-1)/24/60;jad0(tn)=jad0(tn)+1e-5;		%	don't ask
jadson=jazx(1,1:bx);
jadscs=jaz(1,1:d);
jadmns=jazy(1,1:by);

for ii=2:length(jadscs);
   if jadscs(ii)<0,jadscs(ii)=jadscs(ii-1)+.000001;end;
end;

%%%%%%%%**********   sonic loop
clear jazx1;
jazxp=interp1(jadson,jazx',jad0);
jazx1=jazxp';

df=false;
if df
    %fill in beginning empty spots
    jk=1;
    while jad0(jk)<jadson(1),
        jazx1(:,jk)=jazx(:,1);
        jk=jk+1;
    end;

    for i=jk:tn;
        t=find(jadson>jad0(i));
        if ~isempty(t),
            j=t(1);
        else
            j=i;
        end;
        dt=jadson(j)-jadson(j-1);
        Dt=jad0(i)-jadson(j-1);
        jazx1(:,i)=jazx(:,j-1)*(1-Dt/dt)+jazx(:,j)*Dt/dt;
        if j==bx;
            break;
        end;
    end;

    jk=tn;
    while jad0(jk)>=jadson(bx),
        jazx1(:,jk)=jazx(:,bx);
        jk=jk-1;
    end;
end; % if df --not done

% means loop
clear jazy1;
%fill in beginning empty spots 
jk=1;
while jad0(jk)<jadmns(1),
   jazy1(:,jk)=jazy(:,1);
   jk=jk+1;
end;

for i=jk:tn;
   t=find(jadmns>jad0(i));
   if length(t)>0,
      j=t(1);
   else
      j=2;
   end;
   dt=jadmns(j)-jadmns(j-1);Dt=jad0(i)-jadmns(j-1);
   jazy1(:,i)=jazy(:,j-1)*(1-Dt/dt)+jazy(:,j)*Dt/dt;
   if j==by;break;end;
end;
jk=tn;
while jad0(jk)>=jadmns(by),
   jazy1(:,jk)=jazy(:,by);
   jk=jk-1;
end;
%jazy1=interp1(jadmns,jazy',jad0)';

%%%%%%%%%%%%%****************   ship data loop
%first setup to interpolate fricking angles
cog=jaz(4,:);
sog=jaz(5,:);
lrg=jaz(11,:);
imu=jaz(12,:);
imd=jaz(13,:);

iic=find(sog<0 | sog>10);cog(iic)=NaN;sog(iic)=NaN;lrg(iic)=NaN;
iid=find(imu<0 | imu>50);imu(iid)=NaN;imd(iid)=NaN;

rdcon=pi/180;
jaz(4,:)=sog.*cos(cog*rdcon);			%sogx
jaz(5,:)=sog.*sin(cog*rdcon);			%sogy
jaz(12,:)=imu.*cos(imd*rdcon);		%imux
jaz(13,:)=imu.*sin(imd*rdcon);		%imuy

lrgc=cos(lrg*rdcon);
lrgs=sin(lrg*rdcon);

clear jaz1;
												%fill in beginning empty spots 
jk=1;
while jad0(jk)<jadscs(1),
   jaz1(:,jk)=jaz(:,1);
   lrgcm(jk)=lrgc(1);
   lrgsm(jk)=lrgs(1);
   jk=jk+1;
end;

for i=jk:tn;
   t=find(jadscs>jad0(i));
   if~isempty(t)
      j=t(1);
      dt=jadscs(j)-jadscs(j-1);Dt=jad0(i)-jadscs(j-1);
      jaz1(:,i)=jaz(:,j-1)*(1-Dt/dt)+jaz(:,j)*Dt/dt;
      lrgcm(i)=lrgc(j-1)*(1-Dt/dt)+lrgc(j)*Dt/dt;
      lrgsm(i)=lrgs(j-1)*(1-Dt/dt)+lrgs(j)*Dt/dt;
   else
      break;
   end;
end;

jk=tn;
while jk>0 & jad0(jk)>=jadscs(d),
   jaz1(:,jk)=jaz(:,d);
   lrgcm(jk)=lrgc(d);
   lrgsm(jk)=lrgs(d);
   jk=jk-1;
end;

% do 30 minute averages
jax=1440;
jk=48;

ode=jaz1(15,:);%if no ODEC, use SOG or NaN

sodx=ode.*lrgcm;
sody=ode.*lrgsm;
uu=jazx1(2,:)/.95;				%.95 is empirical distortion factor
vv=jazx1(3,:)/1.15;				%1.15 is other empirical distortion factor

un=uu.*lrgcm+vv.*lrgsm;
ue=uu.*lrgsm-vv.*lrgcm;

clear jaz5 jazy5 jazx5 lrgc5 lrgs5 sodx5 sody5 un5 ue5 vrsp jad5;

for j=1:jk
   jl=30*(j-1)+1;
   jj=min(jax,30*j);
   jazx5(:,j)=mean(jazx1(:,jl:jj)')';
   jazy5(:,j)=mean(jazy1(:,jl:jj)')';
   jaz5(:,j)=mean(jaz1(:,jl:jj)')';
   lrgc5(j)=mean(lrgcm(jl:jj)')';
   lrgs5(j)=mean(lrgsm(jl:jj)')';
   sodx5(j)=mean(sodx(jl:jj)')';
   sody5(j)=mean(sody(jl:jj)')';
   un5(j)=mean(un(jl:jj)')';
   ue5(j)=mean(ue(jl:jj)')';
   vrsp(j)=var(jaz1(5,jl:jj))+var(jaz1(6,jl:jj));
   jad5(j)=jad0(jl);
end;
   
% then uncomponentize angles
U=jazx5(2,:);										
V=jazx5(3,:);										
W=jazx5(4,:);
reldir=mod(atan2(-V,U)/rdcon+180+360,360);		%signs convert jeff coord to conventional
relsp=sqrt(U.*U+V.*V);
sogx=jaz5(4,:);
sogy=jaz5(5,:);
imux=jaz5(12,:);%imet true wind
imuy=jaz5(13,:);%imet true wind

sogm5=sqrt(sogx.*sogx+sogy.*sogy);					%1-min sog
cogm5=mod(atan2(sogy,sogx)/rdcon+360,360);
imum5=sqrt(imux.*imux+imuy.*imuy);					%1-min imet speed, m/s
imdm5=mod(atan2(imuy,imux)/rdcon+360,360);
lrgm5=mod(atan2(lrgs5,lrgc5)/rdcon+360,360);		%laser ring gyro
odec5=sqrt(sodx5.*sodx5+sody5.*sody5);

% do some true wind calcs
un1=un5+sogx;
ue1=ue5+sogy;
s1=sqrt(un1.*un1+ue1.*ue1);
s3=sqrt(un1.*un1+ue1.*ue1+W.*W);
s_shp=sqrt(imux.*imux+imuy.*imuy);

dir1=mod(atan2(-ue1,-un1)*180/pi+360,360);%- sign converts to to from
dir_shp=mod(atan2(imuy,imux)*180/pi+360,360);

% sea temps
ts_snk=jazy5(6,:)+tsnake_adjustment; %correct tsnake according to tsnake_adjustment set in read_parameters
ts_tsg=jaz5(6,:);
sal_tsg=jaz5(7,:);
ts_tsg(ts_tsg<0)=NaN;
ts_snk(ts_snk<0)=NaN;
ts=ts_snk;					%water temperature if snake is in
% ts(isnan(ts_snk))=ts_tsg(isnan(ts_snk)); % use thermosalinograph intake T if snake is missing
zt=zsnk*.05;				%depth of water temp measurement
% zt(isnan(ts_snk))=5.0; 
qs=qsea(ts);

% air temps
ta_etl=jazy5(7,:)+ta_etl_adjustment; % correction for etl from read_parameters
ta_im=jaz5(9,:)+ta_im_adjustment;
ta_etl(ta_etl>50 | ta_etl<0)=NaN; % reject out of reasonable range
ta_im(ta_im>50 | ta_im<0)=NaN;

% pressure
pressm=jazy5(12,:);

% relative humidity
rh_etl=jazy5(8,:);  %adjustment ???
rh_im=jaz5(8,:); %9-->8 18 Jan 2008 SPdeS
rh_im(rh_im<0)=NaN;  % Also reject RH>100% ??? SPdeS
rh_etl(rh_etl<0)=NaN;

qa_etl=qair3([pressm' ta_etl' min(rh_etl',100)])'; % min(  ,100) rejects supersaturation
qa_im=qair3([pressm' ta_im' min(rh_im',100)])';

% solar fluxes
psp1=jazy5(3,:);
psp_im=jaz5(10,:);
psp=psp1;
ii=find(psp<-50);psp(ii)=psp_im(ii);

% IR fluxes
pir=jazy5(2,:);
pir_im=jaz5(16,:);%imet ir added 12/04

org=jazy5(9,:);
lat=jaz5(2,:);
lon=jaz5(3,:);

%shp_spd=odec5;
shp_spd=sogm5;
shp_hed=lrgm5;
sig_sp=sqrt(vrsp);

% % % Choose Ta and RH for flux calculation.
% % % VOCALS 2008: It looks like there are problems with PSD RH, if not Tair
% % % also. WHOI had recently calibrated the ship's IMET Rotronics TRH sensor.
% % % As of end of VOCALS RHB leg 1, use ship's IMET Ta and RH.
% % % Awaiting word from postcalibration comparison. SPdeS
% % ta=ta_im; % VOCALS leg 1 PSD T was 0.5 C warmer; don't know which is better. Use IMET for consistency with RH.
% % qa=qa_im; % VOCALS leg 1 PSD RH was 10% low and responded slow, use ship's IMET recently calibrated by WHOI.

ta=ta_etl; 
qa=qa_etl;

for k=1:jk
   %   u    0  ts  ta      qs     qa   rs      rl   rain     zi  p   zu   zt   zq
   x=[s1(k) 0 ts(k) ta(k) qs(k) qa(k) psp(k) pir(k) org(k) 600 press(k) zus zts zqs lat(k) 1 0 1 1];
   y=cor30b(x);
   hsb(k)=y(1);
   hlb(k)=y(2);
   taub(k)=y(3);
   tcol(k)=y(11);
   qcol(k)=y(12);
   rf(k)=y(14);
   usr(k)=y(8);
   zo(k)=y(4);
end;

Rns=psp*.955;
Rnl=.97*(pir-5.67e-8*(ts_snk+273.15).^4);
hnet=Rns+Rnl-hsb-hlb-rf;
hm=nanmean(hnet);sm=nanmean(hsb);lm=nanmean(hlb);Sm=nanmean(Rns);Lm=nanmean(Rnl);

if saveit;
    np=jk;
    eval([ 'prt_flux_30_', currentMname])
end;

if clearit
    clear ztsg zt zsnk zo y x w vv vrsp uu usr un5 un1 un ue5 ue1 ue ts_tsg ts_snk ts tn tm_s tm_ot tm_n tm_in tcol taub ta_im ta_etl ta t sogy sogx
    clear sogm5 sog sody5 sody sodx5 sodx sm skip_it sig_sp shp_spd shp_hed sal_tsg s3 s1 rh_im rh_etl rf relsp reldir rdcon qs qcol qa_im qa_etl qa psp_im psp1 psp prt_it prt_f pir1 pir org odec5 ode npz
    clear np n minn m lrgsm lrgs5 lrgs lrgm5 lrgcm lrgc5 lrgc lrg lon lm lat k jl jk jj jg jf jdx jd jazz jazy5 jazy1 jazy jazx5 jazx1 jazx jaz5 jaz1 jaz jax jaq jadson jadscs jadmns jad5
    clear jad0 j imuy imux imum5 imu imdm5 imd ii i hsb hr hnet hm hlb g frac flist3 flist2 flist f e dt dir1 d cogm5 cog c by bx ay ax ans V U Sm Rns
    clear Rnl Lm Dt
end;
