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