disp('nau_correct_GOMECC_07.m')

Licrate=10;%data rate in Hz
Lenlic=3600*Licrate;
Nlenlic=floor(Lenlic/6);

G=sqrt(mean(accx)^2+mean(accy)^2+mean(accz)^2);

ratex=detrend(ratex);
ratey=detrend(ratey);
ratez=detrend(ratez);
maccx=mean(accx);
maccy=mean(accy);
maccz=mean(accz);

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

wp=1/40/(Licrate/2);
ws=1/20/(Licrate/2);   %Edson coeff...
[n,wn]=buttord(wp,ws,1,Licrate);  

[bhigh,ahigh]=butter(n,wn,'high');

clear Fc wn wp ws n

sdg=['load ', way_data, '\fxyz.prn fxyz -ascii'];
disp(sdg)
eval(sdg)

fx=fxyz(1);
fy=fxyz(2);
fz=fxyz(3);
px=fxyz(4);
py=fxyz(5);
pz=fxyz(6);

if fx==0 | isempty(fx),
   fx=maccy/G;
end;
if fy==0 | isempty(fy),
   fy=maccx/G;
end;
if fz==0 | isempty(fz),
   fz=mean(heading)/180*pi;
end;

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

thetaslowh = atan2(-accx,G) - filtfilt(bhigh,ahigh,atan2(-accx,G));
theta = thetaslowh + filtfilt(bhigh,ahigh,(cumsum(ratey) -.5*ratey -.5*ratey(1))/Licrate);

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

phislowh = atan2(accy,G) - filtfilt(bhigh,ahigh,atan2(accy,G));
phi = phislowh + filtfilt(bhigh,ahigh,(cumsum(ratex) - .5*ratex - .5*ratex(1))/Licrate);

heading=heading/180*pi;                     %convert to radians
psirate=ratez*tau_c-heading;                              %+ to - 09/15/98
%psi=filter(bf,af,psirate,mean(heading));
% psi=filter(bf,af,psirate,fz);

psislowh = (-heading - filtfilt(bhigh,ahigh,-heading));
psi = psislowh + filtfilt(bhigh,ahigh,(cumsum(ratez) -.5*ratez - .5*ratez(1))/Licrate);

fx=filtic(bf,af,[phi(length(phi)) phi(length(phi)-1)],[phirate(length(phirate)) phirate(length(phirate)-1)]);
fy=filtic(bf,af,[theta(length(theta)) theta(length(theta)-1)],[thetarate(length(thetarate)) thetarate(length(thetarate)-1)]);
fz=filtic(bf,af,[psi(length(psi)) psi(length(psi)-1)],[psirate(length(psirate)) psirate(length(psirate)-1)]);

%Compute Omega X R for sonics

Rz= 0.7;
Rx= 0;
Ry= 0;

OXRx=ratey*Rz-ratez*Ry;     
OXRy=-ratex*Rz+ratez*Rx;     
OXRz=+ratex*Ry-ratey*Rx;      

OM=[OXRx OXRy OXRz];
OM=OM';

%TOM=0.*OM';

clear maccx maccy maccz ahigh bhigh
clear thetarate phirate psirate OXRx OXRy OXRz

%clear OXRy OXRz Rx Ry Rz;

Vel=[U              %
     V              %
     W];            %

Plat=[accx
      accy
      accz];
  
Platacc=Plat;
Velrot=Vel;

%*****************************************************************
%            Start Edson correction
%*****************************************************************

deg_rate=[ratex
   ratey
   ratez];

if clearit
clear ratex ratey ratez
end;

R=[Rx
   Ry
   Rz];

[bhi,ahi,bf2,af2,tau2,spinup2]=Filtcoef(Licrate,1/30);

[euler,dr]=angles(ahi,bhi,Licrate,Plat,deg_rate,heading);

roll=euler(1,:);  roll=resample(roll,1,Licrate);   roll=roll';
pitch=euler(2,:); pitch=resample(pitch,1,Licrate); pitch=pitch';

tilts=[roll pitch];

%sdg=['save d:\gasex2\ti_',dfx,'.dat tilts -ascii'];
%disp(sdg)
%eval(sdg)

%break

[uvwplat,xyzplat]=accels(bhi,ahi,Licrate,Plat,euler);clear xyzplat
[uvw]=sonic(Vel,dr,euler,uvwplat,R);

AAAA=uvw;

Jcovtest=cov(uvw(3,:),uvwplat(3,:));
Jcovtest=Jcovtest(1,2);
lagflag=0;
for lag=1:16;
   Jcovtestlag=cov(uvw(3,1:length(uvw)-lag),uvwplat(3,lag+1:length(uvw)));
   Jcovtestlag=Jcovtestlag(1,2);
   if Jcovtestlag < 0 & Jcovtest < 0,
      if abs(Jcovtestlag) > abs(Jcovtest),
         Jcovtest=Jcovtestlag;
         lagflag=lag;
      end;
   elseif Jcovtestlag > 0 & Jcovtest > 0,
      if abs(Jcovtestlag) < abs(Jcovtest),
         Jcovtest=Jcovtestlag;
         lagflag=lag;
      end;
   elseif Jcovtestlag < 0 & Jcovtest > 0,
         Jcovtest=Jcovtestlag;
         lagflag=lag;
   end;
end;
for lag=-1:-1:-16;
   Jcovtestlag=cov(uvw(3,-lag+1:length(uvw)),uvwplat(3,1:length(uvw)+lag));
   Jcovtestlag=Jcovtestlag(1,2);
   if Jcovtestlag < 0 & Jcovtest < 0,
      if abs(Jcovtestlag) > abs(Jcovtest),
         Jcovtest=Jcovtestlag;
         lagflag=lag;
      end;
   elseif Jcovtestlag > 0 & Jcovtest > 0,
      if abs(Jcovtestlag) < abs(Jcovtest),
         Jcovtest=Jcovtestlag;
         lagflag=lag;
      end;
   elseif Jcovtestlag < 0 & Jcovtest > 0,
         Jcovtest=Jcovtestlag;
         lagflag=lag;
   end;
end;

uvw(1,:)=uvw(1,:) + uvwplat(1,:);
uvw(2,:)=uvw(2,:) + uvwplat(2,:);
CCCC=uvwplat;

if lagflag > 0,
   uvw(3,1:length(uvw)-lagflag)=uvw(3,1:length(uvw)-lagflag) + uvwplat(3,lagflag+1:length(uvw));
   uvw(3,length(uvw)-lagflag+1:length(uvw))= uvw(3,length(uvw)-lagflag+1:length(uvw)) + uvwplat(3,1:lagflag);
else
   uvw(3,-lagflag+1:length(uvw)) = uvw(3,-lagflag+1:length(uvw)) + uvwplat(3,1:length(uvw)+lagflag);
   uvw(3,1:-lagflag) = uvw(3,1:-lagflag) + uvwplat(3,length(uvw)+lagflag+1:length(uvw));
end;   
lagflag
DDDD=uvw;
DDDD=DDDD';

%Zap=['save c:\byron\W',dfl(13:18),'.dat DDDD -ascii'];
%eval(Zap)

%*****************************************************************
%            End Edson correction
%*****************************************************************

%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 psi phi R Rx Ry Rz Vel af2 ahi bf2 bhi deg_rate dr %euler

%matrix for a right-handed coordinate system
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 sin_theta cos_phi sin_phi cos_psi sin_psi;


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.*(U+OM(1,:)) + T12.*(V+OM(2,:)) + T13.*(W+OM(3,:));
Velrot(2,:) = T21.*(U+OM(1,:)) + T22.*(V+OM(2,:)) + T23.*(W+OM(3,:));
Velrot(3,:) = T31.*(U+OM(1,:)) + T32.*(V+OM(2,:)) + T33.*(W+OM(3,:));

if clearit
clear Plat T OM ssx ssy ssz accy accz;
end;

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;                 %10 seconds

tau=1/(Fc*2*pi);bc=Fc/(Licrate/2);[bf af]=butter(1,bc);dt=1/Licrate;
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);      
   end;
end;

clear Platacc Acclow kk

%disp('Integrating')

Fc=1/60;           %seconds

tau=1/(Fc*2*pi);bc=Fc/(Licrate/2);[bf af]=butter(1,bc);dt=1/Licrate;
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;

if isempty(px)| isnan(px),
   px=0;
end;
if isempty(py)|| isnan(py),
   py=0;
end;
if isempty(pz)|| isnan(pz),
   pz=0;
end;

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

px=filtic(bf,af,[Platvel(1,length(Platvel)) Platvel(1,length(Platvel)-1)],[Platrate(1,length(Platrate)) Platrate(1,length(Platrate)-1)]);
py=filtic(bf,af,[Platvel(2,length(Platvel)) Platvel(2,length(Platvel)-1)],[Platrate(2,length(Platrate)) Platrate(2,length(Platrate)-1)]);
pz=filtic(bf,af,[Platvel(3,length(Platvel)) Platvel(3,length(Platvel)-1)],[Platrate(3,length(Platrate)) Platrate(3,length(Platrate)-1)]);

if isempty(px)| isnan(px),
   px=0;
end;
if isempty(py)|| isnan(py),
   py=0;
end;
if isempty(pz)|| isnan(pz),
   pz=0;
end;

fxyz=[fx fy fz px py pz];

sdg=['save ', drive_w, '\fxyz.prn fxyz -ascii'];
disp(sdg)
eval(sdg)

if clearit
clear 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 accx;
clear Platrate Platform;
end;

Result=Velrot;

Jcovtest=cov(Velrot(3,:),Platvel(3,:));
Jcovtest=Jcovtest(1,2);
lagflag=0;

for lag=1:16;
   Jcovtestlag=cov(Velrot(3,1:length(Velrot)-lag),Platvel(3,lag+1:length(Velrot)));
   Jcovtestlag=Jcovtestlag(1,2);
   if Jcovtestlag < 0 & Jcovtest < 0,
      if abs(Jcovtestlag) > abs(Jcovtest),
         Jcovtest=Jcovtestlag;
         lagflag=lag;
      end;
   elseif Jcovtestlag > 0 & Jcovtest > 0,
      if abs(Jcovtestlag) < abs(Jcovtest),
         Jcovtest=Jcovtestlag;
         lagflag=lag;
      end;
   elseif Jcovtestlag < 0 & Jcovtest > 0,
         Jcovtest=Jcovtestlag;
         lagflag=lag;
   end;
end;
for lag=-1:-1:-16;
   Jcovtestlag=cov(Velrot(3,-lag+1:length(Velrot)),Platvel(3,1:length(Velrot)+lag));
   Jcovtestlag=Jcovtestlag(1,2);
   if Jcovtestlag < 0 & Jcovtest < 0,
      if abs(Jcovtestlag) > abs(Jcovtest),
         Jcovtest=Jcovtestlag;
         lagflag=lag;
      end;
   elseif Jcovtestlag > 0 & Jcovtest > 0,
      if abs(Jcovtestlag) < abs(Jcovtest),
         Jcovtest=Jcovtestlag;
         lagflag=lag;
      end;
   elseif Jcovtestlag < 0 & Jcovtest > 0,
         Jcovtest=Jcovtestlag;
         lagflag=lag;
   end;
end;

if clearit
clear Jcovtest Jcovtestlag lag
end;

Result(1,:)=Velrot(1,:)+Platvel(1,:);
Result(2,:)=Velrot(2,:)+Platvel(2,:);

if lagflag > 0,
   Result(3,1:length(Velrot)-lagflag)=Velrot(3,1:length(Velrot)-lagflag) + Platvel(3,lagflag+1:length(Velrot));
   Result(3,length(Velrot)-lagflag+1:length(Velrot))= Velrot(3,length(Velrot)-lagflag+1:length(Velrot)) + Platvel(3,1:lagflag);
else
   Result(3,-lagflag+1:length(Velrot)) = Velrot(3,-lagflag+1:length(Velrot)) + Platvel(3,1:length(Velrot)+lagflag);
   Result(3,1:-lagflag) = Velrot(3,1:-lagflag) + Platvel(3,length(Velrot)+lagflag+1:length(Velrot));
end;
lagflag
if clearit
clear lagflag
end;

[Ustream,Vstream,Wstream,alpha,beta]=rotcorr(Result(1,:),Result(2,:),Result(3,:));           

clear Uu Vv Ww ans tm Umean Vmean Wmean 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

Us=0;
zi=600;
P=1013;

%******************************************************************

zu=18.5;zu_ship=12;
zt=15.5;zt_ship=12;
zq=15.5;zq_ship=12;

%**************************************!!!!!!!!!!!!!!!!!!!!!!!!
%abshum=abshum.*1.2;        %Gain Correction Only
%**************************************!!!!!!!!!!!!!!!!!!!!!!!!

%*********************************************************************************
%       Compute products
%*********************************************************************************

%*************************Minute 0-10********************************************

mean_sog_1=mean(sogm(01:10));if mean_sog_1 < 0.2, mean_sog_1=0;end;std_sog_1=std(sogm(01:10));min_sog_1=min(sogm(01:10));max_sog_1=max(sogm(01:10));
mean_sog_2=mean(sogm(11:20));if mean_sog_2 < 0.2, mean_sog_2=0;end;std_sog_2=std(sogm(11:20));min_sog_2=min(sogm(11:20));max_sog_2=max(sogm(11:20));
mean_sog_3=mean(sogm(21:30));if mean_sog_3 < 0.2, mean_sog_3=0;end;std_sog_3=std(sogm(21:30));min_sog_3=min(sogm(21:30));max_sog_3=max(sogm(21:30));
mean_sog_4=mean(sogm(31:40));if mean_sog_4 < 0.2, mean_sog_4=0;end;std_sog_4=std(sogm(31:40));min_sog_4=min(sogm(31:40));max_sog_4=max(sogm(31:40));
mean_sog_5=mean(sogm(41:50));if mean_sog_5 < 0.2, mean_sog_5=0;end;std_sog_5=std(sogm(41:50));min_sog_5=min(sogm(41:50));max_sog_5=max(sogm(41:50));
mean_sog_6=mean(sogm(51:60));if mean_sog_6 < 0.2, mean_sog_6=0;end;std_sog_6=std(sogm(51:60));min_sog_6=min(sogm(51:60));max_sog_6=max(sogm(51:60));

