disp('fastexrt.m')
%clear;
%openmotion;opengyro;openship;opensons;openophir
% 

if hhh<10,
   dfx2=['0',num2str(hhh)];
else,
   dfx2=num2str(hhh);
end;

if ddd<10,
   dfx1=['00',num2str(ddd)];
elseif ddd<100,
   dfx1=['0',num2str(ddd)];
else,
   dfx1=num2str(ddd);
end;

dfx=[dfx1 dfx2];

%sfl=['/home/jeh/fastex/P_',dfx,'.dat'];
sfl=['c:\hare\fastex\P_',dfx,'.dat'];

tgyro(length(tgyro))=tmotion(length(tmotion));
qq=1;shead=(1:length(tmotion)).*0;chead=(1:length(tmotion)).*0;heading=[heading heading(length(heading))];
for pp=1:length(tmotion);
   if tmotion(pp)==tgyro(qq),
      shead(pp)=sin(heading(qq)./180.*pi);
      chead(pp)=cos(heading(qq)./180.*pi);
   elseif tmotion(pp) > tgyro(qq),
      shead(pp)=(sin(heading(qq+1)./180.*pi)-sin(heading(qq)./180.*pi))/(tgyro(qq+1)-tgyro(qq))*(tmotion(pp)-tgyro(qq))+sin(heading(qq)./180.*pi);
      chead(pp)=(cos(heading(qq+1)./180.*pi)-cos(heading(qq)./180.*pi))/(tgyro(qq+1)-tgyro(qq))*(tmotion(pp)-tgyro(qq))+cos(heading(qq)./180.*pi);
      if tmotion(pp) > tgyro(qq+1), qq=qq+1;end;
   end;
end;

clear heading;
heading=atan2(shead,chead).*180./pi;
clear tgyro qq pp;
for ppp=1:length(heading);
   if heading(ppp)<0,
      heading(ppp)=heading(ppp)+360;
   end;
end;
heading=heading';

G=sqrt(mean(accx)^2+mean(accy)^2+mean(accz)^2);
ratex=detrend(ratex);ratey=detrend(ratey);ratez=detrend(ratez);

ssx=(1:length(accx)).*0;ssx(1)=su(1);
ssy=(1:length(accx)).*0;ssy(1)=sv(1);
ssz=(1:length(accx)).*0;ssz(1)=sw(1);
sst=(1:length(accx)).*0;sst(1)=Tsonic(1);
jj=2;

for ii=2:length(accx);
   while tsons(jj)<tmotion(ii) & jj < length(su);
      jj=jj+1;
   end;
   ssx(ii)=(su(jj)-su(jj-1))/(tsons(jj)-tsons(jj-1))*(tmotion(ii)-tsons(jj-1))+su(jj-1);
   ssy(ii)=(sv(jj)-sv(jj-1))/(tsons(jj)-tsons(jj-1))*(tmotion(ii)-tsons(jj-1))+sv(jj-1);
   ssz(ii)=(sw(jj)-sw(jj-1))/(tsons(jj)-tsons(jj-1))*(tmotion(ii)-tsons(jj-1))+sw(jj-1);
   sst(ii)=(Tsonic(jj)-Tsonic(jj-1))/(tsons(jj)-tsons(jj-1))*(tmotion(ii)-tsons(jj-1))+Tsonic(jj-1);
end;

clear su sv sw headr tsons flist;
fclose('all');

maccx=mean(accx);
maccy=mean(accy);
maccz=mean(accz);

%disp('Complementary filtering')

%Fc=1/60;                 %60 seconds
Fc=1/300;                  %10 seconds

tau=1/(Fc*2*pi);bc=Fc/10;[bf af]=butter(1,bc);dt=1/20;
om_c=tan(dt/2/tau);tau_c=dt/2/om_c;
clear Fc;

thetarate=ratey*tau_c+accx/G;                              %pitchrate
theta=filter(bf,af,thetarate,maccx/G);

phirate=ratex*tau_c+(accy/G)./cos(theta);
phi=filter(bf,af,phirate,maccy/G);

