
% -- 
% Jeff Sevadjian
% Monterey Bay Aquarium Research Institute
% 7700 Sandholdt Rd
% Moss Landing, CA 95039, USA
% ph: +1 831-775-1888

% Loads pCO2 data for a WG deployment;
% Queries BOG database for underway data during that deployment;
% GUI to select MBTS lines for both WG and UW
% Saves .mat files as cell structures (one cell for each MBTS line) 
%  in U:/Waveglider/wgXXX_yyyymmdd/

% Scripts for accessing databases (VARS, BOG, etc.) are now in 
%  U:/m/database/

atlas = '//atlas.shore.mbari.org/';

%inputs
 veh = 'wgTiny';
 dep = '20170605';

%load stuff for map
 load('C:/EK500/GIS/GIS.mat');
 clear ML WH;

%define vector to project vehicle positions onto
 [v(1),v(2)] = lonlat2km(C1(1),C1(2),M2(1),M2(2));
 vHat = v./sqrt(v(1)^2 + v(2)^2);
 kmBin = -5:60;
 clear v;

%CO2SYS settings
 PAR1TYPE = 1; %TA
 pHSCALEIN = 1; %total scale
 K1K2CONSTANTS = 4; %Merbach
 KSO4CONSTANTS = 1; %Dickson
 PRESIN = 0;
 PRESOUT = 0;
 SI = 10;
 PO4 = 1;

%--%

%axes param's
 ax.FontSize = 14;
 ax.Box = 'on';
 ax.NextPlot = 'add';
 %map
  ax1.XLim = [-122.4695 -121.722];
  ax1.YLim = [36.55 37];
  ax1.XTick = [round(ax1.XLim(1)*10)/10:0.2:ax1.XLim(2)];
  ax1.YTick = [round(ax1.YLim(1)*10)/10:0.2:ax1.YLim(2)];
  ax1.XGrid = 'on';
  ax1.YGrid = 'on';
  ax1.DataAspectRatio = [1 mean(cos(ax1.YLim*pi/180)) 1];
 %time vs lon
  ax2.XTick = [];

%--%

%load WG data
 mydir = [atlas 'WaveGlider/deployment_data/' ...
    veh '/' dep '/QC/'];
 dn = datenum([1970 1 1]);
 %Position
  wgPos = csvread([mydir 'Position.csv'],1,0);
  wgPos(wgPos==0) = NaN;
  wgPos(:,1) = dn + wgPos(:,1)/86400/1000;
 %CO2 [pH, pCO2water, pCO2air]
  wgCO2 = csvread([mydir 'CO2.csv'],1,0);
  wgCO2(wgCO2==0) = NaN;
  wgCO2(:,1) = dn + wgCO2(:,1)/86400/1000;
  wgCO2(:,5) = interp1(wgPos(:,1),wgPos(:,2),wgCO2(:,1));
  wgCO2(:,6) = interp1(wgPos(:,1),wgPos(:,3),wgCO2(:,1));
 %CTD
  wgCTD = csvread([mydir 'CTD.csv'],1,0);
  wgCTD(wgCTD==0) = NaN;
  wgCTD(:,1) = dn + wgCTD(:,1)/86400/1000;
  wgCTD(:,5) = interp1(wgPos(:,1),wgPos(:,2),wgCTD(:,1));
  wgCTD(:,6) = interp1(wgPos(:,1),wgPos(:,3),wgCTD(:,1));  
 clear mydir dn wgPos;

%add dist along MBTS
 %co2
  [xa,ya] = lonlat2km(C1(1),C1(2),wgCO2(:,5),wgCO2(:,6));
  wgCO2(:,7) = dot([xa ya],repmat(vHat,length(xa),1),2); %km
 %ctd
  [xa,ya] = lonlat2km(C1(1),C1(2),wgCTD(:,5),wgCTD(:,6));
  wgCTD(:,7) = dot([xa ya],repmat(vHat,length(xa),1),2);
 clear xa ya;
 
%more axes param's
 ax2.XLim = [min(wgCO2(:,1)) max(wgCO2(:,1))];
 ax2.YLim = [min(wgCO2(:,5)) max(wgCO2(:,5))] + [-0.01 0.01];