mean_cog_1=mod(360+180.*atan2(mean(sin(cogm(01:10)./180.*pi)),mean(cos(cogm(01:10)./180.*pi)))./pi,360);
mean_cog_2=mod(360+180.*atan2(mean(sin(cogm(11:20)./180.*pi)),mean(cos(cogm(11:20)./180.*pi)))./pi,360);
mean_cog_3=mod(360+180.*atan2(mean(sin(cogm(21:30)./180.*pi)),mean(cos(cogm(21:30)./180.*pi)))./pi,360);
mean_cog_4=mod(360+180.*atan2(mean(sin(cogm(31:40)./180.*pi)),mean(cos(cogm(31:40)./180.*pi)))./pi,360);
mean_cog_5=mod(360+180.*atan2(mean(sin(cogm(41:50)./180.*pi)),mean(cos(cogm(41:50)./180.*pi)))./pi,360);
mean_cog_6=mod(360+180.*atan2(mean(sin(cogm(51:60)./180.*pi)),mean(cos(cogm(51:60)./180.*pi)))./pi,360);

std_cog_1=std(unwrap(cog(0001:0300)./180.*pi)*180/pi);
std_cog_2=std(unwrap(cog(0301:0600)./180.*pi)*180/pi);
std_cog_3=std(unwrap(cog(0601:0900)./180.*pi)*180/pi);
std_cog_4=std(unwrap(cog(0901:1200)./180.*pi)*180/pi);
std_cog_5=std(unwrap(cog(1201:1500)./180.*pi)*180/pi);
std_cog_6=std(unwrap(cog(1501:1800)./180.*pi)*180/pi);

max_cog_1=max(cog(0001:0300));
max_cog_2=max(cog(0301:0600));
max_cog_3=max(cog(0601:0900));
max_cog_4=max(cog(0901:1200));
max_cog_5=max(cog(1201:1500));
max_cog_6=max(cog(1501:1800));

min_cog_1=min(cog(0001:0300));
min_cog_2=min(cog(0301:0600));
min_cog_3=min(cog(0601:0900));
min_cog_4=min(cog(0901:1200));
min_cog_5=min(cog(1201:1500));
min_cog_6=min(cog(1501:1800));

mean_sinsog_1=mean(sogm(01:10).*sin(cogm(01:10)./180.*pi));std_sinsog_1=std(sogm(01:10).*sin(cogm(01:10)./180.*pi));min_sinsog_1=min(sogm(01:10).*sin(cogm(01:10)./180.*pi));max_sinsog_1=max(sogm(01:10).*sin(cogm(01:10)./180.*pi));
mean_sinsog_2=mean(sogm(11:20).*sin(cogm(11:20)./180.*pi));std_sinsog_2=std(sogm(11:20).*sin(cogm(11:20)./180.*pi));min_sinsog_2=min(sogm(11:20).*sin(cogm(11:20)./180.*pi));max_sinsog_2=max(sogm(11:20).*sin(cogm(11:20)./180.*pi));
mean_sinsog_3=mean(sogm(21:30).*sin(cogm(21:30)./180.*pi));std_sinsog_3=std(sogm(21:30).*sin(cogm(21:30)./180.*pi));min_sinsog_3=min(sogm(21:30).*sin(cogm(21:30)./180.*pi));max_sinsog_3=max(sogm(21:30).*sin(cogm(21:30)./180.*pi));
mean_sinsog_4=mean(sogm(31:40).*sin(cogm(31:40)./180.*pi));std_sinsog_4=std(sogm(31:40).*sin(cogm(31:40)./180.*pi));min_sinsog_4=min(sogm(31:40).*sin(cogm(31:40)./180.*pi));max_sinsog_4=max(sogm(31:40).*sin(cogm(31:40)./180.*pi));
mean_sinsog_5=mean(sogm(41:50).*sin(cogm(41:50)./180.*pi));std_sinsog_5=std(sogm(41:50).*sin(cogm(41:50)./180.*pi));min_sinsog_5=min(sogm(41:50).*sin(cogm(41:50)./180.*pi));max_sinsog_5=max(sogm(41:50).*sin(cogm(41:50)./180.*pi));
mean_sinsog_6=mean(sogm(51:60).*sin(cogm(51:60)./180.*pi));std_sinsog_6=std(sogm(51:60).*sin(cogm(51:60)./180.*pi));min_sinsog_6=min(sogm(51:60).*sin(cogm(51:60)./180.*pi));max_sinsog_6=max(sogm(51:60).*sin(cogm(51:60)./180.*pi));

mean_cossog_1=mean(sogm(01:10).*cos(cogm(01:10)./180.*pi));std_cossog_1=std(sogm(01:10).*cos(cogm(01:10)./180.*pi));min_cossog_1=min(sogm(01:10).*cos(cogm(01:10)./180.*pi));max_cossog_1=max(sogm(01:10).*cos(cogm(01:10)./180.*pi));
mean_cossog_2=mean(sogm(11:20).*cos(cogm(11:20)./180.*pi));std_cossog_2=std(sogm(11:20).*cos(cogm(11:20)./180.*pi));min_cossog_2=min(sogm(11:20).*cos(cogm(11:20)./180.*pi));max_cossog_2=max(sogm(11:20).*cos(cogm(11:20)./180.*pi));
mean_cossog_3=mean(sogm(21:30).*cos(cogm(21:30)./180.*pi));std_cossog_3=std(sogm(21:30).*cos(cogm(21:30)./180.*pi));min_cossog_3=min(sogm(21:30).*cos(cogm(21:30)./180.*pi));max_cossog_3=max(sogm(21:30).*cos(cogm(21:30)./180.*pi));
mean_cossog_4=mean(sogm(31:40).*cos(cogm(31:40)./180.*pi));std_cossog_4=std(sogm(31:40).*cos(cogm(31:40)./180.*pi));min_cossog_4=min(sogm(31:40).*cos(cogm(31:40)./180.*pi));max_cossog_4=max(sogm(31:40).*cos(cogm(31:40)./180.*pi));
mean_cossog_5=mean(sogm(41:50).*cos(cogm(41:50)./180.*pi));std_cossog_5=std(sogm(41:50).*cos(cogm(41:50)./180.*pi));min_cossog_5=min(sogm(41:50).*cos(cogm(41:50)./180.*pi));max_cossog_5=max(sogm(41:50).*cos(cogm(41:50)./180.*pi));
mean_cossog_6=mean(sogm(51:60).*cos(cogm(51:60)./180.*pi));std_cossog_6=std(sogm(51:60).*cos(cogm(51:60)./180.*pi));min_cossog_6=min(sogm(51:60).*cos(cogm(51:60)./180.*pi));max_cossog_6=max(sogm(51:60).*cos(cogm(51:60)./180.*pi));

Result(1,0*Nlenlic+1:1*Nlenlic)=Result(1,0*Nlenlic+1:1*Nlenlic)+mean_cossog_1-mean(Platvel(1,0*Nlenlic+1:1*Nlenlic));   %mean_cossog=slow frequency velocity.
Result(1,1*Nlenlic+1:2*Nlenlic)=Result(1,1*Nlenlic+1:2*Nlenlic)+mean_cossog_2-mean(Platvel(1,1*Nlenlic+1:2*Nlenlic));   %mean(Platvel)=high frequency velocity.
Result(1,2*Nlenlic+1:3*Nlenlic)=Result(1,2*Nlenlic+1:3*Nlenlic)+mean_cossog_3-mean(Platvel(1,2*Nlenlic+1:3*Nlenlic));    %see Edson et Al. "direct covariance flux estimates from Mobile platforms ar Sea" p552
Result(1,3*Nlenlic+1:4*Nlenlic)=Result(1,3*Nlenlic+1:4*Nlenlic)+mean_cossog_4-mean(Platvel(1,3*Nlenlic+1:4*Nlenlic));
Result(1,4*Nlenlic+1:5*Nlenlic)=Result(1,4*Nlenlic+1:5*Nlenlic)+mean_cossog_5-mean(Platvel(1,4*Nlenlic+1:5*Nlenlic));
Result(1,5*Nlenlic+1:6*Nlenlic)=Result(1,5*Nlenlic+1:6*Nlenlic)+mean_cossog_6-mean(Platvel(1,5*Nlenlic+1:6*Nlenlic));

% Result(2,0*Nlenlic+1:1*Nlenlic)=Result(2,0*Nlenlic+1:1*Nlenlic)+mean_sinsog_1-mean(Platvel(2,0*Nlenlic+1:1*Nlenlic));
% Result(2,1*Nlenlic+1:2*Nlenlic)=Result(2,1*Nlenlic+1:2*Nlenlic)+mean_sinsog_2-mean(Platvel(2,1*Nlenlic+1:2*Nlenlic));
% Result(2,2*Nlenlic+1:3*Nlenlic)=Result(2,2*Nlenlic+1:3*Nlenlic)+mean_sinsog_3-mean(Platvel(2,2*Nlenlic+1:3*Nlenlic));
% Result(2,3*Nlenlic+1:4*Nlenlic)=Result(2,3*Nlenlic+1:4*Nlenlic)+mean_sinsog_4-mean(Platvel(2,3*Nlenlic+1:4*Nlenlic));     %change sign + by - for mean(Platvel)...???
% Result(2,4*Nlenlic+1:5*Nlenlic)=Result(2,4*Nlenlic+1:5*Nlenlic)+mean_sinsog_5-mean(Platvel(2,4*Nlenlic+1:5*Nlenlic));
% Result(2,5*Nlenlic+1:6*Nlenlic)=Result(2,5*Nlenlic+1:6*Nlenlic)+mean_sinsog_6-mean(Platvel(2,5*Nlenlic+1:6*Nlenlic));
Result(2,0*Nlenlic+1:1*Nlenlic)=Result(2,0*Nlenlic+1:1*Nlenlic)-mean_sinsog_1-mean(Platvel(2,0*Nlenlic+1:1*Nlenlic));
Result(2,1*Nlenlic+1:2*Nlenlic)=Result(2,1*Nlenlic+1:2*Nlenlic)-mean_sinsog_2-mean(Platvel(2,1*Nlenlic+1:2*Nlenlic));
Result(2,2*Nlenlic+1:3*Nlenlic)=Result(2,2*Nlenlic+1:3*Nlenlic)-mean_sinsog_3-mean(Platvel(2,2*Nlenlic+1:3*Nlenlic));
Result(2,3*Nlenlic+1:4*Nlenlic)=Result(2,3*Nlenlic+1:4*Nlenlic)-mean_sinsog_4-mean(Platvel(2,3*Nlenlic+1:4*Nlenlic));     %change sign + by - for mean(Platvel)...???
Result(2,4*Nlenlic+1:5*Nlenlic)=Result(2,4*Nlenlic+1:5*Nlenlic)-mean_sinsog_5-mean(Platvel(2,4*Nlenlic+1:5*Nlenlic));
Result(2,5*Nlenlic+1:6*Nlenlic)=Result(2,5*Nlenlic+1:6*Nlenlic)-mean_sinsog_6-mean(Platvel(2,5*Nlenlic+1:6*Nlenlic));

Result(3,0*Nlenlic+1:1*Nlenlic)=Result(3,0*Nlenlic+1:1*Nlenlic)-mean(Platvel(3,0*Nlenlic+1:1*Nlenlic));
Result(3,1*Nlenlic+1:2*Nlenlic)=Result(3,1*Nlenlic+1:2*Nlenlic)-mean(Platvel(3,1*Nlenlic+1:2*Nlenlic));
Result(3,2*Nlenlic+1:3*Nlenlic)=Result(3,2*Nlenlic+1:3*Nlenlic)-mean(Platvel(3,2*Nlenlic+1:3*Nlenlic));
Result(3,3*Nlenlic+1:4*Nlenlic)=Result(3,3*Nlenlic+1:4*Nlenlic)-mean(Platvel(3,3*Nlenlic+1:4*Nlenlic));
Result(3,4*Nlenlic+1:5*Nlenlic)=Result(3,4*Nlenlic+1:5*Nlenlic)-mean(Platvel(3,4*Nlenlic+1:5*Nlenlic));
Result(3,5*Nlenlic+1:6*Nlenlic)=Result(3,5*Nlenlic+1:6*Nlenlic)-mean(Platvel(3,5*Nlenlic+1:6*Nlenlic));

if ~clearit
%############ Motion correction view #################
if hhh==0 | hhh==12 | hhh==18;
   figure;plot(tmotion(1:length(tmotion)),Result(3,1:length(tmotion))-mean(Result(3,1:length(tmotion))),'g',tmotion(1:length(tmotion)),Platvel(3,1:length(tmotion)),'r',tmotion(1:length(tmotion)),Velrot(3,1:length(tmotion))-mean(Velrot(3,1:length(tmotion))))
   axis([0 100 -2 2]);ylabel('w (m*s-1)');xlabel('time (s)');
   saveas(gcf,[way_rawdata,'\GOMECC\Final\Hare_', jd, hr ,'.fig']);
end;   
%#####################################################
end;

uvw(1,0*Nlenlic+1:1*Nlenlic)=uvw(1,0*Nlenlic+1:1*Nlenlic)+mean_cossog_1-mean(uvwplat(1,0*Nlenlic+1:1*Nlenlic));
uvw(1,1*Nlenlic+1:2*Nlenlic)=uvw(1,1*Nlenlic+1:2*Nlenlic)+mean_cossog_2-mean(uvwplat(1,1*Nlenlic+1:2*Nlenlic));
uvw(1,2*Nlenlic+1:3*Nlenlic)=uvw(1,2*Nlenlic+1:3*Nlenlic)+mean_cossog_3-mean(uvwplat(1,2*Nlenlic+1:3*Nlenlic));
uvw(1,3*Nlenlic+1:4*Nlenlic)=uvw(1,3*Nlenlic+1:4*Nlenlic)+mean_cossog_4-mean(uvwplat(1,3*Nlenlic+1:4*Nlenlic));
uvw(1,4*Nlenlic+1:5*Nlenlic)=uvw(1,4*Nlenlic+1:5*Nlenlic)+mean_cossog_5-mean(uvwplat(1,4*Nlenlic+1:5*Nlenlic));
uvw(1,5*Nlenlic+1:6*Nlenlic)=uvw(1,5*Nlenlic+1:6*Nlenlic)+mean_cossog_6-mean(uvwplat(1,5*Nlenlic+1:6*Nlenlic));

% uvw(2,0*Nlenlic+1:1*Nlenlic)=uvw(2,0*Nlenlic+1:1*Nlenlic)+mean_sinsog_1-mean(uvwplat(2,0*Nlenlic+1:1*Nlenlic));
% uvw(2,1*Nlenlic+1:2*Nlenlic)=uvw(2,1*Nlenlic+1:2*Nlenlic)+mean_sinsog_2-mean(uvwplat(2,1*Nlenlic+1:2*Nlenlic));
% uvw(2,2*Nlenlic+1:3*Nlenlic)=uvw(2,2*Nlenlic+1:3*Nlenlic)+mean_sinsog_3-mean(uvwplat(2,2*Nlenlic+1:3*Nlenlic));
% uvw(2,3*Nlenlic+1:4*Nlenlic)=uvw(2,3*Nlenlic+1:4*Nlenlic)+mean_sinsog_4-mean(uvwplat(2,3*Nlenlic+1:4*Nlenlic));
% uvw(2,4*Nlenlic+1:5*Nlenlic)=uvw(2,4*Nlenlic+1:5*Nlenlic)+mean_sinsog_5-mean(uvwplat(2,4*Nlenlic+1:5*Nlenlic));
% uvw(2,5*Nlenlic+1:6*Nlenlic)=uvw(2,5*Nlenlic+1:6*Nlenlic)+mean_sinsog_6-mean(uvwplat(2,5*Nlenlic+1:6*Nlenlic));