heading=heading/180*pi;
psirate=ratez*tau_c+heading;                              %no turns thru zero allowed
psi=filter(bf,af,psirate,mean(heading));


%Compute Omega X R for sonics
Rz=1;Rx=0;Ry=0;
OXRx=-ratey*Rz+ratez*Ry;
OXRy=-ratex*Rz-ratez*Rx;
OXRz=+ratex*Ry+ratey*Rx;
OM=[OXRx OXRy OXRz];
%TOM=0.*OM';

clear ratex ratey ratez maccx maccy maccz thetarate phirate psirate OXRx;
clear OXRy OXRz Rx Ry Rz;

Plat=[accx'
      accy'
      accz'];

%clear accz accy;

Vel=[ssy              %Matrix defined as [V U W]
     ssx               %
     ssz];             %

Platacc=Plat;
Velrot=Vel;
Plat=Plat';
Vel=Vel';

%disp('Rotating coordinate system');

cos_theta=cos(theta);sin_theta=sin(theta);
cos_phi=cos(phi);sin_phi=sin(phi);
cos_psi=cos(psi);sin_psi=sin(psi);

clear theta;clear psi;clear phi;

%T11=cos_psi.*cos_theta;                           %JEH Masters
%T12=sin_psi.*cos_phi-cos_psi.*sin_theta.*sin_phi;
%T13=-sin_psi.*sin_phi-cos_psi.*sin_theta.*cos_phi;
%T21=-sin_psi.*cos_theta;
%T22=cos_psi.*cos_phi+sin_psi.*sin_theta.*sin_phi;
%T23=-cos_psi.*sin_phi+sin_psi.*sin_theta.*cos_phi;
%T31=sin_theta;
%T32=cos_theta.*sin_phi;
%T33=cos_theta.*cos_phi;

T11=cos_psi.*cos_theta;                           %JBE
T12=sin_psi.*cos_phi+cos_psi.*sin_theta.*sin_phi;
T13=-sin_psi.*sin_phi+cos_psi.*sin_theta.*cos_phi;
T21=-sin_psi.*cos_theta;
T22=cos_psi.*cos_phi-sin_psi.*sin_theta.*sin_phi;
T23=-cos_psi.*sin_phi-sin_psi.*sin_theta.*cos_phi;
T31=-sin_theta;
T32=cos_theta.*sin_phi;
T33=cos_theta.*cos_phi;

clear cos_theta;clear sin_theta;clear cos_phi;
clear sin_phi;clear cos_psi;clear sin_psi;
accx=accx';accy=accy';accz=accz';
%ssx=ssx';ssy=ssy';ssz=ssz';
T11=T11';T12=T12';T13=T13';
T21=T21';T22=T22';T23=T23';
T31=T31';T32=T32';T33=T33';
%for zz=1:length(accx);
%   T=[T11(zz) T12(zz) T13(zz)
%      T21(zz) T22(zz) T23(zz)
%      T31(zz) T32(zz) T33(zz)];
%   Platacc(:,zz)=T*Plat(zz,:)';     %Rotated accelerations
%    Velrot(:,zz)=T*Vel(zz,:)';      %Rotated sonic velocities
%    TOM(:,zz)=T*OM(zz,:)';          %Rotated  Om X R
   Platacc(1,:)=T11.*accx+T12.*accy+T13.*accz;
   Platacc(2,:)=T21.*accx+T22.*accy+T23.*accz;
   Platacc(3,:)=T31.*accx+T32.*accy+T33.*accz;
    Velrot(1,:)=T11.*ssx+T12.*ssy+T13.*ssz;
    Velrot(2,:)=T21.*ssx+T22.*ssy+T23.*ssz;
    Velrot(3,:)=T31.*ssx+T32.*ssy+T33.*ssz;
%end;

clear Plat Vel T OM ssx ssy ssz accy accz;

maccx=mean(Platacc(1,:));
maccy=mean(Platacc(2,:));
maccz=mean(Platacc(3,:));

Platacc(1,:)=Platacc(1,:)-maccx;
Platacc(2,:)=Platacc(2,:)-maccy;
Platacc(3,:)=Platacc(3,:)-maccz;

clear maccx maccy maccz;

Platvel=(1:length(accx)).*0;
Platvel=[Platvel
         Platvel
         Platvel];

Acclow=Platacc.*0;

%Fc=1/60;                 %60 seconds
Fc=1/150;                 %10 seconds
tau=1/(Fc*2*pi);bc=Fc/10;[bf af]=butter(1,bc);dt=1/20;
om_c=tan(dt/2/tau);tau_c=dt/2/om_c;
clear Fc;

Acclow(1,:)=filtfilt(bf,af,Platacc(1,:));
Acclow(2,:)=filtfilt(bf,af,Platacc(2,:));
Acclow(3,:)=filtfilt(bf,af,Platacc(3,:));

Platform=Platacc.*0;

for zz=1:length(accx);
   for kk=1:3;
      Platform(kk,zz)=Platacc(kk,zz)-Acclow(kk,zz);      %Highpassed Rotated Accelerations
   end;
end;

clear Platacc Acclow kk;

%disp('Integrating')

%Fc=1/60;           %seconds
Fc=1/150;                 %10 seconds
tau=1/(Fc*2*pi);bc=Fc/10;[bf af]=butter(1,bc);dt=1/20;
om_c=tan(dt/2/tau);tau_c=dt/2/om_c;
clear Fc;

Platrate(1,:)=Platform(1,:)*tau_c;
Platrate(2,:)=Platform(2,:)*tau_c;
Platrate(3,:)=Platform(3,:)*tau_c;

[Platvel(1,:),Zf]=filter(bf,af,Platrate(1,:),0);
[Platvel(2,:),Zf]=filter(bf,af,Platrate(2,:),0);
[Platvel(3,:),Zf]=filter(bf,af,Platrate(3,:),0);

clear af bf dt ii jj om_c zz tau_c bc G ans T11 T12 T13 T21 T22 T23 T31 T32;
clear T33 Bf0 Bf1 Af0 tau lv dfx accx;
clear Platrate Platform;

Result=Velrot;

Result(1,:)=Velrot(1,:)+Platvel(1,:);
Result(2,:)=Velrot(2,:)+Platvel(2,:);
Result(3,1:length(Result(3,:))-5)=Velrot(3,1:length(Result(3,:))-5)+Platvel(3,6:length(Result(3,:)));

Result(3,length(Result(1,:))-4:length(Result(1,:)))=Velrot(3,length(Result(1,:))-4:length(Result(1,:)))+Platvel(3,length(Result(1,:))-4:length(Result(1,:)));

[A,Fa]=psd(Velrot(1,:),65536,20);[sA,sFa]=specsmoo(A,20);
[B,Fb]=psd(Velrot(2,:),65536,20);[sB,sFb]=specsmoo(B,20);
[C,Fc]=psd(Velrot(3,:),65536,20);[sC,sFc]=specsmoo(C,20);
figure;loglog(sFa,sA.*sFa,'r-',sFb,sB.*sFb,'g-',sFc,sC.*sFc,'b-');title('Velrot')
[A,Fa]=psd(Platvel(1,:),65536,20);[sA,sFa]=specsmoo(A,20);
[B,Fb]=psd(Platvel(2,:),65536,20);[sB,sFb]=specsmoo(B,20);
[C,Fc]=psd(Platvel(3,:),65536,20);[sC,sFc]=specsmoo(C,20);
figure;loglog(sFa,sA.*sFa,'r-',sFb,sB.*sFb,'g-',sFc,sC.*sFc,'b-');title('Platvel')
[A,Fa]=psd(Result(1,:),65536,20);[sA,sFa]=specsmoo(A,20);
[B,Fb]=psd(Result(2,:),65536,20);[sB,sFb]=specsmoo(B,20);
[C,Fc]=psd(Result(3,:),65536,20);[sC,sFc]=specsmoo(C,20);
figure;loglog(sFa,sA.*sFa,'r-',sFb,sB.*sFb,'g-',sFc,sC.*sFc,'b-');title('Result')

%u=Velrot(1,:);
%v=Velrot(2,:);
%w=Velrot(3,:);

%rotcorr;           %converts [u,v,w] to [Ustream, Vstream, Wstream]
%clear Uu Vv Ww ans tm u v w;
%clear Umean Vmean Wmean dfl ss ut uu uv uw vt vv vw wt ww;
%clear yy alpha beta Ub Vb Wb Sb FF Svu1 Spu1 Sru1 Svv1 Spv1 Srv1 Svw1 Spw1 Srw1
%b4a='Before';
%sonicplots;
%

u=Result(1,:);
v=Result(2,:);
w=Result(3,:);

rotcorr;           %converts [u,v,w] to [Ustream, Vstream, Wstream]
clear Uu Vv Ww ans tm u v w;
clear Umean Vmean Wmean dfl ss ut uu uv uw vt vv vw wt ww;
clear yy alpha beta Ub Vb Wb Sb FF Svu1 Spu1 Sru1 Svv1 Spv1 Srv1 Svw1 Spw1 Srw1
UU=[];UW=[];VV=[];VW=[];WW=[];WT=[];WQ=[];
b4a='After '
%sonicplt;

firstq=fix(length(Wstream)/6*1);
seconq=fix(length(Wstream)/6*2);
thirdq=fix(length(Wstream)/6*3);
fourtq=fix(length(Wstream)/6*4);
fifthq=fix(length(Wstream)/6*5);
sixthq=fix(length(Wstream)/6*6);
Us=0;Rs=150;Rl=300;zi=600;P=1013;zu=20;zt=20;zq=20;

if length(ship_Tsea)<60,
   ship_Tsea=[ship_Tsea (1:60-length(ship_Tsea)).*0+ship_Tsea(length(ship_Tsea))];
end;   
if length(Tsea_snake)<60,
   Tsea_snake=[Tsea_snake (1:60-length(Tsea_snake)).*0+Tsea_snake(length(Tsea_snake))];
end;   
if length(ship_Tair)<60,
   ship_Tair=[ship_Tair (1:60-length(ship_Tair)).*0+ship_Tair(length(ship_Tair))];
end;   
if length(rain)<60,
   rain=[rain (1:60-length(rain)).*0+rain(length(rain))];
end;   

wu_1=cov(Ustream(1:firstq),Wstream(1:firstq));
wv_1=cov(Vstream(1:firstq),Wstream(1:firstq));
wt_1=cov(Wstream(1:firstq),sst(1:firstq));
Ubar_1=sqrt(mean(Ustream(1:firstq)).^2+mean(Vstream(1:firstq)).^2+mean(Wstream(1:firstq)).^2);
Tbar_1=mean(ship_Tair(1:10));
rho_1=P*100./(273.15+Tbar_1)./287;
wq_1=cov(Wstream(1:firstq),abshum(1:firstq)./rho_1);
Qbar_1=mean(abshum(1:firstq))./rho_1;
Tsea_1=mean(Tsea_snake(1:10));
Esat=6.112.*exp(17.67.*Tsea_1./(Tsea_1+243.5))*0.98;
Qsea_1=Esat*0.622/(P-0.378*Esat)*1000;
Ss_1=std(shead(1:firstq));Sc_1=std(chead(1:firstq));
Sp_1=std(Platvel(1:firstq));
rain_1=mean(rain(1:10));
x_1=[Ubar_1 Us Tsea_1 Tbar_1 Qsea_1 Qbar_1 Rs Rl rain_1 zi P zu zt zq];
y_1=cool_lkb(x_1);

wu_2=cov(Ustream(firstq+1:seconq),Wstream(firstq+1:seconq));
wv_2=cov(Vstream(firstq+1:seconq),Wstream(firstq+1:seconq));
wt_2=cov(Wstream(firstq+1:seconq),sst(firstq+1:seconq));
Ubar_2=sqrt(mean(Ustream(firstq+1:seconq)).^2+mean(Vstream(firstq+1:seconq)).^2+mean(Wstream(firstq+1:seconq)).^2);
Tbar_2=mean(ship_Tair(11:20));
rho_2=P*100./(273.15+Tbar_2)./287;
wq_2=cov(Wstream(firstq+1:seconq),abshum(firstq+1:seconq)./rho_2);
Qbar_2=mean(abshum(firstq+1:seconq))./rho_2;
Tsea_2=mean(Tsea_snake(11:20));
Esat=6.112.*exp(17.67.*Tsea_2./(Tsea_2+243.5))*0.98;
Qsea_2=Esat*0.622/(P-0.378*Esat)*1000;
Ss_2=std(shead(firstq+1:seconq));Sc_2=std(chead(firstq+1:seconq));
Sp_2=std(Platvel(firstq+1:seconq));
rain_2=mean(rain(11:20));
x_2=[Ubar_2 Us Tsea_2 Tbar_2 Qsea_2 Qbar_2 Rs Rl rain_2 zi P zu zt zq];
y_2=cool_lkb(x_2);

wu_3=cov(Ustream(seconq+1:thirdq),Wstream(seconq+1:thirdq));
wv_3=cov(Vstream(seconq+1:thirdq),Wstream(seconq+1:thirdq));
wt_3=cov(Wstream(seconq+1:thirdq),sst(seconq+1:thirdq));
Ubar_3=sqrt(mean(Ustream(seconq+1:thirdq)).^2+mean(Vstream(seconq+1:thirdq)).^2+mean(Wstream(seconq+1:thirdq)).^2);
Tbar_3=mean(ship_Tair(21:30));
rho_3=P*100./(273.15+Tbar_3)./287;
wq_3=cov(Wstream(seconq+1:thirdq),abshum(seconq+1:thirdq)./rho_3);
Qbar_3=mean(abshum(seconq+1:thirdq))./rho_3;
Tsea_3=mean(Tsea_snake(21:30));
Esat=6.112.*exp(17.67.*Tsea_3./(Tsea_3+243.5))*0.98;
Qsea_3=Esat*0.622/(P-0.378*Esat)*1000;
Ss_3=std(shead(seconq+1:thirdq));Sc_3=std(chead(seconq+1:thirdq));
Sp_3=std(Platvel(seconq+1:thirdq));
rain_3=mean(rain(21:30));
x_3=[Ubar_3 Us Tsea_3 Tbar_3 Qsea_3 Qbar_3 Rs Rl rain_3 zi P zu zt zq];
y_3=cool_lkb(x_3);

wu_4=cov(Ustream(thirdq+1:fourtq),Wstream(thirdq+1:fourtq));
wv_4=cov(Vstream(thirdq+1:fourtq),Wstream(thirdq+1:fourtq));
wt_4=cov(Wstream(thirdq+1:fourtq),sst(thirdq+1:fourtq));
Ubar_4=sqrt(mean(Ustream(thirdq+1:fourtq)).^2+mean(Vstream(thirdq+1:fourtq)).^2+mean(Wstream(thirdq+1:fourtq)).^2);
Tbar_4=mean(ship_Tair(31:40));
rho_4=P*100./(273.15+Tbar_4)./287;
wq_4=cov(Wstream(seconq+1:thirdq),abshum(seconq+1:thirdq)./rho_4);
Qbar_4=mean(abshum(thirdq+1:fourtq))./rho_4;
Tsea_4=mean(Tsea_snake(31:40));
Esat=6.112.*exp(17.67.*Tsea_4./(Tsea_4+243.5))*0.98;
Qsea_4=Esat*0.622/(P-0.378*Esat)*1000;
Ss_4=std(shead(thirdq+1:fourtq));Sc_4=std(chead(thirdq+1:fourtq));
Sp_4=std(Platvel(thirdq+1:fourtq));
rain_4=mean(rain(31:40));
x_4=[Ubar_4 Us Tsea_4 Tbar_4 Qsea_4 Qbar_4 Rs Rl rain_4 zi P zu zt zq];
y_4=cool_lkb(x_4);

wu_5=cov(Ustream(fourtq+1:fifthq),Wstream(fourtq+1:fifthq));
wv_5=cov(Vstream(fourtq+1:fifthq),Wstream(fourtq+1:fifthq));
wt_5=cov(Wstream(fourtq+1:fifthq),sst(fourtq+1:fifthq));
Ubar_5=sqrt(mean(Ustream(fourtq+1:fifthq)).^2+mean(Vstream(fourtq+1:fifthq)).^2+mean(Wstream(fourtq+1:fifthq)).^2);
Tbar_5=mean(ship_Tair(41:50));
rho_5=P*100./(273.15+Tbar_5)./287;
wq_5=cov(Wstream(thirdq+1:fourtq),abshum(thirdq+1:fourtq)./rho_5);
Qbar_5=mean(abshum(fourtq+1:fifthq))./rho_5;
Tsea_5=mean(Tsea_snake(41:50));
Esat=6.112.*exp(17.67.*Tsea_5./(Tsea_5+243.5))*0.98;
Qsea_5=Esat*0.622/(P-0.378*Esat)*1000;
Ss_5=std(shead(fourtq+1:fifthq));Sc_5=std(chead(fourtq+1:fifthq));
Sp_5=std(Platvel(fourtq+1:fifthq));
rain_5=mean(rain(41:50));
x_5=[Ubar_5 Us Tsea_5 Tbar_5 Qsea_5 Qbar_5 Rs Rl rain_5 zi P zu zt zq];
y_5=cool_lkb(x_5);

wu_6=cov(Ustream(fifthq+1:sixthq),Wstream(fifthq+1:sixthq));
wv_6=cov(Vstream(fifthq+1:sixthq),Wstream(fifthq+1:sixthq));
wt_6=cov(Wstream(fifthq+1:sixthq),sst(fifthq+1:sixthq));
Ubar_6=sqrt(mean(Ustream(fifthq+1:sixthq)).^2+mean(Vstream(fifthq+1:sixthq)).^2+mean(Wstream(fifthq+1:sixthq)).^2);
Tbar_6=mean(ship_Tair(51:length(ship_Tair)));
rho_6=P*100./(273.15+Tbar_6)./287;
wq_6=cov(Wstream(fourtq+1:fifthq),abshum(fourtq+1:fifthq)./rho_6);
Qbar_6=mean(abshum(fifthq+1:sixthq))./rho_6;
Tsea_6=mean(Tsea_snake(51:length(Tsea_snake)));
Esat=6.112.*exp(17.67.*Tsea_6./(Tsea_6+243.5))*0.98;
Qsea_6=Esat*0.622/(P-0.378*Esat)*1000;
Ss_6=std(shead(fifthq+1:sixthq));Sc_6=std(chead(fifthq+1:sixthq));
Sp_6=std(Platvel(fifthq+1:sixthq));
rain_6=mean(rain(51:60));
x_6=[Ubar_6 Us Tsea_6 Tbar_6 Qsea_6 Qbar_6 Rs Rl rain_6 zi P zu zt zq];
y_6=cool_lkb(x_6);

y_y=[y_1 Ubar_1 Tbar_1 Tsea_1 Qbar_1 wu_1(1,2) wv_1(1,2) wt_1(1,2) wq_1(1,2) rho_1 stdC_1 stdA_1 medC_1 medA_1 Ss_1 Sc_1 Sp_1 rain_1
     y_2 Ubar_2 Tbar_2 Tsea_2 Qbar_2 wu_2(1,2) wv_2(1,2) wt_2(1,2) wq_2(1,2) rho_2 stdC_2 stdA_2 medC_2 medA_2 Ss_2 Sc_2 Sp_2 rain_2
     y_3 Ubar_3 Tbar_3 Tsea_3 Qbar_3 wu_3(1,2) wv_3(1,2) wt_3(1,2) wq_3(1,2) rho_3 stdC_3 stdA_3 medC_3 medA_3 Ss_3 Sc_3 Sp_3 rain_3
     y_4 Ubar_4 Tbar_4 Tsea_4 Qbar_4 wu_4(1,2) wv_4(1,2) wt_4(1,2) wq_4(1,2) rho_4 stdC_4 stdA_4 medC_4 medA_4 Ss_4 Sc_4 Sp_4 rain_4
     y_5 Ubar_5 Tbar_5 Tsea_5 Qbar_5 wu_5(1,2) wv_5(1,2) wt_5(1,2) wq_5(1,2) rho_5 stdC_5 stdA_5 medC_5 medA_5 Ss_5 Sc_5 Sp_5 rain_5
     y_6 Ubar_6 Tbar_6 Tsea_6 Qbar_6 wu_6(1,2) wv_6(1,2) wt_6(1,2) wq_6(1,2) rho_6 stdC_6 stdA_6 medC_6 medA_6 Ss_6 Sc_6 Sp_6 rain_6];

%clear y_1 y_2 y_3 y_4 y_5 y_6 Ubar_1 Ubar_2 Ubar_3 Ubar_4 Ubar_5 Ubar_6 Tbar_1 Tbar_2 Tbar_3 Tbar_4 Tbar_5 Tbar_6
%clear Tsea_1 Tsea_2 Tsea_3 Tsea_4 Tsea_5 Tsea_6 Qbar_1 Qbar_2 Qbar_3 Qbar_4 Qbar_5 Qbar_6
%clear rho_1 rho_2 rho_3 rho_4 rho_5 rho_6 zi zq zt zu
%   y=[hsb hlb tau zo zot zoq L usr tsr qsr dter dqer tkt RF hl_webb Cd Ch Ce Cdn_10 Chn_10 Cen_10];
%       ********** output data  ***************
%	01 hsb=			sensible heat flux (w/m^2)
%	02 hlb=			latent heat flux (w/m^2)
%	03 tau=			stress (nt/m^2)
%	04 zo=			velocity roughness length (m)
%	05 zot			temperature roughness length (m)
%	06 zoq=			moisture roughness length (m)
%	07 L=			Monin_Obukhov stability length
%	08 usr=			turbulent friction velocity (m/s), including gustiness
%	09 tsr			temperature scaling parameter (K)
%	10 qsr			humidity scaling parameter (g/g)
%	11 dter=		cool skin temperature depression (K)
%	12 dqer=		cool skin humidity depression (g/g)
%	13 tkt=			cool skin thickness (m)
%	14 RF=			rain heat flux(w/m^2)
%	15 hl_webb=		webb addition to latent flux (w/m^2)
%	16 Cd=			velocity drag coefficient at zu, referenced to u
%	17 Ch=			heat transfer coefficient at zt
%	18 Ce=			moisture transfer coefficient at zq
%	19 Cdn_10=		10-m velocity drag coeeficient, including gustiness
%	20 Chn_10=		10-m heat transfer coeeficient, including gustiness
%	21 Cen_10=		10-m humidity transfer coeeficient, including gustiness

%clear Chead Esat FF H_s Mcog Mhead Msog P Platvel Qsea_1 Qsea_2 Qsea_3 Qsea_4 Qsea_5 Qsea_6 Result Rl Rs Shead
%clear Ts Tsonic Twindir Ustream Vstream Ulong Vlong Wlong Velrot Zf abshum ans b4a Us
%clear heading ii mm pp qabshum qx rain sst sx sy tau_u tau_v tb tmotion u uu uw v
%clear vv vw w wq wq_1 wq_2 wq_3 wq_4 wq_5 wq_6 wt wt_1 wt_2 wt_3 wt_4 wt_5 wt_6 wu_1 wu_2 wu_3 wu_4 wu_5 wu_6
%clear wv_1 wv_2 wv_3 wv_4 wv_5 wv_6 ww x_1 x_2 x_3 x_4 x_5 x_6

%ss=['save ',sfl,' y_y -ascii -double'];disp(ss)
%eval(ss);
