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 |
if size(t) ~= size(s) |
20 |
error('T and S must be the same shape and size') |
21 |
end |
22 |
%if size(t,ndims(t)) ~= size(poly3,2) |
23 |
% error('Last dimension of T and S must be the number of levels in P') |
24 |
%end |
25 |
|
26 |
n=size(t); |
27 |
nz=size(poly3,2); |
28 |
|
29 |
t=reshape(t,[prod(size(t))/nz nz]); |
30 |
s=reshape(s,[prod(size(t))/nz nz]); |
31 |
|
32 |
for k=1:nz, |
33 |
tRef=poly3(k).t; |
34 |
sRef=poly3(k).s; |
35 |
dRef=poly3(k).dens; |
36 |
tp=t(:,k)-tRef; |
37 |
sp=s(:,k)-sRef; |
38 |
|
39 |
deltaSig= ... |
40 |
poly3(k).coeffs(1) .*tp ... |
41 |
+poly3(k).coeffs(2) .*sp ... |
42 |
+poly3(k).coeffs(3) .*tp.*tp ... |
43 |
+poly3(k).coeffs(4) .*tp .*sp ... |
44 |
+poly3(k).coeffs(5) .*sp.*sp ... |
45 |
+poly3(k).coeffs(6) .*tp.*tp.*tp ... |
46 |
+poly3(k).coeffs(7) .*tp.*tp .*sp ... |
47 |
+poly3(k).coeffs(8) .*tp .*sp.*sp ... |
48 |
+poly3(k).coeffs(9) .*sp.*sp.*sp ... |
49 |
; |
50 |
dens(:,k)=deltaSig+dRef; |
51 |
end |
52 |
|
53 |
dens=reshape(dens,n); |