uvw(2,0*Nlenlic+1:1*Nlenlic)=uvw(2,0*Nlenlic+1:1*Nlenlic)-mean_sinsog_1-mean(uvwplat(2,0*Nlenlic+1:1*Nlenlic));
uvw(2,1*Nlenlic+1:2*Nlenlic)=uvw(2,1*Nlenlic+1:2*Nlenlic)-mean_sinsog_2-mean(uvwplat(2,1*Nlenlic+1:2*Nlenlic));
uvw(2,2*Nlenlic+1:3*Nlenlic)=uvw(2,2*Nlenlic+1:3*Nlenlic)-mean_sinsog_3-mean(uvwplat(2,2*Nlenlic+1:3*Nlenlic));
uvw(2,3*Nlenlic+1:4*Nlenlic)=uvw(2,3*Nlenlic+1:4*Nlenlic)-mean_sinsog_4-mean(uvwplat(2,3*Nlenlic+1:4*Nlenlic));
uvw(2,4*Nlenlic+1:5*Nlenlic)=uvw(2,4*Nlenlic+1:5*Nlenlic)-mean_sinsog_5-mean(uvwplat(2,4*Nlenlic+1:5*Nlenlic));
uvw(2,5*Nlenlic+1:6*Nlenlic)=uvw(2,5*Nlenlic+1:6*Nlenlic)-mean_sinsog_6-mean(uvwplat(2,5*Nlenlic+1:6*Nlenlic));

uvw(3,0*Nlenlic+1:1*Nlenlic)=uvw(3,0*Nlenlic+1:1*Nlenlic)-mean(uvwplat(3,0*Nlenlic+1:1*Nlenlic));
uvw(3,1*Nlenlic+1:2*Nlenlic)=uvw(3,1*Nlenlic+1:2*Nlenlic)-mean(uvwplat(3,1*Nlenlic+1:2*Nlenlic));
uvw(3,2*Nlenlic+1:3*Nlenlic)=uvw(3,2*Nlenlic+1:3*Nlenlic)-mean(uvwplat(3,2*Nlenlic+1:3*Nlenlic));
uvw(3,3*Nlenlic+1:4*Nlenlic)=uvw(3,3*Nlenlic+1:4*Nlenlic)-mean(uvwplat(3,3*Nlenlic+1:4*Nlenlic));
uvw(3,4*Nlenlic+1:5*Nlenlic)=uvw(3,4*Nlenlic+1:5*Nlenlic)-mean(uvwplat(3,4*Nlenlic+1:5*Nlenlic));
uvw(3,5*Nlenlic+1:6*Nlenlic)=uvw(3,5*Nlenlic+1:6*Nlenlic)-mean(uvwplat(3,5*Nlenlic+1:6*Nlenlic));

if ~clearit
%############ Motion correction view #################
if hhh==0 | hhh==12 | hhh==18,
figure;plot(tmotion,uvw(3,:)-mean(uvw(3,:)),'g',tmotion,uvwplat(3,:),'r',tmotion,AAAA(3,:)-mean(AAAA(3,:)))
axis([0 100 -2 2]);ylabel('w (m*s-1)');xlabel('time (s)')
saveas(gcf,[way_rawdata,'\GOMECC\Final\Edson_', jd, hr ,'.fig']);
end;
%######################################################
end;

[Ustream_1,Vstream_1,Wstream_1,alpha_1,beta_1]=rotcorr(Result(1,0*Nlenlic+1:1*Nlenlic),Result(2,0*Nlenlic+1:1*Nlenlic),Result(3,0*Nlenlic+1:1*Nlenlic));           %converts [u,v,w] to [Ustream, Vstream, Wstream] for ten minutes
[Ustream_2,Vstream_2,Wstream_2,alpha_2,beta_2]=rotcorr(Result(1,1*Nlenlic+1:2*Nlenlic),Result(2,1*Nlenlic+1:2*Nlenlic),Result(3,1*Nlenlic+1:2*Nlenlic));
[Ustream_3,Vstream_3,Wstream_3,alpha_3,beta_3]=rotcorr(Result(1,2*Nlenlic+1:3*Nlenlic),Result(2,2*Nlenlic+1:3*Nlenlic),Result(3,2*Nlenlic+1:3*Nlenlic));
[Ustream_4,Vstream_4,Wstream_4,alpha_4,beta_4]=rotcorr(Result(1,3*Nlenlic+1:4*Nlenlic),Result(2,3*Nlenlic+1:4*Nlenlic),Result(3,3*Nlenlic+1:4*Nlenlic));
[Ustream_5,Vstream_5,Wstream_5,alpha_5,beta_5]=rotcorr(Result(1,4*Nlenlic+1:5*Nlenlic),Result(2,4*Nlenlic+1:5*Nlenlic),Result(3,4*Nlenlic+1:5*Nlenlic));
[Ustream_6,Vstream_6,Wstream_6,alpha_6,beta_6]=rotcorr(Result(1,5*Nlenlic+1:6*Nlenlic),Result(2,5*Nlenlic+1:6*Nlenlic),Result(3,5*Nlenlic+1:6*Nlenlic));


Tsonic_1=Tsonic(0*Nlenlic+1:1*Nlenlic);
Tsonic_2=Tsonic(1*Nlenlic+1:2*Nlenlic);
Tsonic_3=Tsonic(2*Nlenlic+1:3*Nlenlic);
Tsonic_4=Tsonic(3*Nlenlic+1:4*Nlenlic);
Tsonic_5=Tsonic(4*Nlenlic+1:5*Nlenlic);
Tsonic_6=Tsonic(5*Nlenlic+1:6*Nlenlic);

%abshum_1=abshum(00001:12000);
%abshum_2=abshum(12001:24000);
%abshum_3=abshum(24001:36000);
%abshum_4=abshum(36001:48000);
%abshum_5=abshum(48001:60000);
%abshum_6=abshum(60001:72000);

%std_abshum_1=std(abshum_1);
%std_abshum_2=std(abshum_2);
%std_abshum_3=std(abshum_3);
%std_abshum_4=std(abshum_4);
%std_abshum_5=std(abshum_5);
%std_abshum_6=std(abshum_6);

% Licor_CO2=(1:Lenlic).*0;    %Remove when Licor works
% Licor_H2O=(1:Lenlic).*0;    %Remove when Licor works
% Licor_Temp=(1:Lenlic).*0;    %Remove when Licor works
% Licor_Press=(1:Lenlic).*0;    %Remove when Licor works

Licor_CO2_1=Licor_CO2(0*Nlenlic+1:1*Nlenlic);
Licor_CO2_2=Licor_CO2(1*Nlenlic+1:2*Nlenlic);
Licor_CO2_3=Licor_CO2(2*Nlenlic+1:3*Nlenlic);
Licor_CO2_4=Licor_CO2(3*Nlenlic+1:4*Nlenlic);
Licor_CO2_5=Licor_CO2(4*Nlenlic+1:5*Nlenlic);
Licor_CO2_6=Licor_CO2(5*Nlenlic+1:length(Licor_CO2));

Licor_H2O_1=Licor_H2O(0*Nlenlic+1:1*Nlenlic);
Licor_H2O_2=Licor_H2O(1*Nlenlic+1:2*Nlenlic);
Licor_H2O_3=Licor_H2O(2*Nlenlic+1:3*Nlenlic);
Licor_H2O_4=Licor_H2O(3*Nlenlic+1:4*Nlenlic);
Licor_H2O_5=Licor_H2O(4*Nlenlic+1:5*Nlenlic);
Licor_H2O_6=Licor_H2O(5*Nlenlic+1:length(Licor_H2O));

if Licor_Temp(1)>40,
   Licor_Temp(1)=25;
end;
if Licor_Temp(1)<10,
   Licor_Temp(1)=25;
end;

if Licor_Press(1)>120,
   Licor_Press(1)=100.5;
end;
if Licor_Press(1)<90,
   Licor_Press(1)=100.5;
end;

for ii=2:length(Licor_Temp);
    if Licor_Temp(ii)>40,
        Licor_Temp(ii)=Licor_Temp(ii-1);
    end;
    if Licor_Temp(ii)<10,
        Licor_Temp(ii)=Licor_Temp(ii-1);
    end;
end;

for ii=2:length(Licor_Press);
    if Licor_Press(ii)>40,
        Licor_Press(ii)=Licor_Press(ii-1);
    end;
    if Licor_Press(ii)<10,
        Licor_Press(ii)=Licor_Press(ii-1);
    end;
end;

Licor_Temp_1=Licor_Temp(0*Nlenlic+1:1*Nlenlic);
Licor_Temp_2=Licor_Temp(1*Nlenlic+1:2*Nlenlic);
Licor_Temp_3=Licor_Temp(2*Nlenlic+1:3*Nlenlic);
Licor_Temp_4=Licor_Temp(3*Nlenlic+1:4*Nlenlic);
Licor_Temp_5=Licor_Temp(4*Nlenlic+1:5*Nlenlic);
Licor_Temp_6=Licor_Temp(5*Nlenlic+1:length(Licor_Temp));

Licor_Press_1=Licor_Press(0*Nlenlic+1:1*Nlenlic);
Licor_Press_2=Licor_Press(1*Nlenlic+1:2*Nlenlic);
Licor_Press_3=Licor_Press(2*Nlenlic+1:3*Nlenlic);
Licor_Press_4=Licor_Press(3*Nlenlic+1:4*Nlenlic);
Licor_Press_5=Licor_Press(4*Nlenlic+1:5*Nlenlic);
Licor_Press_6=Licor_Press(5*Nlenlic+1:length(Licor_Press));

% speed=[uvw ;uvwplat];
% speed=speed';
%       save([way_rawdata ,'\GOMECC\Final\sonic_corr_',jd, hr,'.txt'], 'speed','-ascii');

[Uj_1,Vj_1,Wj_1,alphaj_1,betaj_1]=rotcorr(uvw(1,0*Nlenlic+1:1*Nlenlic),uvw(2,0*Nlenlic+1:1*Nlenlic),uvw(3,0*Nlenlic+1:1*Nlenlic));
[Uj_2,Vj_2,Wj_2,alphaj_2,betaj_2]=rotcorr(uvw(1,1*Nlenlic+1:2*Nlenlic),uvw(2,1*Nlenlic+1:2*Nlenlic),uvw(3,1*Nlenlic+1:2*Nlenlic));
[Uj_3,Vj_3,Wj_3,alphaj_3,betaj_3]=rotcorr(uvw(1,2*Nlenlic+1:3*Nlenlic),uvw(2,2*Nlenlic+1:3*Nlenlic),uvw(3,2*Nlenlic+1:3*Nlenlic));
[Uj_4,Vj_4,Wj_4,alphaj_4,betaj_4]=rotcorr(uvw(1,3*Nlenlic+1:4*Nlenlic),uvw(2,3*Nlenlic+1:4*Nlenlic),uvw(3,3*Nlenlic+1:4*Nlenlic));
[Uj_5,Vj_5,Wj_5,alphaj_5,betaj_5]=rotcorr(uvw(1,4*Nlenlic+1:5*Nlenlic),uvw(2,4*Nlenlic+1:5*Nlenlic),uvw(3,4*Nlenlic+1:5*Nlenlic));
[Uj_6,Vj_6,Wj_6,alphaj_6,betaj_6]=rotcorr(uvw(1,5*Nlenlic+1:6*Nlenlic),uvw(2,5*Nlenlic+1:6*Nlenlic),uvw(3,5*Nlenlic+1:6*Nlenlic));


Ucat=[Uj_1 Uj_2 Uj_3 Uj_4 Uj_5 Uj_6];
Ucat=Ucat(1:36000);

Vcat=[Vj_1 Vj_2 Vj_3 Vj_4 Vj_5 Vj_6];
Vcat=Vcat(1:36000);

Wcat=[Wj_1 Wj_2 Wj_3 Wj_4 Wj_5 Wj_6];
Wcat=Wcat(1:36000);


   speed=[Ucat; Vcat; Wcat];
   %speed=speed';
   %save([way_rawdata ,'\GOMECC\Final\sonic+O3\sonic_corr_',jd, hr,'.txt'], 'speed','-ascii');
save([way_rawdata ,'\GOMECC\Wcat_',jd, hr,'.txt'], 'speed','-ascii');


Ubarj_1=mean(Uj_1);
Ubarj_2=mean(Uj_2);
Ubarj_3=mean(Uj_3);
Ubarj_4=mean(Uj_4);
Ubarj_5=mean(Uj_5);
Ubarj_6=mean(Uj_6);

Tbar_1=mean(Tvais(1:10));
Tbar_2=mean(Tvais(11:20));
Tbar_3=mean(Tvais(21:30));
Tbar_4=mean(Tvais(31:40));
Tbar_5=mean(Tvais(41:50));
Tbar_6=mean(Tvais(51:60));

Tship_1=mean(imta(1:10));
Tship_2=mean(imta(11:20));
Tship_3=mean(imta(21:30));
Tship_4=mean(imta(31:40));
Tship_5=mean(imta(41:50));
Tship_6=mean(imta(51:60));

rho_1=P*100./(273.15+Tbar_1)./287;
rho_2=P*100./(273.15+Tbar_2)./287;
rho_3=P*100./(273.15+Tbar_3)./287;
rho_4=P*100./(273.15+Tbar_4)./287;
rho_5=P*100./(273.15+Tbar_5)./287;
rho_6=P*100./(273.15+Tbar_6)./287;

wu_1=cov(Ustream_1,Wstream_1);
wu_2=cov(Ustream_2,Wstream_2);
wu_3=cov(Ustream_3,Wstream_3);
wu_4=cov(Ustream_4,Wstream_4);
wu_5=cov(Ustream_5,Wstream_5);
wu_6=cov(Ustream_6,Wstream_6);

wv_1=cov(Vstream_1,Wstream_1);
wv_2=cov(Vstream_2,Wstream_2);
wv_3=cov(Vstream_3,Wstream_3);
wv_4=cov(Vstream_4,Wstream_4);
wv_5=cov(Vstream_5,Wstream_5);
wv_6=cov(Vstream_6,Wstream_6);

