UPP  11.0.0
 All Data Structures Files Functions Variables Pages
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
297 30 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
Definition: MASKS_mod.f:1
subroutine calhel3(LLOW, LUPP, UST, VST, HELI)
Subroutine that computes storm relative helicity.
Definition: CALHEL3.f:56