%'setExtremes' callback fn
 cb = [ 'gx = get(ha(2),''XLim''); ' ...
   'fx = find(wgCO2(:,1) >= gx(1) & wgCO2(:,1) <= gx(2)); ' ...
   'set( pl, ''XData'',wgCO2(fx,5), ''YData'',wgCO2(fx,6) ); ' ...
   'set( ti, ''String'',' ...
   ' [datestr(wgCO2(fx(1),1),''mm/dd/yyyy HH:MM'')' ...
   ' '' - '' datestr(wgCO2(fx(end),1),''mm/dd/yyyy HH:MM'')] ); '];

%index MBTS
 wgInd = [];
 fi = figure('Units','inches','Position',[0 0 11 8.5]);
 ha(1) = tight_subplot(1,1,0,[0.3 0.01],0.2);
 ha(2) = tight_subplot(1,1,0,[0.05 0.75],0.1);

%map
 axes(ha(1));
 set(ha(1),ax,ax1);
 surface(gtLon,gtLat,gtZ,gtIm,...
   'FaceColor','texturemap','EdgeColor','none',...
   'CDataMapping','direct','FaceAlpha',0.5);
 plot(array{1},array{2},'.','Color',[.3 .3 .3]);
 plot(wgCO2(:,5),wgCO2(:,6), '.','Color',[134 159 192]/255);
 pl = plot(wgCO2(:,5),wgCO2(:,6),'.','Color',[11 99 160]/255,...
   'MarkerSize',10);
 plot(C1(1),C1(2),'.','Color',[1 .3 .3],'MarkerSize',10);
 plot(M1(1),M1(2),'.','Color',[1 .3 .3],'MarkerSize',10);
 plot(M2(1),M2(2),'.','Color',[1 .3 .3],'MarkerSize',10);
 ti = title('','FontSize',ax.FontSize+2);

%time vs lon
 axes(ha(2));
 plot(wgCO2(:,1),wgCO2(:,5),'.','Color',[11 99 160]/255);
 set(ha(2),ax,ax2);
 zo = zoom;
 set(zo,'ActionPostCallback', cb,'Motion','horizontal');
 pa = pan;
 set(pa,'Motion','horizontal','ActionPostCallback', cb);
 g = get(ha(2),'Position');
 uicontrol( fi,'Units','normalized',...
   'Style','pushbutton',...
   'String','FLAG',...
   'Position',[g(1)+g(3)-0.08 g(2)+g(4)+0.003 .08 .03],...
   'Callback',['wgInd = [wgInd; wgCO2(fx([1 end]),1)''];' ...
   ' fill( wgCO2(fx([1 end end 1]),1), ax2.YLim([1 1 2 2]),' ...
   ' [.8 .8 .8], ''EdgeColor'',''none'', ''FaceAlpha'',0.5); ']);
 uiwait(fi);
 clear fx gx fi ha pl ti zo pa g;

%re-organize as cell structures
 si = size(wgInd,1);
 wg.co2 = cell(si,1);
 wg.ctd = cell(si,1);
 for n=1:si
  %co2
   f = find( wgCO2(:,1)>=wgInd(n,1) & wgCO2(:,1)<=wgInd(n,2) );
   wg.co2{n} = wgCO2(f,:);
  %ctd
   f = find( wgCTD(:,1)>=wgInd(n,1) & wgCTD(:,1)<=wgInd(n,2) );
   wg.ctd{n} = wgCTD(f,:);
 end
 clear si n f;

%save
 mydir = [atlas 'jsevadjian/Waveglider/' veh '_' dep '/'];
 if ~exist(mydir,'dir')
   mkdir(mydir);
 end 
 save([mydir 'MBTSwg'],'wg');

%--%