wt_1=cov(Wstream_1,Tsonic_1);
wt_2=cov(Wstream_2,Tsonic_2);
wt_3=cov(Wstream_3,Tsonic_3);
wt_4=cov(Wstream_4,Tsonic_4);
wt_5=cov(Wstream_5,Tsonic_5);
wt_6=cov(Wstream_6,Tsonic_6);

Ustr_1=sqrt(mean(U(0*Nlenlic+1:1*Nlenlic)).^2 + mean(V(0*Nlenlic+1:1*Nlenlic)).^2 + mean(W(0*Nlenlic+1:1*Nlenlic)).^2);
Ustr_2=sqrt(mean(U(1*Nlenlic+1:2*Nlenlic)).^2 + mean(V(1*Nlenlic+1:2*Nlenlic)).^2 + mean(W(1*Nlenlic+1:2*Nlenlic)).^2);
Ustr_3=sqrt(mean(U(2*Nlenlic+1:3*Nlenlic)).^2 + mean(V(2*Nlenlic+1:3*Nlenlic)).^2 + mean(W(2*Nlenlic+1:3*Nlenlic)).^2);
Ustr_4=sqrt(mean(U(3*Nlenlic+1:4*Nlenlic)).^2 + mean(V(3*Nlenlic+1:4*Nlenlic)).^2 + mean(W(3*Nlenlic+1:4*Nlenlic)).^2);
Ustr_5=sqrt(mean(U(4*Nlenlic+1:5*Nlenlic)).^2 + mean(V(4*Nlenlic+1:5*Nlenlic)).^2 + mean(W(4*Nlenlic+1:5*Nlenlic)).^2);
Ustr_6=sqrt(mean(U(5*Nlenlic+1:6*Nlenlic)).^2 + mean(V(5*Nlenlic+1:6*Nlenlic)).^2 + mean(W(5*Nlenlic+1:6*Nlenlic)).^2);

wuj_1=cov(Uj_1,Wj_1);
wuj_2=cov(Uj_2,Wj_2);
wuj_3=cov(Uj_3,Wj_3);
wuj_4=cov(Uj_4,Wj_4);
wuj_5=cov(Uj_5,Wj_5);
wuj_6=cov(Uj_6,Wj_6);

wvj_1=cov(Vj_1,Wj_1);
wvj_2=cov(Vj_2,Wj_2);
wvj_3=cov(Vj_3,Wj_3);
wvj_4=cov(Vj_4,Wj_4);
wvj_5=cov(Vj_5,Wj_5);
wvj_6=cov(Vj_6,Wj_6);

wtj_1=cov(Wj_1,Tsonic_1);
wtj_2=cov(Wj_2,Tsonic_2);
wtj_3=cov(Wj_3,Tsonic_3);
wtj_4=cov(Wj_4,Tsonic_4);
wtj_5=cov(Wj_5,Tsonic_5);
wtj_6=cov(Wj_6,Tsonic_6);

Licor_w_1=Wj_1;
Licor_w_2=Wj_2;
Licor_w_3=Wj_3;
Licor_w_4=Wj_4;
Licor_w_5=Wj_5;
Licor_w_6=Wj_6;

I_1=Licor_CO2_1.*0;
I_2=Licor_CO2_2.*0;
I_3=Licor_CO2_3.*0;
I_4=Licor_CO2_4.*0;
I_5=Licor_CO2_5.*0;
I_6=Licor_CO2_6.*0;

for ii=1:length(Licor_CO2_1);
    if Licor_CO2_1(ii)>90000;
        I_1(ii)=1;
    end;
end;
for ii=1:length(Licor_CO2_2);
    if Licor_CO2_2(ii)>90000;
        I_2(ii)=1;
    end;
end;
for ii=1:length(Licor_CO2_3);
    if Licor_CO2_3(ii)>90000;
        I_3(ii)=1;
    end;
end;
for ii=1:length(Licor_CO2_4);
    if Licor_CO2_4(ii)>90000;
        I_4(ii)=1;
    end;
end;
for ii=1:length(Licor_CO2_5);
    if Licor_CO2_5(ii)>90000;
        I_5(ii)=1;
    end;
end;
for ii=1:length(Licor_CO2_6);
    if Licor_CO2_6(ii)>90000;
        I_6(ii)=1;
    end;
end;

Licor_CO2_1=Licor_CO2_1.*8.314.*(mean(Licor_Temp_1)+273.15)./mean(Licor_Press_1);   %convert from mmol/m3 to umol/mol
Licor_CO2_2=Licor_CO2_2.*8.314.*(mean(Licor_Temp_2)+273.15)./mean(Licor_Press_2);   %convert from mmol/m3 to umol/mol
Licor_CO2_3=Licor_CO2_3.*8.314.*(mean(Licor_Temp_3)+273.15)./mean(Licor_Press_3);   %convert from mmol/m3 to umol/mol
Licor_CO2_4=Licor_CO2_4.*8.314.*(mean(Licor_Temp_4)+273.15)./mean(Licor_Press_4);   %convert from mmol/m3 to umol/mol
Licor_CO2_5=Licor_CO2_5.*8.314.*(mean(Licor_Temp_5)+273.15)./mean(Licor_Press_5);   %convert from mmol/m3 to umol/mol
Licor_CO2_6=Licor_CO2_6.*8.314.*(mean(Licor_Temp_6)+273.15)./mean(Licor_Press_6);   %convert from mmol/m3 to umol/mol

Licor_H2O_1=Licor_H2O_1.*18./1000./rho_1;         %convert from mmol/m3 to g/kg
Licor_H2O_2=Licor_H2O_2.*18./1000./rho_2;         %convert from mmol/m3 to g/kg
Licor_H2O_3=Licor_H2O_3.*18./1000./rho_3;         %convert from mmol/m3 to g/kg
Licor_H2O_4=Licor_H2O_4.*18./1000./rho_4;         %convert from mmol/m3 to g/kg
Licor_H2O_5=Licor_H2O_5.*18./1000./rho_5;         %convert from mmol/m3 to g/kg
Licor_H2O_6=Licor_H2O_6.*18./1000./rho_6;         %convert from mmol/m3 to g/kg

Licor_H2O=[Licor_H2O_1 Licor_H2O_2 Licor_H2O_3 Licor_H2O_4 Licor_H2O_5 Licor_H2O_6];
save([way_rawdata, '\GOMECC\Licor_H2O' , jd,hr '.txt'], 'Licor_H2O', '-ascii')

%*  changed to /, cwf 12/04
Licor_CO2_1=Licor_CO2_1./(1-Licor_H2O_1./18.*1000.*rho_1.*8.314.*(mean(Licor_Temp_1)+273.15)./1000./mean(Licor_Press_1)./1000);   %dilution correction
Licor_CO2_2=Licor_CO2_2./(1-Licor_H2O_2./18.*1000.*rho_2.*8.314.*(mean(Licor_Temp_2)+273.15)./1000./mean(Licor_Press_1)./1000);   %dilution correction
Licor_CO2_3=Licor_CO2_3./(1-Licor_H2O_3./18.*1000.*rho_3.*8.314.*(mean(Licor_Temp_3)+273.15)./1000./mean(Licor_Press_1)./1000);   %dilution correction
Licor_CO2_4=Licor_CO2_4./(1-Licor_H2O_4./18.*1000.*rho_4.*8.314.*(mean(Licor_Temp_4)+273.15)./1000./mean(Licor_Press_1)./1000);   %dilution correction
Licor_CO2_5=Licor_CO2_5./(1-Licor_H2O_5./18.*1000.*rho_5.*8.314.*(mean(Licor_Temp_5)+273.15)./1000./mean(Licor_Press_1)./1000);   %dilution correction
Licor_CO2_6=Licor_CO2_6./(1-Licor_H2O_6./18.*1000.*rho_6.*8.314.*(mean(Licor_Temp_6)+273.15)./1000./mean(Licor_Press_1)./1000);   %dilution correction

LicCO2_1=[];LicH2O_1=[];WCO2_1=[];
LicCO2_2=[];LicH2O_2=[];WCO2_2=[];
LicCO2_3=[];LicH2O_3=[];WCO2_3=[];
LicCO2_4=[];LicH2O_4=[];WCO2_4=[];
LicCO2_5=[];LicH2O_5=[];WCO2_5=[];
LicCO2_6=[];LicH2O_6=[];WCO2_6=[];

for ii=1:length(Licor_CO2_1);
    if I_1(ii)==0,
        LicCO2_1=[LicCO2_1 Licor_CO2_1(ii)];
        LicH2O_1=[LicH2O_1 Licor_H2O_1(ii)];
        WCO2_1=[WCO2_1 Licor_w_1(ii)];
    end;
end;
for ii=1:length(Licor_CO2_2);
    if I_2(ii)==0,
        LicCO2_2=[LicCO2_2 Licor_CO2_2(ii)];
        LicH2O_2=[LicH2O_2 Licor_H2O_2(ii)];
        WCO2_2=[WCO2_2 Licor_w_2(ii)];
    end;
end;
for ii=1:length(Licor_CO2_3);
    if I_3(ii)==0,
        LicCO2_3=[LicCO2_3 Licor_CO2_3(ii)];
        LicH2O_3=[LicH2O_3 Licor_H2O_3(ii)];
        WCO2_3=[WCO2_3 Licor_w_3(ii)];
    end;
end;
for ii=1:length(Licor_CO2_4);
    if I_4(ii)==0,
        LicCO2_4=[LicCO2_4 Licor_CO2_4(ii)];
        LicH2O_4=[LicH2O_4 Licor_H2O_4(ii)];
        WCO2_4=[WCO2_4 Licor_w_4(ii)];
    end;
end;
for ii=1:length(Licor_CO2_5);
    if I_5(ii)==0,
        LicCO2_5=[LicCO2_5 Licor_CO2_5(ii)];
        LicH2O_5=[LicH2O_5 Licor_H2O_5(ii)];
        WCO2_5=[WCO2_5 Licor_w_5(ii)];
    end;
end;
for ii=1:length(Licor_CO2_6);
    if I_6(ii)==0,
        LicCO2_6=[LicCO2_6 Licor_CO2_6(ii)];
        LicH2O_6=[LicH2O_6 Licor_H2O_6(ii)];
        WCO2_6=[WCO2_6 Licor_w_6(ii)];
    end;
end;

%LicTemp_1=[LicTemp_1 Licor_Temp_1(ii)];
%LicPress_1=[LicPress_1 Licor_Press_1(ii)];

if length(LicCO2_1)==0,
    LicCO2_1=[1 1];
    LicH2O_1=[1 1];
    WCO2_1=[1 1];
end;
if length(LicCO2_2)==0,
    LicCO2_2=[1 1];
    LicH2O_2=[1 1];
    WCO2_2=[1 1];
end;
if length(LicCO2_3)==0,
    LicCO2_3=[1 1];
    LicH2O_3=[1 1];
    WCO2_3=[1 1];
end;
if length(LicCO2_4)==0,
    LicCO2_4=[1 1];
    LicH2O_4=[1 1];
    WCO2_4=[1 1];
end;
if length(LicCO2_5)==0,
    LicCO2_5=[1 1];
    LicH2O_5=[1 1];
    WCO2_5=[1 1];
end;
if length(LicCO2_6)==0,
    LicCO2_6=[1 1];
    LicH2O_6=[1 1];
    WCO2_6=[1 1];
end;

wlcj_1=cov(WCO2_1,LicCO2_1);wlcj_1=wlcj_1(1,2);
wlcj_2=cov(WCO2_2,LicCO2_2);wlcj_2=wlcj_2(1,2);
wlcj_3=cov(WCO2_3,LicCO2_3);wlcj_3=wlcj_3(1,2);
wlcj_4=cov(WCO2_4,LicCO2_4);wlcj_4=wlcj_4(1,2);
wlcj_5=cov(WCO2_5,LicCO2_5);wlcj_5=wlcj_5(1,2);
wlcj_6=cov(WCO2_6,LicCO2_6);wlcj_6=wlcj_6(1,2);

wlhj_1=cov(WCO2_1,LicH2O_1);wlhj_1=wlhj_1(1,2);
wlhj_2=cov(WCO2_2,LicH2O_2);wlhj_2=wlhj_2(1,2);
wlhj_3=cov(WCO2_3,LicH2O_3);wlhj_3=wlhj_3(1,2);
wlhj_4=cov(WCO2_4,LicH2O_4);wlhj_4=wlhj_4(1,2);
wlhj_5=cov(WCO2_5,LicH2O_5);wlhj_5=wlhj_5(1,2);
wlhj_6=cov(WCO2_6,LicH2O_6);wlhj_6=wlhj_6(1,2);

disp('Licor H2O +')
for iii=0:20;
   A1=cov(WCO2_1(1:length(WCO2_1)-iii),LicH2O_1(1+iii:length(LicH2O_1)));
   A2=cov(WCO2_2(1:length(WCO2_2)-iii),LicH2O_2(1+iii:length(LicH2O_2)));
   A3=cov(WCO2_3(1:length(WCO2_3)-iii),LicH2O_3(1+iii:length(LicH2O_3)));
   A4=cov(WCO2_4(1:length(WCO2_4)-iii),LicH2O_4(1+iii:length(LicH2O_4)));
   A5=cov(WCO2_5(1:length(WCO2_5)-iii),LicH2O_5(1+iii:length(LicH2O_5)));
   A6=cov(WCO2_6(1:length(WCO2_6)-iii),LicH2O_6(1+iii:length(LicH2O_6)));
   if length(A1)<2,
       A1=[1 1
           1 1];
   end;
   if length(A2)<2,
       A2=[1 1
           1 1];
   end;
   if length(A3)<2,
       A3=[1 1
           1 1];
   end;
   if length(A4)<2,
       A4=[1 1
           1 1];
   end;
   if length(A5)<2,
       A5=[1 1
           1 1];
   end;
   if length(A6)<2,
       A6=[1 1
           1 1];
   end;
%   [wlhj_1 wlhj_2 wlhj_3 wlhj_4 wlhj_5 wlhj_6]
   B1=A1(1,2);
   B2=A2(1,2);
   B3=A3(1,2);
   B4=A4(1,2);
   B5=A5(1,2);
   B6=A6(1,2);
   if B1>wlhj_1,
      wlhj_1=B1;
   end;
   if B2>wlhj_2,
      wlhj_2=B2;
   end;
   if B3>wlhj_3,
      wlhj_3=B3;
   end;
   if B4>wlhj_4,
      wlhj_4=B4;
  end;
   if B5>wlhj_5,
      wlhj_5=B5;
   end;
   if B6>wlhj_6,
      wlhj_6=B6;
   end;
