%* compute density of ocean water a longitudinal slice created  
%  by slice_levitus.m
%  draw filled contour transects through the world ocean
%
%* functions
% potential density calculation uses routines created by Phil Morgan, 
% Marine Research division of Australia's Commonwealth Scientific and 
% Industrial Research Organization.  
% sw_adtg.m     sw_pres.m
% sw_dens.m     sw_ptmp.m
% sw_dens0.m    sw_seck.m
% sw_pden.m     sw_smow.m
%
%* variables
% Z:            vertical layer depths in Levitus data set
% maplat:       latitudes of Levitus grid points, defined elsewhere
% salslice:     salinity in p.s.u. in a slice of the Levitus ocean
% tempslice:    temperature in deg C in a slice of the Levitus ocean
% SL:           slice longitude, -ve is west of the prime meridian
%
% rhoslice:    potential density, kg m^-3


%*****
% compute density & potential density
%*****

RP=0;       %select reference pressure, dbar

pslice=sw_pres(-Z'*ones(size(maplat)), ones(size(Z'))*maplat);
rhoslice=sw_dens(salslice,tempslice,pslice);
prhoslice=sw_pden(salslice,tempslice,pslice,RP);



%*********
% plot density in meridional slice through the world ocean
%*********

%* define contour values for filled contours
%  ocean surface density ranges from about 1000 to 1030
cvR=[1022:1:1060];		
cvPR=[1020:0.25:1028];

%* map plot window
windslice=[min(maplat) max(maplat) min(Z) max(Z)];


figure(200)
clf
caxis([min(cvR) max(cvR)])
axis(windslice), hold on
[csR, hR]=contourf(maplat,Z, rhoslice, cvR);
clabel(csR, hR, 'manual')			% use return values to lable contours
xlabel('latitude')
ylabel('depth (m)') 
title('density (ci = 1 kg m^-^1)')
text(100, 0, int2str(SL))
colorbar('horiz')
%shading flat


figure(201)
clf
caxis([min(cvPR) max(cvPR)])
axis(windslice), hold on
[csPR, hPR]=contourf(maplat,Z, prhoslice, cvPR);
clabel(csPR, hPR, 'manual')			% use return values to lable contours
xlabel('latitude')
ylabel('depth (m)') 
title('potential density, sea level pressure reference (ci = 0.25 kg m^-^1)')
text(100, 0, int2str(SL))
colorbar('horiz')
%shading flat