data=importdata('C:\Data\GOMECC\June08\GOMECC07_filter_June08_total.txt');
data2=importdata('C:\Data\GOMECC\June08\GOMECC07_flux_10_June08uf.txt');
% This file inclueds the GOMECC data filtered for reldir, rain, lagtime,
% std(O3>3), delta O3ppb>3, filterchange, calibration, count drop. There
% fore the column O3ppb is for each filter copied and filtered values set
% to NaN

jd=data(:,1);               %julianday
u=data(:,3);                %true wind speed m/s, NaN=u<0.5 m/s
reldir=data(:,4);           %NaN=wrong wind direcrion
rain=data(:,5);             %NaN=rain>5    
lag=data(:,6);              %NaN=lag>70,<40
ppb5_o=data(:,7);             %NaN=O3ppb<5,>80
delta_o=data(:,8);            %NaN=deltaO3ppb>3,<-3
std_o=data(:,9);              %NaN=std>3
ppb5=data(:,10);             %NaN=O3ppb<5,>80
delta=data(:,11);            %NaN=deltaO3ppb>3
std=data(:,12);              %NaN=std>3
filter=data(:,13);          %NaN=filterchange
cali=data(:,14);            %NaN=calibration
miss=data(:,15);            %NaN=missing data
countdrop=data(:,16);       %NaN=count drop
total_o=data(:,17);         %total filter crosscorr without lag
total=data(:,18);           %total filter constant lagtime without lag
O3mean=data2(:,50);   
z=ppb5_o;

jj=find(isnan(lag));
lag(jj)=1;

ii=find(isnan(reldir));
reldir(ii)=2;
z(ii)=NaN;

gg=find(isnan(delta_o));
delta_o(gg)=3;
z(gg)=NaN;
gg2=find(isnan(delta));
delta(gg2)=3;

hh=find(isnan(std_o));
std_o(hh)=4;
z(hh)=NaN;
hh2=find(isnan(std));
std(hh2)=4;

ll=find(isnan(filter));
filter(ll)=5;
z(ll)=NaN;

ff=find(isnan(ppb5_o));
ppb5_o(ff)=6;
ff2=find(isnan(ppb5));
ppb5(ff2)=6;

tt=find(isnan(miss));
miss(tt)=7;
z(tt)=NaN;
 
ss=find(isnan(countdrop));
countdrop(ss)=8;
z(ss)=NaN;

dd=find(isnan(cali));
cali(dd)=9;
z(dd)=NaN;

ij=find(isnan(u));
u(ij)=10;
z(ij)=NaN;

kk=find(isnan(rain));
rain(kk)=11;
z(kk)=NaN;

figure;plot(jd(kk),rain(kk),'y.') ;hold on; plot(jd(ij),u(ij),'c.');plot(jd(dd),cali(dd),'g.');plot(jd(ss),countdrop(ss),'r.');
plot(jd(tt),miss(tt),'k.');plot(jd(ff2),ppb5(ff2),'m.');plot(jd(ll),filter(ll),'y.');plot(jd(hh2),std(hh2),'c.');
plot(jd(gg2),delta(gg2),'g.');plot(jd(ii),reldir(ii),'r.');plot(jd(jj),lag(jj),'m.');
xlabel('day of year');axis([190 220 0 12]);
legend('11:rain (19)','10:wind speed (23)','9:calibration (96)', '8:count drop (95)', '7:missing data (289)', '6:O3<5ppb,>80ppb (205)', '5:filter change (30)','4:std O3 (270)', '3:delta O3 (449)','2:relative wind direction (552)' , '1:lagtime (1480)');
ylabel('filter No.');

l=find(isnan(total));
figure;plot(jd(l),O3ppb(l),'.');hold on;plot(jd,total,'r.');
xlabel('Day of Year'); ylabel('[O3]/ppb');title('blue: eliminated, red: remainig; excluding lag time filter');
axis([190 220 -1 90])

k=find(isnan(total_o));
figure;plot(jd(k),O3ppb(k),'.');hold on;plot(jd,total_o,'r.');
xlabel('Day of Year'); ylabel('[O3]/ppb');title('blue: eliminated, red: remainig; including lag time filter');
axis([190 220 -1 90])

%j=find(isfinite(total_o));
%i=find(isnan(total_o));
%figure;plot(jd(i),slowO3(i),'.');hold on;plot(jd(j),slowO3(j),'r.');
%xlabel('Day of Year'); ylabel('[O3]/ppb');title('slowO3;blue: eliminated, red: remainig; excluding lag time filter');
%axis([190 220 -1 90])

%figure;plot(jd,O3ppb,'.');hold on;plot(jd,slowO3,'.r')
%xlabel('Day of Year'); ylabel('[O3]/ppb');
%title('red:O3 ML8810, blue: O3 FRCI');
%axis([190 220 -1 90])

