UPP (upp-srw-2.2.0)
Loading...
Searching...
No Matches
MDL2AGL.f
Go to the documentation of this file.
1
2!
49 SUBROUTINE mdl2agl
50
51!
52!
53 use vrbls3d, only: zmid, zint, dbz, dbzr, dbzi, dbzc, uh, vh, pmid, t, q, ref_10cm
54 use vrbls2d, only: refd_max, up_heli_max, up_heli_max16, grpl_max, &
55 up_heli_min, up_heli_min16, up_heli_max02, &
56 up_heli_min02, up_heli_max03, up_heli_min03, &
57 rel_vort_max, rel_vort_max01, hail_max2d, hail_maxk1,&
58 hail_maxhailcast,refdm10c_max, rel_vort_maxhy1, &
59 ltg1_max, ltg2_max, ltg3_max, up_heli, up_heli16, &
60 nci_ltg, nca_ltg, nci_wq, nca_wq, nci_refd, nca_refd,&
61 u10, v10, u10h, v10h
62 use masks, only: lmh, lmv
63 use params_mod, only: dbzmin, small, eps, rd
64 use ctlblk_mod, only: spval, lm, modelname, grib, cfld, fld_info, datapd,&
65 ifhr, global, jsta_m, jend_m, mpi_comm_comp, &
66 jsta_2l, jend_2u, im, jm, jsta, jend, imp_physics, &
67 ista, iend, ista_2l, iend_2u, ista_m, iend_m
68 use rqstfld_mod, only: iget, lvls, iavblfld, lvlsxml, id
69 use gridspec_mod, only: gridtype
70!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
71 implicit none
72 include "mpif.h"
73!
74! INCLUDE MODEL DIMENSIONS. SET/DERIVE OTHER PARAMETERS.
75! GAMMA AND RGAMOG ARE USED IN THE EXTRAPOLATION OF VIRTUAL
76! TEMPERATURES BEYOND THE UPPER OF LOWER LIMITS OF DATA.
77!
78 integer,PARAMETER :: LAGL=2,lagl2=1
79!
80! DECLARE VARIABLES.
81!
82 LOGICAL IOOMG,IOALL
83 REAL,dimension(ista_2l:iend_2u,jsta_2l:jend_2u) :: grid1
84 REAL,dimension(ista_2l:iend_2u,jsta_2l:jend_2u) :: UAGL, VAGL, tagl, pagl, qagl
85!
86 INTEGER,dimension(ista_2l:iend_2u,jsta_2l:jend_2u) :: NL1X
87 integer,dimension(jm) :: IHE, IHW
88 INTEGER LXXX,IERR, maxll, minll
89 INTEGER ISTART,ISTOP,JSTART,JSTOP
90!
91!
92!--- Definition of the following 2D (horizontal) dummy variables
93!
94! C1D - total condensate
95! QW1 - cloud water mixing ratio
96! QI1 - cloud ice mixing ratio
97! QR1 - rain mixing ratio
98! QS1 - snow mixing ratio
99! DBZ1 - radar reflectivity
100! DBZR1 - radar reflectivity from rain
101! DBZI1 - radar reflectivity from ice (snow + graupel + sleet)
102! DBZC1 - radar reflectivity from parameterized convection (bogused)
103!
104! REAL C1D(IM,JM),QW1(IM,JM),QI1(IM,JM),QR1(IM,JM)
105! &, QS1(IM,JM) ,DBZ1(IM,JM)
106 REAL,dimension(ista:iend,jsta:jend) :: DBZ1, DBZR1, DBZI1, DBZC1, dbz1log
107 real,dimension(lagl) :: ZAGL
108 real,dimension(lagl2) :: ZAGL2, ZAGL3
109 real PAGLU,PAGLL,TAGLU,TAGLL,QAGLU,QAGLL, pv, rho
110
111 integer I,J,L,II,JJ,LP,LL,LLMH,ie,iw,jn,js,iget1,iget2,iget3,iget4
112 real UAGLL,UAGLU,VAGLL,VAGLU,FACT,ZDUM
113!
114!
115!******************************************************************************
116!
117! START MDL2P.
118!
119! SET TOTAL NUMBER OF POINTS ON OUTPUT GRID.
120!
121!---------------------------------------------------------------
122
123
124 zagl(1) = 4000.
125 zagl(2) = 1000.
126 zagl2(1) = 609.6 ! 2000 ft
127! CRA
128 zagl3(1) = 80.
129! CRA
130
131!
132! *** PART I ***
133!
134! VERTICAL INTERPOLATION OF EVERYTHING ELSE. EXECUTE ONLY
135! IF THERE'S SOMETHING WE WANT.
136!
137 IF (iget(253)>0 .OR. iget(279)>0 .OR. iget(280)>0 .OR. &
138 & iget(281)>0 ) THEN
139!
140!---------------------------------------------------------------------
141!***
142!*** BECAUSE SIGMA LAYERS DO NOT GO UNDERGROUND, DO ALL
143!*** INTERPOLATION ABOVE GROUND NOW.
144!***
145!
146 DO 310 lp=1,lagl
147 iget1 = -1 ; iget2 = -1 ; iget3 = -1 ; iget4 = -1
148 if (iget(253) > 0) iget1 = lvls(lp,iget(253))
149 if (iget(279) > 0) iget2 = lvls(lp,iget(279))
150 if (iget(280) > 0) iget3 = lvls(lp,iget(280))
151 if (iget(281) > 0) iget4 = lvls(lp,iget(281))
152 IF (iget1 > 0 .or. iget2 > 0 .or. iget3 > 0 .or. iget4 > 0) then
153!
154 jj=float(jsta+jend)/2.0
155 ii=float(ista+iend)/3.0
156
157 DO j=jsta,jend
158 DO i=ista,iend
159 dbz1(i,j) = spval
160 dbzr1(i,j) = spval
161 dbzi1(i,j) = spval
162 dbzc1(i,j) = spval
163!
164!*** LOCATE VERTICAL INDEX OF MODEL MIDLAYER JUST BELOW
165!*** THE AGL LEVEL TO WHICH WE ARE INTERPOLATING.
166!
167 llmh = nint(lmh(i,j))
168 nl1x(i,j) = llmh+1
169 DO l=llmh,2,-1
170 zdum = zmid(i,j,l)-zint(i,j,llmh+1)
171 IF(zdum >= zagl(lp)) THEN
172 nl1x(i,j) = l+1
173 exit
174 ENDIF
175 ENDDO
176!
177! IF THE AGL LEVEL IS BELOW THE LOWEST MODEL MIDLAYER
178! BUT STILL ABOVE THE LOWEST MODEL BOTTOM INTERFACE,
179! WE WILL NOT CONSIDER IT UNDERGROUND AND THE INTERPOLATION
180! WILL EXTRAPOLATE TO THAT POINT
181!
182 IF(nl1x(i,j) == (llmh+1) .AND. zagl(lp) > 0.) THEN
183 nl1x(i,j) = lm
184 ENDIF
185!
186! if(NL1X(I,J)==LMP1)print*,'Debug: NL1X=LMP1 AT '
187! 1 ,i,j,lp
188 ENDDO
189 ENDDO
190!
191!mptest IF(NHOLD==0)GO TO 310
192!
193!!$omp parallel do
194!!$omp& private(nn,i,j,ll,fact,qsat,rhl)
195!hc DO 220 NN=1,NHOLD
196!hc I=IHOLD(NN)
197!hc J=JHOLD(NN)
198! DO 220 J=JSTA,JEND
199
200 DO j=jsta,jend
201 DO i=ista,iend
202 ll = nl1x(i,j)
203!---------------------------------------------------------------------
204!*** VERTICAL INTERPOLATION OF GEOPOTENTIAL, TEMPERATURE, SPECIFIC
205!*** HUMIDITY, CLOUD WATER/ICE, OMEGA, WINDS, AND TKE.
206!---------------------------------------------------------------------
207!
208!HC IF(NL1X(I,J)<=LM)THEN
209 llmh = nint(lmh(i,j))
210 IF(nl1x(i,j)<=llmh)THEN
211 IF(zmid(i,j,ll)<spval.and.zmid(i,j,ll-1)<spval)THEN
212!
213!---------------------------------------------------------------------
214! INTERPOLATE LINEARLY IN LOG(P)
215!*** EXTRAPOLATE ABOVE THE TOPMOST MIDLAYER OF THE MODEL
216!*** INTERPOLATION BETWEEN NORMAL LOWER AND UPPER BOUNDS
217!*** EXTRAPOLATE BELOW LOWEST MODEL MIDLAYER (BUT STILL ABOVE GROUND)
218!---------------------------------------------------------------------
219!
220
221! FACT=(ALSL(LP)-ALOG(PMID(I,J,LL)))/
222! & (ALOG(PMID(I,J,LL))-ALOG(PMID(I,J,LL-1)))
223 zdum=zagl(lp)+zint(i,j,nint(lmh(i,j))+1)
224 fact=(zdum-zmid(i,j,ll))/(zmid(i,j,ll)-zmid(i,j,ll-1))
225!
226! KRF: Use arw/nmm output if thompson
227 if (imp_physics==8) then
228 dbz1(i,j)=ref_10cm(i,j,ll)+(ref_10cm(i,j,ll)-ref_10cm(i,j,ll-1))*fact
229 else
230 dbz1(i,j)=dbz(i,j,ll)+(dbz(i,j,ll)-dbz(i,j,ll-1))*fact
231 end if
232 ! DBZ1(I,J) = DBZ(I,J,LL) + (DBZ(I,J,LL)-DBZ(I,J,LL-1))*FACT
233 dbzr1(i,j) = dbzr(i,j,ll) + (dbzr(i,j,ll)-dbzr(i,j,ll-1))*fact
234 dbzi1(i,j) = dbzi(i,j,ll) + (dbzi(i,j,ll)-dbzi(i,j,ll-1))*fact
235 dbzc1(i,j) = dbzc(i,j,ll) + (dbzc(i,j,ll)-dbzc(i,j,ll-1))*fact
236 if(modelname=='RAPR') then
237 if(dbz1(i,j)>0.) then
238 dbz1log(i,j)= 10.*log10(dbz1(i,j))
239 else
240 dbz1log(i,j)= -100.
241 endif
242 endif
243! IF(I==ii.and.j==jj)print*,'Debug AGL RADAR REF',
244! & i,j,ll,zagl(lp),ZINT(I,J,NINT(LMH(I,J))+1)
245! & ,ZMID(I,J,LL-1),ZMID(I,J,LL)
246! & ,DBZ(I,J,LL-1),DBZ(I,J,LL),DBZ1(I,J)
247! & ,DBZR(I,J,LL-1),DBZR(I,J,LL),DBZR1(I,J)
248! & ,DBZI(I,J,LL-1),DBZI(I,J,LL),DBZI1(I,J)
249! & ,DBZC(I,J,LL-1),DBZC(I,J,LL),DBZC1(I,J)
250 if(modelname=='RAPR') then
251 dbz1log(i,j)=max(dbz1log(i,j),dbzmin)
252 else
253 dbz1(i,j)=max(dbz1(i,j),dbzmin)
254 endif
255 dbzr1(i,j) = max(dbzr1(i,j),dbzmin)
256 dbzi1(i,j) = max(dbzi1(i,j),dbzmin)
257 dbzc1(i,j) = max(dbzc1(i,j),dbzmin)
258 ENDIF !end ZMID(I,J,LL)<spval
259!
260! FOR UNDERGROUND AGL LEVELS, ASSUME TEMPERATURE TO CHANGE
261! ADIABATICLY, RH TO BE THE SAME AS THE AVERAGE OF THE 2ND AND 3RD
262! LAYERS FROM THE GOUND, WIND TO BE THE SAME AS THE LOWEST LEVEL ABOVE
263! GOUND
264 ELSE
265 dbz1log(i,j) = dbzmin
266 dbzr1(i,j) = dbzmin
267 dbzi1(i,j) = dbzmin
268 dbzc1(i,j) = dbzmin
269 END IF
270 enddo
271 enddo
272!
273!
274!---------------------------------------------------------------------
275! *** PART II ***
276!---------------------------------------------------------------------
277!
278! OUTPUT SELECTED FIELDS.
279!
280!---------------------------------------------------------------------
281!
282!
283!--- Radar Reflectivity
284 IF((iget(253)>0) )THEN
285 if(modelname=='RAPR') then
286 DO j=jsta,jend
287 DO i=ista,iend
288 grid1(i,j)=dbz1log(i,j)
289 ENDDO
290 ENDDO
291 else
292 DO j=jsta,jend
293 DO i=ista,iend
294 grid1(i,j)=dbz1(i,j)
295 ENDDO
296 ENDDO
297 endif
298 if(grib=='grib2') then
299 cfld=cfld+1
300 fld_info(cfld)%ifld=iavblfld(iget(253))
301 fld_info(cfld)%lvl=lvlsxml(lp,iget(253))
302 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
303 endif
304 END IF
305!--- Radar reflectivity from rain
306 IF((iget(279)>0) )THEN
307 DO j=jsta,jend
308 DO i=ista,iend
309 grid1(i,j)=dbzr1(i,j)
310 ENDDO
311 ENDDO
312 if(grib=='grib2') then
313 cfld=cfld+1
314 fld_info(cfld)%ifld=iavblfld(iget(279))
315 fld_info(cfld)%lvl=lvlsxml(lp,iget(279))
316 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
317 endif
318 END IF
319!--- Radar reflectivity from all ice habits (snow + graupel + sleet, etc.)
320 IF((iget(280)>0) )THEN
321 DO j=jsta,jend
322 DO i=ista,iend
323 grid1(i,j)=dbzi1(i,j)
324 ENDDO
325 ENDDO
326 if(grib=='grib2') then
327 cfld=cfld+1
328 fld_info(cfld)%ifld=iavblfld(iget(280))
329 fld_info(cfld)%lvl=lvlsxml(lp,iget(280))
330 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
331 endif
332 END IF
333!--- Radar reflectivity from parameterized convection
334 IF((iget(281)>0) )THEN
335 DO j=jsta,jend
336 DO i=ista,iend
337 grid1(i,j)=dbzc1(i,j)
338 ENDDO
339 ENDDO
340 if(grib=='grib2') then
341 cfld=cfld+1
342 fld_info(cfld)%ifld=iavblfld(iget(281))
343 fld_info(cfld)%lvl=lvlsxml(lp,iget(281))
344 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
345 endif
346 END IF
347!
348 ENDIF ! FOR LEVEL
349!
350!*** END OF MAIN VERTICAL LOOP
351!
352 310 CONTINUE
353!*** ENDIF FOR IF TEST SEEING IF WE WANT ANY OTHER VARIABLES
354!
355 ENDIF
356! SRD
357 lp=1
358!--- Max Derived Radar Reflectivity
359 IF((iget(421)>0) )THEN
360 DO j=jsta,jend
361 DO i=ista,iend
362 grid1(i,j)=refd_max(i,j)
363 ENDDO
364 ENDDO
365 if(grib=='grib2') then
366 cfld=cfld+1
367 fld_info(cfld)%ifld=iavblfld(iget(421))
368 fld_info(cfld)%lvl=lvlsxml(lp,iget(421))
369 fld_info(cfld)%tinvstat=1
370 if (ifhr > 0) then
371 fld_info(cfld)%tinvstat=1
372 else
373 fld_info(cfld)%tinvstat=0
374 endif
375 fld_info(cfld)%ntrange=1
376 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
377 endif
378 END IF
379
380!--- Max Derived Radar Reflectivity at -10C
381 IF((iget(785)>0) )THEN
382 DO j=jsta,jend
383 DO i=ista,iend
384 grid1(i,j)=refdm10c_max(i,j)
385 ENDDO
386 ENDDO
387 if(grib=='grib2') then
388 cfld=cfld+1
389 fld_info(cfld)%ifld=iavblfld(iget(785))
390 fld_info(cfld)%lvl=lvlsxml(lp,iget(785))
391 if (ifhr > 0) then
392 fld_info(cfld)%tinvstat=1
393 else
394 fld_info(cfld)%tinvstat=0
395 endif
396 fld_info(cfld)%ntrange=1
397 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
398 endif
399 END IF
400
401!--- Max Updraft Helicity
402 IF((iget(420)>0) )THEN
403 DO j=jsta,jend
404 DO i=ista,iend
405 grid1(i,j)=up_heli_max(i,j)
406 ENDDO
407 ENDDO
408 if(grib=='grib2') then
409 cfld=cfld+1
410 fld_info(cfld)%ifld=iavblfld(iget(420))
411 fld_info(cfld)%lvl=lvlsxml(lp,iget(420))
412 if (ifhr > 0) then
413 fld_info(cfld)%tinvstat = 1
414 else
415 fld_info(cfld)%tinvstat = 0
416 endif
417 fld_info(cfld)%ntrange = 1
418 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
419 endif
420 END IF
421
422!--- Max Updraft Helicity 1-6 km
423 IF((iget(700)>0) )THEN
424 DO j=jsta,jend
425 DO i=ista,iend
426 grid1(i,j)=up_heli_max16(i,j)
427 ENDDO
428 ENDDO
429 if(grib=='grib2') then
430 cfld=cfld+1
431 fld_info(cfld)%ifld=iavblfld(iget(700))
432 fld_info(cfld)%lvl=lvlsxml(lp,iget(700))
433 if (ifhr == 0) then
434 fld_info(cfld)%tinvstat = 0
435 else
436 fld_info(cfld)%tinvstat = 1
437 endif
438 fld_info(cfld)%ntrange = 1
439 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
440 endif
441 END IF
442
443!--- Min Updraft Helicity
444 IF((iget(786)>0) )THEN
445 DO j=jsta,jend
446 DO i=ista,iend
447 grid1(i,j)=up_heli_min(i,j)
448 ENDDO
449 ENDDO
450 if(grib=='grib2') then
451 cfld=cfld+1
452 fld_info(cfld)%ifld=iavblfld(iget(786))
453 fld_info(cfld)%lvl=lvlsxml(lp,iget(786))
454 if (ifhr > 0) then
455 fld_info(cfld)%tinvstat = 1
456 else
457 fld_info(cfld)%tinvstat = 0
458 endif
459 fld_info(cfld)%ntrange = 1
460 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
461 endif
462 END IF
463
464!--- Min Updraft Helicity 1-6 km
465 IF((iget(787)>0) )THEN
466 DO j=jsta,jend
467 DO i=ista,iend
468 grid1(i,j)=up_heli_min16(i,j)
469 ENDDO
470 ENDDO
471 if(grib=='grib2') then
472 cfld=cfld+1
473 fld_info(cfld)%ifld=iavblfld(iget(787))
474 fld_info(cfld)%lvl=lvlsxml(lp,iget(787))
475 if (ifhr == 0) then
476 fld_info(cfld)%tinvstat = 0
477 else
478 fld_info(cfld)%tinvstat = 1
479 endif
480 fld_info(cfld)%ntrange = 1
481 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
482 endif
483 END IF
484
485!--- Max Updraft Helicity 0-2 km
486 IF((iget(788)>0) )THEN
487 DO j=jsta,jend
488 DO i=ista,iend
489 grid1(i,j)=up_heli_max02(i,j)
490 ENDDO
491 ENDDO
492 if(grib=='grib2') then
493 cfld=cfld+1
494 fld_info(cfld)%ifld=iavblfld(iget(788))
495 fld_info(cfld)%lvl=lvlsxml(lp,iget(788))
496 if (ifhr > 0) then
497 fld_info(cfld)%tinvstat = 1
498 else
499 fld_info(cfld)%tinvstat = 0
500 endif
501 fld_info(cfld)%ntrange = 1
502 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
503 endif
504 END IF
505!--- Min Updraft Helicity 0-2 km
506 IF((iget(789)>0) )THEN
507 DO j=jsta,jend
508 DO i=ista,iend
509 grid1(i,j)=up_heli_min02(i,j)
510 ENDDO
511 ENDDO
512 if(grib=='grib2') then
513 cfld=cfld+1
514 fld_info(cfld)%ifld=iavblfld(iget(789))
515 fld_info(cfld)%lvl=lvlsxml(lp,iget(789))
516 if (ifhr == 0) then
517 fld_info(cfld)%tinvstat = 0
518 else
519 fld_info(cfld)%tinvstat = 1
520 endif
521 fld_info(cfld)%ntrange = 1
522 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
523 endif
524 END IF
525
526!--- Max Updraft Helicity 0-3 km
527 IF((iget(790)>0) )THEN
528 DO j=jsta,jend
529 DO i=ista,iend
530 grid1(i,j)=up_heli_max03(i,j)
531 ENDDO
532 ENDDO
533 if(grib=='grib2') then
534 cfld=cfld+1
535 fld_info(cfld)%ifld=iavblfld(iget(790))
536 fld_info(cfld)%lvl=lvlsxml(lp,iget(790))
537 if (ifhr > 0) then
538 fld_info(cfld)%tinvstat = 1
539 else
540 fld_info(cfld)%tinvstat = 0
541 endif
542 fld_info(cfld)%ntrange = 1
543 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
544 endif
545 END IF
546
547!--- Min Updraft Helicity 0-3 km
548 IF((iget(791)>0) )THEN
549 DO j=jsta,jend
550 DO i=ista,iend
551 grid1(i,j)=up_heli_min03(i,j)
552 ENDDO
553 ENDDO
554 if(grib=='grib2') then
555 cfld=cfld+1
556 fld_info(cfld)%ifld=iavblfld(iget(791))
557 fld_info(cfld)%lvl=lvlsxml(lp,iget(791))
558 if (ifhr == 0) then
559 fld_info(cfld)%tinvstat = 0
560 else
561 fld_info(cfld)%tinvstat = 1
562 endif
563 fld_info(cfld)%ntrange = 1
564 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
565 endif
566 END IF
567
568!--- Max Relative Vertical Vorticity 0-2 km
569 IF((iget(792)>0) )THEN
570 DO j=jsta,jend
571 DO i=ista,iend
572 grid1(i,j)=rel_vort_max(i,j)
573 ENDDO
574 ENDDO
575 if(grib=='grib2') then
576 cfld=cfld+1
577 fld_info(cfld)%ifld=iavblfld(iget(792))
578 fld_info(cfld)%lvl=lvlsxml(lp,iget(792))
579 if (ifhr > 0) then
580 fld_info(cfld)%tinvstat = 1
581 else
582 fld_info(cfld)%tinvstat = 0
583 endif
584 fld_info(cfld)%ntrange = 1
585 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
586 endif
587 END IF
588
589!--- Max Relative Vertical Vorticity 0-1 km
590 IF((iget(793)>0) )THEN
591 DO j=jsta,jend
592 DO i=ista,iend
593 grid1(i,j)=rel_vort_max01(i,j)
594 ENDDO
595 ENDDO
596 if(grib=='grib2') then
597 cfld=cfld+1
598 fld_info(cfld)%ifld=iavblfld(iget(793))
599 fld_info(cfld)%lvl=lvlsxml(lp,iget(793))
600 if (ifhr > 0) then
601 fld_info(cfld)%tinvstat = 1
602 else
603 fld_info(cfld)%tinvstat = 0
604 endif
605 fld_info(cfld)%ntrange = 1
606 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
607 endif
608 END IF
609!--- Max Relative Vertical Vorticity @ hybrid level 1
610 IF((iget(890)>0) )THEN
611 DO j=jsta,jend
612 DO i=ista,iend
613 grid1(i,j)=rel_vort_maxhy1(i,j)
614 ENDDO
615 ENDDO
616 if(grib=='grib2') then
617 cfld=cfld+1
618 fld_info(cfld)%ifld=iavblfld(iget(890))
619 fld_info(cfld)%lvl=lvlsxml(lp,iget(890))
620 if (ifhr > 0) then
621 fld_info(cfld)%tinvstat = 1
622 else
623 fld_info(cfld)%tinvstat = 0
624 endif
625 fld_info(cfld)%ntrange = 1
626 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
627 endif
628 END IF
629
630!--- Max Hail Diameter in Column
631 IF((iget(794)>0) )THEN
632 DO j=jsta,jend
633 DO i=ista,iend
634 grid1(i,j)=hail_max2d(i,j)
635 ENDDO
636 ENDDO
637 if(grib=='grib2') then
638 cfld=cfld+1
639 fld_info(cfld)%ifld=iavblfld(iget(794))
640 fld_info(cfld)%lvl=lvlsxml(lp,iget(794))
641 if (ifhr == 0) then
642 fld_info(cfld)%tinvstat = 0
643 else
644 fld_info(cfld)%tinvstat = 1
645 endif
646 fld_info(cfld)%ntrange = 1
647 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
648 endif
649 END IF
650
651!--- Max Hail Diameter at k=1
652 IF((iget(795)>0) )THEN
653 DO j=jsta,jend
654 DO i=ista,iend
655 grid1(i,j)=hail_maxk1(i,j)
656 ENDDO
657 ENDDO
658 if(grib=='grib2') then
659 cfld=cfld+1
660 fld_info(cfld)%ifld=iavblfld(iget(795))
661 fld_info(cfld)%lvl=lvlsxml(lp,iget(795))
662 if (ifhr == 0) then
663 fld_info(cfld)%tinvstat = 0
664 else
665 fld_info(cfld)%tinvstat = 1
666 endif
667 fld_info(cfld)%ntrange = 1
668 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
669 endif
670 END IF
671
672!--- Max hail diameter at surface from WRF HAILCAST algorithm (HRRR
673!applications)
674! (J. Kenyon/GSD, added 1 May 2019)
675 IF((iget(728)>0) )THEN
676 DO j=jsta,jend
677 DO i=ista,iend
678 grid1(i,j)=hail_maxhailcast(i,j)/1000.0 ! convert mm to m
679 ENDDO
680 ENDDO
681 if(grib=='grib2') then
682 cfld=cfld+1
683 fld_info(cfld)%ifld=iavblfld(iget(728))
684 fld_info(cfld)%lvl=lvlsxml(lp,iget(728))
685 if (ifhr == 0) then
686 fld_info(cfld)%tinvstat = 0
687 else
688 fld_info(cfld)%tinvstat = 1
689 endif
690 fld_info(cfld)%ntrange = 1
691 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
692 endif
693 END IF
694
695!--- Max Column Integrated Graupel
696 IF((iget(429)>0) )THEN
697 DO j=jsta,jend
698 DO i=ista,iend
699 grid1(i,j)=grpl_max(i,j)
700 ENDDO
701 ENDDO
702 if(grib=='grib2') then
703 cfld=cfld+1
704 fld_info(cfld)%ifld=iavblfld(iget(429))
705 fld_info(cfld)%lvl=lvlsxml(lp,iget(429))
706 if (ifhr == 0) then
707 fld_info(cfld)%tinvstat = 0
708 else
709 fld_info(cfld)%tinvstat = 1
710 endif
711 fld_info(cfld)%ntrange = 1
712 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
713 endif
714 END IF
715
716!--- Max Lightning Threat 1
717 IF((iget(702)>0) )THEN
718 DO j=jsta,jend
719 DO i=ista,iend
720 grid1(i,j)=ltg1_max(i,j)
721 ENDDO
722 ENDDO
723 if(grib=='grib2') then
724 cfld=cfld+1
725 fld_info(cfld)%ifld=iavblfld(iget(702))
726 fld_info(cfld)%lvl=lvlsxml(lp,iget(702))
727 if (ifhr == 0) then
728 fld_info(cfld)%tinvstat = 0
729 else
730 fld_info(cfld)%tinvstat = 1
731 endif
732 fld_info(cfld)%ntrange = 1
733 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
734 endif
735 END IF
736
737!--- Max Lightning Threat 2
738 IF((iget(703)>0) )THEN
739 DO j=jsta,jend
740 DO i=ista,iend
741 grid1(i,j)=ltg2_max(i,j)
742 ENDDO
743 ENDDO
744 if(grib=='grib2') then
745 cfld=cfld+1
746 fld_info(cfld)%ifld=iavblfld(iget(703))
747 fld_info(cfld)%lvl=lvlsxml(lp,iget(703))
748 if (ifhr == 0) then
749 fld_info(cfld)%tinvstat = 0
750 else
751 fld_info(cfld)%tinvstat = 1
752 endif
753 fld_info(cfld)%ntrange = 1
754 datapd(1:iend-ista+1,1:jend-jsta+1,cfld) = grid1(ista:iend,jsta:jend)
755 endif
756 END IF
757
758!--- Max Lightning Threat 3
759 IF((iget(704)>0) )THEN
760 DO j=jsta,jend
761 DO i=ista,iend
762 grid1(i,j)=ltg3_max(i,j)
763 ENDDO
764 ENDDO
765 if(grib=='grib2') then
766 cfld=cfld+1
767 fld_info(cfld)%ifld=iavblfld(iget(704))
768 fld_info(cfld)%lvl=lvlsxml(lp,iget(704))
769 if (ifhr == 0) then
770 fld_info(cfld)%tinvstat = 0
771 else
772 fld_info(cfld)%tinvstat = 1
773 endif
774 fld_info(cfld)%ntrange = 1
775 do j=1,jend-jsta+1
776 jj = jsta+j-1
777 do i=1,iend-ista+1
778 ii = ista+i-1
779 datapd(i,j,cfld) = grid1(ii,jj)
780 enddo
781 enddo
782 endif
783 !datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=GRID1(ista:iend,jsta:jend)
784 END IF
785
786!--- GSD Updraft Helicity
787 IF((iget(727)>0) )THEN
788 DO j=jsta,jend
789 DO i=ista,iend
790 grid1(i,j)=up_heli(i,j)
791 ENDDO
792 ENDDO
793 if(grib=='grib2') then
794 cfld=cfld+1
795 fld_info(cfld)%ifld=iavblfld(iget(727))
796 fld_info(cfld)%lvl=lvlsxml(lp,iget(727))
797 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
798 endif
799 END IF
800
801!--- Updraft Helicity 1-6 km layer
802 IF((iget(701)>0) )THEN
803 DO j=jsta,jend
804 DO i=ista,iend
805 grid1(i,j)=up_heli16(i,j)
806 ENDDO
807 ENDDO
808 if(grib=='grib2') then
809 cfld=cfld+1
810 fld_info(cfld)%ifld=iavblfld(iget(701))
811 fld_info(cfld)%lvl=lvlsxml(lp,iget(701))
812 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
813 endif
814 END IF
815
816!--- Convective Initiation Lightning
817 IF((iget(705)>0) )THEN
818 DO j=jsta,jend
819 DO i=ista,iend
820 grid1(i,j)=nci_ltg(i,j)/60.0
821 ENDDO
822 ENDDO
823 if(grib=='grib2') then
824 cfld=cfld+1
825 fld_info(cfld)%ifld=iavblfld(iget(705))
826 fld_info(cfld)%lvl=lvlsxml(lp,iget(705))
827 if (ifhr == 0) then
828 fld_info(cfld)%tinvstat = 0
829 else
830 fld_info(cfld)%tinvstat = 1
831 endif
832 fld_info(cfld)%ntrange = 1
833 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
834 endif
835 END IF
836
837!--- Convective Activity Lightning
838 IF((iget(706)>0) )THEN
839 DO j=jsta,jend
840 DO i=ista,iend
841 grid1(i,j)=nca_ltg(i,j)/60.0
842 ENDDO
843 ENDDO
844 if(grib=='grib2') then
845 cfld=cfld+1
846 fld_info(cfld)%ifld=iavblfld(iget(706))
847 fld_info(cfld)%lvl=lvlsxml(lp,iget(706))
848 if (ifhr == 0) then
849 fld_info(cfld)%tinvstat = 0
850 else
851 fld_info(cfld)%tinvstat = 1
852 endif
853 fld_info(cfld)%ntrange = 1
854 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
855 endif
856 END IF
857
858!--- Convective Initiation Vertical Hydrometeor Flux
859 IF((iget(707)>0) )THEN
860 DO j=jsta,jend
861 DO i=ista,iend
862 grid1(i,j)=nci_wq(i,j)/60.0
863 ENDDO
864 ENDDO
865 if(grib=='grib2') then
866 cfld=cfld+1
867 fld_info(cfld)%ifld=iavblfld(iget(707))
868 fld_info(cfld)%lvl=lvlsxml(lp,iget(707))
869 if (ifhr == 0) then
870 fld_info(cfld)%tinvstat = 0
871 else
872 fld_info(cfld)%tinvstat = 1
873 endif
874 fld_info(cfld)%ntrange = 1
875 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
876 endif
877 END IF
878
879!--- Convective Activity Vertical Hydrometeor Flux
880 IF((iget(708)>0) )THEN
881 DO j=jsta,jend
882 DO i=ista,iend
883 grid1(i,j)=nca_wq(i,j)/60.0
884 ENDDO
885 ENDDO
886 if(grib=='grib2') then
887 cfld=cfld+1
888 fld_info(cfld)%ifld=iavblfld(iget(708))
889 fld_info(cfld)%lvl=lvlsxml(lp,iget(708))
890 if (ifhr == 0) then
891 fld_info(cfld)%tinvstat = 0
892 else
893 fld_info(cfld)%tinvstat = 1
894 endif
895 fld_info(cfld)%ntrange = 1
896 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
897 endif
898 END IF
899
900!--- Convective Initiation Reflectivity
901 IF((iget(709)>0) )THEN
902 DO j=jsta,jend
903 DO i=ista,iend
904 grid1(i,j)=nci_refd(i,j)/60.0
905 ENDDO
906 ENDDO
907 if(grib=='grib2') then
908 cfld=cfld+1
909 fld_info(cfld)%ifld=iavblfld(iget(709))
910 fld_info(cfld)%lvl=lvlsxml(lp,iget(709))
911 if (ifhr == 0) then
912 fld_info(cfld)%tinvstat = 0
913 else
914 fld_info(cfld)%tinvstat = 1
915 endif
916 fld_info(cfld)%ntrange = 1
917 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
918 endif
919 END IF
920
921!--- Convective Activity Reflectivity
922 IF((iget(710)>0) )THEN
923 DO j=jsta,jend
924 DO i=ista,iend
925 grid1(i,j)=nca_refd(i,j)/60.0
926 ENDDO
927 ENDDO
928 if(grib=='grib2') then
929 cfld=cfld+1
930 fld_info(cfld)%ifld=iavblfld(iget(710))
931 fld_info(cfld)%lvl=lvlsxml(lp,iget(710))
932 if (ifhr == 0) then
933 fld_info(cfld)%tinvstat = 0
934 else
935 fld_info(cfld)%tinvstat = 1
936 endif
937 fld_info(cfld)%ntrange = 1
938 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
939 endif
940 END IF
941!
942! SRD
943
944!
945 IF((iget(259)>0) )THEN
946!
947!---------------------------------------------------------------------
948!***
949!*** BECAUSE SIGMA LAYERS DO NOT GO UNDERGROUND, DO ALL
950!*** INTERPOLATION ABOVE GROUND NOW.
951!***
952!
953 iget2 = -1
954 if (iget(253) > 0 ) iget2 = iavblfld(iget(253))
955 iget2 = iget(253)
956 DO 320 lp=1,lagl2
957 iget1 = -1
958 if (iget(259) > 0 ) iget1 = lvls(lp,iget(259))
959 IF(iget1 > 0 .or. iget2 > 0) THEN
960!
961 jj=(jsta+jend)/2
962 ii=(ista+iend)/2
963 DO j=jsta,jend
964 DO i=ista,iend
965 uagl(i,j) = spval
966 vagl(i,j) = spval
967!
968!*** LOCATE VERTICAL INDEX OF MODEL MIDLAYER JUST BELOW
969!*** THE AGL LEVEL TO WHICH WE ARE INTERPOLATING.
970!
971 llmh=nint(lmh(i,j))
972 nl1x(i,j) = llmh+1
973 DO l=llmh,2,-1
974 zdum=zmid(i,j,l)-zint(i,j,llmh+1)
975 IF(zdum >= zagl2(lp))THEN
976 nl1x(i,j)=l+1
977 exit
978 ENDIF
979 ENDDO
980!
981! IF THE AGL LEVEL IS BELOW THE LOWEST MODEL MIDLAYER
982! BUT STILL ABOVE THE LOWEST MODEL BOTTOM INTERFACE,
983! WE WILL NOT CONSIDER IT UNDERGROUND AND THE INTERPOLATION
984! WILL EXTRAPOLATE TO THAT POINT
985!
986 IF(nl1x(i,j) == (llmh+1) .AND. zagl2(lp) > 0.) THEN
987 nl1x(i,j)=lm
988 ENDIF
989!
990! if(NL1X(I,J)==LMP1)print*,'Debug: NL1X=LMP1 AT '
991! 1 ,i,j,lp
992 ENDDO
993 ENDDO
994!
995!mptest IF(NHOLD==0)GO TO 310
996!
997!!$omp parallel do
998!!$omp& private(nn,i,j,ll,fact,qsat,rhl)
999!hc DO 220 NN=1,NHOLD
1000!hc I=IHOLD(NN)
1001!hc J=JHOLD(NN)
1002! DO 220 J=JSTA,JEND
1003 DO j=jsta,jend
1004 IF(gridtype=='A')THEN
1005 ihw(j)=-1
1006 ihe(j)=1
1007 ELSE IF(gridtype=='E')THEN
1008 ihw(j)=-mod(j,2)
1009 ihe(j)=ihw(j)+1
1010 END IF
1011 ENDDO
1012 IF(global)then
1013 istart=ista
1014 istop=iend
1015 jstart=jsta
1016 jstop=jend
1017 ELSE
1018 istart=ista_m
1019 istop=iend_m
1020 jstart=jsta_m
1021 jstop=jend_m
1022 END IF
1023 IF(gridtype/='A')THEN
1024! MAXLL=maxval(NL1X)
1025 minll=minval(nl1x)
1026! print*,'MINLL before all reduce= ',MINLL
1027 CALL mpi_allreduce(minll,lxxx,1,mpi_integer,mpi_min,mpi_comm_comp,ierr)
1028 minll=lxxx
1029! print*,'exchange wind in MDL2AGL from ',MINLL
1030 DO ll=minll,lm
1031 call exch(uh(ista_2l:iend_2u,jsta_2l:jend_2u,ll))
1032 call exch(vh(ista_2l:iend_2u,jsta_2l:jend_2u,ll))
1033 END DO
1034 END IF
1035 DO 230 j=jstart,jstop
1036 DO 230 i=istart,istop
1037 ll=nl1x(i,j)
1038!---------------------------------------------------------------------
1039!*** VERTICAL INTERPOLATION OF GEOPOTENTIAL, TEMPERATURE, SPECIFIC
1040!*** HUMIDITY, CLOUD WATER/ICE, OMEGA, WINDS, AND TKE.
1041!---------------------------------------------------------------------
1042!
1043!HC IF(NL1X(I,J)<=LM)THEN
1044 llmh = nint(lmh(i,j))
1045 IF(nl1x(i,j)<=llmh)THEN
1046!
1047!---------------------------------------------------------------------
1048! INTERPOLATE LINEARLY IN LOG(P)
1049!*** EXTRAPOLATE ABOVE THE TOPMOST MIDLAYER OF THE MODEL
1050!*** INTERPOLATION BETWEEN NORMAL LOWER AND UPPER BOUNDS
1051!*** EXTRAPOLATE BELOW LOWEST MODEL MIDLAYER (BUT STILL ABOVE GROUND)
1052!---------------------------------------------------------------------
1053!
1054! FACT=(ALSL(LP)-ALOG(PMID(I,J,LL)))/
1055! & (ALOG(PMID(I,J,LL))-ALOG(PMID(I,J,LL-1)))
1056 zdum=zagl2(lp)+zint(i,j,nint(lmh(i,j))+1)
1057 fact=(zdum-zmid(i,j,ll))/(zmid(i,j,ll)-zmid(i,j,ll-1))
1058!
1059 IF(gridtype=='A')THEN
1060 uaglu=uh(i,j,ll-1)
1061 uagll=uh(i,j,ll)
1062
1063 vaglu=vh(i,j,ll-1)
1064 vagll=vh(i,j,ll)
1065 ELSE IF(gridtype=='E')THEN
1066 uaglu=(uh(i+ihe(j),j,ll-1)+uh(i+ihw(j),j,ll-1)+ &
1067 & uh(i,j-1,ll-1)+uh(i,j+1,ll-1))/4.0
1068 uagll=(uh(i+ihe(j),j,ll)+uh(i+ihw(j),j,ll)+ &
1069 & uh(i,j-1,ll)+uh(i,j+1,ll))/4.0
1070
1071 vaglu=(vh(i+ihe(j),j,ll-1)+vh(i+ihw(j),j,ll-1)+ &
1072 & vh(i,j-1,ll-1)+vh(i,j+1,ll-1))/4.0
1073 vagll=(vh(i+ihe(j),j,ll)+vh(i+ihw(j),j,ll)+ &
1074 & vh(i,j-1,ll)+vh(i,j+1,ll))/4.0
1075 ELSE IF(gridtype=='B')THEN
1076 ie=i
1077 iw=i-1
1078 jn=j
1079 js=j-1
1080 uaglu=(uh(ie,j,ll-1)+uh(iw,j,ll-1)+ &
1081 & uh(ie,js,ll-1)+uh(iw,js,ll-1))/4.0
1082 uagll=(uh(ie,j,ll)+uh(iw,j,ll)+ &
1083 & uh(ie,js,ll)+uh(iw,js,ll))/4.0
1084
1085 vaglu=(vh(ie,j,ll-1)+vh(iw,j,ll-1)+ &
1086 & vh(ie,js,ll-1)+vh(iw,js,ll-1))/4.0
1087 vagll=(vh(ie,j,ll)+vh(iw,j,ll)+ &
1088 & vh(ie,js,ll)+vh(iw,js,ll))/4.0
1089 END IF
1090 uagl(i,j)=uagll+(uagll-uaglu)*fact
1091 vagl(i,j)=vagll+(vagll-vaglu)*fact
1092 IF(i==ii.AND.j==jj)print*, &
1093 & 'DEBUG LLWS: I,J,NL1X,UU,UL,VU,VL,ZSFC,ZMIDU,ZMIDL,U,V= ' &
1094 &, i,j,ll,uaglu,uagll,vaglu,vagll,zint(i,j,nint(lmh(i,j))+1)&
1095 &, zmid(i,j,ll-1),zmid(i,j,ll),uagl(i,j),vagl(i,j) &
1096 &, u10(i,j),v10(i,j)
1097!
1098! FOR UNDERGROUND AGL LEVELS, ASSUME TEMPERATURE TO CHANGE
1099! ADIABATICLY, RH TO BE THE SAME AS THE AVERAGE OF THE 2ND AND 3RD
1100! LAYERS FROM THE GOUND, WIND TO BE THE SAME AS THE LOWEST LEVEL ABOVE
1101! GOUND
1102 ELSE
1103 IF(gridtype=='A')THEN
1104 uagl(i,j)=uh(i,j,nint(lmv(i,j)))
1105 vagl(i,j)=vh(i,j,nint(lmv(i,j)))
1106 ELSE IF(gridtype=='E')THEN
1107 uagl(i,j)=(uh(i+ihe(j),j,nint(lmv(i+ihe(j),j))) &
1108 & +uh(i+ihw(j),j,nint(lmv(i+ihw(j),j)))+ &
1109 & uh(i,j-1,nint(lmv(i,j-1)))+uh(i,j+1,nint(lmv(i,j+1))))/4.0
1110 vagl(i,j)=(vh(i+ihe(j),j,nint(lmv(i+ihe(j),j))) &
1111 & +vh(i+ihw(j),j,nint(lmv(i+ihw(j),j)))+ &
1112 & vh(i,j-1,nint(lmv(i,j-1)))+vh(i,j+1,nint(lmv(i,j+1))))/4.0
1113 ELSE IF(gridtype=='B')THEN
1114 ie=i
1115 iw=i-1
1116 jn=j
1117 js=j-1
1118 uagl(i,j)=(uh(ie,j,nint(lmv(ie,j))) &
1119 & +uh(iw,j,nint(lmv(iw,j)))+ &
1120 & uh(ie,js,nint(lmv(ie,js)))+uh(iw,js,nint(lmv(iw,js))))/4.0
1121 vagl(i,j)=(vh(ie,j,nint(lmv(ie,j))) &
1122 & +vh(iw,j,nint(lmv(iw,j)))+ &
1123 & vh(ie,js,nint(lmv(ie,js)))+vh(iw,js,nint(lmv(iw,js))))/4.0
1124 END IF
1125 END IF
1126 230 CONTINUE
1127!
1128!
1129!---------------------------------------------------------------------
1130! *** PART II ***
1131!---------------------------------------------------------------------
1132!
1133! OUTPUT SELECTED FIELDS.
1134!
1135!---------------------------------------------------------------------
1136!
1137!
1138!--- Wind Shear (wind speed difference in knots between sfc and 2000 ft)
1139
1140 DO j=jsta,jend
1141 DO i=ista,iend
1142 IF(abs(uagl(i,j)-spval)>small .AND. &
1143 abs(vagl(i,j)-spval)>small)THEN
1144 IF(gridtype=='B' .OR. gridtype=='E')THEN
1145 grid1(i,j)=sqrt((uagl(i,j)-u10h(i,j))**2+ &
1146 (vagl(i,j)-v10h(i,j))**2)*1.943*zagl2(lp)/ &
1147 (zagl2(lp)-10.)
1148 ELSE
1149 grid1(i,j)=sqrt((uagl(i,j)-u10(i,j))**2+ &
1150 (vagl(i,j)-v10(i,j))**2)*1.943*zagl2(lp)/ &
1151 (zagl2(lp)-10.)
1152 END IF
1153 ELSE
1154 grid1(i,j)=spval
1155 END IF
1156 ENDDO
1157 ENDDO
1158 if(grib=="grib2" )then
1159 cfld=cfld+1
1160 fld_info(cfld)%ifld=iavblfld(iget(259))
1161 fld_info(cfld)%lvl=lvlsxml(lp,iget(259))
1162 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1163 endif
1164!
1165 ENDIF ! FOR LEVEL
1166!
1167!*** END OF MAIN VERTICAL LOOP
1168!
1169 320 CONTINUE
1170!*** ENDIF FOR IF TEST SEEING IF WE WANT ANY OTHER VARIABLES
1171!
1172 ENDIF
1173! CRA
1174 IF (iget(411)>0 .OR. iget(412)>0 .OR. iget(413)>0) THEN
1175!
1176!---------------------------------------------------------------------
1177!***
1178!*** BECAUSE SIGMA LAYERS DO NOT GO UNDERGROUND, DO ALL
1179!*** INTERPOLATION ABOVE GROUND NOW.
1180!***
1181!
1182 DO 330 lp=1,lagl2
1183 iget1 = -1 ; iget2 = -1 ; iget3 = -1
1184 if (iget(411) > 0) iget1 = lvls(lp,iget(411))
1185 if (iget(412) > 0) iget2 = lvls(lp,iget(412))
1186 if (iget(413) > 0) iget3 = lvls(lp,iget(413))
1187 IF (iget1 > 0 .or. iget2 > 0 .or. iget3 > 0) then
1188
1189!
1190 jj = float(jsta+jend)/2.0
1191 ii = float(ista+iend)/3.0
1192 DO j=jsta_2l,jend_2u
1193 DO i=ista_2l,iend_2u
1194!
1195 pagl(i,j) = spval
1196 tagl(i,j) = spval
1197 qagl(i,j) = spval
1198 uagl(i,j) = spval
1199 vagl(i,j) = spval
1200!
1201!*** LOCATE VERTICAL INDEX OF MODEL MIDLAYER JUST BELOW
1202!*** THE AGL LEVEL TO WHICH WE ARE INTERPOLATING.
1203!
1204 llmh = nint(lmh(i,j))
1205 nl1x(i,j) = llmh+1
1206 DO l=llmh,2,-1
1207 zdum = zmid(i,j,l)-zint(i,j,llmh+1)
1208 IF(zdum >= zagl3(lp))THEN
1209 nl1x(i,j) = l+1
1210 exit
1211 ENDIF
1212 ENDDO
1213!
1214! IF THE AGL LEVEL IS BELOW THE LOWEST MODEL MIDLAYER
1215! BUT STILL ABOVE THE LOWEST MODEL BOTTOM INTERFACE,
1216! WE WILL NOT CONSIDER IT UNDERGROUND AND THE INTERPOLATION
1217! WILL EXTRAPOLATE TO THAT POINT
1218!
1219 IF(nl1x(i,j)==(llmh+1) .AND. zagl3(lp)>0.)THEN
1220 nl1x(i,j) = lm
1221 ENDIF
1222!
1223! if(NL1X(I,J)==LMP1)print*,'Debug: NL1X=LMP1 AT '
1224! 1 ,i,j,lp
1225 ENDDO
1226 ENDDO
1227!
1228!mptest IF(NHOLD==0)GO TO 310
1229!
1230!!$omp parallel do
1231!!$omp& private(nn,i,j,ll,fact,qsat,rhl)
1232!chc DO 220 NN=1,NHOLD
1233!chc I=IHOLD(NN)
1234!chc J=JHOLD(NN)
1235! DO 220 J=JSTA,JEND
1236 DO 240 j=jsta_2l,jend_2u
1237 DO 240 i=ista_2l,iend_2u
1238 ll = nl1x(i,j)
1239!---------------------------------------------------------------------
1240!*** VERTICAL INTERPOLATION OF GEOPOTENTIAL, TEMPERATURE, SPECIFIC
1241!*** HUMIDITY, CLOUD WATER/ICE, OMEGA, WINDS, AND TKE.
1242!---------------------------------------------------------------------
1243!
1244!CHC IF(NL1X(I,J)<=LM)THEN
1245 llmh = nint(lmh(i,j))
1246 IF(nl1x(i,j)<=llmh)THEN
1247!
1248!---------------------------------------------------------------------
1249! INTERPOLATE LINEARLY IN LOG(P)
1250!*** EXTRAPOLATE ABOVE THE TOPMOST MIDLAYER OF THE MODEL
1251!*** INTERPOLATION BETWEEN NORMAL LOWER AND UPPER BOUNDS
1252!*** EXTRAPOLATE BELOW LOWEST MODEL MIDLAYER (BUT STILL ABOVE GROUND)
1253!---------------------------------------------------------------------
1254!
1255! FACT=(ALSL(LP)-ALOG(PMID(I,J,LL)))/
1256! & (ALOG(PMID(I,J,LL))-ALOG(PMID(I,J,LL-1)))
1257 zdum=zagl3(lp)+zint(i,j,nint(lmh(i,j))+1)
1258 fact = (zdum-zmid(i,j,ll)) &
1259 / (zmid(i,j,ll)-zmid(i,j,ll-1))
1260!
1261 paglu = log(pmid(i,j,ll-1))
1262 pagll = log(pmid(i,j,ll))
1263
1264 taglu = t(i,j,ll-1)
1265 tagll = t(i,j,ll)
1266
1267 qaglu = q(i,j,ll-1)
1268 qagll = q(i,j,ll)
1269
1270 uaglu = uh(i,j,ll-1)
1271 uagll = uh(i,j,ll)
1272
1273 vaglu = vh(i,j,ll-1)
1274 vagll = vh(i,j,ll)
1275
1276 pagl(i,j) = exp(pagll+(pagll-paglu)*fact)
1277 tagl(i,j) = tagll+(tagll-taglu)*fact
1278 qagl(i,j) = qagll+(qagll-qaglu)*fact
1279 uagl(i,j) = uagll+(uagll-uaglu)*fact
1280 vagl(i,j) = vagll+(vagll-vaglu)*fact
1281!
1282! FOR UNDERGROUND AGL LEVELS, ASSUME TEMPERATURE TO CHANGE
1283! ADIABATICLY, RH TO BE THE SAME AS THE AVERAGE OF THE 2ND AND 3RD
1284! LAYERS FROM THE GOUND, WIND TO BE THE SAME AS THE LOWEST LEVEL ABOVE
1285! GOUND
1286 ELSE
1287 pagl(i,j) = pmid(i,j,nint(lmv(i,j)))
1288 tagl(i,j) = t(i,j,nint(lmv(i,j)))
1289 qagl(i,j) = q(i,j,nint(lmv(i,j)))
1290 uagl(i,j) = uh(i,j,nint(lmv(i,j)))
1291 vagl(i,j) = vh(i,j,nint(lmv(i,j)))
1292 END IF
1293 240 CONTINUE
1294!
1295!
1296!---------------------------------------------------------------------
1297! *** PART II ***
1298!---------------------------------------------------------------------
1299!
1300! OUTPUT SELECTED FIELDS.
1301!
1302!---------------------------------------------------------------------
1303!
1304!
1305!--- Wind Energy Potential -- 0.5 * moist air density * wind speed^3
1306 IF((iget(411)>0) ) THEN
1307 DO j=jsta,jend
1308 DO i=ista,iend
1309 IF(qagl(i,j)<spval.and.pagl(i,j)<spval.and.tagl(i,j)<spval.and.&
1310 uagl(i,j)<spval.and.vagl(i,j)<spval)THEN
1311 qagl(i,j)=qagl(i,j)/1000.0
1312 pv=qagl(i,j)*pagl(i,j)/(eps*(1-qagl(i,j)) + qagl(i,j))
1313 rho=(1/tagl(i,j))*(((pagl(i,j)-pv)/rd) + pv/461.495)
1314 grid1(i,j)=0.5*rho*(sqrt(uagl(i,j)**2+vagl(i,j)**2))**3
1315 ELSE
1316 grid1(i,j)=spval
1317 ENDIF
1318 ENDDO
1319 ENDDO
1320 if(grib=="grib2" )then
1321 cfld=cfld+1
1322 fld_info(cfld)%ifld=iavblfld(iget(411))
1323 fld_info(cfld)%lvl=lvlsxml(lp,iget(411))
1324 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1325 endif
1326 ENDIF
1327!--- U Component of wind
1328 IF((iget(412)>0) ) THEN
1329 DO j=jsta,jend
1330 DO i=ista,iend
1331 grid1(i,j)=uagl(i,j)
1332 ENDDO
1333 ENDDO
1334 if(grib=="grib2" )then
1335 cfld=cfld+1
1336 fld_info(cfld)%ifld=iavblfld(iget(412))
1337 fld_info(cfld)%lvl=lvlsxml(lp,iget(412))
1338 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1339 endif
1340 ENDIF
1341!--- V Component of wind
1342 IF((iget(413)>0) ) THEN
1343 DO j=jsta,jend
1344 DO i=ista,iend
1345 grid1(i,j)=vagl(i,j)
1346 ENDDO
1347 ENDDO
1348 if(grib=="grib2" )then
1349 cfld=cfld+1
1350 fld_info(cfld)%ifld=iavblfld(iget(413))
1351 fld_info(cfld)%lvl=lvlsxml(lp,iget(413))
1352 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1353 endif
1354 ENDIF
1355!
1356 ENDIF ! FOR LEVEL
1357
1358!
1359!*** END OF MAIN VERTICAL LOOP
1360!
1361 330 CONTINUE
1362!*** ENDIF FOR IF TEST SEEING IF WE WANT ANY OTHER VARIABLES
1363!
1364 ENDIF
1365! CRA
1366!
1367! END OF ROUTINE.
1368!
1369 RETURN
1370 END
1371
mdl2agl
subroutine mdl2agl
SUBPROGRAM: MDL2P VERT INTRP OF MODEL LVLS TO AGL HEIGHT PRGRMMR: CHUANG ORG: W/NP22 DATE: 05-05-23
Definition MDL2AGL.f:50