UPP v11.0.0
Loading...
Searching...
No Matches
MDL2P.f
Go to the documentation of this file.
1
35 SUBROUTINE mdl2p(iostatusD3D)
36
37!
38!
39 use vrbls4d, only: dust, smoke
40 use vrbls3d, only: pint, o3, pmid, t, q, uh, vh, wh, omga, q2, cwm, &
41 qqw, qqi, qqr, qqs, qqg, dbz, f_rimef, ttnd, cfr, &
42 rlwtt, rswtt, vdifftt, tcucn, tcucns, &
43 train, vdiffmois, dconvmois, sconvmois,nradtt, &
44 o3vdiff, o3prod, o3tndy, mwpv, unknown, vdiffzacce, &
45 zgdrag, cnvctvmmixing, vdiffmacce, mgdrag, &
46 cnvctummixing, ncnvctcfrac, cnvctumflx, cnvctdetmflx, &
47 cnvctzgdrag, cnvctmgdrag, zmid, zint, pmidv, &
48 cnvctdmflx
49 use vrbls2d, only: t500,t700,w_up_max,w_dn_max,w_mean,pslp,fis,z1000,z700,&
50 z500
51 use masks, only: lmh, sm
52 use physcons_post,only: con_fvirt, con_rog, con_eps, con_epsm1
53 use params_mod, only: h1m12, dbzmin, h1, pq0, a2, a3, a4, rhmin, g, &
54 rgamog, rd, d608, gi, erad, pi, small, h100, &
55 h99999, gamma
56 use ctlblk_mod, only: modelname, lp1, me, jsta, jend, lm, spval, spl, &
57 alsl, jend_m, smflag, grib, cfld, fld_info, datapd,&
58 td3d, ifhr, ifmin, im, jm, nbin_du, jsta_2l, &
59 jend_2u, lsm, d3d_on, gocart_on, ioform, nbin_sm, &
60 imp_physics, ista, iend, ista_m, iend_m, ista_2l, iend_2u
61 use rqstfld_mod, only: iget, lvls, id, iavblfld, lvlsxml
62 use gridspec_mod, only: gridtype, maptype, dxval
63 use upp_physics, only: fpvsnew, calrh, calvor
64!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
65!
66 implicit none
67!
68! INCLUDE MODEL DIMENSIONS. SET/DERIVE OTHER PARAMETERS.
69! GAMMA AND RGAMOG ARE USED IN THE EXTRAPOLATION OF VIRTUAL
70! TEMPERATURES BEYOND THE UPPER OF LOWER LIMITS OF DATA.
71!
72!
73 real,parameter:: gammam=-1*gamma,zshul=75.,tvshul=290.66
74!
75! DECLARE VARIABLES.
76!
77 real,PARAMETER :: CAPA=0.28589641,p1000=1000.e2
78 LOGICAL IOOMG,IOALL
79 real, dimension(im,jm) :: GRID1, GRID2
80 real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u) :: FSL, TSL, QSL, OSL, USL, VSL &
81 &, Q2SL, WSL, CFRSL, O3SL, TDSL &
82 &, EGRID1, EGRID2 &
83 &, FSL_OLD, USL_OLD, VSL_OLD &
84 &, OSL_OLD, OSL995
85! REAL D3DSL(IM,JM,27),DUSTSL(IM,JM,NBIN_DU)
86 REAL, allocatable :: D3DSL(:,:,:), DUSTSL(:,:,:), SMOKESL(:,:,:)
87!
88 integer,intent(in) :: iostatusD3D
89 INTEGER, dimension(ista_2l:iend_2u,jsta_2l:jend_2u) :: NL1X, NL1XF
90 real, dimension(ISTA_2L:IEND_2U,JSTA_2L:JEND_2U,LSM) :: TPRS, QPRS, FPRS
91!
92 INTEGER K, NSMOOTH
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! QG1 - graupel mixing ratio
102! DBZ1 - radar reflectivity
103!
104 REAL, dimension(ista_2l:iend_2u,jsta_2l:jend_2u) :: C1D, QW1, QI1, QR1, QS1, QG1, DBZ1 &
105 , frime, rad, haines
106
107 REAL SDUMMY(IM,2)
108
109! SAVE RH, U,V, for Icing, CAT, LLWS computation
110 REAL SAVRH(ista:iend,jsta:jend)
111!jw
112 integer I,J,L,LP,LL,LLMH,JJB,JJE,II,JJ,LI,IFINCR,ITD3D,istaa,imois,luhi,la
113 real fact,ALPSL,PSFC,QBLO,PNL1,TBLO,TVRL,TVRBLO,FAC,PSLPIJ, &
114 alpth,ahf,pdv,ql,tvu,tvd,gammas,qsat,rhl,zl,tl,pl,es,part,dum1
115 logical log1
116 real dxm, tem, zero
117!
118!******************************************************************************
119!
120! START MDL2P.
121!
122 if(me==0) print*, 'MDL2P SMFLAG=',smflag
123
124 if (modelname == 'GFS') then
125 zero = 0.0
126 else
127 zero = h1m12
128 endif
129 if (d3d_on) then
130 if (.not. allocated(d3dsl)) allocate(d3dsl(im,jm,27))
131!$omp parallel do private(i,j,l)
132 do l=1,27
133 do j=1,jm
134 do i=1,im
135 d3dsl(i,j,l) = spval
136 enddo
137 enddo
138 enddo
139 endif
140 if (gocart_on) then
141 if (.not. allocated(dustsl)) allocate(dustsl(im,jm,nbin_du))
142!$omp parallel do private(i,j,l)
143 do l=1,nbin_du
144 do j=1,jm
145 do i=1,im
146 dustsl(i,j,l) = spval
147 enddo
148 enddo
149 enddo
150 endif
151 if (.not. allocated(smokesl)) allocate(smokesl(im,jm,nbin_sm))
152!$omp parallel do private(i,j,l)
153 do l=1,nbin_sm
154 do j=1,jm
155 do i=1,im
156 smokesl(i,j,l) = spval
157 enddo
158 enddo
159 enddo
160!
161! SET TOTAL NUMBER OF POINTS ON OUTPUT GRID.
162!
163!---------------------------------------------------------------
164!
165! *** PART I ***
166!
167! VERTICAL INTERPOLATION OF EVERYTHING ELSE. EXECUTE ONLY
168! IF THERE'S SOMETHING WE WANT.
169!
170 IF((iget(012) > 0) .OR. (iget(013) > 0) .OR. &
171 (iget(014) > 0) .OR. (iget(015) > 0) .OR. &
172 (iget(016) > 0) .OR. (iget(017) > 0) .OR. &
173 (iget(018) > 0) .OR. (iget(019) > 0) .OR. &
174 (iget(020) > 0) .OR. (iget(030) > 0) .OR. &
175 (iget(021) > 0) .OR. (iget(022) > 0) .OR. &
176 (iget(023) > 0) .OR. (iget(085) > 0) .OR. &
177 (iget(086) > 0) .OR. (iget(284) > 0) .OR. &
178 (iget(153) > 0) .OR. (iget(166) > 0) .OR. &
179 (iget(183) > 0) .OR. (iget(184) > 0) .OR. &
180 (iget(198) > 0) .OR. (iget(251) > 0) .OR. &
181 (iget(257) > 0) .OR. (iget(258) > 0) .OR. &
182 (iget(294) > 0) .OR. (iget(268) > 0) .OR. &
183 (iget(331) > 0) .OR. (iget(326) > 0) .OR. &
184! add D3D fields
185 (iget(354) > 0) .OR. (iget(355) > 0) .OR. &
186 (iget(356) > 0) .OR. (iget(357) > 0) .OR. &
187 (iget(358) > 0) .OR. (iget(359) > 0) .OR. &
188 (iget(360) > 0) .OR. (iget(361) > 0) .OR. &
189 (iget(362) > 0) .OR. (iget(363) > 0) .OR. &
190 (iget(364) > 0) .OR. (iget(365) > 0) .OR. &
191 (iget(366) > 0) .OR. (iget(367) > 0) .OR. &
192 (iget(368) > 0) .OR. (iget(369) > 0) .OR. &
193 (iget(370) > 0) .OR. (iget(371) > 0) .OR. &
194 (iget(372) > 0) .OR. (iget(373) > 0) .OR. &
195 (iget(374) > 0) .OR. (iget(375) > 0) .OR. &
196 (iget(391) > 0) .OR. (iget(392) > 0) .OR. &
197 (iget(393) > 0) .OR. (iget(394) > 0) .OR. &
198 (iget(395) > 0) .OR. (iget(379) > 0) .OR. &
199! ADD DUST FIELDS
200 (iget(438) > 0) .OR. (iget(439) > 0) .OR. &
201 (iget(440) > 0) .OR. (iget(441) > 0) .OR. &
202 (iget(442) > 0) .OR. (iget(455) > 0) .OR. &
203! ADD SMOKE FIELDS
204 (iget(738) > 0) .OR. (modelname == 'RAPR') .OR.&
205! LIFTED INDEX needs 500 mb T
206 (iget(030)>0) .OR. (iget(031)>0) .OR. (iget(075)>0)) THEN
207!
208!---------------------------------------------------------------------
209!***
210!*** BECAUSE SIGMA LAYERS DO NOT GO UNDERGROUND, DO ALL
211!*** INTERPOLATION ABOVE GROUND NOW.
212!***
213!
214! print*,'LSM= ',lsm
215
216 if(gridtype == 'B' .or. gridtype == 'E') &
217 call exch(pint(ista_2l:iend_2u,jsta_2l:jend_2u,lp1))
218
219 DO lp=1,lsm
220
221! if(me == 0) print *,'in LP loop me=',me,'UH=',UH(1:10,JSTA,LP), &
222! 'JSTA_2L=',JSTA_2L,'JEND_2U=',JEND_2U,'JSTA=',JSTA,JEND, &
223! 'PMID(1,1,L)=',(PMID(1,1,LI),LI=1,LM),'SPL(LP)=',SPL(LP)
224
225! if(me ==0) print *,'in mdl2p,LP loop o3=',maxval(o3(1:im,jsta:jend,lm))
226!
227!$omp parallel do private(i,j,l)
228 DO j=jsta_2l,jend_2u
229 DO i=ista_2l,iend_2u
230 tsl(i,j) = spval
231 qsl(i,j) = spval
232 fsl(i,j) = spval
233 osl(i,j) = spval
234 usl(i,j) = spval
235 vsl(i,j) = spval
236! dong initialize wsl
237 wsl(i,j) = spval
238 q2sl(i,j) = spval
239 c1d(i,j) = spval ! Total condensate
240 qw1(i,j) = spval ! Cloud water
241 qi1(i,j) = spval ! Cloud ice
242 qr1(i,j) = spval ! Rain
243 qs1(i,j) = spval ! Snow (precip ice)
244 qg1(i,j) = spval ! Graupel
245 dbz1(i,j) = spval
246 frime(i,j) = spval
247 rad(i,j) = spval
248 o3sl(i,j) = spval
249 cfrsl(i,j) = spval
250!
251!*** LOCATE VERTICAL INDEX OF MODEL MIDLAYER JUST BELOW
252!*** THE PRESSURE LEVEL TO WHICH WE ARE INTERPOLATING.
253!
254 nl1x(i,j) = lp1
255 DO l=2,lm
256 IF(nl1x(i,j) == lp1 .AND. pmid(i,j,l) > spl(lp)) THEN
257 nl1x(i,j) = l
258 ENDIF
259 ENDDO
260!
261! IF THE PRESSURE LEVEL IS BELOW THE LOWEST MODEL MIDLAYER
262! BUT STILL ABOVE THE LOWEST MODEL BOTTOM INTERFACE,
263! WE WILL NOT CONSIDER IT UNDERGROUND AND THE INTERPOLATION
264! WILL EXTRAPOLATE TO THAT POINT
265!
266 IF(nl1x(i,j) == lp1 .AND. pint(i,j,lp1) > spl(lp)) THEN
267 nl1x(i,j) = lm
268 ENDIF
269
270 nl1xf(i,j) = lp1 + 1
271 DO l=2,lp1
272 IF(nl1xf(i,j) == (lp1+1) .AND. pint(i,j,l) > spl(lp)) THEN
273 nl1xf(i,j) = l
274 ENDIF
275 ENDDO
276! if(NL1X(I,J) == LMP1)print*,'Debug: NL1X=LMP1 AT ' 1 ,i,j,lp
277 ENDDO
278 ENDDO ! end of j loop
279
280!
281!mptest IF(NHOLD == 0)GO TO 310
282!
283!!$omp parallel do
284!!$omp& private(nn,i,j,ll,fact,qsat,rhl)
285!hc DO 220 NN=1,NHOLD
286!hc I=IHOLD(NN)
287!hc J=JHOLD(NN)
288! DO 220 J=JSTA,JEND
289
290 ii = (ista+iend)/2
291 jj = (jsta+jend)/2
292
293!$omp parallel do private(i,j,k,l,ll,llmh,la,tvd,tvu,fact,fac,ahf,rhl,tl,pl,ql,zl,es,qsat,part,tvrl,tvrblo,tblo,qblo,gammas,pnl1)
294 DO j=jsta,jend
295 DO i=ista,iend
296!---------------------------------------------------------------------
297!*** VERTICAL INTERPOLATION OF GEOPOTENTIAL, TEMPERATURE, SPECIFIC
298!*** HUMIDITY, CLOUD WATER/ICE, OMEGA, WINDS, AND TKE.
299!---------------------------------------------------------------------
300!
301 ll = nl1x(i,j)
302 llmh = nint(lmh(i,j))
303
304!HC IF(NL1X(I,J)<=LM)THEN
305
306 IF(spl(lp) < pint(i,j,2)) THEN ! Above second interface
307 IF(t(i,j,1) < spval) tsl(i,j) = t(i,j,1)
308 IF(q(i,j,1) < spval) qsl(i,j) = q(i,j,1)
309
310 IF(gridtype == 'A')THEN
311 IF(uh(i,j,1) < spval) usl(i,j) = uh(i,j,1)
312 IF(vh(i,j,1) < spval) vsl(i,j) = vh(i,j,1)
313 END IF
314
315! if ( J == JSTA.and. I == 1.and.me == 0) &
316! print *,'1 USL=',USL(I,J),UH(I,J,1)
317
318 IF(wh(i,j,1) < spval) wsl(i,j) = wh(i,j,1)
319 IF(omga(i,j,1) < spval) osl(i,j) = omga(i,j,1)
320 IF(q2(i,j,1) < spval) q2sl(i,j) = q2(i,j,1)
321 IF(cwm(i,j,1) < spval) c1d(i,j) = cwm(i,j,1)
322 c1d(i,j) = max(c1d(i,j),zero) ! Total condensate
323 IF(qqw(i,j,1) < spval) qw1(i,j) = qqw(i,j,1)
324 qw1(i,j) = max(qw1(i,j),zero) ! Cloud water
325 IF(qqi(i,j,1) < spval) qi1(i,j) = qqi(i,j,1)
326 qi1(i,j) = max(qi1(i,j),zero) ! Cloud ice
327 IF(qqr(i,j,1) < spval) qr1(i,j) = qqr(i,j,1)
328 qr1(i,j) = max(qr1(i,j),zero) ! Rain
329 IF(qqs(i,j,1) < spval) qs1(i,j) = qqs(i,j,1)
330 qs1(i,j) = max(qs1(i,j),zero) ! Snow (precip ice)
331 IF(qqg(i,j,1) < spval) qg1(i,j) = qqg(i,j,1)
332 qg1(i,j) = max(qg1(i,j),zero) ! Graupel (precip ice)
333 IF(dbz(i,j,1) < spval) dbz1(i,j) = dbz(i,j,1)
334 dbz1(i,j) = max(dbz1(i,j),dbzmin)
335 IF(f_rimef(i,j,1) < spval) frime(i,j) = f_rimef(i,j,1)
336 frime(i,j) = max(frime(i,j),h1)
337 IF(ttnd(i,j,1) < spval) rad(i,j) = ttnd(i,j,1)
338 IF(o3(i,j,1) < spval) o3sl(i,j) = o3(i,j,1)
339 IF(cfr(i,j,1) < spval) cfrsl(i,j) = cfr(i,j,1)
340! DUST
341 if (gocart_on) then
342 DO k = 1, nbin_du
343 IF(dust(i,j,1,k) < spval) dustsl(i,j,k) = dust(i,j,1,k)
344 ENDDO
345 endif
346 DO k = 1, nbin_sm
347 IF(smoke(i,j,1,k) < spval) smokesl(i,j,k)=smoke(i,j,1,k)
348 ENDDO
349
350! only interpolate GFS d3d fields when reqested
351! if(iostatusD3D ==0 .and. d3d_on)then
352 if (d3d_on) then
353 IF((iget(354) > 0) .OR. (iget(355) > 0) .OR. &
354 (iget(356) > 0) .OR. (iget(357) > 0) .OR. &
355 (iget(358) > 0) .OR. (iget(359) > 0) .OR. &
356 (iget(360) > 0) .OR. (iget(361) > 0) .OR. &
357 (iget(362) > 0) .OR. (iget(363) > 0) .OR. &
358 (iget(364) > 0) .OR. (iget(365) > 0) .OR. &
359 (iget(366) > 0) .OR. (iget(367) > 0) .OR. &
360 (iget(368) > 0) .OR. (iget(369) > 0) .OR. &
361 (iget(370) > 0) .OR. (iget(371) > 0) .OR. &
362 (iget(372) > 0) .OR. (iget(373) > 0) .OR. &
363 (iget(374) > 0) .OR. (iget(375) > 0) .OR. &
364 (iget(391) > 0) .OR. (iget(392) > 0) .OR. &
365 (iget(393) > 0) .OR. (iget(394) > 0) .OR. &
366 (iget(395) > 0) .OR. (iget(379) > 0)) THEN
367 d3dsl(i,j,1) = rlwtt(i,j,1)
368 d3dsl(i,j,2) = rswtt(i,j,1)
369 d3dsl(i,j,3) = vdifftt(i,j,1)
370 d3dsl(i,j,4) = tcucn(i,j,1)
371 d3dsl(i,j,5) = tcucns(i,j,1)
372 d3dsl(i,j,6) = train(i,j,1)
373 d3dsl(i,j,7) = vdiffmois(i,j,1)
374 d3dsl(i,j,8) = dconvmois(i,j,1)
375 d3dsl(i,j,9) = sconvmois(i,j,1)
376 d3dsl(i,j,10) = nradtt(i,j,1)
377 d3dsl(i,j,11) = o3vdiff(i,j,1)
378 d3dsl(i,j,12) = o3prod(i,j,1)
379 d3dsl(i,j,13) = o3tndy(i,j,1)
380 d3dsl(i,j,14) = mwpv(i,j,1)
381 d3dsl(i,j,15) = unknown(i,j,1)
382 d3dsl(i,j,16) = vdiffzacce(i,j,1)
383 d3dsl(i,j,17) = zgdrag(i,j,1)
384 d3dsl(i,j,18) = cnvctummixing(i,j,1)
385 d3dsl(i,j,19) = vdiffmacce(i,j,1)
386 d3dsl(i,j,20) = mgdrag(i,j,1)
387 d3dsl(i,j,21) = cnvctvmmixing(i,j,1)
388 d3dsl(i,j,22) = ncnvctcfrac(i,j,1)
389 d3dsl(i,j,23) = cnvctumflx(i,j,1)
390 d3dsl(i,j,24) = cnvctdmflx(i,j,1)
391 d3dsl(i,j,25) = cnvctdetmflx(i,j,1)
392 d3dsl(i,j,26) = cnvctzgdrag(i,j,1)
393 d3dsl(i,j,27) = cnvctmgdrag(i,j,1)
394 end if
395 end if
396
397 ELSE IF(ll <= llmh)THEN
398!
399!---------------------------------------------------------------------
400! INTERPOLATE LINEARLY IN LOG(P)
401!*** EXTRAPOLATE ABOVE THE TOPMOST MIDLAYER OF THE MODEL
402!*** INTERPOLATION BETWEEN NORMAL LOWER AND UPPER BOUNDS
403!*** EXTRAPOLATE BELOW LOWEST MODEL MIDLAYER (BUT STILL ABOVE GROUND)
404!---------------------------------------------------------------------
405!
406!KRF: Need ncar and nmm wrf core checks as well?
407 IF (modelname == 'RAPR' .OR. modelname == 'NCAR' .OR. modelname == 'NMM') THEN
408 fact = (alsl(lp)-log(pmid(i,j,ll)))/ &
409 max(1.e-6,(log(pmid(i,j,ll))-log(pmid(i,j,ll-1))))
410 fact = max(-10.0,min(fact, 10.0))
411 ELSEIF (modelname == 'GFS' .OR. modelname == 'FV3R') THEN
412 fact = (alsl(lp)-log(pmid(i,j,ll)))/ &
413 max(1.e-6,(log(pmid(i,j,ll))-log(pmid(i,j,ll-1))))
414 fact = max(-10.0,min(fact, 10.0))
415 IF ( abs(pmid(i,j,ll)-pmid(i,j,ll-1)) < 0.5 ) THEN
416 fact = -1.
417 ENDIF
418 ELSE
419 fact = (alsl(lp)-log(pmid(i,j,ll)))/ &
420 (log(pmid(i,j,ll))-log(pmid(i,j,ll-1)))
421 ENDIF
422 IF(t(i,j,ll) < spval .AND. t(i,j,ll-1) < spval) &
423 tsl(i,j) = t(i,j,ll)+(t(i,j,ll)-t(i,j,ll-1))*fact
424 IF(q(i,j,ll) < spval .AND. q(i,j,ll-1) < spval) &
425 qsl(i,j) = q(i,j,ll)+(q(i,j,ll)-q(i,j,ll-1))*fact
426
427 IF(gridtype=='A')THEN
428 IF(uh(i,j,ll) < spval .AND. uh(i,j,ll-1) < spval) &
429 usl(i,j) = uh(i,j,ll)+(uh(i,j,ll)-uh(i,j,ll-1))*fact
430 IF(vh(i,j,ll) < spval .AND. vh(i,j,ll-1) < spval) &
431 vsl(i,j) = vh(i,j,ll)+(vh(i,j,ll)-vh(i,j,ll-1))*fact
432 END IF
433
434 IF(wh(i,j,ll) < spval .AND. wh(i,j,ll-1) < spval) &
435 wsl(i,j) = wh(i,j,ll)+(wh(i,j,ll)-wh(i,j,ll-1))*fact
436 IF(omga(i,j,ll) < spval .AND. omga(i,j,ll-1) < spval) &
437 osl(i,j) = omga(i,j,ll)+(omga(i,j,ll)-omga(i,j,ll-1))*fact
438 IF(q2(i,j,ll) < spval .AND. q2(i,j,ll-1) < spval) &
439 q2sl(i,j) = q2(i,j,ll)+(q2(i,j,ll)-q2(i,j,ll-1))*fact
440
441! IF(ZMID(I,J,LL) < SPVAL .AND. ZMID(I,J,LL-1) < SPVAL) &
442! & FSL(I,J) = ZMID(I,J,LL)+(ZMID(I,J,LL)-ZMID(I,J,LL-1))*FACT
443! FSL(I,J) = FSL(I,J)*G
444
445 if (modelname == 'GFS') then
446 es = min(fpvsnew(tsl(i,j)), spl(lp))
447 qsat = con_eps*es/(spl(lp)+con_epsm1*es)
448 else
449 qsat = pq0/spl(lp)*exp(a2*(tsl(i,j)-a3)/(tsl(i,j)-a4))
450 endif
451!
452 rhl = max(rhmin, min(1.0, qsl(i,j)/qsat))
453 qsl(i,j) = rhl*qsat
454
455! if(tsl(i,j) > 330. .or. tsl(i,j) < 100.)print*, &
456! 'bad isobaric T Q',i,j,lp,tsl(i,j),qsl(i,j) &
457! ,T(I,J,LL),T(I,J,LL-1),Q(I,J,LL),Q(I,J,LL-1)
458
459 IF(q2sl(i,j) < 0.0) q2sl(i,j) = 0.0
460!
461!HC ADD FERRIER'S HYDROMETEOR
462 IF(cwm(i,j,ll) < spval .AND. cwm(i,j,ll-1) < spval) &
463 c1d(i,j) = cwm(i,j,ll) + (cwm(i,j,ll)-cwm(i,j,ll-1))*fact
464 c1d(i,j) = max(c1d(i,j),zero) ! Total condensate
465
466 IF(qqw(i,j,ll) < spval .AND. qqw(i,j,ll-1) < spval) &
467 qw1(i,j) = qqw(i,j,ll) + (qqw(i,j,ll)-qqw(i,j,ll-1))*fact
468 qw1(i,j) = max(qw1(i,j),zero) ! Cloud water
469
470 IF(qqi(i,j,ll) < spval .AND. qqi(i,j,ll-1) < spval) &
471 qi1(i,j) = qqi(i,j,ll) + (qqi(i,j,ll)-qqi(i,j,ll-1))*fact
472 qi1(i,j) = max(qi1(i,j),zero) ! Cloud ice
473
474 IF(qqr(i,j,ll) < spval .AND. qqr(i,j,ll-1) < spval) &
475 qr1(i,j) = qqr(i,j,ll) + (qqr(i,j,ll)-qqr(i,j,ll-1))*fact
476 qr1(i,j) = max(qr1(i,j),zero) ! Rain
477
478 IF(qqs(i,j,ll) < spval .AND. qqs(i,j,ll-1) < spval) &
479 qs1(i,j) = qqs(i,j,ll) + (qqs(i,j,ll)-qqs(i,j,ll-1))*fact
480 qs1(i,j) = max(qs1(i,j),zero) ! Snow (precip ice)
481
482 IF(qqg(i,j,ll) < spval .AND. qqg(i,j,ll-1) < spval) &
483 qg1(i,j) = qqg(i,j,ll) + (qqg(i,j,ll)-qqg(i,j,ll-1))*fact
484 qg1(i,j) = max(qg1(i,j),zero) ! GRAUPEL (precip ice)
485
486 IF(dbz(i,j,ll) < spval .AND. dbz(i,j,ll-1) < spval) &
487 dbz1(i,j) = dbz(i,j,ll) + (dbz(i,j,ll)-dbz(i,j,ll-1))*fact
488 dbz1(i,j) = max(dbz1(i,j),dbzmin)
489
490 IF(f_rimef(i,j,ll) < spval .AND. f_rimef(i,j,ll-1) < spval) &
491 frime(i,j) = f_rimef(i,j,ll) + (f_rimef(i,j,ll) - f_rimef(i,j,ll-1))*fact
492 frime(i,j)=max(frime(i,j),h1)
493
494 IF(ttnd(i,j,ll) < spval .AND. ttnd(i,j,ll-1) < spval) &
495 rad(i,j) = ttnd(i,j,ll) + (ttnd(i,j,ll)-ttnd(i,j,ll-1))*fact
496
497 IF(o3(i,j,ll) < spval .AND. o3(i,j,ll-1) < spval) &
498 o3sl(i,j) = o3(i,j,ll) + (o3(i,j,ll)-o3(i,j,ll-1))*fact
499
500 IF(cfr(i,j,ll) < spval .AND. cfr(i,j,ll-1) < spval) &
501 cfrsl(i,j) = cfr(i,j,ll) + (cfr(i,j,ll)-cfr(i,j,ll-1))*fact
502! DUST
503 if (gocart_on) then
504 DO k = 1, nbin_du
505 IF(dust(i,j,ll,k) < spval .AND. dust(i,j,ll-1,k) < spval) &
506 dustsl(i,j,k) = dust(i,j,ll,k) + (dust(i,j,ll,k)-dust(i,j,ll-1,k))*fact
507 ENDDO
508 endif
509 DO k = 1, nbin_sm
510 IF(smoke(i,j,ll,k) < spval .AND. smoke(i,j,ll-1,k) < spval) &
511 smokesl(i,j,k)=smoke(i,j,ll,k)+(smoke(i,j,ll,k)-smoke(i,j,ll-1,k))*fact
512 ENDDO
513
514! only interpolate GFS d3d fields when == ested
515! if(iostatusD3D==0)then
516 if (d3d_on) then
517 IF((iget(354) > 0) .OR. (iget(355) > 0) .OR. &
518 (iget(356) > 0) .OR. (iget(357) > 0) .OR. &
519 (iget(358) > 0) .OR. (iget(359) > 0) .OR. &
520 (iget(360) > 0) .OR. (iget(361) > 0) .OR. &
521 (iget(362) > 0) .OR. (iget(363) > 0) .OR. &
522 (iget(364) > 0) .OR. (iget(365) > 0) .OR. &
523 (iget(366) > 0) .OR. (iget(367) > 0) .OR. &
524 (iget(368) > 0) .OR. (iget(369) > 0) .OR. &
525 (iget(370) > 0) .OR. (iget(371) > 0) .OR. &
526 (iget(372) > 0) .OR. (iget(373) > 0) .OR. &
527 (iget(374) > 0) .OR. (iget(375) > 0) .OR. &
528 (iget(391) > 0) .OR. (iget(392) > 0) .OR. &
529 (iget(393) > 0) .OR. (iget(394) > 0) .OR. &
530 (iget(395) > 0) .OR. (iget(379) > 0))THEN
531 d3dsl(i,j,1) = rlwtt(i,j,ll)+(rlwtt(i,j,ll) &
532 - rlwtt(i,j,ll-1))*fact
533 d3dsl(i,j,2) = rswtt(i,j,ll)+(rswtt(i,j,ll) &
534 - rswtt(i,j,ll-1))*fact
535 d3dsl(i,j,3) = vdifftt(i,j,ll)+(vdifftt(i,j,ll) &
536 - vdifftt(i,j,ll-1))*fact
537 d3dsl(i,j,4) = tcucn(i,j,ll)+(tcucn(i,j,ll) &
538 - tcucn(i,j,ll-1))*fact
539 d3dsl(i,j,5) = tcucns(i,j,ll)+(tcucns(i,j,ll) &
540 - tcucns(i,j,ll-1))*fact
541 d3dsl(i,j,6) = train(i,j,ll)+(train(i,j,ll) &
542 - train(i,j,ll-1))*fact
543 d3dsl(i,j,7) = vdiffmois(i,j,ll)+ &
544 (vdiffmois(i,j,ll)-vdiffmois(i,j,ll-1))*fact
545 d3dsl(i,j,8) = dconvmois(i,j,ll)+ &
546 (dconvmois(i,j,ll)-dconvmois(i,j,ll-1))*fact
547 d3dsl(i,j,9) = sconvmois(i,j,ll)+ &
548 (sconvmois(i,j,ll)-sconvmois(i,j,ll-1))*fact
549 d3dsl(i,j,10) = nradtt(i,j,ll)+ &
550 (nradtt(i,j,ll)-nradtt(i,j,ll-1))*fact
551 d3dsl(i,j,11) = o3vdiff(i,j,ll)+ &
552 (o3vdiff(i,j,ll)-o3vdiff(i,j,ll-1))*fact
553 d3dsl(i,j,12) = o3prod(i,j,ll)+ &
554 (o3prod(i,j,ll)-o3prod(i,j,ll-1))*fact
555 d3dsl(i,j,13) = o3tndy(i,j,ll)+ &
556 (o3tndy(i,j,ll)-o3tndy(i,j,ll-1))*fact
557 d3dsl(i,j,14) = mwpv(i,j,ll)+ &
558 (mwpv(i,j,ll)-mwpv(i,j,ll-1))*fact
559 d3dsl(i,j,15) = unknown(i,j,ll)+ &
560 (unknown(i,j,ll)-unknown(i,j,ll-1))*fact
561 d3dsl(i,j,16) = vdiffzacce(i,j,ll)+ &
562 (vdiffzacce(i,j,ll)-vdiffzacce(i,j,ll-1))*fact
563 d3dsl(i,j,17) = zgdrag(i,j,ll)+ &
564 (zgdrag(i,j,ll)-zgdrag(i,j,ll-1))*fact
565 d3dsl(i,j,18) = cnvctummixing(i,j,ll)+ &
566 (cnvctummixing(i,j,ll)-cnvctummixing(i,j,ll-1))*fact
567 d3dsl(i,j,19) = vdiffmacce(i,j,ll)+ &
568 (vdiffmacce(i,j,ll)-vdiffmacce(i,j,ll-1))*fact
569 d3dsl(i,j,20) = mgdrag(i,j,ll)+ &
570 (mgdrag(i,j,ll)-mgdrag(i,j,ll-1))*fact
571 d3dsl(i,j,21) = cnvctvmmixing(i,j,ll)+ &
572 (cnvctvmmixing(i,j,ll)-cnvctvmmixing(i,j,ll-1))*fact
573 d3dsl(i,j,22) = ncnvctcfrac(i,j,ll)+ &
574 (ncnvctcfrac(i,j,ll)-ncnvctcfrac(i,j,ll-1))*fact
575 d3dsl(i,j,23) = cnvctumflx(i,j,ll)+ &
576 (cnvctumflx(i,j,ll)-cnvctumflx(i,j,ll-1))*fact
577 d3dsl(i,j,24) = cnvctdmflx(i,j,ll)+ &
578 (cnvctdmflx(i,j,ll)-cnvctdmflx(i,j,ll-1))*fact
579 d3dsl(i,j,25) = cnvctdetmflx(i,j,ll)+ &
580 (cnvctdetmflx(i,j,ll)-cnvctdetmflx(i,j,ll-1))*fact
581 d3dsl(i,j,26) = cnvctzgdrag(i,j,ll)+ &
582 (cnvctzgdrag(i,j,ll)-cnvctzgdrag(i,j,ll-1))*fact
583 d3dsl(i,j,27) = cnvctmgdrag(i,j,ll)+ &
584 (cnvctmgdrag(i,j,ll)-cnvctmgdrag(i,j,ll-1))*fact
585 end if
586 end if ! if d3d_on test
587
588! FOR UNDERGROUND PRESSURE LEVELS, ASSUME TEMPERATURE TO CHANGE
589! ADIABATICLY, RH TO BE THE SAME AS THE AVERAGE OF THE 2ND AND 3RD
590! LAYERS FROM THE GOUND, WIND TO BE THE SAME AS THE LOWEST LEVEL ABOVE
591! GOUND
592 ELSE ! underground
593 IF(modelname == 'GFS')THEN ! GFS deduce T and H using Shuell
594 tvu = t(i,j,lm) * (1.+con_fvirt*q(i,j,lm))
595 if(zmid(i,j,lm) > zshul) then
596 tvd = tvu + gamma*zmid(i,j,lm)
597 if(tvd > tvshul) then
598 if(tvu > tvshul) then
599 tvd = tvshul - 5.e-3*(tvu-tvshul)*(tvu-tvshul)
600 else
601 tvd = tvshul
602 endif
603 endif
604 gammas = (tvu-tvd)/zmid(i,j,lm)
605 else
606 gammas = 0.
607 endif
608 part = con_rog*(alsl(lp)-log(pmid(i,j,lm)))
609 fsl(i,j) = zmid(i,j,lm) - tvu*part/(1.+0.5*gammas*part)
610! tp(k) = t(1)+gammas*(hp(k)-h(1))
611 tsl(i,j) = t(i,j,lm) - gamma*(fsl(i,j)-zmid(i,j,lm))
612 fsl(i,j) = fsl(i,j)*g ! just use NAM G for now since FSL will be divided by GI later
613!
614! Compute RH at lowest model layer because Iredell and Chuang decided to compute
615! underground GFS Q to maintain RH
616 es = min(fpvsnew(t(i,j,lm)), pmid(i,j,lm))
617 qsat = con_eps*es/(pmid(i,j,lm)+con_epsm1*es)
618 rhl = q(i,j,lm)/qsat
619! compute saturation water vapor at isobaric level
620 es = min(fpvsnew(tsl(i,j)), spl(lp))
621 qsat = con_eps*es/(spl(lp)+con_epsm1*es)
622! Q at isobaric level is computed by maintaining constant RH
623 qsl(i,j) = rhl*qsat
624
625 ELSE
626 pl = pint(i,j,lm-1)
627 zl = zint(i,j,lm-1)
628 tl = 0.5*(t(i,j,lm-2)+t(i,j,lm-1))
629 ql = 0.5*(q(i,j,lm-2)+q(i,j,lm-1))
630! TMT0=TL-A0
631! TMT15=MIN(TMT0,-15.)
632! AI=0.008855
633! BI=1.
634! IF(TMT0 < -20.)THEN
635! AI=0.007225
636! BI=0.9674
637! ENDIF
638
639 qsat = pq0/pl*exp(a2*(tl-a3)/(tl-a4))
640 rhl = ql/qsat
641!
642 IF(rhl > 1.)THEN
643 rhl = 1.
644 ql = rhl*qsat
645 ENDIF
646!
647 IF(rhl < rhmin)THEN
648 rhl = rhmin
649 ql = rhl*qsat
650 ENDIF
651!
652 tvrl = tl*(1.+0.608*ql)
653 tvrblo = tvrl*(spl(lp)/pl)**rgamog
654 tblo = tvrblo/(1.+0.608*ql)
655!
656! TMT0=TBLO-A3
657! TMT15=MIN(TMT0,-15.)
658! AI=0.008855
659! BI=1.
660! IF(TMT0 < -20.)THEN
661! AI=0.007225
662! BI=0.9674
663! ENDIF
664
665 qsat = pq0/spl(lp)*exp(a2*(tblo-a3)/(tblo-a4))
666 tsl(i,j) = tblo
667 qblo = rhl*qsat
668 qsl(i,j) = max(1.e-12,qblo)
669 END IF ! endif loop for deducing T and H differently for GFS
670
671! if(tsl(i,j) > 330. .or. tsl(i,j) < 100.)print*, &
672! 'bad isobaric T Q',i,j,lp,tsl(i,j),qsl(i,j),tl,ql,pl
673
674 IF(gridtype == 'A')THEN
675 usl(i,j) = uh(i,j,llmh)
676 vsl(i,j) = vh(i,j,llmh)
677 END IF
678! if ( J == JSTA.and. I == 1.and.me == 0) &
679! & print *,'3 USL=',USL(I,J),UH(I,J,LLMH),LLMH
680 wsl(i,j) = wh(i,j,llmh)
681 osl(i,j) = omga(i,j,llmh)
682 q2sl(i,j) = max(0.0,0.5*(q2(i,j,llmh-1)+q2(i,j,llmh)))
683 pnl1 = pint(i,j,ll)
684 fac = 0.0
685 ahf = 0.0
686
687! FSL(I,J)=(PNL1-SPL(LP))/(SPL(LP)+PNL1)
688! 1 *(TSL(I,J))*(1.+0.608*QSL(I,J))
689! 2 *RD*2.+ZINT(I,J,NL1X(I,J))*G
690
691! FSL(I,J)=FPRS(I,J,LP-1)-RD*(TPRS(I,J,LP-1)
692! 1 *(H1+D608*QPRS(I,J,LP-1))
693! 2 +TSL(I,J)*(H1+D608*QSL(I,J)))
694! 3 *LOG(SPL(LP)/SPL(LP-1))/2.0
695
696! if(abs(SPL(LP)-97500.0) < 0.01)then
697! if(gdlat(i,j) > 35.0.and.gdlat(i,j)<=37.0 .and. &
698! gdlon(i,j) > -100.0 .and. gdlon(i,j) < -96.0)print*, &
699! 'Debug:I,J,FPRS(LP-1),TPRS(LP-1),TSL,SPL(LP),SPL(LP-1)=' &
700! ,i,j,FPRS(I,J,LP-1),TPRS(I,J,LP-1),TSL(I,J),SPL(LP) &
701! ,SPL(LP-1)
702! if(gdlat(i,j) > 35.0.and.gdlat(i,j)<=37.0 .and.
703! 1 gdlon(i,j) > -100.0 .and. gdlon(i,j) < -96.0)print*,
704! 2 'Debug:I,J,PNL1,TSL,NL1X,ZINT,FSL= ',I,J,PNL1,TSL(I,J)
705! 3 ,NL1X(I,J),ZINT(I,J,NL1X(I,J)),FSL(I,J)/G
706! end if
707! if(lp == lsm)print*,'Debug:undergound T,Q,U,V,FSL='
708! 1,TSL(I,J),QSL(I,J),USL(I,J),VSL(I,J),FSL(I,J)
709!
710!--- Set hydrometeor fields to zero below ground
711 c1d(i,j) = 0.
712 qw1(i,j) = 0.
713 qi1(i,j) = 0.
714 qr1(i,j) = 0.
715 qs1(i,j) = 0.
716 qg1(i,j) = 0.
717 dbz1(i,j) = dbzmin
718 frime(i,j) = 1.
719 rad(i,j) = 0.
720 o3sl(i,j) = o3(i,j,llmh)
721 cfrsl(i,j) = 0.
722 END IF
723! Compute heights by interpolating from heights on interface for NAM but
724! hydrostaticJ integration for GFS
725
726 IF(modelname == 'GFS') then
727 l=ll
728 IF(spl(lp) < pmid(i,j,1)) THEN ! above model domain
729 tvd = t(i,j,1)*(1+con_fvirt*q(i,j,1))
730 fsl(i,j) = zmid(i,j,1)-con_rog*tvd *(alsl(lp)-log(pmid(i,j,1)))
731 fsl(i,j) = fsl(i,j)*g
732 ELSE IF(l <= llmh)THEN
733 tvd = t(i,j,l)*(1+con_fvirt*q(i,j,l))
734 tvu = tsl(i,j)*(1+con_fvirt*qsl(i,j))
735 fsl(i,j) = zmid(i,j,l)-con_rog*0.5*(tvd+tvu) &
736 * (alsl(lp)-log(pmid(i,j,l)))
737 fsl(i,j) = fsl(i,j)*g
738 END IF
739 ELSE
740 la = nl1xf(i,j)
741 IF(nl1xf(i,j)<=(llmh+1)) THEN
742 fact = (alsl(lp)-log(pint(i,j,la)))/ &
743 (log(pint(i,j,la))-log(pint(i,j,la-1)))
744 IF(zint(i,j,la) < spval .AND. zint(i,j,la-1) < spval) &
745 fsl(i,j) = zint(i,j,la)+(zint(i,j,la)-zint(i,j,la-1))*fact
746 fsl(i,j) = fsl(i,j)*g
747 ELSE
748 fsl(i,j) = fprs(i,j,lp-1)-rd*(tprs(i,j,lp-1) &
749 * (h1+d608*qprs(i,j,lp-1)) &
750 + tsl(i,j)*(h1+d608*qsl(i,j))) &
751 * log(spl(lp)/spl(lp-1))/2.0
752 END IF
753 END IF
754
755 enddo ! End of i loop
756 enddo ! End of J loop
757
758!
759!*** FILL THE 3-D-IN-PRESSURE ARRAYS FOR THE MEMBRANE SLP REDUCTION
760!
761!$omp parallel do private(i,j)
762 DO j=jsta,jend
763 DO i=ista,iend
764 tprs(i,j,lp) = tsl(i,j)
765 qprs(i,j,lp) = qsl(i,j)
766 fprs(i,j,lp) = fsl(i,j)
767 ENDDO
768 ENDDO
769!
770! VERTICAL INTERPOLATION FOR WIND FOR E GRID
771!
772 IF(gridtype == 'E')THEN
773 DO j=jsta,jend
774! DO I=2,IM-MOD(J,2)
775 DO i=ista_m,iend-mod(j,2)
776! IF(i == im/2 .and. j == (jsta+jend)/2)then
777! do l=1,lm
778! print*,'PMIDV=',PMIDV(i,j,l)
779! end do
780! end if
781!
782!*** LOCATE VERTICAL INDEX OF MODEL MIDLAYER FOR V POINT JUST BELOW
783!*** THE PRESSURE LEVEL TO WHICH WE ARE INTERPOLATING.
784!
785 nl1x(i,j) = lp1
786 DO l=2,lm
787
788! IF(J == 1 .AND. I < IM)THEN !SOUTHERN BC
789! PDV=0.5*(PMID(I,J,L)+PMID(I+1,J,L))
790! ELSE IF(J == JM .AND. I < IM)THEN !NORTHERN BC
791! PDV=0.5*(PMID(I,J,L)+PMID(I+1,J,L))
792! ELSE IF(I == 1 .AND. MOD(J,2) == 0) THEN !WESTERN EVEN BC
793! PDV=0.5*(PMID(I,J-1,L)+PMID(I,J+1,L))
794! ELSE IF(I == IM .AND. MOD(J,2) == 0) THEN !EASTERN EVEN BC
795! PDV=0.5*(PMID(I,J-1,L)+PMID(I,J+1,L))
796! ELSE IF (MOD(J,2) < 1) THEN
797! PDV=0.25*(PMID(I,J,L)+PMID(I-1,J,L)
798! & +PMID(I,J+1,L)+PMID(I,J-1,L))
799! ELSE
800! PDV=0.25*(PMID(I,J,L)+PMID(I+1,J,L)
801! & +PMID(I,J+1,L)+PMID(I,J-1,L))
802! END IF
803! JJB=JSTA
804! IF(MOD(JSTA,2) == 0)JJB=JSTA+1
805! JJE=JEND
806! IF(MOD(JEND,2) == 0)JJE=JEND-1
807! DO J=JJB,JJE,2 !chc
808! PDV(IM,J)=PDV(IM-1,J)
809! END DO
810
811 IF(nl1x(i,j) == lp1.AND.pmidv(i,j,l) > spl(lp))THEN
812 nl1x(i,j) = l
813! IF(i == im/2 .and. j == jm/2)print*, &
814! 'Wind Debug:LP,NL1X',LP,NL1X(I,J)
815 ENDIF
816 ENDDO
817!
818! IF THE PRESSURE LEVEL IS BELOW THE LOWEST MODEL MIDLAYER
819! BUT STILL ABOVE THE LOWEST MODEL BOTTOM INTERFACE,
820! WE WILL NOT CONSIDER IT UNDERGROUND AND THE INTERPOLATION
821! WILL EXTRAPOLATE TO THAT POINT
822!
823! IF(NL1X(I,J) == LMP1.AND.PINT(I,J,LMP1) > SPL(LP))THEN
824 IF(nl1x(i,j) == lp1)THEN
825 IF(j == jsta .AND. i < iend)THEN !SOUTHERN BC
826 pdv = 0.5*(pint(i,j,lp1)+pint(i+1,j,lp1))
827 ELSE IF(j == jend .AND. i < iend)THEN !NORTHERN BC
828 pdv = 0.5*(pint(i,j,lp1)+pint(i+1,j,lp1))
829 ELSE IF(i == ista .AND. mod(j,2) == 0) THEN !WESTERN EVEN BC
830 pdv = 0.5*(pint(i,j-1,lp1)+pint(i,j+1,lp1))
831 ELSE IF(i == iend .AND. mod(j,2) == 0) THEN !EASTERN EVEN BC
832 pdv = 0.5*(pint(i,j-1,lp1)+pint(i,j+1,lp1))
833 ELSE IF (mod(j,2) < 1) THEN
834 pdv = 0.25*(pint(i,j,lp1)+pint(i-1,j,lp1) &
835 + pint(i,j+1,lp1)+pint(i,j-1,lp1))
836 ELSE
837 pdv = 0.25*(pint(i,j,lp1)+pint(i+1,j,lp1) &
838 + pint(i,j+1,lp1)+pint(i,j-1,lp1))
839 END IF
840 IF(pdv > spl(lp))THEN
841 nl1x(i,j) = lm
842 END IF
843 ENDIF
844!
845 ENDDO
846 ENDDO
847!
848 DO j=jsta,jend
849! DO I=1,IM-MOD(j,2)
850 DO i=ista,iend-mod(j,2)
851 ll = nl1x(i,j)
852!---------------------------------------------------------------------
853!*** VERTICAL INTERPOLATION OF WINDS FOR A-E GRID
854!---------------------------------------------------------------------
855!
856!HC IF(NL1X(I,J)<=LM)THEN
857 llmh = nint(lmh(i,j))
858
859 IF(spl(lp) < pint(i,j,2))THEN ! Above second interface
860 IF(uh(i,j,1) < spval) usl(i,j) = uh(i,j,1)
861 IF(vh(i,j,1) < spval) vsl(i,j) = vh(i,j,1)
862
863 ELSE IF(nl1x(i,j)<=llmh)THEN
864!
865!---------------------------------------------------------------------
866! INTERPOLATE LINEARLY IN LOG(P)
867!*** EXTRAPOLATE ABOVE THE TOPMOST MIDLAYER OF THE MODEL
868!*** INTERPOLATION BETWEEN NORMAL LOWER AND UPPER BOUNDS
869!*** EXTRAPOLATE BELOW LOWEST MODEL MIDLAYER (BUT STILL ABOVE GROUND)
870!---------------------------------------------------------------------
871!
872
873 fact = (alsl(lp)-log(pmidv(i,j,ll)))/ &
874 (log(pmidv(i,j,ll))-log(pmidv(i,j,ll-1)))
875 IF(uh(i,j,ll) < spval .AND. uh(i,j,ll-1) < spval) &
876 usl(i,j) = uh(i,j,ll)+(uh(i,j,ll)-uh(i,j,ll-1))*fact
877 IF(vh(i,j,ll) < spval .AND. vh(i,j,ll-1) < spval) &
878 vsl(i,j) = vh(i,j,ll)+(vh(i,j,ll)-vh(i,j,ll-1))*fact
879! IF(i == im/2 .and. j == jm/2)print*, &
880! 'Wind Debug:LP,NL1X,FACT=',LP,NL1X(I,J),FACT
881!
882! FOR UNDERGROUND PRESSURE LEVELS, ASSUME TEMPERATURE TO CHANGE
883! ADIABATICLY, RH TO BE THE SAME AS THE AVERAGE OF THE 2ND AND 3RD
884! LAYERS FROM THE GOUND, WIND TO BE THE SAME AS THE LOWEST LEVEL ABOVE
885! GOUND
886 ELSE
887 IF(uh(i,j,llmh) < spval) usl(i,j) = uh(i,j,llmh)
888 IF(vh(i,j,llmh) < spval) vsl(i,j) = vh(i,j,llmh)
889 END IF
890 ENDDO ! end of i loop
891 ENDDO ! end of j loop
892
893! if(me == 0) print *,'after 230 me=',me,'USL=',USL(1:10,JSTA)
894 jjb = jsta
895 IF(mod(jsta,2) == 0) jjb = jsta+1
896 jje = jend
897 IF(mod(jend,2) == 0) jje = jend-1
898 DO j=jjb,jje,2 !chc
899 usl(iend,j) = usl(iend-1,j)
900 vsl(iend,j) = vsl(iend-1,j)
901 END DO
902 ELSE IF(gridtype=='B')THEN ! B grid wind interpolation
903 DO j=jsta,jend_m
904! DO I=1,IM-1
905 DO i=ista,iend_m
906!*** LOCATE VERTICAL INDEX OF MODEL MIDLAYER FOR V POINT JUST BELOW
907!*** THE PRESSURE LEVEL TO WHICH WE ARE INTERPOLATING.
908!
909 nl1x(i,j)=lp1
910 DO l=2,lm
911 IF(nl1x(i,j) == lp1.AND.pmidv(i,j,l) > spl(lp))THEN
912 nl1x(i,j) = l
913! IF(i == im/2 .and. j == jm/2)print*, &
914! 'Wind Debug for B grid:LP,NL1X',LP,NL1X(I,J)
915 ENDIF
916 ENDDO
917!
918! IF THE PRESSURE LEVEL IS BELOW THE LOWEST MODEL MIDLAYER
919! BUT STILL ABOVE THE LOWEST MODEL BOTTOM INTERFACE,
920! WE WILL NOT CONSIDER IT UNDERGROUND AND THE INTERPOLATION
921! WILL EXTRAPOLATE TO THAT POINT
922!
923 IF(nl1x(i,j)==lp1)THEN
924 pdv = 0.25*(pint(i,j,lp1)+pint(i+1,j,lp1) &
925 + pint(i,j+1,lp1)+pint(i+1,j+1,lp1))
926 IF(pdv > spl(lp))THEN
927 nl1x(i,j)=lm
928 END IF
929 ENDIF
930!
931 ENDDO
932 ENDDO
933!
934 DO j=jsta,jend_m
935! DO I=1,IM-1
936 DO i=ista,iend_m
937 ll = nl1x(i,j)
938!---------------------------------------------------------------------
939!*** VERTICAL INTERPOLATION OF WINDS FOR A-E GRID
940!---------------------------------------------------------------------
941!
942!HC IF(NL1X(I,J)<=LM)THEN
943 llmh = nint(lmh(i,j))
944
945 IF(spl(lp) < pint(i,j,2))THEN ! Above second interface
946 IF(uh(i,j,1) < spval) usl(i,j) = uh(i,j,1)
947 IF(vh(i,j,1) < spval) vsl(i,j) = vh(i,j,1)
948
949 ELSE IF(nl1x(i,j)<=llmh)THEN
950!
951!---------------------------------------------------------------------
952! INTERPOLATE LINEARLY IN LOG(P)
953!*** EXTRAPOLATE ABOVE THE TOPMOST MIDLAYER OF THE MODEL
954!*** INTERPOLATION BETWEEN NORMAL LOWER AND UPPER BOUNDS
955!*** EXTRAPOLATE BELOW LOWEST MODEL MIDLAYER (BUT STILL ABOVE GROUND)
956!---------------------------------------------------------------------
957!
958
959 fact = (alsl(lp)-log(pmidv(i,j,ll)))/ &
960 (log(pmidv(i,j,ll))-log(pmidv(i,j,ll-1)))
961 IF(uh(i,j,ll) < spval .AND. uh(i,j,ll-1) < spval) &
962 usl(i,j)=uh(i,j,ll)+(uh(i,j,ll)-uh(i,j,ll-1))*fact
963 IF(vh(i,j,ll) < spval .AND. vh(i,j,ll-1) < spval) &
964 vsl(i,j)=vh(i,j,ll)+(vh(i,j,ll)-vh(i,j,ll-1))*fact
965! IF(i == im/2 .and. j == jm/2)print*, &
966! 'Wind Debug:LP,NL1X,FACT=',LP,NL1X(I,J),FACT
967!
968! FOR UNDERGROUND PRESSURE LEVELS, ASSUME TEMPERATURE TO CHANGE
969! ADIABATICLY, RH TO BE THE SAME AS THE AVERAGE OF THE 2ND AND 3RD
970! LAYERS FROM THE GOUND, WIND TO BE THE SAME AS THE LOWEST LEVEL ABOVE
971! GOUND
972 ELSE
973 IF(uh(i,j,llmh) < spval)usl(i,j)=uh(i,j,llmh)
974 IF(vh(i,j,llmh) < spval)vsl(i,j)=vh(i,j,llmh)
975 END IF
976 enddo
977 enddo
978 END IF ! END OF WIND INTERPOLATION FOR NMM
979! if(me == 0) print *,'after 230 if me=',me,'USL=',USL(1:10,JSTA)
980
981
982!
983!---------------------------------------------------------------------
984! LOAD GEOPOTENTIAL AND TEMPERATURE INTO STANDARD LEVEL
985! ARRAYS FOR THE NEXT PASS.
986!---------------------------------------------------------------------
987!
988!*** SAVE 500MB TEMPERATURE FOR LIFTED INDEX.
989!
990 IF(nint(spl(lp)) == 50000)THEN
991!$omp parallel do private(i,j)
992 DO j=jsta,jend
993 DO i=ista,iend
994 t500(i,j) = tsl(i,j)
995 z500(i,j) = fsl(i,j)*gi
996 ENDDO
997 ENDDO
998 ENDIF
999
1000!
1001!*** SAVE 700MB TEMPERATURE FOR LIFTED INDEX.
1002!
1003 IF(nint(spl(lp)) == 70000)THEN
1004!$omp parallel do private(i,j)
1005 DO j=jsta,jend
1006 DO i=ista,iend
1007 t700(i,j) = tsl(i,j)
1008 z700(i,j) = fsl(i,j)*gi
1009 ENDDO
1010 ENDDO
1011 ENDIF
1012
1013!
1014!---------------------------------------------------------------------
1015!*** CALCULATE 1000MB GEOPOTENTIALS CONSISTENT WITH SLP OBTAINED
1016!*** FROM THE MESINGER OR NWS SHUELL SLP REDUCTION.
1017!---------------------------------------------------------------------
1018!
1019!*** FROM MESINGER SLP
1020!
1021!HC MOVE THIS PART TO THE END OF THIS SUBROUTINE AFTER PSLP IS COMPUTED
1022!HC IF(IGET(023) > 0.AND.NINT(SPL(LP)) == 100000)THEN
1023!HC ALPTH=LOG(1.E5)
1024!HC!$omp parallel do private(i,j)
1025!HC DO J=JSTA,JEND
1026!HC DO I=1,IM
1027!HC IF(FSL(I,J) < SPVAL) THEN
1028!HC PSLPIJ=PSLP(I,J)
1029!HC ALPSL=LOG(PSLPIJ)
1030!HC PSFC=PINT(I,J,NINT(LMH(I,J))+1)
1031!HC IF(ABS(PSLPIJ-PSFC) < 5.E2) THEN
1032!HC FSL(I,J)=R*TSL(I,J)*(ALPSL-ALPTH)
1033!HC ELSE
1034!HC FSL(I,J)=FIS(I,J)/(ALPSL-LOG(PSFC))*
1035!HC 1 (ALPSL-ALPTH)
1036!HC ENDIF
1037!HC Z1000(I,J)=FSL(I,J)*GI
1038!HC ELSE
1039!HC Z1000(I,J)=SPVAL
1040!HC ENDIF
1041!HC ENDDO
1042!HC ENDDO
1043!
1044!*** FROM NWS SHUELL SLP. NGMSLP2 COMPUTES 1000MB GEOPOTENTIAL.
1045!
1046!HC ELSEIF(IGET(023)<=0.AND.LP == LSM)THEN
1047!HC IF(IGET(023)<=0.AND.LP == LSM)THEN
1048!!$omp parallel do private(i,j)
1049!HC DO J=JSTA,JEND
1050!HC DO I=1,IM
1051!HC IF(Z1000(I,J) < SPVAL) THEN
1052!HC FSL(I,J)=Z1000(I,J)*G
1053!HC ELSE
1054!HC FSL(I,J)=SPVAL
1055!HC ENDIF
1056!HC ENDDO
1057!HC ENDDO
1058!HC ENDIF
1059!---------------------------------------------------------------------
1060!---------------------------------------------------------------------
1061! *** PART II ***
1062!---------------------------------------------------------------------
1063!---------------------------------------------------------------------
1064!
1065! INTERPOLATE/OUTPUT SELECTED FIELDS.
1066!
1067!---------------------------------------------------------------------
1068!
1069!*** OUTPUT GEOPOTENTIAL (SCALE BY GI)
1070!
1071 IF(iget(012) > 0)THEN
1072 IF(lvls(lp,iget(012)) > 0)THEN
1073 IF((iget(023) > 0 .OR. iget(445) > 0) .AND. nint(spl(lp)) == 100000) THEN
1074 ! GO TO 222
1075 ELSE
1076!$omp parallel do private(i,j)
1077 DO j=jsta,jend
1078 DO i=ista,iend
1079 IF(fsl(i,j) < spval) THEN
1080 grid1(i,j) = fsl(i,j)*gi
1081 ELSE
1082 grid1(i,j) = spval
1083 ENDIF
1084 ENDDO
1085 ENDDO
1086
1087 IF (smflag) THEN
1088!tgs - smoothing of geopotential heights
1089 if(maptype == 6) then
1090 if(grib=='grib2') then
1091 dxm = (dxval / 360.)*(erad*2.*pi)/1.d6 ! [mm]
1092 endif
1093 else
1094 dxm = dxval
1095 endif
1096 if(grib == 'grib2')then
1097 dxm=dxm/1000.0
1098 endif
1099! print *,'dxm=',dxm
1100 nsmooth = nint(5.*(13500./dxm))
1101 call allgetherv(grid1)
1102 do k=1,nsmooth
1103 CALL smooth(grid1,sdummy,im,jm,0.5)
1104 end do
1105 ENDIF
1106 if(grib == 'grib2')then
1107 cfld = cfld + 1
1108 fld_info(cfld)%ifld=iavblfld(iget(012))
1109 fld_info(cfld)%lvl=lvlsxml(lp,iget(012))
1110!$omp parallel do private(i,j,ii,jj)
1111 do j=1,jend-jsta+1
1112 jj = jsta+j-1
1113 do i=1,iend-ista+1
1114 ii=ista+i-1
1115 datapd(i,j,cfld) = grid1(ii,jj)
1116 enddo
1117 enddo
1118 endif
1119 END IF
1120 ENDIF
1121 ENDIF
1122!222 CONTINUE
1123!
1124!*** TEMPERATURE
1125!
1126 IF(iget(013) > 0) THEN
1127 IF(lvls(lp,iget(013)) > 0)THEN
1128!$omp parallel do private(i,j)
1129 DO j=jsta,jend
1130 DO i=ista,iend
1131 grid1(i,j) = tsl(i,j)
1132 ENDDO
1133 ENDDO
1134
1135 IF (smflag) THEN
1136 nsmooth = nint(3.*(13500./dxm))
1137 call allgetherv(grid1)
1138 do k=1,nsmooth
1139 CALL smooth(grid1,sdummy,im,jm,0.5)
1140 end do
1141 ENDIF
1142
1143 if(grib == 'grib2')then
1144 cfld = cfld + 1
1145 fld_info(cfld)%ifld = iavblfld(iget(013))
1146 fld_info(cfld)%lvl = lvlsxml(lp,iget(013))
1147!$omp parallel do private(i,j,ii,jj)
1148 do j=1,jend-jsta+1
1149 jj = jsta+j-1
1150 do i=1,iend-ista+1
1151 ii=ista+i-1
1152 datapd(i,j,cfld) = grid1(ii,jj)
1153 enddo
1154 enddo
1155 endif
1156 ENDIF
1157 ENDIF
1158
1159!*** virtual TEMPERATURE
1160!
1161 IF(iget(910)>0) THEN
1162 IF(lvls(lp,iget(910))>0)THEN
1163!$omp parallel do private(i,j)
1164 DO j=jsta,jend
1165 DO i=ista,iend
1166 IF(tsl(i,j) < spval .AND. qsl(i,j) < spval) THEN
1167 grid1(i,j) = tsl(i,j)*(1.+0.608*qsl(i,j))
1168 ELSE
1169 grid1(i,j) = spval
1170 ENDIF
1171 ENDDO
1172 ENDDO
1173
1174 IF (smflag) THEN
1175 nsmooth = nint(3.*(13500./dxm))
1176 call allgetherv(grid1)
1177 do k=1,nsmooth
1178 CALL smooth(grid1,sdummy,im,jm,0.5)
1179 end do
1180 ENDIF
1181
1182 if(grib=='grib2')then
1183 cfld=cfld+1
1184 fld_info(cfld)%ifld = iavblfld(iget(910))
1185 fld_info(cfld)%lvl = lvlsxml(lp,iget(910))
1186!$omp parallel do private(i,j,ii,jj)
1187 do j=1,jend-jsta+1
1188 jj = jsta+j-1
1189 do i=1,iend-ista+1
1190 ii=ista+i-1
1191 datapd(i,j,cfld) = grid1(ii,jj)
1192 enddo
1193 enddo
1194 endif
1195 ENDIF
1196 ENDIF
1197
1198!
1199!*** POTENTIAL TEMPERATURE.
1200!
1201 IF(iget(014) > 0)THEN
1202 IF(lvls(lp,iget(014)) > 0)THEN
1203
1204 tem = (p1000/spl(lp)) ** capa
1205!$omp parallel do private(i,j)
1206 DO j=jsta,jend
1207 DO i=ista,iend
1208 IF(tsl(i,j) < spval) THEN
1209 grid1(i,j) = tsl(i,j) * tem
1210 ELSE
1211 grid1(i,j) = spval
1212 ENDIF
1213 ENDDO
1214 ENDDO
1215!!$omp parallel do private(i,j)
1216! DO J=JSTA,JEND
1217! DO I=1,IM
1218! EGRID2(I,J) = SPL(LP)
1219! ENDDO
1220! ENDDO
1221!
1222! CALL CALPOT(EGRID2,TSL,EGRID1)
1223!!$omp parallel do private(i,j)
1224! DO J=JSTA,JEND
1225! DO I=1,IM
1226! GRID1(I,J) = EGRID1(I,J)
1227! ENDDO
1228! ENDDO
1229
1230 if(grib == 'grib2')then
1231 cfld = cfld + 1
1232 fld_info(cfld)%ifld=iavblfld(iget(014))
1233 fld_info(cfld)%lvl=lvlsxml(lp,iget(014))
1234!$omp parallel do private(i,j,ii,jj)
1235 do j=1,jend-jsta+1
1236 jj = jsta+j-1
1237 do i=1,iend-ista+1
1238 ii=ista+i-1
1239 datapd(i,j,cfld) = grid1(ii,jj)
1240 enddo
1241 enddo
1242 endif
1243 ENDIF
1244 ENDIF
1245!
1246!*** RELATIVE HUMIDITY.
1247!
1248
1249 IF(iget(017) > 0 .OR. iget(257) > 0)THEN
1250! if ( me == 0) print *,'IGET(17)=',IGET(017),'LP=',LP,IGET(257), &
1251! 'LVLS=',LVLS(1,4)
1252 log1=.false.
1253 IF(iget(017) > 0.) then
1254 if(lvls(lp,iget(017)) > 0 ) log1=.true.
1255 endif
1256 IF(iget(257) > 0) then
1257 if(lvls(lp,iget(257)) > 0 ) log1=.true.
1258 endif
1259 if ( log1 ) then
1260!$omp parallel do private(i,j)
1261 DO j=jsta,jend
1262 DO i=ista,iend
1263 egrid2(i,j) = spl(lp)
1264 ENDDO
1265 ENDDO
1266!
1267 CALL calrh(egrid2(ista:iend,jsta:jend),tsl(ista:iend,jsta:jend),qsl(ista:iend,jsta:jend),egrid1(ista:iend,jsta:jend))
1268
1269!$omp parallel do private(i,j)
1270 DO j=jsta,jend
1271 DO i=ista,iend
1272 IF(egrid1(i,j) < spval) THEN
1273 grid1(i,j) = egrid1(i,j)*100.
1274 ELSE
1275 grid1(i,j) = egrid1(i,j)
1276 ENDIF
1277 ENDDO
1278 ENDDO
1279
1280 IF (smflag) THEN
1281 nsmooth=nint(2.*(13500./dxm))
1282 call allgetherv(grid1)
1283 do k=1,nsmooth
1284 CALL smooth(grid1,sdummy,im,jm,0.5)
1285 end do
1286 ENDIF
1287 if(grib == 'grib2')then
1288 cfld = cfld + 1
1289 fld_info(cfld)%ifld=iavblfld(iget(017))
1290 fld_info(cfld)%lvl=lvlsxml(lp,iget(017))
1291!$omp parallel do private(i,j,ii,jj)
1292 do j=1,jend-jsta+1
1293 jj = jsta+j-1
1294 do i=1,iend-ista+1
1295 ii=ista+i-1
1296 datapd(i,j,cfld) = grid1(ii,jj)
1297 enddo
1298 enddo
1299 endif
1300
1301!$omp parallel do private(i,j)
1302 DO j=jsta,jend
1303 DO i=ista,iend
1304 savrh(i,j) = grid1(i,j)
1305 ENDDO
1306 ENDDO
1307
1308 ENDIF
1309 ENDIF
1310!
1311!*** CLOUD FRACTION.
1312!
1313 IF(iget(331) > 0)THEN
1314 IF(lvls(lp,iget(331)) > 0)THEN
1315!$omp parallel do private(i,j)
1316 DO j=jsta,jend
1317 DO i=ista,iend
1318 grid1(i,j) = spval
1319 cfrsl(i,j) = min(max(0.0,cfrsl(i,j)),1.0)
1320 IF(abs(cfrsl(i,j)-spval) > small) &
1321 grid1(i,j) = cfrsl(i,j)*h100
1322 ENDDO
1323 ENDDO
1324 if(grib == 'grib2')then
1325 cfld = cfld + 1
1326 fld_info(cfld)%ifld = iavblfld(iget(331))
1327 fld_info(cfld)%lvl = lvlsxml(lp,iget(331))
1328!$omp parallel do private(i,j,ii,jj)
1329 do j=1,jend-jsta+1
1330 jj = jsta+j-1
1331 do i=1,iend-ista+1
1332 ii=ista+i-1
1333 datapd(i,j,cfld) = grid1(ii,jj)
1334 enddo
1335 enddo
1336 endif
1337 ENDIF
1338 ENDIF
1339!
1340!*** DEWPOINT TEMPERATURE.
1341!
1342 IF(iget(015) > 0)THEN
1343 IF(lvls(lp,iget(015)) > 0)THEN
1344!$omp parallel do private(i,j)
1345 DO j=jsta,jend
1346 DO i=ista,iend
1347 egrid2(i,j) = spl(lp)
1348 ENDDO
1349 ENDDO
1350!
1351 CALL caldwp(egrid2(ista:iend,jsta:jend),qsl(ista:iend,jsta:jend),egrid1(ista:iend,jsta:jend),tsl(ista:iend,jsta:jend))
1352!$omp parallel do private(i,j)
1353 DO j=jsta,jend
1354 DO i=ista,iend
1355 IF(tsl(i,j) < spval) THEN
1356 grid1(i,j) = egrid1(i,j)
1357 ELSE
1358 grid1(i,j) = spval
1359 ENDIF
1360 ENDDO
1361 ENDDO
1362 if(grib == 'grib2')then
1363 cfld = cfld + 1
1364 fld_info(cfld)%ifld=iavblfld(iget(015))
1365 fld_info(cfld)%lvl=lvlsxml(lp,iget(015))
1366!$omp parallel do private(i,j,ii,jj)
1367 do j=1,jend-jsta+1
1368 jj = jsta+j-1
1369 do i=1,iend-ista+1
1370 ii=ista+i-1
1371 datapd(i,j,cfld) = grid1(ii,jj)
1372 enddo
1373 enddo
1374 endif
1375 ENDIF
1376 ENDIF
1377!
1378!*** SPECIFIC HUMIDITY.
1379!
1380 IF(iget(016) > 0)THEN
1381 IF(lvls(lp,iget(016)) > 0)THEN
1382!$omp parallel do private(i,j)
1383 DO j=jsta,jend
1384 DO i=ista,iend
1385 grid1(i,j) = qsl(i,j)
1386 ENDDO
1387 ENDDO
1388 CALL bound(grid1,zero,h99999)
1389 if(grib == 'grib2')then
1390 cfld = cfld + 1
1391 fld_info(cfld)%ifld=iavblfld(iget(016))
1392 fld_info(cfld)%lvl=lvlsxml(lp,iget(016))
1393!$omp parallel do private(i,j,ii,jj)
1394 do j=1,jend-jsta+1
1395 jj = jsta+j-1
1396 do i=1,iend-ista+1
1397 ii=ista+i-1
1398 datapd(i,j,cfld) = grid1(ii,jj)
1399 enddo
1400 enddo
1401 endif
1402 ENDIF
1403 ENDIF
1404!
1405!*** OMEGA
1406!
1407 IF(iget(020) > 0)THEN
1408 IF(lvls(lp,iget(020)) > 0)THEN
1409!$omp parallel do private(i,j)
1410 DO j=jsta,jend
1411 DO i=ista,iend
1412 grid1(i,j) = osl(i,j)
1413 ENDDO
1414 ENDDO
1415
1416 IF (smflag .or. ioform == 'binarympiio' ) THEN
1417 call allgetherv(grid1)
1418 if (ioform == 'binarympiio') then
1419! nsmooth = max(2, min(30,nint(jm/94.0)))
1420! do k=1,5
1421 CALL smoothc(grid1,sdummy,im,jm,0.5)
1422 CALL smoothc(grid1,sdummy,im,jm,-0.5)
1423! enddo
1424 else
1425 nsmooth = nint(3.*(13500./dxm))
1426! endif
1427 do k=1,nsmooth
1428 CALL smooth(grid1,sdummy,im,jm,0.5)
1429 end do
1430 endif
1431 ENDIF
1432
1433 if(grib == 'grib2')then
1434 cfld = cfld + 1
1435 fld_info(cfld)%ifld=iavblfld(iget(020))
1436 fld_info(cfld)%lvl=lvlsxml(lp,iget(020))
1437!$omp parallel do private(i,j,ii,jj)
1438 do j=1,jend-jsta+1
1439 jj = jsta+j-1
1440 do i=1,iend-ista+1
1441 ii=ista+i-1
1442 datapd(i,j,cfld) = grid1(ii,jj)
1443 enddo
1444 enddo
1445 endif
1446 ENDIF
1447 ENDIF
1448!
1449!*** W
1450!
1451 IF(iget(284) > 0)THEN
1452 IF(lvls(lp,iget(284)) > 0)THEN
1453!$omp parallel do private(i,j)
1454 DO j=jsta,jend
1455 DO i=ista,iend
1456 grid1(i,j) = wsl(i,j)
1457 ENDDO
1458 ENDDO
1459 if(grib == 'grib2')then
1460 cfld = cfld + 1
1461 fld_info(cfld)%ifld=iavblfld(iget(284))
1462 fld_info(cfld)%lvl=lvlsxml(lp,iget(284))
1463!$omp parallel do private(i,j,ii,jj)
1464 do j=1,jend-jsta+1
1465 jj = jsta+j-1
1466 do i=1,iend-ista+1
1467 ii=ista+i-1
1468 datapd(i,j,cfld) = grid1(ii,jj)
1469 enddo
1470 enddo
1471 endif
1472 ENDIF
1473 ENDIF
1474!
1475!*** MOISTURE CONVERGENCE
1476!
1477 IF(iget(085) > 0)THEN
1478 IF(lvls(lp,iget(085)) > 0)THEN
1479 CALL calmcvg(qsl(ista_2l,jsta_2l),usl(ista_2l,jsta_2l),vsl(ista_2l,jsta_2l),egrid1(ista_2l,jsta_2l))
1480! if(me == 0) print *,'after calmcvgme=',me,'USL=',USL(1:10,JSTA)
1481!$omp parallel do private(i,j)
1482 DO j=jsta,jend
1483 DO i=ista,iend
1484 grid1(i,j) = egrid1(i,j)
1485 ENDDO
1486 ENDDO
1487!MEB NOT SURE IF I STILL NEED THIS
1488! CONVERT TO DIVERGENCE FOR GRIB UNITS
1489!
1490! CALL SCLFLD(GRID1(ista:iend,jsta:jend),-1.0,IM,JM)
1491!MEB NOT SURE IF I STILL NEED THIS
1492 if(grib == 'grib2')then
1493 cfld = cfld + 1
1494 fld_info(cfld)%ifld=iavblfld(iget(085))
1495 fld_info(cfld)%lvl=lvlsxml(lp,iget(085))
1496!$omp parallel do private(i,j,ii,jj)
1497 do j=1,jend-jsta+1
1498 jj = jsta+j-1
1499 do i=1,iend-ista+1
1500 ii=ista+i-1
1501 datapd(i,j,cfld) = grid1(ii,jj)
1502 enddo
1503 enddo
1504! if(me==0) print *,'in mdl2p,mconv, lp=',fld_info(cfld)%lvl,'lp=',lp
1505 endif
1506 ENDIF
1507 ENDIF
1508!
1509!*** U AND/OR V WIND
1510!
1511 IF(iget(018) > 0.OR.iget(019) > 0)THEN
1512 log1=.false.
1513 IF(iget(018) > 0.) then
1514 if(lvls(lp,iget(018)) > 0 ) log1=.true.
1515 endif
1516 IF(iget(019) > 0) then
1517 if(lvls(lp,iget(019)) > 0 ) log1=.true.
1518 endif
1519 if ( log1 ) then
1520!$omp parallel do private(i,j)
1521 DO j=jsta,jend
1522 DO i=ista,iend
1523 grid1(i,j) = usl(i,j)
1524 grid2(i,j) = vsl(i,j)
1525 ENDDO
1526 ENDDO
1527
1528 IF (smflag) THEN
1529 nsmooth=nint(5.*(13500./dxm))
1530 call allgetherv(grid1)
1531 do k=1,nsmooth
1532 CALL smooth(grid1,sdummy,im,jm,0.5)
1533 end do
1534 nsmooth=nint(5.*(13500./dxm))
1535 call allgetherv(grid2)
1536 do k=1,nsmooth
1537 CALL smooth(grid2,sdummy,im,jm,0.5)
1538 end do
1539 ENDIF
1540
1541 if(grib == 'grib2')then
1542 cfld = cfld + 1
1543 fld_info(cfld)%ifld=iavblfld(iget(018))
1544 fld_info(cfld)%lvl=lvlsxml(lp,iget(018))
1545!$omp parallel do private(i,j,ii,jj)
1546 do j=1,jend-jsta+1
1547 jj = jsta+j-1
1548 do i=1,iend-ista+1
1549 ii=ista+i-1
1550 datapd(i,j,cfld) = grid1(ii,jj)
1551 enddo
1552 enddo
1553
1554 cfld = cfld + 1
1555 fld_info(cfld)%ifld=iavblfld(iget(019))
1556 fld_info(cfld)%lvl=lvlsxml(lp,iget(019))
1557!$omp parallel do private(i,j,ii,jj)
1558 do j=1,jend-jsta+1
1559 jj = jsta+j-1
1560 do i=1,iend-ista+1
1561 ii=ista+i-1
1562 datapd(i,j,cfld) = grid2(ii,jj)
1563 enddo
1564 enddo
1565 endif
1566 ENDIF
1567 ENDIF
1568!
1569!*** ABSOLUTE VORTICITY
1570!
1571 IF (iget(021) > 0) THEN
1572 IF (lvls(lp,iget(021)) > 0) THEN
1573 CALL calvor(usl,vsl,egrid1)
1574! print *,'me=',me,'EGRID1=',EGRID1(1:10,JSTA)
1575!$omp parallel do private(i,j)
1576 DO j=jsta,jend
1577 DO i=ista,iend
1578 grid1(i,j) = egrid1(i,j)
1579 ENDDO
1580 ENDDO
1581
1582 IF (smflag .or. ioform == 'binarympiio' ) THEN
1583 call allgetherv(grid1)
1584 if (ioform == 'binarympiio') then
1585! nsmooth = max(2, min(30,nint(jm/94.0)))
1586! do k=1,5
1587 CALL smoothc(grid1,sdummy,im,jm,0.5)
1588 CALL smoothc(grid1,sdummy,im,jm,-0.5)
1589! enddo
1590 else
1591 nsmooth = nint(4.*(13500./dxm))
1592! endif
1593 do k=1,nsmooth
1594 CALL smooth(grid1,sdummy,im,jm,0.5)
1595 end do
1596 endif
1597 ENDIF
1598
1599 if(grib == 'grib2')then
1600 cfld = cfld + 1
1601 fld_info(cfld)%ifld=iavblfld(iget(021))
1602 fld_info(cfld)%lvl=lvlsxml(lp,iget(021))
1603!$omp parallel do private(i,j,ii,jj)
1604 do j=1,jend-jsta+1
1605 jj = jsta+j-1
1606 do i=1,iend-ista+1
1607 ii=ista+i-1
1608 datapd(i,j,cfld) = grid1(ii,jj)
1609 enddo
1610 enddo
1611 endif
1612 ENDIF
1613 ENDIF
1614!
1615! GEOSTROPHIC STREAMFUNCTION.
1616 IF (iget(086) > 0) THEN
1617 IF (lvls(lp,iget(086)) > 0) THEN
1618!$omp parallel do private(i,j)
1619 DO j=jsta,jend
1620 DO i=ista,iend
1621 IF(fsl(i,j)<spval)THEN
1622 egrid2(i,j) = fsl(i,j)*gi
1623 ENDIF
1624 ENDDO
1625 ENDDO
1626 CALL calstrm(egrid2(ista:iend,jsta:jend),egrid1(ista:iend,jsta:jend))
1627!$omp parallel do private(i,j)
1628 DO j=jsta,jend
1629 DO i=ista,iend
1630 IF(fsl(i,j) < spval) THEN
1631 grid1(i,j) = egrid1(i,j)
1632 ELSE
1633 grid1(i,j) = spval
1634 ENDIF
1635 ENDDO
1636 ENDDO
1637 if(grib == 'grib2')then
1638 cfld = cfld + 1
1639 fld_info(cfld)%ifld=iavblfld(iget(086))
1640 fld_info(cfld)%lvl=lvlsxml(lp,iget(086))
1641!$omp parallel do private(i,j,ii,jj)
1642 do j=1,jend-jsta+1
1643 jj = jsta+j-1
1644 do i=1,iend-ista+1
1645 ii=ista+i-1
1646 datapd(i,j,cfld) = grid1(ii,jj)
1647 enddo
1648 enddo
1649 endif
1650 ENDIF
1651 ENDIF
1652!
1653!*** TURBULENT KINETIC ENERGY
1654!
1655 IF (iget(022) > 0) THEN
1656 IF (lvls(lp,iget(022)) > 0) THEN
1657!$omp parallel do private(i,j)
1658 DO j=jsta,jend
1659 DO i=ista,iend
1660 grid1(i,j) = q2sl(i,j)
1661 ENDDO
1662 ENDDO
1663 if(grib == 'grib2')then
1664 cfld = cfld + 1
1665 fld_info(cfld)%ifld=iavblfld(iget(022))
1666 fld_info(cfld)%lvl=lvlsxml(lp,iget(022))
1667!$omp parallel do private(i,j,ii,jj)
1668 do j=1,jend-jsta+1
1669 jj = jsta+j-1
1670 do i=1,iend-ista+1
1671 ii=ista+i-1
1672 datapd(i,j,cfld) = grid1(ii,jj)
1673 enddo
1674 enddo
1675 endif
1676 ENDIF
1677 ENDIF
1678!
1679!*** CLOUD WATER
1680!
1681 IF (iget(153) > 0) THEN
1682 IF (lvls(lp,iget(153)) > 0) THEN
1683 IF(imp_physics==99 .or. imp_physics==98)then
1684! GFS does not seperate cloud water from ice, hoping to do that in Feb 08 implementation
1685!$omp parallel do private(i,j)
1686 DO j=jsta,jend
1687 DO i=ista,iend
1688 IF(qw1(i,j) < spval .AND. qi1(i,j) < spval) THEN
1689 grid1(i,j) = qw1(i,j) + qi1(i,j)
1690 qi1(i,j) = spval
1691 ELSE
1692 grid1(i,j) = spval
1693 ENDIF
1694 ENDDO
1695 ENDDO
1696 ELSE
1697!$omp parallel do private(i,j)
1698 DO j=jsta,jend
1699 DO i=ista,iend
1700 grid1(i,j) = qw1(i,j)
1701 ENDDO
1702 ENDDO
1703 END IF
1704 if(grib == 'grib2')then
1705 cfld = cfld + 1
1706 fld_info(cfld)%ifld=iavblfld(iget(153))
1707 fld_info(cfld)%lvl=lvlsxml(lp,iget(153))
1708!$omp parallel do private(i,j,ii,jj)
1709 do j=1,jend-jsta+1
1710 jj = jsta+j-1
1711 do i=1,iend-ista+1
1712 ii=ista+i-1
1713 datapd(i,j,cfld) = grid1(ii,jj)
1714 enddo
1715 enddo
1716 endif
1717 ENDIF
1718 ENDIF
1719!
1720!*** CLOUD ICE
1721!
1722 IF (iget(166) > 0) THEN
1723 IF (lvls(lp,iget(166)) > 0) THEN
1724!$omp parallel do private(i,j)
1725 DO j=jsta,jend
1726 DO i=ista,iend
1727 grid1(i,j) = qi1(i,j)
1728 ENDDO
1729 ENDDO
1730 if(grib == 'grib2')then
1731 cfld = cfld + 1
1732 fld_info(cfld)%ifld=iavblfld(iget(166))
1733 fld_info(cfld)%lvl=lvlsxml(lp,iget(166))
1734!$omp parallel do private(i,j,ii,jj)
1735 do j=1,jend-jsta+1
1736 jj = jsta+j-1
1737 do i=1,iend-ista+1
1738 ii=ista+i-1
1739 datapd(i,j,cfld) = grid1(ii,jj)
1740 enddo
1741 enddo
1742 endif
1743 ENDIF
1744 ENDIF
1745!
1746!--- RAIN
1747 IF (iget(183) > 0) THEN
1748 IF (lvls(lp,iget(183)) > 0) THEN
1749!$omp parallel do private(i,j)
1750 DO j=jsta,jend
1751 DO i=ista,iend
1752 grid1(i,j) = qr1(i,j)
1753 ENDDO
1754 ENDDO
1755 if(grib == 'grib2')then
1756 cfld = cfld + 1
1757 fld_info(cfld)%ifld=iavblfld(iget(183))
1758 fld_info(cfld)%lvl=lvlsxml(lp,iget(183))
1759!$omp parallel do private(i,j,ii,jj)
1760 do j=1,jend-jsta+1
1761 jj = jsta+j-1
1762 do i=1,iend-ista+1
1763 ii=ista+i-1
1764 datapd(i,j,cfld) = grid1(ii,jj)
1765 enddo
1766 enddo
1767 endif
1768 ENDIF
1769 ENDIF
1770!
1771!--- SNOW
1772 IF (iget(184) > 0) THEN
1773 IF (lvls(lp,iget(184)) > 0) THEN
1774!$omp parallel do private(i,j)
1775 DO j=jsta,jend
1776 DO i=ista,iend
1777 grid1(i,j) = qs1(i,j)
1778 ENDDO
1779 ENDDO
1780 if(grib == 'grib2')then
1781 cfld = cfld + 1
1782 fld_info(cfld)%ifld=iavblfld(iget(184))
1783 fld_info(cfld)%lvl=lvlsxml(lp,iget(184))
1784!$omp parallel do private(i,j,ii,jj)
1785 do j=1,jend-jsta+1
1786 jj = jsta+j-1
1787 do i=1,iend-ista+1
1788 ii=ista+i-1
1789 datapd(i,j,cfld) = grid1(ii,jj)
1790 enddo
1791 enddo
1792 endif
1793 ENDIF
1794 ENDIF
1795!
1796!--- GRAUPEL
1797 IF (iget(416) > 0) THEN
1798 IF (lvls(lp,iget(416)) > 0) THEN
1799!$omp parallel do private(i,j)
1800 DO j=jsta,jend
1801 DO i=ista,iend
1802 grid1(i,j) = qg1(i,j)
1803 ENDDO
1804 ENDDO
1805 if(grib == 'grib2')then
1806 cfld = cfld + 1
1807 fld_info(cfld)%ifld=iavblfld(iget(416))
1808 fld_info(cfld)%lvl=lvlsxml(lp,iget(416))
1809!$omp parallel do private(i,j,ii,jj)
1810 do j=1,jend-jsta+1
1811 jj = jsta+j-1
1812 do i=1,iend-ista+1
1813 ii=ista+i-1
1814 datapd(i,j,cfld) = grid1(ii,jj)
1815 enddo
1816 enddo
1817 endif
1818 ENDIF
1819 ENDIF
1820
1821!
1822!--- TOTAL CONDENSATE
1823 IF (iget(198) > 0) THEN
1824 IF (lvls(lp,iget(198)) > 0) THEN
1825!$omp parallel do private(i,j)
1826 DO j=jsta,jend
1827 DO i=ista,iend
1828 grid1(i,j) = c1d(i,j)
1829 ENDDO
1830 ENDDO
1831 if(grib == 'grib2')then
1832 cfld = cfld + 1
1833 fld_info(cfld)%ifld=iavblfld(iget(198))
1834 fld_info(cfld)%lvl=lvlsxml(lp,iget(198))
1835!$omp parallel do private(i,j,ii,jj)
1836 do j=1,jend-jsta+1
1837 jj = jsta+j-1
1838 do i=1,iend-ista+1
1839 ii=ista+i-1
1840 datapd(i,j,cfld) = grid1(ii,jj)
1841 enddo
1842 enddo
1843 endif
1844 ENDIF
1845 ENDIF
1846!
1847!--- RIME FACTOR
1848 IF (iget(263) > 0) THEN
1849 IF (lvls(lp,iget(263)) > 0) THEN
1850!$omp parallel do private(i,j)
1851 DO j=jsta,jend
1852 DO i=ista,iend
1853 grid1(i,j) = frime(i,j)
1854 ENDDO
1855 ENDDO
1856 if(grib == 'grib2')then
1857 cfld = cfld + 1
1858 fld_info(cfld)%ifld=iavblfld(iget(263))
1859 fld_info(cfld)%lvl=lvlsxml(lp,iget(263))
1860!$omp parallel do private(i,j,ii,jj)
1861 do j=1,jend-jsta+1
1862 jj = jsta+j-1
1863 do i=1,iend-ista+1
1864 ii=ista+i-1
1865 datapd(i,j,cfld) = grid1(ii,jj)
1866 enddo
1867 enddo
1868 endif
1869 ENDIF
1870 ENDIF
1871!
1872!--- Temperature tendency by all radiation: == ested by AFWA
1873 IF (iget(294) > 0) THEN
1874 IF (lvls(lp,iget(294)) > 0) THEN
1875!$omp parallel do private(i,j)
1876 DO j=jsta,jend
1877 DO i=ista,iend
1878 grid1(i,j) = rad(i,j)
1879 ENDDO
1880 ENDDO
1881 if(grib == 'grib2')then
1882 cfld = cfld + 1
1883 fld_info(cfld)%ifld=iavblfld(iget(294))
1884 fld_info(cfld)%lvl=lvlsxml(lp,iget(294))
1885!$omp parallel do private(i,j,ii,jj)
1886 do j=1,jend-jsta+1
1887 jj = jsta+j-1
1888 do i=1,iend-ista+1
1889 ii=ista+i-1
1890 datapd(i,j,cfld) = grid1(ii,jj)
1891 enddo
1892 enddo
1893 endif
1894 ENDIF
1895 ENDIF
1896!
1897!--- Radar Reflectivity
1898 IF (iget(251) > 0) THEN
1899 IF (lvls(lp,iget(251)) > 0) THEN
1900!$omp parallel do private(i,j)
1901 DO j=jsta,jend
1902 DO i=ista,iend
1903 grid1(i,j) = dbz1(i,j)
1904 ENDDO
1905 ENDDO
1906 if(grib == 'grib2')then
1907 cfld = cfld + 1
1908 fld_info(cfld)%ifld=iavblfld(iget(251))
1909 fld_info(cfld)%lvl=lvlsxml(lp,iget(251))
1910!$omp parallel do private(i,j,ii,jj)
1911 do j=1,jend-jsta+1
1912 jj = jsta+j-1
1913 do i=1,iend-ista+1
1914 ii=ista+i-1
1915 datapd(i,j,cfld) = grid1(ii,jj)
1916 enddo
1917 enddo
1918 endif
1919 ENDIF
1920 ENDIF
1921!
1922!--- IN-FLIGHT ICING CONDITION: ADD BY B. ZHOU
1923 IF(iget(257) > 0)THEN
1924 IF(lvls(lp,iget(257)) > 0)THEN
1925 CALL calicing(tsl(ista:iend,jsta:jend), savrh, osl(ista:iend,jsta:jend), egrid1(ista:iend,jsta:jend))
1926
1927!$omp parallel do private(i,j)
1928 DO j=jsta,jend
1929 DO i=ista,iend
1930 grid1(i,j) = egrid1(i,j)
1931 ENDDO
1932 ENDDO
1933 if(grib == 'grib2')then
1934 cfld = cfld + 1
1935 fld_info(cfld)%ifld=iavblfld(iget(257))
1936 fld_info(cfld)%lvl=lvlsxml(lp,iget(257))
1937!$omp parallel do private(i,j,ii,jj)
1938 do j=1,jend-jsta+1
1939 jj = jsta+j-1
1940 do i=1,iend-ista+1
1941 ii=ista+i-1
1942 datapd(i,j,cfld) = grid1(ii,jj)
1943 enddo
1944 enddo
1945 endif
1946 ENDIF
1947 ENDIF
1948
1949!
1950
1951!--- CLEAR AIR TURBULENCE (CAT): ADD BY B. ZHOU
1952 IF (lp > 1) THEN
1953 IF(iget(258) > 0)THEN
1954 IF(lvls(lp,iget(258)) > 0)THEN
1955!$omp parallel do private(i,j)
1956 DO j=jsta,jend
1957 DO i=ista,iend
1958 IF(fsl(i,j)<spval)THEN
1959 grid1(i,j) = fsl(i,j)*gi
1960 ELSE
1961 grid1(i,j) = spval
1962 ENDIF
1963 egrid1(i,j) = spval
1964 ENDDO
1965 ENDDO
1966 CALL calcat(usl(ista_2l,jsta_2l),vsl(ista_2l,jsta_2l),grid1(ista_2l,jsta_2l) &
1967 ,usl_old(ista_2l,jsta_2l),vsl_old(ista_2l,jsta_2l) &
1968 ,fsl_old(ista_2l,jsta_2l),egrid1(ista_2l,jsta_2l))
1969!$omp parallel do private(i,j)
1970 DO j=jsta,jend
1971 DO i=ista,iend
1972 grid1(i,j) = egrid1(i,j)
1973! IF(GRID1(I,J) > 3. .OR. GRID1(I,J) < 0.)
1974! + print*,'bad CAT',i,j,GRID1(I,J)
1975 ENDDO
1976 ENDDO
1977 if(grib == 'grib2')then
1978 cfld = cfld + 1
1979 fld_info(cfld)%ifld=iavblfld(iget(258))
1980 fld_info(cfld)%lvl=lvlsxml(lp,iget(258))
1981!$omp parallel do private(i,j,ii,jj)
1982 do j=1,jend-jsta+1
1983 jj = jsta+j-1
1984 do i=1,iend-ista+1
1985 ii=ista+i-1
1986 datapd(i,j,cfld) = grid1(ii,jj)
1987 enddo
1988 enddo
1989 endif
1990 end if
1991 end if
1992 end if
1993!
1994
1995!$omp parallel do private(i,j)
1996 DO j=jsta_2l,jend_2u
1997 DO i=ista_2l,iend_2u
1998 usl_old(i,j) = usl(i,j)
1999 vsl_old(i,j) = vsl(i,j)
2000 IF(fsl(i,j)<spval)THEN
2001 fsl_old(i,j) = fsl(i,j)*gi
2002 ELSE
2003 fsl_old(i,j) = spval
2004 ENDIF
2005 ENDDO
2006 ENDDO
2007!
2008!--- OZONE
2009 IF (iget(268) > 0) THEN
2010 IF (lvls(lp,iget(268)) > 0) THEN
2011!$omp parallel do private(i,j)
2012 DO j=jsta,jend
2013 DO i=ista,iend
2014 grid1(i,j) = o3sl(i,j)
2015 ENDDO
2016 ENDDO
2017! print *,'in mdl2p,o3sl=',minval(o3sl(1:im,jsta:jend)), &
2018! minval(o3sl(1:im,jsta:jend))
2019 if(grib == 'grib2')then
2020 cfld = cfld + 1
2021 fld_info(cfld)%ifld=iavblfld(iget(268))
2022 fld_info(cfld)%lvl=lvlsxml(lp,iget(268))
2023!$omp parallel do private(i,j,ii,jj)
2024 do j=1,jend-jsta+1
2025 jj = jsta+j-1
2026 do i=1,iend-ista+1
2027 ii=ista+i-1
2028 datapd(i,j,cfld) = grid1(ii,jj)
2029 enddo
2030 enddo
2031 endif
2032 ENDIF
2033 ENDIF
2034! E. James - 8 Dec 2017: SMOKE from WRF-CHEM
2035!--- SMOKE
2036 IF (iget(738) > 0) THEN
2037 IF (lvls(lp,iget(738)) > 0) THEN
2038!$omp parallel do private(i,j)
2039 DO j=jsta,jend
2040 DO i=ista,iend
2041 IF(smokesl(i,j,1)<spval.and.spl(lp)<spval.and.tsl(i,j)<spval)THEN
2042 grid1(i,j) = (1./rd)*smokesl(i,j,1)*(spl(lp)/tsl(i,j))
2043 ELSE
2044 grid1(i,j) = spval
2045 ENDIF
2046 ENDDO
2047 ENDDO
2048 if(grib == 'grib2')then
2049 cfld = cfld + 1
2050 fld_info(cfld)%ifld=iavblfld(iget(738))
2051 fld_info(cfld)%lvl=lvlsxml(lp,iget(738))
2052!$omp parallel do private(i,j,ii,jj)
2053 do j=1,jend-jsta+1
2054 jj = jsta+j-1
2055 do i=1,iend-ista+1
2056 ii=ista+i-1
2057 datapd(i,j,cfld) = grid1(ii,jj)
2058 enddo
2059 enddo
2060 endif
2061 ENDIF
2062 ENDIF
2063 if (gocart_on) then
2064!--- DUST
2065 IF (iget(438) > 0) THEN
2066 IF (lvls(lp,iget(438)) > 0) THEN
2067!$omp parallel do private(i,j)
2068 DO j=jsta,jend
2069 DO i=ista,iend
2070 grid1(i,j) = dustsl(i,j,1)
2071 ENDDO
2072 ENDDO
2073 if(grib == 'grib2')then
2074 cfld = cfld + 1
2075 fld_info(cfld)%ifld=iavblfld(iget(438))
2076 fld_info(cfld)%lvl=lvlsxml(lp,iget(438))
2077!$omp parallel do private(i,j,ii,jj)
2078 do j=1,jend-jsta+1
2079 jj = jsta+j-1
2080 do i=1,iend-ista+1
2081 ii=ista+i-1
2082 datapd(i,j,cfld) = grid1(ii,jj)
2083 enddo
2084 enddo
2085 endif
2086 ENDIF
2087 ENDIF
2088
2089 IF (iget(439) > 0) THEN
2090 IF (lvls(lp,iget(439)) > 0) THEN
2091!$omp parallel do private(i,j)
2092 DO j=jsta,jend
2093 DO i=ista,iend
2094 grid1(i,j) = dustsl(i,j,2)
2095 ENDDO
2096 ENDDO
2097 if(grib == 'grib2')then
2098 cfld = cfld + 1
2099 fld_info(cfld)%ifld=iavblfld(iget(439))
2100 fld_info(cfld)%lvl=lvlsxml(lp,iget(439))
2101!$omp parallel do private(i,j,ii,jj)
2102 do j=1,jend-jsta+1
2103 jj = jsta+j-1
2104 do i=1,iend-ista+1
2105 ii=ista+i-1
2106 datapd(i,j,cfld) = grid1(ii,jj)
2107 enddo
2108 enddo
2109 endif
2110 ENDIF
2111 ENDIF
2112
2113 IF (iget(440) > 0) THEN
2114 IF (lvls(lp,iget(440)) > 0) THEN
2115!$omp parallel do private(i,j)
2116 DO j=jsta,jend
2117 DO i=ista,iend
2118 grid1(i,j) = dustsl(i,j,3)
2119 ENDDO
2120 ENDDO
2121 if(grib == 'grib2')then
2122 cfld = cfld + 1
2123 fld_info(cfld)%ifld=iavblfld(iget(440))
2124 fld_info(cfld)%lvl=lvlsxml(lp,iget(440))
2125!$omp parallel do private(i,j,ii,jj)
2126 do j=1,jend-jsta+1
2127 jj = jsta+j-1
2128 do i=1,iend-ista+1
2129 ii=ista+i-1
2130 datapd(i,j,cfld) = grid1(ii,jj)
2131 enddo
2132 enddo
2133 endif
2134 ENDIF
2135 ENDIF
2136
2137 IF (iget(441) > 0) THEN
2138 IF (lvls(lp,iget(441)) > 0) THEN
2139!$omp parallel do private(i,j)
2140 DO j=jsta,jend
2141 DO i=ista,iend
2142 grid1(i,j) = dustsl(i,j,4)
2143 ENDDO
2144 ENDDO
2145 if(grib == 'grib2')then
2146 cfld = cfld + 1
2147 fld_info(cfld)%ifld=iavblfld(iget(441))
2148 fld_info(cfld)%lvl=lvlsxml(lp,iget(441))
2149!$omp parallel do private(i,j,ii,jj)
2150 do j=1,jend-jsta+1
2151 jj = jsta+j-1
2152 do i=1,iend-ista+1
2153 ii=ista+i-1
2154 datapd(i,j,cfld) = grid1(ii,jj)
2155 enddo
2156 enddo
2157 endif
2158 ENDIF
2159 ENDIF
2160
2161 IF (iget(442) > 0) THEN
2162 IF (lvls(lp,iget(442)) > 0) THEN
2163!$omp parallel do private(i,j)
2164 DO j=jsta,jend
2165 DO i=ista,iend
2166 grid1(i,j) = dustsl(i,j,5)
2167 ENDDO
2168 ENDDO
2169 if(grib == 'grib2')then
2170 cfld = cfld + 1
2171 fld_info(cfld)%ifld=iavblfld(iget(442))
2172 fld_info(cfld)%lvl=lvlsxml(lp,iget(442))
2173!$omp parallel do private(i,j,ii,jj)
2174 do j=1,jend-jsta+1
2175 jj = jsta+j-1
2176 do i=1,iend-ista+1
2177 ii=ista+i-1
2178 datapd(i,j,cfld) = grid1(ii,jj)
2179 enddo
2180 enddo
2181 endif
2182 ENDIF
2183 ENDIF
2184 endif ! if gocart_on
2185
2186
2187 if(iostatusd3d==0 .and. d3d_on) then
2188!--- longwave tendency
2189 IF (iget(355) > 0) THEN
2190 IF (lvls(lp,iget(355)) > 0) THEN
2191!$omp parallel do private(i,j)
2192 DO j=jsta,jend
2193 DO i=ista,iend
2194 grid1(i,j) = d3dsl(i,j,1)
2195 ENDDO
2196 ENDDO
2197 id(1:25)=0
2198 itd3d = nint(td3d)
2199 if (itd3d /= 0) then
2200 ifincr = mod(ifhr,itd3d)
2201 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2202 else
2203 ifincr = 0
2204 endif
2205 id(18) = 0
2206 id(19) = ifhr
2207 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2208 id(20) = 3
2209 IF (ifincr == 0) THEN
2210 id(18) = ifhr-itd3d
2211 ELSE
2212 id(18) = ifhr-ifincr
2213 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2214 ENDIF
2215 if(grib == 'grib2')then
2216 cfld = cfld + 1
2217 fld_info(cfld)%ifld=iavblfld(iget(355))
2218 fld_info(cfld)%lvl=lvlsxml(lp,iget(355))
2219 if(itd3d==0) then
2220 fld_info(cfld)%ntrange=0
2221 else
2222 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2223 endif
2224 fld_info(cfld)%tinvstat=itd3d
2225!$omp parallel do private(i,j,ii,jj)
2226 do j=1,jend-jsta+1
2227 jj = jsta+j-1
2228 do i=1,iend-ista+1
2229 ii=ista+i-1
2230 datapd(i,j,cfld) = grid1(ii,jj)
2231 enddo
2232 enddo
2233 endif
2234 ENDIF
2235 ENDIF
2236!--- longwave tendency
2237 IF (iget(354) > 0) THEN
2238 IF (lvls(lp,iget(354)) > 0) THEN
2239!$omp parallel do private(i,j)
2240 DO j=jsta,jend
2241 DO i=ista,iend
2242 grid1(i,j) = d3dsl(i,j,2)
2243 ENDDO
2244 ENDDO
2245 id(1:25)=0
2246 itd3d = nint(td3d)
2247 if (itd3d /= 0) then
2248 ifincr = mod(ifhr,itd3d)
2249 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2250 else
2251 ifincr = 0
2252 endif
2253 id(18) = 0
2254 id(19) = ifhr
2255 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2256 id(20) = 3
2257 IF (ifincr == 0) THEN
2258 id(18) = ifhr-itd3d
2259 ELSE
2260 id(18) = ifhr-ifincr
2261 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2262 ENDIF
2263 if(grib == 'grib2')then
2264 cfld = cfld + 1
2265 fld_info(cfld)%ifld=iavblfld(iget(354))
2266 fld_info(cfld)%lvl=lvlsxml(lp,iget(354))
2267 if(itd3d==0) then
2268 fld_info(cfld)%ntrange=0
2269 else
2270 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2271 endif
2272 fld_info(cfld)%tinvstat=itd3d
2273!$omp parallel do private(i,j,ii,jj)
2274 do j=1,jend-jsta+1
2275 jj = jsta+j-1
2276 do i=1,iend-ista+1
2277 ii=ista+i-1
2278 datapd(i,j,cfld) = grid1(ii,jj)
2279 enddo
2280 enddo
2281 endif
2282 ENDIF
2283 ENDIF
2284!--- longwave tendency
2285 IF (iget(356) > 0) THEN
2286 IF (lvls(lp,iget(356)) > 0) THEN
2287!$omp parallel do private(i,j)
2288 DO j=jsta,jend
2289 DO i=ista,iend
2290 grid1(i,j) = d3dsl(i,j,3)
2291 ENDDO
2292 ENDDO
2293 id(1:25)=0
2294 itd3d = nint(td3d)
2295 if (itd3d /= 0) then
2296 ifincr = mod(ifhr,itd3d)
2297 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2298 else
2299 ifincr = 0
2300 endif
2301 id(18) = 0
2302 id(19) = ifhr
2303 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2304 id(20) = 3
2305 IF (ifincr == 0) THEN
2306 id(18) = ifhr-itd3d
2307 ELSE
2308 id(18) = ifhr-ifincr
2309 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2310 ENDIF
2311 if(grib == 'grib2')then
2312 cfld = cfld + 1
2313 fld_info(cfld)%ifld=iavblfld(iget(356))
2314 fld_info(cfld)%lvl=lvlsxml(lp,iget(356))
2315 if(itd3d==0) then
2316 fld_info(cfld)%ntrange=0
2317 else
2318 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2319 endif
2320 fld_info(cfld)%tinvstat=itd3d
2321!$omp parallel do private(i,j,ii,jj)
2322 do j=1,jend-jsta+1
2323 jj = jsta+j-1
2324 do i=1,iend-ista+1
2325 ii=ista+i-1
2326 datapd(i,j,cfld) = grid1(ii,jj)
2327 enddo
2328 enddo
2329 endif
2330 ENDIF
2331 ENDIF
2332!--- longwave tendency
2333 IF (iget(357) > 0) THEN
2334 IF (lvls(lp,iget(357)) > 0) THEN
2335!$omp parallel do private(i,j)
2336 DO j=jsta,jend
2337 DO i=ista,iend
2338 grid1(i,j) = d3dsl(i,j,4)
2339 ENDDO
2340 ENDDO
2341 id(1:25)=0
2342 itd3d = nint(td3d)
2343 if (itd3d /= 0) then
2344 ifincr = mod(ifhr,itd3d)
2345 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2346 else
2347 ifincr = 0
2348 endif
2349 id(18) = 0
2350 id(19) = ifhr
2351 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2352 id(20) = 3
2353 IF (ifincr == 0) THEN
2354 id(18) = ifhr-itd3d
2355 ELSE
2356 id(18) = ifhr-ifincr
2357 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2358 ENDIF
2359 if(grib == 'grib2')then
2360 cfld = cfld + 1
2361 fld_info(cfld)%ifld=iavblfld(iget(357))
2362 fld_info(cfld)%lvl=lvlsxml(lp,iget(357))
2363 if(itd3d==0) then
2364 fld_info(cfld)%ntrange=0
2365 else
2366 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2367 endif
2368 fld_info(cfld)%tinvstat=itd3d
2369!$omp parallel do private(i,j,ii,jj)
2370 do j=1,jend-jsta+1
2371 jj = jsta+j-1
2372 do i=1,iend-ista+1
2373 ii=ista+i-1
2374 datapd(i,j,cfld) = grid1(ii,jj)
2375 enddo
2376 enddo
2377 endif
2378 ENDIF
2379 ENDIF
2380!--- longwave tendency
2381 IF (iget(358) > 0) THEN
2382 IF (lvls(lp,iget(358)) > 0) THEN
2383!$omp parallel do private(i,j)
2384 DO j=jsta,jend
2385 DO i=ista,iend
2386 grid1(i,j) = d3dsl(i,j,5)
2387 ENDDO
2388 ENDDO
2389 id(1:25)=0
2390 itd3d = nint(td3d)
2391 if (itd3d /= 0) then
2392 ifincr = mod(ifhr,itd3d)
2393 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2394 else
2395 ifincr = 0
2396 endif
2397 id(18) = 0
2398 id(19) = ifhr
2399 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2400 id(20) = 3
2401 IF (ifincr == 0) THEN
2402 id(18) = ifhr-itd3d
2403 ELSE
2404 id(18) = ifhr-ifincr
2405 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2406 ENDIF
2407 if(grib == 'grib2')then
2408 cfld = cfld + 1
2409 fld_info(cfld)%ifld=iavblfld(iget(358))
2410 fld_info(cfld)%lvl=lvlsxml(lp,iget(358))
2411 if(itd3d==0) then
2412 fld_info(cfld)%ntrange=0
2413 else
2414 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2415 endif
2416 fld_info(cfld)%tinvstat=itd3d
2417!$omp parallel do private(i,j,ii,jj)
2418 do j=1,jend-jsta+1
2419 jj = jsta+j-1
2420 do i=1,iend-ista+1
2421 ii=ista+i-1
2422 datapd(i,j,cfld) = grid1(ii,jj)
2423 enddo
2424 enddo
2425 endif
2426 ENDIF
2427 ENDIF
2428!--- longwave tendency
2429 IF (iget(359) > 0) THEN
2430 IF (lvls(lp,iget(359)) > 0) THEN
2431!$omp parallel do private(i,j)
2432 DO j=jsta,jend
2433 DO i=ista,iend
2434 grid1(i,j) = d3dsl(i,j,6)
2435 ENDDO
2436 ENDDO
2437 id(1:25)=0
2438 itd3d = nint(td3d)
2439 if (itd3d /= 0) then
2440 ifincr = mod(ifhr,itd3d)
2441 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2442 else
2443 ifincr = 0
2444 endif
2445 id(18) = 0
2446 id(19) = ifhr
2447 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2448 id(20) = 3
2449 IF (ifincr == 0) THEN
2450 id(18) = ifhr-itd3d
2451 ELSE
2452 id(18) = ifhr-ifincr
2453 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2454 ENDIF
2455 if(grib == 'grib2')then
2456 cfld = cfld + 1
2457 fld_info(cfld)%ifld=iavblfld(iget(359))
2458 fld_info(cfld)%lvl=lvlsxml(lp,iget(359))
2459 if(itd3d==0) then
2460 fld_info(cfld)%ntrange=0
2461 else
2462 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2463 endif
2464 fld_info(cfld)%tinvstat=itd3d
2465!$omp parallel do private(i,j,ii,jj)
2466 do j=1,jend-jsta+1
2467 jj = jsta+j-1
2468 do i=1,iend-ista+1
2469 ii=ista+i-1
2470 datapd(i,j,cfld) = grid1(ii,jj)
2471 enddo
2472 enddo
2473 endif
2474 ENDIF
2475 ENDIF
2476!--- longwave tendency
2477 IF (iget(360) > 0) THEN
2478 IF (lvls(lp,iget(360)) > 0) THEN
2479!$omp parallel do private(i,j)
2480 DO j=jsta,jend
2481 DO i=ista,iend
2482 grid1(i,j) = d3dsl(i,j,7)
2483 ENDDO
2484 ENDDO
2485 id(1:25)=0
2486 itd3d = nint(td3d)
2487 if (itd3d /= 0) then
2488 ifincr = mod(ifhr,itd3d)
2489 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2490 else
2491 ifincr = 0
2492 endif
2493 id(18) = 0
2494 id(19) = ifhr
2495 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2496 id(20) = 3
2497 IF (ifincr == 0) THEN
2498 id(18) = ifhr-itd3d
2499 ELSE
2500 id(18) = ifhr-ifincr
2501 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2502 ENDIF
2503 if(grib == 'grib2')then
2504 cfld = cfld + 1
2505 fld_info(cfld)%ifld=iavblfld(iget(360))
2506 fld_info(cfld)%lvl=lvlsxml(lp,iget(360))
2507 if(itd3d==0) then
2508 fld_info(cfld)%ntrange=0
2509 else
2510 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2511 endif
2512 fld_info(cfld)%tinvstat=itd3d
2513!$omp parallel do private(i,j,ii,jj)
2514 do j=1,jend-jsta+1
2515 jj = jsta+j-1
2516 do i=1,iend-ista+1
2517 ii=ista+i-1
2518 datapd(i,j,cfld) = grid1(ii,jj)
2519 enddo
2520 enddo
2521 endif
2522 ENDIF
2523 ENDIF
2524!--- longwave tendency
2525 IF (iget(361) > 0) THEN
2526 IF (lvls(lp,iget(361)) > 0) THEN
2527!$omp parallel do private(i,j)
2528 DO j=jsta,jend
2529 DO i=ista,iend
2530 grid1(i,j) = d3dsl(i,j,8)
2531 ENDDO
2532 ENDDO
2533 id(1:25)=0
2534 itd3d = nint(td3d)
2535 if (itd3d /= 0) then
2536 ifincr = mod(ifhr,itd3d)
2537 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2538 else
2539 ifincr = 0
2540 endif
2541 id(18) = 0
2542 id(19) = ifhr
2543 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2544 id(20) = 3
2545 IF (ifincr == 0) THEN
2546 id(18) = ifhr-itd3d
2547 ELSE
2548 id(18) = ifhr-ifincr
2549 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2550 ENDIF
2551 if(grib == 'grib2')then
2552 cfld = cfld + 1
2553 fld_info(cfld)%ifld=iavblfld(iget(361))
2554 fld_info(cfld)%lvl=lvlsxml(lp,iget(361))
2555 if(itd3d==0) then
2556 fld_info(cfld)%ntrange=0
2557 else
2558 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2559 endif
2560 fld_info(cfld)%tinvstat=itd3d
2561!$omp parallel do private(i,j,ii,jj)
2562 do j=1,jend-jsta+1
2563 jj = jsta+j-1
2564 do i=1,iend-ista+1
2565 ii=ista+i-1
2566 datapd(i,j,cfld) = grid1(ii,jj)
2567 enddo
2568 enddo
2569 endif
2570 ENDIF
2571 ENDIF
2572!--- longwave tendency
2573 IF (iget(362) > 0) THEN
2574 IF (lvls(lp,iget(362)) > 0) THEN
2575!$omp parallel do private(i,j)
2576 DO j=jsta,jend
2577 DO i=ista,iend
2578 grid1(i,j) = d3dsl(i,j,9)
2579 ENDDO
2580 ENDDO
2581 id(1:25)=0
2582 itd3d = nint(td3d)
2583 if (itd3d /= 0) then
2584 ifincr = mod(ifhr,itd3d)
2585 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2586 else
2587 ifincr = 0
2588 endif
2589 id(18) = 0
2590 id(19) = ifhr
2591 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2592 id(20) = 3
2593 IF (ifincr == 0) THEN
2594 id(18) = ifhr-itd3d
2595 ELSE
2596 id(18) = ifhr-ifincr
2597 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2598 ENDIF
2599 if(grib == 'grib2')then
2600 cfld = cfld + 1
2601 fld_info(cfld)%ifld=iavblfld(iget(362))
2602 fld_info(cfld)%lvl=lvlsxml(lp,iget(362))
2603 if(itd3d==0) then
2604 fld_info(cfld)%ntrange=0
2605 else
2606 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2607 endif
2608 fld_info(cfld)%tinvstat=itd3d
2609!$omp parallel do private(i,j,ii,jj)
2610 do j=1,jend-jsta+1
2611 jj = jsta+j-1
2612 do i=1,iend-ista+1
2613 ii=ista+i-1
2614 datapd(i,j,cfld) = grid1(ii,jj)
2615 enddo
2616 enddo
2617 endif
2618 ENDIF
2619 ENDIF
2620!--- longwave tendency
2621 IF (iget(363) > 0) THEN
2622 IF (lvls(lp,iget(363)) > 0) THEN
2623!$omp parallel do private(i,j)
2624 DO j=jsta,jend
2625 DO i=ista,iend
2626 grid1(i,j) = d3dsl(i,j,10)
2627 ENDDO
2628 ENDDO
2629 id(1:25)=0
2630 itd3d = nint(td3d)
2631 if (itd3d /= 0) then
2632 ifincr = mod(ifhr,itd3d)
2633 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2634 else
2635 ifincr = 0
2636 endif
2637 id(02)=133 ! Table 133
2638 id(18) = 0
2639 id(19) = ifhr
2640 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2641 id(20) = 3
2642 IF (ifincr == 0) THEN
2643 id(18) = ifhr-itd3d
2644 ELSE
2645 id(18) = ifhr-ifincr
2646 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2647 ENDIF
2648 if(grib == 'grib2')then
2649 cfld = cfld + 1
2650 fld_info(cfld)%ifld=iavblfld(iget(363))
2651 fld_info(cfld)%lvl=lvlsxml(lp,iget(363))
2652 if(itd3d==0) then
2653 fld_info(cfld)%ntrange=0
2654 else
2655 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2656 endif
2657 fld_info(cfld)%tinvstat=itd3d
2658!$omp parallel do private(i,j,ii,jj)
2659 do j=1,jend-jsta+1
2660 jj = jsta+j-1
2661 do i=1,iend-ista+1
2662 ii=ista+i-1
2663 datapd(i,j,cfld) = grid1(ii,jj)
2664 enddo
2665 enddo
2666 endif
2667 ENDIF
2668 ENDIF
2669!--- longwave tendency
2670 IF (iget(364) > 0) THEN
2671 IF (lvls(lp,iget(364)) > 0) THEN
2672!$omp parallel do private(i,j)
2673 DO j=jsta,jend
2674 DO i=ista,iend
2675 grid1(i,j) = d3dsl(i,j,11)
2676 ENDDO
2677 ENDDO
2678 id(1:25)=0
2679 itd3d = nint(td3d)
2680 if (itd3d /= 0) then
2681 ifincr = mod(ifhr,itd3d)
2682 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2683 else
2684 ifincr = 0
2685 endif
2686 id(02)=133 ! Table 133
2687 id(18) = 0
2688 id(19) = ifhr
2689 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2690 id(20) = 3
2691 IF (ifincr == 0) THEN
2692 id(18) = ifhr-itd3d
2693 ELSE
2694 id(18) = ifhr-ifincr
2695 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2696 ENDIF
2697 if(grib == 'grib2')then
2698 cfld = cfld + 1
2699 fld_info(cfld)%ifld=iavblfld(iget(364))
2700 fld_info(cfld)%lvl=lvlsxml(lp,iget(364))
2701 if(itd3d==0) then
2702 fld_info(cfld)%ntrange=0
2703 else
2704 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2705 endif
2706 fld_info(cfld)%tinvstat=itd3d
2707!$omp parallel do private(i,j,ii,jj)
2708 do j=1,jend-jsta+1
2709 jj = jsta+j-1
2710 do i=1,iend-ista+1
2711 ii=ista+i-1
2712 datapd(i,j,cfld) = grid1(ii,jj)
2713 enddo
2714 enddo
2715 endif
2716 ENDIF
2717 ENDIF
2718!--- longwave tendency
2719 IF (iget(365) > 0) THEN
2720 IF (lvls(lp,iget(365)) > 0) THEN
2721!$omp parallel do private(i,j)
2722 DO j=jsta,jend
2723 DO i=ista,iend
2724 grid1(i,j) = d3dsl(i,j,12)
2725 ENDDO
2726 ENDDO
2727 id(1:25)=0
2728 itd3d = nint(td3d)
2729 if (itd3d /= 0) then
2730 ifincr = mod(ifhr,itd3d)
2731 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2732 else
2733 ifincr = 0
2734 endif
2735 id(02)=133 ! Table 133
2736 id(18) = 0
2737 id(19) = ifhr
2738 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2739 id(20) = 3
2740 IF (ifincr == 0) THEN
2741 id(18) = ifhr-itd3d
2742 ELSE
2743 id(18) = ifhr-ifincr
2744 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2745 ENDIF
2746 if(grib == 'grib2')then
2747 cfld = cfld + 1
2748 fld_info(cfld)%ifld=iavblfld(iget(365))
2749 fld_info(cfld)%lvl=lvlsxml(lp,iget(365))
2750 if(itd3d==0) then
2751 fld_info(cfld)%ntrange=0
2752 else
2753 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2754 endif
2755 fld_info(cfld)%tinvstat=itd3d
2756!$omp parallel do private(i,j,ii,jj)
2757 do j=1,jend-jsta+1
2758 jj = jsta+j-1
2759 do i=1,iend-ista+1
2760 ii=ista+i-1
2761 datapd(i,j,cfld) = grid1(ii,jj)
2762 enddo
2763 enddo
2764 endif
2765 ENDIF
2766 ENDIF
2767!--- longwave tendency
2768 IF (iget(366) > 0) THEN
2769 IF (lvls(lp,iget(366)) > 0) THEN
2770!$omp parallel do private(i,j)
2771 DO j=jsta,jend
2772 DO i=ista,iend
2773 grid1(i,j) = d3dsl(i,j,13)
2774 ENDDO
2775 ENDDO
2776 id(1:25)=0
2777 itd3d = nint(td3d)
2778 if (itd3d /= 0) then
2779 ifincr = mod(ifhr,itd3d)
2780 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2781 else
2782 ifincr = 0
2783 endif
2784 id(02)=133 ! Table 133
2785 id(18) = 0
2786 id(19) = ifhr
2787 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2788 id(20) = 3
2789 IF (ifincr == 0) THEN
2790 id(18) = ifhr-itd3d
2791 ELSE
2792 id(18) = ifhr-ifincr
2793 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2794 ENDIF
2795 if(grib == 'grib2')then
2796 cfld = cfld + 1
2797 fld_info(cfld)%ifld=iavblfld(iget(366))
2798 fld_info(cfld)%lvl=lvlsxml(lp,iget(366))
2799 if(itd3d==0) then
2800 fld_info(cfld)%ntrange=0
2801 else
2802 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2803 endif
2804 fld_info(cfld)%tinvstat=itd3d
2805!$omp parallel do private(i,j,ii,jj)
2806 do j=1,jend-jsta+1
2807 jj = jsta+j-1
2808 do i=1,iend-ista+1
2809 ii=ista+i-1
2810 datapd(i,j,cfld) = grid1(ii,jj)
2811 enddo
2812 enddo
2813 endif
2814 ENDIF
2815 ENDIF
2816!--- longwave tendency
2817 IF (iget(367) > 0) THEN
2818 IF (lvls(lp,iget(367)) > 0) THEN
2819!$omp parallel do private(i,j)
2820 DO j=jsta,jend
2821 DO i=ista,iend
2822 grid1(i,j) = d3dsl(i,j,14)
2823 ENDDO
2824 ENDDO
2825 id(1:25)=0
2826 itd3d = nint(td3d)
2827 if (itd3d /= 0) then
2828 ifincr = mod(ifhr,itd3d)
2829 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2830 else
2831 ifincr = 0
2832 endif
2833 id(02)=133 ! Table 133
2834 id(18) = 0
2835 id(19) = ifhr
2836 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2837 id(20) = 3
2838 IF (ifincr == 0) THEN
2839 id(18) = ifhr-itd3d
2840 ELSE
2841 id(18) = ifhr-ifincr
2842 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2843 ENDIF
2844 if(grib == 'grib2')then
2845 cfld = cfld + 1
2846 fld_info(cfld)%ifld=iavblfld(iget(367))
2847 fld_info(cfld)%lvl=lvlsxml(lp,iget(367))
2848 if(itd3d==0) then
2849 fld_info(cfld)%ntrange=0
2850 else
2851 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2852 endif
2853 fld_info(cfld)%tinvstat=itd3d
2854!$omp parallel do private(i,j,ii,jj)
2855 do j=1,jend-jsta+1
2856 jj = jsta+j-1
2857 do i=1,iend-ista+1
2858 ii=ista+i-1
2859 datapd(i,j,cfld) = grid1(ii,jj)
2860 enddo
2861 enddo
2862 endif
2863 ENDIF
2864 ENDIF
2865!--- longwave tendency
2866 IF (iget(368) > 0) THEN
2867 IF (lvls(lp,iget(368)) > 0) THEN
2868!$omp parallel do private(i,j)
2869 DO j=jsta,jend
2870 DO i=ista,iend
2871 grid1(i,j) = d3dsl(i,j,15)
2872 ENDDO
2873 ENDDO
2874 id(1:25)=0
2875 itd3d = nint(td3d)
2876 if (itd3d /= 0) then
2877 ifincr = mod(ifhr,itd3d)
2878 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2879 else
2880 ifincr = 0
2881 endif
2882 id(02)=133 ! Table 133
2883 id(18) = 0
2884 id(19) = ifhr
2885 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2886 id(20) = 3
2887 IF (ifincr == 0) THEN
2888 id(18) = ifhr-itd3d
2889 ELSE
2890 id(18) = ifhr-ifincr
2891 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2892 ENDIF
2893 if(grib == 'grib2')then
2894 cfld = cfld + 1
2895 fld_info(cfld)%ifld=iavblfld(iget(368))
2896 fld_info(cfld)%lvl=lvlsxml(lp,iget(368))
2897 if(itd3d==0) then
2898 fld_info(cfld)%ntrange=0
2899 else
2900 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2901 endif
2902 fld_info(cfld)%tinvstat=itd3d
2903!$omp parallel do private(i,j,ii,jj)
2904 do j=1,jend-jsta+1
2905 jj = jsta+j-1
2906 do i=1,iend-ista+1
2907 ii=ista+i-1
2908 datapd(i,j,cfld) = grid1(ii,jj)
2909 enddo
2910 enddo
2911 endif
2912 ENDIF
2913 ENDIF
2914!--- longwave tendency
2915 IF (iget(369) > 0) THEN
2916 IF (lvls(lp,iget(369)) > 0) THEN
2917!$omp parallel do private(i,j)
2918 DO j=jsta,jend
2919 DO i=ista,iend
2920 grid1(i,j) = d3dsl(i,j,16)
2921 ENDDO
2922 ENDDO
2923 id(1:25)=0
2924 itd3d = nint(td3d)
2925 if (itd3d /= 0) then
2926 ifincr = mod(ifhr,itd3d)
2927 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2928 else
2929 ifincr = 0
2930 endif
2931 id(18) = 0
2932 id(19) = ifhr
2933 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2934 id(20) = 3
2935 IF (ifincr == 0) THEN
2936 id(18) = ifhr-itd3d
2937 ELSE
2938 id(18) = ifhr-ifincr
2939 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2940 ENDIF
2941 if(grib == 'grib2')then
2942 cfld = cfld + 1
2943 fld_info(cfld)%ifld=iavblfld(iget(369))
2944 fld_info(cfld)%lvl=lvlsxml(lp,iget(369))
2945 if(itd3d==0) then
2946 fld_info(cfld)%ntrange=0
2947 else
2948 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2949 endif
2950 fld_info(cfld)%tinvstat=itd3d
2951!$omp parallel do private(i,j,ii,jj)
2952 do j=1,jend-jsta+1
2953 jj = jsta+j-1
2954 do i=1,iend-ista+1
2955 ii=ista+i-1
2956 datapd(i,j,cfld) = grid1(ii,jj)
2957 enddo
2958 enddo
2959 endif
2960 ENDIF
2961 ENDIF
2962!--- longwave tendency
2963 IF (iget(370) > 0) THEN
2964 IF (lvls(lp,iget(370)) > 0) THEN
2965!$omp parallel do private(i,j)
2966 DO j=jsta,jend
2967 DO i=ista,iend
2968 grid1(i,j) = d3dsl(i,j,17)
2969 ENDDO
2970 ENDDO
2971 id(1:25)=0
2972 itd3d = nint(td3d)
2973 if (itd3d /= 0) then
2974 ifincr = mod(ifhr,itd3d)
2975 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
2976 else
2977 ifincr = 0
2978 endif
2979 id(02)=133 ! Table 133
2980 id(18) = 0
2981 id(19) = ifhr
2982 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
2983 id(20) = 3
2984 IF (ifincr == 0) THEN
2985 id(18) = ifhr-itd3d
2986 ELSE
2987 id(18) = ifhr-ifincr
2988 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
2989 ENDIF
2990 if(grib == 'grib2')then
2991 cfld = cfld + 1
2992 fld_info(cfld)%ifld=iavblfld(iget(370))
2993 fld_info(cfld)%lvl=lvlsxml(lp,iget(370))
2994 if(itd3d==0) then
2995 fld_info(cfld)%ntrange=0
2996 else
2997 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
2998 endif
2999 fld_info(cfld)%tinvstat=itd3d
3000!$omp parallel do private(i,j,ii,jj)
3001 do j=1,jend-jsta+1
3002 jj = jsta+j-1
3003 do i=1,iend-ista+1
3004 ii=ista+i-1
3005 datapd(i,j,cfld) = grid1(ii,jj)
3006 enddo
3007 enddo
3008 endif
3009 ENDIF
3010 ENDIF
3011!--- longwave tendency
3012 IF (iget(371) > 0) THEN
3013 IF (lvls(lp,iget(371)) > 0) THEN
3014!$omp parallel do private(i,j)
3015 DO j=jsta,jend
3016 DO i=ista,iend
3017 grid1(i,j) = d3dsl(i,j,18)
3018 ENDDO
3019 ENDDO
3020 id(1:25)=0
3021 itd3d = nint(td3d)
3022 if (itd3d /= 0) then
3023 ifincr = mod(ifhr,itd3d)
3024 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
3025 else
3026 ifincr = 0
3027 endif
3028 id(02)=133 ! Table 133
3029 id(18) = 0
3030 id(19) = ifhr
3031 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
3032 id(20) = 3
3033 IF (ifincr == 0) THEN
3034 id(18) = ifhr-itd3d
3035 ELSE
3036 id(18) = ifhr-ifincr
3037 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
3038 ENDIF
3039 if(grib == 'grib2')then
3040 cfld = cfld + 1
3041 fld_info(cfld)%ifld=iavblfld(iget(371))
3042 fld_info(cfld)%lvl=lvlsxml(lp,iget(371))
3043 if(itd3d==0) then
3044 fld_info(cfld)%ntrange=0
3045 else
3046 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
3047 endif
3048 fld_info(cfld)%tinvstat=itd3d
3049!$omp parallel do private(i,j,ii,jj)
3050 do j=1,jend-jsta+1
3051 jj = jsta+j-1
3052 do i=1,iend-ista+1
3053 ii=ista+i-1
3054 datapd(i,j,cfld) = grid1(ii,jj)
3055 enddo
3056 enddo
3057 endif
3058 ENDIF
3059 ENDIF
3060!--- longwave tendency
3061 IF (iget(372) > 0) THEN
3062 IF (lvls(lp,iget(372)) > 0) THEN
3063!$omp parallel do private(i,j)
3064 DO j=jsta,jend
3065 DO i=ista,iend
3066 grid1(i,j) = d3dsl(i,j,19)
3067 ENDDO
3068 ENDDO
3069 id(1:25)=0
3070 itd3d = nint(td3d)
3071 if (itd3d /= 0) then
3072 ifincr = mod(ifhr,itd3d)
3073 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
3074 else
3075 ifincr = 0
3076 endif
3077 id(18) = 0
3078 id(19) = ifhr
3079 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
3080 id(20) = 3
3081 IF (ifincr == 0) THEN
3082 id(18) = ifhr-itd3d
3083 ELSE
3084 id(18) = ifhr-ifincr
3085 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
3086 ENDIF
3087 if(grib == 'grib2')then
3088 cfld = cfld + 1
3089 fld_info(cfld)%ifld=iavblfld(iget(372))
3090 fld_info(cfld)%lvl=lvlsxml(lp,iget(372))
3091 if(itd3d==0) then
3092 fld_info(cfld)%ntrange=0
3093 else
3094 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
3095 endif
3096 fld_info(cfld)%tinvstat=itd3d
3097!$omp parallel do private(i,j,ii,jj)
3098 do j=1,jend-jsta+1
3099 jj = jsta+j-1
3100 do i=1,iend-ista+1
3101 ii=ista+i-1
3102 datapd(i,j,cfld) = grid1(ii,jj)
3103 enddo
3104 enddo
3105 endif
3106 ENDIF
3107 ENDIF
3108!--- longwave tendency
3109 IF (iget(373) > 0) THEN
3110 IF (lvls(lp,iget(373)) > 0) THEN
3111!$omp parallel do private(i,j)
3112 DO j=jsta,jend
3113 DO i=ista,iend
3114 grid1(i,j) = d3dsl(i,j,20)
3115 ENDDO
3116 ENDDO
3117 id(1:25)=0
3118 itd3d = nint(td3d)
3119 if (itd3d /= 0) then
3120 ifincr = mod(ifhr,itd3d)
3121 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
3122 else
3123 ifincr = 0
3124 endif
3125 id(02)=133 ! Table 133
3126 id(18) = 0
3127 id(19) = ifhr
3128 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
3129 id(20) = 3
3130 IF (ifincr == 0) THEN
3131 id(18) = ifhr-itd3d
3132 ELSE
3133 id(18) = ifhr-ifincr
3134 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
3135 ENDIF
3136 if(grib == 'grib2')then
3137 cfld = cfld + 1
3138 fld_info(cfld)%ifld=iavblfld(iget(373))
3139 fld_info(cfld)%lvl=lvlsxml(lp,iget(373))
3140 if(itd3d==0) then
3141 fld_info(cfld)%ntrange=0
3142 else
3143 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
3144 endif
3145 fld_info(cfld)%tinvstat=itd3d
3146!$omp parallel do private(i,j,ii,jj)
3147 do j=1,jend-jsta+1
3148 jj = jsta+j-1
3149 do i=1,iend-ista+1
3150 ii=ista+i-1
3151 datapd(i,j,cfld) = grid1(ii,jj)
3152 enddo
3153 enddo
3154 endif
3155 ENDIF
3156 ENDIF
3157!--- longwave tendency
3158 IF (iget(374) > 0) THEN
3159 IF (lvls(lp,iget(374)) > 0) THEN
3160!$omp parallel do private(i,j)
3161 DO j=jsta,jend
3162 DO i=ista,iend
3163 grid1(i,j) = d3dsl(i,j,21)
3164 ENDDO
3165 ENDDO
3166 id(1:25)=0
3167 itd3d = nint(td3d)
3168 if (itd3d /= 0) then
3169 ifincr = mod(ifhr,itd3d)
3170 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
3171 else
3172 ifincr = 0
3173 endif
3174 id(02)=133 ! Table 133
3175 id(18) = 0
3176 id(19) = ifhr
3177 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
3178 id(20) = 3
3179 IF (ifincr == 0) THEN
3180 id(18) = ifhr-itd3d
3181 ELSE
3182 id(18) = ifhr-ifincr
3183 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
3184 ENDIF
3185 if(grib == 'grib2')then
3186 cfld = cfld + 1
3187 fld_info(cfld)%ifld=iavblfld(iget(374))
3188 fld_info(cfld)%lvl=lvlsxml(lp,iget(374))
3189 if(itd3d==0) then
3190 fld_info(cfld)%ntrange=0
3191 else
3192 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
3193 endif
3194 fld_info(cfld)%tinvstat=itd3d
3195!$omp parallel do private(i,j,ii,jj)
3196 do j=1,jend-jsta+1
3197 jj = jsta+j-1
3198 do i=1,iend-ista+1
3199 ii=ista+i-1
3200 datapd(i,j,cfld) = grid1(ii,jj)
3201 enddo
3202 enddo
3203 endif
3204 ENDIF
3205 ENDIF
3206!--- longwave tendency
3207 IF (iget(375) > 0) THEN
3208 IF (lvls(lp,iget(375)) > 0) THEN
3209!$omp parallel do private(i,j)
3210 DO j=jsta,jend
3211 DO i=ista,iend
3212 grid1(i,j) = d3dsl(i,j,22)
3213 ENDDO
3214 ENDDO
3215 id(1:25)=0
3216 itd3d = nint(td3d)
3217 if (itd3d /= 0) then
3218 ifincr = mod(ifhr,itd3d)
3219 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
3220 else
3221 ifincr = 0
3222 endif
3223 id(18) = 0
3224 id(19) = ifhr
3225 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
3226 id(20) = 3
3227 IF (ifincr == 0) THEN
3228 id(18) = ifhr-itd3d
3229 ELSE
3230 id(18) = ifhr-ifincr
3231 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
3232 ENDIF
3233 if(grib == 'grib2') then
3234 cfld = cfld + 1
3235 fld_info(cfld)%ifld=iavblfld(iget(375))
3236 fld_info(cfld)%lvl=lvlsxml(lp,iget(375))
3237 if(itd3d==0) then
3238 fld_info(cfld)%ntrange=0
3239 else
3240 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
3241 endif
3242 fld_info(cfld)%tinvstat=itd3d
3243!$omp parallel do private(i,j,ii,jj)
3244 do j=1,jend-jsta+1
3245 jj = jsta+j-1
3246 do i=1,iend-ista+1
3247 ii=ista+i-1
3248 datapd(i,j,cfld) = grid1(ii,jj)
3249 enddo
3250 enddo
3251 endif
3252 ENDIF
3253 ENDIF
3254!--- total diabatic heating
3255 IF (iget(379) > 0) THEN
3256 IF (lvls(lp,iget(379)) > 0) THEN
3257!$omp parallel do private(i,j)
3258 DO j=jsta,jend
3259 DO i=ista,iend
3260 IF(d3dsl(i,j,1)/=spval)THEN
3261 grid1(i,j) = d3dsl(i,j,1) + d3dsl(i,j,2) &
3262 + d3dsl(i,j,3) + d3dsl(i,j,4) &
3263 + d3dsl(i,j,5) + d3dsl(i,j,6)
3264 ELSE
3265 grid1(i,j) = spval
3266 END IF
3267 ENDDO
3268 ENDDO
3269 id(1:25)=0
3270 itd3d = nint(td3d)
3271 if (itd3d /= 0) then
3272 ifincr = mod(ifhr,itd3d)
3273 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
3274 else
3275 ifincr = 0
3276 endif
3277 id(18) = 0
3278 id(19) = ifhr
3279 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
3280 id(20) = 3
3281 IF (ifincr == 0) THEN
3282 id(18) = ifhr-itd3d
3283 ELSE
3284 id(18) = ifhr-ifincr
3285 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
3286 ENDIF
3287 if(grib == 'grib2') then
3288 cfld = cfld + 1
3289 fld_info(cfld)%ifld=iavblfld(iget(379))
3290 fld_info(cfld)%lvl=lvlsxml(lp,iget(379))
3291 if(itd3d==0) then
3292 fld_info(cfld)%ntrange=0
3293 else
3294 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
3295 endif
3296 fld_info(cfld)%tinvstat=itd3d
3297!$omp parallel do private(i,j,ii,jj)
3298 do j=1,jend-jsta+1
3299 jj = jsta+j-1
3300 do i=1,iend-ista+1
3301 ii=ista+i-1
3302 datapd(i,j,cfld) = grid1(ii,jj)
3303 enddo
3304 enddo
3305 endif
3306 ENDIF
3307 ENDIF
3308!--- convective updraft
3309 IF (iget(391) > 0) THEN
3310 IF (lvls(lp,iget(391)) > 0) THEN
3311!$omp parallel do private(i,j)
3312 DO j=jsta,jend
3313 DO i=ista,iend
3314 grid1(i,j) = d3dsl(i,j,23)
3315 ENDDO
3316 ENDDO
3317 id(1:25)=0
3318 itd3d = nint(td3d)
3319 if (itd3d /= 0) then
3320 ifincr = mod(ifhr,itd3d)
3321 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
3322 else
3323 ifincr = 0
3324 endif
3325 id(02)=133 ! Table 133
3326 id(18) = 0
3327 id(19) = ifhr
3328 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
3329 id(20) = 3
3330 IF (ifincr == 0) THEN
3331 id(18) = ifhr-itd3d
3332 ELSE
3333 id(18) = ifhr-ifincr
3334 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
3335 ENDIF
3336 if(grib == 'grib2') then
3337 cfld = cfld + 1
3338 fld_info(cfld)%ifld=iavblfld(iget(391))
3339 fld_info(cfld)%lvl=lvlsxml(lp,iget(391))
3340 if(itd3d==0) then
3341 fld_info(cfld)%ntrange=0
3342 else
3343 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
3344 endif
3345 fld_info(cfld)%tinvstat=itd3d
3346!$omp parallel do private(i,j,ii,jj)
3347 do j=1,jend-jsta+1
3348 jj = jsta+j-1
3349 do i=1,iend-ista+1
3350 ii=ista+i-1
3351 datapd(i,j,cfld) = grid1(ii,jj)
3352 enddo
3353 enddo
3354 endif
3355 ENDIF
3356 ENDIF
3357!--- convective downdraft
3358 IF (iget(392) > 0) THEN
3359 IF (lvls(lp,iget(392)) > 0) THEN
3360!$omp parallel do private(i,j)
3361 DO j=jsta,jend
3362 DO i=ista,iend
3363 grid1(i,j) = d3dsl(i,j,24)
3364 ENDDO
3365 ENDDO
3366 id(1:25)=0
3367 itd3d = nint(td3d)
3368 if (itd3d /= 0) then
3369 ifincr = mod(ifhr,itd3d)
3370 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
3371 else
3372 ifincr = 0
3373 endif
3374 id(02)=133 ! Table 133
3375 id(18) = 0
3376 id(19) = ifhr
3377 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
3378 id(20) = 3
3379 IF (ifincr == 0) THEN
3380 id(18) = ifhr-itd3d
3381 ELSE
3382 id(18) = ifhr-ifincr
3383 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
3384 ENDIF
3385 if(grib == 'grib2') then
3386 cfld = cfld + 1
3387 fld_info(cfld)%ifld=iavblfld(iget(392))
3388 fld_info(cfld)%lvl=lvlsxml(lp,iget(392))
3389 if(itd3d==0) then
3390 fld_info(cfld)%ntrange=0
3391 else
3392 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
3393 endif
3394 fld_info(cfld)%tinvstat=itd3d
3395!$omp parallel do private(i,j,ii,jj)
3396 do j=1,jend-jsta+1
3397 jj = jsta+j-1
3398 do i=1,iend-ista+1
3399 ii=ista+i-1
3400 datapd(i,j,cfld) = grid1(ii,jj)
3401 enddo
3402 enddo
3403 endif
3404 ENDIF
3405 ENDIF
3406!--- convective detraintment
3407 IF (iget(393) > 0) THEN
3408 IF (lvls(lp,iget(393)) > 0) THEN
3409!$omp parallel do private(i,j)
3410 DO j=jsta,jend
3411 DO i=ista,iend
3412 grid1(i,j) = d3dsl(i,j,25)
3413 ENDDO
3414 ENDDO
3415 id(1:25)=0
3416 itd3d = nint(td3d)
3417 if (itd3d /= 0) then
3418 ifincr = mod(ifhr,itd3d)
3419 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
3420 else
3421 ifincr = 0
3422 endif
3423 id(02)=133 ! Table 133
3424 id(18) = 0
3425 id(19) = ifhr
3426 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
3427 id(20) = 3
3428 IF (ifincr == 0) THEN
3429 id(18) = ifhr-itd3d
3430 ELSE
3431 id(18) = ifhr-ifincr
3432 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
3433 ENDIF
3434 if(grib == 'grib2') then
3435 cfld = cfld + 1
3436 fld_info(cfld)%ifld=iavblfld(iget(393))
3437 fld_info(cfld)%lvl=lvlsxml(lp,iget(393))
3438 if(itd3d==0) then
3439 fld_info(cfld)%ntrange=0
3440 else
3441 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
3442 endif
3443 fld_info(cfld)%tinvstat=itd3d
3444!$omp parallel do private(i,j,ii,jj)
3445 do j=1,jend-jsta+1
3446 jj = jsta+j-1
3447 do i=1,iend-ista+1
3448 ii=ista+i-1
3449 datapd(i,j,cfld) = grid1(ii,jj)
3450 enddo
3451 enddo
3452 endif
3453 ENDIF
3454 ENDIF
3455!--- convective gravity drag zonal acce
3456 IF (iget(394) > 0) THEN
3457 IF (lvls(lp,iget(394)) > 0) THEN
3458!$omp parallel do private(i,j)
3459 DO j=jsta,jend
3460 DO i=ista,iend
3461 grid1(i,j) = d3dsl(i,j,26)
3462 ENDDO
3463 ENDDO
3464 id(1:25)=0
3465 itd3d = nint(td3d)
3466 if (itd3d /= 0) then
3467 ifincr = mod(ifhr,itd3d)
3468 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
3469 else
3470 ifincr = 0
3471 endif
3472 id(02)=133 ! Table 133
3473 id(18) = 0
3474 id(19) = ifhr
3475 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
3476 id(20) = 3
3477 IF (ifincr == 0) THEN
3478 id(18) = ifhr-itd3d
3479 ELSE
3480 id(18) = ifhr-ifincr
3481 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
3482 ENDIF
3483 if(grib == 'grib2') then
3484 cfld = cfld + 1
3485 fld_info(cfld)%ifld=iavblfld(iget(394))
3486 fld_info(cfld)%lvl=lvlsxml(lp,iget(394))
3487 if(itd3d==0) then
3488 fld_info(cfld)%ntrange=0
3489 else
3490 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
3491 endif
3492 fld_info(cfld)%tinvstat=itd3d
3493!$omp parallel do private(i,j,ii,jj)
3494 do j=1,jend-jsta+1
3495 jj = jsta+j-1
3496 do i=1,iend-ista+1
3497 ii=ista+i-1
3498 datapd(i,j,cfld) = grid1(ii,jj)
3499 enddo
3500 enddo
3501 endif
3502 ENDIF
3503 ENDIF
3504!--- convective gravity drag meridional acce
3505 IF (iget(395) > 0) THEN
3506 IF (lvls(lp,iget(395)) > 0) THEN
3507!$omp parallel do private(i,j)
3508 DO j=jsta,jend
3509 DO i=ista,iend
3510 grid1(i,j) = d3dsl(i,j,27)
3511 ENDDO
3512 ENDDO
3513 id(1:25)=0
3514 itd3d = nint(td3d)
3515 if (itd3d /= 0) then
3516 ifincr = mod(ifhr,itd3d)
3517 IF(ifmin >= 1)ifincr= mod(ifhr*60+ifmin,itd3d*60)
3518 else
3519 ifincr = 0
3520 endif
3521 id(02)=133 ! Table 133
3522 id(18) = 0
3523 id(19) = ifhr
3524 IF(ifmin >= 1)id(19)=ifhr*60+ifmin
3525 id(20) = 3
3526 IF (ifincr == 0) THEN
3527 id(18) = ifhr-itd3d
3528 ELSE
3529 id(18) = ifhr-ifincr
3530 IF(ifmin >= 1)id(18)=ifhr*60+ifmin-ifincr
3531 ENDIF
3532 if(grib == 'grib2') then
3533 cfld = cfld + 1
3534 fld_info(cfld)%ifld=iavblfld(iget(395))
3535 fld_info(cfld)%lvl=lvlsxml(lp,iget(395))
3536 if(itd3d==0) then
3537 fld_info(cfld)%ntrange=0
3538 else
3539 fld_info(cfld)%ntrange=(ifhr-id(18))/itd3d
3540 endif
3541 fld_info(cfld)%tinvstat=itd3d
3542!$omp parallel do private(i,j,ii,jj)
3543 do j=1,jend-jsta+1
3544 jj = jsta+j-1
3545 do i=1,iend-ista+1
3546 ii=ista+i-1
3547 datapd(i,j,cfld) = grid1(ii,jj)
3548 enddo
3549 enddo
3550 endif
3551 ENDIF
3552 ENDIF
3553 END IF ! end of d3d output
3554
3555! CHUANG: COMPUTE HAINES INDEX
3556 IF (iget(455) > 0) THEN
3557 ii=(ista+iend)/2+100
3558 jj=(jsta+jend)/2-100
3559 IF(abs(spl(lp)-50000.)<small) luhi=lp
3560 IF(abs(spl(lp)-70000.)<small) THEN ! high evevation
3561! HAINES=SPVAL
3562! print*,'computing dew point for Haine Index at ',SPL(LP)
3563!$omp parallel do private(i,j)
3564 DO j=jsta,jend
3565 DO i=ista,iend
3566 haines(i,j) = spval
3567 egrid2(i,j) = spl(lp)
3568 ENDDO
3569 ENDDO
3570 CALL caldwp(egrid2(ista:iend,jsta:jend),qsl(ista:iend,jsta:jend),tdsl(ista:iend,jsta:jend),tsl(ista:iend,jsta:jend))
3571
3572!$omp parallel do private(i,j,dum1,istaa,imois)
3573 DO j=jsta,jend
3574 DO i=ista,iend
3575 IF(sm(i,j) < 1.0 .AND. zint(i,j,lm+1) < fsl(i,j)*gi) THEN
3576 dum1 = tsl(i,j)-tprs(i,j,luhi)
3577 IF(dum1 <= 17.)THEN
3578 istaa = 1
3579 ELSE IF(dum1 > 17. .AND. dum1 <= 21.) THEN
3580 istaa = 2
3581 ELSE
3582 istaa = 3
3583 END IF
3584 dum1 = tsl(i,j)-tdsl(i,j)
3585 IF(dum1 <= 14.) THEN
3586 imois = 1
3587 ELSE IF(dum1>14. .AND. dum1<=20.) THEN
3588 imois = 2
3589 ELSE
3590 imois = 3
3591 END IF
3592 IF(tsl(i,j)<spval.and.tprs(i,j,luhi)<spval.and.tdsl(i,j)<spval)THEN
3593 haines(i,j) = istaa + imois
3594 ELSE
3595 haines(i,j) = spval
3596 ENDIF
3597! if(i==570 .and. j==574)print*,'high hainesindex:',i,j,luhi,tsl(i,j) &
3598! ,tprs(i,j,luhi),tdsl(i,j),ista,imois,spl(luhi),spl(lp),haines(i,j)
3599 END IF
3600 END DO
3601 END DO
3602
3603 luhi = lp
3604 ENDIF
3605
3606 IF(abs(spl(lp)-85000.)<small)THEN ! mid evevation
3607! print*,'computing dew point for Haine Index at ',SPL(LP)
3608!$omp parallel do private(i,j)
3609 DO j=jsta,jend
3610 DO i=ista,iend
3611 egrid2(i,j) = spl(lp)
3612 ENDDO
3613 ENDDO
3614 CALL caldwp(egrid2(ista:iend,jsta:jend),qsl(ista:iend,jsta:jend),tdsl(ista:iend,jsta:jend),tsl(ista:iend,jsta:jend))
3615
3616!$omp parallel do private(i,j,dum1,istaa,imois)
3617 DO j=jsta,jend
3618 DO i=ista,iend
3619 IF(sm(i,j) < 1.0 .AND. zint(i,j,lm+1) < fsl(i,j)*gi) THEN
3620 dum1 = tsl(i,j)-tprs(i,j,luhi)
3621 IF(dum1 <=5. ) THEN
3622 istaa = 1
3623 ELSE IF(dum1 > 5. .AND. dum1 <= 10.) THEN
3624 istaa = 2
3625 ELSE
3626 istaa = 3
3627 END IF
3628 dum1 = tsl(i,j)-tdsl(i,j)
3629 IF(dum1 <= 5.) THEN
3630 imois = 1
3631 ELSE IF(dum1 > 5. .AND. dum1 <= 12.) THEN
3632 imois = 2
3633 ELSE
3634 imois = 3
3635 END IF
3636! if(i==570 .and. j==574)print*,'mid haines index:',i,j,luhi,tsl(i,j) &
3637! ,tprs(i,j,luhi),tdsl(i,j),ista,imois,spl(luhi),spl(lp),haines(i,j)
3638 IF(tsl(i,j)<spval.and.tprs(i,j,luhi)<spval.and.tdsl(i,j)<spval)THEN
3639 haines(i,j) = istaa + imois
3640 ELSE
3641 haines(i,j) = spval
3642 ENDIF
3643 END IF
3644 END DO
3645 END DO
3646
3647 luhi = lp
3648 ENDIF
3649
3650 IF(abs(spl(lp)-95000.)<small)THEN ! LOW evevation
3651! print*,'computing dew point for Haine Index at ',SPL(LP)
3652!$omp parallel do private(i,j)
3653 DO j=jsta,jend
3654 DO i=ista,iend
3655 egrid2(i,j)=spl(lp)
3656 ENDDO
3657 ENDDO
3658 CALL caldwp(egrid2(ista:iend,jsta:jend),qsl(ista:iend,jsta:jend),tdsl(ista:iend,jsta:jend),tsl(ista:iend,jsta:jend))
3659
3660!$omp parallel do private(i,j,dum1,istaa,imois)
3661 DO j=jsta,jend
3662 DO i=ista,iend
3663 IF(sm(i,j) < 1.0 .AND. zint(i,j,lm+1) < fsl(i,j)*gi) THEN
3664 dum1 = tsl(i,j)-tprs(i,j,luhi)
3665 IF(dum1 <= 3.)THEN
3666 istaa = 1
3667 ELSE IF(dum1 > 3. .AND. dum1 <=7. ) THEN
3668 istaa = 2
3669 ELSE
3670 istaa = 3
3671 END IF
3672 dum1 = tsl(i,j)-tdsl(i,j)
3673 IF(dum1 <=5. ) THEN
3674 imois = 1
3675 ELSE IF(dum1 > 5. .AND. dum1 <= 9.)THEN
3676 imois = 2
3677 ELSE
3678 imois = 3
3679 END IF
3680! if(i==570 .and. j==574)print*,'low haines index:',i,j,luhi,tsl(i,j) &
3681! ,tprs(i,j,luhi),tdsl(i,j),ista,imois,spl(luhi),spl(lp),haines(i,j)
3682 IF(tsl(i,j)<spval.and.tprs(i,j,luhi)<spval.and.tdsl(i,j)<spval)THEN
3683 haines(i,j) = istaa + imois
3684 ELSE
3685 haines(i,j) = spval
3686 ENDIF
3687 END IF
3688 END DO
3689 END DO
3690
3691 if(grib == 'grib2') then
3692 cfld = cfld + 1
3693 fld_info(cfld)%ifld=iavblfld(iget(455))
3694!$omp parallel do private(i,j,ii,jj)
3695 do j=1,jend-jsta+1
3696 jj = jsta+j-1
3697 do i=1,iend-ista+1
3698 ii=ista+i-1
3699 datapd(i,j,cfld) = haines(ii,jj)
3700 enddo
3701 enddo
3702 endif
3703
3704 ENDIF
3705
3706 ENDIF
3707!
3708 ENDDO !*** END OF MAIN VERTICAL LOOP DO LP=1,LSM
3709!*** ENDIF FOR IF TEST SEEING IF WE WANT ANY OTHER VARIABLES
3710 ENDIF
3711! SRD
3712!
3713! MAX VERTICAL VELOCITY UPDRAFT
3714!
3715 IF (iget(423) > 0) THEN
3716! LP=22 ! 400 MB
3717 lp=46 ! 1000 MB
3718!$omp parallel do private(i,j)
3719 DO j=jsta,jend
3720 DO i=ista,iend
3721 grid1(i,j) = w_up_max(i,j)
3722! print *,' writing w_up_max, i,j, = ', w_up_max(i,j)
3723 ENDDO
3724 ENDDO
3725 if(grib == 'grib2')then
3726 cfld = cfld + 1
3727 fld_info(cfld)%ifld = iavblfld(iget(423))
3728 fld_info(cfld)%lvl = lvlsxml(lp,iget(423))
3729 if (ifhr > 0) then
3730 fld_info(cfld)%tinvstat=1
3731 else
3732 fld_info(cfld)%tinvstat=0
3733 endif
3734 fld_info(cfld)%ntrange=1
3735!$omp parallel do private(i,j,ii,jj)
3736 do j=1,jend-jsta+1
3737 jj = jsta+j-1
3738 do i=1,iend-ista+1
3739 ii=ista+i-1
3740 datapd(i,j,cfld) = grid1(ii,jj)
3741 enddo
3742 enddo
3743 endif
3744 ENDIF
3745!
3746! MAX VERTICAL VELOCITY DOWNDRAFT
3747!
3748 IF (iget(424) > 0) THEN
3749 lp = 46 ! 1000 MB
3750!$omp parallel do private(i,j)
3751 DO j=jsta,jend
3752 DO i=ista,iend
3753 grid1(i,j) = w_dn_max(i,j)
3754 ENDDO
3755 ENDDO
3756 if(grib == 'grib2')then
3757 cfld = cfld + 1
3758 fld_info(cfld)%ifld=iavblfld(iget(424))
3759 fld_info(cfld)%lvl=lvlsxml(lp,iget(424))
3760 if (ifhr > 0) then
3761 fld_info(cfld)%tinvstat=1
3762 else
3763 fld_info(cfld)%tinvstat=0
3764 endif
3765 fld_info(cfld)%ntrange=1
3766!$omp parallel do private(i,j,ii,jj)
3767 do j=1,jend-jsta+1
3768 jj = jsta+j-1
3769 do i=1,iend-ista+1
3770 ii=ista+i-1
3771 datapd(i,j,cfld) = grid1(ii,jj)
3772 enddo
3773 enddo
3774 endif
3775 ENDIF
3776!
3777! MEAN VERTICAL VELOCITY
3778!
3779! This hourly mean field is, in fact, based on the column
3780! mean bounded by sigma levels, but I chose to keep the
3781! output here with the other diagnostic vertical
3782! velocity fields
3783!
3784 IF (iget(425) > 0) THEN
3785 lp = 46 ! 1000 MB
3786!$omp parallel do private(i,j)
3787 DO j=jsta,jend
3788 DO i=ista,iend
3789 grid1(i,j) = w_mean(i,j)
3790 ENDDO
3791 ENDDO
3792 if(grib == 'grib2')then
3793 cfld = cfld + 1
3794 fld_info(cfld)%ifld = iavblfld(iget(425))
3795 fld_info(cfld)%lvl = lvlsxml(lp,iget(425))
3796 if (ifhr == 0) then
3797 fld_info(cfld)%tinvstat = 0
3798 else
3799 fld_info(cfld)%tinvstat = 1
3800 endif
3801 fld_info(cfld)%ntrange = 1
3802!$omp parallel do private(i,j,ii,jj)
3803 do j=1,jend-jsta+1
3804 jj = jsta+j-1
3805 do i=1,iend-ista+1
3806 ii=ista+i-1
3807 datapd(i,j,cfld) = grid1(ii,jj)
3808 enddo
3809 enddo
3810 endif
3811 ENDIF
3812! SRD
3813
3814!
3815! CALL MEMBRANE SLP REDUCTION IF REQESTED IN CONTROL FILE
3816!
3817! OUTPUT MEMBRANCE SLP
3818 IF(iget(023) > 0)THEN
3819 IF(gridtype == 'A'.OR. gridtype == 'B') then
3820 if(me==0)print*,'CALLING MEMSLP for A or B grid'
3821 CALL memslp(tprs,qprs,fprs)
3822 if(me==0)print*,'aft CALLING MEMSLP for A or B grid,pslp=', &
3823 maxval(pslp(ista:iend,jsta:jend)),minval(pslp(ista:iend,jsta:jend)),pslp((ista+iend)/2,(jsta+jend)/2)
3824 ELSE IF (gridtype == 'E')THEN
3825 if(me==0)print*,'CALLING MEMSLP_NMM for E grid'
3826! CALL MEMSLP_NMM(TPRS,QPRS,FPRS)
3827 ELSE
3828 print*,'unknow grid type-> WONT DERIVE MESINGER SLP'
3829 END IF
3830!$omp parallel do private(i,j)
3831 DO j=jsta,jend
3832 DO i=ista,iend
3833 grid1(i,j) = pslp(i,j)
3834 ENDDO
3835 ENDDO
3836! print *,'inmdl2p,pslp=',maxval(pslp(1:im,jsta:jend)),minval(pslp(1:im,jsta:jend))
3837! print *,'inmdl2p,point pslp=',pslp(im/2,(jsta+jend)/2),pslp(1,jsta),'cfld=',cfld
3838 if(grib == 'grib2')then
3839 cfld = cfld + 1
3840 fld_info(cfld)%ifld = iavblfld(iget(023))
3841!$omp parallel do private(i,j,ii,jj)
3842 do j=1,jend-jsta+1
3843 jj = jsta+j-1
3844 do i=1,iend-ista+1
3845 ii=ista+i-1
3846 datapd(i,j,cfld) = grid1(ii,jj)
3847 enddo
3848 enddo
3849 endif
3850 ENDIF
3851
3852! OUTPUT of MAPS SLP
3853 IF(iget(445) > 0)THEN
3854 if(me==0)print*,'CALLING MAPS SLP'
3855 CALL mapsslp(tprs)
3856!$omp parallel do private(i,j)
3857 DO j=jsta,jend
3858 DO i=ista,iend
3859 grid1(i,j) = pslp(i,j)
3860 ENDDO
3861 ENDDO
3862 if(grib == 'grib2') then
3863 cfld = cfld + 1
3864 fld_info(cfld)%ifld = iavblfld(iget(445))
3865!$omp parallel do private(i,j,ii,jj)
3866 do j=1,jend-jsta+1
3867 jj = jsta+j-1
3868 do i=1,iend-ista+1
3869 ii=ista+i-1
3870 datapd(i,j,cfld) = grid1(ii,jj)
3871 enddo
3872 enddo
3873 endif
3874 ENDIF
3875
3876
3877! ADJUST 1000 MB HEIGHT TO MEMBEANCE SLP
3878 IF(iget(023) > 0.OR.iget(445) > 0)THEN
3879 IF(iget(012) > 0)THEN
3880! dong add missing value to 1000 mb hgt
3881 grid1= spval
3882 DO lp=lsm,1,-1
3883 IF(abs(spl(lp)-1.0e5) <= 1.0e-5)THEN
3884 IF(lvls(lp,iget(012)) > 0)THEN
3885 alpth = log(1.e5)
3886 IF(modelname == 'GFS')THEN
3887! GFS does not want to adjust 1000 mb H to membrane SLP
3888! because MOS can't adjust to the much lower H
3889!$omp parallel do private(i,j)
3890 DO j=jsta,jend
3891 DO i=ista,iend
3892 IF(fsl(i,j)<spval)THEN
3893 grid1(i,j) = fsl(i,j)*gi
3894 ELSE
3895 grid1(i,j) = spval
3896 ENDIF
3897 ENDDO
3898 ENDDO
3899 ELSE
3900!$omp parallel do private(i,j,PSLPIJ,ALPSL,PSFC)
3901 DO j=jsta,jend
3902 DO i=ista,iend
3903 IF(pslp(i,j) < spval) THEN
3904 pslpij = pslp(i,j)
3905 alpsl = log(pslpij)
3906 psfc = pint(i,j,nint(lmh(i,j))+1)
3907 IF(abs(pslpij-psfc) < 5.e2) THEN
3908 grid1(i,j) = rd*tprs(i,j,lp)*(alpsl-alpth)
3909 ELSE
3910 grid1(i,j) = fis(i,j)/(alpsl-log(psfc))*(alpsl-alpth)
3911 ENDIF
3912 z1000(i,j) = grid1(i,j)*gi
3913 grid1(i,j) = z1000(i,j)
3914 ELSE
3915 grid1(i,j) = spval
3916 END IF
3917 ENDDO
3918 ENDDO
3919 END IF
3920
3921 IF (smflag) THEN
3922!tgs - smoothing of geopotential heights
3923 nsmooth = nint(5.*(13500./dxm))
3924 call allgetherv(grid1)
3925 do k=1,nsmooth
3926 CALL smooth(grid1,sdummy,im,jm,0.5)
3927 end do
3928 ENDIF
3929
3930 if(grib == 'grib2') then
3931 cfld = cfld + 1
3932 fld_info(cfld)%ifld = iavblfld(iget(012))
3933 fld_info(cfld)%lvl = lvlsxml(lp,iget(012))
3934!$omp parallel do private(i,j,ii,jj)
3935 do j=1,jend-jsta+1
3936 jj = jsta+j-1
3937 do i=1,iend-ista+1
3938 ii=ista+i-1
3939 datapd(i,j,cfld) = grid1(ii,jj)
3940 enddo
3941 enddo
3942 endif
3943 exit
3944 ENDIF
3945 ENDIF
3946 END DO
3947 ENDIF
3948 ENDIF
3949!
3950if(allocated(d3dsl)) deallocate(d3dsl)
3951if(allocated(dustsl)) deallocate(dustsl)
3952if(allocated(smokesl)) deallocate(smokesl)
3953! END OF ROUTINE.
3954!
3955 RETURN
3956 END
calicing
subroutine calicing(t1, rh, omga, icing)
Computes In-Flight Icing.
Definition AVIATION.f:161
calcat
subroutine calcat(u, v, h, u_old, v_old, h_old, cat)
Computes Clear Air Turbulence Index.
Definition AVIATION.f:265
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
caldwp
subroutine caldwp(p1d, q1d, tdwp, t1d)
Computes dewpoint from P, T, and Q.
Definition CALDWP.f:28
smooth
subroutine smooth(field, hold, ix, iy, smth)
Definition SMOOTH.f:36
smoothc
subroutine smoothc(field, hold, ix, iy, smth)
smoothc() smooths a meteorological field using Shapiro smoother.
Definition SMOOTH.f:128
upp_physics
calcape() computes CAPE/CINS and other storm related variables.
Definition UPP_PHYSICS.f:27
upp_physics::fpvsnew
elemental real function, public fpvsnew(t)
Definition UPP_PHYSICS.f:346