| 1 |
function [dens] = dens_poly3(poly3,t,s) |
| 2 |
% D=DENS_POLY3(P,T,S) |
| 3 |
% |
| 4 |
% Calculates in-situ density as approximated by the POLY3 method |
| 5 |
% used in the MITgcm. |
| 6 |
% P - coefficients read from file 'POLY3.COEFFS' using INI_POLY3 |
| 7 |
% T - potential temperature |
| 8 |
% S - salinity |
| 9 |
% |
| 10 |
% eg. |
| 11 |
% >> P=ini_poly3; |
| 12 |
% >> T=rdmds('T',100); |
| 13 |
% >> S=rdmds('S',100); |
| 14 |
% >> D=dens_poly3(P,T,S); |
| 15 |
% |
| 16 |
% or to work within a single model level |
| 17 |
% >> D=dens_poly3(P(3),T(:,:,3),S(:,:,3)); |
| 18 |
|
| 19 |
% $Header: $ |
| 20 |
% $Name: $ |
| 21 |
|
| 22 |
if size(t) ~= size(s) |
| 23 |
error('T and S must be the same shape and size') |
| 24 |
end |
| 25 |
%if size(t,ndims(t)) ~= size(poly3,2) |
| 26 |
% error('Last dimension of T and S must be the number of levels in P') |
| 27 |
%end |
| 28 |
|
| 29 |
n=size(t); |
| 30 |
nz=size(poly3,2); |
| 31 |
|
| 32 |
t=reshape(t,[prod(size(t))/nz nz]); |
| 33 |
s=reshape(s,[prod(size(t))/nz nz]); |
| 34 |
|
| 35 |
for k=1:nz, |
| 36 |
tRef=poly3(k).t; |
| 37 |
sRef=poly3(k).s; |
| 38 |
dRef=poly3(k).dens; |
| 39 |
tp=t(:,k)-tRef; |
| 40 |
sp=s(:,k)-sRef; |
| 41 |
|
| 42 |
deltaSig= ... |
| 43 |
poly3(k).coeffs(1) .*tp ... |
| 44 |
+poly3(k).coeffs(2) .*sp ... |
| 45 |
+poly3(k).coeffs(3) .*tp.*tp ... |
| 46 |
+poly3(k).coeffs(4) .*tp .*sp ... |
| 47 |
+poly3(k).coeffs(5) .*sp.*sp ... |
| 48 |
+poly3(k).coeffs(6) .*tp.*tp.*tp ... |
| 49 |
+poly3(k).coeffs(7) .*tp.*tp .*sp ... |
| 50 |
+poly3(k).coeffs(8) .*tp .*sp.*sp ... |
| 51 |
+poly3(k).coeffs(9) .*sp.*sp.*sp ... |
| 52 |
; |
| 53 |
dens(:,k)=deltaSig+dRef; |
| 54 |
end |
| 55 |
|
| 56 |
dens=reshape(dens,n); |
| 57 |
dens( find(t==0 & s==0) )=0; |
Parent Directory
| 
Revision Log
| 
Revision Graph