/[MITgcm]/MITgcm_contrib/gael/profilesMatlabProcessing/profiles_misc/sw_ptmp.m
ViewVC logotype

Contents of /MITgcm_contrib/gael/profilesMatlabProcessing/profiles_misc/sw_ptmp.m

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.1 - (show annotations) (download)
Thu Oct 21 19:40:16 2010 UTC (14 years, 8 months ago) by gforget
Branch: MAIN
CVS Tags: checkpoint65x, checkpoint65r, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, HEAD
up-to-date matlab codes to process in-situ
data sets to MITgcm/pkg/profiles input format.

1 function PT = sw_ptmp(S,T,P,PR)
2
3 % SW_PTMP Potential temperature
4 %===========================================================================
5 % SW_PTMP $Revision: 1.3 $ $Date: 1994/10/10 05:45:13 $
6 % Copyright (C) CSIRO, Phil Morgan 1992.
7 %
8 % USAGE: ptmp = sw_ptmp(S,T,P,PR)
9 %
10 % DESCRIPTION:
11 % Calculates potential temperature as per UNESCO 1983 report.
12 %
13 % INPUT: (all must have same dimensions)
14 % S = salinity [psu (PSS-78) ]
15 % T = temperature [degree C (IPTS-68)]
16 % P = pressure [db]
17 % PR = Reference pressure [db]
18 % (P & PR may have dims 1x1, mx1, 1xn or mxn for S(mxn) )
19 %
20 % OUTPUT:
21 % ptmp = Potential temperature relative to PR [degree C (IPTS-68)]
22 %
23 % AUTHOR: Phil Morgan 92-04-06 (morgan@ml.csiro.au)
24 %
25 % DISCLAIMER:
26 % This software is provided "as is" without warranty of any kind.
27 % See the file sw_copy.m for conditions of use and licence.
28 %
29 % REFERENCES:
30 % Fofonoff, P. and Millard, R.C. Jr
31 % Unesco 1983. Algorithms for computation of fundamental properties of
32 % seawater, 1983. _Unesco Tech. Pap. in Mar. Sci._, No. 44, 53 pp.
33 % Eqn.(31) p.39
34 %
35 % Bryden, H. 1973.
36 % "New Polynomials for thermal expansion, adiabatic temperature gradient
37 % and potential temperature of sea water."
38 % DEEP-SEA RES., 1973, Vol20,401-408.
39 %=========================================================================
40
41 % CALLER: general purpose
42 % CALLEE: sw_adtg.m
43
44 %-------------
45 % CHECK INPUTS
46 %-------------
47 if nargin ~= 4
48 error('sw_ptmp.m: Must pass 4 parameters ')
49 end %if
50
51 % CHECK S,T,P dimensions and verify consistent
52 [ms,ns] = size(S);
53 [mt,nt] = size(T);
54 [mp,np] = size(P);
55
56 [mpr,npr] = size(PR);
57
58
59 % CHECK THAT S & T HAVE SAME SHAPE
60 if (ms~=mt) | (ns~=nt)
61 error('check_stp: S & T must have same dimensions')
62 end %if
63
64 % CHECK OPTIONAL SHAPES FOR P
65 if mp==1 & np==1 % P is a scalar. Fill to size of S
66 P = P(1)*ones(ms,ns);
67 elseif np==ns & mp==1 % P is row vector with same cols as S
68 P = P( ones(1,ms), : ); % Copy down each column.
69 elseif mp==ms & np==1 % P is column vector
70 P = P( :, ones(1,ns) ); % Copy across each row
71 elseif mp==ms & np==ns % PR is a matrix size(S)
72 % shape ok
73 else
74 error('check_stp: P has wrong dimensions')
75 end %if
76 [mp,np] = size(P);
77
78
79 % CHECK OPTIONAL SHAPES FOR PR
80 if mpr==1 & npr==1 % PR is a scalar. Fill to size of S
81 PR = PR(1)*ones(ms,ns);
82 elseif npr==ns & mpr==1 % PR is row vector with same cols as S
83 PR = PR( ones(1,ms), : ); % Copy down each column.
84 elseif mpr==ms & npr==1 % P is column vector
85 PR = PR( :, ones(1,ns) ); % Copy across each row
86 elseif mpr==ms & npr==ns % PR is a matrix size(S)
87 % shape ok
88 else
89 error('check_stp: PR has wrong dimensions')
90 end %if
91 [mpr,npr] = size(PR);
92
93
94 % IF ALL ROW VECTORS ARE PASSED THEN LET US PRESERVE SHAPE ON RETURN.
95 Transpose = 0;
96 if mp == 1 % row vector
97 P = P(:);
98 T = T(:);
99 S = S(:);
100
101 PR = PR(:);
102
103 Transpose = 1;
104 end %if
105 %***check_stp
106
107 %------
108 % BEGIN
109 %------
110
111 % theta1
112 del_P = PR - P;
113 del_th = del_P.*sw_adtg(S,T,P);
114 th = T + 0.5*del_th;
115 q = del_th;
116
117 % theta2
118 del_th = del_P.*sw_adtg(S,th,P+0.5*del_P);
119 th = th + (1 - 1/sqrt(2))*(del_th - q);
120 q = (2-sqrt(2))*del_th + (-2+3/sqrt(2))*q;
121
122 % theta3
123 del_th = del_P.*sw_adtg(S,th,P+0.5*del_P);
124 th = th + (1 + 1/sqrt(2))*(del_th - q);
125 q = (2 + sqrt(2))*del_th + (-2-3/sqrt(2))*q;
126
127 % theta4
128 del_th = del_P.*sw_adtg(S,th,P+del_P);
129 PT = th + (del_th - 2*q)/6;
130
131 if Transpose
132 PT = PT';
133 end %if
134
135 return
136 %=========================================================================

  ViewVC Help
Powered by ViewVC 1.1.22