% HiWinGS2013_CO2specreader.m - loads spectra from concatenated project
% data file HiWinGS_2013_CO2spec_UH_v2_vCO2decorr1.txt  10-min spectra
% produced by HiWinGS_2013_CO2_process_vCO2decorr.m script and concatenated
% by run_HiWinGS_2013_CO2.m.

% file structure as follows:
%{
columns     Varible Name                Description
1           Decimal_DOY                 Decimal day-of-year
2:49        Fpic                        Frequency array for all except ogives
50:97       SCO2_7200PSD                LI-7200, wet mf, spectrum
98:145      CwCO2_7200PSD               LI-7200, wet mf, W cospectrum
146:193     SCO2_7500PSD                LI-7500, dry mf, spectrum (Webb corrected)
194:241     CwCO2_7500PSD               LI-7500, dry mf, cospectrum (Webb corrected)
242:289     SH2O_7200PSD                LI-7200, g/kg dry air, spectrum
290:337     CwH2O_7200PSD               LI-7200, g/kg dry air, W cospectrum
338:385     SH2O_7500PSD                LI-7500, g/kg wet air, spectrum
386:433     CwH2O_7500PSD               LI-7500, g/kg wet air, W cospectrum
434:481     CCO2H2O_7500PSD             LI-7500, CO2:H2O cospectrum
482:529     CCO2H2O_7200PSD             LI-7200, CO2:H2O cospectrum
530:577     Cqc7500PSD                  LI-7500, cross-correlation corrected qc cospectrum, Cqc=CCO2H2O-muw.*SH2O;
578:625     Scc7500PSD                  LI-7500, cross-correlation corrected c spectrum, Scc=SCO2-muw.*[CCO2H2O+Cqc];
626:673     muw7500PSD                  LI-7500, cross-correlation mu spectrum
674:721     Cqc7200PSD                  LI-7200, cross-correlation corrected qc cospectrum, Cqc=CCO2H2O-muw.*SH2O;
722:769     Scc7200PSD                  LI-7200, cross-correlation corrected c spectrum, Scc=SCO2-muw.*[CCO2H2O+Cqc];
770:817     muw7200PSD                  LI-7200, cross-correlation mu spectrum
818:865     C_H2O7200PSD_H2O7500PSD     LI-7200 : LI7500 water vapor cospectrum - relative response test
866:913     SCO2d_7200PSD               LI-7200, CO2 dry mf output, spectrum
914:961     CwCO2d_7200PSD              LI-7200, CO2 dry mf output, W cospectrum
962:1009    CwCO2c7500PSD               LI-7500, cross-correlation corrected cospectrum
1010:1057   CwCO2c7200PSD               LI-7200, cross-correlation corrected cospectrum
1058:1105   Srhoc7200PSDPKT             LI-7200, PKT corrected CO2 dry mf spectrum
1106:1153   Cwrhoc7200PSDPKT            LI-7200, PKT corrected CO2 dry mf W cospectrum
1154:1201   Srhoc7500PSDPKT             LI-7500, PKT corrected CO2 dry mf spectrum
1202:1249   Cwrhoc7500PSDPKT            LI-7500, PKT corrected CO2 dry mf W cospectrum
1250:1297   CCO2H2Or_7200PSD            LI-7200, cross-correlation mu spectrum
1298:1347   CwCO2d_7200PSD_ogive        LI-7200, dry air W cospectrum ogive
1348:1397   FwCO2d_7200PSD_ogive        LI-7200, frequency array, ogive
1398:1445   CwCO2d_7200PSD_cumtrapz     LI-7200, dry air integrated W cospectrum
1446:1493   Cwtson                      Sonic T W cospectrum
1494:1543   Cwt_ogive                   Sonic T W cospectrum ogive
1544:1593   Fwt_ogive                   Sonic T, frequency array, ogive
1594:1641   Cwt_cumtrapz                Sonic T, integrated W cospectrum
1642:1689   SCO2_7200PSDrecomp          LI-7200, recomputed CO2 dry mf, spectrum
1690:1737   CwCO2_7200PSDrecomp         LI-7200, recomputed CO2 dry mf, W cospectrum
1738:1785   SH2O_7200PSDrecomp          LI-7200, recomputed H2O g/kg dry air, spectrum
1786:1833   CwH2O_7200PSDrecomp         LI-7200, recomputed H2O g/kg dry air, W cospectrum
1834:1881   CTH2O_7500PSD               LI-7500, H2O:Tsonic cospectrum, cross-correlation correction
1882:1929   CTH2O_7200PSD               LI-7200, H2O:Tsonic cospectrum, cross-correlation correction
1930:1977   STson_PSD                   Sonic T spectrum
1978:2025   RqT7500PSD                  LI-7500, Rqt(f), Edson 2011
2026:2073   newmu7500PSD                LI-7500, Mu(f), Edson 2011
2074:2121   CwCO2_7500PSDdecorr         LI-7500, Corrected CO2 cospectrum, Edson 2011
2122:2171   CwCO2d_7200PSDi             LI-7200, Raw uncorrected cospectrum, dry mf output, equiv to CwCO2d_7200PSD above
2172:2221   CwCO2d_7200PSDichop         LI-7200, despiked cospectrum, dry mf output
2222:2271   F_chop                      Frequency array for above two cospetra
%}


