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