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