datasetv='_v2';
dcnames='_vCO2decorr1';     %means CO2 is decorrelated WRT motion
decorr=0;
dcname=[];
data_path='/Volumes/User/bblomquist/Documents/Cruises/HiWinGS_2013/Analysis_BWB/Final_CO2_Flux_Files/';
if decorr==1
    dcname='_vdecorr';      %means W wsa decorrelated as well
%    dcname='_vCO2decorr_noWebbon7500mast';
end;

SA=load([data_path,'HiWinGS_2013_CO2spec_UH',datasetv,dcname,dcnames,'.txt']);

% split out the various spectra
DOY                     =   SA(:,1);
fco2                    =   SA(:,2:49);
SCO2_7200PSD            =   SA(:,50:97);
CwCO2_7200PSD           =   SA(:,98:145);
SCO2_7500PSD            =   SA(:,146:193);
CwCO2_7500PSD           =   SA(:,194:241);
SH2O_7200PSD            =   SA(:,242:289);
CwH2O_7200PSD           =   SA(:,290:337);
SH2O_7500PSD            =   SA(:,338:385);
CwH2O_7500PSD           =   SA(:,386:433);
CCO2H2O_7500PSD         =   SA(:,434:481);
CCO2H2O_7200PSD         =   SA(:,482:529);
Cqc7500PSD              =   SA(:,530:577);
Scc7500PSD              =   SA(:,578:625);
muw7500PSD              =   SA(:,626:673);
Cqc7200PSD              =   SA(:,674:721);
Scc7200PSD              =   SA(:,722:769);
muw7200PSD              =   SA(:,770:817);
C_H2O7200PSD_H2O7500PSD =   SA(:,818:865);
SCO2d_7200PSD           =   SA(:,866:913);
CwCO2d_7200PSD          =   SA(:,914:961); 
CwCO2c7500PSD           =   SA(:,962:1009);
CwCO2c7200PSD           =   SA(:,1010:1057);
Srhoc7200PSDPKT         =   SA(:,1058:1105); 
Cwrhoc7200PSDPKT        =   SA(:,1106:1153);
Srhoc7500PSDPKT         =   SA(:,1154:1201);
Cwrhoc7500PSDPKT        =   SA(:,1202:1249);
CCO2H2Or_7200PSD        =   SA(:,1250:1297);
CwCO2d_7200PSD_ogive    =   SA(:,1298:1347); %50 bins
FwCO2d_7200PSD_ogive    =   SA(:,1348:1397); %50 bins
CwCO2d_7200PSD_cumtrapz =   SA(:,1398:1445); 
Cwtson                  =   SA(:,1446:1493);
Cwt_ogive               =   SA(:,1494:1543); %50 bins
Fwt_ogive               =   SA(:,1544:1593); %50 bins
Cwt_cumtrapz            =   SA(:,1594:1641);
SCO2_7200PSDrecomp      =   SA(:,1642:1689);
CwCO2_7200PSDrecomp     =   SA(:,1690:1737);
SH2O_7200PSDrecomp      =   SA(:,1738:1785);
CwH2O_7200PSDrecomp     =   SA(:,1786:1833);
CTH2O_7500PSD           =   SA(:,1834:1881);
CTH2O_7200PSD           =   SA(:,1882:1929);
STson_PSD               =   SA(:,1930:1977);
RqT7500PSD              =   SA(:,1978:2025);
newmu7500PSD            =   SA(:,2026:2073);
CwCO2_7500PSDdecorr     =   SA(:,2074:2121);
CwCO2d_7200PSDi         =   SA(:,2122:2171); %50 bins
CwCO2d_7200PSDichop     =   SA(:,2172:2221); %50 bins
F_chop                  =   SA(:,2222:2271); %50 bins

% load 10-min mean QC met data



[c, ia, ib] = intersect(jdy,DOY);  %indexes common to both AB and AC
%QC criteria
jj=find((reldir(ia)>270 | reldir(ia)<90) & sig_h(ia)<5 & sig_u(ia)<1.4 &sp2(ia)<1.2 & org(ia)<5 & sqrt(ww(ia))./ugw(ia)<2+.0015*U(ia).^2 & sqrt(vv(ia))./ugu(ia)<2);
SteadyStateCO2criteria;

