UPP (develop)
Loading...
Searching...
No Matches
MISCLN.f
Go to the documentation of this file.
1
55!-----------------------------------------------------------------------------------------------------
57 SUBROUTINE miscln
58
59!
60!
61 use vrbls3d, only: pmid, uh, vh, t, zmid, zint, pint, alpint, q, omga
62 use vrbls3d, only: catedr,mwt,gtg, cit
63 use vrbls2d, only: pblh, cprate, fis, t500, t700, z500, z700,&
64 teql,ieql, cape,cin
65 use masks, only: lmh
66 use params_mod, only: d00, d50, h99999, h100, h1, h1m12, pq0, a2, a3, a4, &
67 rhmin, rgamog, tfrz, small, g
68 use ctlblk_mod, only: grib, cfld, fld_info, datapd, im, jsta, jend, jm, jsta_m, jend_m, &
69 nbnd, nbin_du, lm, htfd, spval, pthresh, nfd, petabnd, me,&
70 jsta_2l, jend_2u, modelname, submodelname, &
71 ista, iend, ista_m, iend_m, ista_2l, iend_2u, &
72 ifi_flight_levels, gtg_on
73 use rqstfld_mod, only: iget, lvls, id, iavblfld, lvlsxml
74 use grib2_module, only: pset
75 use upp_physics, only: fpvsnew,calrh_pw,calcape,calcape2,tvirtual
76 use gridspec_mod, only: gridtype
77!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
78 implicit none
79!
80! SET LOCAL PARAMETERS. MAKE SURE NFD AND NBND AGREE
81! WITH THE VALUES SET IN SUBROUTINES FDLVL AND BNDLYR,
82! RESPECTIVELY.
83 real,PARAMETER :: C2K=273.15
84 real,parameter :: con_rd =2.8705e+2 ! gas constant air (J/kg/K)
85 real,parameter :: con_rv =4.6150e+2 ! gas constant H2O
86 real,parameter :: con_eps =con_rd/con_rv
87 real,parameter :: con_epsm1 =con_rd/con_rv-1
88 real,parameter :: cpthresh =0.000004
89 real,PARAMETER :: D1000=1000
90 real,PARAMETER :: D1500=1500
91 real,PARAMETER :: D2000=2000
92 real,PARAMETER :: HCONST=42000000.
93 real,PARAMETER :: K2C=273.16
94 REAL,PARAMETER :: DM9999=-9999.0
95
96!
97! DECLARE VARIABLES.
98!
99 LOGICAL NORTH, FIELD1,FIELD2, NEED_IFI
100 LOGICAL, dimension(ISTA:IEND,JSTA:JEND) :: DONE, DONE1
101
102 INTEGER, allocatable :: LVLBND(:,:,:),LB2(:,:)
103! INTEGER LVLBND(IM,JM,NBND),LB2(IM,JM),LPBL(IM,JM)
104
105 real,dimension(im,jm) :: GRID1, GRID2
106 real,dimension(ista:iend,jsta:jend) :: P1D, T1D, Q1D, U1D, V1D, SHR1D, Z1D, &
107 rh1d, egrid1, egrid2, egrid3, egrid4, &
108 egrid5, egrid6, egrid7, egrid8, &
109 mlcape,mlcin,mllcl,mucape,mucin, &
110 freezelvl,muq1d,slcl,the,maxthe
111 integer,dimension(ista:iend,jsta:jend) :: MAXTHEPOS
112 real, dimension(:,:,:),allocatable :: OMGBND, PWTBND, QCNVBND, &
113 pbnd, tbnd, qbnd, &
114 ubnd, vbnd, rhbnd, &
115 wbnd, t7d, q7d, &
116 u7d, v6d, p7d, &
117 icingfd,midcal
118
119 real, dimension(:,:),allocatable :: QM8510, RH4710, RH8498, &
120 rh4796, rh1847, ust, vst, &
121 rh3310, rh6610, rh3366, &
122 pw3310, rh4410, rh7294, &
123 rh4472, &
124 t78483, t89671, p78483, p89671
125
126 REAL, dimension(:,:,:),allocatable :: HELI
127 real, dimension(:,:), allocatable :: USHR1, VSHR1, USHR6, VSHR6, &
128 maxwp, maxwz, maxwu, maxwv, &
129 maxwt
130 INTEGER,dimension(:,:),allocatable :: LLOW,LUPP,LLOW_ZINT,IEQL_ZINT, &
131 z_midcal
132 REAL, dimension(:,:),allocatable :: CANGLE,ESHR,UVECT,VVECT,&
133 effust,effvst,fshr,htsfc,&
134 esrh,z_temp
135!
136 integer I,J,ii,jj,L,ITYPE,ISVALUE,LBND,ILVL,IFD,ITYPEFDLVL(NFD), &
137 iget1, iget2, iget3, llmh,imax,jmax,lmax
138 real DPBND,PKL1,PKU1,FAC1,FAC2,PL,TL,QL,QSAT,RHL,TVRL,TVRBLO, &
139 es1,es2,qs1,qs2,rh1,rh2,zsf,depth(2),work1,work2,work3, &
140 scintmp,mucapetmp,mucintmp,mllcltmp,eshrtmp,mlcapetmp,stp,&
141 fshrtmp,mlcintmp,slcltmp,lapse,ship
142
143 integer IE,IW,JN,JS,IVE(JM),IVW(JM),JVN,JVS
144 integer ISTART,ISTOP,JSTART,JSTOP
145 real dummy(ista:iend,jsta:jend)
146 integer idummy(ista:iend,jsta:jend)
147! NEW VARIABLES USED FOR EFFECTIVE LAYER
148 INTEGER,dimension(:,:),allocatable :: EL_BASE, EL_TOPS
149 LOGICAL,dimension(:,:),allocatable :: FOUND_BASE, FOUND_TOPS
150 INTEGER,dimension(:,:),allocatable :: L_THETAE_MAX
151 INTEGER,dimension(:,:),allocatable :: CAPE9, CINS9
152 CHARACTER(LEN=5) :: IM_CH, JSTA_CH, JEND_CH, ME_CH
153 CHARACTER(LEN=60) :: EFFL_FNAME
154 CHARACTER(LEN=60) :: EFFL_FNAME2
155 INTEGER :: IREC, IUNIT
156 INTEGER :: IREC2, IUNIT2
157 LOGICAL :: debugprint
158 INTEGER :: LLL
159 INTEGER :: LLCL_PAR, LEQL_PAR
160 REAL :: LMASK, PSFC, CAPE_PAR, CINS_PAR, LPAR0
161 REAL, DIMENSION(4) :: PARCEL0
162 REAL, DIMENSION(:), ALLOCATABLE :: TPAR_B, TPAR_T
163 REAL, DIMENSION(:), ALLOCATABLE :: TPAR_TMP
164 REAL, DIMENSION(:), ALLOCATABLE :: P_AMB, T_AMB, Q_AMB, ZINT_AMB
165 REAL, DIMENSION(:,:,:), ALLOCATABLE :: TPAR_BASE, TPAR_TOPS
166
167! Variables introduced to allow FD levels from control file - Y Mao
168 integer :: N,NFDCTL
169 REAL, allocatable :: HTFDCTL(:)
170 integer, allocatable :: ITYPEFDLVLCTL(:)
171 real, allocatable :: QIN(:,:,:,:), QFD(:,:,:,:)
172 character, allocatable :: QTYPE(:)
173
174 integer, parameter :: NFDMAX=10 ! Max number of fields with the same HTFDCTL
175 integer :: IDS(NFDMAX) ! All field IDs with the same HTFDCTL
176 integer :: nFDS ! How many fields with the same HTFDCTL in the control file
177 integer :: iID ! which field with HTFDCTL
178!
179!****************************************************************************
180! START MISCLN HERE.
181!
182 debugprint = .false.
183
184 need_ifi = iget(1007)>0 .or. iget(1008)>0 .or. iget(1009)>0 .or. iget(1010)>0
185
186 allocate(ushr1(ista_2l:iend_2u,jsta_2l:jend_2u),vshr1(ista_2l:iend_2u,jsta_2l:jend_2u), &
187 ushr6(ista_2l:iend_2u,jsta_2l:jend_2u),vshr6(ista_2l:iend_2u,jsta_2l:jend_2u))
188 allocate(ust(ista_2l:iend_2u,jsta_2l:jend_2u),vst(ista_2l:iend_2u,jsta_2l:jend_2u), &
189 heli(ista_2l:iend_2u,jsta_2l:jend_2u,2),fshr(ista_2l:iend_2u,jsta_2l:jend_2u))
190!
191! HELICITY AND STORM MOTION.
192 iget1 = iget(162)
193 iget2 = -1
194 iget3 = -1
195 if (iget1 > 0) then
196 iget2 = lvls(1,iget1)
197 iget3 = lvls(2,iget1)
198 endif
199 IF (iget1 > 0 .OR. iget(163) > 0 .OR. iget(164) > 0) THEN
200 depth(1) = 3000.0
201 depth(2) = 1000.0
202 CALL calhel(depth,ust,vst,heli,ushr1,vshr1,ushr6,vshr6)
203 IF (iget2 > 0) then
204!$omp parallel do private(i,j)
205 DO j=jsta,jend
206 DO i=ista,iend
207 grid1(i,j) = heli(i,j,1)
208 ENDDO
209 ENDDO
210 if(grib=='grib2') then
211 cfld=cfld+1
212 fld_info(cfld)%ifld=iavblfld(iget1)
213 fld_info(cfld)%lvl=lvlsxml(1,iget1)
214!$omp parallel do private(i,j,ii,jj)
215 do j=1,jend-jsta+1
216 jj = jsta+j-1
217 do i=1,iend-ista+1
218 ii = ista+i-1
219 datapd(i,j,cfld) = grid1(ii,jj)
220 enddo
221 enddo
222 endif
223 ENDIF
224
225 IF (iget3 > 0) then
226!$omp parallel do private(i,j)
227 DO j=jsta,jend
228 DO i=ista,iend
229 grid1(i,j) = heli(i,j,2)
230 ENDDO
231 ENDDO
232 if(grib=='grib2') then
233 cfld=cfld+1
234 fld_info(cfld)%ifld=iavblfld(iget1)
235 fld_info(cfld)%lvl=lvlsxml(2,iget1)
236!$omp parallel do private(i,j,ii,jj)
237 do j=1,jend-jsta+1
238 jj = jsta+j-1
239 do i=1,iend-ista+1
240 ii = ista+i-1
241 datapd(i,j,cfld) = grid1(ii,jj)
242 enddo
243 enddo
244 endif
245 ENDIF
246
247 IF (iget(163) > 0) THEN
248!$omp parallel do private(i,j)
249 DO j=jsta,jend
250 DO i=ista,iend
251 grid1(i,j) = ust(i,j)
252 ENDDO
253 ENDDO
254 if(grib=='grib2') then
255 cfld=cfld+1
256 fld_info(cfld)%ifld=iavblfld(iget(163))
257!$omp parallel do private(i,j,ii,jj)
258 do j=1,jend-jsta+1
259 jj = jsta+j-1
260 do i=1,iend-ista+1
261 ii = ista+i-1
262 datapd(i,j,cfld) = grid1(ii,jj)
263 enddo
264 enddo
265 endif
266 ENDIF
267 IF (iget(164) > 0) THEN
268!$omp parallel do private(i,j)
269 DO j=jsta,jend
270 DO i=ista,iend
271 grid1(i,j) = vst(i,j)
272 ENDDO
273 ENDDO
274 if(grib=='grib2') then
275 cfld=cfld+1
276 fld_info(cfld)%ifld=iavblfld(iget(164))
277!$omp parallel do private(i,j,ii,jj)
278 do j=1,jend-jsta+1
279 jj = jsta+j-1
280 do i=1,iend-ista+1
281 ii = ista+i-1
282 datapd(i,j,cfld) = grid1(ii,jj)
283 enddo
284 enddo
285 endif
286 ENDIF
287 ENDIF
288
289! UPDRAFT HELICITY
290
291 if (iget(427) > 0) THEN
292 CALL calupdhel(grid1(ista_2l:iend_2u,jsta_2l:jend_2u))
293 if(grib=='grib2') then
294 cfld=cfld+1
295 fld_info(cfld)%ifld=iavblfld(iget(427))
296!$omp parallel do private(i,j,ii,jj)
297 do j=1,jend-jsta+1
298 jj = jsta+j-1
299 do i=1,iend-ista+1
300 ii = ista+i-1
301 datapd(i,j,cfld) = grid1(ii,jj)
302 enddo
303 enddo
304 endif
305
306 ENDIF
307
308! CRA 0-1 KM AND 0-6 KM SHEAR
309
310 IF(iget(430) > 0 .OR. iget(431) > 0 .OR. iget(432) > 0 &
311 .OR. iget(433) > 0) THEN
312
313 depth = 6000.0
314 CALL calhel(depth,ust,vst,heli,ushr1,vshr1,ushr6,vshr6)
315! 0-6 km shear magnitude
316!$omp parallel do private(i,j)
317 DO j=jsta,jend
318 DO i=ista,iend
319 fshr(i,j) = sqrt(ushr6(i,j)**2+vshr6(i,j)**2)
320 ENDDO
321 ENDDO
322 IF(iget(430) > 0) THEN
323!$omp parallel do private(i,j)
324 DO j=jsta,jend
325 DO i=ista,iend
326 grid1(i,j) = ushr1(i,j)
327 ENDDO
328 ENDDO
329 if(grib=='grib2') then
330 cfld=cfld+1
331 fld_info(cfld)%ifld=iavblfld(iget(430))
332!$omp parallel do private(i,j,ii,jj)
333 do j=1,jend-jsta+1
334 jj = jsta+j-1
335 do i=1,iend-ista+1
336 ii = ista+i-1
337 datapd(i,j,cfld) = grid1(ii,jj)
338 enddo
339 enddo
340 endif
341 ENDIF
342 IF(iget(431) > 0) THEN
343!$omp parallel do private(i,j)
344 DO j=jsta,jend
345 DO i=ista,iend
346 grid1(i,j) = vshr1(i,j)
347 ENDDO
348 ENDDO
349 if(grib=='grib2') then
350 cfld=cfld+1
351 fld_info(cfld)%ifld=iavblfld(iget(431))
352!$omp parallel do private(i,j,ii,jj)
353 do j=1,jend-jsta+1
354 jj = jsta+j-1
355 do i=1,iend-ista+1
356 ii = ista+i-1
357 datapd(i,j,cfld) = grid1(ii,jj)
358 enddo
359 enddo
360 endif
361 ENDIF
362 IF(iget(432) > 0) THEN
363!$omp parallel do private(i,j)
364 DO j=jsta,jend
365 DO i=ista,iend
366 grid1(i,j) = ushr6(i,j)
367 ENDDO
368 ENDDO
369 if(grib=='grib2') then
370 cfld=cfld+1
371 fld_info(cfld)%ifld=iavblfld(iget(432))
372!$omp parallel do private(i,j,ii,jj)
373 do j=1,jend-jsta+1
374 jj = jsta+j-1
375 do i=1,iend-ista+1
376 ii = ista+i-1
377 datapd(i,j,cfld) = grid1(ii,jj)
378 enddo
379 enddo
380 endif
381 ENDIF
382 IF(iget(433) > 0) THEN
383!$omp parallel do private(i,j)
384 DO j=jsta,jend
385 DO i=ista,iend
386 grid1(i,j) = vshr6(i,j)
387 ENDDO
388 ENDDO
389 if(grib=='grib2') then
390 cfld=cfld+1
391 fld_info(cfld)%ifld=iavblfld(iget(433))
392!$omp parallel do private(i,j,ii,jj)
393 do j=1,jend-jsta+1
394 jj = jsta+j-1
395 do i=1,iend-ista+1
396 ii = ista+i-1
397 datapd(i,j,cfld) = grid1(ii,jj)
398 enddo
399 enddo
400 endif
401 ENDIF
402 ENDIF
403!
404 if (allocated(ushr1)) deallocate(ushr1)
405 if (allocated(vshr1)) deallocate(vshr1)
406 if (allocated(ushr6)) deallocate(ushr6)
407 if (allocated(vshr6)) deallocate(vshr6)
408 if (allocated(ust)) deallocate(ust)
409 if (allocated(vst)) deallocate(vst)
410 if (allocated(heli)) deallocate(heli)
411! CRA
412!
413!
414! ***BLOCK 1: TROPOPAUSE P, Z, T, U, V, AND WIND SHEAR.
415!
416 IF ((iget(054)>0).OR.(iget(055)>0).OR. &
417 (iget(056)>0).OR.(iget(057)>0).OR. &
418 (iget(177)>0).OR. &
419 (iget(058)>0).OR.(iget(108)>0) ) THEN
420! Chuang: Use GFS algorithm per Iredell's and DiMego's decision on unification
421!$omp parallel do private(i,j)
422 DO j=jsta,jend
423 DO i=ista,iend
424
425 if(pmid(i,j,1)<spval) then
426! INPUT
427 CALL tpause(lm,pmid(i,j,1:lm),uh(i,j,1:lm) &
428! INPUT
429 ,vh(i,j,1:lm),t(i,j,1:lm),zmid(i,j,1:lm) &
430! OUTPUT
431 ,p1d(i,j),u1d(i,j),v1d(i,j),t1d(i,j) &
432! OUTPUT
433 ,z1d(i,j),shr1d(i,j)) ! OUTPUT
434 else
435 p1d(i,j) = spval
436 u1d(i,j) = spval
437 v1d(i,j) = spval
438 t1d(i,j) = spval
439 z1d(i,j) = spval
440 shr1d(i,j) = spval
441 endif
442
443 END DO
444 END DO
445!
446! TROPOPAUSE PRESSURE.
447 IF (iget(054) > 0) THEN
448!$omp parallel do private(i,j)
449 DO j=jsta,jend
450 DO i=ista,iend
451 grid1(i,j) = p1d(i,j)
452 ENDDO
453 ENDDO
454 if(grib=='grib2') then
455 cfld=cfld+1
456 fld_info(cfld)%ifld=iavblfld(iget(054))
457!$omp parallel do private(i,j,ii,jj)
458 do j=1,jend-jsta+1
459 jj = jsta+j-1
460 do i=1,iend-ista+1
461 ii = ista+i-1
462 datapd(i,j,cfld) = grid1(ii,jj)
463 enddo
464 enddo
465 endif
466 ENDIF
467
468! ICAO HEIGHT OF TROPOPAUSE
469 IF (iget(399)>0) THEN
470 CALL icaoheight(p1d, grid1(ista:iend,jsta:jend))
471! print*,'sample TROPOPAUSE ICAO HEIGHTS',GRID1(im/2,(jsta+jend)/2)
472 if(grib=='grib2') then
473 cfld=cfld+1
474 fld_info(cfld)%ifld=iavblfld(iget(399))
475!$omp parallel do private(i,j,ii,jj)
476 do j=1,jend-jsta+1
477 jj = jsta+j-1
478 do i=1,iend-ista+1
479 ii = ista+i-1
480 datapd(i,j,cfld) = grid1(ii,jj)
481 enddo
482 enddo
483 endif
484 ENDIF
485
486! TROPOPAUSE HEIGHT.
487 IF (iget(177) > 0) THEN
488!$omp parallel do private(i,j)
489 DO j=jsta,jend
490 DO i=ista,iend
491 grid1(i,j) = z1d(i,j)
492 ENDDO
493 ENDDO
494 if(grib=='grib2') then
495 cfld=cfld+1
496 fld_info(cfld)%ifld=iavblfld(iget(177))
497!$omp parallel do private(i,j,ii,jj)
498 do j=1,jend-jsta+1
499 jj = jsta+j-1
500 do i=1,iend-ista+1
501 ii = ista+i-1
502 datapd(i,j,cfld) = grid1(ii,jj)
503 enddo
504 enddo
505 endif
506 ENDIF
507!
508! TROPOPAUSE TEMPERATURE.
509 IF (iget(055) > 0) THEN
510!$omp parallel do private(i,j)
511 DO j=jsta,jend
512 DO i=ista,iend
513 grid1(i,j) = t1d(i,j)
514 ENDDO
515 ENDDO
516 if(grib=='grib2') then
517 cfld=cfld+1
518 fld_info(cfld)%ifld=iavblfld(iget(055))
519!$omp parallel do private(i,j,ii,jj)
520 do j=1,jend-jsta+1
521 jj = jsta+j-1
522 do i=1,iend-ista+1
523 ii = ista+i-1
524 datapd(i,j,cfld) = grid1(ii,jj)
525 enddo
526 enddo
527 endif
528 ENDIF
529!
530! TROPOPAUSE POTENTIAL TEMPERATURE.
531 IF (iget(108) > 0) THEN
532 CALL calpot(p1d,t1d,grid1(ista:iend,jsta:jend))
533 if(grib=='grib2') then
534 cfld=cfld+1
535 fld_info(cfld)%ifld=iavblfld(iget(108))
536!$omp parallel do private(i,j,ii,jj)
537 do j=1,jend-jsta+1
538 jj = jsta+j-1
539 do i=1,iend-ista+1
540 ii = ista+i-1
541 datapd(i,j,cfld) = grid1(ii,jj)
542 enddo
543 enddo
544 endif
545 ENDIF
546!
547! TROPOPAUSE U WIND AND/OR V WIND.
548 IF ((iget(056) > 0).OR.(iget(057) > 0)) THEN
549!$omp parallel do private(i,j)
550 DO j=jsta,jend
551 DO i=ista,iend
552 grid1(i,j)=u1d(i,j)
553 grid2(i,j)=v1d(i,j)
554 ENDDO
555 ENDDO
556 if(grib=='grib2') then
557 if(iget(056)>0) then
558 cfld=cfld+1
559 fld_info(cfld)%ifld=iavblfld(iget(056))
560!$omp parallel do private(i,j,ii,jj)
561 do j=1,jend-jsta+1
562 jj = jsta+j-1
563 do i=1,iend-ista+1
564 ii = ista+i-1
565 datapd(i,j,cfld) = grid1(ii,jj)
566 enddo
567 enddo
568 endif
569 if(iget(057)>0) then
570 cfld=cfld+1
571 fld_info(cfld)%ifld=iavblfld(iget(057))
572!$omp parallel do private(i,j,ii,jj)
573 do j=1,jend-jsta+1
574 jj = jsta+j-1
575 do i=1,iend-ista+1
576 ii = ista+i-1
577 datapd(i,j,cfld) = grid2(ii,jj)
578 enddo
579 enddo
580 endif
581 endif
582 ENDIF
583!
584! TROPOPAUSE WIND SHEAR.
585 IF (iget(058) > 0) THEN
586!$omp parallel do private(i,j)
587 DO j=jsta,jend
588 DO i=ista,iend
589 grid1(i,j) = shr1d(i,j)
590 ENDDO
591 ENDDO
592 if(grib=='grib2') then
593 cfld=cfld+1
594 fld_info(cfld)%ifld=iavblfld(iget(058))
595!$omp parallel do private(i,j,ii,jj)
596 do j=1,jend-jsta+1
597 jj = jsta+j-1
598 do i=1,iend-ista+1
599 ii = ista+i-1
600 datapd(i,j,cfld) = grid1(ii,jj)
601 enddo
602 enddo
603 endif
604 ENDIF
605 ENDIF
606!
607!
608! ***BLOCK 2: MAX WIND LEVEL P, Z, U, AND V
609!
610! MAX WIND LEVEL CALCULATIONS
611 IF ((iget(173)>0) .OR. (iget(174)>0) .OR. &
612 (iget(175)>0) .OR. (iget(176)>0)) THEN
613
614 allocate(maxwp(ista:iend,jsta:jend), maxwz(ista:iend,jsta:jend), &
615 maxwu(ista:iend,jsta:jend), maxwv(ista:iend,jsta:jend),maxwt(ista:iend,jsta:jend))
616!$omp parallel do private(i,j)
617 DO j=jsta,jend
618 DO i=ista,iend
619 maxwp(i,j)=spval
620 maxwz(i,j)=spval
621 maxwu(i,j)=spval
622 maxwv(i,j)=spval
623 ENDDO
624 ENDDO
625
626! CALL CALMXW(MAXWP,MAXWZ,MAXWU,MAXWV,MAXWT)
627! Chuang: Use GFS algorithm per Iredell's and DiMego's decision on unification
628!$omp parallel do private(i,j)
629 DO j=jsta,jend
630 loopi:DO i=ista,iend
631 DO l=1,lm
632 IF (abs(pmid(i,j,l)-spval)<=small .OR. &
633 abs(uh(i,j,l)-spval)<=small .OR. &
634 abs(uh(i,j,l)-spval)<=small .OR. &
635 abs(vh(i,j,l)-spval)<=small .OR. &
636 abs(t(i,j,l)-spval)<=small .OR. &
637 abs(zmid(i,j,l)-spval)<=small) cycle loopi
638 ENDDO
639! INPUT
640 CALL mxwind(lm,pmid(i,j,1:lm),uh(i,j,1:lm) &
641! INPUT
642 ,vh(i,j,1:lm),t(i,j,1:lm),zmid(i,j,1:lm) &
643! OUTPUT
644 ,maxwp(i,j),maxwu(i,j),maxwv(i,j) &
645! OUTPUT
646 ,maxwt(i,j),maxwz(i,j))
647 ENDDO loopi
648 END DO
649! PRESSURE OF MAX WIND LEVEL
650 IF (iget(173) > 0) THEN
651!$omp parallel do private(i,j)
652 DO j=jsta,jend
653 DO i=ista,iend
654 grid1(i,j) = maxwp(i,j)
655 ENDDO
656 ENDDO
657 if(grib=='grib2') then
658 cfld=cfld+1
659 fld_info(cfld)%ifld=iavblfld(iget(173))
660!$omp parallel do private(i,j,ii,jj)
661 do j=1,jend-jsta+1
662 jj = jsta+j-1
663 do i=1,iend-ista+1
664 ii = ista+i-1
665 datapd(i,j,cfld) = grid1(ii,jj)
666 enddo
667 enddo
668 endif
669 ENDIF
670! ICAO HEIGHT OF MAX WIND LEVEL
671 IF (iget(398)>0) THEN
672 CALL icaoheight(maxwp, grid1(ista:iend,jsta:jend))
673! print*,'sample MAX WIND ICAO HEIGHTS',GRID1(im/2,(jsta+jend)/2)
674 if(grib=='grib2') then
675 cfld=cfld+1
676 fld_info(cfld)%ifld=iavblfld(iget(398))
677!$omp parallel do private(i,j,ii,jj)
678 do j=1,jend-jsta+1
679 jj = jsta+j-1
680 do i=1,iend-ista+1
681 ii = ista+i-1
682 datapd(i,j,cfld) = grid1(ii,jj)
683 enddo
684 enddo
685 endif
686 ENDIF
687! HEIGHT OF MAX WIND LEVEL
688 IF (iget(174) > 0) THEN
689!$omp parallel do private(i,j)
690 DO j=jsta,jend
691 DO i=ista,iend
692 grid1(i,j) = maxwz(i,j)
693 ENDDO
694 ENDDO
695 if(grib=='grib2') then
696 cfld=cfld+1
697 fld_info(cfld)%ifld=iavblfld(iget(174))
698!$omp parallel do private(i,j,ii,jj)
699 do j=1,jend-jsta+1
700 jj = jsta+j-1
701 do i=1,iend-ista+1
702 ii = ista+i-1
703 datapd(i,j,cfld) = grid1(ii,jj)
704 enddo
705 enddo
706 endif
707 ENDIF
708
709! MAX WIND LEVEL U WIND AND/OR V WIND.
710 IF ((iget(175) > 0).OR.(iget(176) > 0)) THEN
711!$omp parallel do private(i,j)
712 DO j=jsta,jend
713 DO i=ista,iend
714 grid1(i,j) = maxwu(i,j)
715 grid2(i,j) = maxwv(i,j)
716 ENDDO
717 ENDDO
718 if(grib=='grib2') then
719 cfld=cfld+1
720 fld_info(cfld)%ifld=iavblfld(iget(175))
721!$omp parallel do private(i,j,ii,jj)
722 do j=1,jend-jsta+1
723 jj = jsta+j-1
724 do i=1,iend-ista+1
725 ii = ista+i-1
726 datapd(i,j,cfld) = grid1(ii,jj)
727 enddo
728 enddo
729 cfld=cfld+1
730 fld_info(cfld)%ifld=iavblfld(iget(176))
731!$omp parallel do private(i,j,ii,jj)
732 do j=1,jend-jsta+1
733 jj = jsta+j-1
734 do i=1,iend-ista+1
735 ii = ista+i-1
736 datapd(i,j,cfld) = grid2(ii,jj)
737 enddo
738 enddo
739 endif
740 ENDIF
741! TEMPERATURE OF MAX WIND LEVEL
742 IF (iget(314) > 0) THEN
743!$omp parallel do private(i,j)
744 DO j=jsta,jend
745 DO i=ista,iend
746 grid1(i,j)=maxwt(i,j)
747 ENDDO
748 ENDDO
749 if(grib=='grib2') then
750 cfld=cfld+1
751 fld_info(cfld)%ifld=iavblfld(iget(314))
752!$omp parallel do private(i,j,ii,jj)
753 do j=1,jend-jsta+1
754 jj = jsta+j-1
755 do i=1,iend-ista+1
756 ii = ista+i-1
757 datapd(i,j,cfld) = grid1(ii,jj)
758 enddo
759 enddo
760 endif
761 ENDIF
762 deallocate(maxwp,maxwz,maxwu,maxwv,maxwt)
763 ENDIF
764
765!
766! ***BLOCK 3-1: FD LEVEL (selected) T, Q, U, AND V.
767!
768 IF ( (iget(059)>0.or.iget(586)>0).OR.iget(911)>0.OR. &
769 (iget(060)>0.or.iget(576)>0).OR. &
770 (iget(061)>0.or.iget(577)>0).OR. &
771 (iget(451)>0.or.iget(578)>0).OR.iget(580)>0 ) THEN
772
773 ALLOCATE(t7d(ista:iend,jsta:jend,nfd), q7d(ista:iend,jsta:jend,nfd), &
774 u7d(ista:iend,jsta:jend,nfd), v6d(ista:iend,jsta:jend,nfd), &
775 p7d(ista:iend,jsta:jend,nfd), icingfd(ista:iend,jsta:jend,nfd))
776
777!
778! DETERMINE WHETHER TO DO MSL OR AGL FD LEVELS
779!
780 itypefdlvl=1
781 DO ifd = 1,nfd
782 IF (iget(059)>0) THEN
783 IF (lvls(ifd,iget(059))>1) itypefdlvl(ifd)=2
784 ENDIF
785 IF (iget(911)>0) THEN
786 IF (lvls(ifd,iget(911))>1) itypefdlvl(ifd)=2
787 ENDIF
788!for grib2, spec hgt only
789 IF (iget(586)>0) THEN
790 IF(lvls(ifd,iget(586))>0) itypefdlvl(ifd)=2
791 ENDIF
792 IF (iget(060)>0) THEN
793 IF (lvls(ifd,iget(060))>1) itypefdlvl(ifd)=2
794 ENDIF
795 IF (iget(576)>0) THEN
796 IF(lvls(ifd,iget(576))>0) itypefdlvl(ifd)=2
797 ENDIF
798 IF (iget(061)>0) THEN
799 IF (lvls(ifd,iget(061))>1) itypefdlvl(ifd)=2
800 ENDIF
801 IF (iget(577)>0) then
802 if(lvls(ifd,iget(577))>0) itypefdlvl(ifd)=2
803 ENDIF
804 IF (iget(451)>0) THEN
805 IF (lvls(ifd,iget(451))>1) itypefdlvl(ifd)=2
806 ENDIF
807 IF (iget(578)>0) then
808 if(lvls(ifd,iget(578))>0) itypefdlvl(ifd)=2
809 ENDIF
810
811 IF (iget(580)>0) then
812 if(lvls(ifd,iget(580))>1) itypefdlvl(ifd)=2
813 ENDIF
814 IF (iget(587)>0) then
815 if(lvls(ifd,iget(587))>0) itypefdlvl(ifd)=2
816 ENDIF
817
818 ENDDO
819! print *,'call FDLVL with ITYPEFDLVL: ', ITYPEFDLVL,'for tmp,lvls=',LVLS(1:15,iget(59)), &
820! 'grib2tmp lvs=',LVLS(1:15,iget(586))
821
822 CALL fdlvl(itypefdlvl,t7d,q7d,u7d,v6d,p7d,icingfd)
823!
824 loop_10: DO ifd = 1,nfd
825!
826! FD LEVEL TEMPERATURE.
827 iget1 = iget(059)
828 iget2 = iget(586)
829 if (iget1 > 0) then
830 work1 = lvls(ifd,iget1)
831 else
832 work1 = 0.0
833 endif
834 if (iget2 > 0) then
835 work2 = lvls(ifd,iget2)
836 else
837 work2 = 0.0
838 endif
839 IF (iget1 > 0 .or. iget2 > 0) THEN
840 IF (work1 > 0 .or. work2 > 0) THEN
841
842!$omp parallel do private(i,j)
843 DO j=jsta,jend
844 DO i=ista,iend
845 grid1(i,j) = t7d(i,j,ifd)
846 ENDDO
847 ENDDO
848 IF(work1 > 0) then
849 if(grib == 'grib2') then
850 cfld = cfld + 1
851 fld_info(cfld)%ifld = iavblfld(iget1)
852 fld_info(cfld)%lvl = lvlsxml(ifd,iget1)
853!$omp parallel do private(i,j,ii,jj)
854 do j=1,jend-jsta+1
855 jj = jsta+j-1
856 do i=1,iend-ista+1
857 ii = ista+i-1
858 datapd(i,j,cfld) = grid1(ii,jj)
859 enddo
860 enddo
861 endif
862 ENDIF
863 IF (work2 > 0) THEN
864 if(grib == 'grib2') then
865 cfld = cfld + 1
866 fld_info(cfld)%ifld = iavblfld(iget2)
867 fld_info(cfld)%lvl = lvlsxml(ifd,iget2)
868!$omp parallel do private(i,j,ii,jj)
869 do j=1,jend-jsta+1
870 jj = jsta+j-1
871 do i=1,iend-ista+1
872 ii = ista+i-1
873 datapd(i,j,cfld) = grid1(ii,jj)
874 enddo
875 enddo
876 endif
877 ENDIF
878 ENDIF
879 ENDIF
880
881! FD LEVEL VIRTUAL TEMPERATURE.
882 IF (iget(911)>0) THEN
883 IF (lvls(ifd,iget(911))>0) THEN
884 DO j=jsta,jend
885 DO i=ista,iend
886 if ( t7d(i,j,ifd) > 600 ) then
887 grid1(i,j)=spval
888 else
889 grid1(i,j)=t7d(i,j,ifd)*(1.+0.608*q7d(i,j,ifd))
890 endif
891 !print *, "grid value ",T7D(I,J,IFD),Q7D(I,J,IFD),T7D(I,J,IFD)*(1.+0.608*Q7D(I,J,IFD)),GRID1(I,J)
892 ENDDO
893 ENDDO
894 IF(lvls(ifd,iget(911))>0) then
895 if(grib=='grib2') then
896 cfld=cfld+1
897 fld_info(cfld)%ifld=iavblfld(iget(911))
898 fld_info(cfld)%lvl=lvlsxml(ifd,iget(911))
899 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
900 endif
901 ENDIF
902 ENDIF
903 ENDIF
904
905!
906! FD LEVEL SPEC HUMIDITY.
907 iget1 = iget(451)
908 iget2 = iget(578)
909 if (iget1 > 0) then
910 work1 = lvls(ifd,iget1)
911 else
912 work1 = 0.0
913 endif
914 if (iget2 > 0) then
915 work2 = lvls(ifd,iget2)
916 else
917 work2 = 0.0
918 endif
919 IF (iget1 > 0 .or. iget2 > 0) THEN
920 IF (work1 > 0 .or. work2 > 0)THEN
921!$omp parallel do private(i,j)
922 DO j=jsta,jend
923 DO i=ista,iend
924 grid1(i,j) = q7d(i,j,ifd)
925 ENDDO
926 ENDDO
927 if(work1 > 0) then
928 if(grib == 'grib2') then
929 cfld = cfld + 1
930 fld_info(cfld)%ifld = iavblfld(iget1)
931 fld_info(cfld)%lvl = lvlsxml(ifd,iget1)
932!$omp parallel do private(i,j,ii,jj)
933 do j=1,jend-jsta+1
934 jj = jsta+j-1
935 do i=1,iend-ista+1
936 ii = ista+i-1
937 datapd(i,j,cfld) = grid1(ii,jj)
938 enddo
939 enddo
940 endif
941 endif
942 if(work2 > 0) then
943 if(grib == 'grib2') then
944 cfld = cfld + 1
945 fld_info(cfld)%ifld = iavblfld(iget2)
946 fld_info(cfld)%lvl = lvlsxml(ifd,iget2)
947!$omp parallel do private(i,j,ii,jj)
948 do j=1,jend-jsta+1
949 jj = jsta+j-1
950 do i=1,iend-ista+1
951 ii = ista+i-1
952 datapd(i,j,cfld) = grid1(ii,jj)
953 enddo
954 enddo
955 endif
956 endif
957 ENDIF
958 ENDIF
959!
960! FD LEVEL PRESSURE
961 iget1 = iget(482)
962 iget2 = iget(579)
963 if (iget1 > 0) then
964 work1 = lvls(ifd,iget1)
965 else
966 work1 = 0.0
967 endif
968 if (iget2 > 0) then
969 work2 = lvls(ifd,iget2)
970 else
971 work2 = 0.0
972 endif
973 IF (iget1 > 0 .or. iget2 > 0) THEN
974 IF (work1 > 0 .or. work2 > 0) THEN
975!$omp parallel do private(i,j)
976 DO j=jsta,jend
977 DO i=ista,iend
978 grid1(i,j) = p7d(i,j,ifd)
979 ENDDO
980 ENDDO
981 if(work1 > 0) then
982 if(grib == 'grib2') then
983 cfld = cfld + 1
984 fld_info(cfld)%ifld = iavblfld(iget1)
985 fld_info(cfld)%lvl = lvlsxml(ifd,iget1)
986!$omp parallel do private(i,j,ii,jj)
987 do j=1,jend-jsta+1
988 jj = jsta+j-1
989 do i=1,iend-ista+1
990 ii = ista+i-1
991 datapd(i,j,cfld) = grid1(ii,jj)
992 enddo
993 enddo
994 endif
995 endif
996 if(work2 > 0) then
997 if(grib == 'grib2') then
998 cfld = cfld + 1
999 fld_info(cfld)%ifld = iavblfld(iget2)
1000 fld_info(cfld)%lvl = lvlsxml(ifd,iget2)
1001!$omp parallel do private(i,j,ii,jj)
1002 do j=1,jend-jsta+1
1003 jj = jsta+j-1
1004 do i=1,iend-ista+1
1005 ii = ista+i-1
1006 datapd(i,j,cfld) = grid1(ii,jj)
1007 enddo
1008 enddo
1009 endif
1010 endif
1011 ENDIF
1012 ENDIF
1013!
1014! FD LEVEL ICING
1015 iget1 = iget(580)
1016 iget2 = iget(587)
1017 if (iget1 > 0) then
1018 work1 = lvls(ifd,iget1)
1019 else
1020 work1 = 0.0
1021 endif
1022 if (iget2 > 0) then
1023 work2 = lvls(ifd,iget2)
1024 else
1025 work2 = 0.0
1026 endif
1027 IF (iget1 > 0 .or. iget2 > 0) THEN
1028 IF (work1 > 0 .or. work2 > 0) THEN
1029!$omp parallel do private(i,j)
1030 DO j=jsta,jend
1031 DO i=ista,iend
1032 grid1(i,j) = icingfd(i,j,ifd)
1033 ENDDO
1034 ENDDO
1035 if(iget1 > 0) then
1036 if(grib == 'grib2') then
1037 cfld = cfld + 1
1038 fld_info(cfld)%ifld = iavblfld(iget1)
1039 fld_info(cfld)%lvl = lvlsxml(ifd,iget1)
1040!$omp parallel do private(i,j,ii,jj)
1041 do j=1,jend-jsta+1
1042 jj = jsta+j-1
1043 do i=1,iend-ista+1
1044 ii = ista+i-1
1045 datapd(i,j,cfld) = grid1(ii,jj)
1046 enddo
1047 enddo
1048 endif
1049 endif
1050 if(work2 > 0) then
1051 if(grib == 'grib2') then
1052 cfld = cfld + 1
1053 fld_info(cfld)%ifld = iavblfld(iget2)
1054 fld_info(cfld)%lvl = lvlsxml(ifd,iget2)
1055!$omp parallel do private(i,j,ii,jj)
1056 do j=1,jend-jsta+1
1057 jj = jsta+j-1
1058 do i=1,iend-ista+1
1059 ii = ista+i-1
1060 datapd(i,j,cfld) = grid1(ii,jj)
1061 enddo
1062 enddo
1063 endif
1064 endif
1065
1066 ENDIF
1067 ENDIF
1068!
1069
1070!
1071!
1072! FD LEVEL U WIND AND/OR V WIND.
1073 IF ((iget(060)>0).OR.(iget(061)>0)) THEN
1074!$omp parallel do private(i,j)
1075 DO j=jsta,jend
1076 DO i=ista,iend
1077 grid1(i,j)=u7d(i,j,ifd)
1078 grid2(i,j)=v6d(i,j,ifd)
1079 ENDDO
1080 ENDDO
1081 IF (iget(060)>0) THEN
1082 IF (lvls(ifd,iget(060))>0) then
1083 if(grib=='grib2') then
1084 cfld=cfld+1
1085 fld_info(cfld)%ifld=iavblfld(iget(060))
1086 fld_info(cfld)%lvl=lvlsxml(ifd,iget(060))
1087!$omp parallel do private(i,j,ii,jj)
1088 do j=1,jend-jsta+1
1089 jj = jsta+j-1
1090 do i=1,iend-ista+1
1091 ii = ista+i-1
1092 datapd(i,j,cfld) = grid1(ii,jj)
1093 enddo
1094 enddo
1095 endif
1096 ENDIF
1097 ENDIF
1098 IF (iget(061)>0) THEN
1099 IF (lvls(ifd,iget(061))>0) THEN
1100 if(grib=='grib2') then
1101 cfld=cfld+1
1102 fld_info(cfld)%ifld=iavblfld(iget(061))
1103 fld_info(cfld)%lvl=lvlsxml(ifd,iget(061))
1104!$omp parallel do private(i,j,ii,jj)
1105 do j=1,jend-jsta+1
1106 jj = jsta+j-1
1107 do i=1,iend-ista+1
1108 ii = ista+i-1
1109 datapd(i,j,cfld) = grid2(ii,jj)
1110 enddo
1111 enddo
1112 endif
1113 ENDIF
1114 ENDIF
1115 ENDIF
1116!
1117! FD LEVEL U WIND AND/OR V WIND.
1118 IF ((iget(576)>0).OR.(iget(577)>0)) THEN
1119!$omp parallel do private(i,j)
1120 DO j=jsta,jend
1121 DO i=ista,iend
1122 grid1(i,j) = u7d(i,j,ifd)
1123 grid2(i,j) = v6d(i,j,ifd)
1124 ENDDO
1125 ENDDO
1126 IF (iget(576)>0) THEN
1127 IF (lvls(ifd,iget(576))>0) then
1128 if(grib=='grib2') then
1129 cfld=cfld+1
1130 fld_info(cfld)%ifld=iavblfld(iget(576))
1131 fld_info(cfld)%lvl=lvlsxml(ifd,iget(576))
1132!$omp parallel do private(i,j,ii,jj)
1133 do j=1,jend-jsta+1
1134 jj = jsta+j-1
1135 do i=1,iend-ista+1
1136 ii = ista+i-1
1137 datapd(i,j,cfld) = grid1(ii,jj)
1138 enddo
1139 enddo
1140 endif
1141 ENDIF
1142 ENDIF
1143 IF (iget(577)>0) THEN
1144 IF (lvls(ifd,iget(577))>0) THEN
1145 if(grib=='grib2') then
1146 cfld=cfld+1
1147 fld_info(cfld)%ifld=iavblfld(iget(577))
1148 fld_info(cfld)%lvl=lvlsxml(ifd,iget(577))
1149!$omp parallel do private(i,j,ii,jj)
1150 do j=1,jend-jsta+1
1151 jj = jsta+j-1
1152 do i=1,iend-ista+1
1153 ii = ista+i-1
1154 datapd(i,j,cfld) = grid2(ii,jj)
1155 enddo
1156 enddo
1157 endif
1158 ENDIF
1159 ENDIF
1160 ENDIF
1161
1162 END DO loop_10
1163 DEALLOCATE(t7d,q7d,u7d,v6d,p7d,icingfd)
1164 ENDIF
1165
1166!
1167! ***BLOCK 3-2: FD LEVEL (from control file) GTG
1168!
1169 IF(gtg_on .and. (iget(467)>0.or.iget(468)>0.or.iget(469)>0.or.iget(477)>0)) THEN
1170 ! MASS FIELDS INTERPOLATION
1171 if(allocated(qin)) deallocate(qin)
1172 if(allocated(qtype)) deallocate(qtype)
1173 ALLOCATE(qin(ista:iend,jsta:jend,lm,nfdmax))
1174 ALLOCATE(qtype(nfdmax))
1175
1176! INITIALIZE INPUTS
1177 nfds = 0
1178 IF(iget(467) > 0) THEN
1179 nfds = nfds + 1
1180 ids(nfds) = 467
1181 qin(ista:iend,jsta:jend,1:lm,nfds)=gtg(ista:iend,jsta:jend,1:lm)
1182 qtype(nfds)="O"
1183 end if
1184 IF(iget(468) > 0) THEN
1185 nfds = nfds + 1
1186 ids(nfds) = 468
1187 qin(ista:iend,jsta:jend,1:lm,nfds)=catedr(ista:iend,jsta:jend,1:lm)
1188 qtype(nfds)="O"
1189 end if
1190 IF(iget(469) > 0) THEN
1191 nfds = nfds + 1
1192 ids(nfds) = 469
1193 qin(ista:iend,jsta:jend,1:lm,nfds)=mwt(ista:iend,jsta:jend,1:lm)
1194 qtype(nfds)="O"
1195 end if
1196
1197 IF(iget(477) > 0) THEN
1198 nfds = nfds + 1
1199 ids(nfds) = 477
1200 qin(ista:iend,jsta:jend,1:lm,nfds)=cit(ista:iend,jsta:jend,1:lm)
1201 qtype(nfds)="O"
1202 end if
1203
1204
1205! FOR Regional GTG, ALL LEVLES OF DIFFERENT VARIABLES ARE THE SAME
1206 iid=467
1207 n = iavblfld(iget(iid))
1208 nfdctl=size(pset%param(n)%level)
1209 if(allocated(itypefdlvlctl)) deallocate(itypefdlvlctl)
1210 allocate(itypefdlvlctl(nfdctl))
1211 DO ifd = 1,nfdctl
1212 itypefdlvlctl(ifd)=lvls(ifd,iget(iid))
1213 ENDDO
1214 if(allocated(htfdctl)) deallocate(htfdctl)
1215 allocate(htfdctl(nfdctl))
1216 htfdctl=pset%param(n)%level
1217
1218 if(allocated(qfd)) deallocate(qfd)
1219 ALLOCATE(qfd(ista:iend,jsta:jend,nfdctl,nfds))
1220 qfd=spval
1221
1222 ! pset%param(N)%level will not be used in FDLVL_MASS() because QTYPE='O'
1223 call fdlvl_mass(itypefdlvlctl,nfdctl,pset%param(n)%level,htfdctl,nfds,qin,qtype,qfd)
1224
1225! Adjust values before output
1226 DO n=1,nfds
1227 iid=ids(n)
1228 if(iid==467 .or. iid==468 .or. iid==469 .or. iid==477) then
1229 DO ifd = 1,nfdctl
1230 DO j=jsta,jend
1231 DO i=ista,iend
1232 if(qfd(i,j,ifd,n) < spval) then
1233 qfd(i,j,ifd,n)=max(0.0,qfd(i,j,ifd,n))
1234 qfd(i,j,ifd,n)=min(1.0,qfd(i,j,ifd,n))
1235 endif
1236 ENDDO
1237 ENDDO
1238 ENDDO
1239 endif
1240 ENDDO
1241
1242! Output
1243 DO n=1,nfds
1244 iid=ids(n)
1245
1246! For regional GTG, output the max value of EDPARM(ID=467) in the whole vertical column
1247! to MXEDPRM(ID=476)
1248 if (iid == 467 .and. iget(476) > 0) then
1249 egrid1 = spval
1250 DO ifd = 1,nfdctl
1251 DO j=jsta,jend
1252 DO i=ista,iend
1253 work1=qfd(i,j,ifd,n)
1254 if(egrid1(i,j)>=spval) then
1255 egrid1(i,j)=work1
1256 elseif(work1<spval) then
1257 if(egrid1(i,j)<work1) egrid1(i,j)=work1
1258 endif
1259 ENDDO
1260 ENDDO
1261 ENDDO
1262 DO j=jsta,jend
1263 DO i=ista,iend
1264 grid1(i,j)=egrid1(i,j)
1265 ENDDO
1266 ENDDO
1267 if(grib=='grib2') then
1268 cfld=cfld+1
1269 fld_info(cfld)%ifld=iavblfld(iget(476)) ! MXEDPRM ID
1270!$omp parallel do private(i,j,ii,jj)
1271 do j=1,jend-jsta+1
1272 jj = jsta+j-1
1273 do i=1,iend-ista+1
1274 ii = ista+i-1
1275 datapd(i,j,cfld) = grid1(ii,jj)
1276 enddo
1277 enddo
1278 endif
1279 end if
1280
1281 DO ifd = 1,nfdctl
1282 IF (lvls(ifd,iget(iid)) > 0) THEN
1283!$omp parallel do private(i,j)
1284 DO j=jsta,jend
1285 DO i=ista,iend
1286 grid1(i,j)=qfd(i,j,ifd,n)
1287 ENDDO
1288 ENDDO
1289 if(grib=='grib2') then
1290 cfld=cfld+1
1291 fld_info(cfld)%ifld=iavblfld(iget(iid))
1292 fld_info(cfld)%lvl=lvlsxml(ifd,iget(iid))
1293!$omp parallel do private(i,j,ii,jj)
1294 do j=1,jend-jsta+1
1295 jj = jsta+j-1
1296 do i=1,iend-ista+1
1297 ii = ista+i-1
1298 datapd(i,j,cfld) = grid1(ii,jj)
1299 enddo
1300 enddo
1301 endif
1302 ENDIF
1303 ENDDO
1304 ENDDO
1305
1306 DEALLOCATE(qin,qfd)
1307 DEALLOCATE(qtype)
1308 if(allocated(itypefdlvlctl)) deallocate(itypefdlvlctl)
1309 if(allocated(htfdctl)) deallocate(htfdctl)
1310
1311 ENDIF
1312!
1313!
1314! ***BLOCK 4: FREEZING LEVEL Z, RH AND P.
1315!
1316 IF ( (iget(062)>0).OR.(iget(063)>0) ) THEN
1317 CALL frzlvl(z1d,rh1d,p1d)
1318!
1319! FREEZING LEVEL HEIGHT.
1320 IF (iget(062)>0) THEN
1321!$omp parallel do private(i,j)
1322 DO j=jsta,jend
1323 DO i=ista,iend
1324 grid1(i,j)=z1d(i,j)
1325 IF (submodelname == 'RTMA') THEN
1326 freezelvl(i,j)=grid1(i,j)
1327 ENDIF
1328 ENDDO
1329 ENDDO
1330 CALL bound (grid1,d00,h99999)
1331 if(grib=='grib2') then
1332 cfld=cfld+1
1333 fld_info(cfld)%ifld=iavblfld(iget(062))
1334!$omp parallel do private(i,j,ii,jj)
1335 do j=1,jend-jsta+1
1336 jj = jsta+j-1
1337 do i=1,iend-ista+1
1338 ii = ista+i-1
1339 datapd(i,j,cfld) = grid1(ii,jj)
1340 enddo
1341 enddo
1342 endif
1343 ENDIF
1344!
1345! FREEZING LEVEL RELATIVE HUMIDITY.
1346 IF (iget(063)>0) THEN
1347!$omp parallel do private(i,j)
1348 DO j=jsta,jend
1349 DO i=ista,iend
1350 grid1(i,j) = rh1d(i,j)
1351 ENDDO
1352 ENDDO
1353 CALL sclfld(grid1(ista:iend,jsta:jend),h100,im,jm)
1354 CALL bound(grid1,h1,h100)
1355 if(grib=='grib2') then
1356 cfld=cfld+1
1357 fld_info(cfld)%ifld=iavblfld(iget(063))
1358!$omp parallel do private(i,j,ii,jj)
1359 do j=1,jend-jsta+1
1360 jj = jsta+j-1
1361 do i=1,iend-ista+1
1362 ii = ista+i-1
1363 datapd(i,j,cfld) = grid1(ii,jj)
1364 enddo
1365 enddo
1366 endif
1367 ENDIF
1368!
1369! FREEZING LEVEL PRESSURE
1370 IF (iget(753)>0) THEN
1371!$omp parallel do private(i,j)
1372 DO j=jsta,jend
1373 DO i=ista,iend
1374 grid1(i,j) = p1d(i,j)
1375 ENDDO
1376 ENDDO
1377 if(grib=='grib2') then
1378 cfld=cfld+1
1379 fld_info(cfld)%ifld=iavblfld(iget(753))
1380!$omp parallel do private(i,j,ii,jj)
1381 do j=1,jend-jsta+1
1382 jj = jsta+j-1
1383 do i=1,iend-ista+1
1384 ii = ista+i-1
1385 datapd(i,j,cfld) = grid1(ii,jj)
1386 enddo
1387 enddo
1388 endif
1389
1390 ENDIF
1391 ENDIF
1392
1393 IF (iget(165)>0 .OR. iget(350)>0.OR. iget(756)>0) THEN
1394 CALL frzlvl2(tfrz,z1d,rh1d,p1d)
1395!
1396! HIGHEST FREEZING LEVEL HEIGHT.
1397 IF (iget(165)>0)THEN
1398!$omp parallel do private(i,j)
1399 DO j=jsta,jend
1400 DO i=ista,iend
1401 grid1(i,j)=z1d(i,j)
1402 ENDDO
1403 ENDDO
1404 CALL bound (grid1,d00,h99999)
1405 if(grib=='grib2') then
1406 cfld=cfld+1
1407 fld_info(cfld)%ifld=iavblfld(iget(165))
1408!$omp parallel do private(i,j,ii,jj)
1409 do j=1,jend-jsta+1
1410 jj = jsta+j-1
1411 do i=1,iend-ista+1
1412 ii = ista+i-1
1413 datapd(i,j,cfld) = grid1(ii,jj)
1414 enddo
1415 enddo
1416 endif
1417 END IF
1418
1419! HIGHEST FREEZING LEVEL RELATIVE HUMIDITY
1420 IF (iget(350)>0)THEN
1421 grid1=spval
1422!$omp parallel do private(i,j)
1423 DO j=jsta,jend
1424 DO i=ista,iend
1425 IF(rh1d(i,j) < spval) grid1(i,j)=rh1d(i,j)*100.
1426 ENDDO
1427 ENDDO
1428 CALL bound (grid1,h1,h100)
1429 if(grib=='grib2') then
1430 cfld=cfld+1
1431 fld_info(cfld)%ifld=iavblfld(iget(350))
1432!$omp parallel do private(i,j,ii,jj)
1433 do j=1,jend-jsta+1
1434 jj = jsta+j-1
1435 do i=1,iend-ista+1
1436 ii = ista+i-1
1437 datapd(i,j,cfld) = grid1(ii,jj)
1438 enddo
1439 enddo
1440 endif
1441 END IF
1442
1443! HIGHEST FREEZING LEVEL PRESSURE
1444 IF (iget(756)>0) THEN
1445!$omp parallel do private(i,j)
1446 DO j=jsta,jend
1447 DO i=ista,iend
1448 grid1(i,j) = p1d(i,j)
1449 ENDDO
1450 ENDDO
1451 if(grib=='grib2') then
1452 cfld=cfld+1
1453 fld_info(cfld)%ifld=iavblfld(iget(756))
1454!$omp parallel do private(i,j,ii,jj)
1455 do j=1,jend-jsta+1
1456 jj = jsta+j-1
1457 do i=1,iend-ista+1
1458 ii = ista+i-1
1459 datapd(i,j,cfld) = grid1(ii,jj)
1460 enddo
1461 enddo
1462 endif
1463 ENDIF
1464
1465 ENDIF
1466!
1467! HIGHEST -10 C ISOTHERM VALUES
1468!
1469 IF (iget(776)>0 .OR. iget(777)>0.OR. iget(778)>0) THEN
1470 CALL frzlvl2(263.15,z1d,rh1d,p1d)
1471!
1472! HIGHEST -10C ISOTHERM LEVEL HEIGHT.
1473 IF (iget(776)>0)THEN
1474!$omp parallel do private(i,j)
1475 DO j=jsta,jend
1476 DO i=ista,iend
1477 grid1(i,j)=z1d(i,j)
1478 ENDDO
1479 ENDDO
1480 CALL bound (grid1,d00,h99999)
1481 if(grib=='grib2') then
1482 cfld=cfld+1
1483 fld_info(cfld)%ifld=iavblfld(iget(776))
1484!$omp parallel do private(i,j,ii,jj)
1485 do j=1,jend-jsta+1
1486 jj = jsta+j-1
1487 do i=1,iend-ista+1
1488 ii = ista+i-1
1489 datapd(i,j,cfld) = grid1(ii,jj)
1490 enddo
1491 enddo
1492 endif
1493 END IF
1494
1495! HIGHEST -10C ISOTHERM RELATIVE HUMIDITY
1496 IF (iget(777)>0)THEN
1497 grid1=spval
1498!$omp parallel do private(i,j)
1499 DO j=jsta,jend
1500 DO i=ista,iend
1501 IF(rh1d(i,j) < spval) grid1(i,j)=rh1d(i,j)*100.
1502 ENDDO
1503 ENDDO
1504 CALL bound (grid1,h1,h100)
1505 if(grib=='grib2') then
1506 cfld=cfld+1
1507 fld_info(cfld)%ifld=iavblfld(iget(777))
1508!$omp parallel do private(i,j,ii,jj)
1509 do j=1,jend-jsta+1
1510 jj = jsta+j-1
1511 do i=1,iend-ista+1
1512 ii = ista+i-1
1513 datapd(i,j,cfld) = grid1(ii,jj)
1514 enddo
1515 enddo
1516 endif
1517 END IF
1518
1519! HIGHEST -10C ISOTHERM LEVEL PRESSURE
1520 IF (iget(778)>0) THEN
1521!$omp parallel do private(i,j)
1522 DO j=jsta,jend
1523 DO i=ista,iend
1524 grid1(i,j)=p1d(i,j)
1525 ENDDO
1526 ENDDO
1527 if(grib=='grib2') then
1528 cfld=cfld+1
1529 fld_info(cfld)%ifld=iavblfld(iget(778))
1530!$omp parallel do private(i,j,ii,jj)
1531 do j=1,jend-jsta+1
1532 jj = jsta+j-1
1533 do i=1,iend-ista+1
1534 ii = ista+i-1
1535 datapd(i,j,cfld) = grid1(ii,jj)
1536 enddo
1537 enddo
1538 endif
1539 ENDIF
1540
1541 ENDIF
1542!
1543! HIGHEST -20 C ISOTHERM VALUES
1544!
1545 IF (iget(779)>0 .OR. iget(780)>0.OR. iget(781)>0) THEN
1546 CALL frzlvl2(253.15,z1d,rh1d,p1d)
1547!
1548! HIGHEST -20C ISOTHERM LEVEL HEIGHT.
1549 IF (iget(779)>0)THEN
1550!$omp parallel do private(i,j)
1551 DO j=jsta,jend
1552 DO i=ista,iend
1553 grid1(i,j)=z1d(i,j)
1554 ENDDO
1555 ENDDO
1556 CALL bound (grid1,d00,h99999)
1557 if(grib=='grib2') then
1558 cfld=cfld+1
1559 fld_info(cfld)%ifld=iavblfld(iget(779))
1560!$omp parallel do private(i,j,ii,jj)
1561 do j=1,jend-jsta+1
1562 jj = jsta+j-1
1563 do i=1,iend-ista+1
1564 ii = ista+i-1
1565 datapd(i,j,cfld) = grid1(ii,jj)
1566 enddo
1567 enddo
1568 endif
1569 END IF
1570
1571! HIGHEST -20C ISOTHERM RELATIVE HUMIDITY
1572 IF (iget(780)>0)THEN
1573 grid1=spval
1574!$omp parallel do private(i,j)
1575 DO j=jsta,jend
1576 DO i=ista,iend
1577 IF(rh1d(i,j) < spval) grid1(i,j)=rh1d(i,j)*100.
1578 ENDDO
1579 ENDDO
1580 CALL bound (grid1,h1,h100)
1581 if(grib=='grib2') then
1582 cfld=cfld+1
1583 fld_info(cfld)%ifld=iavblfld(iget(780))
1584!$omp parallel do private(i,j,ii,jj)
1585 do j=1,jend-jsta+1
1586 jj = jsta+j-1
1587 do i=1,iend-ista+1
1588 ii = ista+i-1
1589 datapd(i,j,cfld) = grid1(ii,jj)
1590 enddo
1591 enddo
1592 endif
1593 END IF
1594
1595! HIGHEST -20C ISOTHERM LEVEL PRESSURE
1596 IF (iget(781)>0) THEN
1597!$omp parallel do private(i,j)
1598 DO j=jsta,jend
1599 DO i=ista,iend
1600 grid1(i,j)=p1d(i,j)
1601 ENDDO
1602 ENDDO
1603 if(grib=='grib2') then
1604 cfld=cfld+1
1605 fld_info(cfld)%ifld=iavblfld(iget(781))
1606!$omp parallel do private(i,j,ii,jj)
1607 do j=1,jend-jsta+1
1608 jj = jsta+j-1
1609 do i=1,iend-ista+1
1610 ii = ista+i-1
1611 datapd(i,j,cfld) = grid1(ii,jj)
1612 enddo
1613 enddo
1614 endif
1615 ENDIF
1616
1617 ENDIF
1618!
1619 allocate(pbnd(ista:iend,jsta:jend,nbnd), tbnd(ista:iend,jsta:jend,nbnd), &
1620 qbnd(ista:iend,jsta:jend,nbnd), ubnd(ista:iend,jsta:jend,nbnd), &
1621 vbnd(ista:iend,jsta:jend,nbnd), rhbnd(ista:iend,jsta:jend,nbnd), &
1622 wbnd(ista:iend,jsta:jend,nbnd))
1623
1624!
1625! ***BLOCK 5: BOUNDARY LAYER FIELDS.
1626!
1627 IF ( (iget(067)>0).OR.(iget(068)>0).OR. &
1628 (iget(069)>0).OR.(iget(070)>0).OR. &
1629 (iget(071)>0).OR.(iget(072)>0).OR. &
1630 (iget(073)>0).OR.(iget(074)>0).OR. &
1631 (iget(088)>0).OR.(iget(089)>0).OR. &
1632 (iget(090)>0).OR.(iget(075)>0).OR. &
1633 (iget(109)>0).OR.(iget(110)>0).OR. &
1634 (iget(031)>0).OR.(iget(032)>0).OR. &
1635 (iget(573)>0).OR. need_ifi .OR. &
1636 (iget(107)>0).OR.(iget(091)>0).OR. &
1637 (iget(092)>0).OR.(iget(093)>0).OR. &
1638 (iget(094)>0).OR.(iget(095)>0).OR. &
1639 (iget(096)>0).OR.(iget(097)>0).OR. &
1640 (iget(098)>0).OR.(iget(221)>0) ) THEN
1641!
1642 call allocate_cape_arrays
1643! COMPUTE ETA BOUNDARY LAYER FIELDS.
1644 CALL bndlyr(pbnd,tbnd,qbnd,rhbnd,ubnd,vbnd, &
1645 wbnd,omgbnd,pwtbnd,qcnvbnd,lvlbnd)
1646
1647!$omp parallel do private(i,j)
1648 DO j=jsta,jend
1649 DO i=ista,iend
1650 egrid2(i,j) = spval
1651 ENDDO
1652 ENDDO
1653
1654!
1655! LOOP OVER NBND BOUNDARY LAYERS.
1656 boundary_layer_loop: DO lbnd = 1,nbnd
1657!
1658! BOUNDARY LAYER PRESSURE.
1659 IF (iget(067)>0) THEN
1660 IF (lvls(lbnd,iget(067))>0) THEN
1661!$omp parallel do private(i,j)
1662 DO j=jsta,jend
1663 DO i=ista,iend
1664 grid1(i,j) = pbnd(i,j,lbnd)
1665 ENDDO
1666 ENDDO
1667 if(grib=='grib2') then
1668 cfld=cfld+1
1669 fld_info(cfld)%ifld=iavblfld(iget(067))
1670 fld_info(cfld)%lvl=lvlsxml(lbnd,iget(067))
1671!$omp parallel do private(i,j,ii,jj)
1672 do j=1,jend-jsta+1
1673 jj = jsta+j-1
1674 do i=1,iend-ista+1
1675 ii = ista+i-1
1676 datapd(i,j,cfld) = grid1(ii,jj)
1677 enddo
1678 enddo
1679 endif
1680 ENDIF
1681 ENDIF
1682!
1683! BOUNDARY LAYER TEMPERATURE.
1684 IF (iget(068)>0) THEN
1685 IF (lvls(lbnd,iget(068))>0) THEN
1686!$omp parallel do private(i,j)
1687 DO j=jsta,jend
1688 DO i=ista,iend
1689 grid1(i,j)=tbnd(i,j,lbnd)
1690 ENDDO
1691 ENDDO
1692 if(grib=='grib2') then
1693 cfld=cfld+1
1694 fld_info(cfld)%ifld=iavblfld(iget(068))
1695 fld_info(cfld)%lvl=lvlsxml(lbnd,iget(068))
1696!$omp parallel do private(i,j,ii,jj)
1697 do j=1,jend-jsta+1
1698 jj = jsta+j-1
1699 do i=1,iend-ista+1
1700 ii = ista+i-1
1701 datapd(i,j,cfld) = grid1(ii,jj)
1702 enddo
1703 enddo
1704 endif
1705 ENDIF
1706 ENDIF
1707!
1708! BOUNDARY LAYER POTENTIAL TEMPERATURE.
1709 IF (iget(069)>0) THEN
1710 IF (lvls(lbnd,iget(069))>0) THEN
1711 CALL calpot(pbnd(ista,jsta,lbnd),tbnd(ista,jsta,lbnd),grid1(ista:iend,jsta:jend))
1712 if(grib=='grib2') then
1713 cfld=cfld+1
1714 fld_info(cfld)%ifld=iavblfld(iget(069))
1715 fld_info(cfld)%lvl=lvlsxml(ifd,iget(069))
1716!$omp parallel do private(i,j,ii,jj)
1717 do j=1,jend-jsta+1
1718 jj = jsta+j-1
1719 do i=1,iend-ista+1
1720 ii = ista+i-1
1721 datapd(i,j,cfld) = grid1(ii,jj)
1722 enddo
1723 enddo
1724 endif
1725 ENDIF
1726 ENDIF
1727!
1728! BOUNDARY LAYER RELATIVE HUMIDITY.
1729 IF (iget(072)>0) THEN
1730 IF (lvls(lbnd,iget(072))>0) THEN
1731!$omp parallel do private(i,j)
1732 DO j=jsta,jend
1733 DO i=ista,iend
1734 grid1(i,j)=rhbnd(i,j,lbnd)
1735 ENDDO
1736 ENDDO
1737 CALL sclfld(grid1(ista:iend,jsta:jend),h100,im,jm)
1738 CALL bound(grid1,h1,h100)
1739 if(grib=='grib2') then
1740 cfld=cfld+1
1741 fld_info(cfld)%lvl=lvlsxml(lbnd,iget(072))
1742 fld_info(cfld)%ifld=iavblfld(iget(072))
1743!$omp parallel do private(i,j,ii,jj)
1744 do j=1,jend-jsta+1
1745 jj = jsta+j-1
1746 do i=1,iend-ista+1
1747 ii = ista+i-1
1748 datapd(i,j,cfld) = grid1(ii,jj)
1749 enddo
1750 enddo
1751 endif
1752 ENDIF
1753 ENDIF
1754!
1755! BOUNDARY LAYER DEWPOINT TEMPERATURE.
1756 IF (iget(070)>0) THEN
1757 IF (lvls(lbnd,iget(070))>0) THEN
1758 CALL caldwp(pbnd(ista:iend,jsta:jend,lbnd), qbnd(ista:iend,jsta:jend,lbnd), &
1759 grid1(ista:iend,jsta:jend), tbnd(ista:iend,jsta:jend,lbnd))
1760 if(grib=='grib2') then
1761 cfld=cfld+1
1762 fld_info(cfld)%ifld=iavblfld(iget(070))
1763 fld_info(cfld)%lvl=lvlsxml(lbnd,iget(070))
1764!$omp parallel do private(i,j,ii,jj)
1765 do j=1,jend-jsta+1
1766 jj = jsta+j-1
1767 do i=1,iend-ista+1
1768 ii = ista+i-1
1769 datapd(i,j,cfld) = grid1(ii,jj)
1770 enddo
1771 enddo
1772 endif
1773 ENDIF
1774 ENDIF
1775!
1776! BOUNDARY LAYER SPECIFIC HUMIDITY.
1777 IF (iget(071)>0) THEN
1778 IF (lvls(lbnd,iget(071))>0) THEN
1779!$omp parallel do private(i,j)
1780 DO j=jsta,jend
1781 DO i=ista,iend
1782 grid1(i,j)=qbnd(i,j,lbnd)
1783 ENDDO
1784 ENDDO
1785 CALL bound(grid1,h1m12,h99999)
1786 if(grib=='grib2') then
1787 cfld=cfld+1
1788 fld_info(cfld)%ifld=iavblfld(iget(071))
1789 fld_info(cfld)%lvl=lvlsxml(lbnd,iget(071))
1790!$omp parallel do private(i,j,ii,jj)
1791 do j=1,jend-jsta+1
1792 jj = jsta+j-1
1793 do i=1,iend-ista+1
1794 ii = ista+i-1
1795 datapd(i,j,cfld) = grid1(ii,jj)
1796 enddo
1797 enddo
1798 endif
1799 ENDIF
1800 ENDIF
1801!
1802! BOUNDARY LAYER MOISTURE CONVERGENCE.
1803 IF (iget(088)>0) THEN
1804 IF (lvls(lbnd,iget(088))>0) THEN
1805!$omp parallel do private(i,j)
1806 DO j=jsta,jend
1807 DO i=ista,iend
1808 grid1(i,j) = qcnvbnd(i,j,lbnd)
1809 ENDDO
1810 ENDDO
1811 if(grib=='grib2') then
1812 cfld=cfld+1
1813 fld_info(cfld)%ifld=iavblfld(iget(088))
1814 fld_info(cfld)%lvl=lvlsxml(lbnd,iget(088))
1815!$omp parallel do private(i,j,ii,jj)
1816 do j=1,jend-jsta+1
1817 jj = jsta+j-1
1818 do i=1,iend-ista+1
1819 ii = ista+i-1
1820 datapd(i,j,cfld) = grid1(ii,jj)
1821 enddo
1822 enddo
1823 endif
1824 ENDIF
1825 ENDIF
1826!
1827! BOUNDARY LAYER U WIND AND/OR V WIND.
1828!
1829 field1=.false.
1830 field2=.false.
1831!
1832 IF(iget(073)>0)THEN
1833 IF(lvls(lbnd,iget(073))>0)field1=.true.
1834 ENDIF
1835 IF(iget(074)>0)THEN
1836 IF(lvls(lbnd,iget(074))>0)field2=.true.
1837 ENDIF
1838!
1839 IF(field1.OR.field2)THEN
1840!$omp parallel do private(i,j)
1841 DO j=jsta,jend
1842 DO i=ista,iend
1843 grid1(i,j) = ubnd(i,j,lbnd)
1844 grid2(i,j) = vbnd(i,j,lbnd)
1845 ENDDO
1846 ENDDO
1847!
1848 IF (iget(073)>0) THEN
1849 IF (lvls(lbnd,iget(073))>0) then
1850 if(grib=='grib2') then
1851 cfld=cfld+1
1852 fld_info(cfld)%ifld=iavblfld(iget(073))
1853 fld_info(cfld)%lvl=lvlsxml(lbnd,iget(073))
1854!$omp parallel do private(i,j,ii,jj)
1855 do j=1,jend-jsta+1
1856 jj = jsta+j-1
1857 do i=1,iend-ista+1
1858 ii = ista+i-1
1859 datapd(i,j,cfld) = grid1(ii,jj)
1860 enddo
1861 enddo
1862 endif
1863 ENDIF
1864 ENDIF
1865 IF (iget(074)>0) THEN
1866 IF (lvls(lbnd,iget(074))>0) THEN
1867 if(grib=='grib2') then
1868 cfld=cfld+1
1869 fld_info(cfld)%ifld=iavblfld(iget(074))
1870 fld_info(cfld)%lvl=lvlsxml(lbnd,iget(074))
1871!$omp parallel do private(i,j,ii,jj)
1872 do j=1,jend-jsta+1
1873 jj = jsta+j-1
1874 do i=1,iend-ista+1
1875 ii = ista+i-1
1876 datapd(i,j,cfld) = grid2(ii,jj)
1877 enddo
1878 enddo
1879 endif
1880 ENDIF
1881 ENDIF
1882 ENDIF
1883!
1884! BOUNDARY LAYER OMEGA.
1885 IF (iget(090)>0) THEN
1886 IF (lvls(lbnd,iget(090))>0) THEN
1887!$omp parallel do private(i,j)
1888 DO j=jsta,jend
1889 DO i=ista,iend
1890 grid1(i,j) = omgbnd(i,j,lbnd)
1891 ENDDO
1892 ENDDO
1893 if(grib=='grib2') then
1894 cfld=cfld+1
1895 fld_info(cfld)%ifld=iavblfld(iget(090))
1896 fld_info(cfld)%lvl=lvlsxml(lbnd,iget(090))
1897!$omp parallel do private(i,j,ii,jj)
1898 do j=1,jend-jsta+1
1899 jj = jsta+j-1
1900 do i=1,iend-ista+1
1901 ii = ista+i-1
1902 datapd(i,j,cfld) = grid1(ii,jj)
1903 enddo
1904 enddo
1905 endiF
1906 ENDIF
1907 ENDIF
1908!
1909! BOUNDARY LAYER PRECIPITBLE WATER.
1910 IF (iget(089)>0) THEN
1911 IF (lvls(lbnd,iget(089))>0) THEN
1912!$omp parallel do private(i,j)
1913 DO j=jsta,jend
1914 DO i=ista,iend
1915 grid1(i,j) = pwtbnd(i,j,lbnd)
1916 ENDDO
1917 ENDDO
1918 CALL bound(grid1,d00,h99999)
1919 if(grib=='grib2') then
1920 cfld=cfld+1
1921 fld_info(cfld)%ifld=iavblfld(iget(089))
1922 fld_info(cfld)%lvl=lvlsxml(lbnd,iget(089))
1923!$omp parallel do private(i,j,ii,jj)
1924 do j=1,jend-jsta+1
1925 jj = jsta+j-1
1926 do i=1,iend-ista+1
1927 ii = ista+i-1
1928 datapd(i,j,cfld) = grid1(ii,jj)
1929 enddo
1930 enddo
1931 endif
1932 ENDIF
1933 ENDIF
1934!
1935! BOUNDARY LAYER LIFTED INDEX.
1936 IF (iget(075)>0 .OR. iget(031)>0 .OR. iget(573)>0) THEN
1937 CALL otlft(pbnd(ista,jsta,lbnd),tbnd(ista,jsta,lbnd), &
1938 qbnd(ista,jsta,lbnd),grid1(ista:iend,jsta:jend))
1939 IF(iget(075)>0)THEN
1940 IF (lvls(lbnd,iget(075))>0) THEN
1941 if(grib=='grib2') then
1942 cfld=cfld+1
1943 fld_info(cfld)%ifld=iavblfld(iget(075))
1944 fld_info(cfld)%lvl=lvlsxml(lbnd,iget(075))
1945!$omp parallel do private(i,j,ii,jj)
1946 do j=1,jend-jsta+1
1947 jj = jsta+j-1
1948 do i=1,iend-ista+1
1949 ii = ista+i-1
1950 datapd(i,j,cfld) = grid1(ii,jj)
1951 enddo
1952 enddo
1953 endif
1954 ENDIF
1955 END IF
1956 IF(iget(031)>0 .or. iget(573)>0)THEN
1957!$omp parallel do private(i,j)
1958 DO j=jsta,jend
1959 DO i=ista,iend
1960 egrid2(i,j) = min(egrid2(i,j),grid1(i,j))
1961 END DO
1962 END DO
1963 END IF
1964 ENDIF
1965!
1966! END OF ETA BOUNDARY LAYER LOOP.
1967 END DO boundary_layer_loop
1968 deallocate(omgbnd,pwtbnd,qcnvbnd)
1969!
1970! BEST LIFTED INDEX FROM BOUNDARY LAYER FIELDS.
1971!
1972 IF (iget(031)>0 .OR. iget(573)>0 ) THEN
1973! DO J=JSTA,JEND
1974! DO I=ISTA,IEND
1975! EGRID1(I,J) = H99999
1976! EGRID2(I,J) = H99999
1977! ENDDO
1978! ENDDO
1979!
1980! DO 50 LBND = 1,NBND
1981! CALL OTLFT(PBND(1,1,LBND),TBND(1,1,LBND), &
1982! QBND(1,1,LBND),EGRID2)
1983! DO J=JSTA,JEND
1984! DO I=ISTA,IEND
1985! EGRID1(I,J)=AMIN1(EGRID1(I,J),EGRID2(I,J))
1986! ENDDO
1987! ENDDO
1988! 50 CONTINUE
1989!$omp parallel do private(i,j)
1990 DO j=jsta,jend
1991 DO i=ista,iend
1992 grid1(i,j)=egrid2(i,j)
1993 ENDDO
1994 ENDDO
1995! print*,'writting out best lifted index'
1996
1997 if (iget(031)>0) then
1998 if(grib=='grib2') then
1999 cfld=cfld+1
2000 fld_info(cfld)%ifld=iavblfld(iget(031))
2001 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
2002 endif
2003 endif
2004
2005 if(iget(573)> 0 ) THEN
2006 if(grib=='grib2') then
2007 cfld=cfld+1
2008 fld_info(cfld)%ifld=iavblfld(iget(573))
2009!$omp parallel do private(i,j,ii,jj)
2010 do j=1,jend-jsta+1
2011 jj = jsta+j-1
2012 do i=1,iend-ista+1
2013 ii = ista+i-1
2014 datapd(i,j,cfld) = grid1(ii,jj)
2015 enddo
2016 enddo
2017 endif
2018 endif
2019
2020 END IF
2021!
2022! BEST BOUNDARY LAYER CAPE AND CINS.
2023!
2024 field1=.false.
2025 field2=.false.
2026!
2027 IF(iget(032)>0)THEN
2028 IF(lvls(2,iget(032))>0)field1=.true.
2029 ENDIF
2030 IF(iget(107)>0)THEN
2031 IF(lvls(2,iget(107))>0)field2=.true.
2032 ENDIF
2033!
2034 IF(iget(566)>0)THEN
2035 field1=.true.
2036 ENDIF
2037 IF(iget(567)>0)THEN
2038 field2=.true.
2039 ENDIF
2040!
2041 !if(grib=="grib2") print *,'in MISCLN.f,iget(566)=', &
2042 ! iget(566), 'iget(567)=',iget(567),'LVLSXML(1,IGET(566)=', &
2043 ! LVLSXML(1,IGET(566)),'LVLSXML(1,IGET(567)=', &
2044 ! LVLSXML(1,IGET(567)),'field1=',field1,'field2=',field2
2045!
2046 IF(field1.OR.field2.OR.need_ifi)THEN
2047 itype = 2
2048 call allocate_cape_arrays
2049!
2050!$omp parallel do private(i,j)
2051 DO j=jsta,jend
2052 DO i=ista,iend
2053 egrid1(i,j) = -h99999
2054 egrid2(i,j) = -h99999
2055 ENDDO
2056 ENDDO
2057!
2058 loop_80: DO lbnd = 1,nbnd
2059 CALL calthte(pbnd(ista,jsta,lbnd),tbnd(ista,jsta,lbnd), &
2060 qbnd(ista,jsta,lbnd),egrid1)
2061!$omp parallel do private(i,j)
2062 DO j=jsta,jend
2063 DO i=ista,iend
2064 IF (egrid1(i,j) > egrid2(i,j)) THEN
2065 egrid2(i,j) = egrid1(i,j)
2066 lb2(i,j) = lvlbnd(i,j,lbnd)
2067 p1d(i,j) = pbnd(i,j,lbnd)
2068 t1d(i,j) = tbnd(i,j,lbnd)
2069 q1d(i,j) = qbnd(i,j,lbnd)
2070 ENDIF
2071 ENDDO
2072 ENDDO
2073 ENDDO loop_80
2074!
2075 dpbnd = 0.
2076 CALL calcape(itype,dpbnd,p1d,t1d,q1d,lb2,egrid1, &
2077 egrid2,egrid3,egrid4,egrid5)
2078
2079!
2080 IF(iget(566)>0 .or. need_ifi) THEN
2081 grid1=spval
2082!$omp parallel do private(i,j)
2083 DO j=jsta,jend
2084 DO i=ista,iend
2085 IF(t1d(i,j) < spval) grid1(i,j) = egrid1(i,j)
2086 ENDDO
2087 ENDDO
2088 CALL bound(grid1,d00,h99999)
2089!$omp parallel do private(i,j)
2090 DO j=jsta,jend
2091 DO i=ista,iend
2092 cape(i,j) = grid1(i,j)
2093 ENDDO
2094 ENDDO
2095 ENDIF
2096
2097 IF (iget(566)>0) THEN
2098 if(grib=='grib2') then
2099 cfld=cfld+1
2100 fld_info(cfld)%ifld=iavblfld(iget(566))
2101 fld_info(cfld)%lvl=lvlsxml(1,iget(566))
2102!$omp parallel do private(i,j,ii,jj)
2103 do j=1,jend-jsta+1
2104 jj = jsta+j-1
2105 do i=1,iend-ista+1
2106 ii = ista+i-1
2107 datapd(i,j,cfld) = cape(ii,jj)
2108 enddo
2109 enddo
2110 endif
2111 ENDIF
2112!
2113 IF (iget(567) > 0 .or. need_ifi) THEN
2114! dong add missing value for CIN
2115 grid1=spval
2116!$omp parallel do private(i,j)
2117 DO j=jsta,jend
2118 DO i=ista,iend
2119 IF(t1d(i,j) < spval) grid1(i,j) = - egrid2(i,j)
2120 ENDDO
2121 ENDDO
2122!
2123 CALL bound(grid1,d00,h99999)
2124!
2125!$omp parallel do private(i,j)
2126 DO j=jsta,jend
2127 DO i=ista,iend
2128 IF(t1d(i,j) < spval) grid1(i,j) = - grid1(i,j)
2129 cin(i,j) = grid1(i,j)
2130 ENDDO
2131 ENDDO
2132 ENDIF
2133
2134 IF(iget(567) > 0) THEN
2135 if(grib=='grib2') then
2136 cfld=cfld+1
2137 fld_info(cfld)%ifld=iavblfld(iget(567))
2138 fld_info(cfld)%lvl=lvlsxml(1,iget(567))
2139!$omp parallel do private(i,j,ii,jj)
2140 do j=1,jend-jsta+1
2141 jj = jsta+j-1
2142 do i=1,iend-ista+1
2143 ii = ista+i-1
2144 datapd(i,j,cfld) = cin(ii,jj)
2145 enddo
2146 enddo
2147 endif
2148 ENDIF
2149 ENDIF
2150!
2151
2152! PBL HEIGHT
2153 IF(iget(221) > 0) THEN
2154!$omp parallel do private(i,j)
2155 DO j=jsta,jend
2156 DO i=ista,iend
2157 grid1(i,j) = pblh(i,j)
2158 ENDDO
2159 ENDDO
2160 if(grib=='grib2') then
2161 cfld=cfld+1
2162 fld_info(cfld)%ifld=iavblfld(iget(221))
2163!$omp parallel do private(i,j,ii,jj)
2164 do j=1,jend-jsta+1
2165 jj = jsta+j-1
2166 do i=1,iend-ista+1
2167 ii = ista+i-1
2168 datapd(i,j,cfld) = grid1(ii,jj)
2169 enddo
2170 enddo
2171 endif
2172 END IF
2173! BOUNDARY LAYER LIFTING CONDENSATION PRESSURE AND HEIGHT.
2174! EGRID1 IS LCL PRESSURE. EGRID2 IS LCL HEIGHT.
2175!
2176 IF ( (iget(109)>0).OR.(iget(110)>0) ) THEN
2177 CALL callcl(pbnd(ista,jsta,1),tbnd(ista,jsta,1), &
2178 qbnd(ista,jsta,1),egrid1,egrid2)
2179 IF (iget(109)>0) THEN
2180 grid1=spval
2181!$omp parallel do private(i,j)
2182 DO j=jsta,jend
2183 DO i=ista,iend
2184 IF(tbnd(i,j,1) < spval) grid1(i,j) = egrid2(i,j)
2185 ENDDO
2186 ENDDO
2187 if(grib=='grib2') then
2188 cfld=cfld+1
2189 fld_info(cfld)%ifld=iavblfld(iget(109))
2190!$omp parallel do private(i,j,ii,jj)
2191 do j=1,jend-jsta+1
2192 jj = jsta+j-1
2193 do i=1,iend-ista+1
2194 ii = ista+i-1
2195 datapd(i,j,cfld) = grid1(ii,jj)
2196 enddo
2197 enddo
2198 endif
2199 ENDIF
2200 IF (iget(110)>0) THEN
2201 grid1=spval
2202!$omp parallel do private(i,j)
2203 DO j=jsta,jend
2204 DO i=ista,iend
2205 IF(tbnd(i,j,1) < spval) grid1(i,j) = egrid1(i,j)
2206 ENDDO
2207 ENDDO
2208 if(grib=='grib2') then
2209 cfld=cfld+1
2210 fld_info(cfld)%ifld=iavblfld(iget(110))
2211!$omp parallel do private(i,j,ii,jj)
2212 do j=1,jend-jsta+1
2213 jj = jsta+j-1
2214 do i=1,iend-ista+1
2215 ii = ista+i-1
2216 datapd(i,j,cfld) = grid1(ii,jj)
2217 enddo
2218 enddo
2219 endif
2220 ENDIF
2221 ENDIF
2222!
2223! NGM BOUNDARY LAYER FIELDS.
2224!
2225 IF ( (iget(091)>0).OR.(iget(092)>0).OR. &
2226 (iget(093)>0).OR.(iget(094)>0).OR. &
2227 (iget(095)>0).OR.(iget(095)>0).OR. &
2228 (iget(096)>0).OR.(iget(097)>0).OR. &
2229 (iget(098)>0) ) THEN
2230
2231 allocate(t78483(ista:iend,jsta:jend), t89671(ista:iend,jsta:jend), &
2232 p78483(ista:iend,jsta:jend), p89671(ista:iend,jsta:jend))
2233!
2234! COMPUTE SIGMA 0.89671 AND 0.78483 TEMPERATURES
2235! INTERPOLATE LINEAR IN LOG P
2236 IF (iget(097)>0.OR.iget(098)>0) THEN
2237!$omp parallel do private(i,j)
2238 DO j=jsta,jend
2239 DO i=ista,iend
2240 p78483(i,j) = log(pint(i,j,nint(lmh(i,j)))*0.78483)
2241 p89671(i,j) = log(pint(i,j,nint(lmh(i,j)))*0.89671)
2242 ENDDO
2243 ENDDO
2244 done =.false.
2245 done1=.false.
2246!!$omp parallel do private(fac1,fac2,pkl1,pku1,t78483,t89671)
2247 DO l=2,lm
2248 DO j=jsta,jend
2249 DO i=ista,iend
2250 pkl1=0.5*(alpint(i,j,l)+alpint(i,j,l+1))
2251 pku1=0.5*(alpint(i,j,l)+alpint(i,j,l-1))
2252! IF(I==1 .AND. J==1)PRINT*,'L,P89671,PKL1,PKU1= ', &
2253! L,P89671(I,J), PKL1, PKU1
2254 IF(p78483(i,j) < pkl1.AND.p78483(i,j) > pku1) THEN
2255 fac1 = (pkl1-p78483(i,j))/(pkl1-pku1)
2256 fac2 = (p78483(i,j)-pku1)/(pkl1-pku1)
2257 t78483(i,j) = t(i,j,l)*fac2 + t(i,j,l-1)*fac1
2258 done1(i,j)=.true.
2259 ENDIF
2260 IF(p89671(i,j) < pkl1.AND.p89671(i,j) > pku1)THEN
2261 fac1 = (pkl1-p89671(i,j))/(pkl1-pku1)
2262 fac2 = (p89671(i,j)-pku1)/(pkl1-pku1)
2263 t89671(i,j) = t(i,j,l)*fac2 + t(i,j,l-1)*fac1
2264 done(i,j) = .true.
2265! IF(I==1 .AND. J==1)PRINT*,'done(1,1)= ',done(1,1)
2266 ENDIF
2267 ENDDO
2268 ENDDO
2269 ENDDO
2270! print*,'done(1,1)= ',done(1,1)
2271!$omp parallel do private(i,j,pl,tl,ql,qsat,rhl)
2272 DO j=jsta,jend
2273 DO i=ista,iend
2274 IF(.NOT. done(i,j)) THEN
2275 pl = pint(i,j,lm-1)
2276 tl = 0.5*(t(i,j,lm-2)+t(i,j,lm-1))
2277 ql = 0.5*(q(i,j,lm-2)+q(i,j,lm-1))
2278 qsat = pq0/pl *exp(a2*(tl-a3)/(tl-a4))
2279!
2280 rhl = ql/qsat
2281!
2282 IF(rhl > 1.)THEN
2283 rhl = 1.
2284 ql = rhl*qsat
2285 ENDIF
2286!
2287 IF(rhl < rhmin)THEN
2288 rhl = rhmin
2289 ql = rhl*qsat
2290 ENDIF
2291!
2292! TVRL = TL*(1.+0.608*QL)
2293! TVRBLO = TVRL*(P89671(I,J)/PL)**RGAMOG
2294! T89671(I,J) = TVRBLO/(1.+0.608*QL)
2295
2296 t89671(i,j) = tl * (p89671(i,j)/pl)**rgamog
2297!
2298
2299! PKL1=0.5*(ALPINT(I,J,LM)+ALPINT(I,J,LM+1))
2300! PKU1=0.5*(ALPINT(I,J,LM-1)+ALPINT(I,J,LM))
2301! T89671(I,J)=T(I,J,LM)+(T(I,J,LM)-T(I,J,LM-1))*
2302! + (P89671(I,J)-PKL1)/(PKL1-PKU1)
2303
2304! print*,'Debug T89671= ',i,j
2305! + ,P89671(I,J), PKL1, PKU1
2306! + ,T89671(I,J),T(I,J,LM-1),T(I,J,LM)
2307 END IF
2308
2309 IF(.NOT. done1(i,j))THEN
2310 pl = pint(i,j,lm-1)
2311 tl = 0.5*(t(i,j,lm-2)+t(i,j,lm-1))
2312 ql = 0.5*(q(i,j,lm-2)+q(i,j,lm-1))
2313 qsat = pq0/pl *exp(a2*(tl-a3)/(tl-a4))
2314!
2315 rhl = ql/qsat
2316!
2317 IF(rhl>1.)THEN
2318 rhl = 1.
2319 ql = rhl*qsat
2320 ENDIF
2321!
2322 IF(rhl<rhmin)THEN
2323 rhl = rhmin
2324 ql = rhl*qsat
2325 ENDIF
2326!
2327! TVRL = TL*(1.+0.608*QL)
2328! TVRBLO = TVRL*(P78483(I,J)/PL)**RGAMOG
2329! T78483(I,J) =TVRBLO/(1.+0.608*QL)
2330
2331 t78483(i,j) = tl * (p78483(i,j)/pl)**rgamog
2332!
2333 END IF
2334
2335 END DO
2336 END DO
2337!
2338! SIGMA 0.89671 TEMPERATURE
2339 IF (iget(097) > 0) THEN
2340 grid1=spval
2341!$omp parallel do private(i,j)
2342 DO j=jsta,jend
2343 DO i=ista,iend
2344 IF(t(i,j,lm) < spval) grid1(i,j) = t89671(i,j)
2345! IF(T89671(I,J)>350.)PRINT*,'LARGE T89671 ', &
2346! I,J,T89671(I,J)
2347 ENDDO
2348 ENDDO
2349 if(grib=='grib2') then
2350 cfld=cfld+1
2351 fld_info(cfld)%ifld=iavblfld(iget(097))
2352 fld_info(cfld)%lvl=lvlsxml(1,iget(097))
2353!$omp parallel do private(i,j,ii,jj)
2354 do j=1,jend-jsta+1
2355 jj = jsta+j-1
2356 do i=1,iend-ista+1
2357 ii = ista+i-1
2358 datapd(i,j,cfld) = grid1(ii,jj)
2359 enddo
2360 enddo
2361 endif
2362 ENDIF
2363!
2364! SIGMA 0.78483 TEMPERATURE
2365 IF (iget(098)>0) THEN
2366 grid1=spval
2367!$omp parallel do private(i,j)
2368 DO j=jsta,jend
2369 DO i=ista,iend
2370 IF(t(i,j,lm) < spval) grid1(i,j) = t78483(i,j)
2371 ENDDO
2372 ENDDO
2373 if(grib=='grib2') then
2374 cfld=cfld+1
2375 fld_info(cfld)%ifld=iavblfld(iget(098))
2376 fld_info(cfld)%lvl=lvlsxml(1,iget(098))
2377!$omp parallel do private(i,j,ii,jj)
2378 do j=1,jend-jsta+1
2379 jj = jsta+j-1
2380 do i=1,iend-ista+1
2381 ii = ista+i-1
2382 datapd(i,j,cfld) = grid1(ii,jj)
2383 enddo
2384 enddo
2385 endif
2386 ENDIF
2387 deallocate(t78483, t89671, p78483, p89671)
2388 ENDIF
2389!
2390! NGM SIGMA LAYER 0.98230 FIELDS. THESE FIELDS ARE
2391! THE FIRST ETA LAYER BOUNDARY LAYER FIELDS.
2392!
2393!
2394 IF ( (iget(091)>0).OR.(iget(092)>0).OR. &
2395 (iget(093)>0).OR.(iget(094)>0).OR. &
2396 (iget(095)>0).OR.(iget(095)>0).OR. &
2397 (iget(096)>0) ) THEN
2398!
2399!
2400! PRESSURE.
2401 IF (iget(091)>0) THEN
2402!$omp parallel do private(i,j)
2403 DO j=jsta,jend
2404 DO i=ista,iend
2405 grid1(i,j) = pbnd(i,j,1)
2406 ENDDO
2407 ENDDO
2408 if(grib=='grib2') then
2409 cfld=cfld+1
2410 fld_info(cfld)%ifld=iavblfld(iget(091))
2411!$omp parallel do private(i,j,ii,jj)
2412 do j=1,jend-jsta+1
2413 jj = jsta+j-1
2414 do i=1,iend-ista+1
2415 ii = ista+i-1
2416 datapd(i,j,cfld) = grid1(ii,jj)
2417 enddo
2418 enddo
2419 endif
2420 ENDIF
2421!
2422! TEMPERATURE.
2423 IF (iget(092)>0) THEN
2424!$omp parallel do private(i,j)
2425 DO j=jsta,jend
2426 DO i=ista,iend
2427 grid1(i,j) = tbnd(i,j,1)
2428 ENDDO
2429 ENDDO
2430 if(grib=='grib2') then
2431 cfld=cfld+1
2432 fld_info(cfld)%ifld=iavblfld(iget(092))
2433 fld_info(cfld)%lvl=lvlsxml(1,iget(092))
2434!$omp parallel do private(i,j,ii,jj)
2435 do j=1,jend-jsta+1
2436 jj = jsta+j-1
2437 do i=1,iend-ista+1
2438 ii = ista+i-1
2439 datapd(i,j,cfld) = grid1(ii,jj)
2440 enddo
2441 enddo
2442 endif
2443 ENDIF
2444!
2445! SPECIFIC HUMIDITY.
2446 IF (iget(093)>0) THEN
2447!$omp parallel do private(i,j)
2448 DO j=jsta,jend
2449 DO i=ista,iend
2450 grid1(i,j) = qbnd(i,j,1)
2451 ENDDO
2452 ENDDO
2453 CALL bound(grid1,h1m12,h99999)
2454 if(grib=='grib2') then
2455 cfld=cfld+1
2456 fld_info(cfld)%ifld=iavblfld(iget(093))
2457 fld_info(cfld)%lvl=lvlsxml(1,iget(093))
2458!$omp parallel do private(i,j,ii,jj)
2459 do j=1,jend-jsta+1
2460 jj = jsta+j-1
2461 do i=1,iend-ista+1
2462 ii = ista+i-1
2463 datapd(i,j,cfld) = grid1(ii,jj)
2464 enddo
2465 enddo
2466 endif
2467 ENDIF
2468!
2469! RELATIVE HUMIDITY.
2470 IF (iget(094)>0) THEN
2471!$omp parallel do private(i,j)
2472 DO j=jsta,jend
2473 DO i=ista,iend
2474 grid1(i,j) = rhbnd(i,j,1)
2475 ENDDO
2476 ENDDO
2477 CALL sclfld(grid1(ista:iend,jsta:jend),h100,im,jm)
2478 CALL bound(grid1,h1,h100)
2479 if(grib=='grib2') then
2480 cfld=cfld+1
2481 fld_info(cfld)%ifld=iavblfld(iget(094))
2482 fld_info(cfld)%lvl=lvlsxml(1,iget(094))
2483!$omp parallel do private(i,j,ii,jj)
2484 do j=1,jend-jsta+1
2485 jj = jsta+j-1
2486 do i=1,iend-ista+1
2487 ii = ista+i-1
2488 datapd(i,j,cfld) = grid1(ii,jj)
2489 enddo
2490 enddo
2491 endif
2492 ENDIF
2493!
2494! U AND/OR V WIND.
2495 IF ((iget(095)>0).OR.(iget(096)>0)) THEN
2496!$omp parallel do private(i,j)
2497 DO j=jsta,jend
2498 DO i=ista,iend
2499 grid1(i,j) = ubnd(i,j,1)
2500 grid2(i,j) = vbnd(i,j,1)
2501 ENDDO
2502 ENDDO
2503 IF (iget(095)>0) then
2504 if(grib=='grib2') then
2505 cfld=cfld+1
2506 fld_info(cfld)%ifld=iavblfld(iget(095))
2507 fld_info(cfld)%lvl=lvlsxml(1,iget(095))
2508!$omp parallel do private(i,j,ii,jj)
2509 do j=1,jend-jsta+1
2510 jj = jsta+j-1
2511 do i=1,iend-ista+1
2512 ii = ista+i-1
2513 datapd(i,j,cfld) = grid1(ii,jj)
2514 enddo
2515 enddo
2516 endif
2517 ENDIF
2518 IF (iget(096)>0) then
2519 if(grib=='grib2') then
2520 cfld=cfld+1
2521 fld_info(cfld)%ifld=iavblfld(iget(096))
2522 fld_info(cfld)%lvl=lvlsxml(1,iget(096))
2523!$omp parallel do private(i,j,ii,jj)
2524 do j=1,jend-jsta+1
2525 jj = jsta+j-1
2526 do i=1,iend-ista+1
2527 ii = ista+i-1
2528 datapd(i,j,cfld) = grid2(ii,jj)
2529 enddo
2530 enddo
2531 endif
2532 ENDIF
2533 ENDIF
2534 ENDIF
2535 ENDIF
2536!
2537! ENDIF FOR BOUNDARY LAYER BLOCK.
2538!
2539 ENDIF
2540!
2541!
2542! ***BLOCK 6: MISCELLANEOUS LAYER MEAN LFM AND NGM FIELDS.
2543!
2544 IF ( (iget(066)>0).OR.(iget(081)>0).OR. &
2545 (iget(082)>0).OR.(iget(104)>0).OR. &
2546 (iget(099)>0).OR.(iget(100)>0).OR. &
2547 (iget(101)>0).OR.(iget(102)>0).OR. &
2548 (iget(103)>0) ) THEN
2549!
2550! LFM "MEAN" RELATIVE HUMIDITIES AND PRECIPITABLE WATER.
2551!
2552 IF ( (iget(066)>0).OR.(iget(081)>0).OR. &
2553 (iget(082)>0).OR.(iget(104)>0) ) THEN
2554 allocate(rh3310(ista:iend,jsta:jend),rh6610(ista:iend,jsta:jend), &
2555 rh3366(ista:iend,jsta:jend),pw3310(ista:iend,jsta:jend))
2556 CALL lfmfld(rh3310,rh6610,rh3366,pw3310)
2557!
2558! SIGMA 0.33-1.00 MEAN RELATIVE HUMIIDITY.
2559 IF (iget(066)>0) THEN
2560!$omp parallel do private(i,j)
2561 DO j=jsta,jend
2562 DO i=ista,iend
2563 grid1(i,j) = rh3310(i,j)
2564 ENDDO
2565 ENDDO
2566 CALL sclfld(grid1(ista:iend,jsta:jend),h100,im,jm)
2567 CALL bound(grid1,h1,h100)
2568 if(grib=='grib2') then
2569 cfld=cfld+1
2570 fld_info(cfld)%ifld=iavblfld(iget(066))
2571 fld_info(cfld)%lvl=lvlsxml(1,iget(066))
2572!$omp parallel do private(i,j,ii,jj)
2573 do j=1,jend-jsta+1
2574 jj = jsta+j-1
2575 do i=1,iend-ista+1
2576 ii = ista+i-1
2577 datapd(i,j,cfld) = grid1(ii,jj)
2578 enddo
2579 enddo
2580! print *,'in miscln,RH0.33-1.0,cfld=',cfld,'fld=', &
2581! IAVBLFLD(IGET(066)),'lvl=',fld_info(cfld)%lvl
2582 endif
2583 ENDIF
2584!
2585! SIGMA 0.66-1.00 MEAN RELATIVE HUMIIDITY.
2586 IF (iget(081)>0) THEN
2587!$omp parallel do private(i,j)
2588 DO j=jsta,jend
2589 DO i=ista,iend
2590 grid1(i,j) = rh6610(i,j)
2591 ENDDO
2592 ENDDO
2593 CALL sclfld(grid1(ista:iend,jsta:jend),h100,im,jm)
2594 CALL bound(grid1,h1,h100)
2595 if(grib=='grib2') then
2596 cfld=cfld+1
2597 fld_info(cfld)%ifld=iavblfld(iget(081))
2598 fld_info(cfld)%lvl=lvlsxml(1,iget(081))
2599!$omp parallel do private(i,j,ii,jj)
2600 do j=1,jend-jsta+1
2601 jj = jsta+j-1
2602 do i=1,iend-ista+1
2603 ii = ista+i-1
2604 datapd(i,j,cfld) = grid1(ii,jj)
2605 enddo
2606 enddo
2607 endif
2608 ENDIF
2609!
2610! SIGMA 0.33-0.66 MEAN RELATIVE HUMIIDITY.
2611 IF (iget(082)>0) THEN
2612!$omp parallel do private(i,j)
2613 DO j=jsta,jend
2614 DO i=ista,iend
2615 grid1(i,j) = rh3366(i,j)
2616 ENDDO
2617 ENDDO
2618 CALL sclfld(grid1(ista:iend,jsta:jend),h100,im,jm)
2619 CALL bound(grid1,h1,h100)
2620 if(grib=='grib2') then
2621 cfld=cfld+1
2622 fld_info(cfld)%ifld=iavblfld(iget(082))
2623 fld_info(cfld)%lvl=lvlsxml(1,iget(082))
2624!$omp parallel do private(i,j,ii,jj)
2625 do j=1,jend-jsta+1
2626 jj = jsta+j-1
2627 do i=1,iend-ista+1
2628 ii = ista+i-1
2629 datapd(i,j,cfld) = grid1(ii,jj)
2630 enddo
2631 enddo
2632 endif
2633 ENDIF
2634!
2635! SIGMA 0.33-1.00 PRECIPITABLE WATER.
2636 IF (iget(104)>0) THEN
2637!$omp parallel do private(i,j)
2638 DO j=jsta,jend
2639 DO i=ista,iend
2640 grid1(i,j) = pw3310(i,j)
2641 ENDDO
2642 ENDDO
2643 CALL bound(grid1,d00,h99999)
2644 if(grib=='grib2') then
2645 cfld=cfld+1
2646 fld_info(cfld)%ifld=iavblfld(iget(104))
2647 fld_info(cfld)%lvl=lvlsxml(1,iget(104))
2648!$omp parallel do private(i,j,ii,jj)
2649 do j=1,jend-jsta+1
2650 jj = jsta+j-1
2651 do i=1,iend-ista+1
2652 ii = ista+i-1
2653 datapd(i,j,cfld) = grid1(ii,jj)
2654 enddo
2655 enddo
2656 endif
2657 ENDIF
2658 deallocate(rh3310,rh6610,rh3366,pw3310)
2659 ENDIF
2660!
2661! VARIOUS LAYER MEAN NGM SIGMA FIELDS.
2662!
2663 IF ( (iget(099)>0).OR.(iget(100)>0).OR. &
2664 (iget(101)>0).OR.(iget(102)>0).OR. &
2665 (iget(103)>0) ) THEN
2666 allocate(rh4710(ista_2l:iend_2u,jsta_2l:jend_2u),rh4796(ista_2l:iend_2u,jsta_2l:jend_2u), &
2667 rh1847(ista_2l:iend_2u,jsta_2l:jend_2u))
2668 allocate(rh8498(ista_2l:iend_2u,jsta_2l:jend_2u),qm8510(ista_2l:iend_2u,jsta_2l:jend_2u))
2669
2670 CALL ngmfld(rh4710,rh4796,rh1847,rh8498,qm8510)
2671!
2672! SIGMA 0.47191-1.00000 RELATIVE HUMIDITY.
2673 IF (iget(099)>0) THEN
2674!$omp parallel do private(i,j)
2675 DO j=jsta,jend
2676 DO i=ista,iend
2677 grid1(i,j) = rh4710(i,j)
2678 ENDDO
2679 ENDDO
2680 CALL sclfld(grid1(ista:iend,jsta:jend),h100,im,jm)
2681 CALL bound(grid1,h1,h100)
2682 if(grib=='grib2') then
2683 cfld=cfld+1
2684 fld_info(cfld)%ifld=iavblfld(iget(099))
2685 fld_info(cfld)%lvl=lvlsxml(1,iget(099))
2686!$omp parallel do private(i,j,ii,jj)
2687 do j=1,jend-jsta+1
2688 jj = jsta+j-1
2689 do i=1,iend-ista+1
2690 ii = ista+i-1
2691 datapd(i,j,cfld) = grid1(ii,jj)
2692 enddo
2693 enddo
2694 endif
2695 ENDIF
2696!
2697! SIGMA 0.47191-0.96470 RELATIVE HUMIDITY.
2698 IF (iget(100)>0) THEN
2699!$omp parallel do private(i,j)
2700 DO j=jsta,jend
2701 DO i=ista,iend
2702 grid1(i,j) = rh4796(i,j)
2703 ENDDO
2704 ENDDO
2705 CALL sclfld(grid1(ista:iend,jsta:jend),h100,im,jm)
2706 CALL bound(grid1,h1,h100)
2707 if(grib=='grib2') then
2708 cfld=cfld+1
2709 fld_info(cfld)%ifld=iavblfld(iget(100))
2710 fld_info(cfld)%lvl=lvlsxml(1,iget(100))
2711!$omp parallel do private(i,j,ii,jj)
2712 do j=1,jend-jsta+1
2713 jj = jsta+j-1
2714 do i=1,iend-ista+1
2715 ii = ista+i-1
2716 datapd(i,j,cfld) = grid1(ii,jj)
2717 enddo
2718 enddo
2719 endif
2720 ENDIF
2721!
2722! SIGMA 0.18019-0.47191 RELATIVE HUMIDITY.
2723 IF (iget(101)>0) THEN
2724!$omp parallel do private(i,j)
2725 DO j=jsta,jend
2726 DO i=ista,iend
2727 grid1(i,j) = rh1847(i,j)
2728 ENDDO
2729 ENDDO
2730 CALL sclfld(grid1(ista:iend,jsta:jend),h100,im,jm)
2731 CALL bound(grid1,h1,h100)
2732 if(grib=='grib2') then
2733 cfld=cfld+1
2734 fld_info(cfld)%ifld=iavblfld(iget(101))
2735 fld_info(cfld)%lvl=lvlsxml(1,iget(101))
2736!$omp parallel do private(i,j,ii,jj)
2737 do j=1,jend-jsta+1
2738 jj = jsta+j-1
2739 do i=1,iend-ista+1
2740 ii = ista+i-1
2741 datapd(i,j,cfld) = grid1(ii,jj)
2742 enddo
2743 enddo
2744 endif
2745 ENDIF
2746!
2747! SIGMA 0.84368-0.98230 RELATIVE HUMIDITY.
2748 IF (iget(102)>0) THEN
2749!$omp parallel do private(i,j)
2750 DO j=jsta,jend
2751 DO i=ista,iend
2752 grid1(i,j) = rh8498(i,j)
2753 ENDDO
2754 ENDDO
2755 CALL sclfld(grid1(ista:iend,jsta:jend),h100,im,jm)
2756 CALL bound(grid1,h1,h100)
2757 if(grib=='grib2') then
2758 cfld=cfld+1
2759 fld_info(cfld)%ifld=iavblfld(iget(102))
2760 fld_info(cfld)%lvl=lvlsxml(1,iget(102))
2761!$omp parallel do private(i,j,ii,jj)
2762 do j=1,jend-jsta+1
2763 jj = jsta+j-1
2764 do i=1,iend-ista+1
2765 ii = ista+i-1
2766 datapd(i,j,cfld) = grid1(ii,jj)
2767 enddo
2768 enddo
2769 endif
2770 ENDIF
2771!
2772! SIGMA 0.85000-1.00000 MOISTURE CONVERGENCE.
2773 IF (iget(103)>0) THEN
2774 grid1=spval
2775! CONVERT TO DIVERGENCE FOR GRIB
2776!$omp parallel do private(i,j)
2777 DO j=jsta,jend
2778 DO i=ista,iend
2779 IF(qm8510(i,j) < spval) grid1(i,j) = -1.0*qm8510(i,j)
2780 ENDDO
2781 ENDDO
2782 if(grib=='grib2') then
2783 cfld=cfld+1
2784 fld_info(cfld)%ifld=iavblfld(iget(103))
2785 fld_info(cfld)%lvl=lvlsxml(1,iget(103))
2786!$omp parallel do private(i,j,ii,jj)
2787 do j=1,jend-jsta+1
2788 jj = jsta+j-1
2789 do i=1,iend-ista+1
2790 ii = ista+i-1
2791 datapd(i,j,cfld) = grid1(ii,jj)
2792 enddo
2793 enddo
2794 endif
2795 ENDIF
2796 deallocate(rh4710,rh4796,rh1847)
2797 deallocate(rh8498,qm8510)
2798 ENDIF
2799 ENDIF
2800
2801 IF ( (iget(318)>0).OR.(iget(319)>0).OR. &
2802 (iget(320)>0))THEN
2803 allocate(rh4410(ista:iend,jsta:jend),rh7294(ista:iend,jsta:jend), &
2804 rh4472(ista:iend,jsta:jend),rh3310(ista:iend,jsta:jend))
2805 CALL lfmfld_gfs(rh4410,rh7294,rh4472,rh3310)
2806!
2807! SIGMA 0.44-1.00 MEAN RELATIVE HUMIIDITY.
2808 IF (iget(318)>0) THEN
2809 grid1=spval
2810!$omp parallel do private(i,j)
2811 DO j=jsta,jend
2812 DO i=ista,iend
2813 IF(rh4410(i,j) < spval) grid1(i,j) = rh4410(i,j)*100.
2814 ENDDO
2815 ENDDO
2816 CALL bound(grid1,d00,h100)
2817 if(grib=='grib2') then
2818 cfld=cfld+1
2819 fld_info(cfld)%ifld=iavblfld(iget(318))
2820 fld_info(cfld)%lvl=lvlsxml(1,iget(318))
2821!$omp parallel do private(i,j,ii,jj)
2822 do j=1,jend-jsta+1
2823 jj = jsta+j-1
2824 do i=1,iend-ista+1
2825 ii = ista+i-1
2826 datapd(i,j,cfld) = grid1(ii,jj)
2827 enddo
2828 enddo
2829 endif
2830 ENDIF
2831!
2832! SIGMA 0.72-0.94 MEAN RELATIVE HUMIIDITY.
2833 IF (iget(319)>0) THEN
2834 grid1=spval
2835!$omp parallel do private(i,j)
2836 DO j=jsta,jend
2837 DO i=ista,iend
2838 IF(rh7294(i,j) < spval) grid1(i,j) = rh7294(i,j)*100.
2839 ENDDO
2840 ENDDO
2841 CALL bound(grid1,d00,h100)
2842 if(grib=='grib2') then
2843 cfld=cfld+1
2844 fld_info(cfld)%ifld=iavblfld(iget(319))
2845 fld_info(cfld)%lvl=lvlsxml(1,iget(319))
2846!$omp parallel do private(i,j,ii,jj)
2847 do j=1,jend-jsta+1
2848 jj = jsta+j-1
2849 do i=1,iend-ista+1
2850 ii = ista+i-1
2851 datapd(i,j,cfld) = grid1(ii,jj)
2852 enddo
2853 enddo
2854 endif
2855 ENDIF
2856!
2857! SIGMA 0.44-0.72 MEAN RELATIVE HUMIIDITY.
2858 IF (iget(320)>0) THEN
2859 grid1=spval
2860!$omp parallel do private(i,j)
2861 DO j=jsta,jend
2862 DO i=ista,iend
2863 IF(rh4472(i,j) < spval) grid1(i,j)=rh4472(i,j)*100.
2864 ENDDO
2865 ENDDO
2866 CALL bound(grid1,d00,h100)
2867 if(grib=='grib2') then
2868 cfld=cfld+1
2869 fld_info(cfld)%ifld=iavblfld(iget(320))
2870 fld_info(cfld)%lvl=lvlsxml(1,iget(320))
2871!$omp parallel do private(i,j,ii,jj)
2872 do j=1,jend-jsta+1
2873 jj = jsta+j-1
2874 do i=1,iend-ista+1
2875 ii = ista+i-1
2876 datapd(i,j,cfld) = grid1(ii,jj)
2877 enddo
2878 enddo
2879 endif
2880 ENDIF
2881 deallocate(rh4410,rh7294,rh4472,rh3310)
2882 END IF
2883
2884! GFS computes sigma=0.9950 T, THETA, U, V from lowest two model level fields
2885 IF ( (iget(321)>0).OR.(iget(322)>0).OR. &
2886 (iget(323)>0).OR.(iget(324)>0).OR. &
2887 (iget(325)>0).OR.(iget(326)>0)) THEN
2888!$omp parallel do private(i,j)
2889 DO j=jsta,jend
2890 DO i=ista,iend
2891 egrid2(i,j) = 0.995*pint(i,j,lm+1)
2892 egrid1(i,j) = log(pmid(i,j,lm)/egrid2(i,j)) &
2893 / log(pmid(i,j,lm)/pmid(i,j,lm-1))
2894
2895 IF (modelname == 'GFS' .OR. modelname == 'FV3R') THEN
2896 egrid1(i,j) = log(pmid(i,j,lm)/egrid2(i,j)) &
2897 / max(1.e-6,log(pmid(i,j,lm)/pmid(i,j,lm-1)))
2898 egrid1(i,j) =max(-10.0,min(egrid1(i,j), 10.0))
2899 IF ( abs(pmid(i,j,lm)-pmid(i,j,lm-1)) < 0.5 ) THEN
2900 egrid1(i,j) = -1.
2901 ENDIF
2902 ENDIF
2903
2904 END DO
2905 END DO
2906! Temperature
2907 IF (iget(321)>0) THEN
2908 grid1=spval
2909!$omp parallel do private(i,j)
2910 DO j=jsta,jend
2911 DO i=ista,iend
2912 IF(t(i,j,lm)<spval.and.t(i,j,lm-1)<spval.and.egrid1(i,j)<spval)&
2913 grid1(i,j) = t(i,j,lm)+(t(i,j,lm-1)-t(i,j,lm)) &
2914 * egrid1(i,j)
2915 ENDDO
2916 ENDDO
2917 if(grib=='grib2') then
2918 cfld=cfld+1
2919 fld_info(cfld)%ifld=iavblfld(iget(321))
2920 fld_info(cfld)%lvl=lvlsxml(1,iget(321))
2921!$omp parallel do private(i,j,ii,jj)
2922 do j=1,jend-jsta+1
2923 jj = jsta+j-1
2924 do i=1,iend-ista+1
2925 ii = ista+i-1
2926 datapd(i,j,cfld) = grid1(ii,jj)
2927 enddo
2928 enddo
2929 endif
2930! print *,'in miscln,iget(321,temp,sigmadata=',maxval(GRID1(1:im,jsta:jend)), &
2931! minval(GRID1(1:im,jsta:jend)),'grib=',grib
2932 ENDIF
2933! Potential Temperature
2934 IF (iget(322)>0) THEN
2935 grid2=spval
2936!$omp parallel do private(i,j)
2937 DO j=jsta,jend
2938 DO i=ista,iend
2939 IF(t(i,j,lm)<spval.and.t(i,j,lm-1)<spval.and.egrid1(i,j)<spval)&
2940 grid2(i,j) = t(i,j,lm)+(t(i,j,lm-1)-t(i,j,lm)) &
2941 * egrid1(i,j)
2942 ENDDO
2943 ENDDO
2944 CALL calpot(egrid2,grid2(ista:iend,jsta:jend),grid1(ista:iend,jsta:jend))
2945 if(grib=='grib2') then
2946 cfld=cfld+1
2947 fld_info(cfld)%ifld=iavblfld(iget(322))
2948 fld_info(cfld)%lvl=lvlsxml(1,iget(322))
2949!$omp parallel do private(i,j,ii,jj)
2950 do j=1,jend-jsta+1
2951 jj = jsta+j-1
2952 do i=1,iend-ista+1
2953 ii = ista+i-1
2954 datapd(i,j,cfld) = grid1(ii,jj)
2955 enddo
2956 enddo
2957 endif
2958 ENDIF
2959! RH
2960 IF (iget(323)>0) THEN
2961 grid1=spval
2962!$omp parallel do private(i,j,es1,qs1,rh1,es2,qs2,rh2)
2963 DO j=jsta,jend
2964 DO i=ista,iend
2965 IF(pmid(i,j,lm)<spval.and.pmid(i,j,lm-1)<spval.and.&
2966 q(i,j,lm)<spval.and.q(i,j,lm-1)<spval) THEN
2967 es1 = min(pmid(i,j,lm),fpvsnew(t(i,j,lm)))
2968 qs1 = con_eps*es1/(pmid(i,j,lm)+con_epsm1*es1)
2969 rh1 = q(i,j,lm)/qs1
2970 es2 = min(pmid(i,j,lm-1),fpvsnew(t(i,j,lm-1)))
2971 qs2 = con_eps*es2/(pmid(i,j,lm-1)+con_epsm1*es2)
2972 rh2 = q(i,j,lm-1)/qs2
2973 grid1(i,j) = (rh1+(rh2-rh1)*egrid1(i,j))*100.
2974 ENDIF
2975 ENDDO
2976 ENDDO
2977 CALL bound(grid1,d00,h100)
2978 if(grib=='grib2') then
2979 cfld=cfld+1
2980 fld_info(cfld)%ifld=iavblfld(iget(323))
2981 fld_info(cfld)%lvl=lvlsxml(1,iget(323))
2982!$omp parallel do private(i,j,ii,jj)
2983 do j=1,jend-jsta+1
2984 jj = jsta+j-1
2985 do i=1,iend-ista+1
2986 ii = ista+i-1
2987 datapd(i,j,cfld) = grid1(ii,jj)
2988 enddo
2989 enddo
2990 endif
2991 ENDIF
2992! U
2993 IF (iget(324)>0) THEN
2994 grid1=spval
2995!$omp parallel do private(i,j)
2996 DO j=jsta,jend
2997 DO i=ista,iend
2998 IF(uh(i,j,lm)<spval.and.uh(i,j,lm-1)<spval.and.egrid1(i,j)<spval)&
2999 grid1(i,j) = uh(i,j,lm)+(uh(i,j,lm-1)-uh(i,j,lm)) &
3000 * egrid1(i,j)
3001 ENDDO
3002 ENDDO
3003 if(grib=='grib2') then
3004 cfld=cfld+1
3005 fld_info(cfld)%ifld=iavblfld(iget(324))
3006 fld_info(cfld)%lvl=lvlsxml(1,iget(324))
3007!$omp parallel do private(i,j,ii,jj)
3008 do j=1,jend-jsta+1
3009 jj = jsta+j-1
3010 do i=1,iend-ista+1
3011 ii = ista+i-1
3012 datapd(i,j,cfld) = grid1(ii,jj)
3013 enddo
3014 enddo
3015 endif
3016 ENDIF
3017! V
3018 IF (iget(325)>0) THEN
3019 grid1=spval
3020!$omp parallel do private(i,j)
3021 DO j=jsta,jend
3022 DO i=ista,iend
3023 IF(vh(i,j,lm)<spval.and.vh(i,j,lm-1)<spval.and.egrid1(i,j)<spval)&
3024 grid1(i,j) = vh(i,j,lm)+(vh(i,j,lm-1)-vh(i,j,lm)) &
3025 * egrid1(i,j)
3026 ENDDO
3027 ENDDO
3028 if(grib=='grib2') then
3029 cfld=cfld+1
3030 fld_info(cfld)%ifld=iavblfld(iget(325))
3031 fld_info(cfld)%lvl=lvlsxml(1,iget(325))
3032!$omp parallel do private(i,j,ii,jj)
3033 do j=1,jend-jsta+1
3034 jj = jsta+j-1
3035 do i=1,iend-ista+1
3036 ii = ista+i-1
3037 datapd(i,j,cfld) = grid1(ii,jj)
3038 enddo
3039 enddo
3040 endif
3041 ENDIF
3042! OMGA
3043 IF (iget(326)>0) THEN
3044 grid1=spval
3045!$omp parallel do private(i,j)
3046 DO j=jsta,jend
3047 DO i=ista,iend
3048 IF(omga(i,j,lm)<spval.and.omga(i,j,lm-1)<spval.and.egrid1(i,j)<spval)&
3049 grid1(i,j) = omga(i,j,lm)+(omga(i,j,lm-1)-omga(i,j,lm))&
3050 * egrid1(i,j)
3051 ENDDO
3052 ENDDO
3053 if(grib=='grib2') then
3054 cfld=cfld+1
3055 fld_info(cfld)%ifld=iavblfld(iget(326))
3056 fld_info(cfld)%lvl=lvlsxml(1,iget(326))
3057!$omp parallel do private(i,j,ii,jj)
3058 do j=1,jend-jsta+1
3059 jj = jsta+j-1
3060 do i=1,iend-ista+1
3061 ii = ista+i-1
3062 datapd(i,j,cfld) = grid1(ii,jj)
3063 enddo
3064 enddo
3065 endif
3066 ENDIF
3067 END IF
3068
3069! MIXED LAYER CAPE AND CINS
3070!
3071 field1=.false.
3072 field2=.false.
3073!
3074 IF(iget(032)>0)THEN
3075 IF(lvls(3,iget(032))>0)field1=.true.
3076 ENDIF
3077 IF(iget(107)>0)THEN
3078 IF(lvls(3,iget(107))>0)field2=.true.
3079 ENDIF
3080!
3081 IF(iget(582)>0)THEN
3082 field1=.true.
3083 ENDIF
3084 IF(iget(583)>0)THEN
3085 field2=.true.
3086 ENDIF
3087
3088 IF(field1.OR.field2)THEN
3089 itype = 2
3090!
3091!$omp parallel do private(i,j)
3092 DO j=jsta,jend
3093 DO i=ista,iend
3094 egrid1(i,j) = -h99999
3095 egrid2(i,j) = -h99999
3096 lb2(i,j) = (lvlbnd(i,j,1) + lvlbnd(i,j,2) + &
3097 lvlbnd(i,j,3))/3
3098 p1d(i,j) = (pbnd(i,j,1) + pbnd(i,j,2) + pbnd(i,j,3))/3
3099 t1d(i,j) = (tbnd(i,j,1) + tbnd(i,j,2) + tbnd(i,j,3))/3
3100 q1d(i,j) = (qbnd(i,j,1) + qbnd(i,j,2) + qbnd(i,j,3))/3
3101 ENDDO
3102 ENDDO
3103!
3104 dpbnd = 0.
3105 CALL calcape(itype,dpbnd,p1d,t1d,q1d,lb2,egrid1, &
3106 egrid2,egrid3,egrid4,egrid5)
3107
3108 IF (iget(582)>0) THEN
3109! dong add missing value for cape
3110 grid1=spval
3111!$omp parallel do private(i,j)
3112 DO j=jsta,jend
3113 DO i=ista,iend
3114 IF(t1d(i,j) < spval) THEN
3115 grid1(i,j) = egrid1(i,j)
3116 IF (submodelname == 'RTMA')mlcape(i,j)=grid1(i,j)
3117 ENDIF
3118 ENDDO
3119 ENDDO
3120 CALL bound(grid1,d00,h99999)
3121 if(grib=='grib2') then
3122 cfld=cfld+1
3123 fld_info(cfld)%ifld=iavblfld(iget(582))
3124 fld_info(cfld)%lvl=lvlsxml(1,iget(582))
3125!$omp parallel do private(i,j,ii,jj)
3126 do j=1,jend-jsta+1
3127 jj = jsta+j-1
3128 do i=1,iend-ista+1
3129 ii = ista+i-1
3130 datapd(i,j,cfld) = grid1(ii,jj)
3131 enddo
3132 enddo
3133 endif
3134 ENDIF
3135 IF (iget(583)>0) THEN
3136! dong add missing value for cape
3137 grid1=spval
3138!$omp parallel do private(i,j)
3139 DO j=jsta,jend
3140 DO i=ista,iend
3141 IF(t1d(i,j) < spval) grid1(i,j) = - egrid2(i,j)
3142 ENDDO
3143 ENDDO
3144!
3145 CALL bound(grid1,d00,h99999)
3146!
3147!$omp parallel do private(i,j)
3148 DO j=jsta,jend
3149 DO i=ista,iend
3150 IF(t1d(i,j) < spval) THEN
3151 grid1(i,j) = - grid1(i,j)
3152 IF (submodelname == 'RTMA') mlcin(i,j)=grid1(i,j)
3153 ENDIF
3154 ENDDO
3155 ENDDO
3156!
3157 if(grib=='grib2') then
3158 cfld=cfld+1
3159 fld_info(cfld)%ifld=iavblfld(iget(583))
3160 fld_info(cfld)%lvl=lvlsxml(1,iget(583))
3161!$omp parallel do private(i,j,ii,jj)
3162 do j=1,jend-jsta+1
3163 jj = jsta+j-1
3164 do i=1,iend-ista+1
3165 ii = ista+i-1
3166 datapd(i,j,cfld) = grid1(ii,jj)
3167 enddo
3168 enddo
3169 endif
3170 ENDIF
3171 ENDIF
3172
3173! MIXED LAYER LIFTING CONDENSATION PRESSURE AND HEIGHT.
3174! EGRID1 IS LCL PRESSURE. EGRID2 IS LCL HEIGHT.
3175!
3176 IF ( (iget(109)>0).OR.(iget(110)>0) ) THEN
3177 CALL callcl(p1d,t1d,q1d,egrid1,egrid2)
3178 IF (iget(109)>0) THEN
3179 grid1=spval
3180 DO j=jsta,jend
3181 DO i=ista,iend
3182 IF(t1d(i,j) < spval) grid1(i,j)=egrid2(i,j)
3183 IF (submodelname == 'RTMA') mllcl(i,j) = grid1(i,j)
3184 ENDDO
3185 ENDDO
3186!
3187! ID(1:25) = 0
3188!
3189! CALL GRIBIT(IGET(109),1,
3190! X GRID1,IM,JM)
3191 ENDIF
3192!
3193! IF (IGET(110)>0) THEN
3194! DO J=JSTA,JEND
3195! DO I=ISTA,IEND
3196! GRID1(I,J)=EGRID1(I,J)
3197! ENDDO
3198! ENDDO
3199!
3200! ID(1:25) = 0
3201!
3202! CALL GRIBIT(IGET(110),1,
3203! X GRID1,IM,JM)
3204! ENDIF
3205 ENDIF
3206!
3207! MOST UNSTABLE CAPE-LOWEST 300 MB
3208!
3209 field1=.false.
3210 field2=.false.
3211!
3212 IF(iget(032)>0)THEN
3213 IF(lvls(4,iget(032))>0)field1=.true.
3214 ENDIF
3215
3216 IF(iget(107)>0)THEN
3217 IF(lvls(4,iget(107))>0)field2=.true.
3218 ENDIF
3219!
3220 IF(iget(584)>0)THEN
3221 field1=.true.
3222 ENDIF
3223 IF(iget(585)>0)THEN
3224 field2=.true.
3225 ENDIF
3226!
3227 IF(field1.OR.field2.OR.need_ifi)THEN
3228 itype = 1
3229 call allocate_cape_arrays
3230!
3231!$omp parallel do private(i,j)
3232 DO j=jsta,jend
3233 DO i=ista,iend
3234 egrid1(i,j) = -h99999
3235 egrid2(i,j) = -h99999
3236 ENDDO
3237 ENDDO
3238
3239 dpbnd = 300.e2
3240 CALL calcape(itype,dpbnd,p1d,t1d,q1d,lb2,egrid1, &
3241 egrid2,egrid3,egrid4,egrid5)
3242 IF (iget(584)>0 .or. need_ifi) THEN
3243! dong add missing value to cin
3244 grid1 = spval
3245!$omp parallel do private(i,j)
3246 DO j=jsta,jend
3247 DO i=ista,iend
3248 IF(t1d(i,j) < spval) THEN
3249 grid1(i,j) = egrid1(i,j)
3250 IF (submodelname == 'RTMA') mucape(i,j)=grid1(i,j)
3251 ENDIF
3252 ENDDO
3253 ENDDO
3254 CALL bound(grid1,d00,h99999)
3255! IF (SUBMODELNAME == 'RTMA') THEN
3256! CALL BOUND(MUCAPE,D00,H99999)
3257! ENDIF
3258!$omp parallel do private(i,j)
3259 DO j=jsta,jend
3260 DO i=ista,iend
3261 cape(i,j) = grid1(i,j)
3262 ENDDO
3263 ENDDO
3264 if(iget(584)>0 .and. grib=='grib2') then
3265 cfld=cfld+1
3266 fld_info(cfld)%ifld=iavblfld(iget(584))
3267 fld_info(cfld)%lvl=lvlsxml(1,iget(584))
3268!$omp parallel do private(i,j,ii,jj)
3269 do j=1,jend-jsta+1
3270 jj = jsta+j-1
3271 do i=1,iend-ista+1
3272 ii = ista+i-1
3273 datapd(i,j,cfld) = grid1(ii,jj)
3274 enddo
3275 enddo
3276 endif
3277
3278 ENDIF
3279
3280 IF (iget(585)>0 .or. need_ifi) THEN
3281! dong add missing value to cin
3282 grid1 = spval
3283!$omp parallel do private(i,j)
3284 DO j=jsta,jend
3285 DO i=ista,iend
3286 IF(t1d(i,j) < spval) grid1(i,j) = - egrid2(i,j)
3287 ENDDO
3288 ENDDO
3289 CALL bound(grid1,d00,h99999)
3290 DO j=jsta,jend
3291 DO i=ista,iend
3292 IF(t1d(i,j) < spval) THEN
3293 grid1(i,j) = - grid1(i,j)
3294 IF (submodelname == 'RTMA')THEN
3295 mucape(i,j) = grid1(i,j)
3296 muq1d(i,j) = q1d(i,j)
3297 ENDIF
3298 ENDIF
3299 ENDDO
3300 ENDDO
3301
3302!$omp parallel do private(i,j)
3303 DO j=jsta,jend
3304 DO i=ista,iend
3305 cin(i,j) = grid1(i,j)
3306 ENDDO
3307 ENDDO
3308
3309 if(iget(585)>0 .and. grib=='grib2') then
3310 cfld=cfld+1
3311 fld_info(cfld)%ifld=iavblfld(iget(585))
3312 fld_info(cfld)%lvl=lvlsxml(1,iget(585))
3313!$omp parallel do private(i,j,ii,jj)
3314 do j=1,jend-jsta+1
3315 jj = jsta+j-1
3316 do i=1,iend-ista+1
3317 ii = ista+i-1
3318 datapd(i,j,cfld) = grid1(ii,jj)
3319 enddo
3320 enddo
3321 endif
3322
3323 ENDIF
3324
3325! EQUILLIBRIUM HEIGHT
3326 IF (iget(443)>0) THEN
3327 grid1 = spval
3328!$omp parallel do private(i,j)
3329 DO j=jsta,jend
3330 DO i=ista,iend
3331 IF(t1d(i,j) < spval) grid1(i,j) = egrid4(i,j)
3332 ENDDO
3333 ENDDO
3334 if(grib=='grib2') then
3335 cfld=cfld+1
3336 fld_info(cfld)%ifld=iavblfld(iget(443))
3337 fld_info(cfld)%lvl=lvlsxml(1,iget(443))
3338!$omp parallel do private(i,j,ii,jj)
3339 do j=1,jend-jsta+1
3340 jj = jsta+j-1
3341 do i=1,iend-ista+1
3342 ii = ista+i-1
3343 datapd(i,j,cfld) = grid1(ii,jj)
3344 enddo
3345 enddo
3346 endif
3347 ENDIF
3348!Equilibrium Temperature
3349 IF (iget(982)>0) THEN
3350 DO j=jsta,jend
3351 DO i=ista,iend
3352 grid1(i,j) = teql(i,j)
3353 ENDDO
3354 ENDDO
3355 if(grib=='grib2') then
3356 cfld=cfld+1
3357 fld_info(cfld)%ifld=iavblfld(iget(982))
3358 fld_info(cfld)%lvl=lvlsxml(1,iget(982))
3359!$omp parallel do private(i,j,ii,jj)
3360 do j=1,jend-jsta+1
3361 jj = jsta+j-1
3362 do i=1,iend-ista+1
3363 ii = ista+i-1
3364 datapd(i,j,cfld) = grid1(ii,jj)
3365 enddo
3366 enddo
3367 endif
3368 ENDIF
3369
3370
3371! PRESSURE OF LEVEL FROM WHICH 300 MB MOST UNSTABLE CAPE
3372! PARCEL WAS LIFTED (eq. PRESSURE OF LEVEL OF HIGHEST THETA-E)
3373 IF (iget(246)>0) THEN
3374 grid1 = spval
3375!$omp parallel do private(i,j)
3376 DO j=jsta,jend
3377 DO i=ista,iend
3378 IF(t1d(i,j) < spval) grid1(i,j) = egrid3(i,j)
3379 ENDDO
3380 ENDDO
3381 CALL bound(grid1,d00,h99999)
3382! print *,'in miscln,PLPL=',maxval(grid1(1:im,jsta:jend)), &
3383! minval(grid1(1:im,jsta:jend))
3384 if(grib=='grib2') then
3385 cfld=cfld+1
3386 fld_info(cfld)%ifld=iavblfld(iget(246))
3387 fld_info(cfld)%lvl=lvlsxml(1,iget(246))
3388!$omp parallel do private(i,j,ii,jj)
3389 do j=1,jend-jsta+1
3390 jj = jsta+j-1
3391 do i=1,iend-ista+1
3392 ii = ista+i-1
3393 datapd(i,j,cfld) = grid1(ii,jj)
3394 enddo
3395 enddo
3396 endif
3397 ENDIF ! 246
3398
3399! GENERAL THUNDER PARAMETER ??? 458 ???
3400 IF (iget(444)>0) THEN
3401 grid1 = spval
3402!$omp parallel do private(i,j)
3403 DO j=jsta,jend
3404 DO i=ista,iend
3405 IF(cprate(i,j) < spval) THEN
3406 IF (cprate(i,j) > pthresh) THEN
3407 grid1(i,j) = egrid5(i,j)
3408 ELSE
3409 grid1(i,j) = 0
3410 ENDIF
3411 ENDIF
3412 ENDDO
3413 ENDDO
3414 CALL bound(grid1,d00,h99999)
3415 if(grib=='grib2') then
3416 cfld=cfld+1
3417 fld_info(cfld)%ifld=iavblfld(iget(444))
3418 fld_info(cfld)%lvl=lvlsxml(1,iget(444))
3419!$omp parallel do private(i,j,ii,jj)
3420 do j=1,jend-jsta+1
3421 jj = jsta+j-1
3422 do i=1,iend-ista+1
3423 ii = ista+i-1
3424 datapd(i,j,cfld) = grid1(ii,jj)
3425 enddo
3426 enddo
3427 endif
3428 ENDIF
3429 ENDIF
3430
3431
3432 IF (submodelname == 'RTMA')THEN
3433
3434!
3435! --- Effective (inflow) Layer (EL)
3436!
3437
3438 ALLOCATE(el_base(ista_2l:iend_2u,jsta_2l:jend_2u))
3439 ALLOCATE(el_tops(ista_2l:iend_2u,jsta_2l:jend_2u))
3440 ALLOCATE(found_base(ista_2l:iend_2u,jsta_2l:jend_2u))
3441 ALLOCATE(found_tops(ista_2l:iend_2u,jsta_2l:jend_2u))
3442!$omp parallel do private(i,j)
3443 DO j=jsta,jend
3444 DO i=ista,iend
3445 el_base(i,j) = lm
3446 el_tops(i,j) = lm
3447 found_base(i,j) = .false.
3448 found_tops(i,j) = .false.
3449 ENDDO
3450 ENDDO
3451
3452!
3453
3454 itype = 2
3455 dpbnd = 0.
3456
3457 DO l = lm, 1, -1
3458
3459! SET AIR PARCELS FOR LEVEL L
3460!$omp parallel do private(i,j)
3461 DO j=jsta,jend
3462 DO i=ista,iend
3463 egrid1(i,j) = -h99999
3464 egrid2(i,j) = -h99999
3465 idummy(i,j) = 0
3466 p1d(i,j) = pmid(i,j,l)
3467 t1d(i,j) = t(i,j,l)
3468 q1d(i,j) = q(i,j,l)
3469 ENDDO
3470 ENDDO
3471
3472!--- CALCULATE CAPE/CIN FOR ALL AIR PARCELS on LEVEL L
3473 IF (debugprint) WRITE(1000+me,'(A,I3)') &
3474 ' CALCULATING CAPE/CINS ON LEVEL:',l
3475 CALL calcape(itype,dpbnd,p1d,t1d,q1d,idummy,egrid1, &
3476 egrid2,egrid3,egrid4,egrid5)
3477
3478!--- CHECK CAPE/CIN OF EACH AIR PARCELS WITH EL CRITERIA
3479!$omp parallel do private(i,j)
3480 DO j=jsta,jend
3481 DO i=ista,iend
3482 IF ( .NOT. found_base(i,j) ) THEN
3483 IF ( egrid1(i,j) >= 100. .AND. egrid2(i,j) >= -250. ) THEN
3484 el_base(i,j) = l
3485 found_base(i,j) = .true.
3486 ELSE
3487 el_base(i,j) = lm
3488 found_base(i,j) = .false.
3489 END IF
3490 ELSE
3491 IF ( .NOT. found_tops(i,j) ) THEN
3492 IF ( egrid1(i,j) < 100. .OR. egrid2(i,j) < -250. ) THEN
3493 el_tops(i,j) = l + 1
3494 found_tops(i,j) = .true.
3495 ELSE
3496 el_tops(i,j) = lm
3497 found_tops(i,j) = .false.
3498 END IF
3499 END IF
3500 END IF
3501 ENDDO
3502 ENDDO
3503
3504 END DO ! L
3505
3506
3507 IF (ALLOCATED(found_base)) DEALLOCATE(found_base)
3508 IF (ALLOCATED(found_tops)) DEALLOCATE(found_tops)
3509
3510 IF (debugprint) THEN
3511 WRITE(im_ch,'(I5.5)') im
3512 WRITE(jsta_ch,'(I5.5)') jsta
3513 WRITE(jend_ch,'(I5.5)') jend
3514 effl_fname="EFFL_NEW_"//im_ch//"_"//jsta_ch//"_"//jend_ch &
3515 //".dat"
3516 effl_fname2="EFFL_NEW_LVLS_"//im_ch//"_"//jsta_ch//"_"//jend_ch &
3517 //".dat"
3518 iunit=10000+jsta
3519 iunit2=20000+jsta
3520 irec=0
3521 irec2=0
3522 OPEN(iunit,file=trim(adjustl(effl_fname)),form='FORMATTED')
3523
3524
3525 DO j=jsta,jend
3526 DO i=ista,iend
3527 irec = irec + 1
3528 irec2 = irec2 + 1
3529 WRITE(iunit,'(1x,I6,2x,I6,2(2x,I6,2x,F12.3))') i, j, &
3530 el_base(i,j),pmid(i,j,el_base(i,j)), &
3531 el_tops(i,j),pmid(i,j,el_tops(i,j))
3532 END DO
3533 ENDDO
3534 CLOSE(iunit)
3535 ENDIF
3536
3537 IF(ALLOCATED(tpar_base)) DEALLOCATE(tpar_base)
3538 IF(ALLOCATED(tpar_tops)) DEALLOCATE(tpar_tops)
3539
3540 ENDIF
3541!
3542! EXPAND HRRR CAPE/CIN RELATED VARIABLES
3543!
3544! CAPE AND CINS 0-3KM, FOLLOW ML PROCEDURE WITH HEIGHT 0-3KM
3545!
3546 field1=.false.
3547 field2=.false.
3548!
3549 IF(iget(032)>0)THEN
3550 IF(lvls(3,iget(032))>0)field1=.true.
3551 ENDIF
3552 IF(iget(107)>0)THEN
3553 IF(lvls(3,iget(107))>0)field2=.true.
3554 ENDIF
3555!
3556 IF(iget(950)>0)THEN
3557 field1=.true.
3558 ENDIF
3559 IF(iget(951)>0)THEN
3560 field2=.true.
3561 ENDIF
3562 IF(modelname == "FV3R" .and. submodelname == "RTMA") THEN
3563 field1=.true.
3564 field2=.true.
3565 ENDIF
3566!
3567! IF(FIELD1)ITYPE=2
3568! IF(FIELD2)ITYPE=2
3569
3570 IF(field1.OR.field2)THEN
3571 itype = 2
3572
3573 call allocate_cape_arrays
3574
3575!
3576!$omp parallel do private(i,j)
3577 DO j=jsta,jend
3578 DO i=ista,iend
3579 egrid1(i,j) = -h99999
3580 egrid2(i,j) = -h99999
3581 egrid3(i,j) = -h99999
3582 egrid4(i,j) = -h99999
3583 egrid5(i,j) = -h99999
3584 egrid6(i,j) = -h99999
3585 egrid7(i,j) = -h99999
3586 egrid8(i,j) = -h99999
3587! ENDDO
3588! ENDDO
3589! DO J=JSTA,JEND
3590! DO I=ISTA,IEND
3591 lb2(i,j) = (lvlbnd(i,j,1) + lvlbnd(i,j,2) + &
3592 lvlbnd(i,j,3))/3
3593 p1d(i,j) = (pbnd(i,j,1) + pbnd(i,j,2) + pbnd(i,j,3))/3
3594 t1d(i,j) = (tbnd(i,j,1) + tbnd(i,j,2) + tbnd(i,j,3))/3
3595 q1d(i,j) = (qbnd(i,j,1) + qbnd(i,j,2) + qbnd(i,j,3))/3
3596 ENDDO
3597 ENDDO
3598
3599 dpbnd = 0.
3600 CALL calcape2(itype,dpbnd,p1d,t1d,q1d,lb2, &
3601 egrid1,egrid2,egrid3,egrid4,egrid5, &
3602 egrid6,egrid7,egrid8)
3603
3604! CAPE1, CINS2, LFC3, ESRHL4,ESRHH5,
3605! DCAPE6,DGLD7, ESP8)
3606!
3607 IF (iget(950)>0) THEN
3608! dong add missing value for cape
3609 grid1=spval
3610!$omp parallel do private(i,j)
3611 DO j=jsta,jend
3612 DO i=ista,iend
3613 IF(t1d(i,j) < spval) grid1(i,j) = egrid1(i,j)
3614 ENDDO
3615 ENDDO
3616 CALL bound(grid1,d00,h99999)
3617 if(grib=='grib2') then
3618 cfld=cfld+1
3619 fld_info(cfld)%ifld=iavblfld(iget(950))
3620 fld_info(cfld)%lvl=lvlsxml(1,iget(950))
3621!$omp parallel do private(i,j,ii,jj)
3622 do j=1,jend-jsta+1
3623 jj = jsta+j-1
3624 do i=1,iend-ista+1
3625 ii = ista+i-1
3626 datapd(i,j,cfld) = grid1(ii,jj)
3627 enddo
3628 enddo
3629 endif
3630 ENDIF !950
3631!
3632 IF (iget(951)>0) THEN
3633! dong add missing value for cape
3634 grid1=spval
3635!$omp parallel do private(i,j)
3636 DO j=jsta,jend
3637 DO i=ista,iend
3638 IF(t1d(i,j) < spval) grid1(i,j) = - egrid2(i,j)
3639 ENDDO
3640 ENDDO
3641!
3642 CALL bound(grid1,d00,h99999)
3643!
3644!$omp parallel do private(i,j)
3645 DO j=jsta,jend
3646 DO i=ista,iend
3647 IF(t1d(i,j) < spval) grid1(i,j) = - grid1(i,j)
3648 ENDDO
3649 ENDDO
3650!
3651 if(grib=='grib2') then
3652 cfld=cfld+1
3653 fld_info(cfld)%ifld=iavblfld(iget(951))
3654 fld_info(cfld)%lvl=lvlsxml(1,iget(951))
3655!$omp parallel do private(i,j,ii,jj)
3656 do j=1,jend-jsta+1
3657 jj = jsta+j-1
3658 do i=1,iend-ista+1
3659 ii = ista+i-1
3660 datapd(i,j,cfld) = grid1(ii,jj)
3661 enddo
3662 enddo
3663 endif
3664 ENDIF !951
3665
3666! LFC HEIGHT
3667
3668 IF (iget(952)>0) THEN
3669 grid1=spval
3670!$omp parallel do private(i,j)
3671 DO j=jsta,jend
3672 DO i=ista,iend
3673 IF(t1d(i,j) < spval) grid1(i,j) = egrid3(i,j)
3674 ENDDO
3675 ENDDO
3676 CALL bound(grid1,d00,h99999)
3677 if(grib=='grib2') then
3678 cfld=cfld+1
3679 fld_info(cfld)%ifld=iavblfld(iget(952))
3680 fld_info(cfld)%lvl=lvlsxml(1,iget(952))
3681!$omp parallel do private(i,j,ii,jj)
3682 do j=1,jend-jsta+1
3683 jj = jsta+j-1
3684 do i=1,iend-ista+1
3685 ii = ista+i-1
3686 datapd(i,j,cfld) = grid1(ii,jj)
3687 enddo
3688 enddo
3689 endif
3690 ENDIF !952
3691
3692
3693! EFFECTIVE STORM RELATIVE HELICITY AND STORM MOTION.
3694
3695 allocate(ust(ista_2l:iend_2u,jsta_2l:jend_2u),vst(ista_2l:iend_2u,jsta_2l:jend_2u), &
3696 heli(ista_2l:iend_2u,jsta_2l:jend_2u,2))
3697 allocate(llow(ista_2l:iend_2u,jsta_2l:jend_2u),lupp(ista_2l:iend_2u,jsta_2l:jend_2u), &
3698 cangle(ista_2l:iend_2u,jsta_2l:jend_2u))
3699 allocate(llow_zint(ista_2l:iend_2u,jsta_2l:jend_2u), &
3700 ieql_zint(ista_2l:iend_2u,jsta_2l:jend_2u),z_temp(ista_2l:iend_2u,jsta_2l:jend_2u))
3701 allocate(midcal(ista_2l:iend_2u,jsta_2l:jend_2u,1:lm))
3702 allocate(z_midcal(ista_2l:iend_2u,jsta_2l:jend_2u))
3703 iget1 = iget(953)
3704 iget2 = -1
3705 iget3 = -1
3706 if (iget1 > 0) then
3707 iget2 = lvls(1,iget1)
3708 iget3 = lvls(2,iget1)
3709 endif
3710 if(me==0) write(*,*) '953 ',iget1,iget2,iget3
3711 IF (iget1 > 0 .OR. iget(162) > 0 .OR. iget(953) > 0) THEN
3712 depth(1) = 3000.0
3713 depth(2) = 1000.0
3714 IF (submodelname == 'RTMA') THEN
3715!--- IF USSING EL BASE & TOP COMPUTED BY NEW SCHEME FOR THE
3716!RELATED VARIABLES
3717!$omp parallel do private(i,j)
3718 DO j=jsta,jend
3719 DO i=ista,iend
3720 llow(i,j) = el_base(i,j)
3721 lupp(i,j) = el_tops(i,j)
3722 ENDDO
3723 ENDDO
3724 ELSE
3725!$omp parallel do private(i,j)
3726 DO j=jsta,jend
3727 DO i=ista,iend
3728 llow(i,j) = int(egrid4(i,j))
3729 lupp(i,j) = int(egrid5(i,j))
3730 ENDDO
3731 ENDDO
3732 ENDIF
3733!--- OUTPUT EL BASE & TOP COMPUTED BY OLD SCHEME
3734 IF (debugprint) THEN
3735 WRITE(im_ch,'(I5.5)') im
3736 WRITE(jsta_ch,'(I5.5)') jsta
3737 WRITE(jend_ch,'(I5.5)') jend
3738 effl_fname="EFFL_OLD_"//im_ch//"_"//jsta_ch//"_"//jend_ch &
3739 //".dat"
3740 iunit=10000+jsta
3741 irec=0
3742 OPEN(iunit,file=trim(adjustl(effl_fname)),form='FORMATTED')
3743 DO j=jsta,jend
3744 DO i=ista,iend
3745 irec = irec + 1
3746 WRITE(iunit,'(1x,I6,2x,I6,2(2x,I6,2x,F12.3))') i, j, &
3747 llow(i,j),pmid(i,j,llow(i,j)), &
3748 lupp(i,j),pmid(i,j,lupp(i,j))
3749 END DO
3750 ENDDO
3751 CLOSE(iunit)
3752 ENDIF
3753
3754
3755 CALL calhel2(llow,lupp,depth,ust,vst,heli,cangle)
3756!
3757 IF (iget2 > 0) then
3758!$omp parallel do private(i,j)
3759 DO j=jsta,jend
3760 DO i=ista,iend
3761 grid1(i,j) = heli(i,j,1)
3762 ENDDO
3763 ENDDO
3764 if(grib=='grib2') then
3765 cfld=cfld+1
3766 fld_info(cfld)%ifld=iavblfld(iget1)
3767 fld_info(cfld)%lvl=lvlsxml(1,iget1)
3768!$omp parallel do private(i,j,ii,jj)
3769 do j=1,jend-jsta+1
3770 jj = jsta+j-1
3771 do i=1,iend-ista+1
3772 ii = ista+i-1
3773 datapd(i,j,cfld) = grid1(ii,jj)
3774 enddo
3775 enddo
3776 endif
3777 ENDIF
3778
3779 ENDIF !953
3780
3781
3782 IF (submodelname == 'RTMA') THEN !Start RTMA block
3783
3784!EL field allocation
3785
3786 allocate(eshr(ista_2l:iend_2u,jsta_2l:jend_2u),uvect(ista_2l:iend_2u,jsta_2l:jend_2u),&
3787 vvect(ista_2l:iend_2u,jsta_2l:jend_2u),htsfc(ista_2l:iend_2u,jsta_2l:jend_2u))
3788 allocate(effust(ista_2l:iend_2u,jsta_2l:jend_2u),effvst(ista_2l:iend_2u,jsta_2l:jend_2u),&
3789 esrh(ista_2l:iend_2u,jsta_2l:jend_2u))
3790
3791!
3792 DO j=jsta,jend
3793 DO i=ista,iend
3794 maxthe(i,j)=-h99999
3795 the(i,j)=-h99999
3796 maxthepos(i,j)=0
3797 muq1d(i,j) = 0.
3798 ENDDO
3799 ENDDO
3800 DO l=lm,1,-1
3801
3802 DO j=jsta,jend
3803 DO i=ista,iend
3804 egrid1(i,j) = -h99999
3805 p1d(i,j)=pmid(i,j,l)
3806 t1d(i,j)=t(i,j,l)
3807 q1d(i,j)=q(i,j,l)
3808 ENDDO
3809 ENDDO
3810 CALL calthte(p1d,t1d,q1d,egrid1)
3811 DO j=jsta,jend
3812 DO i=ista,iend
3813 the(i,j)=egrid1(i,j)
3814 IF(the(i,j)>=maxthe(i,j))THEN
3815 maxthe(i,j)=the(i,j)
3816 maxthepos(i,j)=l
3817 muq1d(i,j) = q(i,j,l) ! save the Q of air parcel with max theta-e (MU Parcel)
3818 ENDIF
3819 ENDDO
3820 ENDDO
3821 DO j=jsta,jend
3822 DO i=ista,iend
3823 llow(i,j) = el_base(i,j)
3824 llow_zint(i,j)=zint(i,j,llow(i,j))
3825 ieql_zint(i,j)=zint(i,j,ieql(i,j))
3826 z_temp(i,j)=llow_zint(i,j)+d50*(ieql_zint(i,j)-llow_zint(i,j))
3827 midcal(i,j,l)=abs(zint(i,j,l)-z_temp(i,j))
3828 ENDDO
3829 ENDDO
3830 ENDDO
3831 z_midcal=minloc(midcal,dim=3)
3832
3833!
3834!get surface height
3835 IF(gridtype == 'E')THEN
3836 jvn = 1
3837 jvs = -1
3838 do j=jsta,jend
3839 ive(j) = mod(j,2)
3840 ivw(j) = ive(j)-1
3841 enddo
3842 istart = ista_m
3843 istop = iend_m
3844 jstart = jsta_m
3845 jstop = jend_m
3846 ELSE IF(gridtype == 'B')THEN
3847 jvn = 1
3848 jvs = 0
3849 do j=jsta,jend
3850 ive(j)=1
3851 ivw(j)=0
3852 enddo
3853 istart = ista_m
3854 istop = iend_m
3855 jstart = jsta_m
3856 jstop = jend_m
3857 ELSE
3858 jvn = 0
3859 jvs = 0
3860 do j=jsta,jend
3861 ive(j) = 0
3862 ivw(j) = 0
3863 enddo
3864 istart = ista
3865 istop = iend
3866 jstart = jsta
3867 jstop = jend
3868 END IF
3869
3870 IF(gridtype /= 'A') CALL exch(fis(ista_2l:iend_2u,jsta_2l:jend_2u))
3871 DO j=jstart,jstop
3872 DO i=istart,istop
3873 ie = i+ive(j)
3874 iw = i+ivw(j)
3875 jn = j+jvn
3876 js = j+jvs
3877 IF (gridtype=='B')THEN
3878 htsfc(i,j)=(0.25/g)*(fis(iw,j)+fis(ie,j)+fis(i,jn)+fis(ie,jn))
3879 ELSE
3880 htsfc(i,j)=(0.25/g)*(fis(iw,j)+fis(ie,j)+fis(i,jn)+fis(i,js))
3881 ENDIF
3882 ENDDO
3883 ENDDO
3884
3885!Height of effbot
3886 IF (iget(979)>0) THEN
3887 grid1=spval
3888 DO j=jsta,jend
3889 DO i=ista,iend
3890 IF(zint(i,j,llow(i,j))<spval.and.htsfc(i,j)<spval)&
3891 grid1(i,j) = zint(i,j,llow(i,j)) - htsfc(i,j)
3892 ENDDO
3893 ENDDO
3894 if(grib=='grib2') then
3895 cfld=cfld+1
3896 fld_info(cfld)%ifld=iavblfld(iget(979))
3897 fld_info(cfld)%lvl=lvlsxml(1,iget(979))
3898!$omp parallel do private(i,j,ii,jj)
3899 do j=1,jend-jsta+1
3900 jj = jsta+j-1
3901 do i=1,iend-ista+1
3902 ii = ista+i-1
3903 datapd(i,j,cfld) = grid1(ii,jj)
3904 enddo
3905 enddo
3906 endif
3907 ENDIF
3908!Height of effbot
3909 IF (iget(980)>0) THEN
3910 grid1=spval
3911 DO j=jsta,jend
3912 DO i=ista,iend
3913 IF(zint(i,j,lupp(i,j))<spval.and.htsfc(i,j)<spval)&
3914 grid1(i,j) = zint(i,j,lupp(i,j)) - htsfc(i,j)
3915 ENDDO
3916 ENDDO
3917 if(grib=='grib2') then
3918 cfld=cfld+1
3919 fld_info(cfld)%ifld=iavblfld(iget(980))
3920 fld_info(cfld)%lvl=lvlsxml(1,iget(980))
3921!$omp parallel do private(i,j,ii,jj)
3922 do j=1,jend-jsta+1
3923 jj = jsta+j-1
3924 do i=1,iend-ista+1
3925 ii = ista+i-1
3926 datapd(i,j,cfld) = grid1(ii,jj)
3927 enddo
3928 enddo
3929 endif
3930 ENDIF
3931
3932!U inflow based to 50% EL shear vector
3933
3934 IF (iget(983)>0) THEN
3935 grid1=spval
3936 DO j=jsta,jend
3937 DO i=ista,iend
3938 IF(llow(i,j)<spval.and.lupp(i,j)<spval) THEN
3939 uvect(i,j)=uh(i,j,z_midcal(i,j))-uh(i,j,llow(i,j))
3940 grid1(i,j)=uvect(i,j)
3941 ENDIF
3942 ENDDO
3943 ENDDO
3944 if(grib=='grib2') then
3945 cfld=cfld+1
3946 fld_info(cfld)%ifld=iavblfld(iget(983))
3947 fld_info(cfld)%lvl=lvlsxml(1,iget(983))
3948!$omp parallel do private(i,j,ii,jj)
3949 do j=1,jend-jsta+1
3950 jj = jsta+j-1
3951 do i=1,iend-ista+1
3952 ii = ista+i-1
3953 datapd(i,j,cfld) = grid1(ii,jj)
3954 enddo
3955 enddo
3956 endif
3957 ENDIF
3958
3959!V inflow based to 50% EL shear vector
3960 IF (iget(984)>0) THEN
3961 grid1=spval
3962 DO j=jsta,jend
3963 DO i=ista,iend
3964 IF(llow(i,j)<spval.and.lupp(i,j)<spval) THEN
3965 vvect(i,j)=vh(i,j,z_midcal(i,j))-vh(i,j,llow(i,j))
3966 grid1(i,j)=vvect(i,j)
3967 ENDIF
3968 ENDDO
3969 ENDDO
3970 if(grib=='grib2') then
3971 cfld=cfld+1
3972 fld_info(cfld)%ifld=iavblfld(iget(984))
3973 fld_info(cfld)%lvl=lvlsxml(1,iget(984))
3974!$omp parallel do private(i,j,ii,jj)
3975 do j=1,jend-jsta+1
3976 jj = jsta+j-1
3977 do i=1,iend-ista+1
3978 ii = ista+i-1
3979 datapd(i,j,cfld) = grid1(ii,jj)
3980 enddo
3981 enddo
3982 endif
3983 ENDIF
3984
3985!Inflow based (ESFC) to (50%) EL shear magnitude
3986 IF (iget(985)>0) THEN
3987 grid1=spval
3988 DO j=jsta,jend
3989 DO i=ista,iend
3990 IF(uvect(i,j)<spval.and.vvect(i,j)<spval) THEN
3991 eshr(i,j)=sqrt((uvect(i,j)**2)+(vvect(i,j))**2)
3992 !effshear
3993 !calc
3994 grid1(i,j)=eshr(i,j) !Effective
3995 ENDIF
3996 ENDDO
3997 ENDDO
3998 if(grib=='grib2') then
3999 cfld=cfld+1
4000 fld_info(cfld)%ifld=iavblfld(iget(985))
4001 fld_info(cfld)%lvl=lvlsxml(1,iget(985))
4002!$omp parallel do private(i,j,ii,jj)
4003 do j=1,jend-jsta+1
4004 jj = jsta+j-1
4005 do i=1,iend-ista+1
4006 ii = ista+i-1
4007 datapd(i,j,cfld) = grid1(ii,jj)
4008 enddo
4009 enddo
4010 endif
4011 ENDIF
4012
4013! Effective Helicity
4014!$omp parallel do private(i,j)
4015 DO j=jsta,jend
4016 DO i=ista,iend
4017 llow(i,j) = el_base(i,j)
4018 lupp(i,j) = el_tops(i,j)
4019 ENDDO
4020 ENDDO
4021
4022 CALL calhel3(llow,lupp,effust,effvst,esrh)
4023
4024!U Bunkers Effective right motion
4025!
4026
4027 IF (iget(986)>0) THEN
4028 grid1=spval
4029 DO j=jsta,jend
4030 DO i=ista,iend
4031 IF(llow(i,j)<spval.and.lupp(i,j)<spval)&
4032 grid1(i,j)=effust(i,j)
4033 ENDDO
4034 ENDDO
4035 if(grib=='grib2') then
4036 cfld=cfld+1
4037 fld_info(cfld)%ifld=iavblfld(iget(986))
4038 fld_info(cfld)%lvl=lvlsxml(1,iget(986))
4039! $omp parallel do private(i,j,ii,jj)
4040 do j=1,jend-jsta+1
4041 jj = jsta+j-1
4042 do i=1,iend-ista+1
4043 ii = ista+i-1
4044 datapd(i,j,cfld) = grid1(ii,jj)
4045 enddo
4046 enddo
4047 endif
4048 ENDIF
4049
4050!V Bunkers Effective right motion
4051 IF (iget(987)>0) THEN
4052 grid1=spval
4053 DO j=jsta,jend
4054 DO i=ista,iend
4055 IF(llow(i,j)<spval.and.lupp(i,j)<spval)&
4056 grid1(i,j)=effvst(i,j)
4057 ENDDO
4058 ENDDO
4059 if(grib=='grib2') then
4060 cfld=cfld+1
4061 fld_info(cfld)%ifld=iavblfld(iget(987))
4062 fld_info(cfld)%lvl=lvlsxml(1,iget(987))
4063! $omp parallel do private(i,j,ii,jj)
4064 do j=1,jend-jsta+1
4065 jj = jsta+j-1
4066 do i=1,iend-ista+1
4067 ii = ista+i-1
4068 datapd(i,j,cfld) = grid1(ii,jj)
4069 enddo
4070 enddo
4071 endif
4072 ENDIF
4073
4074!Effective layer helicity
4075 IF (iget(988)>0) THEN
4076 grid1=spval
4077 DO j=jsta,jend
4078 DO i=ista,iend
4079 IF(llow(i,j)<spval.and.lupp(i,j)<spval)&
4080 grid1(i,j)=esrh(i,j)
4081 ENDDO
4082 ENDDO
4083 if(grib=='grib2') then
4084 cfld=cfld+1
4085 fld_info(cfld)%ifld=iavblfld(iget(988))
4086 fld_info(cfld)%lvl=lvlsxml(1,iget(988))
4087! $omp parallel do private(i,j,ii,jj)
4088 do j=1,jend-jsta+1
4089 jj = jsta+j-1
4090 do i=1,iend-ista+1
4091 ii = ista+i-1
4092 datapd(i,j,cfld) = grid1(ii,jj)
4093 enddo
4094 enddo
4095 endif
4096 ENDIF
4097
4098!Effective Layer Tornado Parameter
4099 IF (iget(989)>0) THEN
4100 DO j=jsta,jend
4101 DO i=ista,iend
4102 IF (mllcl(i,j)>d2000) THEN
4103 mllcltmp=d00
4104 ELSEIF (mllcl(i,j)<d1000) THEN
4105 mllcltmp=1.0
4106 ELSE
4107 mllcltmp=((d2000-mllcl(i,j))/d1000)
4108 ENDIF
4109 IF (eshr(i,j)<12.5) THEN
4110 eshrtmp=d00
4111 ELSEIF (eshr(i,j)>30.0) THEN
4112 eshrtmp=1.5
4113 ELSE
4114 eshrtmp=(eshr(i,j)/20.)
4115 ENDIF
4116 IF (mlcin(i,j)>-50.) THEN
4117 mlcintmp=1.0
4118 ELSEIF (mlcin(i,j)<-200.) THEN
4119 mlcintmp=d00
4120 ELSE
4121 mlcintmp=(200.+mlcin(i,j))/150.
4122 ENDIF
4123 stp=(mlcape(i,j)/d1500)*mllcltmp*(esrh(i,j)/150.)*&
4124 eshrtmp*mlcintmp
4125 grid1(i,j) = spval
4126 IF(llow(i,j)<spval.and.lupp(i,j)<spval) THEN
4127 IF (stp>0) THEN
4128 grid1(i,j)=stp
4129 ELSE
4130 grid1(i,j)=d00
4131 ENDIF
4132 ENDIF
4133 ENDDO
4134 ENDDO
4135 if(grib=='grib2') then
4136 cfld=cfld+1
4137 fld_info(cfld)%ifld=iavblfld(iget(989))
4138 fld_info(cfld)%lvl=lvlsxml(1,iget(989))
4139! $omp parallel do private(i,j,ii,jj)
4140 do j=1,jend-jsta+1
4141 jj = jsta+j-1
4142 do i=1,iend-ista+1
4143 ii = ista+i-1
4144 datapd(i,j,cfld) = grid1(ii,jj)
4145 enddo
4146 enddo
4147 endif
4148 ENDIF
4149
4150!Fixed Layer Tornado Parameter
4151 IF (iget(990)>0) THEN
4152 DO j=jsta,jend
4153 DO i=ista,iend
4154 llmh = nint(lmh(i,j))
4155 p1d(i,j) = pmid(i,j,llmh)
4156 t1d(i,j) = t(i,j,llmh)
4157 q1d(i,j) = q(i,j,llmh)
4158 ENDDO
4159 ENDDO
4160 CALL callcl(p1d,t1d,q1d,egrid1,egrid2)
4161 DO j=jsta,jend
4162 DO i=ista,iend
4163 slcl(i,j)=egrid2(i,j)
4164 ENDDO
4165 ENDDO
4166 itype = 1
4167 dpbnd = 10.e2
4168 dummy = 0.
4169 idummy = 0
4170 CALL calcape(itype,dpbnd,dummy,dummy,dummy,&
4171 idummy,egrid1,egrid2,&
4172 egrid3,dummy,dummy)
4173
4174 DO j=jsta,jend
4175 DO i=ista,iend
4176 IF (slcl(i,j)>d2000) THEN
4177 slcltmp=d00
4178 ELSEIF (slcl(i,j)<=d1000) THEN
4179 slcltmp=1.0
4180 ELSE
4181 slcltmp=((d2000-slcl(i,j))/d1000)
4182 ENDIF
4183 IF (fshr(i,j)<12.5) THEN
4184 fshrtmp=d00
4185 ELSEIF (fshr(i,j)>30.0) THEN
4186 fshrtmp=1.5
4187 ELSE
4188 fshrtmp=(fshr(i,j)/20.)
4189 ENDIF
4190 IF (egrid2(i,j)>-50.) THEN
4191 scintmp=1.0
4192 ELSEIF (egrid2(i,j)<-200.) THEN
4193 scintmp=d00
4194 ELSE
4195 scintmp=((200.+egrid2(i,j)/150.))
4196 ENDIF
4197 stp=(egrid1(i,j)/d1500)*slcltmp*(heli(i,j,2)/150.)*&
4198 fshrtmp*scintmp
4199 grid1(i,j) = spval
4200 IF(t1d(i,j) < spval) THEN
4201 IF (stp>0) THEN
4202 grid1(i,j)=stp
4203 ELSE
4204 grid1(i,j)=d00
4205 ENDIF
4206 ENDIF
4207 ENDDO
4208 ENDDO
4209 if(grib=='grib2') then
4210 cfld=cfld+1
4211 fld_info(cfld)%ifld=iavblfld(iget(990))
4212 fld_info(cfld)%lvl=lvlsxml(1,iget(990))
4213! $omp parallel do private(i,j,ii,jj)
4214 do j=1,jend-jsta+1
4215 jj = jsta+j-1
4216 do i=1,iend-ista+1
4217 ii = ista+i-1
4218 datapd(i,j,cfld) = grid1(ii,jj)
4219 enddo
4220 enddo
4221 endif
4222 ENDIF
4223
4224!Effective Layer Supercell Parameter
4225 IF (iget(991)>0) THEN
4226 DO j=jsta,jend
4227 DO i=ista,iend
4228 IF (eshr(i,j)<10.) THEN
4229 eshrtmp=d00
4230 ELSEIF (eshr(i,j)>20.0) THEN
4231 eshrtmp=1
4232 ELSE
4233 eshrtmp=(eshr(i,j)/20.)
4234 ENDIF
4235 IF (mucin(i,j)>-40.) THEN
4236 mucintmp=1.0
4237 ELSE
4238 mucintmp=(-40./mucin(i,j))
4239 ENDIF
4240 stp=(mucape(i,j)/d1000)*(esrh(i,j)/50.)*&
4241 eshrtmp*mucintmp
4242 grid1(i,j) = spval
4243 IF(t1d(i,j) < spval) THEN
4244 IF (stp>0) THEN
4245 grid1(i,j)=stp
4246 ELSE
4247 grid1(i,j)=d00
4248 ENDIF
4249 ENDIF
4250 ENDDO
4251 ENDDO
4252 if(grib=='grib2') then
4253 cfld=cfld+1
4254 fld_info(cfld)%ifld=iavblfld(iget(991))
4255 fld_info(cfld)%lvl=lvlsxml(1,iget(991))
4256! $omp parallel do private(i,j,ii,jj)
4257 do j=1,jend-jsta+1
4258 jj = jsta+j-1
4259 do i=1,iend-ista+1
4260 ii = ista+i-1
4261 datapd(i,j,cfld) = grid1(ii,jj)
4262 enddo
4263 enddo
4264 endif
4265 ENDIF
4266
4267!Mixed Layer (100 mb) Virtual LFC
4268
4269 IF (iget(992)>0) THEN
4270!$omp parallel do private(i,j)
4271 DO j=jsta,jend
4272 DO i=ista,iend
4273 egrid1(i,j) = -h99999
4274 egrid2(i,j) = -h99999
4275 egrid3(i,j) = -h99999
4276 egrid4(i,j) = -h99999
4277 egrid5(i,j) = -h99999
4278 egrid6(i,j) = -h99999
4279 egrid7(i,j) = -h99999
4280 egrid8(i,j) = -h99999
4281 lb2(i,j) = (lvlbnd(i,j,1) + lvlbnd(i,j,2) + &
4282 lvlbnd(i,j,3))/3
4283 p1d(i,j) = (pbnd(i,j,1) + pbnd(i,j,2) + pbnd(i,j,3))/3
4284 t1d(i,j) = (tvirtual(tbnd(i,j,1),qbnd(i,j,1)) + &
4285 tvirtual(tbnd(i,j,2),qbnd(i,j,2)) + &
4286 tvirtual(tbnd(i,j,3),qbnd(i,j,3)))/3
4287 q1d(i,j) = (qbnd(i,j,1) + qbnd(i,j,2) + qbnd(i,j,3))/3
4288 ENDDO
4289 ENDDO
4290
4291 dpbnd = 0.
4292 itype = 2
4293! EGRID3 is Virtual LFC
4294 CALL calcape2(itype,dpbnd,p1d,t1d,q1d,lb2, &
4295 egrid1,egrid2,egrid3,egrid4,egrid5, &
4296 egrid6,egrid7,egrid8)
4297
4298 grid1=spval
4299 DO j=jsta,jend
4300 DO i=ista,iend
4301 IF(t1d(i,j) < spval) grid1(i,j) = egrid3(i,j)
4302 ENDDO
4303 ENDDO
4304 CALL bound(grid1,d00,h99999)
4305 if(grib=='grib2') then
4306 cfld=cfld+1
4307 fld_info(cfld)%ifld=iavblfld(iget(992))
4308 fld_info(cfld)%lvl=lvlsxml(1,iget(992))
4309!$omp parallel do private(i,j,ii,jj)
4310 do j=1,jend-jsta+1
4311 jj = jsta+j-1
4312 do i=1,iend-ista+1
4313 ii = ista+i-1
4314 datapd(i,j,cfld) = grid1(ii,jj)
4315 enddo
4316 enddo
4317 endif
4318 ENDIF !992
4319
4320
4321 IF (iget(763)>0) THEN
4322!$omp parallel do private(i,j)
4323! EGRID3 is Virtual LFC
4324 DO j=jsta,jend
4325 DO i=ista,iend
4326 grid1(i,j) = q1d(i,j)
4327 ENDDO
4328 ENDDO
4329 if(grib=='grib2') then
4330 cfld=cfld+1
4331 fld_info(cfld)%ifld=iavblfld(iget(763))
4332 fld_info(cfld)%lvl=lvlsxml(1,iget(763))
4333!$omp parallel do private(i,j,ii,jj)
4334 do j=1,jend-jsta+1
4335 jj = jsta+j-1
4336 do i=1,iend-ista+1
4337 ii = ista+i-1
4338 datapd(i,j,cfld) = grid1(ii,jj)
4339 enddo
4340 enddo
4341 endif
4342 ENDIF
4343
4344!Hail parameter
4345 IF (iget(993)>0) THEN
4346 grid1=spval
4347 DO j=jsta,jend
4348 DO i=ista,iend
4349 lapse=-((t700(i,j)-t500(i,j))/((z700(i,j)-z500(i,j))))
4350 ship=(mucape(i,j)*d1000*muq1d(i,j)*lapse*(t500(i,j)-k2c)*fshr(i,j))/hconst
4351 IF (mucape(i,j)<1300.)THEN
4352 ship=ship*(mucape(i,j)/1300.)
4353 ENDIF
4354 IF (lapse < 5.8)THEN
4355 ship=ship*(lapse/5.8)
4356 ENDIF
4357 IF (freezelvl(i,j) < 2400.)THEN
4358 ship=ship*(freezelvl(i,j)/2400.)
4359 ENDIF
4360 grid1(i,j)=ship
4361 ENDDO
4362 ENDDO
4363 if(grib=='grib2') then
4364 cfld=cfld+1
4365 fld_info(cfld)%ifld=iavblfld(iget(993))
4366 fld_info(cfld)%lvl=lvlsxml(1,iget(993))
4367! $omp parallel do private(i,j,ii,jj)
4368 do j=1,jend-jsta+1
4369 jj = jsta+j-1
4370 do i=1,iend-ista+1
4371 ii = ista+i-1
4372 datapd(i,j,cfld) = grid1(ii,jj)
4373 enddo
4374 enddo
4375 endif
4376 ENDIF
4377
4378 ENDIF !END RTMA BLOCK
4379
4380
4381! Critical Angle
4382
4383 IF (iget(957)>0) THEN
4384 grid1=spval
4385!$omp parallel do private(i,j)
4386 DO j=jsta,jend
4387 DO i=ista,iend
4388 IF(t1d(i,j) < spval ) grid1(i,j) = cangle(i,j)
4389 ! IF(EGRID1(I,J)<100. .OR. EGRID2(I,J)>-250.) THEN
4390 ! GRID1(I,J) = 0.
4391 ! ENDIF
4392 ENDDO
4393 ENDDO
4394 if(grib=='grib2') then
4395 cfld=cfld+1
4396 fld_info(cfld)%ifld=iavblfld(iget(957))
4397 fld_info(cfld)%lvl=lvlsxml(1,iget(957))
4398!$omp parallel do private(i,j,ii,jj)
4399 do j=1,jend-jsta+1
4400 jj = jsta+j-1
4401 do i=1,iend-ista+1
4402 ii = ista+i-1
4403 datapd(i,j,cfld) = grid1(ii,jj)
4404 enddo
4405 enddo
4406 endif
4407 ENDIF !957
4408
4409! Dendritic Layer Depth, -17C < T < -12C
4410
4411 IF (iget(955)>0) THEN
4412 grid1=spval
4413!$omp parallel do private(i,j)
4414 DO j=jsta,jend
4415 DO i=ista,iend
4416 IF(t1d(i,j) < spval ) grid1(i,j) = egrid7(i,j)
4417 ENDDO
4418 ENDDO
4419 CALL bound(grid1,d00,h99999)
4420 if(grib=='grib2') then
4421 cfld=cfld+1
4422 fld_info(cfld)%ifld=iavblfld(iget(955))
4423 fld_info(cfld)%lvl=lvlsxml(1,iget(955))
4424!$omp parallel do private(i,j,ii,jj)
4425 do j=1,jend-jsta+1
4426 jj = jsta+j-1
4427 do i=1,iend-ista+1
4428 ii = ista+i-1
4429 datapd(i,j,cfld) = grid1(ii,jj)
4430 enddo
4431 enddo
4432 endif
4433 ENDIF !955
4434
4435! Enhanced Stretching Potential
4436
4437 IF (iget(956)>0) THEN
4438 grid1=spval
4439!$omp parallel do private(i,j)
4440 DO j=jsta,jend
4441 DO i=ista,iend
4442 IF(t1d(i,j) < spval ) grid1(i,j) = egrid8(i,j)
4443 ENDDO
4444 ENDDO
4445 CALL bound(grid1,d00,h99999)
4446 if(grib=='grib2') then
4447 cfld=cfld+1
4448 fld_info(cfld)%ifld=iavblfld(iget(956))
4449 fld_info(cfld)%lvl=lvlsxml(1,iget(956))
4450!$omp parallel do private(i,j,ii,jj)
4451 do j=1,jend-jsta+1
4452 jj = jsta+j-1
4453 do i=1,iend-ista+1
4454 ii = ista+i-1
4455 datapd(i,j,cfld) = grid1(ii,jj)
4456 enddo
4457 enddo
4458 endif
4459 ENDIF !956
4460
4461! Downdraft CAPE
4462
4463! ITYPE = 1
4464! DO J=JSTA,JEND
4465! DO I=ISTA,IEND
4466! LB2(I,J) = (LVLBND(I,J,1) + LVLBND(I,J,2) + &
4467! LVLBND(I,J,3))/3
4468! P1D(I,J) = (PBND(I,J,1) + PBND(I,J,2) + PBND(I,J,3))/3
4469! T1D(I,J) = (TBND(I,J,1) + TBND(I,J,2) + TBND(I,J,3))/3
4470! Q1D(I,J) = (QBND(I,J,1) + QBND(I,J,2) + QBND(I,J,3))/3
4471! ENDDO
4472! ENDDO
4473
4474! DPBND = 400.E2
4475! CALL CALCAPE2(ITYPE,DPBND,P1D,T1D,Q1D,LB2, &
4476! EGRID1,EGRID2,EGRID3,EGRID4,EGRID5, &
4477! EGRID6,EGRID7,EGRID8)
4478
4479 IF (iget(954)>0) THEN
4480 grid1 = spval
4481!$omp parallel do private(i,j)
4482 DO j=jsta,jend
4483 DO i=ista,iend
4484 IF(t1d(i,j) < spval) grid1(i,j) = -egrid6(i,j)
4485 ENDDO
4486 ENDDO
4487 CALL bound(grid1,d00,h99999)
4488 if(grib=='grib2') then
4489 cfld=cfld+1
4490 fld_info(cfld)%ifld=iavblfld(iget(954))
4491 fld_info(cfld)%lvl=lvlsxml(1,iget(954))
4492!$omp parallel do private(i,j,ii,jj)
4493 do j=1,jend-jsta+1
4494 jj = jsta+j-1
4495 do i=1,iend-ista+1
4496 ii = ista+i-1
4497 datapd(i,j,cfld) = grid1(ii,jj)
4498 enddo
4499 enddo
4500 endif
4501
4502 ENDIF !954
4503
4504 if (allocated(ushr1)) deallocate(ushr1)
4505 if (allocated(vshr1)) deallocate(vshr1)
4506 if (allocated(ushr6)) deallocate(ushr6)
4507 if (allocated(vshr6)) deallocate(vshr6)
4508 if (allocated(ust)) deallocate(ust)
4509 if (allocated(vst)) deallocate(vst)
4510 if (allocated(heli)) deallocate(heli)
4511 if (allocated(llow)) deallocate(llow)
4512 if (allocated(lupp)) deallocate(lupp)
4513 if (allocated(cangle))deallocate(cangle)
4514 if (allocated(effust))deallocate(effust)
4515 if (allocated(effvst))deallocate(effvst)
4516 if (allocated(eshr)) deallocate(eshr)
4517 if (allocated(uvect)) deallocate(uvect)
4518 if (allocated(vvect)) deallocate(vvect)
4519 if (allocated(esrh)) deallocate(esrh)
4520 if (allocated(htsfc)) deallocate(htsfc)
4521 if (allocated(fshr)) deallocate(fshr)
4522 if (allocated(llow_zint)) deallocate(llow_zint)
4523 if (allocated(ieql_zint)) deallocate(ieql_zint)
4524 if (allocated(z_temp)) deallocate(z_temp)
4525 if (allocated(midcal)) deallocate(midcal)
4526 if (allocated(z_midcal)) deallocate(z_midcal)
4527 if (allocated(el_base)) deallocate(el_base)
4528 if (allocated(el_tops)) deallocate(el_tops)
4529
4530 ENDIF
4531
4532 if (allocated(pbnd)) deallocate(pbnd)
4533 if (allocated(tbnd)) deallocate(tbnd)
4534 if (allocated(qbnd)) deallocate(qbnd)
4535 if (allocated(ubnd)) deallocate(ubnd)
4536 if (allocated(vbnd)) deallocate(vbnd)
4537 if (allocated(rhbnd)) deallocate(rhbnd)
4538 if (allocated(wbnd)) deallocate(wbnd)
4539 if (allocated(lvlbnd)) deallocate(lvlbnd)
4540 if (allocated(lb2)) deallocate(lb2)
4541
4542!
4543!
4544! RELATIVE HUMIDITY WITH RESPECT TO PRECIPITABLE WATER
4545 IF (iget(749)>0) THEN
4546 CALL calrh_pw(grid1(ista:iend,jsta:jend))
4547 if(grib=='grib2') then
4548 cfld=cfld+1
4549 fld_info(cfld)%ifld=iavblfld(iget(749))
4550!$omp parallel do private(i,j,ii,jj)
4551 do j=1,jend-jsta+1
4552 jj = jsta+j-1
4553 do i=1,iend-ista+1
4554 ii = ista+i-1
4555 datapd(i,j,cfld) = grid1(ii,jj)
4556 enddo
4557 enddo
4558 endif
4559 ENDIF
4560
4561
4562! END OF ROUTINE.
4563!
4564 RETURN
4565 CONTAINS
4567 subroutine allocate_cape_arrays
4568 if(.not.allocated(omgbnd)) allocate(omgbnd(ista:iend,jsta:jend,nbnd))
4569 if(.not.allocated(pwtbnd)) allocate(pwtbnd(ista:iend,jsta:jend,nbnd))
4570 if(.not.allocated(qcnvbnd)) allocate(qcnvbnd(ista:iend,jsta:jend,nbnd))
4571 if(.not.allocated(lvlbnd)) allocate(lvlbnd(ista:iend,jsta:jend,nbnd))
4572 if(.not.allocated(lb2)) allocate(lb2(ista:iend,jsta:jend))
4573 end subroutine allocate_cape_arrays
4574
4575 END SUBROUTINE miscln
bndlyr
subroutine bndlyr(pbnd, tbnd, qbnd, rhbnd, ubnd, vbnd, wbnd, omgbnd, pwtbnd, qcnvbnd, lvlbnd)
Computes boundary layer fields.
Definition BNDLYR.f:57
bound
subroutine bound(fld, fmin, fmax)
Clips data in passed array.
Definition BOUND.f:37
caldwp
subroutine caldwp(p1d, q1d, tdwp, t1d)
Computes dewpoint from P, T, and Q.
Definition CALDWP.f:28
calhel2
subroutine calhel2(llow, lupp, depth, ust, vst, heli, cangle)
This routine computes estimated storm motion and storm-relative environmental helicity.
Definition CALHEL2.f:49
calhel3
subroutine calhel3(llow, lupp, ust, vst, heli)
This routine computes estimated storm motion and storm-relative environmental helicity.
Definition CALHEL3.f:50
calhel
subroutine calhel(depth, ust, vst, heli, ushr1, vshr1, ushr6, vshr6)
Computes storm relative helicity and storm motion.
Definition CALHEL.f:47
callcl
subroutine callcl(p1d, t1d, q1d, plcl, zlcl)
Subroutine that computes the lifting condensation level (LCL) height (above ground level) and pressur...
Definition CALLCL.f:42
calpot
subroutine calpot(p1d, t1d, theta)
Subroutine that computes potential temperature.
Definition CALPOT.f:29
calthte
subroutine calthte(p1d, t1d, q1d, thte)
Subroutine that computes Theta-E.
Definition CALTHTE.f:32
calupdhel
subroutine calupdhel(updhel)
Subroutine that computes the updraft helicity.
Definition CALUPDHEL.f:25
exch
subroutine exch(a)
exch() Subroutine that exchanges one halo row.
Definition EXCH.f:27
fdlvl
subroutine fdlvl(itype, tfd, qfd, ufd, vfd, pfd, icingfd)
fdlvl() Subroutine that computes T, Q, U, V, P, ICING on the flight levels (FD).
Definition FDLVL.f:60
fdlvl_mass
subroutine fdlvl_mass(itype, nfd, ptfd, htfd, nin, qin, qtype, qfd)
Computes FD level for mass variables.
Definition FDLVL.f:799
frzlvl2
subroutine frzlvl2(isotherm, zfrz, rhfrz, pfrzl)
FRZLVL2 computes FRZING LVL, Z and RH.
Definition FRZLVL2.f:51
frzlvl
subroutine frzlvl(zfrz, rhfrz, pfrzl)
FRZLVL() Subroutine that computes FRZING LVL, Z and RH.
Definition FRZLVL.f:49
mxwind
subroutine mxwind(km, p, u, v, t, h, pmw, umw, vmw, tmw, hmw)
mxwind() computes maximum wind level fields.
Definition GFSPOST.F:460
tpause
subroutine tpause(km, p, u, v, t, h, ptp, utp, vtp, ttp, htp, shrtp)
tpause() computes tropopause level fields.
Definition GFSPOST.F:373
lfmfld
subroutine lfmfld(rh3310, rh6610, rh3366, pw3310)
LFMFLD() computes layer mean LFM fields.
Definition LFMFLD.f:53
lfmfld_gfs
subroutine lfmfld_gfs(rh4410, rh7294, rh4472, rh3310)
LFMFLD_GFS() computes layer mean LFM fields.
Definition LFMFLD_GFS.f:55
miscln
subroutine miscln
MISCLN posts miscellaneous fields.
Definition MISCLN.f:58
ngmfld
subroutine ngmfld(rh4710, rh4796, rh1847, rh8498, qm8510)
ngmfld() computes layer mean NGM fields
Definition NGMFLD.f:60
otlft
subroutine otlft(pbnd, tbnd, qbnd, slindx)
Computes lifted index at 500mb.
Definition OTLFT.f:32
sclfld
subroutine sclfld(fld, scale, imo, jmo)
sclfld() scale array element by constant.
Definition SCLFLD.f:34