1 | ; Name: bn2 |
---|
2 | ; ----- |
---|
3 | ; |
---|
4 | ; Purpose : |
---|
5 | ; -------- |
---|
6 | ; Compute the local Brunt-Vaisala frequency |
---|
7 | ; |
---|
8 | ; Method : |
---|
9 | ; ------- |
---|
10 | ; neos = 0 : UNESCO sea water properties |
---|
11 | ; The brunt-vaisala frequency is computed using the polynomial |
---|
12 | ; polynomial expression of McDougall (1987): |
---|
13 | ; N^2 = g * beta * ( alpha/beta*dk[ tn ] - dk[ sn ] )/e3w |
---|
14 | ; |
---|
15 | ; neos = 1 : linear equation of state (temperature only) |
---|
16 | ; N^2 = g * ralpha * dk[ tn ]/e3w |
---|
17 | ; |
---|
18 | ; neos = 2 : linear equation of state (temperature & salinity) |
---|
19 | ; N^2 = g * (ralpha * dk[ tn ] - rbeta * dk[ sn ] ) / e3w |
---|
20 | ; |
---|
21 | ; The use of potential density to compute N^2 introduces e r r o r |
---|
22 | ; in the sign of N^2 at great depths. We recommand the use of |
---|
23 | ; neos = 0, except for academical studies. |
---|
24 | ; |
---|
25 | ; N.B. N^2 is set to zero at the first level (JK=1) in inidtr |
---|
26 | ; and is never used at this level. |
---|
27 | ; |
---|
28 | ; Input : potential temperature tn, salinity sn |
---|
29 | ; ------ |
---|
30 | ; |
---|
31 | ; Output : brunt-vaisala frequency |
---|
32 | ; ------- |
---|
33 | ; |
---|
34 | ; References : |
---|
35 | ; ----------- |
---|
36 | ; McDougall, T. J., J. Phys. Oceanogr., 17, 1950-1964, 1987. |
---|
37 | ; |
---|
38 | ; Modifications : |
---|
39 | ; -------------- |
---|
40 | ; Original : 94-07 (G. Madec, M. Imbard) |
---|
41 | ; Additions : 97-07 (G. Madec) introduction of statement functions |
---|
42 | ; Idl version : 98-12 (M. Levy) |
---|
43 | ;---------------------------------------------------------------------- |
---|
44 | FUNCTION bn2, tn, sn, neos |
---|
45 | @common |
---|
46 | case n_params() of |
---|
47 | 0:begin |
---|
48 | print, 'rentrer nom des champs de temperature et/ou salinite' |
---|
49 | goto, sortie |
---|
50 | end |
---|
51 | 1:begin |
---|
52 | print, 'pas de salinite : valeur de neos imposee a 1' |
---|
53 | neos = 1 |
---|
54 | end |
---|
55 | 2:begin |
---|
56 | print, 'calcul Unesco par default' |
---|
57 | neos = 0 |
---|
58 | END |
---|
59 | 3: |
---|
60 | endcase |
---|
61 | domdef |
---|
62 | ; coordonnees z : on s'assure qu'on a autant de points W que de points |
---|
63 | ; T sur le domaine defini par domdef. |
---|
64 | ; on va se servir dans les calculs de e3w (e3) et gdepw (gdep) |
---|
65 | IF (nzt NE nzw ) THEN BEGIN |
---|
66 | print, 'le domaine doit etre choisi avec autant de points' |
---|
67 | print, 'sur niveaux T que sur les niveaux W pour calculer N2' |
---|
68 | z = -1 |
---|
69 | GOTO, sortie |
---|
70 | endif |
---|
71 | IF (nzt LE 1 ) THEN BEGIN |
---|
72 | print, 'il faut un minimum de deux points sur la verticale' |
---|
73 | print, ' pour calculer N2' |
---|
74 | z = -1 |
---|
75 | GOTO, sortie |
---|
76 | endif |
---|
77 | cdgrid = 'W' |
---|
78 | g = 9.81 |
---|
79 | grille, mask, glam, gphi, gdep, nx, ny,nz, premierx, premiery, premierz, dernierx, derniery, dernierz ,e1,e2,e3, TRI = tri |
---|
80 | zgde3w = reform( (fltarr(nx*ny)+1)#(g/e3), nx, ny, nz) |
---|
81 | ; |
---|
82 | varname='Brunt-Vaisala frequency' |
---|
83 | varunit='s-1' |
---|
84 | ;vargrid='W' |
---|
85 | ; |
---|
86 | ; mask tn and sn |
---|
87 | tn=tn*tmask |
---|
88 | sn=sn*tmask |
---|
89 | ; |
---|
90 | ; Interior points ( 2=< jk =< jpkm1 ) |
---|
91 | ; ---------------- |
---|
92 | ; |
---|
93 | case 1 of |
---|
94 | (neos EQ 0): begin |
---|
95 | ; ... UNESCO seawater equation of state |
---|
96 | ; ... temperature, salinity anomalie |
---|
97 | zzt = 0.5*( tn + shift( tn, 0, 0, 1 )) |
---|
98 | zzs = 0.5*( sn + shift( sn, 0, 0, 1) ) - 35.0 |
---|
99 | zzp = reform( (fltarr(nx*ny)+1)#gdep, nx , ny, nz) |
---|
100 | ; ... ratio alpha/beta and beta |
---|
101 | zalbet = fsalbt( zzt, zzs, zzp ) |
---|
102 | zbeta = fsbeta( zzt, zzs, zzp ) |
---|
103 | ; ... N^2 = g/e3w * beta * ( alpha/beta * dk[tn] - dk[sn]) |
---|
104 | z = zgde3w * zbeta * mask $ |
---|
105 | * ( zalbet * ( shift(tn, 0, 0, 1) - tn ) $ |
---|
106 | - ( shift(sn, 0, 0, 1) - sn ) ) |
---|
107 | end |
---|
108 | ( neos EQ 1 ): begin |
---|
109 | ; ... First Linear density formulation (function of temperature only) |
---|
110 | ; |
---|
111 | ralpha = 1.176e-4 |
---|
112 | z = zgde3w * ralpha * mask* (shift(tn, 0, 0, 1) - tn) |
---|
113 | |
---|
114 | end |
---|
115 | ( neos EQ 2 ): begin |
---|
116 | ; ... Second linear density formulation (function of temp. and salinity) |
---|
117 | ; |
---|
118 | ralpha = 2.e-4 |
---|
119 | rbeta = 0.001 |
---|
120 | z = zgde3w * mask $ |
---|
121 | * ( ralpha * (shift(tn, 0, 0, 1) - tn) $ |
---|
122 | - rbeta * (shift(sn, 0, 0, 1) - sn) ) |
---|
123 | end |
---|
124 | (neos NE 0) AND (neos NE 1) AND (neos NE 2): BEGIN |
---|
125 | print, 'neos doit etre egal a 0, 1 ou 2' |
---|
126 | GOTO, sortie |
---|
127 | end |
---|
128 | endcase |
---|
129 | ; |
---|
130 | ; |
---|
131 | ; first and last levels |
---|
132 | ; ------------------------ |
---|
133 | ; bn^2=0. at first vertical point and jk=jpk |
---|
134 | z(*, *, 0) = 0 |
---|
135 | IF (gdep(nz-1) eq gdepw(jpk-1)) THEN z(*, *, nz-1) = 0 |
---|
136 | sortie: |
---|
137 | return, z |
---|
138 | END |
---|