% SH2O_7500PSD((meanH2O_7500PSD(jj(ii7500CO2PSD))>23 | meanH2O_7500PSD(jj(ii7500CO2PSD)) <14),:)=NaN;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
figure;
set(gcf,'defaultaxesfonts',12)
set(gca,'fontname','time','fontsize',12);
loglog(fco2(1,:),fco2(1,:).*nanmedian(SH2O_7500PSD(jj(ii7500CO2PSD),:)),'b--','LineWidth',2)
hold on;loglog(fco2(1,:),fco2(1,:).*nanmedian(SH2O_7200PSD(jj(ii7200CO2PSD),:)),'b')
grid
% fco2(1,:),fco2(1,:).^(-4/3)*nanmedian(SCO2_7200PSD(jj,36)),'--','LineWidth',2);
title(['HIWINGS2013 - S_{H2O} averaged']);xlabel('freqency (Hz)');ylabel('f*S_{H2O}(f)');
legend('7500PSD','7200PSD');%ylim([1e-5 1e-2])

% jj=find((jdy(ia)>311 ) & (reldir(ia)>270 | reldir(ia)<90) & sig_h(ia)<5 & sig_u(ia)<0.6 &sp2(ia)<0.8 & org(ia)<5 & sqrt(ww(ia))./ugw(ia)<.8+.0015*U(ia).^2 & sqrt(vv(ia))./ugu(ia)<1.5);
% toto=find(DOY(jj)>341.25);
% hold on;loglog(fco2(1,:),fco2(1,:).*nanmedian(SH2O_7200PSD(jj(toto),:))*1.2,'m','LineWidth',2)


Hf=1./(1+(2.*pi.*0.4.*fco2(1,:)).^2);

%%%%%%%%%%%%%%%%%%%%%%%%   different SCO2 plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

figure;
set(gcf,'defaultaxesfonts',12)
set(gca,'fontname','time','fontsize',12);
loglog(fco2(1,:),nanmedian(SCO2_7500PSD(jj(ii7500CO2PSD),:)),'b','LineWidth',2)
hold on;loglog(fco2(1,:),nanmedian(SCO2d_7200PSD(jj(ii7200CO2PSD),:)),'m','LineWidth',2)
hold on;loglog(fco2(1,:),fco2(1,:).^(-5/3)*nanmedian(SCO2_7200PSD(jj,36)),'--','LineWidth',2);
% hold on;loglog(fco2(1,:),nanmedian(SCO2_7200UC(jj,36)),'k--')
grid
title(['HIWINGS2013 - S_{CO2} averaged']);xlabel('freqency (Hz)');ylabel('S_{CO2}(f)');
% hold on;loglog(fco2(1,:),nanmedian(SH2O_7200UC(jj,:)),'g--',fco2(1,:),nanmedian(SH2O_7500UC(jj,:)),'b--','LineWidth',2)
hold on;loglog(fco2(1,:),nanmedian(SCO2d_7200PSD(jj(ii7200CO2PSD),:))./Hf,'k')
legend('7500PSD','7200PSD','-5/3','7200PSD./H(f)');%ylim([1e-5 1e-2])

figure;grid;
for ll=1:length(ii7200CO2PSD)
%     loglog(fco2(1,:),fco2(1,:).*(SCO2d_7200PSD(jj(ii7200CO2PSD(ll)),:)));
    semilogx(fco2(1,:),fco2(1,:).*(CwCO2d_7200PSD(jj(ii7200CO2PSD(ll)),:)));
    hold on;
end;
% for ll=1:length(ii7200CO2UC)
%        loglog(fco2(1,:),fco2(1,:).*(SCO2_7200UC(jj(ii7200CO2UC(ll)),:)),'r');
%     hold on;
% end;
%%%%%%%%%%%%%%%%%%%%%%%%%   CwCO2  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
figure;
semilogx(fco2(1,:),fco2(1,:).*nanmedian(CwCO2_7500PSD(jj(ii7500CO2PSD),:)),'b--','LineWidth',2);%,fco2(1,:),fco2(1,:).*nanmedian(CwCO2_7500UC(jj(ii7500CO2UC),:)),'r--','LineWidth',2)
hold on;semilogx(fco2(1,:),fco2(1,:).*nanmedian(CwCO2_7200PSD(jj(ii7200CO2PSD),:)),'m','LineWidth',2)
% hold on;loglog(fco2(1,:),fco2(1,:).*nanmedian(SH2O_7200UC(jj,:)),'g','LineWidth',2)
grid
% fco2(1,:),fco2(1,:).^(-4/3)*nanmedian(SCO2_7200PSD(jj,36)),'--','LineWidth',2);
title(['HIWINGS2013  - C_{wCO2} averaged']);xlabel('freqency (Hz)');ylabel('f*C_{wCO2}(f)');
legend('7500PSD','7200PSD');%ylim([1e-5 1e-2])
hold on;semilogx(fco2(1,:),fco2(1,:).*nanmedian(CwCO2_7200PSD(jj(ii7200CO2PSD),:))./sqrt(Hf),'k')



