% loads & plots sp files

%{
10-min spec output arrays for spectra and cospectra: 144x481

column designations as follows, all spectra are 48 bin

COLUMNS     VARIABLE        DESCRIPTION
1           doy             decimal day of year
2:49        Fxx             frequency bins
50:97       Puu             U spectrum
98:145      Pvv             V spectrum
146:193     Pww             W spectrum
194:241     Ptt             T spectrum
242:289     Pqq             q spectrum, LI7500
290:337     Cuw             wu cospectrum
338:385     Cvw             wv cospectrum
386:433     Ctw             wt cospectrum
434:481     Cqw             wq cospectrum
%}

% insert path to desired sp file here...
% and edit file name as required
sp_path = '/Users/byronb/Raw_Data_Files.noindex/CalWater_2015/RHB/flux/Processed/sp';
da_path = '/Users/byronb/Raw_Data_Files.noindex/CalWater_2015/RHB/flux/Processed/da';
cr_path = '/Users/byronb/Raw_Data_Files.noindex/CalWater_2015/RHB/flux/Processed/cr';

sp_name = ['sp_v2_032_decorr'];
da_name = ['da_v2_032_decorr'];
cr_name = ['cr_v2_032_decorr'];

fin_sp = fullfile(sp_path,[sp_name,'.txt']);
fin_da = fullfile(da_path,[da_name,'.txt']);
fin_cr = fullfile(cr_path,[cr_name,'.txt']);
load (fin_sp);
load (fin_da);
load (fin_cr);

SP = sp_v2_032_decorr;
DA = da_v2_032_decorr;
CR = cr_v2_032_decorr;

% apply selected criteria from DA file
wtv = -DA(:,74).*(DA(:,79)+.61.*(DA(:,17)+273.16).*DA(:,80));
ws = (9.83./(DA(:,17)+273.16).*abs(wtv)*600).^.333;  % convective velocity
ug = sqrt(DA(:,74).^2+ws.^2);
ugw = sqrt(1.5*DA(:,74).^2+ws.^2);
ugu = sqrt(8*DA(:,74).^2+1.*ws.^2);
gg = sqrt(1+1.6*ws.^2./DA(:,16).^2);
kk = find((DA(:,56)<90 | DA(:,56)>-90) & DA(:,113)<5 &...
           DA(:,98)<0.6 & DA(:,130)<0.8 &...
           (sqrt(DA(:,127))./ugw < .7+.0015*DA(:,16).^2) &...
           (sqrt(DA(:,126))./ugu<1.4));
hh = intersect(kk, find(CR(:,17) < 60));

% define variables
Fx = SP(1,2:49);
Su = mean(SP(kk,50:97),1);
Sv = mean(SP(kk,98:145),1);
Sw = mean(SP(kk,146:193),1);
St = mean(SP(kk,194:241),1);
Sq = mean(SP(hh,242:289),1);
Cuw = mean(SP(kk,290:337),1);
Cvw = mean(SP(kk,338:385),1);
Ctw = mean(SP(kk,386:433),1);
Cqw = mean(SP(hh,434:481),1);

figure(1);
loglog(Fx,Fx.*Su,'b-',Fx,Fx.*Sv,'r-');
title('U,V spectra, Modified Tilt Corr, decorr');
legend('U','V','location','best');
grid;

figure(2);
loglog(Fx,Fx.*Sw,'b-');
title('W spectrum, Modified Tilt Corr, decorr');
grid;

figure(3);
semilogx(Fx,Fx.*Cuw,'b-',Fx,Fx.*Ctw,'r-',Fx,Fx.*Cqw,'g-');
title('Cospectra, Modified Tilt Corr, decorr');
legend('Cuw','Ctw','Cqw','location','best');
grid;

figure(4);
loglog(Fx,Fx.*St,'b-');
title('T spectrum, Modified Tilt Corr, decorr');
grid;
