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