%load underway data from BOG
 tBeg = datestr(wgCO2(1,1),'mm/dd/yyyy');
 tEnd = datestr(wgCO2(end,1),'mm/dd/yyyy');
 sql = [ 'SELECT date_time,dec_lat,dec_long,tmp,sal,pco2,tco2,sio4 '...
   'FROM underway ' ...
   'WHERE date_time>''' tBeg ''' AND ' ...
         'date_time<''' tEnd ''' AND ' ...
         'pco2>0 ' ...
   'ORDER BY date_time' ];
 c = database('jdbc:jtds:sqlserver://solstice.shore.mbari.org:1433/BOG',...
   'everyone','guest');
 UW = jdbcquery(c, sql);
 close(c);
 clear tBeg tEnd sql c;

if ~isempty(UW.date_time)
 
  %add dist along MBTS
   [xa,ya] = lonlat2km(C1(1),C1(2),UW.dec_long,UW.dec_lat);
   UW.dist_along = dot([xa ya],repmat(vHat,length(xa),1),2);
   clear xa ya;

  %add pH from pCO2w
   PAR2TYPE = 4; %pCO2
   SAL = UW.sal;
   PAR1 = 2150.0 + 44.0.*(SAL-31.25); %TA (GF)
   PAR2 = UW.pco2;
   TEMPIN = UW.tmp;
   TEMPOUT = TEMPIN;
   [RESULT,~,~] = CO2SYS(PAR1,PAR2,PAR1TYPE,PAR2TYPE,...
     SAL,TEMPIN,TEMPOUT,PRESIN,PRESOUT,SI,PO4,pHSCALEIN,...
     K1K2CONSTANTS,KSO4CONSTANTS);
   UW.pH = RESULT(:,18);
   clear PAR2TYPE SAL PAR1 PAR2 TEMPIN TEMPOUT RESULT;

  %more axes param's
   ax2.XLim = [min(UW.date_time) max(UW.date_time)];
   ax2.YLim = [min(UW.dec_long) max(UW.dec_long)] + [-0.01 0.01];
  %'setExtremes' callback fn
   cb = [ 'gx = get(ha(2),''XLim''); ' ...
     'fx = find(UW.date_time >= gx(1) & UW.date_time <= gx(2)); ' ...
     'set( pl, ''XData'',UW.dec_long(fx), ''YData'',UW.dec_lat(fx) ); ' ...
     'set( ti, ''String'',' ...
     ' [datestr(UW.date_time(fx(1)),''mm/dd/yyyy HH:MM'')' ...
     ' '' - '' datestr(UW.date_time(fx(end)),''mm/dd/yyyy HH:MM'')] ); '];

  %index MBTS
   uwInd = [];
   fi = figure('Units','inches','Position',[0 0 11 8.5]);
   ha(1) = tight_subplot(1,1,0,[0.3 0.01],0.2);
   ha(2) = tight_subplot(1,1,0,[0.05 0.75],0.1);

  %map
   axes(ha(1));
   set(ha(1),ax,ax1);
   surface(gtLon,gtLat,gtZ,gtIm,...
     'FaceColor','texturemap','EdgeColor','none',...
     'CDataMapping','direct','FaceAlpha',0.5);
   plot(array{1},array{2},'.','Color',[.3 .3 .3]);
   plot(UW.dec_long,UW.dec_lat, '.','Color',[134 159 192]/255);
   pl = plot(UW.dec_long,UW.dec_lat,'.','Color',[11 99 160]/255,...
     'MarkerSize',10);
   plot(C1(1),C1(2),'.','Color',[1 .3 .3],'MarkerSize',10);
   plot(M1(1),M1(2),'.','Color',[1 .3 .3],'MarkerSize',10);
   plot(M2(1),M2(2),'.','Color',[1 .3 .3],'MarkerSize',10);
   ti = title('','FontSize',ax.FontSize+2);

  %time vs lon
   axes(ha(2));
   plot(UW.date_time,UW.dec_long,'.','Color',[11 99 160]/255);
   set(ha(2),ax,ax2);
   zo = zoom;
   set(zo,'ActionPostCallback', cb,'Motion','horizontal');
   pa = pan;
   set(pa,'Motion','horizontal','ActionPostCallback', cb);
   g = get(ha(2),'Position');
   uicontrol( fi,'Units','normalized',...
     'Style','pushbutton',...
     'String','FLAG',...
     'Position',[g(1)+g(3)-0.08 g(2)+g(4)+0.003 .08 .03],...
     'Callback',['uwInd = [uwInd; UW.date_time(fx([1 end]))''];' ...
     ' fill( UW.date_time(fx([1 end end 1])), ax2.YLim([1 1 2 2]),' ...
     ' [.8 .8 .8], ''EdgeColor'',''none'', ''FaceAlpha'',0.5); ']);
   uiwait(fi);
   clear fx gx fi ha pl ti zo pa g;

  %re-organize as cell structures
   si = size(uwInd,1);
   uw = cell(si,1);
   for n=1:si
     f = find( UW.date_time>=uwInd(n,1) & UW.date_time(:,1)<=uwInd(n,2) );
     uw{n} = [UW.date_time(f) UW.dec_long(f) UW.dec_lat(f) UW.tmp(f) ...
       UW.sal(f) UW.pco2(f) UW.tco2(f) UW.sio4(f) UW.dist_along(f) ...
       UW.pH(f)];
   end
   clear si n f;

  %save
   save([atlas 'jsevadjian/Waveglider/' veh '_' dep '/MBTSuw'],'uw');

end