end;

disp('Licor H2O -')
for iii=-20:0;
   A1=cov(WCO2_1(1-iii:length(WCO2_1)),LicH2O_1(1:length(LicH2O_1)+iii));
   A2=cov(WCO2_2(1-iii:length(WCO2_2)),LicH2O_2(1:length(LicH2O_2)+iii));
   A3=cov(WCO2_3(1-iii:length(WCO2_3)),LicH2O_3(1:length(LicH2O_3)+iii));
   A4=cov(WCO2_4(1-iii:length(WCO2_4)),LicH2O_4(1:length(LicH2O_4)+iii));
   A5=cov(WCO2_5(1-iii:length(WCO2_5)),LicH2O_5(1:length(LicH2O_5)+iii));
   A6=cov(WCO2_6(1-iii:length(WCO2_6)),LicH2O_6(1:length(LicH2O_6)+iii));
   if length(A1)<2,
       A1=[1 1
           1 1];
   end;
   if length(A2)<2,
       A2=[1 1
           1 1];
   end;
   if length(A3)<2,
       A3=[1 1
           1 1];
   end;
   if length(A4)<2,
       A4=[1 1
           1 1];
   end;
   if length(A5)<2,
       A5=[1 1
           1 1];
   end;
   if length(A6)<2,
       A6=[1 1
           1 1];
   end;
   B1=A1(1,2);
   B2=A2(1,2);
   B3=A3(1,2);
   B4=A4(1,2);
   B5=A5(1,2);
   B6=A6(1,2);
%   [wlhj_1 wlhj_2 wlhj_3 wlhj_4 wlhj_5 wlhj_6]
   if B1>wlhj_1,
      wlhj_1=B1;
   end;
   if B2>wlhj_2,
      wlhj_2=B2;
   end;
   if B3>wlhj_3,
      wlhj_3=B3;
   end;
   if B4>wlhj_4,
      wlhj_4=B4;
  end;
   if B5>wlhj_5,
      wlhj_5=B5;
   end;
   if B6>wlhj_6,
      wlhj_6=B6;
   end;
end;

wlcj_1=0;
wlcj_2=0;
wlcj_3=0;
wlcj_4=0;
wlcj_5=0;
wlcj_6=0;

disp('Licor CO2 +')
for iii=0:20;
   A1=cov(WCO2_1(1:length(WCO2_1)-iii),LicCO2_1(1+iii:length(LicCO2_1)));
   A2=cov(WCO2_2(1:length(WCO2_2)-iii),LicCO2_2(1+iii:length(LicCO2_2)));
   A3=cov(WCO2_3(1:length(WCO2_3)-iii),LicCO2_3(1+iii:length(LicCO2_3)));
   A4=cov(WCO2_4(1:length(WCO2_4)-iii),LicCO2_4(1+iii:length(LicCO2_4)));
   A5=cov(WCO2_5(1:length(WCO2_5)-iii),LicCO2_5(1+iii:length(LicCO2_5)));
   A6=cov(WCO2_6(1:length(WCO2_6)-iii),LicCO2_6(1+iii:length(LicCO2_6)));
   if length(A1)<2,
       A1=[1 1
           1 1];
   end;
   if length(A2)<2,
       A2=[1 1
           1 1];
   end;
   if length(A3)<2,
       A3=[1 1
           1 1];
   end;
   if length(A4)<2,
       A4=[1 1
           1 1];
   end;
   if length(A5)<2,
       A5=[1 1
           1 1];
   end;
   if length(A6)<2,
       A6=[1 1
           1 1];
   end;
   B1=A1(1,2);
   B2=A2(1,2);
   B3=A3(1,2);
   B4=A4(1,2);
   B5=A5(1,2);
   B6=A6(1,2);
%   [wlcj_1 wlcj_2 wlcj_3 wlcj_4 wlcj_5 wlcj_6]
   if B1>wlcj_1,
      wlcj_1=B1;
   end;
   if B2>wlcj_2,
      wlcj_2=B2;
   end;
   if B3>wlcj_3,
      wlcj_3=B3;
   end;
   if B4>wlcj_4,
      wlcj_4=B4;
   end;
   if B5>wlcj_5,
      wlcj_5=B5;
   end;
   if B6>wlcj_6,
      wlcj_6=B6;
   end;
end;

disp('Licor CO2 -')
for iii=-20:0;
   A1=cov(WCO2_1(1-iii:length(WCO2_1)),LicCO2_1(1:length(LicCO2_1)+iii));
   A2=cov(WCO2_2(1-iii:length(WCO2_2)),LicCO2_2(1:length(LicCO2_2)+iii));
   A3=cov(WCO2_3(1-iii:length(WCO2_3)),LicCO2_3(1:length(LicCO2_3)+iii));
   A4=cov(WCO2_4(1-iii:length(WCO2_4)),LicCO2_4(1:length(LicCO2_4)+iii));
   A5=cov(WCO2_5(1-iii:length(WCO2_5)),LicCO2_5(1:length(LicCO2_5)+iii));
   A6=cov(WCO2_6(1-iii:length(WCO2_6)),LicCO2_6(1:length(LicCO2_6)+iii));
   if length(A1)<2,
       A1=[1 1
           1 1];
   end;
   if length(A2)<2,
       A2=[1 1
           1 1];
   end;
   if length(A3)<2,
       A3=[1 1
           1 1];
   end;
   if length(A4)<2,
       A4=[1 1
           1 1];
   end;
   if length(A5)<2,
       A5=[1 1
           1 1];
   end;
   if length(A6)<2,
       A6=[1 1
           1 1];
   end;
   B1=A1(1,2);
   B2=A2(1,2);
   B3=A3(1,2);
   B4=A4(1,2);
   B5=A5(1,2);
   B6=A6(1,2);
%   [wlcj_1 wlcj_2 wlcj_3 wlcj_4 wlcj_5 wlcj_6]
   if B1>wlcj_1,
      wlcj_1=B1;
   end;
   if B2>wlcj_2,
      wlcj_2=B2;
   end;
   if B3>wlcj_3,
      wlcj_3=B3;
   end;
   if B4>wlcj_4,
      wlcj_4=B4;
   end;
   if B5>wlcj_5,
      wlcj_5=B5;
   end;
   if B6>wlcj_6,
      wlcj_6=B6;
   end;
end;

wohj_1=0;
wohj_2=0;
wohj_3=0;
wohj_4=0;
wohj_5=0;
wohj_6=0;

wlcj_1=wlcj_1.*28.96.*365.24.*24.*60.*60.*1e-6;   %units are mol/m2/year
wlcj_2=wlcj_2.*28.96.*365.24.*24.*60.*60.*1e-6;   %units are mol/m2/year
wlcj_3=wlcj_3.*28.96.*365.24.*24.*60.*60.*1e-6;   %units are mol/m2/year
wlcj_4=wlcj_4.*28.96.*365.24.*24.*60.*60.*1e-6;   %units are mol/m2/year
wlcj_5=wlcj_5.*28.96.*365.24.*24.*60.*60.*1e-6;   %units are mol/m2/year
wlcj_6=wlcj_6.*28.96.*365.24.*24.*60.*60.*1e-6;   %units are mol/m2/year

wahj_1=0;
wahj_2=0;
wahj_3=0;
wahj_4=0;
wahj_5=0;
wahj_6=0;



mean_wind_1=mean(imum(01:10));std_wind_1=std(imum(01:10));max_wind_1=max(imum(01:10));min_wind_1=min(imum(01:10));
mean_wind_2=mean(imum(11:20));std_wind_2=std(imum(11:20));max_wind_2=max(imum(11:20));min_wind_2=min(imum(11:20));
mean_wind_3=mean(imum(21:30));std_wind_3=std(imum(21:30));max_wind_3=max(imum(21:30));min_wind_3=min(imum(21:30));
mean_wind_4=mean(imum(31:40));std_wind_4=std(imum(31:40));max_wind_4=max(imum(31:40));min_wind_4=min(imum(31:40));
mean_wind_5=mean(imum(41:50));std_wind_5=std(imum(41:50));max_wind_5=max(imum(41:50));min_wind_5=min(imum(41:50));
mean_wind_6=mean(imum(51:60));std_wind_6=std(imum(51:60));max_wind_6=max(imum(51:60));min_wind_6=min(imum(51:60));

mean_dirC_1=median(cos(imdm(01:10)./180.*pi));
mean_dirC_2=median(cos(imdm(11:20)./180.*pi));
mean_dirC_3=median(cos(imdm(21:30)./180.*pi));
mean_dirC_4=median(cos(imdm(31:40)./180.*pi));
mean_dirC_5=median(cos(imdm(41:50)./180.*pi));
mean_dirC_6=median(cos(imdm(51:60)./180.*pi));

mean_dirS_1=median(sin(imdm(01:10)./180.*pi));
mean_dirS_2=median(sin(imdm(11:20)./180.*pi));
mean_dirS_3=median(sin(imdm(21:30)./180.*pi));
mean_dirS_4=median(sin(imdm(31:40)./180.*pi));
mean_dirS_5=median(sin(imdm(41:50)./180.*pi));
mean_dirS_6=median(sin(imdm(51:60)./180.*pi));

mean_dir_1=mod(atan2(mean_dirS_1,mean_dirC_1).*180./pi+360,360);		%
mean_dir_2=mod(atan2(mean_dirS_2,mean_dirC_2).*180./pi+360,360);		%
mean_dir_3=mod(atan2(mean_dirS_3,mean_dirC_3).*180./pi+360,360);		%
mean_dir_4=mod(atan2(mean_dirS_4,mean_dirC_4).*180./pi+360,360);		%
mean_dir_5=mod(atan2(mean_dirS_5,mean_dirC_5).*180./pi+360,360);		%reported IMET true wind direction
mean_dir_6=mod(atan2(mean_dirS_6,mean_dirC_6).*180./pi+360,360);		%

if clearit
clear mean_dirC_1 mean_dirC_2 mean_dirC_3 mean_dirC_4 mean_dirC_5 mean_dirC_6 
clear mean_dirS_1 mean_dirS_2 mean_dirS_3 mean_dirS_4 mean_dirS_5 mean_dirS_6 
end;

std_dir_1=std(unwrap(imdm(01:10)));max_dir_1=max(unwrap(imdm(01:10)));min_dir_1=min(unwrap(imdm(01:10)));
std_dir_2=std(unwrap(imdm(11:20)));max_dir_2=max(unwrap(imdm(11:20)));min_dir_2=min(unwrap(imdm(11:20)));
std_dir_3=std(unwrap(imdm(21:30)));max_dir_3=max(unwrap(imdm(21:30)));min_dir_3=min(unwrap(imdm(21:30)));
std_dir_4=std(unwrap(imdm(31:40)));max_dir_4=max(unwrap(imdm(31:40)));min_dir_4=min(unwrap(imdm(31:40)));
std_dir_5=std(unwrap(imdm(41:50)));max_dir_5=max(unwrap(imdm(41:50)));min_dir_5=min(unwrap(imdm(41:50)));
std_dir_6=std(unwrap(imdm(51:60)));max_dir_6=max(unwrap(imdm(51:60)));min_dir_6=min(unwrap(imdm(51:60)));

Reldir_1=atan2(mean(V(0*Nlenlic+1:1*Nlenlic)),-mean(U(0*Nlenlic+1:1*Nlenlic)))*180/pi;
Reldir_2=atan2(mean(V(1*Nlenlic+1:2*Nlenlic)),-mean(U(1*Nlenlic+1:2*Nlenlic)))*180/pi;
Reldir_3=atan2(mean(V(2*Nlenlic+1:3*Nlenlic)),-mean(U(2*Nlenlic+1:3*Nlenlic)))*180/pi;
Reldir_4=atan2(mean(V(3*Nlenlic+1:4*Nlenlic)),-mean(U(3*Nlenlic+1:4*Nlenlic)))*180/pi;
Reldir_5=atan2(mean(V(4*Nlenlic+1:5*Nlenlic)),-mean(U(4*Nlenlic+1:5*Nlenlic)))*180/pi;
Reldir_6=atan2(mean(V(5*Nlenlic+1:6*Nlenlic)),-mean(U(5*Nlenlic+1:6*Nlenlic)))*180/pi;

U=U./0.95;              %Flow Distortion correction
V=V./1.15;

%SS=[(1:900).*0+sodecx(1) sodecx (1:900).*0+sodecx(length(sodecx))];         %resample sodecx
%TT=[(1:900).*0+sodecy(1) sodecy (1:900).*0+sodecy(length(sodecy))];         %resample sodecy

%YY=resample(SS,144000,3600);YY=YY(36001:108000);         %resample sodecx
%ZZ=resample(TT,144000,3600);ZZ=ZZ(36001:108000);         %resample sodecy

%SS=[(1:900).*0+sog(1) sog (1:900).*0+sog(length(sog))];         %resample sog
%sinTT=[(1:900).*0+sin(cog(1)./180.*pi) sin(cog./180.*pi) (1:900).*0+sin(cog(length(cog)))./180.*pi];         %resample sin(cog)
%cosTT=[(1:900).*0+cos(cog(1)./180.*pi) cos(cog./180.*pi) (1:900).*0+cos(cog(length(cog)))./180.*pi];         %resample cos(cog)

%YY=resample(SS,144000,3600);YY=YY(36001:108000);                %resample sog
%sinZZ=resample(sinTT,144000,3600);sinZZ=sinZZ(36001:108000);                %resample cog
%cosZZ=resample(cosTT,144000,3600);cosZZ=cosZZ(36001:108000);                %resample cog

%ZZ=atan2(sinZZ,cosZZ);

clear SS TT

