function [S,T,TSvol]=TSvolumetric(theta,salt,gamma,lon,lat,depth,gammarange,thetarange,saltrange); % Use [S,T,TSvol]=TSvolumetric(theta,salt,gamma,lon,lat,depth,gammarange,thetarange,saltrange) % % This routine calculates the volumetric theta-S diagram for waters with neutral densities between gammarange(1) and gammarange(2). The routine will calculate the volume span by theta-S values that fall between the elements in the thetarange and saltrange values % % theta is the 3D potential temperature field % salt is the 3D salinity field % gamma is the 3D neutral density field % lon is the 1D longitude field % lat is the 1D latitude field % depth is the 1D pressure field % gammarange is a two element vector that gives the maximum and minimum neutral density to be considered % thetarange is a 1D vector with the potential temperature bins % saltrange is a 1D vector with the temperature bins % % Example % [S,T,TSvol]=TSvolumetric(t,s,g,lon,lat,press,[28.1,28.5],[-1:.1:4],[34.4:.01:35]); disp('Setting up parameters'); % Compute the lateral area span by each point in the data set [X,Y]=meshgrid(lon,lat); R=sw_dist([-90 90]/pi,[0 0],'km')*1e3; ll=size(X); dlat=mean(diff(lat(:)))*pi/180; dlon=mean(diff(lon(:)))*pi/180; xx=X*pi/180;yy=Y*pi/180; al(2:ll(1),:)=2*pi*(sin((yy(1:ll(1)-1,:)+yy(2:ll(1),:))/2)-sin(lat(1)*pi/180)); al(1,:)=2*pi*(sin((yy(1,:)-dlat+yy(1,:))/2)-sin(lat(1)*pi/180)); al(ll(1)+1,:)=2*pi*(sin((yy(end,:)+dlat+yy(end,:))/2)-sin(lat(1)*pi/180)); al=diff(al)/ll(2)*R^2; % Compute the height span by each point in the data set press=[0;diff(depth(:))]; % Compute the volume span by each point in the data set for i=1:length(press) volume(:,:,i)=al*press(i); end % Compute the range of neutral densities g1=gammarange(1); g2=gammarange(2); kk=find((gamma>g1) & (gamma