disp('red_epic_micflux')

%get cloud data






%get flux data
x=load('d:\cwf\data\epic2001\proc_hp\flux5_epic2_nohd.txt');%10-min averaged mailbox microwave data

jdy=x(:,1);%julian day at beginning of time average
   U=x(:,2);%true wind,ETL sonic (m/s)
   dir=x(:,3);%true wind direction, ETL sonic (deg)
   tsnk=x(:,4);%sea snake temperature,ETL, 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, ETL (C)   
   qse=x(:,8);%sea surface specific humidity,from snake (g/kg)
   qa=x(:,9);%air specific humidity, ETL (g/kg)nic anemometer(W/m^2)
   rs=x(:,10);%downward solar flux, ETL units (W/m^2)
   rl=x(:,11);%downward IR flux,  ETL units (W/m^2)
   org1=x(:,12);%rainrate, ETL 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, ETL (m/s)
   reldir=x(:,16);%relative wind direction (from),clockwise rel ship's bow, ETL sonic (deg)
   latx=x(:,17);%latitude, deg  (SCS pcode)
   lonx=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)
   org2=x(:,25);
   org3=x(:,26);
   org4=x(:,27);
   org1_c=x(:,28);
   org2_c=x(:,29);
   org=org3;
   
   zu_etl=18.5;%height of mean wind sensor, 17.7 m
   zt_etl=15.5;%height of mean air temperature sensor, 15.5 m
   zq_etl=15.5;%height of mean air humidity sensor, 15.5 m
   zq_ass=6.7;
   %*****   ships imet and scs data
   sog=x(:,13);%speed over ground, SCS gps, (m/s)
   cog=x(:,14);%%course over ground, SCS gps, (m/s)
   
   lat=latdecode(latx); 
   lon=-latdecode(lonx);%negative for W long
   clear yy;
for i=1:length(x);  
   xx=[U(i) 0 tsnk(i) ta(i) qse(i) qa(i) rs(i) rl(i) org(i) 600  1010 zu_etl zt_etl zq_etl 10 1] ;		%set data for basic flux alogithm
   y=cor26(xx);
   yy(i,:)=y;
end;
zou=yy(:,4);
zoq=yy(:,6);
cdn10=yy(:,19);
cen10=yy(:,21);
hlb=yy(:,2);%account for old algorithm
taub=yy(:,3);%account for old algorithm
hsb=yy(:,1);%account for old algorithm
dter=yy(:,11);
dqer=yy(:,12)*1e3;
L=yy(:,7);
qsr=yy(:,10);



   
   %*******************  set constants  ****************
tdk=273.16;
grav=9.83;%grv(lat)
Rgas=287.1;
cpa=1004.67;
rhoair=1010*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);

z=load('d:\cwf\data\epic2001\micro\microceilo_10m_nohd.txt');%10-min averaged mailbox microwave data
zz=interp1(z(:,1),z,jdy);
jdz=zz(:,1);%1    Julian Day (2001), time at the beginning of 10-min ave
cilof=zz(:,2);%2    Cloud fraction (0 to 1)
zb1=zz(:,3);%3    Cloud base height, 16% of cumulative distribution in 10 min
zbm=zz(:,4);%4    Cloud base height, 50% of cumulative distribution in 10 min (median)
zb2=zz(:,5);%5    Cloud base height, 84% of cumulative distribution in 10 min
iv=zz(:,6);%6    Mailbox column water vapor (cm)
il=zz(:,7);%7    Mailbox column water liquid (mm)
tb23=zz(:,8);%8    Mailbox 23 Ghz brightness temperature (K)
tb31=zz(:,9);%9    Mailbox 31 Ghz brightness temperature (K)
ivh=zz(:,10);%10   MMCR column water vapor (cm)
ilh=zz(:,11);%11   MMCR column water liquid (mm)
tb20=zz(:,12);%12   MMCR 21 Ghz brightness temperature (K)
tb30=zz(:,13);%13   MMCR 31 Ghz brightness temperature (K)
tb90=zz(:,14);%14   MMCR 90 Ghz brightness temperature (K)

ii=find(tb31>35);iv(ii)=NaN;il(ii)=NaN;




ii=find(isfinite(tsnk));Le=(2.501-.00237*median(tsnk(ii)))*1e6;
rlclr=(.52+.13/60*abs(lat)+(.082-.03/60*abs(lat)).*sqrt(qa)).*(5.67e-8*(ta+tdk).^4); 
js=283;je=297;
stdt=283;endt=floor(jdy(length(jdy)))+1;
epic01_Rsclear_5;%get clear sky solar flux data
rsclr=Rscl;
rlnet=0.97*(rl-5.67e-8*(tsnk+tdk).^4);
rlnetc=0.97*(rlclr-5.67e-8*(tsnk+tdk).^4);
hnet=0.95*rs+rlnet-hsb-hlb-hrain;


figure;plot(jdy,rsclr,jdy,rs,'.');set(gca,'fontsize',12);xlabel('Julian Day (2001)');ylabel('Solar Flux (W/m^2)');      
 j=294;axis([j j+1 0 1200])
plot(jdy,rsclr,jdy,rs,'.');set(gca,'fontsize',12);xlabel('Julian Day (2001)');ylabel('Solar Flux (W/m^2)');  

figure;subplot(2,1,1);plot(jdz,il*1000,jdz,ilh/2*1000);set(gca,'fontsize',12);xlabel('Julian Day (2001)');ylabel(' LWP (gm/m^2)');      
j=294;axis([j+.5 j+1 0 500])
ii=find(rsclr>50);subplot(2,1,2);plot(jdy(ii),rs(ii)./rsclr(ii),'.',[280 299],[1 1]);set(gca,'fontsize',12);xlabel('Julian Day (2001)');ylabel('Solar Transmission Coeff.');      
j=294;axis([j+.5 j+1 0 1.2])

lwp=il;
ij=find(isnan(il));lwp(ij)=ilh(ij)/2;
ik=find(~isnan(il) & ~isnan(ilh));lwp(ik)=(il(ik)+ilh(ik)/2)/2;
ilx=find(il(ik)-ilh(ik)/2>.3*lwp);lwp(ik(ilx))=min(il(ik(ilx)),ilh(ik(ilx))/2);

ii=find(rsclr>50);plot(il(ii)*1e3,rs(ii)./rsclr(ii),'.');set(gca,'fontsize',12);xlabel('LWP (gm/m^2)');ylabel('Solar Transmission Coeff.');      