Un_1=mean(U(0*Nlenlic+1:1*Nlenlic).*cos(heading(0*Nlenlic+1:1*Nlenlic))+V(0*Nlenlic+1:1*Nlenlic).*sin(heading(0*Nlenlic+1:1*Nlenlic)))+mean(sog(0001:0300)./1.94.*cos(cog(0001:0300)./180.*pi));
Un_2=mean(U(1*Nlenlic+1:2*Nlenlic).*cos(heading(1*Nlenlic+1:2*Nlenlic))+V(1*Nlenlic+1:2*Nlenlic).*sin(heading(1*Nlenlic+1:2*Nlenlic)))+mean(sog(0301:0600)./1.94.*cos(cog(0301:0600)./180.*pi));
Un_3=mean(U(2*Nlenlic+1:3*Nlenlic).*cos(heading(2*Nlenlic+1:3*Nlenlic))+V(2*Nlenlic+1:3*Nlenlic).*sin(heading(2*Nlenlic+1:3*Nlenlic)))+mean(sog(0601:0900)./1.94.*cos(cog(0601:0900)./180.*pi));
Un_4=mean(U(3*Nlenlic+1:4*Nlenlic).*cos(heading(3*Nlenlic+1:4*Nlenlic))+V(3*Nlenlic+1:4*Nlenlic).*sin(heading(3*Nlenlic+1:4*Nlenlic)))+mean(sog(0901:1200)./1.94.*cos(cog(0901:1200)./180.*pi));
Un_5=mean(U(4*Nlenlic+1:5*Nlenlic).*cos(heading(4*Nlenlic+1:5*Nlenlic))+V(4*Nlenlic+1:5*Nlenlic).*sin(heading(4*Nlenlic+1:5*Nlenlic)))+mean(sog(1201:1500)./1.94.*cos(cog(1201:1500)./180.*pi));
Un_6=mean(U(5*Nlenlic+1:6*Nlenlic).*cos(heading(5*Nlenlic+1:6*Nlenlic))+V(5*Nlenlic+1:6*Nlenlic).*sin(heading(5*Nlenlic+1:6*Nlenlic)))+mean(sog(1501:1800)./1.94.*cos(cog(1501:1800)./180.*pi));

Ue_1=mean(U(0*Nlenlic+1:1*Nlenlic).*sin(heading(0*Nlenlic+1:1*Nlenlic))-V(0*Nlenlic+1:1*Nlenlic).*cos(heading(0*Nlenlic+1:1*Nlenlic)))+mean(sog(0001:0300)./1.94.*sin(cog(0001:0300)./180.*pi));
Ue_2=mean(U(1*Nlenlic+1:2*Nlenlic).*sin(heading(1*Nlenlic+1:2*Nlenlic))-V(1*Nlenlic+1:2*Nlenlic).*cos(heading(1*Nlenlic+1:2*Nlenlic)))+mean(sog(0301:0600)./1.94.*sin(cog(0301:0600)./180.*pi));
Ue_3=mean(U(2*Nlenlic+1:3*Nlenlic).*sin(heading(2*Nlenlic+1:3*Nlenlic))-V(2*Nlenlic+1:3*Nlenlic).*cos(heading(2*Nlenlic+1:3*Nlenlic)))+mean(sog(0601:0900)./1.94.*sin(cog(0601:0900)./180.*pi));
Ue_4=mean(U(3*Nlenlic+1:4*Nlenlic).*sin(heading(3*Nlenlic+1:4*Nlenlic))-V(3*Nlenlic+1:4*Nlenlic).*cos(heading(3*Nlenlic+1:4*Nlenlic)))+mean(sog(0901:1200)./1.94.*sin(cog(0901:1200)./180.*pi));
Ue_5=mean(U(4*Nlenlic+1:5*Nlenlic).*sin(heading(4*Nlenlic+1:5*Nlenlic))-V(4*Nlenlic+1:5*Nlenlic).*cos(heading(4*Nlenlic+1:5*Nlenlic)))+mean(sog(1201:1500)./1.94.*sin(cog(1201:1500)./180.*pi));
Ue_6=mean(U(5*Nlenlic+1:6*Nlenlic).*sin(heading(5*Nlenlic+1:6*Nlenlic))-V(5*Nlenlic+1:6*Nlenlic).*cos(heading(5*Nlenlic+1:6*Nlenlic)))+mean(sog(1501:1800)./1.94.*sin(cog(1501:1800)./180.*pi));

Ubar_1=sqrt(Un_1.^2 + Ue_1.^2);
Ubar_2=sqrt(Un_2.^2 + Ue_2.^2);
Ubar_3=sqrt(Un_3.^2 + Ue_3.^2);
Ubar_4=sqrt(Un_4.^2 + Ue_4.^2);
Ubar_5=sqrt(Un_5.^2 + Ue_5.^2);
Ubar_6=sqrt(Un_6.^2 + Ue_6.^2);

Twindir_1=atan2(-Ue_1,-Un_1).*180/pi;if Reldir_1<0,Reldir_1=Reldir_1+360;end;if Twindir_1<0,Twindir_1=Twindir_1+360;end;
Twindir_2=atan2(-Ue_2,-Un_2).*180/pi;if Reldir_2<0,Reldir_2=Reldir_2+360;end;if Twindir_2<0,Twindir_2=Twindir_2+360;end;
Twindir_3=atan2(-Ue_3,-Un_3).*180/pi;if Reldir_3<0,Reldir_3=Reldir_3+360;end;if Twindir_3<0,Twindir_3=Twindir_3+360;end;
Twindir_4=atan2(-Ue_4,-Un_4).*180/pi;if Reldir_4<0,Reldir_4=Reldir_4+360;end;if Twindir_4<0,Twindir_4=Twindir_4+360;end;
Twindir_5=atan2(-Ue_5,-Un_5).*180/pi;if Reldir_5<0,Reldir_5=Reldir_5+360;end;if Twindir_5<0,Twindir_5=Twindir_5+360;end;
Twindir_6=atan2(-Ue_6,-Un_6).*180/pi;if Reldir_6<0,Reldir_6=Reldir_6+360;end;if Twindir_6<0,Twindir_6=Twindir_6+360;end;

%Licor_H2O=resample(Licor_H2O,72000,18000);
%Licor_CO2=resample(Licor_CO2,72000,18000);

%wa_1=cov(Wstream_1,abshum_1./rho_1);   % m/s * g/kg
%wa_2=cov(Wstream_2,abshum_2./rho_2);   % m/s * g/kg
%wa_3=cov(Wstream_3,abshum_3./rho_3);   % m/s * g/kg
%wa_4=cov(Wstream_4,abshum_4./rho_4);   % m/s * g/kg
%wa_5=cov(Wstream_5,abshum_5./rho_5);   % m/s * g/kg
%wa_6=cov(Wstream_6,abshum_6./rho_6);   % m/s * g/kg

wocj_1=cov(Licor_w_1,Licor_CO2_1);   % m/s * g/kg
wocj_2=cov(Licor_w_2,Licor_CO2_2);   % m/s * g/kg
wocj_3=cov(Licor_w_3,Licor_CO2_3);   % m/s * g/kg
wocj_4=cov(Licor_w_4,Licor_CO2_4);   % m/s * g/kg
wocj_5=cov(Licor_w_5,Licor_CO2_5);   % m/s * g/kg
wocj_6=cov(Licor_w_6,Licor_CO2_6);   % m/s * g/kg

wocj_1=wocj_1(1,2);
wocj_2=wocj_2(1,2);
wocj_3=wocj_3(1,2);
wocj_4=wocj_4(1,2);
wocj_5=wocj_5(1,2);
wocj_6=wocj_6(1,2);

ship_RH_1=mean(imrh(1:10));
ship_RH_2=mean(imrh(11:20));
ship_RH_3=mean(imrh(21:30));
ship_RH_4=mean(imrh(31:40));
ship_RH_5=mean(imrh(41:50));
ship_RH_6=mean(imrh(51:60));

RHv_1=mean(Rhvais(1:10));
RHv_2=mean(Rhvais(11:20));
RHv_3=mean(Rhvais(21:30));
RHv_4=mean(Rhvais(31:40));
RHv_5=mean(Rhvais(41:50));
RHv_6=mean(Rhvais(51:60));

Qship_1=0.622*ship_RH_1/100*6.112*exp((17.67*Tship_1)/(Tship_1+243.5))/1013.25*1000;
Qship_2=0.622*ship_RH_2/100*6.112*exp((17.67*Tship_2)/(Tship_2+243.5))/1013.25*1000;
Qship_3=0.622*ship_RH_3/100*6.112*exp((17.67*Tship_3)/(Tship_3+243.5))/1013.25*1000;
Qship_4=0.622*ship_RH_4/100*6.112*exp((17.67*Tship_4)/(Tship_4+243.5))/1013.25*1000;
Qship_5=0.622*ship_RH_5/100*6.112*exp((17.67*Tship_5)/(Tship_5+243.5))/1013.25*1000;
Qship_6=0.622*ship_RH_6/100*6.112*exp((17.67*Tship_6)/(Tship_6+243.5))/1013.25*1000;

Qbarv_1=0.622*RHv_1/100*6.112*exp((17.67*Tbar_1)/(Tbar_1+243.5))/1013.25*1000;
Qbarv_2=0.622*RHv_2/100*6.112*exp((17.67*Tbar_2)/(Tbar_2+243.5))/1013.25*1000;
Qbarv_3=0.622*RHv_3/100*6.112*exp((17.67*Tbar_3)/(Tbar_3+243.5))/1013.25*1000;
Qbarv_4=0.622*RHv_4/100*6.112*exp((17.67*Tbar_4)/(Tbar_4+243.5))/1013.25*1000;
Qbarv_5=0.622*RHv_5/100*6.112*exp((17.67*Tbar_5)/(Tbar_5+243.5))/1013.25*1000;
Qbarv_6=0.622*RHv_6/100*6.112*exp((17.67*Tbar_6)/(Tbar_6+243.5))/1013.25*1000;

Tseasnake_1=mean(Tsea(1:10));
Tseasnake_2=mean(Tsea(11:20));
Tseasnake_3=mean(Tsea(21:30));
Tseasnake_4=mean(Tsea(31:40));
Tseasnake_5=mean(Tsea(41:50));
Tseasnake_6=mean(Tsea(51:60));

% Esatsnake_1=6.112.*exp(17.67.*Tseasnake_1./(Tseasnake_1+243.5))*0.98;
% Esatsnake_2=6.112.*exp(17.67.*Tseasnake_2./(Tseasnake_2+243.5))*0.98;
% Esatsnake_3=6.112.*exp(17.67.*Tseasnake_3./(Tseasnake_3+243.5))*0.98;
% Esatsnake_4=6.112.*exp(17.67.*Tseasnake_4./(Tseasnake_4+243.5))*0.98;
% Esatsnake_5=6.112.*exp(17.67.*Tseasnake_5./(Tseasnake_5+243.5))*0.98;
% Esatsnake_6=6.112.*exp(17.67.*Tseasnake_6./(Tseasnake_6+243.5))*0.98;
% 
% Qseasnake_1=Esatsnake_1*0.622/(P-0.378*Esatsnake_1)*1000;
% Qseasnake_2=Esatsnake_2*0.622/(P-0.378*Esatsnake_2)*1000;
% Qseasnake_3=Esatsnake_3*0.622/(P-0.378*Esatsnake_3)*1000;
% Qseasnake_4=Esatsnake_4*0.622/(P-0.378*Esatsnake_4)*1000;
% Qseasnake_5=Esatsnake_5*0.622/(P-0.378*Esatsnake_5)*1000;
% Qseasnake_6=Esatsnake_6*0.622/(P-0.378*Esatsnake_6)*1000;
Qseasnake_1=qsea(Tseasnake_1);
Qseasnake_2=qsea(Tseasnake_2);
Qseasnake_3=qsea(Tseasnake_3);
Qseasnake_4=qsea(Tseasnake_4);
Qseasnake_5=qsea(Tseasnake_5);
Qseasnake_6=qsea(Tseasnake_6);

Tsea_1=mean(tsg(1:10));
Tsea_2=mean(tsg(11:20));
Tsea_3=mean(tsg(21:30));
Tsea_4=mean(tsg(31:40));
Tsea_5=mean(tsg(41:50));
Tsea_6=mean(tsg(51:60));

% Esat_1=6.112.*exp(17.67.*Tsea_1./(Tsea_1+243.5))*0.98;
% Esat_2=6.112.*exp(17.67.*Tsea_2./(Tsea_2+243.5))*0.98;
% Esat_3=6.112.*exp(17.67.*Tsea_3./(Tsea_3+243.5))*0.98;
% Esat_4=6.112.*exp(17.67.*Tsea_4./(Tsea_4+243.5))*0.98;
% Esat_5=6.112.*exp(17.67.*Tsea_5./(Tsea_5+243.5))*0.98;
% Esat_6=6.112.*exp(17.67.*Tsea_6./(Tsea_6+243.5))*0.98;
% 
% Qsea_1=Esat_1*0.622/(P-0.378*Esat_1)*1000;
% Qsea_2=Esat_2*0.622/(P-0.378*Esat_2)*1000;
% Qsea_3=Esat_3*0.622/(P-0.378*Esat_3)*1000;
% Qsea_4=Esat_4*0.622/(P-0.378*Esat_4)*1000;
% Qsea_5=Esat_5*0.622/(P-0.378*Esat_5)*1000;
% Qsea_6=Esat_6*0.622/(P-0.378*Esat_6)*1000;
Qsea_1=qsea(Tsea_1);
Qsea_2=qsea(Tsea_2);
Qsea_3=qsea(Tsea_3);
Qsea_4=qsea(Tsea_4);
Qsea_5=qsea(Tsea_5);
Qsea_6=qsea(Tsea_6);

Ssin_1=std(shead(0*Nlenlic+1:1*Nlenlic));Scos_1=std(chead(0*Nlenlic+1:1*Nlenlic));
Ssin_2=std(shead(1*Nlenlic+1:2*Nlenlic));Scos_2=std(chead(1*Nlenlic+1:2*Nlenlic));
Ssin_3=std(shead(2*Nlenlic+1:3*Nlenlic));Scos_3=std(chead(2*Nlenlic+1:3*Nlenlic));
Ssin_4=std(shead(3*Nlenlic+1:4*Nlenlic));Scos_4=std(chead(3*Nlenlic+1:4*Nlenlic));
Ssin_5=std(shead(4*Nlenlic+1:5*Nlenlic));Scos_5=std(chead(4*Nlenlic+1:5*Nlenlic));
Ssin_6=std(shead(5*Nlenlic+1:6*Nlenlic));Scos_6=std(chead(5*Nlenlic+1:6*Nlenlic));

Sp1_1=std(Platvel(1,0*Nlenlic+1:1*Nlenlic));Sp2_1=std(Platvel(2,0*Nlenlic+1:1*Nlenlic));Sp3_1=std(Platvel(3,0*Nlenlic+1:1*Nlenlic));
Sp1_2=std(Platvel(1,1*Nlenlic+1:2*Nlenlic));Sp2_2=std(Platvel(2,1*Nlenlic+1:2*Nlenlic));Sp3_2=std(Platvel(3,1*Nlenlic+1:2*Nlenlic));
Sp1_3=std(Platvel(1,2*Nlenlic+1:3*Nlenlic));Sp2_3=std(Platvel(2,2*Nlenlic+1:3*Nlenlic));Sp3_3=std(Platvel(3,2*Nlenlic+1:3*Nlenlic));
Sp1_4=std(Platvel(1,3*Nlenlic+1:4*Nlenlic));Sp2_4=std(Platvel(2,3*Nlenlic+1:4*Nlenlic));Sp3_4=std(Platvel(3,3*Nlenlic+1:4*Nlenlic));
Sp1_5=std(Platvel(1,4*Nlenlic+1:5*Nlenlic));Sp2_5=std(Platvel(2,4*Nlenlic+1:5*Nlenlic));Sp3_5=std(Platvel(3,4*Nlenlic+1:5*Nlenlic));
Sp1_6=std(Platvel(1,5*Nlenlic+1:6*Nlenlic));Sp2_6=std(Platvel(2,5*Nlenlic+1:6*Nlenlic));Sp3_6=std(Platvel(3,5*Nlenlic+1:6*Nlenlic));