%Fluxes and deposition velocity
jday=data2(:,1);
lagtime=data2(:,43);            % lagtime between w and O3 calculated by crosscovriance, 10 min lagtime
O3w_o=data2(:,44);              % w'O3' calculated with crosscovariance, 10 min lagtime, no webb correction
O3w_h=data2(:,45);              % w'O3' calculated with crosscovariance, flux at hourly mean lagtime, no webb correction
O3w_d=data2(:,46);              % w'O3' calculated with crosscovariance, flux at dayly mean lagtime, no webb correction
O3mean_o=data2(:,47);           % mean O3 concentration cross (ppb)
O3dev_o=data2(:,48);            % stdev ozone cross (ppb)
O3w=data2(:,49);                % O3'w' with constant lagtime of 5.1 s, and function cov, no webb correction
O3mean=data2(:,50);             % mean O3 concentration constant lagtime (ppb)
O3dev=data2(:,51);              % stdev ozone constant lagtime (ppb)
Vd_o=data2(:,52);               % deposition velocity 10 min lagtime (cm/s)
Vd_h=data2(:,53);               % deposition velocity hourly lagtime (cm/s)
Vd_d=data2(:,54);               % deposition velocity dayly lagtime (cm/s)
Vd=data2(:,55);                 % deposition velocity 5.1s lagtime (cm/s)
slowO3=data2(:,56);
lat=data2(:,23);
lon=data2(:,24);

map;
%excluding lag time filter
%const
Vd=Vd(1:1100);
p=find(isfinite(total(1:1100)));

Vdmean=mean(Vd(p))
Vdmedian=median(Vd(p))
ij=find(Vd(p)>0);
Vdmean2=mean(Vd(p(ij)))
Vdmedian2=median(Vd(p(ij)))

figure; plot(jday(p),Vd(p),'.');title('Vd constant lag time excluding lag time filter');
xlabel('Day of Year');ylabel('Vd / cm*s-1');

%Histogram
x=-5:0.005:5;
figure;hist(Vd(p),x);title('Vd constant lag time excluding lag time filter');
xlabel('Vd / cm*s-1'); ylabel('Number of Data Points');
hold on; plot([mean(Vd(p)) mean(Vd(p))],[0 190], 'r')

%cross
Vd_o=Vd_o(1:1100);
s=z(1:1100);
p=find(isfinite(s));
Vdmean_o=mean(Vd_o(p))
Vdmedian_o=median(Vd_o(p))
ij=find(Vd_o(p)>0);
Vdmean_o2=mean(Vd_o(p(ij)))
Vdmedian_o2=median(Vd_o(p(ij)))

figure; plot(jday(p),Vd_o(p),'.');
xlabel('Day of Year');ylabel('Vd / cm*s-1');title('Vd cross correl including lag time filter');


%Histogram
x=-5:0.005:5;
figure;hist(Vd_o(p),x);
xlabel('Vd / cm*s-1'); ylabel('Number of Data Points');title('Vd cross correl including lag time filter');
hold on; plot([mean(Vd_o(p)) mean(Vd_o(p))],[0 190], 'r')

% With lagtime filter
Vd=data2(:,55); 
lag=data(:,6);   
i=find(isnan(lag));
total(i)=NaN;
g=find(isfinite(total(1:1100)));

Vdmean3=mean(Vd(g))
Vdmedian3=median(Vd(g))
ij=find(Vd(g)>0);
Vdmean4=mean(Vd(g(ij)))
Vdmedian4=median(Vd(g(ij)))

figure; plot(jday(g),Vd(g),'.');
xlabel('Day of Year');ylabel('Vd / cm*s-1');title('Vd const lag time including lag time filter');

%Histogramm data filtered
x=-5:0.005:5;
figure;hist(Vd(g),x);
xlabel('Vd / cm*s-1'); ylabel('Number of Data Points');title('Vd constant lag time including lag time filter');
hold on; plot([mean(Vd(p)) mean(Vd(p))],[0 190], 'r')

%cross incl lag time filter
g=find(isfinite(total_o(1:1100)));
Vdmean_o3=mean(Vd_o(g))
Vdmedian_o3=median(Vd_o(g))
ij=find(Vd_o(g)>0);
Vdmean_o4=mean(Vd_o(g(ij)))
Vdmedian_o4=median(Vd_o(g(ij)))

figure; plot(jday(g),Vd_o(g),'.');
xlabel('Day of Year');ylabel('Vd / cm*s-1');title('Vd cross correl including lag time filter');

%Histogramm data filtered
x=-5:0.005:5;
figure;hist(Vd_o(g),x);
xlabel('Vd / cm*s-1'); ylabel('Number of Data Points');title('Vd cross correl including lag time filter');
hold on; plot([mean(Vd_o(p)) mean(Vd_o(p))],[0 190], 'r')