yy=[];dta=10;
x1=load('c:\data\cwf\Dropbox\DATA\dynamo\files\DYNAMO_2011_PSD_flux_30min_leg1.txt');%read file with 10min 
jdy=x1(:,1);[jdav1,y1,yvar1]=timave_2_nonan(jdy,dta,x1);
x2=load('c:\data\cwf\Dropbox\DATA\dynamo\files\DYNAMO_2011_PSD_flux_30min_leg2.txt');%read file with 10min 
jdy=x2(:,1);[jdav2,y2,yvar1]=timave_2_nonan(jdy,dta,x2);
x3=load('c:\data\cwf\Dropbox\DATA\dynamo\files\DYNAMO_2011_PSD_flux_30min_leg3.txt');%read file with 10min 
jdy=x3(:,1);[jdav3,y3,yvar1]=timave_2_nonan(jdy,dta,x3);
x4=load('c:\data\cwf\Dropbox\DATA\dynamo\files\DYNAMO_2011_PSD_flux_30min_leg4.txt');%read file with 10min 
jdy=x4(:,1);ii=find(jdy>200);[jdav4,y4,yvar1]=timave_2_nonan(jdy(ii),dta,x4(ii,:));
yy=[y1' y2' y3' y4']';
xx=[];
ii=find(~isnan(yy(:,2)) );xx=yy(ii,:);ii=find(~isnan(xx(:,7)) );x=xx(ii,:);
jdy=x(:,1);%julian day at beginning of time average
U=x(:,2);%true wind,PSD sonic (m/s)
dir=x(:,3);%true wind direction, PSD sonic (deg)
tsnk=x(:,4);%sea snake temperature,PSD, 0.05 m depth (C)
tsg=x(:,5);%tsg water temperature, 5 m depth, (C)
sal=x(:,6);%tsg salinity (psu)
ta=x(:,7);%air temperature, PSD (C)
qse=x(:,8);%sea surface specific humidity,from snake (g/kg)
qa=x(:,9);%air specific humidity, PSD (g/kg)
rs=x(:,10)*1.00;%downward solar flux, PSD units (W/m^2)
rl=x(:,11);%downward IR flux,  PSD units (W/m^2)
org=x(:,12);%rainrate, PSD STI optical rain gauge, uncorrected (mm/hr)
ushp=x(:,13);%doppler log, SCS (m/s)
head=x(:,14);%ship heading, deg clockwise rel north, SCS laser ring gyro (deg)
urel=x(:,15);%relative wind speed, PSD (m/s)
reldir=x(:,16);%relative wind direction (from),clockwise rel ship's bow, PSD sonic (deg)
lat=x(:,17);%latitude, deg  (SCS pcode)
lon=x(:,18);%longitude, deg (SCS pcode)
zts=x(:,19);%depth for bulk flux Ts reference, =0.05 when snake is used
sig_u=x(:,20);%std dev of ship speed, m/s (>.2 indicates maneuver during half hour average)
taub=x(:,21);%bulk wind stress along mean wind, (N/m^2)
hsb=x(:,22);%bulk sensible heat flux, (W/m^2)
hlb=x(:,23);%bulk latent heat flux, W/m^2 (includes Webb et al. correction)
hrain=x(:,24);%rain heat flux, as per Gosnell et al 1995, JGR, 18437-18442, (W/m^2)
ta_im=x(:,25);%ship IMET air temp @21 m (degree C)
qa_im=x(:,26);%ship IMET air specific humidity (g/kg)
s_ship=x(:,27);%ship truye wind speed  (m/s)
dir_ship=x(:,28);%ship true wind dir (degrees)
psp_ship=x(:,29);%ship IMET pyranometer (solar flux, W/m^2)
pir_ship=x(:,30);%ship IMET pirgeometer (IR flux, W/m^2)
barpress=x(:,31);%IMET barometric pressure (hPa)
RH=x(:,32);%PSD RHvais, %
rh_ship=relhum5([ta_im qa_im barpress]);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  clean up PSD data
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
ii=find(tsg>35 | tsg<25 );tsg(ii)=tsg(ii-1);
ii=find(tsnk>35 | tsnk<25 | abs(tsnk-tsg)>5 | isnan(tsnk));tsnk(ii)=tsg(ii);
ii=find(jdy>274);tsnk(ii(1):ii(1)+11)=tsg(ii(1):ii(1)+11);
ii=find(jdy>311);tsnk(ii(1):ii(1)+40)=tsg(ii(1):ii(1)+40);
ii=find(jdy>350);tsnk(ii(1):ii(1)+49)=tsg(ii(1):ii(1)+49);
ii=find(jdy>248 & jdy<256 & abs(tsnk-tsg)>1 );tsnk(ii)=tsg(ii);
ii=find(ta<22);ta(ii)=ta_im(ii);ii=find(isnan(ta));ta(ii)=ta_im(ii);
ii=find(rl<300);rl(ii)=400;
RH=RH/.95;ii=find(jdy>243.5 & jdy<245.4 | RH>120);RH(ii)=rh_ship(ii);RHx=RH;
ii=find(isnan(RH));RH(ii)=rh_ship(ii);
qa=qair3([barpress ta RH]);
ii=find(isnan(rl));rl(ii)=pir_ship(ii);
qse = qsat26sea(tsnk,barpress);
tdw=tdew([barpress qa]);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%   set up for warm layer calc
zu_etl=18*ones(length(jdy),1);
zt_etl=15*ones(length(jdy),1);
zq_etl=15*ones(length(jdy),1);
Lon=lon;Lat=lat;
Ts=tsg;   ts_depth=4;%Use tsg for water T, depth of snake is 4 m
%Ts=tsnk;   ts_depth=.05;%Use sea snake for water T, depth of snake is 0.05 m
tsg_depth=4;
tdk=273.15;
zts=tsg_depth*ones(length(jdy),1);
Rgas=287.1;
rhoa=barpress*100./(Rgas*(ta+tdk).*(1+0.61.*qa/1000));
%Hsb=hsb;Hlb=hlb;Taub=taub;usb=sqrt(taub./rhoa);

