
% Finally getting to this, today 11/30/2015, b/c
%  need some info for WG paper Table 1.

%depl info
 veh='SV2'; %SV2 or SV3
 dep='20130123'; %yyyymmdd of deployment date
 pCO2stnd=428.87; %from WG_cal_config.xlsx

if ismac
  pre='/Volumes/';
else
  pre='//atlas.shore.mbari.org/';
end

%working directory
 mydir=[pre 'WaveGlider/deployment_data/' veh '/' dep '/'];
 myfile=[veh '_' dep];

%add path to m-files/functions called from here
% i.e., codes which load M1 data for comparison
 addpath([pre 'WaveGlider/code/']);

%load WG mat data from loadSVx.m
 load([mydir myfile '.mat']);

%Get an idea of 'where' and 'when' there are data from each payload.
 load([pre 'jsevadjian/GIS/coastline_from_etopo1_Out.mat']);
 M1loc=[-122.02203 36.75060]; %a recent fix
 C1=[-121.847 36.797];
 M2=[-122.378 36.697];
 figure('Name','When&Where','Units','inches','Position',[0 0 11 8.5]);
 %map
  ha1=tight_subplot(1,1,0,[0.65 0.05],0.375);
  plot(coast(:,1),coast(:,2),'k.'); hold on;
  scatter(glider.lon,glider.lat,3,glider.time,'filled');
  plot(M1loc(1),M1loc(2),'k.','MarkerSize',10);
  plot(M2(1),M2(2),'k.','MarkerSize',10);
  plot(C1(1),C1(2),'k.','MarkerSize',10);
  colormap(jet);
  set(ha1,'XLim',[min(glider.lon) max(glider.lon)]+20*[-1 1]);
  set(ha1,'YLim',[min(glider.lat) max(glider.lat)]+20*[-1 1]);
 %time-ser
  w=whos('glider*');
  ha2=tight_subplot(length(w),1,0.01,[0.05 0.4],0.2);
  for n=1:length(w)
    axes(ha2(n));
    e=eval([w(n).name]);
    plot(e.time,e.lat,'.-');
    ylabel({w(n).name; 'lat'});
    xl(n,:)=[min(e.time) max(e.time)];
  end
  set(ha2,'XLim',[min(xl(:,1)) max(xl(:,2))],'XTick',[]);
  datetick('x','keeplimits');
  clear M2 C1 ha1 w ha2 n e xl coast;

%Do not load M1 data (G. of M.)
  
%QC

