% Extract2mat
%
%    extract the grid stations at standard depths to put in a .mat file
%      save pt and s as a 2d grid (depth, cast)
%      also save grid location, lat, lon, bottom depth, cast number and ctdcastlist
%
nz=33;
stddep=[0 10 20 30 50 75 100 125 150 200 250 300 400 500 600 700 800 900 ...
1000 1100 1200 1300 1400 1500 1750 2000 2500 3000 3500 4000 4500 5000 5500]';

%SetDefaults;

DepthResolution='StdDepth';
Line='grid';
FixGrid=0;

%CruiseName='NBP0103';o2units=1;
CruiseName='NBP0104';o2units=44.636;
%    convert oxygen to umol/kg for cruise NBP0104

CRUISEINFO=0;
SetupPlots;

[n,ncast] = size(ctdcastlist);
pt=NaN*ones(nz,ncast);
t=NaN*ones(nz,ncast);
s=NaN*ones(nz,ncast);
o2=NaN*ones(nz,ncast);
xgrid=zeros(1,ncast);
ygrid=zeros(1,ncast);
lat=zeros(1,ncast);
lon=zeros(1,ncast);
depth=zeros(1,ncast);
for i = 1:ncast
  cast=ctdcastlist(i);
  file= [DataDir '/' CruiseName '/Ctd/' DepthResolution '/ctd' castnum(cast) '.dat'];
   [nlines,ndepths,desc,p,T,S,PT,r0,r1,bv,t0,t1,s0,s1,O2,par,tr,fl] = readctd(file);
  xgrid(i)=ctdxgrid(cast);
  ygrid(i)=ctdygrid(cast);
  lat(i)=ctdlat(cast);
  lon(i)=ctdlon(cast);
  depth(i)=ctdtotdepth(cast);
  t(1:ndepths,i)=T;
  pt(1:ndepths,i)=PT;
  s(1:ndepths,i)=S;
  o2(1:ndepths,i)=O2;
end

o2=o2units*o2;

save([CruiseName '.gridstns.mat'],'stddep','ctdcastlist','xgrid','ygrid',...
     'lat','lon','depth','t','pt','s','o2');


