function [dens]=density(p,t,s)
%
% function [dens]=density(p,t,s)
% p	m*1 column vector of standard pressures in decibars
% t	m*n matrix of temperatures (NaN used as flags) in celsius
% s	m*n matrix of salinities (NaN used as flags) in parts per
%                  thousand
% dens  m*n matrix of sigma, the density anomoly, i.e., in situ
%                  density in kg/m^3 - 1000.

%     Copyright J. V. Mansbridge, CSIRO, march 16 1992


[m,n]=size(t);
dens=zeros(m,n);
p_matrix = p*ones(1, n);
dens=sigma(p_matrix,t,s);
%for j=1:n
%	dens(:,j)=sigma(p,t(:,j),s(:,j));
%end
