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