My Project
 All Classes Files Functions Variables Macros
PVtheta.F
Go to the documentation of this file.
1  SUBROUTINE pvtheta(ilon,ilev,pucov,pvcov,pteta,
2  $ ztfi,zplay,zplev,
3  $ nbteta,theta,pvteta)
4  IMPLICIT none
5 
6 c=======================================================================
7 c
8 c Auteur: I. Musat
9 c -------
10 c
11 c Objet:
12 c ------
13 c
14 c *******************************************************************
15 c Calcul de la vorticite potentielle PVteta sur des iso-theta selon
16 c la methodologie du NCEP/NCAR :
17 c 1) on calcule la stabilite statique N**2=g/T*(dT/dz+g/cp) sur les
18 c niveaux du modele => N2
19 c 2) on interpole les vents, la temperature et le N**2 sur des isentropes
20 c (en fait sur des iso-theta) lineairement en log(theta) =>
21 c ucovteta, vcovteta, N2teta
22 c 3) on calcule la vorticite absolue sur des iso-theta => vorateta
23 c 4) on calcule la densite rho sur des iso-theta => rhoteta
24 c
25 c rhoteta = (T/theta)**(cp/R)*p0/(R*T)
26 c
27 c 5) on calcule la vorticite potentielle sur des iso-theta => PVteta
28 c
29 c PVteta = (vorateta * N2 * theta)/(g * rhoteta) ! en PVU
30 c
31 c NB: 1PVU=10**(-6) K*m**2/(s * kg)
32 c
33 c PVteta = vorateta * N2/(g**2 * rhoteta) ! en 1/(Pa*s)
34 c
35 c
36 c *******************************************************************
37 c
38 c
39 c Variables d'entree : ilon,ilev,pucov,pvcov,pteta,ztfi,zplay,zplev,nbteta,theta
40 c -> sur la grille dynamique
41 c Variable de sortie : PVteta
42 c -> sur la grille physique
43 c=======================================================================
44 
45 #include "dimensions.h"
46 #include "paramet.h"
47 c
48 c variables Input
49 c
50  INTEGER ilon, ilev
51  REAL pvcov(iip1,jjm,ilev)
52  REAL pucov(iip1,jjp1,ilev)
53  REAL pteta(iip1,jjp1,ilev)
54  REAL ztfi(ilon,ilev)
55  REAL zplay(ilon,ilev), zplev(ilon,ilev+1)
56  INTEGER nbteta
57  REAL theta(nbteta)
58 c
59 c variable Output
60 c
61  REAL pvteta(ilon,nbteta)
62 c
63 c variables locales
64 c
65  INTEGER i, j, l, ig0
66  REAL ssum
67  REAL teta(ilon, ilev)
68  REAL ptetau(ip1jmp1, ilev), ptetav(ip1jm, ilev)
69  REAL ucovteta(ip1jmp1,ilev), vcovteta(ip1jm,ilev)
70  REAL n2(ilon,ilev-1), n2teta(ilon,nbteta)
71  REAL ztfiteta(ilon,nbteta)
72  REAL rhoteta(ilon,nbteta)
73  REAL vorateta(iip1,jjm,nbteta)
74  REAL voratetafi(ilon,nbteta), vorpol(iim)
75 c
76 #include "comgeom2.h"
77 #include "comconst.h"
78 #include "comvert.h"
79 c
80 c projection teta sur la grille physique
81 c
82  DO l=1,llm
83  teta(1,l) = pteta(1,1,l)
84  ig0 = 2
85  DO j = 2, jjm
86  DO i = 1, iim
87  teta(ig0,l) = pteta(i,j,l)
88  ig0 = ig0 + 1
89  ENDDO
90  ENDDO
91  teta(ig0,l) = pteta(1,jjp1,l)
92  ENDDO
93 c
94 c calcul pteta sur les grilles U et V
95 c
96  DO l=1, llm
97  DO j=1, jjp1
98  DO i=1, iip1
99  ig0=i+(j-1)*iip1
100  ptetau(ig0,l)=pteta(i,j,l)
101  ENDDO !i
102  ENDDO !j
103  DO j=1, jjm
104  DO i=1, iip1
105  ig0=i+(j-1)*iip1
106  ptetav(ig0,l)=0.5*(pteta(i,j,l)+pteta(i,j+1,l))
107  ENDDO !i
108  ENDDO !j
109  ENDDO !l
110 c
111 c projection pucov, pvcov sur une surface de theta constante
112 c
113  DO l=1, nbteta
114 cIM 1rout CALL tetaleveli1j1(ip1jmp1,llm,.true.,ptetau,theta(l),
115  CALL tetalevel(ip1jmp1,llm,.true.,ptetau,theta(l),
116  . pucov,ucovteta(:,l))
117 cIM 1rout CALL tetaleveli1j(ip1jm,llm,.true.,ptetav,theta(l),
118  CALL tetalevel(ip1jm,llm,.true.,ptetav,theta(l),
119  . pvcov,vcovteta(:,l))
120  ENDDO !l
121 c
122 c calcul vorticite absolue sur une iso-theta : vorateta
123 c
124  CALL tourabs(nbteta,vcovteta,ucovteta,vorateta)
125 c
126 c projection vorateta sur la grille physique => voratetafi
127 c
128  DO l=1,nbteta
129  DO j=2,jjm
130  ig0=1+(j-2)*iim
131  DO i=1,iim
132  voratetafi(ig0+i+1,l) = vorateta( i ,j-1,l) * alpha4(i+1,j) +
133  $ vorateta(i+1,j-1,l) * alpha1(i+1,j) +
134  $ vorateta(i ,j ,l) * alpha3(i+1,j) +
135  $ vorateta(i+1,j ,l) * alpha2(i+1,j)
136  ENDDO
137  voratetafi(ig0 +1,l) = voratetafi(ig0 +1+ iim,l)
138  ENDDO
139  ENDDO
140 c
141  DO l=1,nbteta
142  DO i=1,iim
143  vorpol(i) = vorateta(i,1,l)*aire(i,1)
144  ENDDO
145  voratetafi(1,l)= ssum(iim,vorpol,1)/apoln
146  ENDDO
147 c
148  DO l=1,nbteta
149  DO i=1,iim
150  vorpol(i) = vorateta(i,jjm,l)* aire(i,jjm +1)
151  ENDDO
152  voratetafi(ilon,l)= ssum(iim,vorpol,1)/apols
153  ENDDO
154 c
155 c calcul N**2 sur la grille physique => N2
156 c
157  DO l=1, llm-1
158  DO i=1, ilon
159  n2(i,l) = (g**2 * zplay(i,l) *
160  $ (ztfi(i,l+1)-ztfi(i,l)) )/
161  $ (r*ztfi(i,l)*ztfi(i,l)*
162  $ (zplev(i,l)-zplev(i,l+1)) )+
163  $ (g**2)/(ztfi(i,l)*cpp)
164  ENDDO !i
165  ENDDO !l
166 c
167 c calcul N2 sur une iso-theta => N2teta
168 c
169  DO l=1, nbteta
170  CALL tetalevel(ilon,llm-1,.true.,teta,theta(l),
171  $ n2,n2teta(:,l))
172  CALL tetalevel(ilon,llm,.true.,teta,theta(l),
173  $ ztfi,ztfiteta(:,l))
174  ENDDO !l=1, nbteta
175 c
176 c calcul rho et PV sur une iso-theta : rhoteta, PVteta
177 c
178  DO l=1, nbteta
179  DO i=1, ilon
180  rhoteta(i,l)=(ztfiteta(i,l)/theta(l))**(cpp/r)*
181  $ (preff/(r*ztfiteta(i,l)))
182 c
183 c PVteta en PVU
184 c
185  pvteta(i,l)=(theta(l)*g*voratetafi(i,l)*n2teta(i,l))/
186  $ (g**2*rhoteta(i,l))
187 c
188 c PVteta en 1/(Pa*s)
189 c
190  pvteta(i,l)=(voratetafi(i,l)*n2teta(i,l))/
191  $ (g**2*rhoteta(i,l))
192  ENDDO !i
193  ENDDO !l
194 c
195  RETURN
196  END