read_DYNAMO_ceilo_all;
read_DYNAMO_radar_all;
read_DYNAMO_micro_all_pc;


tdk=273.16;
grav=grv(lat);%9.82;
Rgas=287.1;
cpa=1004.67;
rhoair=barpress.*100./(Rgas*(ta+tdk).*(1+0.61*qa/1000));
visa=1.326e-5*(1+6.542e-3*ta+8.301e-6*ta.*ta-4.84e-9*ta.*ta.*ta);

ii=find(sigvap>0.5);vap(ii)=qa(ii)./3.3;
ii=find(isfinite(tsnk));Le=(2.501-.00237*median(tsnk(ii)))*1e6;
rlclr1=(.5+.13/60*abs(lat)+(.091-.03/60*abs(lat)).*sqrt(qa)).*(5.67e-8*(ta+tdk).^4); 
%rlclr2=(.5+.13/60*abs(lat)+(.091-.03/60*abs(lat)).*sqrt(qa) ).*(5.67e-8*(ta+tdk).^4)-8.5+2.8*vap; 
rlclr2=(.52+.13/60*abs(lat)+(.082-.03/60.*abs(lat)).*sqrt(qa)+0.0084*vap-0.0451).*(5.67e-8*(ta+tdk).^4);
%rlclr3=(.5+.13/60*abs(lat)+(.091-.03/60*abs(lat)).*sqrt(qa) -.0188+.0063*vap).*(5.67e-8*(ta+tdk).^4); 
rlclr=rlclr2;

jd=300;figure;
subplot(3,1,1);semilogy(jdy,Tot_rangegates,'-x');xlim([jd jd+1]);ylabel('Total Radar Rangate Returns');
subplot(3,1,2);plot(jdy,cloudhgt,'-o',jdy,cloud_top*1e3,'-o');xlim([jd jd+1]);ylabel('Cloud top and base heights (m)');legend('Ceilo Cloud Base','Radar Cloud Top');
subplot(3,1,3);plot(jdy,org,'-d');xlim([jd jd+1]);ylabel('Rainrate (mm/hr)');
xlabel('Year Day (2011)');