rain_1=mean(org(1:10));
rain_2=mean(org(11:20));
rain_3=mean(org(21:30));
rain_4=mean(org(31:40));
rain_5=mean(org(41:50));
rain_6=mean(org(51:60));

ship_rain_1=mean(imrain(1:10));
ship_rain_2=mean(imrain(11:20));
ship_rain_3=mean(imrain(21:30));
ship_rain_4=mean(imrain(31:40));
ship_rain_5=mean(imrain(41:50));
ship_rain_6=mean(imrain(51:60));

imir_1=mean(imir(1:10));
imir_2=mean(imir(11:20));
imir_3=mean(imir(21:30));
imir_4=mean(imir(31:40));
imir_5=mean(imir(41:50));
imir_6=mean(imir(51:60));

std_imir_1=std(imir(1:10));
std_imir_2=std(imir(11:20));
std_imir_3=std(imir(21:30));
std_imir_4=std(imir(31:40));
std_imir_5=std(imir(41:50));
std_imir_6=std(imir(51:60));

Rl_1=mean(pir(1:10));Rs_1=mean(psp(1:10));
Rl_2=mean(pir(11:20));Rs_2=mean(psp(11:20));
Rl_3=mean(pir(21:30));Rs_3=mean(psp(21:30));
Rl_4=mean(pir(31:40));Rs_4=mean(psp(31:40));
Rl_5=mean(pir(41:50));Rs_5=mean(psp(41:50));
Rl_6=mean(pir(51:60));Rs_6=mean(psp(51:60));

std_Rl_1=std(pir(1:10));std_Rs_1=std(psp(1:10));
std_Rl_2=std(pir(11:20));std_Rs_2=std(psp(11:20));
std_Rl_3=std(pir(21:30));std_Rs_3=std(psp(21:30));
std_Rl_4=std(pir(31:40));std_Rs_4=std(psp(31:40));
std_Rl_5=std(pir(41:50));std_Rs_5=std(psp(41:50));
std_Rl_6=std(pir(51:60));std_Rs_6=std(psp(51:60));

Rs_ship_1=mean(imsol(1:10));
Rs_ship_2=mean(imsol(11:20));
Rs_ship_3=mean(imsol(21:30));
Rs_ship_4=mean(imsol(31:40));
Rs_ship_5=mean(imsol(41:50));
Rs_ship_6=mean(imsol(51:60));


mean_gpscog_1=mod(360+180.*atan2(mean(sin(gpshead(0001:0600)./180.*pi)),mean(cos(gpshead(0001:0600)./180.*pi)))./pi,360);
mean_gpscog_2=mod(360+180.*atan2(mean(sin(gpshead(0601:1200)./180.*pi)),mean(cos(gpshead(0601:1200)./180.*pi)))./pi,360);
mean_gpscog_3=mod(360+180.*atan2(mean(sin(gpshead(1201:1800)./180.*pi)),mean(cos(gpshead(1201:1800)./180.*pi)))./pi,360);
mean_gpscog_4=mod(360+180.*atan2(mean(sin(gpshead(1801:2400)./180.*pi)),mean(cos(gpshead(1801:2400)./180.*pi)))./pi,360);
mean_gpscog_5=mod(360+180.*atan2(mean(sin(gpshead(2401:3000)./180.*pi)),mean(cos(gpshead(2401:3000)./180.*pi)))./pi,360);
mean_gpscog_6=mod(360+180.*atan2(mean(sin(gpshead(3001:3600)./180.*pi)),mean(cos(gpshead(3001:3600)./180.*pi)))./pi,360);

%mean_gpscog_1=atan2(mean(sin(gpshead(0001:0600)/180*pi)),mean(cos(gpshead(0001:0600)/180*pi)))*180/pi;

std_gpscog_1=std(unwrap(gpshead(0001:0600)./180.*pi)*180/pi);
std_gpscog_2=std(unwrap(gpshead(0601:1200)./180.*pi)*180/pi);
std_gpscog_3=std(unwrap(gpshead(1201:1800)./180.*pi)*180/pi);
std_gpscog_4=std(unwrap(gpshead(1801:2400)./180.*pi)*180/pi);
std_gpscog_5=std(unwrap(gpshead(2401:3000)./180.*pi)*180/pi);
std_gpscog_6=std(unwrap(gpshead(3001:3600)./180.*pi)*180/pi);

max_gpscog_1=max(gpshead(0001:0600));
max_gpscog_2=max(gpshead(0601:1200));
max_gpscog_3=max(gpshead(1201:1800));
max_gpscog_4=max(gpshead(1801:2400));
max_gpscog_5=max(gpshead(2401:3000));
max_gpscog_6=max(gpshead(3001:3600));

min_gpscog_1=min(gpshead(0001:0600));
min_gpscog_2=min(gpshead(0601:1200));
min_gpscog_3=min(gpshead(1201:1800));
min_gpscog_4=min(gpshead(1801:2400));
min_gpscog_5=min(gpshead(2401:3000));
min_gpscog_6=min(gpshead(3001:3600));

%std_gpscog_1=std(gpshead(1:600));min_gpscog_1=min(gpshead(1:600));max_gpscog_1=max(gpshead(1:600));if mean_gpscog_1<0, mean_gpscog_1=mean_gpscog_1+360;end;if min_gpscog_1<0, min_gpscog_1=min_gpscog_1+360;end;if max_gpscog_1<0, max_gpscog_1=max_gpscog_1+360;end;

gpsspeed=gpsspeed./1.94;

mean_gpsspeed_1=mean(gpsspeed(0001:0600));min_gpsspeed_1=min(gpsspeed(0001:0600));max_gpsspeed_1=max(gpsspeed(0001:0600));std_gpsspeed_1=std(gpsspeed(0001:0600));
mean_gpsspeed_2=mean(gpsspeed(0601:1200));min_gpsspeed_2=min(gpsspeed(0601:1200));max_gpsspeed_2=max(gpsspeed(0601:1200));std_gpsspeed_2=std(gpsspeed(0601:1200));
mean_gpsspeed_3=mean(gpsspeed(1201:1800));min_gpsspeed_3=min(gpsspeed(1201:1800));max_gpsspeed_3=max(gpsspeed(1201:1800));std_gpsspeed_3=std(gpsspeed(1201:1800));
mean_gpsspeed_4=mean(gpsspeed(1801:2400));min_gpsspeed_4=min(gpsspeed(1801:2400));max_gpsspeed_4=max(gpsspeed(1801:2400));std_gpsspeed_4=std(gpsspeed(1801:2400));
mean_gpsspeed_5=mean(gpsspeed(2401:3000));min_gpsspeed_5=min(gpsspeed(2401:3000));max_gpsspeed_5=max(gpsspeed(2401:3000));std_gpsspeed_5=std(gpsspeed(2401:3000));
mean_gpsspeed_6=mean(gpsspeed(3001:3600));min_gpsspeed_6=min(gpsspeed(3001:3600));max_gpsspeed_6=max(gpsspeed(3001:3600));std_gpsspeed_6=std(gpsspeed(3001:3600));

Lat_1=gpslat(01);Lon_1=gpslon(01);
Lat_2=gpslat(601);Lon_2=gpslon(601);
Lat_3=gpslat(1201);Lon_3=gpslon(1201);
Lat_4=gpslat(1801);Lon_4=gpslon(1801);
Lat_5=gpslat(2401);Lon_5=gpslon(2401);
Lat_6=gpslat(3001);Lon_6=gpslon(3001);

mean_gyro_1=mod(360+180.*atan2(mean(sin(heading(0*Nlenlic+1:1*Nlenlic))),mean(cos(heading(0*Nlenlic+1:1*Nlenlic))))./pi,360);
mean_gyro_2=mod(360+180.*atan2(mean(sin(heading(1*Nlenlic+1:2*Nlenlic))),mean(cos(heading(1*Nlenlic+1:2*Nlenlic))))./pi,360);
mean_gyro_3=mod(360+180.*atan2(mean(sin(heading(2*Nlenlic+1:3*Nlenlic))),mean(cos(heading(2*Nlenlic+1:3*Nlenlic))))./pi,360);
mean_gyro_4=mod(360+180.*atan2(mean(sin(heading(3*Nlenlic+1:4*Nlenlic))),mean(cos(heading(3*Nlenlic+1:4*Nlenlic))))./pi,360);
mean_gyro_5=mod(360+180.*atan2(mean(sin(heading(4*Nlenlic+1:5*Nlenlic))),mean(cos(heading(4*Nlenlic+1:5*Nlenlic))))./pi,360);
mean_gyro_6=mod(360+180.*atan2(mean(sin(heading(5*Nlenlic+1:6*Nlenlic))),mean(cos(heading(5*Nlenlic+1:6*Nlenlic))))./pi,360);

std_gyro_1=std(unwrap(heading(0*Nlenlic+1:1*Nlenlic)).*180./pi);
std_gyro_2=std(unwrap(heading(1*Nlenlic+1:2*Nlenlic)).*180./pi);
std_gyro_3=std(unwrap(heading(2*Nlenlic+1:3*Nlenlic)).*180./pi);
std_gyro_4=std(unwrap(heading(3*Nlenlic+1:4*Nlenlic)).*180./pi);
std_gyro_5=std(unwrap(heading(4*Nlenlic+1:5*Nlenlic)).*180./pi);
std_gyro_6=std(unwrap(heading(5*Nlenlic+1:6*Nlenlic)).*180./pi);

min_gyro_1=min(heading(0*Nlenlic+1:1*Nlenlic)).*180./pi;
min_gyro_2=min(heading(1*Nlenlic+1:2*Nlenlic)).*180./pi;
min_gyro_3=min(heading(2*Nlenlic+1:3*Nlenlic)).*180./pi;
min_gyro_4=min(heading(3*Nlenlic+1:4*Nlenlic)).*180./pi;
min_gyro_5=min(heading(4*Nlenlic+1:5*Nlenlic)).*180./pi;
min_gyro_6=min(heading(5*Nlenlic+1:6*Nlenlic)).*180./pi;

max_gyro_1=max(heading(0*Nlenlic+1:1*Nlenlic)).*180./pi;
max_gyro_2=max(heading(1*Nlenlic+1:2*Nlenlic)).*180./pi;
max_gyro_3=max(heading(2*Nlenlic+1:3*Nlenlic)).*180./pi;
max_gyro_4=max(heading(3*Nlenlic+1:4*Nlenlic)).*180./pi;
max_gyro_5=max(heading(4*Nlenlic+1:5*Nlenlic)).*180./pi;
max_gyro_6=max(heading(5*Nlenlic+1:6*Nlenlic)).*180./pi;

mean_Tsonic_1=mean(Tsonic_1);
mean_Tsonic_2=mean(Tsonic_2);
mean_Tsonic_3=mean(Tsonic_3);
mean_Tsonic_4=mean(Tsonic_4);
mean_Tsonic_5=mean(Tsonic_5);
mean_Tsonic_6=mean(Tsonic_6);

xship_1=[mean_wind_1 Us Tsea_1 Tship_1 Qsea_1 Qship_1 Rs_ship_1 imir_1 ship_rain_1 zi P zu_ship zt_ship zq_ship];
xship_2=[mean_wind_2 Us Tsea_2 Tship_2 Qsea_2 Qship_2 Rs_ship_2 imir_2 ship_rain_2 zi P zu_ship zt_ship zq_ship];
xship_3=[mean_wind_3 Us Tsea_3 Tship_3 Qsea_3 Qship_3 Rs_ship_3 imir_3 ship_rain_3 zi P zu_ship zt_ship zq_ship];
xship_4=[mean_wind_4 Us Tsea_4 Tship_4 Qsea_4 Qship_4 Rs_ship_4 imir_4 ship_rain_4 zi P zu_ship zt_ship zq_ship];
xship_5=[mean_wind_5 Us Tsea_5 Tship_5 Qsea_5 Qship_5 Rs_ship_5 imir_5 ship_rain_5 zi P zu_ship zt_ship zq_ship];
xship_6=[mean_wind_6 Us Tsea_6 Tship_6 Qsea_6 Qship_6 Rs_ship_6 imir_6 ship_rain_6 zi P zu_ship zt_ship zq_ship];

yship_1=coul_lkb(xship_1);
yship_2=coul_lkb(xship_2);
yship_3=coul_lkb(xship_3);
yship_4=coul_lkb(xship_4);
yship_5=coul_lkb(xship_5);
yship_6=coul_lkb(xship_6);

xsnake_1=[Ubar_1 Us Tseasnake_1 Tbar_1 Qseasnake_1 Qbarv_1 Rs_1 Rl_1 rain_1 zi P zu-0.8 zt zq];
xsnake_2=[Ubar_2 Us Tseasnake_2 Tbar_2 Qseasnake_2 Qbarv_2 Rs_2 Rl_2 rain_2 zi P zu-0.8 zt zq];
xsnake_3=[Ubar_3 Us Tseasnake_3 Tbar_3 Qseasnake_3 Qbarv_3 Rs_3 Rl_3 rain_3 zi P zu-0.8 zt zq];
xsnake_4=[Ubar_4 Us Tseasnake_4 Tbar_4 Qseasnake_4 Qbarv_4 Rs_4 Rl_4 rain_4 zi P zu-0.8 zt zq];
xsnake_5=[Ubar_5 Us Tseasnake_5 Tbar_5 Qseasnake_5 Qbarv_5 Rs_5 Rl_5 rain_5 zi P zu-0.8 zt zq];
xsnake_6=[Ubar_6 Us Tseasnake_6 Tbar_6 Qseasnake_6 Qbarv_6 Rs_6 Rl_6 rain_6 zi P zu-0.8 zt zq];

ysnake_1=coul_lkb(xsnake_1);
ysnake_2=coul_lkb(xsnake_2);
ysnake_3=coul_lkb(xsnake_3);
ysnake_4=coul_lkb(xsnake_4);
ysnake_5=coul_lkb(xsnake_5);
ysnake_6=coul_lkb(xsnake_6);

