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