UPP v11.0.0
Loading...
Searching...
No Matches
MDL2SIGMA.f
Go to the documentation of this file.
1
2!
53 SUBROUTINE mdl2sigma
54
55!
56!
57 use vrbls3d, only: pint, t, q, zint, alpint, pmid, exch_h, uh, &
58 vh, omga, q2, cwm, qqw, qqi, qqr, qqs, cfr, &
59 f_rimef, pmidv
60! use vrbls2d, only:
61 use masks, only: lmh
62 use params_mod, only: d50 , pq0, a2, a3, a4, h1, d01, d608, rgamog,&
63 h1m12, d00, h2, rd, g, gi, h99999
64 use ctlblk_mod, only: jsta_2l, jend_2u, spval, lp1, jsta, jend, lm, &
65 grib, cfld, datapd, fld_info, me, jend_m, im, &
66 jm, im_jm, ista, iend, ista_2l, iend_2u, ista_m, iend_m
67 use rqstfld_mod, only: iget, lvls, id, iavblfld, lvlsxml
68 use gridspec_mod, only :gridtype
69!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
70 implicit none
71!
72! INCLUDE MODEL DIMENSIONS. SET/DERIVE OTHER PARAMETERS.
73! GAMMA AND RGAMOG ARE USED IN THE EXTRAPOLATION OF VIRTUAL
74! TEMPERATURES BEYOND THE UPPER OF LOWER LIMITS OF DATA.
75!
76 integer,PARAMETER :: LSIG=22
77 real,PARAMETER :: PTSIGO=1.0e4
78!
79! DECLARE VARIABLES.
80!
81 LOGICAL READTHK
82 LOGICAL IOOMG,IOALL
83 LOGICAL DONEFSL1,TSLDONE
84 real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u) :: FSL, TSL, QSL, OSL, USL, VSL, Q2SL, &
85 fsl1, cfrsig, egrid1, egrid2
86 REAL GRID1(IM,JM)
87 real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u) :: grid2
88
89 REAL SIGO(LSIG+1),DSIGO(LSIG),ASIGO(LSIG)
90!
91 INTEGER,dimension(ista_2l:iend_2u,jsta_2l:jend_2u) :: NL1X,NL1XF
92!
93!
94!--- Definition of the following 2D (horizontal) dummy variables
95!
96! C1D - total condensate
97! QW1 - cloud water mixing ratio
98! QI1 - cloud ice mixing ratio
99! QR1 - rain mixing ratio
100! QS1 - snow mixing ratio
101!
102 real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u) :: C1D, QW1, QI1, QR1, QS1, QG1, AKH
103!
104 integer I,J,L,LL,LP,LLMH,II,JJ,JJB,JJE,NHOLD
105 real PFSIGO,APFSIGO,PSIGO,APSIGO,PNL1,PU,ZU,TU,QU,QSAT, &
106 rhu,tvru,tvrabv,tabv,qabv,b,ahf,fac,pl,zl,tl,ql, &
107 rhl,tmt0,ai,bi,tvrl,tvrblo,tblo,qblo,fact, &
108 px,bf,facf,ahff,dpsig,tv,pdv,denom,denomf,pnl1f,dum
109!
110!
111!******************************************************************************
112!
113! START MDL2P.
114!
115! SET TOTAL NUMBER OF POINTS ON OUTPUT GRID.
116!
117!---------------------------------------------------------------
118!
119! *** PART I ***
120!
121! VERTICAL INTERPOLATION OF EVERYTHING ELSE. EXECUTE ONLY
122! IF THERE'S SOMETHING WE WANT.
123!
124 IF((iget(205)>0).OR.(iget(206)>0).OR. &
125 (iget(207)>0).OR.(iget(208)>0).OR. &
126 (iget(209)>0).OR.(iget(210)>0).OR. &
127 (iget(216)>0).OR.(iget(217)>0).OR. &
128 (iget(211)>0).OR.(iget(212)>0).OR. &
129 (iget(213)>0).OR.(iget(214)>0).OR. &
130 (iget(215)>0).OR.(iget(222)>0).OR. &
131 (iget(243)>0) ) THEN !!Air Quality (Plee Oct2003)
132!
133!---------------------------------------------------------------------
134!
135!--- VERTICAL INTERPOLATION OF GEOPOTENTIAL, SPECIFIC HUMIDITY, TEMPERATURE,
136! OMEGA, TKE, & CLOUD FIELDS. START AT THE UPPERMOST TARGET SIGMA LEVEL.
137!
138 readthk=.false.
139 IF(readthk)THEN ! EITHER READ DSG THICKNESS
140 READ(41)dsigo !DSIGO FROM TOP TO BOTTOM
141!
142 sigo(1)=0.0
143 DO l=2,lsig+1
144 sigo(l)=sigo(l-1)+dsigo(lsig-l+2)
145 END DO
146 sigo(lsig+1)=1.0
147 DO l=1,lsig
148 asigo(l)=0.5*(sigo(l)+sigo(l+1))
149 END DO
150 ELSE ! SPECIFY SIGO
151 asigo( 1)= 0.0530
152 asigo( 2)= 0.1580
153 asigo( 3)= 0.2605
154 asigo( 4)= 0.3595
155 asigo( 5)= 0.4550
156 asigo( 6)= 0.5470
157 asigo( 7)= 0.6180
158 asigo( 8)= 0.6690
159 asigo( 9)= 0.7185
160 asigo(10)= 0.7585
161 asigo(11)= 0.7890
162 asigo(12)= 0.8190
163 asigo(13)= 0.8480
164 asigo(14)= 0.8755
165 asigo(15)= 0.9015
166 asigo(16)= 0.9260
167 asigo(17)= 0.9490
168 asigo(18)= 0.9650
169 asigo(19)= 0.9745
170 asigo(20)= 0.9835
171 asigo(21)= 0.9915
172 asigo(22)= 0.9975
173!
174 sigo( 1)= 0.0
175 sigo( 2)= 0.1060
176 sigo( 3)= 0.2100
177 sigo( 4)= 0.3110
178 sigo( 5)= 0.4080
179 sigo( 6)= 0.5020
180 sigo( 7)= 0.5920
181 sigo( 8)= 0.6440
182 sigo( 9)= 0.6940
183 sigo(10)= 0.7430
184 sigo(11)= 0.7740
185 sigo(12)= 0.8040
186 sigo(13)= 0.8340
187 sigo(14)= 0.8620
188 sigo(15)= 0.8890
189 sigo(16)= 0.9140
190 sigo(17)= 0.9380
191 sigo(18)= 0.9600
192 sigo(19)= 0.9700
193 sigo(20)= 0.9790
194 sigo(21)= 0.9880
195 sigo(22)= 0.9950
196 sigo(23)= 1.0
197 END IF
198! OBTAIN GEOPOTENTIAL AT 1ST LEVEL
199 DO j=jsta_2l,jend_2u
200 DO i=ista_2l,iend_2u
201 fsl(i,j)=spval
202 akh(i,j)=spval
203 nl1xf(i,j)=lp1
204 DO l=1,lp1
205 IF(nl1xf(i,j)==lp1.AND.pint(i,j,l)>ptsigo)THEN
206 nl1xf(i,j)=l
207 ENDIF
208 ENDDO
209 END DO
210 END DO
211 DO 167 j=jsta,jend
212 DO 167 i=ista_2l,iend_2u
213 donefsl1=.false.
214 pfsigo=ptsigo
215 apfsigo=log(pfsigo)
216 pnl1=pint(i,j,nl1xf(i,j))
217 ll=nl1xf(i,j)
218 llmh = nint(lmh(i,j))
219 IF(nl1xf(i,j)==1 .AND. t(i,j,1)<spval &
220 .AND. t(i,j,2)<spval .AND. q(i,j,1)<spval &
221 .AND. q(i,j,2)<spval)THEN
222 pu=pint(i,j,2)
223 zu=zint(i,j,2)
224 tu=d50*(t(i,j,1)+t(i,j,2))
225 qu=d50*(q(i,j,1)+q(i,j,2))
226 qsat=pq0/pu*exp(a2*(tu-a3)/(tu-a4))
227 rhu =qu/qsat
228 IF(rhu>h1)THEN
229 rhu=h1
230 qu =rhu*qsat
231 ENDIF
232 IF(rhu<d01)THEN
233 rhu=d01
234 qu =rhu*qsat
235 ENDIF
236!
237 tvru=tu*(h1+d608*qu)
238 tvrabv=tvru*(pfsigo/pu)**rgamog
239 tabv=tvrabv/(h1+d608*qu)
240 qsat=pq0/pfsigo*exp(a2*(tabv-a3)/(tabv-a4))
241 qabv =rhu*qsat
242 qabv =max(h1m12,qabv)
243!== B =TABV
244 b =tvrabv !Marina Tsidulko Dec22, 2003
245 ahf =d00
246 fac =d00
247 donefsl1=.true.
248 ELSEIF(nl1xf(i,j)==lp1 .AND. t(i,j,lm-1)<spval &
249 .AND. t(i,j,lm-2)<spval .AND. q(i,j,lm-1)<spval &
250 .AND. q(i,j,lm-2)<spval)THEN
251 pl=pint(i,j,lm-1)
252 zl=zint(i,j,lm-1)
253 tl=d50*(t(i,j,lm-2)+t(i,j,lm-1))
254 ql=d50*(q(i,j,lm-2)+q(i,j,lm-1))
255!
256!--- RH w/r/t water for all conditions (Ferrier, 25 Jan 02)
257!
258 qsat=pq0/pl*exp(a2*(tl-a3)/(tl-a4))
259 rhl=ql/qsat
260 IF(rhl>h1)THEN
261 rhl=h1
262 ql =rhl*qsat
263 ENDIF
264 IF(rhl<d01)THEN
265 rhl=d01
266 ql =rhl*qsat
267 ENDIF
268!
269 tvrl =tl*(h1+d608*ql)
270 tvrblo=tvrl*(pfsigo/pl)**rgamog
271 tblo =tvrblo/(h1+d608*ql)
272 qsat=pq0/pfsigo*exp(a2*(tblo-a3)/(tblo-a4))
273 qblo =rhl*qsat
274 qblo =max(h1m12,qblo)
275!== B =TBLO
276 b =tvrblo !Marina Tsidulko Dec22, 2003
277 ahf =d00
278 fac =d00
279 donefsl1=.true.
280 ELSEIF(t(i,j,nl1xf(i,j))<spval &
281 & .AND. q(i,j,nl1xf(i,j))<spval)THEN
282!== B =T(I,J,NL1XF(I,J)) !Marina Tsidulko Dec22, 2003
283 b =t(i,j,nl1xf(i,j))*(h1+d608*q(i,j,nl1xf(i,j)))
284 denom=(alpint(i,j,nl1xf(i,j)+1)-alpint(i,j,nl1xf(i,j)-1))
285!== AHF =(B-T(I,J,NL1XF(I,J)-1))/DENOM
286 ahf =(b-t(i,j,nl1xf(i,j)-1)*(h1+d608*q(i,j,nl1xf(i,j)-1))) &
287 /denom !Marina Tsidulko Dec22, 2003
288 fac =h2*log(pmid(i,j,nl1xf(i,j)))
289 donefsl1=.true.
290 END IF
291!
292 if(donefsl1)fsl1(i,j)=(pnl1-pfsigo)/(pfsigo+pnl1) &
293 *((apfsigo+alpint(i,j,nl1xf(i,j))-fac)*ahf+b)*rd*h2 &
294 +zint(i,j,nl1xf(i,j))*g
295! COMPUTE EXCHANGE COEFFICIENT ON FIRST INTERFACTE
296 IF(nl1xf(i,j)<=2 .OR. nl1xf(i,j)>(llmh+1))THEN
297 akh(i,j)=0.0
298 ELSE
299 fact=(apfsigo-log(pint(i,j,ll)))/ &
300 & (log(pint(i,j,ll))-log(pint(i,j,ll-1)))
301! EXCH_H is on the bottom of model interfaces
302 IF(exch_h(i,j,ll-2)<spval .AND. exch_h(i,j,ll-1)<spval) &
303 & akh(i,j)=exch_h(i,j,ll-1)+(exch_h(i,j,ll-1) &
304 & -exch_h(i,j,ll-2))*fact
305 END IF
306 167 CONTINUE
307! OUTPUT FIRST LAYER GEOPOTENTIAL
308! GEOPOTENTIAL (SCALE BY GI)
309 IF (iget(205)>0) THEN
310 IF (lvls(1,iget(205))>0) THEN
311!$omp parallel do
312 DO j=jsta,jend
313 DO i=ista,iend
314 IF(fsl1(i,j)<spval) THEN
315 grid1(i,j)=fsl1(i,j)*gi
316 ELSE
317 grid1(i,j)=spval
318 ENDIF
319 ENDDO
320 ENDDO
321 if(grib=='grib2')then
322 cfld=cfld+1
323 fld_info(cfld)%ifld=iavblfld(iget(205))
324 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
325 endif
326 ENDIF
327 ENDIF
328
329! OUTPUT FIRST INTERFACE KH Heat Diffusivity
330 IF (iget(243)>0) THEN !!Air Quality (Plee Oct2003) ^^^^^
331 IF (lvls(1,iget(243))>0) THEN
332!$omp parallel do
333 DO j=jsta,jend
334 DO i=ista,iend
335 grid1(i,j)=akh(i,j)
336 ENDDO
337 ENDDO
338 if(grib=="grib2" )then
339 cfld=cfld+1
340 fld_info(cfld)%ifld=iavblfld(iget(243))
341 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
342 endif
343
344 if(me==0)print*,'output Heat Diffusivity'
345 ENDIF
346 ENDIF
347
348!***
349!*** BECAUSE SIGMA LAYERS DO NOT GO UNDERGROUND, DO ALL
350!*** INTERPOLATION ABOVE GROUND NOW.
351!***
352!
353 DO 310 lp=1,lsig
354 nhold=0
355!
356 DO j=jsta_2l,jend_2u
357 DO i=ista_2l,iend_2u
358
359!
360 tsl(i,j)=spval
361 qsl(i,j)=spval
362 fsl(i,j)=spval
363 osl(i,j)=spval
364 usl(i,j)=spval
365 vsl(i,j)=spval
366 q2sl(i,j)=spval
367 c1d(i,j)=spval ! Total condensate
368 qw1(i,j)=spval ! Cloud water
369 qi1(i,j)=spval ! Cloud ice
370 qr1(i,j)=spval ! Rain
371 qs1(i,j)=spval ! Snow (precip ice)
372 qg1(i,j)=spval
373 cfrsig(i,j)=spval
374!
375!*** LOCATE VERTICAL INDEX OF MODEL MIDLAYER JUST BELOW
376!*** THE PRESSURE LEVEL TO WHICH WE ARE INTERPOLATING.
377!
378 nl1x(i,j)=lp1
379 DO l=2,lm
380 llmh = nint(lmh(i,j))
381 psigo=ptsigo+asigo(lp)*(pint(i,j,llmh+1)-ptsigo)
382 IF(nl1x(i,j)==lp1.AND.pmid(i,j,l)>psigo)THEN
383 nl1x(i,j)=l
384 ENDIF
385 ENDDO
386!
387! IF THE PRESSURE LEVEL IS BELOW THE LOWEST MODEL MIDLAYER
388! BUT STILL ABOVE THE LOWEST MODEL BOTTOM INTERFACE,
389! WE WILL NOT CONSIDER IT UNDERGROUND AND THE INTERPOLATION
390! WILL EXTRAPOLATE TO THAT POINT
391!
392 IF(nl1x(i,j)==lp1.AND.pint(i,j,llmh+1)>=psigo)THEN
393 nl1x(i,j)=lm
394 ENDIF
395!
396! if(NL1X(I,J)==LP1)print*,'Debug: NL1X=LP1 AT '
397! 1 ,i,j,lp
398 ENDDO
399 ENDDO
400!
401!mptest IF(NHOLD==0)GO TO 310
402!
403!$omp parallel do private(i,j,ll,llmh,psigo,apsigo,fact,dum,pl, &
404!$omp & zl,tl,ql,ai,bi,qsat,rhl,tvrl,tvrblo,tblo,tmt0, &
405!$omp & qblo,pnl1,fac,ahf)
406!hc DO 220 NN=1,NHOLD
407!hc I=IHOLD(NN)
408!hc J=JHOLD(NN)
409 DO 220 j=jsta,jend ! Moorthi on Nov 26 2014
410! DO 220 J=JSTA_2L,JEND_2U
411 DO 220 i=ista,iend
412 ll=nl1x(i,j)
413!---------------------------------------------------------------------
414!*** VERTICAL INTERPOLATION OF GEOPOTENTIAL, TEMPERATURE, SPECIFIC
415!*** HUMIDITY, CLOUD WATER/ICE, OMEGA, WINDS, AND TKE.
416!---------------------------------------------------------------------
417!
418!HC IF(NL1X(I,J)<=LM)THEN
419 llmh = nint(lmh(i,j))
420 psigo=ptsigo+asigo(lp)*(pint(i,j,llmh+1)-ptsigo)
421 apsigo=log(psigo)
422 IF(nl1x(i,j)<=llmh)THEN
423!
424!---------------------------------------------------------------------
425! INTERPOLATE LINEARLY IN LOG(P)
426!*** EXTRAPOLATE ABOVE THE TOPMOST MIDLAYER OF THE MODEL
427!*** INTERPOLATION BETWEEN NORMAL LOWER AND UPPER BOUNDS
428!*** EXTRAPOLATE BELOW LOWEST MODEL MIDLAYER (BUT STILL ABOVE GROUND)
429!---------------------------------------------------------------------
430!
431
432 fact=(apsigo-log(pmid(i,j,ll)))/ &
433 & (log(pmid(i,j,ll))-log(pmid(i,j,ll-1)))
434 tsl(i,j)=t(i,j,ll)+(t(i,j,ll)-t(i,j,ll-1))*fact
435 IF(q(i,j,ll)<spval .AND. q(i,j,ll-1)<spval) &
436 & qsl(i,j)=q(i,j,ll)+(q(i,j,ll)-q(i,j,ll-1))*fact
437 IF(gridtype=='A')THEN
438 IF(uh(i,j,ll)<spval .AND. uh(i,j,ll-1)<spval) &
439 & usl(i,j)=uh(i,j,ll)+(uh(i,j,ll)-uh(i,j,ll-1))*fact
440 IF(vh(i,j,ll)<spval .AND. vh(i,j,ll-1)<spval) &
441 & vsl(i,j)=vh(i,j,ll)+(vh(i,j,ll)-vh(i,j,ll-1))*fact
442 END IF
443 IF(omga(i,j,ll)<spval .AND. omga(i,j,ll-1)<spval) &
444 & osl(i,j)=omga(i,j,ll)+(omga(i,j,ll)-omga(i,j,ll-1))*fact
445 IF(q2(i,j,ll)<spval .AND. q2(i,j,ll-1)<spval) &
446 & q2sl(i,j)=q2(i,j,ll)+(q2(i,j,ll)-q2(i,j,ll-1))*fact
447! FSL(I,J)=ZMID(I,J,LL)+(ZMID(I,J,LL)-ZMID(I,J,LL-1))*FACT
448! FSL(I,J)=FSL(I,J)*G
449!hc QSAT=PQ0/PSIGO
450!hc 1 *EXP(A2*(TSL(I,J)-A3)/(TSL(I,J)-A4))
451!
452!hc RHL=QSL(I,J)/QSAT
453!
454!hc IF(RHL>1.) QSL(I,J)=QSAT
455!hc IF(RHL<0.01) QSL(I,J)=0.01*QSAT
456 IF(q2sl(i,j)<0.0) q2sl(i,j)=0.0
457!
458!HC ADD FERRIER'S HYDROMETEOR
459 IF(cwm(i,j,ll)<spval .AND. cwm(i,j,ll-1)<spval) &
460 & c1d(i,j)=cwm(i,j,ll)+(cwm(i,j,ll)-cwm(i,j,ll-1))*fact
461 c1d(i,j)=max(c1d(i,j),h1m12) ! Total condensate
462 IF(qqw(i,j,ll)<spval .AND. qqw(i,j,ll-1)<spval) &
463 & qw1(i,j)=qqw(i,j,ll)+(qqw(i,j,ll)-qqw(i,j,ll-1))*fact
464 qw1(i,j)=max(qw1(i,j),h1m12) ! Cloud water
465 IF(qqi(i,j,ll)<spval .AND. qqi(i,j,ll-1)<spval) &
466 & qi1(i,j)=qqi(i,j,ll)+(qqi(i,j,ll)-qqi(i,j,ll-1))*fact
467 qi1(i,j)=max(qi1(i,j),h1m12) ! Cloud ice
468 IF(qqr(i,j,ll)<spval .AND. qqr(i,j,ll-1)<spval) &
469 & qr1(i,j)=qqr(i,j,ll)+(qqr(i,j,ll)-qqr(i,j,ll-1))*fact
470 qr1(i,j)=max(qr1(i,j),h1m12) ! Rain
471 IF(qqs(i,j,ll)<spval .AND. qqs(i,j,ll-1)<spval) &
472 & qs1(i,j)=qqs(i,j,ll)+(qqs(i,j,ll)-qqs(i,j,ll-1))*fact
473 qs1(i,j)=max(qs1(i,j),h1m12) ! Snow (precip ice)
474 IF(cfr(i,j,ll)<spval .AND. cfr(i,j,ll-1)<spval) &
475 & cfrsig(i,j)=cfr(i,j,ll)+(cfr(i,j,ll)-cfr(i,j,ll-1))*fact
476 cfrsig(i,j)=max(cfrsig(i,j),h1m12)
477 IF(qqs(i,j,ll)<spval .AND. qqs(i,j,ll-1)<spval)THEN
478 dum=f_rimef(i,j,ll)+(f_rimef(i,j,ll)-f_rimef(i,j,ll-1))*fact
479 IF(dum <= 5.0)THEN
480 qg1(i,j)=h1m12
481 ELSE
482 qg1(i,j)=qs1(i,j)
483 qs1(i,j)=h1m12
484 END IF
485 END IF
486! FOR UNDERGROUND PRESSURE LEVELS, ASSUME TEMPERATURE TO CHANGE
487! ADIABATICLY, RH TO BE THE SAME AS THE AVERAGE OF THE 2ND AND 3RD
488! LAYERS FROM THE GOUND, WIND TO BE THE SAME AS THE LOWEST LEVEL ABOVE
489! GOUND
490 ELSE
491! ii=91
492! jj=13
493! if(i==ii.and.j==jj)print*,'Debug: underg extra at i,j,lp' &
494! &, i,j,lp
495 pl=pint(i,j,lm-1)
496 zl=zint(i,j,lm-1)
497 tl=0.5*(t(i,j,lm-2)+t(i,j,lm-1))
498 ql=0.5*(q(i,j,lm-2)+q(i,j,lm-1))
499 ai=0.008855
500 bi=1.
501 tmt0=tl-a3
502 IF(tmt0<-20.)THEN
503 ai=0.007225
504 bi=0.9674
505 ENDIF
506 qsat=pq0/pl*exp(a2*(tl-a3)/(tl-a4))
507!
508 rhl=ql/qsat
509!
510 IF(rhl>1.)THEN
511 rhl=1.
512 ql =rhl*qsat
513 ENDIF
514!
515 IF(rhl<0.01)THEN
516 rhl=0.01
517 ql =rhl*qsat
518 ENDIF
519!
520 tvrl =tl*(1.+0.608*ql)
521 tvrblo=tvrl*(psigo/pl)**rgamog
522 tblo =tvrblo/(1.+0.608*ql)
523!
524 tmt0=tblo-a3
525 ai=0.008855
526 bi=1.
527 IF(tmt0<-20.)THEN
528 ai=0.007225
529 bi=0.9674
530 ENDIF
531 qsat=pq0/psigo*exp(a2*(tblo-a3)/(tblo-a4))
532!
533! TSL(I,J)=TBLO
534 qblo = rhl*qsat
535 qsl(i,j) = max(1.e-12,qblo)
536 IF(gridtype=='A')THEN
537 usl(i,j) = uh(i,j,llmh)
538 vsl(i,j) = vh(i,j,llmh)
539 END IF
540 osl(i,j) = omga(i,j,llmh)
541 q2sl(i,j) = max(0.0,0.5*(q2(i,j,llmh-1)+q2(i,j,llmh)))
542 pnl1 = pint(i,j,nl1x(i,j))
543 fac = 0.
544 ahf = 0.0
545!
546!--- Set hydrometeor fields to zero below ground
547 c1d(i,j)=0.
548 qw1(i,j)=0.
549 qi1(i,j)=0.
550 qr1(i,j)=0.
551 qs1(i,j)=0.
552 qg1(i,j)=0.
553 cfrsig(i,j)=0.
554 END IF
555 220 CONTINUE
556!
557! OBTAIN GEOPOTENTIAL AND KH ON INTERFACES
558 DO j=jsta_2l,jend_2u
559 DO i=ista_2l,iend_2u
560 fsl(i,j)=spval
561 akh(i,j)=spval
562 nl1xf(i,j)=lp1
563 llmh = nint(lmh(i,j))
564 psigo=ptsigo+sigo(lp+1)*(pint(i,j,llmh+1)-ptsigo)
565 DO l=1,lp1
566 IF(nl1xf(i,j)==lp1.AND.pint(i,j,l)>psigo)THEN
567 nl1xf(i,j)=l
568 ENDIF
569 ENDDO
570 END DO
571 END DO
572!
573! DO J=JSTA_2L,JEND_2U
574 DO j=jsta,jend ! Moorthi on 26 Nov 2014
575 DO i=ista,iend
576 donefsl1=.false.
577 tsldone=.false.
578 llmh = nint(lmh(i,j))
579 pfsigo=ptsigo+sigo(lp+1)*(pint(i,j,llmh+1)-ptsigo)
580 psigo=ptsigo+asigo(lp)*(pint(i,j,llmh+1)-ptsigo)
581 apfsigo=log(pfsigo)
582 pnl1f=pint(i,j,nl1xf(i,j))
583 ll=nl1xf(i,j)
584 IF(nl1xf(i,j)==1 .AND. t(i,j,1)<spval &
585 & .AND. t(i,j,2)<spval .AND. q(i,j,1)<spval &
586 & .AND. q(i,j,2)<spval)THEN
587 pu=pint(i,j,2)
588 zu=zint(i,j,2)
589 tu=d50*(t(i,j,1)+t(i,j,2))
590 qu=d50*(q(i,j,1)+q(i,j,2))
591!
592!--- RH w/r/t water for all conditions
593!
594 qsat=pq0/pu*exp(a2*(tu-a3)/(tu-a4))
595 rhu =qu/qsat
596 IF(rhu>h1)THEN
597 rhu=h1
598 qu =rhu*qsat
599 ENDIF
600 IF(rhu<d01)THEN
601 rhu=d01
602 qu =rhu*qsat
603 ENDIF
604!
605 tvru=tu*(h1+d608*qu)
606! TVRABV=TVRU*(SPL(L)/PU)**RGAMOG
607!== TVRABV=TVRU*(PFSIGO/PU)**RGAMOG
608 px=(pfsigo+pnl1f)*0.5 !Marina Tsidulko Dec22, 2003
609 tvrabv=tvru*(px/pu)**rgamog!Marina Tsidulko Dec22, 2003
610 tabv=tvrabv/(h1+d608*qu)
611! ADD ADDITIONAL BLOCK FOR INTERPOLATING HEIGHT ONTO INTERFACES
612!== BF =TABV
613 bf =tvrabv !Marina Tsidulko Dec22, 2003
614 facf =d00
615 ahff =d00
616 donefsl1=.true.
617!
618!
619 ELSEIF(nl1xf(i,j)==lp1 .AND. t(i,j,lm-1)<spval &
620 & .AND. t(i,j,lm-2)<spval .AND. q(i,j,lm-1)<spval &
621 & .AND. q(i,j,lm-2)<spval)THEN
622!
623! EXTRAPOLATION AT LOWER BOUND. THE LOWER BOUND IS
624! LM IF OLDRD=.FALSE. IF OLDRD=.TRUE. THE LOWER
625! BOUND IS THE FIRST ATMOSPHERIC ETA LAYER.
626!
627 pl=pint(i,j,lm-1)
628 zl=zint(i,j,lm-1)
629 tl=d50*(t(i,j,lm-2)+t(i,j,lm-1))
630 ql=d50*(q(i,j,lm-2)+q(i,j,lm-1))
631!
632!--- RH w/r/t water for all conditions (Ferrier, 25 Jan 02)
633!
634 qsat=pq0/pl*exp(a2*(tl-a3)/(tl-a4))
635 rhl=ql/qsat
636 IF(rhl>h1)THEN
637 rhl=h1
638 ql =rhl*qsat
639 ENDIF
640 IF(rhl<d01)THEN
641 rhl=d01
642 ql =rhl*qsat
643 ENDIF
644!
645 tvrl =tl*(h1+d608*ql)
646!== TVRBLO=TVRL*(PFSIGO/PL)**RGAMOG
647 px=(pfsigo+pnl1f)*0.5 !Marina Tsidulko Dec22, 2003
648 tvrblo=tvrl*(px/pl)**rgamog !Marina Tsidulko Dec22, 2003
649 tblo =tvrblo/(h1+d608*ql)
650! ADD ADDITIONAL BLOCK FOR INTERPOLATING HEIGHT ONTO INTERFACES
651!== BF =TBLO
652 bf =tvrblo !Marina Tsidulko Dec22, 2003
653 facf =d00
654 ahff =d00
655 donefsl1=.true.
656 tsldone=.true.
657!
658!
659 ELSEIF(t(i,j,nl1xf(i,j))<spval &
660 & .AND. q(i,j,nl1xf(i,j))<spval)THEN
661!
662! INTERPOLATION BETWEEN LOWER AND UPPER BOUNDS.
663!
664! ADD ADDITIONAL BLOCK FOR INTERPOLATING HEIGHT ONTO INTERFACES
665! if(NL1XF(I,J)==LMp1)
666! + print*,'Debug: using bad temp at',i,j
667!== BF =T(I,J,NL1XF(I,J)) !Marina Tsidulko Dec22, 2003
668 bf =t(i,j,nl1xf(i,j))*(h1+d608*q(i,j,nl1xf(i,j)))
669!HC FACF =H2*LOG(PT+PDSL(I,J)*AETA(NL1XF(I,J)))
670 facf =h2*log(pmid(i,j,nl1xf(i,j)))
671 denomf=(alpint(i,j,nl1xf(i,j)+1)-alpint(i,j,nl1xf(i,j)-1))
672!== AHFF=(BF-T(I,J,NL1XF(I,J)-1))/DENOMF
673 ahff=(bf-t(i,j,nl1xf(i,j)-1)*(h1+d608*q(i,j,nl1xf(i,j)-1))) &
674 /denomf !Marina Tsidulko Dec22, 2003
675!
676 donefsl1=.true.
677 ENDIF
678!
679 IF(donefsl1)then
680 fsl(i,j)=(pnl1f-pfsigo)/(pfsigo+pnl1f) &
681 *((apfsigo+alpint(i,j,nl1xf(i,j))-facf)*ahff+bf)*rd*h2 &
682 +zint(i,j,nl1xf(i,j))*g
683! OBTAIN TEMPERATURE AT MID-SIGMA LAYER BASED ON HYDROSTATIC
684 dpsig=(sigo(lp+1)-sigo(lp))*(pint(i,j,llmh+1)-ptsigo)
685! IF(J==jj.and.i==ii)print*,'Debug:L,OLD T= ',
686! + L,TSL(I,J)
687 IF(.NOT.tsldone) THEN
688 tsl(i,j)=(fsl1(i,j)-fsl(i,j))*psigo/(rd*dpsig)
689 ENDIF
690 fsl1(i,j)=fsl(i,j)
691 IF(.NOT.tsldone) THEN
692 tv=tsl(i,j)
693 tsl(i,j)=tv/(h1+d608*qsl(i,j))
694 ENDIF
695 qsat=pq0/psigo *exp(a2*(tsl(i,j)-a3)/(tsl(i,j)-a4))
696!
697 rhl=qsl(i,j)/qsat
698!
699 IF(rhl>1.) qsl(i,j)=qsat
700 IF(rhl<0.01) qsl(i,j)=0.01*qsat
701 END IF
702! IF(J==jj.and.i==ii)print*,'Debug:L,T,Q,Q2,FSL=',
703! + L,TSL(I,J),QSL(I,J),Q2SL(I,J),FSL(I,J)
704
705! COMPUTE EXCHANGE COEFFICIENT ON INTERFACES
706 IF(nl1xf(i,j)<=2 .OR. nl1xf(i,j)>(llmh+1))THEN
707 akh(i,j)=0.0
708 ELSE
709 fact=(apfsigo-log(pint(i,j,ll)))/ &
710 & (log(pint(i,j,ll))-log(pint(i,j,ll-1)))
711! EXCH_H is on the bottom of model interfaces
712 IF(exch_h(i,j,ll-2)<spval .AND. &
713 & exch_h(i,j,ll-1)<spval) &
714 & akh(i,j)=exch_h(i,j,ll-1)+(exch_h(i,j,ll-1) &
715 & -exch_h(i,j,ll-2))*fact
716 END IF
717
718! END OF HYDROSTATIC TEMPERATURE DERIVATION
719 END DO
720 END DO
721!
722! VERTICAL INTERPOLATION FOR WIND FOR E and B GRIDS
723!
724 if(gridtype=='B' .or. gridtype=='E') &
725 call exch(pint(ista_2l:iend_2u,jsta_2l:jend_2u,lp1))
726 IF(gridtype=='E')THEN
727 DO j=jsta,jend
728! DO I=1,IM-MOD(J,2)
729 DO i=ista,iend-mod(j,2) !Jesse 20211014
730!
731!*** LOCATE VERTICAL INDEX OF MODEL MIDLAYER FOR V POINT JUST BELOW
732!*** THE PRESSURE LEVEL TO WHICH WE ARE INTERPOLATING.
733!
734 llmh = nint(lmh(i,j))
735
736!Jesse 20211014
737! IF(J == 1 .AND. I < IM)THEN !SOUTHERN BC
738! PDV=0.5*(PINT(I,J,LLMH+1)+PINT(I+1,J,LLMH+1))
739! ELSE IF(J==JM .AND. I<IM)THEN !NORTHERN BC
740! PDV=0.5*(PINT(I,J,LLMH+1)+PINT(I+1,J,LLMH+1))
741! ELSE IF(I == 1 .AND. MOD(J,2) == 0) THEN !WESTERN EVEN BC
742! PDV=0.5*(PINT(I,J-1,LLMH+1)+PINT(I,J+1,LLMH+1))
743! ELSE IF(I == IM .AND. MOD(J,2) == 0) THEN !EASTERN EVEN BC
744! PDV=0.5*(PINT(I,J-1,LLMH+1)+PINT(I,J+1,LLMH+1))
745! ELSE IF (MOD(J,2) < 1) THEN
746! PDV=0.25*(PINT(I,J,LLMH+1)+PINT(I-1,J,LLMH+1) &
747! & +PINT(I,J+1,LLMH+1)+PINT(I,J-1,LLMH+1))
748! ELSE
749! PDV=0.25*(PINT(I,J,LLMH+1)+PINT(I+1,J,LLMH+1) &
750! & +PINT(I,J+1,LLMH+1)+PINT(I,J-1,LLMH+1))
751! END IF
752
753 IF(j == 1 .AND. i < iend)THEN !SOUTHERN BC
754 pdv=0.5*(pint(i,j,llmh+1)+pint(i+1,j,llmh+1))
755 ELSE IF(j==jm .AND. i<iend)THEN !NORTHERN BC
756 pdv=0.5*(pint(i,j,llmh+1)+pint(i+1,j,llmh+1))
757 ELSE IF(i == ista .AND. mod(j,2) == 0) THEN !WESTERN EVEN BC
758 pdv=0.5*(pint(i,j-1,llmh+1)+pint(i,j+1,llmh+1))
759 ELSE IF(i == iend .AND. mod(j,2) == 0) THEN !EASTERN EVEN BC
760 pdv=0.5*(pint(i,j-1,llmh+1)+pint(i,j+1,llmh+1))
761 ELSE IF (mod(j,2) < 1) THEN
762 pdv=0.25*(pint(i,j,llmh+1)+pint(i-1,j,llmh+1) &
763 & +pint(i,j+1,llmh+1)+pint(i,j-1,llmh+1))
764 ELSE
765 pdv=0.25*(pint(i,j,llmh+1)+pint(i+1,j,llmh+1) &
766 & +pint(i,j+1,llmh+1)+pint(i,j-1,llmh+1))
767 END IF
768
769 psigo=ptsigo+asigo(lp)*(pdv-ptsigo)
770 nl1x(i,j)=lp1
771 DO l=2,lm
772 IF(nl1x(i,j)==lp1.AND.pmidv(i,j,l)>psigo)THEN
773 nl1x(i,j)=l
774 ENDIF
775 ENDDO
776!
777! IF THE PRESSURE LEVEL IS BELOW THE LOWEST MODEL MIDLAYER
778! BUT STILL ABOVE THE LOWEST MODEL BOTTOM INTERFACE,
779! WE WILL NOT CONSIDER IT UNDERGROUND AND THE INTERPOLATION
780! WILL EXTRAPOLATE TO THAT POINT
781!
782 IF(nl1x(i,j)==lp1.AND. pdv>psigo)THEN
783 nl1x(i,j)=lm
784 ENDIF
785!
786 ENDDO
787 ENDDO
788!
789 DO 230 j=jsta,jend
790! DO 230 I=1,IM-MOD(j,2)
791 DO 230 i=ista,iend-mod(j,2) !Jesse 20211014
792
793 llmh = nint(lmh(i,j))
794
795!Jesse 20211014
796! IF(J == 1 .AND. I < IM)THEN !SOUTHERN BC
797! PDV=0.5*(PINT(I,J,LLMH+1)+PINT(I+1,J,LLMH+1))
798! ELSE IF(J==JM .AND. I<IM)THEN !NORTHERN BC
799! PDV=0.5*(PINT(I,J,LLMH+1)+PINT(I+1,J,LLMH+1))
800! ELSE IF(I == 1 .AND. MOD(J,2) == 0) THEN !WESTERN EVEN BC
801! PDV=0.5*(PINT(I,J-1,LLMH+1)+PINT(I,J+1,LLMH+1))
802! ELSE IF(I == IM .AND. MOD(J,2) == 0) THEN !EASTERN EVEN BC
803! PDV=0.5*(PINT(I,J-1,LLMH+1)+PINT(I,J+1,LLMH+1))
804! ELSE IF (MOD(J,2) < 1) THEN
805! PDV=0.25*(PINT(I,J,LLMH+1)+PINT(I-1,J,LLMH+1) &
806! & +PINT(I,J+1,LLMH+1)+PINT(I,J-1,LLMH+1))
807! ELSE
808! PDV=0.25*(PINT(I,J,LLMH+1)+PINT(I+1,J,LLMH+1) &
809! & +PINT(I,J+1,LLMH+1)+PINT(I,J-1,LLMH+1))
810! END IF
811
812 IF(j == 1 .AND. i < iend)THEN !SOUTHERN BC
813 pdv=0.5*(pint(i,j,llmh+1)+pint(i+1,j,llmh+1))
814 ELSE IF(j==jm .AND. i<iend)THEN !NORTHERN BC
815 pdv=0.5*(pint(i,j,llmh+1)+pint(i+1,j,llmh+1))
816 ELSE IF(i == ista .AND. mod(j,2) == 0) THEN !WESTERN EVEN BC
817 pdv=0.5*(pint(i,j-1,llmh+1)+pint(i,j+1,llmh+1))
818 ELSE IF(i == iend .AND. mod(j,2) == 0) THEN !EASTERN EVEN BC
819 pdv=0.5*(pint(i,j-1,llmh+1)+pint(i,j+1,llmh+1))
820 ELSE IF (mod(j,2) < 1) THEN
821 pdv=0.25*(pint(i,j,llmh+1)+pint(i-1,j,llmh+1) &
822 & +pint(i,j+1,llmh+1)+pint(i,j-1,llmh+1))
823 ELSE
824 pdv=0.25*(pint(i,j,llmh+1)+pint(i+1,j,llmh+1) &
825 & +pint(i,j+1,llmh+1)+pint(i,j-1,llmh+1))
826 END IF
827
828 psigo=ptsigo+asigo(lp)*(pdv-ptsigo)
829 apsigo=log(psigo)
830 ll=nl1x(i,j)
831!---------------------------------------------------------------------
832!*** VERTICAL INTERPOLATION OF WINDS FOR A-E GRID
833!---------------------------------------------------------------------
834!
835!HC IF(NL1X(I,J)<=LM)THEN
836 llmh = nint(lmh(i,j))
837 IF(nl1x(i,j)<=llmh)THEN
838!
839!---------------------------------------------------------------------
840! INTERPOLATE LINEARLY IN LOG(P)
841!*** EXTRAPOLATE ABOVE THE TOPMOST MIDLAYER OF THE MODEL
842!*** INTERPOLATION BETWEEN NORMAL LOWER AND UPPER BOUNDS
843!*** EXTRAPOLATE BELOW LOWEST MODEL MIDLAYER (BUT STILL ABOVE GROUND)
844!---------------------------------------------------------------------
845!
846
847 fact=(apsigo-log(pmidv(i,j,ll)))/ &
848 & (log(pmidv(i,j,ll))-log(pmidv(i,j,ll-1)))
849 IF(uh(i,j,ll)<spval .AND. uh(i,j,ll-1)<spval) &
850 & usl(i,j)=uh(i,j,ll)+(uh(i,j,ll)-uh(i,j,ll-1))*fact
851 IF(vh(i,j,ll)<spval .AND. vh(i,j,ll-1)<spval) &
852 & vsl(i,j)=vh(i,j,ll)+(vh(i,j,ll)-vh(i,j,ll-1))*fact
853!
854! FOR UNDERGROUND PRESSURE LEVELS, ASSUME TEMPERATURE TO CHANGE
855! ADIABATICLY, RH TO BE THE SAME AS THE AVERAGE OF THE 2ND AND 3RD
856! LAYERS FROM THE GOUND, WIND TO BE THE SAME AS THE LOWEST LEVEL ABOVE
857! GOUND
858 ELSE
859 IF(uh(i,j,llmh)<spval)usl(i,j)=uh(i,j,llmh)
860 IF(vh(i,j,llmh)<spval)vsl(i,j)=vh(i,j,llmh)
861 END IF
862 230 CONTINUE
863 jjb=jsta
864 IF(mod(jsta,2)==0)jjb=jsta+1
865 jje=jend
866 IF(mod(jend,2)==0)jje=jend-1
867 DO j=jjb,jje,2 !chc
868 usl(iend,j)=usl(iend-1,j)
869 vsl(iend,j)=vsl(iend-1,j)
870 END DO
871
872 ELSE IF (gridtype=='B')THEN
873 DO j=jsta,jend_m
874 DO i=ista,iend_m
875!
876!*** LOCATE VERTICAL INDEX OF MODEL MIDLAYER FOR V POINT JUST BELOW
877!*** THE PRESSURE LEVEL TO WHICH WE ARE INTERPOLATING.
878!
879 pdv=0.25*(pint(i,j,lp1)+pint(i+1,j,lp1) &
880 +pint(i,j+1,lp1)+pint(i+1,j+1,lp1))
881
882 psigo=ptsigo+asigo(lp)*(pdv-ptsigo)
883 nl1x(i,j)=lp1
884 DO l=2,lm
885 IF(nl1x(i,j)==lp1.AND.pmidv(i,j,l)>psigo)THEN
886 nl1x(i,j)=l
887 ENDIF
888 ENDDO
889!
890! IF THE PRESSURE LEVEL IS BELOW THE LOWEST MODEL MIDLAYER
891! BUT STILL ABOVE THE LOWEST MODEL BOTTOM INTERFACE,
892! WE WILL NOT CONSIDER IT UNDERGROUND AND THE INTERPOLATION
893! WILL EXTRAPOLATE TO THAT POINT
894!
895 IF(nl1x(i,j)==lp1.AND. pdv>psigo)THEN
896 nl1x(i,j)=lm
897 ENDIF
898!
899 ENDDO
900 ENDDO
901!
902 DO 231 j=jsta,jend_m
903 DO 231 i=ista,iend_m
904 pdv=0.25*(pint(i,j,lp1)+pint(i+1,j,lp1) &
905 +pint(i,j+1,lp1)+pint(i+1,j+1,lp1))
906 psigo=ptsigo+asigo(lp)*(pdv-ptsigo)
907 apsigo=log(psigo)
908 ll=nl1x(i,j)
909!---------------------------------------------------------------------
910!*** VERTICAL INTERPOLATION OF WINDS FOR A-E GRID
911!---------------------------------------------------------------------
912!
913!HC IF(NL1X(I,J)<=LM)THEN
914 llmh = nint(lmh(i,j))
915 IF(nl1x(i,j)<=llmh)THEN
916!
917!---------------------------------------------------------------------
918! INTERPOLATE LINEARLY IN LOG(P)
919!*** EXTRAPOLATE ABOVE THE TOPMOST MIDLAYER OF THE MODEL
920!*** INTERPOLATION BETWEEN NORMAL LOWER AND UPPER BOUNDS
921!*** EXTRAPOLATE BELOW LOWEST MODEL MIDLAYER (BUT STILL ABOVE GROUND)
922!---------------------------------------------------------------------
923!
924
925 fact=(apsigo-log(pmidv(i,j,ll)))/ &
926 & (log(pmidv(i,j,ll))-log(pmidv(i,j,ll-1)))
927 IF(uh(i,j,ll)<spval .AND. uh(i,j,ll-1)<spval) &
928 & usl(i,j)=uh(i,j,ll)+(uh(i,j,ll)-uh(i,j,ll-1))*fact
929 IF(vh(i,j,ll)<spval .AND. vh(i,j,ll-1)<spval) &
930 & vsl(i,j)=vh(i,j,ll)+(vh(i,j,ll)-vh(i,j,ll-1))*fact
931!
932! FOR UNDERGROUND PRESSURE LEVELS, ASSUME TEMPERATURE TO CHANGE
933! ADIABATICLY, RH TO BE THE SAME AS THE AVERAGE OF THE 2ND AND 3RD
934! LAYERS FROM THE GOUND, WIND TO BE THE SAME AS THE LOWEST LEVEL ABOVE
935! GOUND
936 ELSE
937 IF(uh(i,j,llmh)<spval)usl(i,j)=uh(i,j,llmh)
938 IF(vh(i,j,llmh)<spval)vsl(i,j)=vh(i,j,llmh)
939 END IF
940 231 CONTINUE
941
942
943
944 END IF ! END OF WIND INTERPOLATION FOR NMM
945
946
947!
948!---------------------------------------------------------------------
949! LOAD GEOPOTENTIAL AND TEMPERATURE INTO STANDARD LEVEL
950! ARRAYS FOR THE NEXT PASS.
951!---------------------------------------------------------------------
952!
953!---------------------------------------------------------------------
954!*** CALCULATE 1000MB GEOPOTENTIALS CONSISTENT WITH SLP OBTAINED
955!*** FROM THE MESINGER OR NWS SHUELL SLP REDUCTION.
956!---------------------------------------------------------------------
957!
958!---------------------------------------------------------------------
959!---------------------------------------------------------------------
960! *** PART II ***
961!---------------------------------------------------------------------
962!---------------------------------------------------------------------
963!
964! INTERPOLATE/OUTPUT SELECTED FIELDS.
965!
966!---------------------------------------------------------------------
967!
968!*** OUTPUT GEOPOTENTIAL (SCALE BY GI)
969!
970 IF(iget(205)>0)THEN
971 IF(lvls(lp+1,iget(205))>0)THEN
972!$omp parallel do
973 DO j=jsta,jend
974 DO i=ista,iend
975 IF(fsl(i,j)<spval) THEN
976 grid1(i,j)=fsl(i,j)*gi
977 ELSE
978 grid1(i,j)=spval
979 ENDIF
980 ENDDO
981 ENDDO
982 if(grib=="grib2" )then
983 cfld=cfld+1
984 fld_info(cfld)%ifld=iavblfld(iget(205))
985 fld_info(cfld)%lvl=lvlsxml(lp+1,iget(205))
986 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
987 endif
988 ENDIF
989 ENDIF
990!
991!*** OUTPUT EXCHANGE COEEFICIENT
992!
993! OUTPUT FIRST INTERFACE KH Heat Diffusivity
994 IF (iget(243)>0) THEN !!Air Quality (Plee Oct2003) ^^^^^
995 IF (lvls(lp+1,iget(243))>0) THEN
996!$omp parallel do
997 DO j=jsta,jend
998 DO i=ista,iend
999 grid1(i,j)=akh(i,j)
1000 IF(lp==(lsig+1))grid1(i,j)=0.0 !! NO SLIP ASSUMTION FOR CMAQ
1001 ENDDO
1002 ENDDO
1003 if(grib=="grib2" )then
1004 cfld=cfld+1
1005 fld_info(cfld)%ifld=iavblfld(iget(243))
1006 fld_info(cfld)%lvl=lvlsxml(lp+1,iget(243))
1007 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1008 endif
1009 if(me==0)print*,'output Heat Diffusivity'
1010 ENDIF
1011 ENDIF
1012!
1013!*** TEMPERATURE
1014!
1015 IF(iget(206)>0) THEN
1016 IF(lvls(lp,iget(206))>0) THEN
1017 DO j=jsta,jend
1018 DO i=ista,iend
1019 grid1(i,j)=tsl(i,j)
1020 ENDDO
1021 ENDDO
1022 if(grib=="grib2" )then
1023 cfld=cfld+1
1024 fld_info(cfld)%ifld=iavblfld(iget(206))
1025 fld_info(cfld)%lvl=lvlsxml(lp,iget(206))
1026 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1027 endif
1028 ENDIF
1029 ENDIF
1030!
1031!*** PRESSURE
1032!
1033 IF(iget(216)>0)THEN
1034 IF(lvls(lp,iget(216))>0)THEN
1035!$omp parallel do
1036 DO j=jsta,jend
1037 DO i=ista,iend
1038 llmh = nint(lmh(i,j))
1039 grid1(i,j)=ptsigo+asigo(lp)*(pint(i,j,llmh+1)-ptsigo)
1040 ENDDO
1041 ENDDO
1042 if(grib=="grib2" )then
1043 cfld=cfld+1
1044 fld_info(cfld)%ifld=iavblfld(iget(216))
1045 fld_info(cfld)%lvl=lvlsxml(lp,iget(216))
1046 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1047 endif
1048 ENDIF
1049 ENDIF
1050!
1051!*** SPECIFIC HUMIDITY.
1052!
1053 IF(iget(207)>0)THEN
1054 IF(lvls(lp,iget(207))>0)THEN
1055 DO j=jsta,jend
1056 DO i=ista,iend
1057 grid1(i,j)=qsl(i,j)
1058 ENDDO
1059 ENDDO
1060 CALL bound(grid1,h1m12,h99999)
1061 if(grib=="grib2" )then
1062 cfld=cfld+1
1063 fld_info(cfld)%ifld=iavblfld(iget(207))
1064 fld_info(cfld)%lvl=lvlsxml(lp,iget(207))
1065 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1066 endif
1067 ENDIF
1068 ENDIF
1069!
1070!*** OMEGA
1071!
1072 IF(iget(210)>0)THEN
1073 IF(lvls(lp,iget(210))>0)THEN
1074 DO j=jsta,jend
1075 DO i=ista,iend
1076 grid1(i,j)=osl(i,j)
1077 ENDDO
1078 ENDDO
1079 if(grib=="grib2" )then
1080 cfld=cfld+1
1081 fld_info(cfld)%ifld=iavblfld(iget(210))
1082 fld_info(cfld)%lvl=lvlsxml(lp,iget(210))
1083 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1084 endif
1085 ENDIF
1086 ENDIF
1087!
1088!*** U AND/OR V WIND
1089!
1090 IF(iget(208)>0.OR.iget(209)>0)THEN
1091 IF(lvls(lp,iget(208))>0.OR.lvls(lp,iget(209))>0) then
1092 DO j=jsta,jend
1093 DO i=ista,iend
1094 grid1(i,j)=usl(i,j)
1095 grid2(i,j)=vsl(i,j)
1096 ENDDO
1097 ENDDO
1098 if(grib=="grib2" )then
1099 cfld=cfld+1
1100 fld_info(cfld)%ifld=iavblfld(iget(208))
1101 fld_info(cfld)%lvl=lvlsxml(lp,iget(208))
1102 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1103 cfld=cfld+1
1104 fld_info(cfld)%ifld=iavblfld(iget(209))
1105 fld_info(cfld)%lvl=lvlsxml(lp,iget(209))
1106 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid2(ista:iend,jsta:jend)
1107 endif
1108 ENDIF
1109 ENDIF
1110!
1111!*** TURBULENT KINETIC ENERGY
1112!
1113 IF (iget(217)>0) THEN
1114 IF (lvls(lp,iget(217))>0) THEN
1115 DO j=jsta,jend
1116 DO i=ista,iend
1117 grid1(i,j)=q2sl(i,j)
1118 ENDDO
1119 ENDDO
1120 if(grib=="grib2" )then
1121 cfld=cfld+1
1122 fld_info(cfld)%ifld=iavblfld(iget(217))
1123 fld_info(cfld)%lvl=lvlsxml(lp,iget(217))
1124 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1125 endif
1126 ENDIF
1127 ENDIF
1128!
1129!*** CLOUD WATER
1130!
1131 IF (iget(211)>0) THEN
1132 IF (lvls(lp,iget(211))>0) THEN
1133 DO j=jsta,jend
1134 DO i=ista,iend
1135 grid1(i,j)=qw1(i,j)
1136 ENDDO
1137 ENDDO
1138 if(grib=="grib2" )then
1139 cfld=cfld+1
1140 fld_info(cfld)%ifld=iavblfld(iget(211))
1141 fld_info(cfld)%lvl=lvlsxml(lp,iget(211))
1142 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1143 endif
1144 ENDIF
1145 ENDIF
1146!
1147!*** CLOUD ICE
1148!
1149 IF (iget(212)>0) THEN
1150 IF (lvls(lp,iget(212))>0) THEN
1151 DO j=jsta,jend
1152 DO i=ista,iend
1153 grid1(i,j)=qi1(i,j)
1154 ENDDO
1155 ENDDO
1156 if(grib=="grib2" )then
1157 cfld=cfld+1
1158 fld_info(cfld)%ifld=iavblfld(iget(212))
1159 fld_info(cfld)%lvl=lvlsxml(lp,iget(212))
1160 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1161 endif
1162 ENDIF
1163 ENDIF
1164!
1165!--- RAIN
1166 IF (iget(213)>0) THEN
1167 IF (lvls(lp,iget(213))>0) THEN
1168 DO j=jsta,jend
1169 DO i=ista,iend
1170 grid1(i,j)=qr1(i,j)
1171 ENDDO
1172 ENDDO
1173 if(grib=="grib2" )then
1174 cfld=cfld+1
1175 fld_info(cfld)%ifld=iavblfld(iget(213))
1176 fld_info(cfld)%lvl=lvlsxml(lp,iget(213))
1177 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1178 endif
1179 ENDIF
1180 ENDIF
1181!
1182!--- SNOW
1183 IF (iget(214)>0) THEN
1184 IF (lvls(lp,iget(214))>0) THEN
1185 DO j=jsta,jend
1186 DO i=ista,iend
1187 grid1(i,j)=qs1(i,j)
1188 ENDDO
1189 ENDDO
1190 if(grib=="grib2" )then
1191 cfld=cfld+1
1192 fld_info(cfld)%ifld=iavblfld(iget(214))
1193 fld_info(cfld)%lvl=lvlsxml(lp,iget(214))
1194 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1195 endif
1196 ENDIF
1197 ENDIF
1198!
1199!--- GRAUPEL
1200 IF (iget(255)>0) THEN
1201 IF (lvls(lp,iget(255))>0) THEN
1202 DO j=jsta,jend
1203 DO i=ista,iend
1204 grid1(i,j)=qg1(i,j)
1205 ENDDO
1206 ENDDO
1207 if(grib=="grib2" )then
1208 cfld=cfld+1
1209 fld_info(cfld)%ifld=iavblfld(iget(255))
1210 fld_info(cfld)%lvl=lvlsxml(lp,iget(255))
1211 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1212 endif
1213 ENDIF
1214 ENDIF
1215!
1216!--- TOTAL CONDENSATE
1217 IF (iget(215)>0) THEN
1218 IF (lvls(lp,iget(215))>0) THEN
1219 DO j=jsta,jend
1220 DO i=ista,iend
1221 grid1(i,j)=c1d(i,j)
1222 ENDDO
1223 ENDDO
1224 if(grib=="grib2" )then
1225 cfld=cfld+1
1226 fld_info(cfld)%ifld=iavblfld(iget(215))
1227 fld_info(cfld)%lvl=lvlsxml(lp,iget(215))
1228 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1229 endif
1230 ENDIF
1231 ENDIF
1232!
1233! TOTAL CLOUD COVER
1234 IF (iget(222)>0) THEN
1235 IF (lvls(lp,iget(222))>0) THEN
1236 DO j=jsta,jend
1237 DO i=ista,iend
1238 grid1(i,j)=cfrsig(i,j)
1239 ENDDO
1240 ENDDO
1241 if(grib=="grib2" )then
1242 cfld=cfld+1
1243 fld_info(cfld)%ifld=iavblfld(iget(222))
1244 fld_info(cfld)%lvl=lvlsxml(lp,iget(222))
1245 datapd(1:iend-ista+1,1:jend-jsta+1,cfld)=grid1(ista:iend,jsta:jend)
1246 endif
1247 ENDIF
1248 ENDIF
1249!*** END OF MAIN VERTICAL LOOP
1250!
1251 310 CONTINUE
1252!*** ENDIF FOR IF TEST SEEING IF WE WANT ANY OTHER VARIABLES
1253!
1254 ENDIF
1255!
1256!
1257!
1258! END OF ROUTINE.
1259!
1260 RETURN
1261 END
bound
subroutine bound(fld, fmin, fmax)
This routine bounds data in the passed array FLD (im x jm elements long) and clips data values such t...
Definition BOUND.f:35
mdl2sigma
subroutine mdl2sigma
SUBPROGRAM: MDL2P VERT INTRP OF MODEL LVLS TO PRESSURE PRGRMMR: BLACK ORG: W/NP22 DATE: 99-09-23
Definition MDL2SIGMA.f:54