%CTD and Sub:
 figure('Name','CTD','Units','inches','Position',[0 0 11 8.5]);
 ha=tight_subplot(5,1,0.03,0.15,0.1);
 set(ha,'Box','on','NextPlot','add');
 axes(ha(1));
   plot(gliderCTD.time,gliderCTD.press,'k.-');
   ylabel('press');
 axes(ha(2));
   plot(gliderCTD.time,gliderCTD.temp,'.-');
   plot(glider.time,glider.tempSub,'c.-');
   le=legend('CTD','Sub');
   set(le,'Color','none');
   ylabel('T');
   clear le;
 axes(ha(3));
   plot(gliderCTD.time,gliderCTD.sal,'.-');
   ylabel('S');
 axes(ha(4));
   plot(gliderCTD.time,gliderCTD.o2umolkg,'.-');
   ylabel({'O_2'; '(umol kg^-^1)'});
 axes(ha(5));
   plot(gliderCTD.time,gliderCTD.o2sol,'.-');
   ylabel('O_2 sol');
 linkaxes(ha,'x');
 %index times in-water
  tIn=gliderCTD.time([9 37856]);
  fc=[]; fg=[];
  for n=1:size(tIn,1)
    fc=[fc; find(gliderCTD.time>=tIn(n,1) & gliderCTD.time<=tIn(n,2))'];
    fg=[fg; find(glider.time>=tIn(n,1) & glider.time<=tIn(n,2))'];
  end
  sc=setdiff([1:length(gliderCTD.time)],fc); %out of w
  sg=setdiff([1:length(glider.time)],fg); 
  clear fc fg n;
 %chop out-of-water data
  %glider
   fn=fieldnames(glider);
   for n=1:length(fn)
     e=eval(['glider.' fn{n}]);
     e(sg)=[];
     glider=setfield(glider,fn{n},e);
   end
   clear sg fn n e;
  %CTD
   fn=fieldnames(gliderCTD);
   for n=1:length(fn)
     e=eval(['gliderCTD.' fn{n}]);
     e(sc)=[];
     gliderCTD=setfield(gliderCTD,fn{n},e);
   end
   clear sc fn n e;
  %check
   for n=1:length(ha)
     axes(ha(n));
     for m=1:size(tIn,1)
       yl=get(ha(n),'YLim');
       plot(tIn(m,1)*[1 1],yl,'k-');
       plot(tIn(m,2)*[1 1],yl,'k-');
       set(ha(n),'YLim',yl);
     end
   end
   clear n m yl;
 %sliding variance filter on salinity
  %inputs
   win=100; %window
   nst=5; %number of stdev's above median
  %run
   sm=slidingmed(gliderCTD.sal,win);
   ss=slidingstdev(gliderCTD.sal,win);
   fs=find(abs(gliderCTD.sal-sm)>nst*ss);
   clear win nst sm ss;
 %check
  axes(ha(3));
    plot(gliderCTD.time(fs),gliderCTD.sal(fs),'r.');
  axes(ha(4));
    plot(gliderCTD.time(fs),gliderCTD.o2umolkg(fs),'r.');
  axes(ha(5));
    plot(gliderCTD.time(fs),gliderCTD.o2sol(fs),'r.');
 %NaN bad data
  gliderCTD.sal(fs)=NaN;
  gliderCTD.sigt(fs)=NaN;
  gliderCTD.o2umolkg(fs)=NaN;
  gliderCTD.o2sol(fs)=NaN;
  clear fs;

 %T vs S
  f=[];
  x=gliderCTD.sal;
  y=gliderCTD.temp;
  figure('Name','TS','Units','inches','Position',[0 0 11 8.5]);
  plot(x,y,'.');
  xlabel('S');
  ylabel('T');
  hold on;
  uicontrol('Units','inches','Style','pushbutton',...
    'Position',[5 8 1 0.3],'String','NaN',...
    'Callback',...
    ['gx=get(gca,''XLim'');' ...
     'gy=get(gca,''YLim'');' ...
     'f=[f find(x>=gx(1) & x<=gx(2) & y>=gy(1) & y<=gy(2))];' ...
     'plot(x(f),y(f),''r.'');' ]);
 %check
  axes(ha(3));
    plot(gliderCTD.time(f),gliderCTD.sal(f),'m.');
  axes(ha(4));
    plot(gliderCTD.time(f),gliderCTD.o2umolkg(f),'m.');
  axes(ha(5));
    plot(gliderCTD.time(f),gliderCTD.o2sol(f),'m.');
 %NaN bad data
  gliderCTD.sal(f)=NaN;
  gliderCTD.sigt(f)=NaN;
  gliderCTD.o2umolkg(f)=NaN;
  gliderCTD.o2sol(f)=NaN;
  clear f x y gx gy;

 %T vs O2
%  Note, O2 in MB is basically never < 62 umol/kg (2.0 mg/l), 
%   and is rarely below 143 umol/kg, according to Booth et al.
  f=[];
  x=gliderCTD.temp;
  y=gliderCTD.o2umolkg;
  figure('Name','TvsO2','Units','inches','Position',[0 0 11 8.5]);
  xlim([min(x) max(x)]);
  ylim([min(y) max(y)]);
  hold on;
  plot(x,y,'.');
  uicontrol('Units','inches','Style','pushbutton',...
    'Position',[5 8 1 0.3],'String','NaN',...
    'Callback',...
    ['gx=get(gca,''XLim'');' ...
     'gy=get(gca,''YLim'');' ...
     'f=[f find(x>=gx(1) & x<=gx(2) & y>=gy(1) & y<=gy(2))];' ...
     'plot(x(f),y(f),''r.'');' ]);
 %check
  axes(ha(4));
    plot(gliderCTD.time(f),gliderCTD.o2umolkg(f),'k.');
 %NaN bad data
  gliderCTD.o2umolkg(f)=NaN;
  clear f x y i gx gy;

 %O2 by itself
  f=[];
  x1=gliderCTD.time;
  y1=gliderCTD.o2umolkg;
  figure('Name','O2vsO2','Units','inches','Position',[0 0 11 8.5]);
  hold on;
  plot(x1,y1,'b.-');
  legend('WG');
  uicontrol('Units','inches','Style','pushbutton',...
    'Position',[5 8 1 0.3],'String','NaN',...
    'Callback',...
    ['gx=get(gca,''XLim'');' ...
     'gy=get(gca,''YLim'');' ...
     'f=[f find(x1>=gx(1) & x1<=gx(2) & y1>=gy(1) & y1<=gy(2))];' ...
     'plot(x1(f),y1(f),''r.'');' ]);
 %check
  f=unique(f);
  axes(ha(4));
    plot(gliderCTD.time(f),gliderCTD.o2umolkg(f),'y.');
 %NaN bad data
  gliderCTD.o2umolkg(f)=NaN;
  clear f x1 x2 y1 y2 gx gy;

 %get rid of variables not continuing forth with
  gliderCTD=rmfield(gliderCTD,{'cond','o2','o2Hz','o2sol','sid'});
  glider=rmfield(glider,'tempSub');


%lon,lat:
 figure('Name','LonLat','Units','inches','Position',[0 0 11 8.5]);
 haLL=tight_subplot(2,1,0.02,0.3,0.1);
 set(haLL,'Box','on','NextPlot','add');
 set(haLL(1:end-1),'XTickLabel',[]);
 axes(haLL(1));
   plot(glider.time,glider.lon,'.-');
 axes(haLL(2));
   plot(glider.time,glider.lat,'.-');
 linkaxes(haLL,'x');
 %index bad data
  fo=[];
  fa=[];
 %check
  axes(haLL(1));
   plot(glider.time(fo),glider.lon(fo),'r.');
  axes(haLL(2));
   plot(glider.time(fa),glider.lat(fa),'r.');
 %NaN bad data
  glider.lon(fo)=NaN;
  glider.lat(fa)=NaN;
  clear fo fa haLL;


%CO2: *doesn't use all CO2 QC measurements, can do manually
 figure('Name','CO2','Units','inches','Position',[0 0 11 8.5]);
 ha=tight_subplot(6,1,0.02,0.15,0.13);
 set(ha,'Box','on','NextPlot','add');
 set(ha(1:end-1),'XTickLabel',[]);
 axes(ha(1));
   plot(gliderCO2.time,gliderCO2.EquilPumpOnPress,'.-',...
     'Color',[.3 .3 1]);
   plot(gliderCO2.time,gliderCO2.EquilPumpOffPress,'.-',...
     'Color',[0 0  .8]);
   plot(gliderCO2.time,gliderCO2.AirPumpOnPress,'.-',...
     'Color',[.3 1 .3]);
   plot(gliderCO2.time,gliderCO2.AirPumpOffPress,'.-',...
     'Color',[0 .8  0]);
   plot(gliderCO2.time,gliderCO2.ZeroPumpOnPress,'.-',...
     'Color',[1 .3 .3]);
   plot(gliderCO2.time,gliderCO2.ZeroPumpOffPress,'.-',...
     'Color',[.8 0 0]);
   plot(gliderCO2.time,gliderCO2.StandardFlowOnPress,'.-',...
     'Color',[.7 .7 .7]);
   plot(gliderCO2.time,gliderCO2.StandardFlowOffPress,'.-',...
     'Color',[.2 .2 .2]);
   g=get(ha(1),'Position');
   le=legend('Equil ON','Equil OFF','Air ON','Air OFF',...
     'Zero ON','Zero OFF','Stnd ON','Stnd OFF');
   set(le,'Color','w','FontSize',5,'Location','EastOutside');
   set(ha(1),'Position',g);
   ylabel({'Pressure'; '(kPa)'});
 axes(ha(2));
   plot(gliderCO2.time,gliderCO2.ZeroPumpOffCO2,'.-');
   ylabel({'Zero'; '(uatm)'});
 axes(ha(3));
   plot(gliderCO2.time,gliderCO2.StandardFlowOffCO2,'.-');
   plot(gliderCO2.time([1 end]),pCO2stnd*[1 1],'k');
   ylabel({'Standard'; '(uatm)'});
 axes(ha(4));
   plot(gliderCO2.time,gliderCO2.pCO2water,'.-');
   ylabel({'pCO2 water'; '(uatm)'});
 axes(ha(5));
   plot(gliderCO2.time,gliderCO2.pCO2air,'.-');
   ylabel({'pCO2 air'; '(uatm)'});
 axes(ha(6));
   plot(gliderCO2.pHtime,gliderCO2.pH,'b.');
   ylabel('pH');
 linkaxes(ha,'x');
 %index times in-water
  f=[];
  for n=1:size(tIn,1)
    f=[f; find(gliderCO2.time>=tIn(n,1) & gliderCO2.time<=tIn(n,2))'];
  end
  s=setdiff([1:length(gliderCO2.time)],f); %out of w
  clear f n;
 %chop out-of-water data
  fn=fieldnames(gliderCO2);
  for n=1:length(fn)
    e=eval(['gliderCO2.' fn{n}]);
    e(:,s)=[];
    gliderCO2=setfield(gliderCO2,fn{n},e);
  end
  clear ft fn n e;
  %check
   for n=1:length(ha)
     axes(ha(n));
     for m=1:size(tIn,1)
       plot(tIn(m,1)*[1 1],get(ha(n),'YLim'),'k-');
       plot(tIn(m,2)*[1 1],get(ha(n),'YLim'),'k-');
     end
   end
   clear n m;
 %index press>90 kPa
  fe=find(gliderCO2.EquilPumpOffPress<90);
  fa=find(gliderCO2.AirPumpOffPress<90);
 %check
  axes(ha(4));
    plot(gliderCO2.time(fe),gliderCO2.pCO2water(fe),'r.');
  axes(ha(5));
    plot(gliderCO2.time(fa),gliderCO2.pCO2air(fa),'r.');
 %NaN bad data
  gliderCO2.pCO2water(fe)=NaN;
  gliderCO2.pCO2air(fa)=NaN;
  clear fe fa;

 %pH and est pH vs time
  f=[];
  x1=gliderCO2.pHtime;
  x2=gliderCO2.time;
  y1=gliderCO2.pH;
  y2=gliderCO2.pHest;
  ne=numel(gliderCO2.pH);
  figure('Name','pHvspH','Units','inches','Position',[0 0 11 8.5]);
  hold on;
  plot(reshape(x1,1,ne),reshape(y1,1,ne),'b.-');
  plot(x2,y2,'c.-');
  legend('pH meas','pH est');
  uicontrol('Units','inches','Style','pushbutton',...
    'Position',[5 8 1 0.3],'String','NaN',...
    'Callback',...
    ['gx=get(gca,''XLim'');' ...
     'gy=get(gca,''YLim'');' ...
     'f=[f; find(x1>=gx(1) & x1<=gx(2) & y1>=gy(1) & y1<=gy(2))];' ...
     'plot(x1(f),y1(f),''r.'');' ]);
 %check
  f=unique(f);
  axes(ha(6));
    plot(gliderCO2.pHtime(f),gliderCO2.pH(f),'k.');
 %NaN bad data
  gliderCO2.pH(f)=NaN;
  clear f x1 x2 y1 y2 ne gx gy;

 %pCO2 air & water vs time
  fa=[];
  fw=[];
  x1=gliderCO2.time;
  x2=gliderCO2.time;
  y1=gliderCO2.pCO2air;
  y2=gliderCO2.pCO2water;
  figure('Name','pCO2vspCO2','Units','inches','Position',[0 0 11 8.5]);
  hac = tight_subplot(2,1,0.03,0.2,0.15);
  set(hac,'NextPlot','add','Box','on');
 %air
  axes(hac(1));
  plot(x1,y1,'b.-');
  legend('WG pCO2air');
 %water
  axes(hac(2));
  plot(x2,y2,'b.-');
  legend('WG pCO2water');
  linkaxes(hac,'x');
 %manual edits only for now
  fw=[60:88];
  fa=[];
 %check
  fa=unique(fa);
  fw=unique(fw);
  axes(ha(4));
    plot(gliderCO2.time(fw),gliderCO2.pCO2water(fw),'k.');
  axes(ha(5));
    plot(gliderCO2.time(fa),gliderCO2.pCO2air(fa),'k.');
 %NaN bad data
  gliderCO2.pCO2air(fa)=NaN;
  gliderCO2.pCO2water(fw)=NaN;
  clear fa fw x1 x2 x3 x4 y1 y2 y3 y4;

 %get rid of variables not continuing forth with
  gliderCO2new.time=gliderCO2.time;
  gliderCO2new.pHtime=gliderCO2.pHtime;
  gliderCO2new.pH=gliderCO2.pH;
  gliderCO2new.pCO2air=gliderCO2.pCO2air;
  gliderCO2new.pCO2water=gliderCO2.pCO2water;
  gliderCO2=gliderCO2new;
  clear gliderCO2new;

%Wx:
 figure('Name','Wx','Units','inches','Position',[0 0 11 8.5]);
 ha=tight_subplot(6,1,0.03,0.1,0.1);
 set(ha,'Box','on','NextPlot','add');
 set(ha(1:end-1),'XTickLabel',[]);
 axes(ha(1));
   plot(gliderWx.time,gliderWx.lon,'.-');
   ylabel('lon');
 axes(ha(2));
   plot(gliderWx.time,gliderWx.lat,'.-');
   ylabel('lat');
 axes(ha(3));
   plot(gliderWx.time,gliderWx.wdirAvg,'.');
   ylabel('wind dir');
 axes(ha(4));
   plot(gliderWx.time,gliderWx.wspdAvg,'.-');
   plot(gliderWx.time,gliderWx.wspdMax,'c.-');
   le=legend('avg','max');
   set(le,'Color','none');
   ylabel({'wind speed'; '(m s^-^1)'});
 axes(ha(5));
   plot(gliderWx.time,gliderWx.airPress,'.-');
   ylabel({'pressure'; '(hPa)'});
 axes(ha(6));
   plot(gliderWx.time,gliderWx.airTemp,'.-');
   ylabel({'air temp'; '(C)'});
 linkaxes(ha,'x');
 %index times in-water
  f=[];
  for n=1:size(tIn,1)
    f=[f; find(gliderWx.time>=tIn(n,1) & gliderWx.time<=tIn(n,2))'];
  end
  s=setdiff([1:length(gliderWx.time)],f); %out of w
  clear f n;
 %chop out-of-water data
  fn=fieldnames(gliderWx);
  for n=1:length(fn)
    e=eval(['gliderWx.' fn{n}]);
    e(:,s)=[];
    gliderWx=setfield(gliderWx,fn{n},e);
  end
  clear ft fn n e;
  %check
   for n=1:length(ha)
     axes(ha(n));
     for m=1:size(tIn,1)
       plot(tIn(m,1)*[1 1],get(ha(n),'YLim'),'k-');
       plot(tIn(m,2)*[1 1],get(ha(n),'YLim'),'k-');
     end
   end
   clear n m;
 %index other bad data
  fw=find(gliderWx.wdirAvg>360);
  fm=find(gliderWx.wspdMax>gliderWx.wspdAvg+20);
  ft=[];
  fp=[];
 %check
  axes(ha(1));
    plot(gliderWx.time(fw),gliderWx.lon(fw),'r.');
  axes(ha(2));
    plot(gliderWx.time(fw),gliderWx.lat(fw),'r.');
  axes(ha(4));
    plot(gliderWx.time(fw),gliderWx.wspdAvg(fw),'r.');
    plot(gliderWx.time(fw),gliderWx.wspdMax(fw),'r.');
    plot(gliderWx.time(fm),gliderWx.wspdMax(fm),'m.');
  axes(ha(5));
    plot(gliderWx.time(fp),gliderWx.airPress(fp),'r.');
   axes(ha(6));
    plot(gliderWx.time(ft),gliderWx.airTemp(ft),'r.');
 %NaN bad data
  gliderWx.lon(fw)=NaN;
  gliderWx.lat(fw)=NaN;
  gliderWx.wspdAvg(fw)=NaN;
  gliderWx.wspdMax(fw)=NaN;
  gliderWx.wspdMax(fm)=NaN;
  gliderWx.airPress(fp)=NaN;
  gliderWx.airTemp(ft)=NaN;
  clear fw fm fp ft;


%No BB2FB data this depl


%Final plot of all data to keep
 figure('Name','Final','Units','normalized','Position',[0 0 1 1]);
 ha=tight_subplot(17,1,0.005,[0.05 0.02],0.1);
 set(ha,'Box','on','NextPlot','add');
 set(ha(1:end-1),'XTickLabel',[]);
 axes(ha(1)); plot(glider.time,glider.lon,'.-'); ylabel('lon');
 axes(ha(2)); plot(glider.time,glider.lat,'.-'); ylabel('lat');
 axes(ha(8)); plot(gliderCO2.time,gliderCO2.pCO2water,'.-'); ylabel('pCO2_w');
 axes(ha(9)); plot(gliderCO2.time,gliderCO2.pCO2air,'.-'); ylabel('pCO2_a');
 axes(ha(10)); plot(gliderCO2.pHtime,gliderCO2.pH,'b.'); ylabel('pH');
 axes(ha(11)); plot(gliderCTD.time,gliderCTD.press,'.-'); ylabel('p');
 axes(ha(12)); plot(gliderCTD.time,gliderCTD.temp,'.-'); ylabel('T');
 axes(ha(13)); plot(gliderCTD.time,gliderCTD.sal,'.-'); ylabel('S');
 axes(ha(14)); plot(gliderCTD.time,gliderCTD.o2umolkg,'.-'); ylabel('O2');
 axes(ha(15)); plot(gliderWx.time,gliderWx.wspdAvg,'.-'); ylabel('Wspd_a');
 axes(ha(16)); plot(gliderWx.time,gliderWx.wspdMax,'.-'); ylabel('Wspd_m');
 axes(ha(17)); plot(gliderWx.time,gliderWx.wdirAvg,'.-'); ylabel('Wdir');
 set(ha(17),'XTick',unique(floor(glider.time)),'XTickLabel',...
   datestr(unique(floor(glider.time)),'mm/dd'));
 linkaxes(ha,'x');


%save
 w=whos('glider*');
 save([mydir myfile 'QC.mat'],w.name);
 clear w;
