1 |
function [n2] = dens_poly3(poly3,t,s,dzc) |
2 |
% D=DENS_POLY3(P,T,S,dzc) |
3 |
% |
4 |
% Calculates Brunt-Vasala frequecy 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 |
% >> N2=n2_poly3(P,T,S); |
15 |
% |
16 |
% or to work within a single model level |
17 |
% >> N2=n2_poly3(P,T(:,:,3),S(:,:,3),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 |
n2=0*t; |
33 |
for k=2:nz, |
34 |
d1=dens_poly3(poly3(k-1),t(:,k-1),s(:,k-1)); |
35 |
d2=dens_poly3(poly3(k-1),t(:,k ),s(:,k )); |
36 |
dRdz=(d1-d2)/dzc(k); |
37 |
dRdz( find( t(:,k)==0 ) )=0; |
38 |
n2(:,k)=-9.81*dRdz; |
39 |
end |
40 |
|
41 |
n2=reshape(n2,n); |