UPP (develop)
Loading...
Searching...
No Matches
AVIATION.f
Go to the documentation of this file.
1
3!
65 SUBROUTINE calllws(U,V,H,LLWS)
66
67!
68 USE vrbls2d, only: fis, u10, v10
69 use params_mod, only: gi
70 use ctlblk_mod, only: jsta, jend, im, jm, lsm, spval, ista, iend
71!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
72 implicit none
73!
74! DECLARE VARIABLES.
75!
76 REAL,DIMENSION(IM,JM,LSM),INTENT(IN) :: U,V,H
77 REAL,DIMENSION(IM,JM),INTENT(INOUT) :: LLWS
78 REAL :: Z1,Z2,HZ1,DH,U2,V2,W2,RT
79 INTEGER :: K1,K2
80 integer I,J,LP
81
82!***************************************************************
83!
84!
85
86 DO 100 j=jsta,jend
87 DO i=ista,iend
88
89 z1 = 10.0 + fis(i,j)*gi !Height of 10m levels geographic height (from sea level)
90
91 IF(z1<h(i,j,lsm)) THEN !First search location of 10m wind level
92 k1 = lsm + 1 !to see it is in which pressure levels
93 ELSE
94 DO lp = lsm,2,-1 !If not found, keep searching upward
95 IF(z1>=h(i,j,lp).AND.z1<h(i,j,lp-1)) THEN
96 k1 = lp
97 END IF
98 END DO
99 END IF
100
101 hz1 = h(i,j,k1-1) - z1 !Distance between K1-1 and 10m level
102
103 dh = 0.0
104
105 IF((hz1+10)>609.6) THEN !Then, search 2000ft(609.6m) location
106 u2= u10(i,j) + (u(i,j,k1-1)-u10(i,j))*599.6/hz1 !found it between K1-1 and K1, then linear
107 v2= v10(i,j) + (v(i,j,k1-1)-v10(i,j))*599.6/hz1 !interpolate to get wind at 2000ft U2,V2
108 z2= fis(i,j)*gi + 609.6
109 ELSE !otherwise, keep on search upward
110 DO lp = k1-1,2,-1
111 dh=dh+(h(i,j,lp-1) - h(i,j,lp))
112 IF((dh+hz1+10)>609.6) THEN !found the 2000ft level
113 z2=fis(i,j)*gi+609.6
114 rt=(z2-h(i,j,lp))/(h(i,j,lp-1)-h(i,j,lp))
115 u2=u(i,j,lp)+rt*(u(i,j,lp-1)-u(i,j,lp))
116 v2=v(i,j,lp)+rt*(v(i,j,lp-1)-v(i,j,lp))
117 k2=lp
118 exit
119 END IF
120 END DO
121 END IF
122
123!computer vector difference
124 llws(i,j) = spval
125 if(u10(i,j)<spval.and.v10(i,j)<spval) &
126 llws(i,j)=sqrt((u2-u10(i,j))**2+(v2-v10(i,j))**2)/ &
127 609.6 * 1.943*609.6 !unit: knot/2000ft
128 ENDDO
129
130100 CONTINUE
131
132 RETURN
133 END
134
160 SUBROUTINE calicing (T1,RH,OMGA, ICING)
161!$$$ SUBPROGRAM DOCUMENTATION BLOCK
162! . . .
163! SUBPROGRAM: CALICING COMPUTES In-Flight Icing
164! PRGRMMR: Binbin Zhou /NCEP/EMC DATE: 2005-08-16
165!
166! ABSTRACT:
167! This program computes the in-flight icing condition
168! with the T-RH-OMGA algorithm provided by S. Silberberg of
169! NCEP/AWC (improved new version)
170!
171! According to S. Silberberg, Icing happens in following
172! situation:
173! (1) -22C < T < 0C to
174! (2) RH > 70 %
175! (3) Ascent air, OMGA < 0
176! (4) Equivalent Potential Vorticity (EPV) < 0
177! (5) Cloud water if SLD (supercooled large droplet)
178!
179! Current version dosn't consider SLD, so cloud water
180! is not used. EPV computation is not available for current
181! NCEP/EMC models(NAM, WRF, RSM), so EPV is also not
182! used
183!
184! USAGE: CALL CALICING(T1,RH,OMGA,ICING)
185! INPUT ARGUMENT LIST:
186! T1 - TEMPERATURE (K)
187! RH - RELATIVE HUMIDITY (DECIMAL FORM)
188! OMGA - Vertical velocity (Pa/sec)
189!
190! OUTPUT ARGUMENT LIST:
191! ICING - ICING CONDITION (1 or 0)
192!
193! OUTPUT FILES:
194! NONE
195!
196! SUBPROGRAMS CALLED:
197! UTILITIES:
198! LIBRARY:
199! NONE
200!
201! ATTRIBUTES:
202! LANGUAGE: FORTRAN 90/77
203! MACHINE : BLUE AT NCEP
204!$$$
205!
206 use ctlblk_mod, only: jsta, jend, im, spval, ista, iend
207!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
208 implicit none
209!
210! DECLARE VARIABLES.
211!
212 REAL, DIMENSION(ista:iend,jsta:jend), INTENT(IN) :: T1,RH,OMGA
213 REAL, DIMENSION(ista:iend,jsta:jend), INTENT(INOUT) :: ICING
214 integer I,J
215!***************************************************************
216!
217!
218 DO j=jsta,jend
219 DO i=ista,iend
220 IF(omga(i,j)<spval.AND.t1(i,j)<spval.AND.rh(i,j)<spval) THEN
221 IF(omga(i,j) < 0.0 .AND. &
222 (t1(i,j) <= 273.0 .AND. t1(i,j) >= 251.0) &
223 .AND. rh(i,j) >= 70.0) THEN
224
225 icing(i,j) = 1.0
226 ELSE
227 icing(i,j) = 0.0
228 END IF
229 ELSE
230 icing(i,j) = spval
231 ENDIF
232 ENDDO
233 ENDDO
234
235 RETURN
236 END
237
264 SUBROUTINE calcat(U,V,H,U_OLD,V_OLD,H_OLD,CAT)
265 use masks, only: dx, dy
266 use ctlblk_mod, only: spval, jsta_2l, jend_2u, jsta_m, jend_m, &
267 im, jm, ista_2l, iend_2u, ista_m, iend_m, ista, iend
268 use gridspec_mod, only: gridtype
269!
270!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
271 implicit none
272
273!
274! DECLARE VARIABLES.
275!
276 REAL,DIMENSION(ista_2l:iend_2u,jsta_2l:jend_2u),INTENT(IN) :: U,V,H, &
277 u_old,v_old,h_old
278! INTEGER,INTENT(IN) :: L
279 REAL,DIMENSION(ista_2l:iend_2u,jsta_2l:jend_2u),INTENT(INOUT) :: CAT
280
281 REAL DSH, DST, DEF, CVG, VWS, TRBINDX
282 INTEGER IHE(JM),IHW(JM)
283 integer I,J
284 integer ISTART,ISTOP,JSTART,JSTOP
285 real VWS1,VWS2,VWS3,VWS4
286
287!***************************************************************
288!
289!
290 cat=spval
291 DO j=jsta_2l,jend_2u
292 IF(gridtype == 'A')THEN
293 ihw(j)=-1
294 ihe(j)=1
295 istart=ista_m
296 istop=iend_m
297 jstart=jsta_m
298 jstop=jend_m
299 ELSE IF(gridtype=='E')THEN
300 ihw(j)=-mod(j,2)
301 ihe(j)=ihw(j)+1
302 istart=ista_m
303 istop=iend_m
304 jstart=jsta_m
305 jstop=jend_m
306 ELSE IF(gridtype=='B')THEN
307 ihw(j)=-1
308 ihe(j)=0
309 istart=ista_m
310 istop=iend_m
311 jstart=jsta_m
312 jstop=jend_m
313 ELSE
314 print*,'no gridtype specified, exit calcat comp'
315 return
316 END IF
317 ENDDO
318
319 call exch(u)
320 call exch(v)
321 call exch(u_old)
322 call exch(v_old)
323 call exch(h)
324 call exch(h_old)
325
326 DO 100 j=jstart,jstop
327 DO i=istart,istop
328!
329 IF(gridtype=='B')THEN
330 IF(u(i,j)<spval.and.u(i,j-1)<spval.and.u(i-1,j)<spval.and.u(i-1,j-1)<spval.and.&
331 v(i,j)<spval.and.v(i,j-1)<spval.and.v(i-1,j)<spval.and.v(i-1,j-1)<spval)THEN
332!dsh=dv/dx+du/dy
333 dsh=(0.5*(v(i,j)+v(i,j-1))-0.5*(v(i-1,j)+v(i-1,j-1)))*10000./dx(i,j) &
334 +(0.5*(u(i,j)+u(i-1,j))-0.5*(u(i,j-1)+u(i-1,j-1)))*10000./dy(i,j)
335!dst=du/dx-dv/dy
336 dst =(0.5*(u(i,j)+u(i,j-1))-0.5*(u(i-1,j)+u(i-1,j-1)))*10000./dx(i,j) &
337 -(0.5*(v(i,j)+v(i-1,j))-0.5*(v(i,j-1)+v(i-1,j-1)))*10000./dy(i,j)
338 def = sqrt(dsh*dsh + dst*dst)
339
340!cvg=-(du/dx+dv/dy)
341 cvg = -((0.5*(u(i,j)+u(i,j-1))-0.5*(u(i-1,j)+u(i-1,j-1)))*10000./dx(i,j) &
342 +(0.5*(v(i,j)+v(i-1,j))-0.5*(v(i,j-1)+v(i-1,j-1)))*10000./dy(i,j))
343 ELSE
344 def = spval
345 cvg = spval
346 ENDIF
347 ELSE
348 IF(u(i,j+1)<spval.and.u(i,j-1)<spval.and.u(i+ihe(j),j)<spval.and.u(i+ihw(j),j)<spval.and.&
349 v(i,j+1)<spval.and.v(i,j-1)<spval.and.v(i+ihe(j),j)<spval.and.v(i+ihw(j),j)<spval)THEN
350!dsh=dv/dx+du/dy
351 dsh = (v(i+ihe(j),j) - v(i+ihw(j),j))*10000./(2*dx(i,j)) &
352 + (u(i,j+1) - u(i,j-1))*10000./(2*dy(i,j))
353
354!dst=du/dx-dv/dy
355 dst = (u(i+ihe(j),j) - u(i+ihw(j),j))*10000./(2*dx(i,j)) &
356 - (v(i,j+1) - v(i,j-1))*10000./(2*dy(i,j))
357
358 def = sqrt(dsh*dsh + dst*dst)
359
360!cvg=-(du/dx+dv/dy)
361 cvg = -( (u(i+ihe(j),j) - u(i+ihw(j),j))*10000./(2*dx(i,j)) &
362 +(v(i,j+1) - v(i,j-1))*10000./(2*dy(i,j)) )
363 ELSE
364 def = spval
365 cvg = spval
366 ENDIF
367 END IF
368
369 IF(gridtype == 'A')THEN
370!vws=d|U|/dz
371 IF(u_old(i,j)<spval.and.u(i,j)<spval.and.&
372 v_old(i,j)<spval.and.v(i,j)<spval.and.&
373 h_old(i,j)<spval.and.h(i,j)<spval)THEN
374 vws = ( sqrt(u_old(i,j)**2+v_old(i,j)**2 ) - &
375 sqrt(u(i,j)**2+v(i,j)**2 ) ) * &
376 1000.0/(h_old(i,j) - h(i,j))
377 ELSE
378 vws = spval
379 ENDIF
380 else IF(gridtype == 'E')THEN
381!vws=d|U|/dz
382 IF(u_old(i+ihe(j),j)<spval.and.u(i+ihe(j),j)<spval.and.&
383 v_old(i+ihe(j),j)<spval.and.v(i+ihe(j),j)<spval)THEN
384
385 vws1 = ( sqrt(u_old(i+ihe(j),j)**2+v_old(i+ihe(j),j)**2 ) -&
386 sqrt(u(i+ihe(j),j)**2+v(i+ihe(j),j)**2 ) )
387 ELSE
388 vws1 = spval
389 ENDIF
390!vws=d|U|/dz
391 IF(u_old(i+ihw(j),j)<spval.and.u(i+ihw(j),j)<spval.and.&
392 v_old(i+ihw(j),j)<spval.and.v(i+ihw(j),j)<spval)THEN
393 vws2 = ( sqrt(u_old(i+ihw(j),j)**2+v_old(i+ihw(j),j)**2 ) -&
394 sqrt(u(i+ihw(j),j)**2+v(i+ihw(j),j)**2 ) )
395 ELSE
396 vws2 = spval
397 ENDIF
398!vws=d|U|/dz
399 IF(u_old(i,j-1)<spval.and.u(i,j-1)<spval.and.&
400 v_old(i,j-1)<spval.and.v(i,j-1)<spval)THEN
401 vws3 = ( sqrt(u_old(i,j-1)**2+v_old(i,j-1)**2 ) - &
402 sqrt(u(i,j-1)**2+v(i,j-1)**2 ) )
403 ELSE
404 vws3 = spval
405 ENDIF
406!vws=d|U|/dz
407 IF(u_old(i,j+1)<spval.and.u(i,j+1)<spval.and.&
408 v_old(i,j+1)<spval.and.v(i,j+1)<spval)THEN
409 vws4 = ( sqrt(u_old(i,j+1)**2+v_old(i,j+1)**2 ) - &
410 sqrt(u(i,j+1)**2+v(i,j+1)**2 ) )
411 ELSE
412 vws4 = spval
413 ENDIF
414
415 IF(vws1<spval.and.vws2<spval.and.vws3<spval.and.vws4<spval.and.&
416 h_old(i,j)<spval.and.h(i,j)<spval)THEN
417 vws=1000.0*(vws1+vws2+vws3+vws4)/4.0/(h_old(i,j) - h(i,j))
418 ELSE
419 vws = spval
420 ENDIF
421 ELSE IF(gridtype == 'B')THEN
422 IF(u_old(i+ihe(j),j)<spval.and.u(i+ihe(j),j)<spval.and.&
423 v_old(i+ihe(j),j)<spval.and.v(i+ihe(j),j)<spval)THEN
424 vws1 = ( sqrt(u_old(i+ihe(j),j)**2+v_old(i+ihe(j),j)**2 ) -&
425 sqrt(u(i+ihe(j),j)**2+v(i+ihe(j),j)**2 ) )
426 ELSE
427 vws1 = spval
428 ENDIF
429!vws=d|U|/dz
430 IF(u_old(i+ihw(j),j)<spval.and.u(i+ihw(j),j)<spval.and.&
431 v_old(i+ihw(j),j)<spval.and.v(i+ihw(j),j)<spval)THEN
432 vws2 = ( sqrt(u_old(i+ihw(j),j)**2+v_old(i+ihw(j),j)**2 ) -&
433 sqrt(u(i+ihw(j),j)**2+v(i+ihw(j),j)**2 ) )
434 ELSE
435 vws2 = spval
436 ENDIF
437!vws=d|U|/dz
438 IF(u_old(i,j-1)<spval.and.u(i,j-1)<spval.and.&
439 v_old(i,j-1)<spval.and.v(i,j-1)<spval)THEN
440 vws3 = ( sqrt(u_old(i,j-1)**2+v_old(i,j-1)**2 ) - &
441 sqrt(u(i,j-1)**2+v(i,j-1)**2 ) )
442 ELSE
443 vws3 = spval
444 ENDIF
445!vws=d|U|/dz
446 IF(u_old(i-1,j-1)<spval.and.u(i-1,j-1)<spval.and.&
447 v_old(i-1,j-1)<spval.and.v(i-1,j-1)<spval)THEN
448 vws4 = ( sqrt(u_old(i-1,j-1)**2+v_old(i-1,j-1)**2 ) - &
449 sqrt(u(i-1,j-1)**2+v(i-1,j-1)**2 ) )
450 ELSE
451 vws4 = spval
452 ENDIF
453
454 IF(vws1<spval.and.vws2<spval.and.vws3<spval.and.vws4<spval.and.&
455 h_old(i,j)<spval.and.h(i,j)<spval)THEN
456 vws=1000.0*(vws1+vws2+vws3+vws4)/4.0/(h_old(i,j) - h(i,j))
457 ELSE
458 vws=spval
459 ENDIF
460 END IF
461
462 IF(vws<spval.and.def<spval.and.cvg<spval)THEN
463 trbindx = abs(vws)*(def + abs(cvg))
464
465 IF(trbindx<=4.) THEN
466 cat(i,j) = 0.0
467 ELSE IF(trbindx<=8.) THEN
468 cat(i,j)=1.0
469 ELSE IF(trbindx<=12.) THEN
470 cat(i,j)=2.0
471 ELSE
472 cat(i,j)=3.0
473 END IF
474 ELSE
475 cat(i,j)=spval
476 ENDIF
477 ENDDO
478
479100 CONTINUE
480
481 RETURN
482 END
483
496 SUBROUTINE calceiling (CLDZ,TCLD,CEILING)
497 USE vrbls2d, only: fis
498 use params_mod, only: small, gi
499 use ctlblk_mod, only: jsta, jend, spval, im, modelname, ista, iend
500!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
501 implicit none
502!
503! DECLARE VARIABLES.
504!
505 REAL, DIMENSION(ista:iend,jsta:jend), INTENT(IN) :: CLDZ, TCLD
506 REAL, DIMENSION(ista:iend,jsta:jend), INTENT(INOUT) :: CEILING
507 integer I,J
508!***************************************************************
509!
510!
511 DO j=jsta,jend
512 DO i=ista,iend
513 IF(abs(tcld(i,j)-spval) <= small) THEN
514 ceiling(i,j)=spval
515 ELSE IF(tcld(i,j) >= 50.) THEN
516 if(modelname == 'RAPR')then
517 ceiling(i,j) = cldz(i,j) - fis(i,j)*gi
518 else
519 ceiling(i,j) = cldz(i,j) ! for RAP/HRRR - FIS(I,J)*GI
520 endif
521 ELSE
522 ceiling(i,j) = 20000.0
523 END IF
524
525 IF(ceiling(i,j) < 0.0) ceiling(i,j)=20000.0
526
527 ENDDO
528 ENDDO
529
530 RETURN
531 END
532
550 SUBROUTINE calfltcnd (CEILING,FLTCND)
551 use vrbls2d, only: vis
552 use ctlblk_mod, only: jsta, jend, im, spval, ista, iend
553!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
554 implicit none
555!
556! DECLARE VARIABLES.
557!
558 REAL, DIMENSION(ista:iend,jsta:jend), INTENT(IN) :: CEILING
559 REAL, DIMENSION(ista:iend,jsta:jend), INTENT(INOUT) :: FLTCND
560 REAL CEIL,VISI
561 integer I,J
562!
563!***************************************************************
564!
565!
566 DO j=jsta,jend
567 DO i=ista,iend
568
569 IF(ceiling(i,j)<spval.and.vis(i,j)<spval)THEN
570 ceil = ceiling(i,j) * 3.2808 !from m -> feet
571 visi = vis(i,j) / 1609.0 !from m -> miles
572
573 IF(ceil<500.0 .OR. visi<1.0 ) THEN
574 fltcnd(i,j) = 1.0
575
576 ELSE IF( (ceil>=500.AND.ceil<1000.0) .OR. &
577 (visi>=1.0.AND.visi<3.0) ) THEN
578 fltcnd(i,j) = 2.0
579
580 ELSE IF( (ceil>=1000.AND.ceil<=3000.0) .OR. &
581 (visi>=3.0.AND.visi<=5.0) ) THEN
582 fltcnd(i,j) = 3.0
583
584 ELSE IF( ceil>3000.0 .OR. visi>5.0) THEN
585 fltcnd(i,j) = 4.0
586
587 END IF
588 ELSE
589 fltcnd(i,j) = spval
590 ENDIF
591 ENDDO
592 ENDDO
593
594 RETURN
595 END
subroutine calceiling(cldz, tcld, ceiling)
Computes ceiling.
Definition AVIATION.f:497
subroutine calfltcnd(ceiling, fltcnd)
Computes Ceiling.
Definition AVIATION.f:551
subroutine calllws(u, v, h, llws)
Computes Low Level Wind Shear (0-2000feet)
Definition AVIATION.f:66
subroutine calicing(t1, rh, omga, icing)
Computes In-Flight Icing.
Definition AVIATION.f:161
subroutine calcat(u, v, h, u_old, v_old, h_old, cat)
Computes Clear Air Turbulence Index.
Definition AVIATION.f:265
subroutine exch(a)
exch() Subroutine that exchanges one halo row.
Definition EXCH.f:27