UPP (upp-srw-2.2.0)
Loading...
Searching...
No Matches
CALHEL.f
Go to the documentation of this file.
1
3!
43 SUBROUTINE calhel(DEPTH,UST,VST,HELI,USHR1,VSHR1,USHR6,VSHR6)
44
45!
46 use vrbls3d, only: zmid, uh, vh, u, v, zint
47 use vrbls2d, only: fis, u10, v10
48 use masks, only: lmv, gdlat
49 use params_mod, only: g
50 use lookup_mod, only: itb,jtb,itbq,jtbq
51 use ctlblk_mod, only: jsta, jend, jsta_m, jend_m, jsta_2l, jend_2u, &
52 lm, im, jm, me, spval, &
53 ista, iend, ista_m, iend_m, ista_2l, iend_2u
54 use gridspec_mod, only: gridtype
55!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
56 implicit none
57!
58 real,PARAMETER :: P150=15000.0,p300=30000.0,s15=15.0
59 real,PARAMETER :: D3000=3000.0,pi6=0.5235987756,pi9=0.34906585
60 real,PARAMETER :: D5500=5500.0,d6000=6000.0,d7000=7000.0
61 real,PARAMETER :: D500=500.0
62! CRA
63 real,PARAMETER :: D1000=1000.0
64 real,PARAMETER :: D1500=1500.0
65! CRA
66
67!
68! DECLARE VARIABLES
69!
70 real,intent(in) :: DEPTH(2)
71 REAL,dimension(ista_2l:iend_2u,jsta_2l:jend_2u), intent(out) :: UST,VST
72 REAL,dimension(ista_2l:iend_2u,jsta_2l:jend_2u,2),intent(out) :: HELI
73!
74 real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u) :: HTSFC, UST6, VST6, UST5, VST5, &
75 ust1, vst1, ushr1, vshr1, &
76 ushr6, vshr6, u1, v1, u2, v2, &
77 hgt1, hgt2, umean, vmean
78! REAL HTSFC(IM,JM)
79!
80! REAL UST6(IM,JM),VST6(IM,JM)
81! REAL UST5(IM,JM),VST5(IM,JM)
82! REAL UST1(IM,JM),VST1(IM,JM)
83! CRA
84! REAL USHR1(IM,JM),VSHR1(IM,JM),USHR6(IM,JM),VSHR6(IM,JM)
85! REAL U1(IM,JM),V1(IM,JM),U2(IM,JM),V2(IM,JM)
86! REAL HGT1(IM,JM),HGT2(IM,JM),UMEAN(IM,JM),VMEAN(IM,JM)
87! CRA
88
89 integer, dimension(ista_2l:iend_2u,jsta_2l:jend_2u) :: COUNT6, COUNT5, COUNT1, L1, L2
90! INTEGER COUNT6(IM,JM),COUNT5(IM,JM),COUNT1(IM,JM)
91! CRA
92! INTEGER L1(IM,JM),L2(IM,JM)
93! CRA
94
95 INTEGER IVE(JM),IVW(JM)
96 integer I,J,IW,IE,JS,JN,JVN,JVS,L,N,lv
97 integer ISTART,ISTOP,JSTART,JSTOP
98 real Z2,DZABV,UMEAN5,VMEAN5,UMEAN1,VMEAN1,UMEAN6,VMEAN6, &
99 denom,z1,z3,dz,dz1,dz2,du1,du2,dv1,dv2
100
101!
102!****************************************************************
103! START CALHEL HERE
104!
105! INITIALIZE ARRAYS.
106!
107!$omp parallel do private(i,j)
108 DO j=jsta,jend
109 DO i=ista,iend
110 ust(i,j) = 0.0
111 vst(i,j) = 0.0
112 heli(i,j,1) = 0.0
113 heli(i,j,2) = 0.0
114 ust1(i,j) = 0.0
115 vst1(i,j) = 0.0
116 ust5(i,j) = 0.0
117 vst5(i,j) = 0.0
118 ust6(i,j) = 0.0
119 vst6(i,j) = 0.0
120 count6(i,j) = 0
121 count5(i,j) = 0
122 count1(i,j) = 0
123! CRA
124 ushr1(i,j) = 0.0
125 vshr1(i,j) = 0.0
126 ushr6(i,j) = 0.0
127 vshr6(i,j) = 0.0
128 u1(i,j) = 0.0
129 u2(i,j) = 0.0
130 v1(i,j) = 0.0
131 v2(i,j) = 0.0
132 umean(i,j) = 0.0
133 vmean(i,j) = 0.0
134 hgt1(i,j) = 0.0
135 hgt2(i,j) = 0.0
136 l1(i,j) = 0
137 l2(i,j) = 0
138! CRA
139
140 ENDDO
141 ENDDO
142 IF(gridtype == 'E')THEN
143 jvn = 1
144 jvs = -1
145 do j=jsta,jend
146 ive(j) = mod(j,2)
147 ivw(j) = ive(j)-1
148 enddo
149 istart = ista_m
150 istop = iend_m
151 jstart = jsta_m
152 jstop = jend_m
153 ELSE IF(gridtype == 'B')THEN
154 jvn = 1
155 jvs = 0
156 do j=jsta,jend
157 ive(j)=1
158 ivw(j)=0
159 enddo
160 istart = ista_m
161 istop = iend_m
162 jstart = jsta_m
163 jstop = jend_m
164 ELSE
165 jvn = 0
166 jvs = 0
167 do j=jsta,jend
168 ive(j) = 0
169 ivw(j) = 0
170 enddo
171 istart = ista
172 istop = iend
173 jstart = jsta
174 jstop = jend
175 END IF
176!
177! LOOP OVER HORIZONTAL GRID.
178!
179! CALL EXCH(RES(1,jsta_2l)
180! CALL EXCH(PD()
181
182! DO L = 1,LP1
183! CALL EXCH(ZINT(1,jsta_2l,L))
184! END DO
185!
186!!$omp parallel do private(htsfc,ie,iw)
187 IF(gridtype /= 'A') CALL exch(fis(ista_2l:iend_2u,jsta_2l:jend_2u))
188 DO l = 1,lm
189 IF(gridtype /= 'A') CALL exch(zmid(ista_2l:iend_2u,jsta_2l:jend_2u,l))
190 DO j=jstart,jstop
191 DO i=istart,istop
192 ie = i+ive(j)
193 iw = i+ivw(j)
194 jn = j+jvn
195 js = j+jvs
196!mp PDSLVK=(PD(IW,J)*RES(IW,J)+PD(IE,J)*RES(IE,J)+
197!mp 1 PD(I,J+1)*RES(I,J+1)+PD(I,J-1)*RES(I,J-1))*0.25
198!mp PSFCK=AETA(LMV(I,J))*PDSLVK+PT
199 IF (gridtype=='B')THEN
200 htsfc(i,j) = (0.25/g)*(fis(iw,j)+fis(ie,j)+fis(i,jn)+fis(ie,jn))
201!
202! COMPUTE MASS WEIGHTED MEAN WIND IN THE 0-6 KM LAYER, THE
203! 0-0.5 KM LAYER, AND THE 5.5-6 KM LAYER
204!
205 z2 = 0.25*(zmid(iw,j,l)+zmid(ie,j,l)+zmid(i,jn,l)+zmid(ie,jn,l))
206 ELSE
207 htsfc(i,j) = (0.25/g)*(fis(iw,j)+fis(ie,j)+fis(i,jn)+fis(i,js))
208!
209! COMPUTE MASS WEIGHTED MEAN WIND IN THE 0-6 KM LAYER, THE
210! 0-0.5 KM LAYER, AND THE 5.5-6 KM LAYER
211!
212 z2 = 0.25*(zmid(iw,j,l)+zmid(ie,j,l)+zmid(i,jn,l)+zmid(i,js,l))
213 END IF
214 dzabv = z2-htsfc(i,j)
215
216 lv = nint(lmv(i,j))
217 IF (dzabv <= d6000 .AND. l <= lv) THEN
218 ust6(i,j) = ust6(i,j) + uh(i,j,l)
219 vst6(i,j) = vst6(i,j) + vh(i,j,l)
220 count6(i,j) = count6(i,j) + 1
221 ENDIF
222
223 IF (dzabv < d6000 .AND. dzabv >= d5500 .AND. l <= lv) THEN
224 ust5(i,j) = ust5(i,j) + uh(i,j,l)
225 vst5(i,j) = vst5(i,j) + vh(i,j,l)
226 count5(i,j) = count5(i,j) + 1
227 ENDIF
228
229 IF (dzabv < d500 .AND. l <= lv) THEN
230 ust1(i,j) = ust1(i,j) + uh(i,j,l)
231 vst1(i,j) = vst1(i,j) + vh(i,j,l)
232 count1(i,j) = count1(i,j) + 1
233 ENDIF
234! CRA
235 IF (dzabv >= d1000 .AND. dzabv <= d1500 .AND. l <= lv) THEN
236 u2(i,j) = u(i,j,l)
237 v2(i,j) = v(i,j,l)
238 hgt2(i,j) = dzabv
239 l2(i,j) = l
240 ENDIF
241
242 IF (dzabv >= d500 .AND. dzabv < d1000 .AND. &
243 l <= lv .AND. l1(i,j) <= l2(i,j)) THEN
244 u1(i,j) = u(i,j,l)
245 v1(i,j) = v(i,j,l)
246 hgt1(i,j) = dzabv
247 l1(i,j) = l
248 ENDIF
249! CRA
250
251 ENDDO
252 ENDDO
253 ENDDO
254!
255! CASE WHERE THERE IS NO LEVEL WITH HEIGHT BETWEEN 5500 AND 6000
256
257 DO j=jstart,jstop
258 DO i=istart,istop
259 IF (count5(i,j) == 0) THEN
260 lloop: DO l=lm,1,-1
261 ie=i+ive(j)
262 iw=i+ivw(j)
263 jn=j+jvn
264 js=j+jvs
265 IF (gridtype=='B')THEN
266 z2 = 0.25*(zmid(iw,j,l)+zmid(ie,j,l)+zmid(i,jn,l)+zmid(ie,jn,l))
267 ELSE
268 z2 = 0.25*(zmid(iw,j,l)+zmid(ie,j,l)+zmid(i,jn,l)+zmid(i,js,l))
269 END IF
270
271 dzabv=z2-htsfc(i,j)
272 IF (dzabv < d7000 .AND. dzabv >= d6000) THEN
273 ust5(i,j) = ust5(i,j) + uh(i,j,l)
274 vst5(i,j) = vst5(i,j) + vh(i,j,l)
275 count5(i,j) = 1
276 exit lloop
277 ENDIF
278 enddo lloop
279 ENDIF
280 ENDDO
281 ENDDO
282
283!
284!$omp parallel do private(i,j,umean6,vmean6,umean5,vmean5,umean1,vmean1,denom)
285
286 DO j=jstart,jstop
287 DO i=istart,istop
288 IF (count6(i,j) > 0 .AND. count1(i,j) > 0 .AND. count5(i,j) > 0) THEN
289 umean5 = ust5(i,j) / count5(i,j)
290 vmean5 = vst5(i,j) / count5(i,j)
291 umean1 = ust1(i,j) / count1(i,j)
292 vmean1 = vst1(i,j) / count1(i,j)
293 umean6 = ust6(i,j) / count6(i,j)
294 vmean6 = vst6(i,j) / count6(i,j)
295
296!
297! COMPUTE STORM MOTION VECTOR
298! IT IS DEFINED AS 7.5 M/S TO THE RIGHT OF THE 0-6 KM MEAN
299! WIND CONSTRAINED ALONG A LINE WHICH IS BOTH PERPENDICULAR
300! TO THE 0-6 KM MEAN VERTICAL WIND SHEAR VECTOR AND PASSES
301! THROUGH THE 0-6 KM MEAN WIND. THE WIND SHEAR VECTOR IS
302! SET AS THE DIFFERENCE BETWEEN THE 5.5-6 KM WIND (THE HEAD
303! OF THE SHEAR VECTOR) AND THE 0-0.5 KM WIND (THE TAIL).
304! THIS IS FOR THE RIGHT-MOVING CASE; WE IGNORE THE LEFT MOVER.
305
306! CRA
307 ushr6(i,j) = umean5 - umean1
308 vshr6(i,j) = vmean5 - vmean1
309
310 denom = ushr6(i,j)*ushr6(i,j)+vshr6(i,j)*vshr6(i,j)
311 IF (denom /= 0.0) THEN
312 IF (gdlat(i,j) >= 0) THEN
313 ust(i,j) = umean6 + (7.5*vshr6(i,j)/sqrt(denom))
314 vst(i,j) = vmean6 - (7.5*ushr6(i,j)/sqrt(denom))
315 ELSE IF (gdlat(i,j) < 0) THEN
316 ust(i,j) = umean6 - (7.5*vshr6(i,j)/sqrt(denom))
317 vst(i,j) = vmean6 + (7.5*ushr6(i,j)/sqrt(denom))
318 ENDIF
319 ELSE
320 ust(i,j) = 0
321 vst(i,j) = 0
322 ENDIF
323 ELSE
324 ust(i,j) = 0.0
325 vst(i,j) = 0.0
326 ushr6(i,j) = 0.0
327 vshr6(i,j) = 0.0
328 ENDIF
329
330 IF(l1(i,j) > 0 .AND. l2(i,j) > 0) THEN
331 umean(i,j) = u1(i,j) + (d1000 - hgt1(i,j))*(u2(i,j) - &
332 u1(i,j))/(hgt2(i,j) - hgt1(i,j))
333 vmean(i,j) = v1(i,j) + (d1000 - hgt1(i,j))*(v2(i,j) - &
334 v1(i,j))/(hgt2(i,j) - hgt1(i,j))
335 ELSE IF(l1(i,j) > 0 .AND. l2(i,j) == 0) THEN
336 umean(i,j) = u1(i,j)
337 vmean(i,j) = v1(i,j)
338 ELSE IF(l1(i,j) == 0 .AND. l2(i,j) > 0) THEN
339 umean(i,j) = u2(i,j)
340 vmean(i,j) = u2(i,j)
341 ELSE
342 umean(i,j) = 0.0
343 vmean(i,j) = 0.0
344 ENDIF
345
346 IF(l1(i,j) > 0 .OR. l2(i,j) > 0) THEN
347 ushr1(i,j) = umean(i,j) - u10(i,j)
348 vshr1(i,j) = vmean(i,j) - v10(i,j)
349 ELSE
350 ushr1(i,j) = 0.0
351 vshr1(i,j) = 0.0
352 ENDIF
353! CRA
354
355!tgs USHR = UMEAN5 - UMEAN1
356! VSHR = VMEAN5 - VMEAN1
357
358! UST(I,J) = UMEAN6 + (7.5*VSHR/SQRT(USHR*USHR+VSHR*VSHR))
359! VST(I,J) = VMEAN6 - (7.5*USHR/SQRT(USHR*USHR+VSHR*VSHR))
360! ELSE
361! UST(I,J) = 0.0
362! VST(I,J) = 0.0
363! ENDIF
364
365 ENDDO
366 ENDDO
367!
368! COMPUTE STORM-RELATIVE HELICITY
369!
370!!$omp parallel do private(i,j,n,l,du1,du2,dv1,dv2,dz,dz1,dz2,dzabv,ie,iw,jn,js,z1,z2,z3)
371 DO n=1,2 ! for dfferent helicity depth
372 DO l = 2,lm-1
373 if(gridtype /= 'A')then
374 call exch(zint(1,jsta_2l,l))
375 call exch(zint(1,jsta_2l,l+1))
376 end if
377 DO j=jstart,jstop
378 DO i=istart,istop
379 iw=i+ivw(j)
380 ie=i+ive(j)
381 jn=j+jvn
382 js=j+jvs
383 IF (gridtype=='B')THEN
384 z2=0.25*(zmid(iw,j,l)+zmid(ie,j,l)+ &
385 zmid(i,jn,l)+zmid(ie,jn,l))
386 ELSE
387 z2=0.25*(zmid(iw,j,l)+zmid(ie,j,l)+ &
388 zmid(i,jn,l)+zmid(i,js,l))
389 END IF
390 dzabv=z2-htsfc(i,j)
391!
392 IF(dzabv < depth(n) .AND. l <= nint(lmv(i,j)))THEN
393 IF (gridtype=='B')THEN
394 z1 = 0.25*(zmid(iw,j,l+1)+zmid(ie,j,l+1)+ &
395 zmid(i,jn,l+1)+zmid(ie,jn,l+1))
396 z3 = 0.25*(zmid(iw,j,l-1)+zmid(ie,j,l-1)+ &
397 zmid(i,jn,l-1)+zmid(ie,jn,l-1))
398 dz = 0.25*((zint(iw,j,l)+zint(ie,j,l)+ &
399 zint(i,jn,l)+zint(ie,jn,l))- &
400 (zint(iw,j,l+1)+zint(ie,j,l+1)+ &
401 zint(i,jn,l+1)+zint(ie,jn,l+1)))
402 ELSE
403 z1 = 0.25*(zmid(iw,j,l+1)+zmid(ie,j,l+1)+ &
404 zmid(i,jn,l+1)+zmid(i,js,l+1))
405 z3 = 0.25*(zmid(iw,j,l-1)+zmid(ie,j,l-1)+ &
406 zmid(i,jn,l-1)+zmid(i,js,l-1))
407 dz = 0.25*((zint(iw,j,l)+zint(ie,j,l)+ &
408 zint(i,js,l)+zint(i,jn,l))- &
409 (zint(iw,j,l+1)+zint(ie,j,l+1)+ &
410 zint(i,js,l+1)+zint(i,jn,l+1)))
411 END IF
412 dz1 = z1-z2
413 dz2 = z2-z3
414 du1 = uh(i,j,l+1)-uh(i,j,l)
415 du2 = uh(i,j,l)-uh(i,j,l-1)
416 dv1 = vh(i,j,l+1)-vh(i,j,l)
417 dv2 = vh(i,j,l)-vh(i,j,l-1)
418 IF( vh(i,j,l) <spval.and.uh(i,j,l) <spval.and. &
419 vh(i,j,l+1)<spval.and.uh(i,j,l+1)<spval.and. &
420 vh(i,j,l-1)<spval.and.uh(i,j,l-1)<spval.and. &
421 vst(i,j) <spval.and.ust(i,j) <spval) &
422 heli(i,j,n) = ((vh(i,j,l)-vst(i,j))* &
423 (dz2*(du1/dz1)+dz1*(du2/dz2)) &
424 - (uh(i,j,l)-ust(i,j))* &
425 (dz2*(dv1/dz1)+dz1*(dv2/dz2))) &
426 *dz/(dz1+dz2)+heli(i,j,n)
427
428! if(i==im/2.and.j==(jsta+jend)/2)print*,'Debug Helicity',depth(N),l,dz1,dz2,du1, &
429! du2,dv1,dv2,ust(i,j),vst(i,j)
430 ENDIF
431 ENDDO
432 ENDDO
433 ENDDO
434 END DO ! end of different helicity depth
435!
436! END OF ROUTINE.
437!
438 RETURN
439 END
calhel
subroutine calhel(depth, ust, vst, heli, ushr1, vshr1, ushr6, vshr6)
This routine computes estimated storm motion and storm-relative environmental helicity.
Definition CALHEL.f:44