%  Tsection.m  - contour potential temperature from ctd on a vertical section.
%
%   need variables set: UsePotTemp RefPress xrange deprange
%                       DataDir CruiseName DepthResolution Line
%                       SectionNx SectionNz ctdcastlist
%                       InfoX InfoZ InfoString
%   sets variables: none
%   calls: readctd sw_ptmp

figtitle=['Temperature for ' Line ' Section'];
%           range of temperature and depth for section plot

%        open an new graphics window.
figure;

%              set interpolation depths
zi=0:SectionNz-1;
zi=zi./(SectionNz-1);
zi=deprange(1) + (deprange(2)-deprange(1))*zi;

Tarray=ones(SectionNz,SectionNx);
Zarray=ones(SectionNz,SectionNx);
Xarray=ones(SectionNz,SectionNx);
bottom=ones(1,SectionNx);
X=ones(1,SectionNx);
Y=ones(1,SectionNx);

[n,stations]=size(ctdcastlist);
Casts=0;
for station=1:stations
%           count number of casts
   Casts=Casts+1;
   cast=ctdcastlist(station);
   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);
   x=ctdxgrid(cast);
   y=ctdygrid(cast);
   bdepth=ctdtotdepth(cast);
%          convert pressure to depth
   Z=-p;
   if UsePotTemp 
%   calculate the potential temperature
      PT=zeros(size(t));
      PT = sw_ptmp(s,t,p,RefPress);
   else
      PT=t;
   end
%              fill array with observations
   bottom(1,station)=-bdepth;
   X(1,station)=x;
   Y(1,station)=y;
   [n,nz]=size(Z);
   ti=interp1(Z,PT,zi,'linear');
   Tarray(1:SectionNz,station)=ti';
   Zarray(1:SectionNz,station)=zi';
end
dx=X(1,1:stations)-X(1,1);
dy=Y(1,1:stations)-Y(1,1);
dx2=dx.*dx;
dy2=dy.*dy;
if Line(1) == 'x'
%             cross shelf arrays
   Xarray(1,1:stations)=Y(1:stations);
elseif Line(1) == 'y'
%             along shelf arrays
   Xarray(1,1:stations)=X(1:stations);
else
   Xarray(1,1:stations)=sqrt( (X(1,1:stations)-X(1,1)).^2 ...
                             +(Y(1,1:stations)-Y(1,1)).^2);
end

for i=2:SectionNz
   Xarray(i,1:stations)=Xarray(1,1:stations);
end

[C,h]=contour(Xarray(:,1:stations),Zarray(:,1:stations),...
              Tarray(:,1:stations),cl);
clabel(C,h);
title(figtitle);
xlabel('distance (km)');
ylabel('depth(m)');
text(InfoX,InfoZ,InfoString);
text(InfoX,1.1*InfoZ,CruiseName);
axis([xrange deprange]);
hold on;
plot(Xarray(1,1:stations),bottom(1:stations));
for station=1:stations
  text(Xarray(1,station),-20,'T');
end
hold off;
%disp('Tsection is done.');

%               print values to a file
fid=fopen('Tsection.asc','w')
for station=1:stations
   fprintf(fid,'%9.4f  %9.4f \n',Xarray(1,station),bottom(station));
   for iz=1:SectionNz
      fprintf(fid,'%9.0f  %9.4f \n',Zarray(iz,station),Tarray(iz,station));
   end
end
fclose(fid)
%    end of Tsection.m



