UPP  V11.0.0
 All Data Structures Files Functions Pages
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
286 30 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
masks
Definition: MASKS_mod.f:1
params_mod
Definition: params.F:1
lookup_mod
Definition: LOOKUP.f:1