%%%%%%% look at cospectra for postive CO2 flux...
jj=find((jdy(ia)>=299.4 &  jdy(ia)<=299.45) & (reldir(ia)>270 | reldir(ia)<90) & sig_h(ia)<5 & sig_u(ia)<1.4 &sp2(ia)<1.2 & org(ia)<5 & sqrt(ww(ia))./ugw(ia)<2+.0015*U(ia).^2 & sqrt(vv(ia))./ugu(ia)<2);
SteadyStateCO2criteria

figure;
%semilogx(fco2(1,:),fco2(1,:).*nanmedian(CwCO2_7500PSD(jj(ii7500CO2PSD),:)),'b--','LineWidth',2);%,fco2(1,:),fco2(1,:).*nanmedian(CwCO2_7500UC(jj(ii7500CO2UC),:)),'r--','LineWidth',2)
semilogx(fco2(1,:),fco2(1,:).*nanmedian(CwCO2_7200PSD(jj(ii7200CO2PSD),:)),'b','LineWidth',2)
% hold on;loglog(fco2(1,:),fco2(1,:).*nanmedian(SH2O_7200UC(jj,:)),'g','LineWidth',2)
grid
% fco2(1,:),fco2(1,:).^(-4/3)*nanmedian(SCO2_7200PSD(jj,36)),'--','LineWidth',2);
title(['HIWINGS2013  - C_{wCO2} averaged']);xlabel('freqency (Hz)');ylabel('f*C_{wCO2}(f)');
% legend('7500PSD','7200PSD');%ylim([1e-5 1e-2])

%%% ouch use spikes below 0.1Hz...

hold on;semilogx(fco2(1,:),fco2(1,:).*nanmedian(CwCO2d_7200PSDi(jj(ii7200CO2PSD),:)),'m','LineWidth',2)
hold on;semilogx(fco2(1,:),fco2(1,:).*nanmedian(CwCO2d_7200PSDichop(jj(ii7200CO2PSD),:)),'g','LineWidth',2)

% jj=find((jdy(ia)>=300.55 &  jdy(ia)<=300.7) & (reldir(ia)>270 | reldir(ia)<90) & sig_h(ia)<5 & sig_u(ia)<1.4 &sp2(ia)<1.2 & org(ia)<5 & sqrt(ww(ia))./ugw(ia)<2+.0015*U(ia).^2 & sqrt(vv(ia))./ugu(ia)<2);
% SteadyStateCO2criteria
% hold on;semilogx(fco2(1,:),fco2(1,:).*nanmedian(CwCO2_7200PSD(jj(ii7200CO2PSD),:)),'g','LineWidth',2)



figure;
set(gcf,'defaultaxesfonts',12)
set(gca,'fontname','time','fontsize',12);
loglog(fco2(1,:),nanmedian(SCO2_7500PSD(jj(ii7500CO2PSD),:)),'b','LineWidth',2)
hold on;loglog(fco2(1,:),nanmedian(SCO2d_7200PSD(jj(ii7200CO2PSD),:)),'m','LineWidth',2)
hold on;loglog(fco2(1,:),fco2(1,:).^(-5/3)*nanmedian(SCO2_7200PSD(jj,36)),'--','LineWidth',2);
% hold on;loglog(fco2(1,:),nanmedian(SCO2_7200UC(jj,36)),'k--')
grid
title(['HIWINGS2013 - S_{CO2} averaged']);xlabel('freqency (Hz)');ylabel('S_{CO2}(f)');
% hold on;loglog(fco2(1,:),nanmedian(SH2O_7200UC(jj,:)),'g--',fco2(1,:),nanmedian(SH2O_7500UC(jj,:)),'b--','LineWidth',2)
hold on;loglog(fco2(1,:),nanmedian(SCO2d_7200PSD(jj(ii7200CO2PSD),:))./Hf,'k')
legend('7500PSD','7200PSD','-5/3','7200PSD./H(f)');%ylim([1e-5 1e-2])

jj=find((jdy(ia)>=300.55 &  jdy(ia)<=300.7) & (reldir(ia)>270 | reldir(ia)<90) & sig_h(ia)<5 & sig_u(ia)<1.4 &sp2(ia)<1.2 & org(ia)<5 & sqrt(ww(ia))./ugw(ia)<2+.0015*U(ia).^2 & sqrt(vv(ia))./ugu(ia)<2);
SteadyStateCO2criteria
hold on;loglog(fco2(1,:),nanmedian(SCO2d_7200PSD(jj(ii7200CO2PSD),:)),'g','LineWidth',2)
