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