% %%%%  Quenching correction 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
qwater=[Licor_H2O_1; Licor_H2O_2;Licor_H2O_3; Licor_H2O_4; Licor_H2O_5; Licor_H2O_6];
mq=nanmedian(qwater(qwater<35 & qwater>=0));
urel=sqrt(U.^2 + V.^2 + W.^2);

qwater(qwater>mq+3)=mq;
qwater(qwater<mq-3)=mq;

%calculating mixing ratio
mr=qwater/1000./(1-qwater/1000);%watch out, Licor H2O is in g/kg (line 916 nau_correct)-->/1000

mr_att=TubingAttenuation(mr,sampleRate,0.4);
mr_att(1:30)=nanmedian(mr_att);

% ozone=(1+5.*mr).*ozone;
% ozone=(1+5.*mr_att).*ozone; % Quenching only
% ozone=(1+6.85.*mr_att).*ozone; %includes density fluctuation corrections (see Lee et al, 2000)

if ddd==205 & hhh>16 | ddd==206 | ddd==207 & hhh<12 
   O3ppb_real=(1+6.85.*mr_att).*(fastO3_real+2299.5)/4696.4;          %calibration Jd 206, for Jd 206, hr 17 to Jd 207 hr 11., incluedes quenching correction
else
   O3ppb_real=(1+6.85.*mr_att).*(fastO3_real-2155.44)/3823.13;       %Calibration average over all calibrations ecxcept Jd 206
end;
   slowO3_real=1.0371*slowO3_real+2.9742;             % incluesd calibration of the ML8810, this formula is also inclueded in the calibration of the FRCI

% save([way_rawdata, '\GOMECC\flux_data\O3ppb_real_interpcorr' , jd,hr '.txt'], 'O3ppb_real', '-ascii')

%Ozone old method / cross correlation

[c1, lags1]=xcov(Wj_1,O3ppb_real(1:6000),100,'unbiased');
    x=1:1:6000;        
    figure;[AX,H1,H2] = plotyy(x,Wj_1,x,O3ppb_real(1:6000),'plot');
    set(get(AX(1),'Ylabel'),'string','W (m*s-1)'); 
    set(get(AX(2),'Ylabel'),'string','O3 (ppb)');
    title([jd,hr]);
    %figure;plot(Wj_1);hold on;plot(O3ppb_real(1:6000),'r'); ylabel('w / m.s and O3 / ppb');
    print([way_rawdata,'\GOMECC\flux_data\lagtime_', jd, hr ,'.jpg'],'-djpeg');
    p = polyfit(lags1,c1,26);
    covEpol = polyval(p,lags1);
    %figure;plot(lags1,c1,'.',lags1,covEpol,'r--')
    ll=find(-100<lags1 & lags1<-35);
    jj=find(covEpol(ll)==min(covEpol(ll)));
    lagtime_1=-lags1(ll(jj));  %lagtime in points
    flux_1_o=min(covEpol(ll(jj)));

[c2, lags2]=xcov(Wj_2,O3ppb_real(6001:12000),100,'unbiased');
    %figure;plot(Wj_2);hold on;plot(O3ppb_real(6001:12000),'r');
    p = polyfit(lags2,c2,26);
    covEpol = polyval(p,lags2);
   % figure;plot(lags2,c2,'.',lags2,covEpol,'r--')
    ll=find(-100<lags2 & lags2<-35);
    jj=find(covEpol(ll)==min(covEpol(ll)));
    lagtime_2=-lags2(ll(jj));  %lagtime in points
    flux_2_o=min(covEpol(ll(jj)));

[c3, lags3]=xcov(Wj_3,O3ppb_real(12001:18000),100,'unbiased');
    %figure;plot(Wj_3);hold on;plot(O3ppb_real(12001:18000),'r');
    p = polyfit(lags3,c3,26);
    covEpol = polyval(p,lags3);
    %figure;plot(lags3,c3,'.',lags3,covEpol,'r--')
    ll=find(-100<lags3 & lags3<-35);
    jj=find(covEpol(ll)==min(covEpol(ll)));
    lagtime_3=-lags3(ll(jj));  %lagtime in points
    flux_3_o=min(covEpol(ll(jj)));

[c4, lags4]=xcov(Wj_4,O3ppb_real(18001:24000),100,'unbiased');
    %figure;plot(Wj_4);hold on;plot(O3ppb_real(18001:24000),'r');
    p = polyfit(lags4,c4,26);
    covEpol = polyval(p,lags4);
    %figure;plot(lags4,c4,'.',lags4,covEpol,'r--')
    ll=find(-100<lags4 & lags4<-35);
    jj=find(covEpol(ll)==min(covEpol(ll)));
    lagtime_4=-lags4(ll(jj));  %lagtime in points
    flux_4_o=min(covEpol(ll(jj)));
    
[c5, lags5]=xcov(Wj_5,O3ppb_real(24001:30000),100,'unbiased');
    %figure;plot(Wj_5);hold on;plot(O3ppb_real(24001:30000),'r');
    p = polyfit(lags5,c5,26);
    covEpol = polyval(p,lags5);
    %figure;plot(lags5,c5,'.',lags5,covEpol,'r--')
    ll=find(-100<lags5 & lags5<-35);
    jj=find(covEpol(ll)==min(covEpol(ll)));
    lagtime_5=-lags5(ll(jj));  %lagtime in points
    flux_5_o=min(covEpol(ll(jj)));

[c6, lags6]=xcov(Wj_6,O3ppb_real(30001:36000),100,'unbiased');
   %figure;plot(Wj_6);hold on;plot(O3ppb_real(30001:36000),'r'); 
    p = polyfit(lags6,c6,26);
    covEpol = polyval(p,lags6);
    %figure;plot(lags6,c6,'.',lags6,covEpol,'r--')
    ll=find(-100<lags6 & lags6<-35);
    jj=find(covEpol(ll)==min(covEpol(ll)));
    lagtime_6=-lags6(ll(jj));  %lagtime in points
    flux_6_o=min(covEpol(ll(jj)));

O3mean_1_o=mean(O3ppb_real(1:6000));
O3mean_2_o=mean(O3ppb_real(6001:12000));
O3mean_3_o=mean(O3ppb_real(12001:18000));
O3mean_4_o=mean(O3ppb_real(18001:24000));
O3mean_5_o=mean(O3ppb_real(24001:30000));
O3mean_6_o=mean(O3ppb_real(30001:36000));

O3stdev_1_o=std(O3ppb_real(1:6000));
O3stdev_2_o=std(O3ppb_real(6001:12000));
O3stdev_3_o=std(O3ppb_real(12001:18000));

O3stdev_4_o=std(O3ppb_real(18001:24000));
O3stdev_5_o=std(O3ppb_real(24001:30000));
O3stdev_6_o=std(O3ppb_real(30001:36000));

slowO3mean_1_o=mean(slowO3_real(1:6000));
slowO3mean_2_o=mean(slowO3_real(6001:12000));
slowO3mean_3_o=mean(slowO3_real(12001:18000));
slowO3mean_4_o=mean(slowO3_real(18001:24000));
slowO3mean_5_o=mean(slowO3_real(24001:30000));
slowO3mean_6_o=mean(slowO3_real(30001:36000));

slowO3stdev_1_o=std(slowO3_real(1:6000));
slowO3stdev_2_o=std(slowO3_real(6001:12000));
slowO3stdev_3_o=std(slowO3_real(12001:18000));
slowO3stdev_4_o=std(slowO3_real(18001:24000));
slowO3stdev_5_o=std(slowO3_real(24001:30000));
slowO3stdev_6_o=std(slowO3_real(30001:36000));




%Ozone new method
O3ppb_real=O3ppb_real(49:35988); 
mr_att=mr_att(49:35988);

O3mean_1=mean(O3ppb_real(1:5990));
O3mean_2=mean(O3ppb_real(5991:11980));
O3mean_3=mean(O3ppb_real(11981:17970));
O3mean_4=mean(O3ppb_real(17971:23960));
O3mean_5=mean(O3ppb_real(23961:29950));
O3mean_6=mean(O3ppb_real(29951:35940));

O3stdev_1=std(O3ppb_real(1:5990));
O3stdev_2=std(O3ppb_real(5991:11980));
O3stdev_3=std(O3ppb_real(11981:17970));
O3stdev_4=std(O3ppb_real(17971:23960));
O3stdev_5=std(O3ppb_real(23961:29950));
O3stdev_6=std(O3ppb_real(29951:35940));

flux_1=cov(Wcat(1:5990),O3ppb_real(1:5990));
flux_1=flux_1(1,2);

flux_2=cov(Wcat(5991:11980),O3ppb_real(5991:11980));
flux_2=flux_2(1,2);

flux_3=cov(Wcat(11981:17970),O3ppb_real(11981:17970));
flux_3=flux_3(1,2);

flux_4=cov(Wcat(17971:23960),O3ppb_real(17971:23960));
flux_4=flux_4(1,2);

flux_5=cov(Wcat(23961:29950),O3ppb_real(23961:29950));
flux_5=flux_5(1,2);

flux_6=cov(Wcat(29951:35940),O3ppb_real(29951:35940));
flux_6=flux_6(1,2);

if hhh==0 | hhh==12 | hhh==18,
    
%            Ustr=[Ustr_1 Ustr_2 Ustr_3 Ustr_4 Ustr_5 Ustr_6];
%            Cw2=mean([Cw2_1_2' Cw2_2_2' Cw2_3_2' Cw2_4_2' Cw2_5_2' Cw2_6_2']);
%            disp(['Calculating spectra for ',num2str(hhh)])
           [Sw,F]=psd2(detrend(O3ppb_real),length(O3ppb_real),10.);
           [Sws,Fsw]=specsmoo(Sw,10.);
           
          %figure; loglog(Fsw(length(Fsw)-31:length(Fsw)),((Fsw(length(Fsw)-31:length(Fsw))).^(-2/3))./(Fsw(length(Fsw)-18).^(-2/3)).*Sws(length(Fsw)-18).*Fsw(length(Fsw)-18),'r');
%            figure;loglog(Fsw,Cw2*(2*pi/mean(Ustr))^-.667*Fsw.^-.66/4,'r');
           figure;loglog(Fsw,Sws.*Fsw,'g');
           
           title(['O3 Spectra:   Day ',jd,hr]);
           xlabel('Frequency in Hz');
           ylabel('Spectrum * Frequency');
           print([way_rawdata,'\GOMECC\Analysis\spectraO3_' jd, hr,  '.jpg'],'-djpeg');
            
           disp(['Calculating cross spectra for ',num2str(hhh)])
           [C_wO3,Fson_wO3]=csd2(detrend(Wcat(1:35940)),detrend(O3ppb_real),length(O3ppb_real),10,hamming(length(O3ppb_real)));
           [C_wO3s,Fson_wO3s]=specsmoo(C_wO3,10);
           figure;semilogx(Fson_wO3s,C_wO3s.*Fson_wO3s,'g');
           
           title(['wO3 Cross Spectra:   Day ',jd,hr]);
           %legend('5/3','O3');       xlabel('Frequency in Hz');
           xlabel('Frequency in Hz');
           ylabel('Spectrum * Frequency');
           print([way_rawdata,'\GOMECC\Analysis\cross_spectraWO3_' jd, hr,  '.jpg'],'-djpeg');
           
        end;%end hhh=0 or 12 or 18           
        
if clearit
clear su sv ship_dir ship_time ship_wind
clear Ustream Vstream Wstream LongIR1 LongIR2 P Platvel RHv Result Solar1 Solar2
clear LongIR1 LongIR2 Platvel RHv Result Solar1 Solar2 Tdry
clear zi zq zt zu shead chead rain ppp px py pz heading fz 
clear firstr firsts firstt firstu
clear seconr secons secont seconu
clear thirdr thirds thirdt thirdu
clear fourtr fourts fourtt fourtu
clear fifthr fifths fiftht fifthu
clear sixthr sixths sixtht sixthu
clear Oabsorb Oclear Quw fx fy fxyz dfx1 dfx2 chead bf af Zf Qbar Qbarv
clear Velrot Tseasnake Tairv Esat ship_Tair ship_Tsea ship_cog ship_gyro ship_hum ship_sog
clear ship_cog ship_gyro RH ship_sog gpslat_sign gpslon_sign gpslat gpslon gpshour
clear gpsmin gpssec gpshead gpstime gpsspeed
clear uvwplat sw
end;

%   xship_1=[mean_wind_1 Us Tsea_1 Tship_1 Qsea_1 Qship_1 Rs_ship_1 Rl_1 ship_rain_1 zi P zu_ship zt_ship zq_ship];

%   1)   mean_wind_1 
%   2)   Us 
%   3)   Tsea_1 
%   4)   Tship_1 
%   5)   Qsea_1 
%   6)   Qship_1 
%   7)   Rs_ship_1 
%   8)   Rl_1 
%   9)   ship_rain_1 
%   10)  zi 
%   11)  P 
%   12)  zu_ship 
%   13)  zt_ship 
%   14)  zq_ship

%   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

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

if clearit
clear tau2 spinup2 ans Rhvais Tsea Tvais U V W cog cogm gpstim lat lon lrgm 
clear imdm imum imta imsol imrh imrain sog sogm tsg tss org pir psp
end;
%clear Uj_1 Uj_2 Uj_3 Uj_4 Uj_5 Uj_6 
%clear Vj_1 Vj_2 Vj_3 Vj_4 Vj_5 Vj_6 
%clear Wj_1 Wj_2 Wj_3 Wj_4 Wj_5 Wj_6 
%clear Ustream_1 Ustream_2 Ustream_3 Ustream_4 Ustream_5 Ustream_6 
%clear Vstream_1 Vstream_2 Vstream_3 Vstream_4 Vstream_5 Vstream_6 
%clear Wstream_1 Wstream_2 Wstream_3 Wstream_4 Wstream_5 Wstream_6 
%clear Tsonic Tsonic_1 Tsonic_2 Tsonic_3 Tsonic_4 Tsonic_5 Tsonic_6 
%clear abshum abshum_1 abshum_2 abshum_3 abshum_4 abshum_5 abshum_6
if hhh==0 | hhh==12 | hhh==18
  [Sw,Fw]=psd2(detrend(AAAA(3,:)),36000,10); [Sws,Fws]=specsmoo(Sw,10);
  [Sw_c,Fw_c]=psd2(detrend(uvw(3,:)),36000,10); [Sws_c,Fws_c]=specsmoo(Sw_c,10);
  figure;loglog(Fws,Fws.*Sws,'--',Fws_c,Fws_c.*Sws_c,'g--')
  xlabel('Frequency in Hz');
  ylabel('Spectrum (wind speed) * Frequency');p
  saveas(gcf,[way_rawdata,'\GOMECC\Final\motion_', jd, hr ,'.fig']);
end;