% plotdat.m - to plot cast data
% this one for EDDIES 07/05 OC415-2

printit=1;

cast

set(0,'defaultaxesfontname','timesroman')
set(0,'defaultaxesfontsize',12)
% keep list of starting points for CTD downcast
% cast
%cast   1     2   3   4   5   6    7    8
st_ind=[2609 180 7910 nan 241 7679 8137 2150];

%           1 2 3 4 5 6 7 8
deepcast = [0 1 0 0 1 1 1 1];

% note: casts w/ index 1 are multi-casts -  use the bcd version

mat_path = '/mnt/oc415-4a/oc415_4/ctd/process/';

%cast = input('Enter cast number: ');
if cast == 1
  minth=18;maxth=30;
  minsal=36.45;maxsal=36.95;
  minsig=23;maxsig=26.5;
elseif cast == 2
  minth=3;maxth=30;
  minsal=35.;maxsal=37;
  minsig=23;maxsig=28;
elseif cast == 7
  minth=0;maxth=30;
  minsal=34.75;maxsal=37;
  minsig=23;maxsig=28;
else
  minth=3;maxth=30;
  minsal=35.;maxsal=37;
  minsig=23;maxsig=28;
end

if cast <10
  fn = ['oc41540',int2str(cast)];
elseif cast <100
  fn = ['oc4154',int2str(cast)];
end
%fn = ['cast0',int2str(cast)];
%dirstr = ['Down  ';'Up    '];

%for ii = 1:size(dirstr,1),
%  eval(['load ',mat_path,fn,'_',lower(dirstr(ii,:))])
  eval(['load ',mat_path,fn])

  readhdr

  % compute quantities
  theta0=sw_ptmp(s0,t0,press,0);
  sigma0=sw_pden(s0,t0,press,0);
    
  % initial plot of pressure vs time
  if 0
    figure
    plot(scan,press)                    
    set(gca,'ydir','reverse')
    orient landscape
    wysiwyg
    xlabel('Scan')
    ylabel('Pressure (db)')
    %title(['Cast: ',fn,' - ',dirstr(ii,:),' cast'])
    title(['Cast: ',fn])
    ax = axis;
    axis([scan(1) scan(length(scan)) -1 1010])
    keyboard
  end
  % this plots all 3 quantites on one plot
  figure(cast)
  clf
  %hold off
  bb = find(max(press) == press);
  bb = st_ind(cast):bb(1);

  cc=find(press(bb)>0 & s0(bb)>20 & theta0(bb)>0);
  bb=bb(cc);

  % here can limit the max press plotted.
  if 0
    % plot either the top 200m or 200-600m
    lim600 = 0;
    if lim600
      aa = find(press(bb) >0 & press(bb) <= 600);
    else
      aa = find(press(bb) >0 & press(bb) <= 200);
    end
    myax = 1;
    %miny = 200;
    miny = 0;
    if lim600
      maxy = 600;
    else
      maxy = 200;
    end
  else
    % plot the entire cast
    myax=0;
    aa=1:length(bb);
    miny = 0;
    maxy = press(length(bb)) +10;
  end
  
  ntics = 10;
  ax1 = axes('Position',[.13 .11 .7750 .775]);
  h1 = line(theta0(bb(aa)),press(bb(aa)));
  set(ax1,'ydir','reverse')
  if myax
    axlim = axis;
    axis([axlim(1:2) miny maxy])
  end
  axis([minth maxth miny maxy])
  grid
  xlabel(['Theta (',setstr(176),' C) (-)'])
  ylabel('Pressure (db)')
  %title(['Cast: ',fn,' - ',dirstr(ii,:),' cast'])
  %keyboard
  hold on
  ax2 = axes('Position',get(ax1,'Position'),'xaxislocation','top', ...
	     'color','none','ydir','reverse');
  h2 = line(s0(bb(aa)),press(bb(aa)),'linestyle','--','parent',ax2);  
  if myax
    axlim = axis;
    axis([axlim(1:2) miny maxy])
  end
  axis([minsal maxsal miny maxy])
  xlabel('Salinity (PSU) (- - -)')
  %title(['Cast: ',fn,' - ',dirstr(ii,:),' cast ',Lat,' ',Lon,' ', ...
	% date,' ',time])
  title(['Cast: ',fn,'  ',Lat,' ',Lon,' ', ...
	 date,' ',time])
  box on
  xlimits = get(ax1,'XLim');
  xinc = (xlimits(2)-xlimits(1))/ntics;
  set(ax1,'XTick',[xlimits(1):xinc:xlimits(2)])
  x2limits = get(ax2,'XLim');
  xinc = (x2limits(2)-x2limits(1))/ntics;
  set(ax2,'XTick',[x2limits(1):xinc:x2limits(2)])
  ax3 = axes('Position',get(ax1,'Position'),'visible','off', ...
	     'color','none','ydir','reverse');
  h3 = line(sigma0(bb(aa))-1000,press(bb(aa)), ...
	    'linestyle','-.','parent',ax3);
  if myax
    axlim = axis;
    axis([axlim(1:2) miny maxy])
  end
  axis([minsig maxsig miny maxy])
  ax3lim = get(ax3,'XLim');
  xinc = (ax3lim(2)-ax3lim(1))/ntics;
  set(ax3,'XTick',[ax3lim(1):xinc:ax3lim(2)])
  axpos = get(ax1,'Position');
  ax4 = axes('position',[axpos(1) .05 axpos(3) .005]);
  set(ax4,'xlim',ax3lim)
  set(ax4,'XTick',[ax3lim(1):xinc:ax3lim(2)],'yticklabel',[])
  txl=xlabel('Sigma (kg/m^3) ( _ . _ . _ )');
  set(txl,'interpreter','none')
  orient tall
  wysiwyg
%end

if printit
  eval(['print -dpng ../process/oc4152',int2str(cast)]);
 % print -Poceanpr2
end
