UPP (upp-srw-2.2.0)
Loading...
Searching...
No Matches
INITPOST_NEMS.f
Go to the documentation of this file.
1
5
20 SUBROUTINE initpost_nems(NREC,nfile)
21
22 use vrbls3d, only: t, q, uh, vh, q2, cwm, f_ice, f_rain, f_rimef, cfr, pint,&
23 pint, alpint, pmid, pmidv, zint, zmid, wh, rlwtt, rswtt,&
24 ttnd, tcucn, train, el_pbl, exch_h, omga, qqni, qqnr, qqw, qqi, &
25 qqr, qqs, qqg, ref_10cm, radius_cloud, radius_ice, radius_snow, &
26 extcof55, aextc55
27 use vrbls2d, only: f, pd, fis, pblh, mixht, ustar, z0, ths, qs, twbs, qwbs, prec,&
28 acprec, cuprec,ancprc, lspa, sno, snoavg, psfcavg, t10avg, t10m, akhsavg, akmsavg,&
29 refd_max, w_up_max, w_dn_max, up_heli_max, si, cldefi, th10, q10, pshltr,&
30 tshltr, qshltr, maxtshltr, mintshltr, maxrhshltr, minrhshltr, akhs, akms, albase,&
31 albedo, czen, cfracl, cfracm, islope, cmc, grnflx, pctsno, soiltb, vegfrc,&
32 acfrcv, acfrst, ssroff, bgroff, czmean, mxsnal, radot, sigt4, tg, sr, cfrach,&
33 rlwin, rlwtoa, alwin, alwout, alwtoa, rswin, rswinc, rswout, aswin,aswout,&
34 aswtoa, sfcshx, sfclhx, subshx, snopcx, sfcuvx, potevp, ncfrcv, ncfrst, u10h,&
35 u10, v10h, v10, u10max, v10max, smstav, smstot, sfcevp, ivgtyp, acsnow, acsnom,&
36 sst, thz0, qz0, uz0, vz0, htop, isltyp, sfcexc, hbot, htopd, htops, cuppt, cprate,&
37 hbotd, hbots, prate_max, fprate_max
38 use soil, only: sldpth, sh2o, smc, stc
39 use masks, only: lmv, lmh, htm, vtm, dx, dy, hbm2, gdlat, gdlon, sm, sice
40 use kinds, only: i_llong
41 use wrf_io_flags_mod, only:
42 use params_mod, only: pi, dtr, g, d608, rd,tfrz
43 use lookup_mod, only: thl, plq, ptbl, ttbl, rdq, rdth, rdp, rdthe, pl, qs0, sqs, sthe, the0,&
44 ttblq, rdpq, rdtheq, stheq, the0q
45 use ctlblk_mod, only: me, mpi_comm_comp, global, icnt, idsp, jsta, ihrst, imin, idat, sdat,&
46 ifhr, ifmin, filename, restrt, imp_physics, isf_surface_physics, icu_physics, jend,&
47 dt, spval, gdsdegr, grib, pdtop, pt, tmaxmin, nsoil, lp1, jend_m, nprec, nphs, avrain,&
48 avcnvc, ardlw, ardsw, asrfc, novegtype, spl, lsm, dtq2, tsrfc, trdlw, trdsw, theat, tclod,&
49 tprec, alsl, lm , im, jm, jsta_2l, jend_2u, ivegsrc, pthresh
50 use gridspec_mod, only: dyval, dxval, cenlat, cenlon, maptype, gridtype, latstart, latlast, latnw,&
51 latse, lonstart, lonlast, lonnw, lonse, latstartv, latlastv, cenlatv, lonstartv,&
52 lonlastv, cenlonv
53 use nemsio_module, only: nemsio_gfile, nemsio_getfilehead, nemsio_close, nemsio_getheadvar
54 use upp_math, only: h2u
55!
56! INCLUDE/SET PARAMETERS.
57 implicit none
58!
59 type(nemsio_gfile),intent(inout) :: nfile
60!
61 include "mpif.h"
62! This version of INITPOST shows how to initialize, open, read from, and
63! close a NetCDF dataset. In order to change it to read an internal (binary)
64! dataset, do a global replacement of _ncd_ with _int_.
65
66 character(len=20) :: VarName
67 character(len=20) :: VcoordName
68 integer :: Status
69 character startdate*19,SysDepInfo*80,cgar*1
70 character startdate2(19)*4
71!
72! NOTE: SOME INTEGER VARIABLES ARE READ INTO DUMMY ( A REAL ). THIS IS OK
73! AS LONG AS REALS AND INTEGERS ARE THE SAME SIZE.
74!
75! ALSO, EXTRACT IS CALLED WITH DUMMY ( A REAL ) EVEN WHEN THE NUMBERS ARE
76! INTEGERS - THIS IS OK AS LONG AS INTEGERS AND REALS ARE THE SAME SIZE.
77 logical, parameter :: debugprint = .false.
78! logical, parameter :: debugprint = .true.
79 logical :: convert_rad_to_deg=.false.
80! logical global
81! CHARACTER BLANK*4
82 integer nfhour ! forecast hour from nems io file
83 INTEGER IDATE(8),JDATE(8)
84!
85! DECLARE VARIABLES.
86!
87 REAL FACT,tsph,tstart
88 REAL RINC(5)
89 REAL ETA1(LM+1), ETA2(LM+1)
90 REAL GARB
91 REAL DUM1D (LM+1)
92 REAL DUMMY ( IM, JM )
93 REAL DUMMY2 ( IM, JM )
94 real, allocatable :: fi(:,:,:)
95
96 integer ibuf(im,jsta_2l:jend_2u)
97 real buf(im,jsta_2l:jend_2u)
98 character*8,allocatable:: recname(:)
99 character*16,allocatable :: reclevtyp(:)
100 integer,allocatable:: reclev(:)
101 real, allocatable:: bufy(:)
102 real, allocatable:: glat1d(:),glon1d(:)
103!jw
104 integer ii,jj,js,je,jev,iyear,imn,iday,itmp,ioutcount,istatus, &
105 nsrfc,nrdlw,nrdsw,nheat,nclod, &
106 iunit,nrec,i,j,l, iret,nframe,impf,jmpf,nframed2, &
107 igdout,ll,n,im1,jm1,iim1,item
108!
109! DATA BLANK/' '/
110!
111!***********************************************************************
112! START INIT HERE.
113!
114 if(me==0)WRITE(6,*)'INITPOST: ENTER INITPOST'
115!
116!
117! STEP 1. READ MODEL OUTPUT FILE
118!
119!***
120! LMH always = LM for sigma-type vert coord
121! LMV always = LM for sigma-type vert coord
122
123 do j = jsta_2l, jend_2u
124 do i = 1, im
125 lmv( i, j ) = lm
126 lmh( i, j ) = lm
127 end do
128 end do
129
130! HTM VTM all 1 for sigma-type vert coord
131
132 do l = 1, lm
133 do j = jsta_2l, jend_2u
134 do i = 1, im
135 htm( i, j, l ) = 1.0
136 vtm( i, j, l ) = 1.0
137 end do
138 end do
139 end do
140
141! The end j row is going to be jend_2u for all variables except for V.
142 js=jsta_2l
143 je=jend_2u
144 IF (jend_2u==jm) THEN
145 jev=jend_2u+1
146 ELSE
147 jev=jend_2u
148 ENDIF
149! sample print point
150 ii=(1+im)/2
151 jj=(1+jm)/2
152! get start date
153 idate=0
154 if (me == 0)then
155 print*,'nrec=',nrec
156 allocate(recname(nrec),reclevtyp(nrec),reclev(nrec))
157
158 call nemsio_getfilehead(nfile,iret=iret &
159 ,idate=idate(1:7),nfhour=nfhour,recname=recname &
160 ,reclevtyp=reclevtyp,reclev=reclev,nframe=nframe)
161! if(iret/=0)print*,'error getting idate,fhour, stopping';stop
162 if (debugprint) then
163 print *,'printing an inventory of NEMS file'
164 do i=1,nrec
165 print *,'recname,reclevtyp,reclev=',trim(recname(i)),' ', &
166 trim(reclevtyp(i)),reclev(i)
167 end do
168 endif
169! print *,'reclevtyp=',(trim(reclevtyp(i)),i=1,nrec)
170! print *,'reclev=',(reclev(i),i=1,nrec)
171 deallocate(recname,reclevtyp,reclev)
172 impf=im+nframe*2
173 jmpf=jm+nframe*2
174! nframed2=nframe/2
175 print*,'nframe,impf,jmpf= ',nframe,impf,jmpf
176 allocate(glat1d(impf*jmpf),glon1d(impf*jmpf) )
177 call nemsio_getfilehead(nfile,dx=glat1d &
178 ,dy=glon1d,iret=iret)
179 if(iret/=0)print*,'did not find dx dy'
180! do j=1,impf*jmpf
181! print*,'dx before scatter= ',j,glat1d(j)
182! end do
183!$omp parallel do private(i,j,item)
184 do j=1,jm
185 item = (j-1)*impf + nframe
186 do i=1,im
187 dummy(i,j) = glat1d(item+i)
188 dummy2(i,j) = glon1d(item+i)
189! dummy(i,j)=glat1d(i-nframe,j-nframe)
190! dummy2(i,j)=glon1d(i-nframe,j-nframe)
191 end do
192 end do
193 deallocate(glat1d,glon1d)
194! latstart=nint(dummy(1,1)*1000.)
195! latlast=nint(dummy(im,jm)*1000.)
196! lonstart=nint(dummy2(1,1)*1000.)
197! lonlast=nint(dummy2(im,jm)*1000.)
198! dyval=nint((dummy(1,2)-dummy(1,1))*1000.)
199! dxval=nint((dummy(2,1)-dummy(1,1))*1000.)
200! cenlat=nint(dummy(ii,jj)*1000.)
201! cenlon=nint(dummy2(ii,jj)*1000.)
202 print*,'idate before broadcast = ',(idate(i),i=1,7)
203 end if
204 call mpi_bcast(idate(1),7,mpi_integer,0,mpi_comm_comp,iret)
205 call mpi_bcast(nfhour,1,mpi_integer,0,mpi_comm_comp,iret)
206 call mpi_bcast(nframe,1,mpi_integer,0,mpi_comm_comp,iret)
207! call mpi_bcast(latstart,1,MPI_INTEGER,0,mpi_comm_comp,iret)
208! call mpi_bcast(latlast,1,MPI_INTEGER,0,mpi_comm_comp,iret)
209! call mpi_bcast(lonstart,1,MPI_INTEGER,0,mpi_comm_comp,iret)
210! call mpi_bcast(lonlast,1,MPI_INTEGER,0,mpi_comm_comp,iret)
211! call mpi_bcast(cenlat,1,MPI_INTEGER,0,mpi_comm_comp,iret)
212! call mpi_bcast(cenlon,1,MPI_INTEGER,0,mpi_comm_comp,iret)
213! print*, 'latstart,latlast A calling bcast=',latstart,latlast
214! print*,'lonstart,lonlast A calling bcast=',lonstart,lonlast
215! print*,'cenlat,cenlon A calling bcast=',cenlat,cenlon
216
217! if(me == 0)then
218! call nemsio_getheadvar(nfile,'global',global,iret)
219! if (iret /= 0) then
220! print*,"global not found in file-Assigned false"
221! global=.FALSE.
222! end if
223! end if
224! call mpi_bcast(global,1,MPI_LOGICAL,0,mpi_comm_comp,iret)
225
226! print*,'Is this a global run ',global
227 IF(.not. global)THEN
228! nframe=0 ! Wang added option to read without halos, so specify nframe=0
229 impf=im+nframe*2
230 jmpf=jm+nframe*2
231! nframed2=nframe/2
232 ELSE
233! nframe=1 !
234 impf=im+1 ! post cut im off because it's the same as i=1 but data from model is till im
235 jmpf=jm
236! nframed2=nframe/2
237 END IF
238
239 if (debugprint) then
240 print*,'impf,jmpf,nframe for reading fields = ',impf,jmpf,nframe
241 print*,'idate after broadcast = ',(idate(i),i=1,7)
242 print*,'nfhour = ',nfhour
243 end if
244 call mpi_scatterv(dummy(1,1),icnt,idsp,mpi_real &
245 ,dx(1,jsta),icnt(me),mpi_real,0,mpi_comm_comp,iret)
246 call mpi_scatterv(dummy2(1,1),icnt,idsp,mpi_real &
247 ,dy(1,jsta),icnt(me),mpi_real,0,mpi_comm_comp,iret)
248
249! print *,'before call EXCH,mype=',me,'max(gdlat)=',maxval(gdlat),'max(gdlon)=', &
250! maxval(gdlon)
251! CALL EXCH(gdlat(1,JSTA_2L))
252! print *,'after call EXCH,mype=',me
253
254 iyear = idate(1)
255 imn = idate(2) ! ask Jun
256 iday = idate(3) ! ask Jun
257 ihrst = idate(4)
258 imin = idate(5)
259 jdate = 0
260 idate = 0
261!
262! read(startdate,15)iyear,imn,iday,ihrst,imin
263 15 format(i4,1x,i2,1x,i2,1x,i2,1x,i2)
264 if (me==0) then
265 print*,'start yr mo day hr min =',iyear,imn,iday,ihrst,imin
266 print*,'processing yr mo day hr min=' &
267 ,idat(3),idat(1),idat(2),idat(4),idat(5)
268 endif
269!
270 idate(1) = iyear
271 idate(2) = imn
272 idate(3) = iday
273 idate(5) = ihrst
274 idate(6) = imin
275 sdat(1) = imn
276 sdat(2) = iday
277 sdat(3) = iyear
278 jdate(1) = idat(3)
279 jdate(2) = idat(1)
280 jdate(3) = idat(2)
281 jdate(5) = idat(4)
282 jdate(6) = idat(5)
283!
284
285! CALL W3DIFDAT(JDATE,IDATE,2,RINC)
286! ifhr=nint(rinc(2))
287!
288 CALL w3difdat(jdate,idate,0,rinc)
289!
290 ifhr=nint(rinc(2)+rinc(1)*24.)
291 ifmin=nint(rinc(3))
292! if(ifhr /= nfhour)print*,'find wrong Model input file';stop
293 if (me==0)print*,' in INITPOST ifhr ifmin fileName=',ifhr,ifmin,trim(filename)
294
295! Getting tstart
296 tstart=0.
297
298! Getiing restart
299
300 restrt=.true. ! set RESTRT as default
301! call ext_int_get_dom_ti_integer(DataHandle,'RESTARTBIN',itmp
302! + ,1,ioutcount,istatus)
303
304! IF(itmp < 1)THEN
305! RESTRT=.FALSE.
306! ELSE
307! RESTRT=.TRUE.
308! END IF
309
310! print*,'status for getting RESTARTBIN= ',istatus
311
312! print*,'Is this a restrt run? ',RESTRT
313
314 IF(tstart > 1.0e-2)THEN
315 ifhr=ifhr+nint(tstart)
316 rinc=0
317 idate=0
318 rinc(2)=-1.0*ifhr
319 call w3movdat(rinc,jdate,idate)
320 sdat(1)=idate(2)
321 sdat(2)=idate(3)
322 sdat(3)=idate(1)
323 ihrst=idate(5)
324 print*,'new forecast hours for restrt run= ',ifhr
325 print*,'new start yr mo day hr min =',sdat(3),sdat(1) &
326 ,sdat(2),ihrst,imin
327 END IF
328
329
330!KRF: Initialize extinction coef for aerosol to zero to avoid failure.
331! These are not in NEMS model output, but new CALVIS_GSD methods uses
332! these fields from ARW, and if not initialized here will cause failure.
333 extcof55=0.
334 aextc55=0.
335
336!Chuang: set default to Ferrier's MP scheme. NPS does not write MP option
337!used in model to nemsio output
338 varname='mp_physi'
339 if(me == 0)then
340 call nemsio_getheadvar(nfile,trim(varname),imp_physics,iret)
341 if (iret /= 0) then
342 print*,varname," not found in file- go to 16 character "
343 varname='mp_physics'
344 call nemsio_getheadvar(nfile,trim(varname),imp_physics,iret)
345 if (iret /= 0) then
346 print*,varname," not found in file-Assigned 1000"
347! imp_physics=1000 ! should this line be uncommented? - Moorthi
348 imp_physics=5
349 end if
350 end if
351 end if
352 call mpi_bcast(imp_physics,1,mpi_integer,0,mpi_comm_comp,iret)
353 if(me==0)print*,'MP_PHYSICS= ',imp_physics
354
355! Initializes constants for Ferrier microphysics
356 if(imp_physics==5 .or. imp_physics==85 .or. imp_physics==95)then
357 CALL microinit(imp_physics)
358 end if
359
360 varname='sf_surface_physi'
361 if(me == 0)then
362 call nemsio_getheadvar(nfile,trim(varname),isf_surface_physics,iret)
363 if (iret /= 0) then
364 print*,varname," not found in file-Assigned 2 for NOAH LSM as default"
365 isf_surface_physics=2
366 end if
367 end if
368 call mpi_bcast(isf_surface_physics,1,mpi_integer,0,mpi_comm_comp,iret)
369 if(me==0) print*,'SF_SURFACE_PHYSICS= ',isf_surface_physics
370
371! IVEGSRC=1 for IGBP and 0 for USGS
372 varname='IVEGSRC'
373 if(me == 0)then
374 call nemsio_getheadvar(nfile,trim(varname),ivegsrc,iret)
375 if (iret /= 0) then
376 print*,varname," not found in file-Assigned 1 for IGBP as default"
377 ivegsrc=1
378 end if
379 end if
380 call mpi_bcast(ivegsrc,1,mpi_integer,0,mpi_comm_comp,iret)
381 if(me==0) print*,'IVEGSRC= ',ivegsrc
382
383! set novegtype based on vegetation classification
384 if(ivegsrc==1)then
385 novegtype=20
386 else if(ivegsrc==0)then
387 novegtype=24
388 end if
389 if(me==0) print*,'novegtype= ',novegtype
390
391 varname='CU_PHYSICS'
392 if(me == 0)then
393 call nemsio_getheadvar(nfile,trim(varname),icu_physics,iret)
394 if (iret /= 0) then
395 print*,varname," not found in file-Assigned 2 for BMJ as default"
396 icu_physics=2
397 end if
398 end if
399 call mpi_bcast(icu_physics,1,mpi_integer,0,mpi_comm_comp,iret)
400 if(me==0) print*,'CU_PHYSICS= ',icu_physics
401
402
403 allocate(bufy(jm))
404 varname='DX'
405! if(me == 0)then
406! call nemsio_getheadvar(nfile,trim(VarName),bufy,iret)
407! if (iret /= 0) then
408! print*,VarName," not found in file-Assigned missing values"
409! dx=spval
410! end if
411! end if
412! call mpi_bcast(bufy,jm,MPI_REAL,0,mpi_comm_comp,iret)
413! do j=jsta,jend
414! do i=1,im
415! dx(i,j)=bufy(j)
416! end do
417! end do
418 if(debugprint)print*,'sample ',varname,' = ',dx(im/2,(jsta+jend)/2)
419
420 varname='DY'
421! if(me == 0)then
422! call nemsio_getheadvar(nfile,trim(VarName),bufy,iret)
423! if (iret /= 0) then
424! print*,VarName," not found in file-Assigned missing values"
425! dx=spval
426! end if
427! end if
428! call mpi_bcast(bufy,jm,MPI_REAL,0,mpi_comm_comp,iret)
429! do j=jsta,jend
430! do i=1,im
431! dy(i,j)=bufy(j)
432! end do
433! end do
434 if(debugprint)print*,'sample ',varname,' = ',dy(im/2,(jsta+jend)/2)
435 deallocate(bufy)
436
437 varname='dt'
438 if(me == 0)then
439 call nemsio_getheadvar(nfile,trim(varname),garb,iret)
440 if (iret /= 0) then
441 print*,varname," not found in file-Assigned missing values"
442 dt=spval
443 else
444 dt=garb
445 end if
446 end if
447 call mpi_bcast(dt,1,mpi_real,0,mpi_comm_comp,iret)
448
449 varname='dphd'
450 if(me == 0)then
451 call nemsio_getheadvar(nfile,trim(varname),garb,iret)
452 if (iret /= 0) then
453 print*,varname," not found in file-Assigned missing values"
454 dyval=spval
455 else
456 dyval=garb*gdsdegr
457 end if
458 end if
459 call mpi_bcast(dyval,1,mpi_real,0,mpi_comm_comp,iret)
460! dyval=106 ! hard wire for AQ domain testing
461
462 varname='dlmd'
463 if(me == 0)then
464 call nemsio_getheadvar(nfile,trim(varname),garb,iret)
465 if (iret /= 0) then
466 print*,varname," not found in file-Assigned missing values"
467 dxval=spval
468 else
469 dxval=garb*gdsdegr
470 end if
471 end if
472 call mpi_bcast(dxval,1,mpi_real,0,mpi_comm_comp,iret)
473! dxval=124 ! hard wire for AQ domain testing
474
475 if(me==0) print*,'DX, DY, DT=',dxval,dyval,dt
476
477 varname='TPH0D'
478 if(me == 0)then
479 call nemsio_getheadvar(nfile,trim(varname),garb,iret)
480 if (iret /= 0) then
481 print*,varname," not found in file-Assigned missing values"
482 cenlat=spval
483 else
484 cenlat=nint(garb*gdsdegr)
485 end if
486 end if
487 call mpi_bcast(cenlat,1,mpi_integer,0,mpi_comm_comp,iret)
488
489 varname='TLM0D'
490 if(me == 0)then
491 call nemsio_getheadvar(nfile,trim(varname),garb,iret)
492 if (iret /= 0) then
493 print*,varname," not found in file-Assigned missing values"
494 cenlon=spval
495 else
496 if(grib=="grib2") then
497 cenlon=nint((garb+360.)*gdsdegr)
498 endif
499 end if
500 end if
501 call mpi_bcast(cenlon,1,mpi_integer,0,mpi_comm_comp,iret)
502
503! VarName='TRUELAT1'
504! call retrieve_index(index,VarName,varname_all,nrecs,iret)
505! if (iret /= 0) then
506! print*,VarName," not found in file"
507! else
508! call mpi_file_read_at(iunit,file_offset(index)+5*4 &
509! ,garb,1,mpi_real4, mpi_status_ignore, ierr)
510! if (ierr /= 0) then
511! print*,"Error reading ", VarName," using MPIIO"
512! else
513! print*,VarName, ' from MPIIO READ= ',garb
514! TRUELAT1=nint(garb*1000.)
515! write(6,*) 'truelat1= ', TRUELAT1
516! end if
517! end if
518
519! VarName='TRUELAT2'
520! call retrieve_index(index,VarName,varname_all,nrecs,iret)
521! if (iret /= 0) then
522! print*,VarName," not found in file"
523! else
524! call mpi_file_read_at(iunit,file_offset(index)+5*4 &
525! ,garb,1,mpi_real4, mpi_status_ignore, ierr)
526! if (ierr /= 0) then
527! print*,"Error reading ", VarName," using MPIIO"
528! else
529! print*,VarName, ' from MPIIO READ= ',garb
530! TRUELAT2=nint(garb*1000.)
531! write(6,*) 'truelat2= ', TRUELAT2
532! end if
533! end if
534
535! VarName='MAP_PROJ'
536! if(me == 0)then
537! call nemsio_getheadvar(nfile,trim(VarName),maptype,iret)
538! if (iret /= 0) then
539! print*,VarName," not found in file-Assigned 1000"
540! maptype=1000
541! end if
542! end if
543! call mpi_bcast(maptype,1,MPI_INTEGER,0,mpi_comm_comp,iret)
544 IF(.not. global)THEN
545 maptype=205 ! for Arakawa-B grid
546 gridtype='B'
547 ELSE
548 maptype=0 ! for global NMMB on latlon grid
549 gridtype='A' ! will put wind on mass point for now to make regular latlon
550 END IF
551 if(me==0) print*,'maptype and gridtype= ',maptype,gridtype
552
553 hbm2=1.0
554
555 varname='glat'
556 vcoordname='sfc'
557 l=1
558 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
559 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
560 ,l,impf,jmpf,nframe,gdlat)
561
562 call collect_loc(gdlat,dummy)
563! decides whether or not to convert to degree
564 if(me==0)then
565 if(maxval(abs(dummy))<pi)then ! convert from radian to degree
566 if(debugprint)print*,'convert from radian to degree'
567 dummy=dummy*180./pi
568 convert_rad_to_deg=.true.
569 end if
570 end if
571! print *,' gdlat1=',gdlat(1,:)
572! print *,' gdlatim=',gdlat(im,:)
573
574 call mpi_bcast(convert_rad_to_deg,1,mpi_logical,0,mpi_comm_comp,iret)
575 if(convert_rad_to_deg)call mpi_scatterv(dummy(1,1),icnt,idsp,mpi_real &
576 ,gdlat(1,jsta),icnt(me),mpi_real,0,mpi_comm_comp,iret)
577 if(debugprint)print*,'sample ',varname,' = ',gdlat(im/2,(jsta+jend)/2)
578 if(debugprint)print*,'max min lat=',maxval(gdlat),minval(gdlat),'im=',im, &
579 'jsta_2l=',jsta_2l,'jend_2u=',jend_2u
580 !call collect_loc(gdlat,dummy)
581 if(me==0.and.debugprint)print*,'after collect lat=',dummy(1,1),dummy(im,jm)
582 if(me==0)then
583 ii=(1+im)/2
584 jj=(1+jm)/2
585 latstart=nint(dummy(1,1)*gdsdegr)
586 latlast=nint(dummy(im,jm)*gdsdegr)
587 latnw=nint(dummy(1,jm)*gdsdegr)
588 latse=nint(dummy(im,1)*gdsdegr)
589! dyval=nint((dummy(1,2)-dummy(1,1))*1000.)
590! dyval=106 ! hard wire for AQ domain testing
591 if(mod(im,2)==0)then
592! cenlat=nint((dummy(ii,jj)+dummy(ii+1,jj)+dummy(ii+1,jj+1)+dummy(ii,jj+1))/4.0*1000.)
593 else
594! cenlat=nint(dummy(ii,jj)*1000.)
595 end if
596 print*,'latstart,latlast B bcast= ',latstart,latlast
597 end if
598 call mpi_bcast(latstart,1,mpi_integer,0,mpi_comm_comp,iret)
599 call mpi_bcast(latlast,1,mpi_integer,0,mpi_comm_comp,iret)
600! call mpi_bcast(dyval,1,MPI_INTEGER,0,mpi_comm_comp,iret)
601! call mpi_bcast(cenlat,1,MPI_INTEGER,0,mpi_comm_comp,iret)
602 if(debugprint) write(6,*) 'latstart,latlast,me A calling bcast=',latstart,latlast,me
603 if(me==0) print*,'dyval, cenlat= ',dyval, cenlat
604
605!$omp parallel do private(i,j)
606 do j=jsta,jend
607 do i=1,im
608 f(i,j) = 1.454441e-4*sin(gdlat(i,j)*dtr) ! 2*omeg*sin(phi)
609 end do
610 end do
611
612 varname='glon'
613 vcoordname='sfc'
614 l=1
615 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
616 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
617 ,l,impf,jmpf,nframe,gdlon)
618 if(convert_rad_to_deg)gdlon=gdlon*180./pi
619 if(global)then
620! do j=jsta,jend
621! do i=1,im
622! if(gdlon(i,j)<0.)gdlon(i,j)=360.+gdlon(i,j)
623! end do
624! end do
625 if(gdlon(1,jsta)>0. .and. gdlon(2,jsta)<0.)then
626 do j=jsta,jend
627 gdlon(1,j) = gdlon(1,j)-360.0
628 end do
629 end if
630 end if
631 if(debugprint)print*,'sample ',varname,' = ',(gdlon(i,(jsta+jend)/2),i=1,im,8)
632 if(debugprint)print*,'max min lon=',maxval(gdlon),minval(gdlon)
633 call collect_loc(gdlon,dummy)
634 if(me==0)then
635 if(grib=='grib2') then
636 if(dummy(1,1)<0) dummy(1,1)=dummy(1,1)+360.
637 if(dummy(im,jm)<0) dummy(im,jm)=dummy(im,jm)+360.
638 endif
639 lonstart=nint(dummy(1,1)*gdsdegr)
640 lonlast=nint(dummy(im,jm)*gdsdegr)
641 lonnw=nint(dummy(1,jm)*gdsdegr)
642 lonse=nint(dummy(im,1)*gdsdegr)
643! dxval=nint((dummy(2,1)-dummy(1,1))*1000.)
644! dxval=124 ! hard wire for AQ domain testing
645 if(mod(im,2)==0)then
646! cenlon=nint((dummy(ii,jj)+dummy(ii+1,jj)+dummy(ii+1,jj+1)+dummy(ii,jj+1))/4.0*1000.)
647 else
648! cenlon=nint(dummy(ii,jj)*1000.)
649 end if
650 end if
651 call mpi_bcast(lonstart,1,mpi_integer,0,mpi_comm_comp,iret)
652 call mpi_bcast(lonlast,1,mpi_integer,0,mpi_comm_comp,iret)
653! call mpi_bcast(dxval,1,MPI_INTEGER,0,mpi_comm_comp,iret)
654! call mpi_bcast(cenlon,1,MPI_INTEGER,0,mpi_comm_comp,iret)
655 write(6,*)'lonstart,lonlast A calling bcast=',lonstart,lonlast
656 print*,'dxval, cenlon= ',dxval, cenlon
657
658 convert_rad_to_deg=.false.
659 varname='vlat'
660 vcoordname='sfc'
661 l=1
662 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
663 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
664 ,l,impf,jmpf,nframe,buf)
665 if(debugprint)print*,'sample ',varname,' = ',buf(im/2,(jsta+jend)/2)
666 if(debugprint)print*,'max min vlat=',maxval(buf),minval(buf)
667 call collect_loc(buf,dummy)
668 if(me==0)then
669 if(maxval(abs(dummy))<pi)then ! convert from radian to degree
670 dummy(1,1)=dummy(1,1)*180./pi
671 dummy(im,jm)=dummy(im,jm)*180./pi
672 convert_rad_to_deg=.true.
673 end if
674 latstartv=nint(dummy(1,1)*gdsdegr)
675 latlastv=nint(dummy(im,jm)*gdsdegr)
676! cenlatv=nint(dummy(ii,jj)*1000.)
677! print*,'latstartv,cenlatv B bcast= ',latstartv,cenlatv
678 end if
679 call mpi_bcast(latstartv,1,mpi_integer,0,mpi_comm_comp,iret)
680 call mpi_bcast(latlastv,1,mpi_integer,0,mpi_comm_comp,iret)
681! call mpi_bcast(cenlatv,1,MPI_INTEGER,0,mpi_comm_comp,iret)
682 cenlatv=cenlat
683 if(debugprint)write(6,*) 'latstartv,cenlatv,latlastv,me A calling bcast=', &
684 latstartv,cenlatv,latlastv,me
685
686 varname='vlon'
687 vcoordname='sfc'
688 l=1
689 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
690 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
691 ,l,impf,jmpf,nframe,buf)
692 if(debugprint)print*,'sample ',varname,' = ',buf(im/2,(jsta+jend)/2)
693 if(debugprint)print*,'max min vlon=',maxval(buf),minval(buf)
694 call collect_loc(buf,dummy)
695 if(me==0)then
696 if(convert_rad_to_deg)then
697 dummy(1,1)=dummy(1,1)*180./pi
698 dummy(im,jm)=dummy(im,jm)*180./pi
699 end if
700 if(grib=='grib2') then
701 if(dummy(1,1)<0) dummy(1,1)=dummy(1,1)+360.
702 endif
703 lonstartv=nint(dummy(1,1)*gdsdegr)
704 lonlastv=nint(dummy(im,jm)*gdsdegr)
705! cenlonv=nint(dummy(ii,jj)*1000.)
706! print*,'lonstartv,cenlonv B bcast= ',lonstartv,cenlonv
707 end if
708 call mpi_bcast(lonstartv,1,mpi_integer,0,mpi_comm_comp,iret)
709 call mpi_bcast(lonlastv,1,mpi_integer,0,mpi_comm_comp,iret)
710! call mpi_bcast(cenlonv,1,MPI_INTEGER,0,mpi_comm_comp,iret)
711 cenlonv=cenlon
712 write(6,*) 'lonstartv,cenlonv,lonlastv,me A calling bcast=', &
713 lonstartv,cenlonv,lonlastv,me
714
715 varname='sm'
716 vcoordname='sfc'
717 l=1
718 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
719 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
720 ,l,impf,jmpf,nframe,sm)
721 if(debugprint)print*,'sample ',varname,' = ',sm(im/2,(jsta+jend)/2)
722
723 varname='sice'
724 vcoordname='sfc'
725 l=1
726 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
727 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
728 ,l,impf,jmpf,nframe,sice)
729 if(debugprint)print*,'sample ',varname,' = ',sice(im/2,(jsta+jend)/2)
730
731
732 varname='dpres'
733 vcoordname='hybrid sig lev'
734 l=1
735 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
736 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
737 ,l,impf,jmpf,nframe,pd)
738 if(debugprint)print*,'sample ',varname,' = ',pd(im/2,(jsta+jend)/2)
739
740! do j = jsta_2l, jend_2u
741! do i = 1, im
742! PD(I,J)=DUMMY2(I,J)
743! enddo
744! enddo
745
746 varname='hgt'
747 vcoordname='sfc'
748 l=1
749 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
750 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
751 ,l,impf,jmpf,nframe,fis)
752 if(debugprint)print*,'sample ',varname,' = ',fis(im/2,(jsta+jend)/2)
753 where(fis /= spval)fis=fis*g ! convert to geopotential
754
755 varname='tmp'
756 vcoordname='mid layer'
757 do l=1,lm
758! ll=lm-l+1
759 ll=l
760 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
761 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
762 ,l,impf,jmpf,nframe,t(1,jsta_2l,ll))
763 if(debugprint)print*,'sample l ',varname,' = ',ll,t(im/2,(jsta+jend)/2,ll)
764! if(debugprint)print*,'i=1 and im ',VarName,' = ',ll,t(1,(jsta+jend)/2,ll), &
765! t(im,(jsta+jend)/2,ll)
766 do i=1,im
767 do j=jsta,jend
768 if(t(i,j,ll)<150.)print*,'abnormal incoming T ',i,j,ll,t(i,j,ll)
769 end do
770 end do
771 end do ! do loop for l
772
773! Chuang: start mods to read variable for other mp schemes
774! model level q
775 if(imp_physics==5 .or. imp_physics==85 .or. imp_physics==95 .or. imp_physics==99)then
776 varname='spfh'
777 vcoordname='mid layer'
778 do l=1,lm
779! ll=lm-l+1
780 ll=l
781 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
782 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
783 ,l,impf,jmpf,nframe,q(1,jsta_2l,ll))
784 if(debugprint)print*,'sample l ',varname,' = ',ll,q(im/2,(jsta+jend)/2,ll)
785 end do ! do loop for l
786 end if
787
788! model level u
789 varname='ugrd'
790 vcoordname='mid layer'
791 do l=1,lm
792! ll=lm-l+1
793 ll=l
794 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
795 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
796 ,l,impf,jmpf,nframe,uh(1,jsta_2l,ll))
797 if(debugprint)print*,'sample l ',varname,' = ',ll,uh(im/2,(jsta+jend)/2,ll)
798! put u on h point for global nmm
799 if(global)then
800 buf(:,:)=uh(:,:,ll)
801 call exch(buf(1,jsta_2l))
802 if(debugprint)print*,'sample l u = ',ll,buf(im/2,(jsta+jend)/2)
803 do j=jsta,jend
804 do i=1,im
805 im1=i-1
806 if(im1<1)im1=im1+im
807 jm1=j-1
808 if(j==1)then
809 ii=i+im/2
810 iim1=ii-1
811 if(iim1<1)iim1=iim1+im
812 if (ii > im) ii = ii - im
813 uh(i,j,ll)=(buf(i,j)+buf(im1,j)+buf(ii,j)+buf(iim1,j))/4.0
814 else
815 uh(i,j,ll)=(buf(i,j)+buf(im1,j)+buf(im1,jm1)+buf(i,jm1))/4.0
816 end if
817 end do
818 end do
819 end if ! end of wind interpolation for global NMM
820 end do ! do loop for l
821
822! model level v
823 varname='vgrd'
824 vcoordname='mid layer'
825 do l=1,lm
826! ll=lm-l+1
827 ll=l
828 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
829 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
830 ,l,impf,jmpf,nframe,vh(1,jsta_2l,ll))
831 if(debugprint)print*,'sample l ',varname,' = ',ll,vh(im/2,(jsta+jend)/2,ll)
832! put u on h point for global nmm
833 if(global)then
834 buf(:,:)=vh(:,:,ll)
835 call exch(buf(1,jsta_2l))
836 if(debugprint)print*,'sample l v = ',ll,buf(im/2,(jsta+jend)/2)
837 do j=jsta,jend
838 do i=1,im
839 im1=i-1
840 if(im1<1)im1=im1+im
841 jm1=j-1
842 if(j==1)then
843 ii=i+im/2
844 iim1=ii-1
845 if(iim1<1)iim1=iim1+im
846 if (ii > im) ii = ii - im
847 vh(i,j,ll)=(buf(i,j)+buf(im1,j)+buf(ii,j)+buf(iim1,j))/4.0
848 else
849 vh(i,j,ll)=(buf(i,j)+buf(im1,j)+buf(im1,jm1)+buf(i,jm1))/4.0
850 end if
851 end do
852 end do
853 end if ! end of wind interpolation for global NMM
854
855 end do ! do loop for l
856
857 varname='sg1'
858 if(me==0)then
859 call nemsio_getheadvar(nfile,trim(varname),eta1,iret)
860 if (iret /= 0) then
861 print*,varname," not found in file-Assigned missing values"
862 eta1=spval
863 end if
864 end if
865 call mpi_bcast(eta1,lm+1,mpi_real,0,mpi_comm_comp,iret)
866
867 varname='sg2'
868 if(me==0)then
869 call nemsio_getheadvar(nfile,trim(varname),eta2,iret)
870 if (iret /= 0) then
871 print*,varname," not found in file-Assigned missing values"
872 eta2=spval
873 end if
874 end if
875 call mpi_bcast(eta2,lm+1,mpi_real,0,mpi_comm_comp,iret)
876
877 open(75,file='ETAPROFILE.txt',form='formatted', &
878 status='unknown')
879 DO l=1,lm+1
880 write(75,1020)l, eta1(lm+2-l), eta2(lm+2-l)
881 END DO
882 1020 format(i3,2e17.10)
883 close (75)
884
885 varname='pdtop'
886 if(me==0)then
887 call nemsio_getheadvar(nfile,trim(varname),pdtop,iret)
888 if (iret /= 0) then
889 print*,varname," not found in file-Assigned missing values"
890 pdtop=spval
891 end if
892 end if
893 call mpi_bcast(pdtop,1,mpi_real,0,mpi_comm_comp,iret)
894
895 varname='pt'
896 if(me==0)then
897 call nemsio_getheadvar(nfile,trim(varname),pt,iret)
898 if (iret /= 0) then
899 print*,varname," not found in file-Assigned missing values"
900 pt=spval
901 end if
902 end if
903 call mpi_bcast(pt,1,mpi_real,0,mpi_comm_comp,iret)
904 if(me==0) print*,'PT, PDTOP= ',pt,pdtop
905
906 varname='pblh'
907 vcoordname='sfc'
908 l=1
909 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
910 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
911 ,l,impf,jmpf,nframe,pblh)
912 if(debugprint)print*,'sample ',varname,' = ',pblh(im/2,(jsta+jend)/2)
913
914 varname='mixht'
915 vcoordname='sfc'
916 l=1
917 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
918 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
919 ,l,impf,jmpf,nframe,mixht)
920 if(debugprint)print*,'sample ',varname,' = ',mixht(im/2,(jsta+jend)/2)
921
922 varname='uustar'
923 vcoordname='sfc'
924 l=1
925 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
926 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
927 ,l,impf,jmpf,nframe,ustar)
928 if(debugprint)print*,'sample ',varname,' = ',ustar(im/2,(jsta+jend)/2)
929
930 varname='zorl'
931 vcoordname='sfc'
932 l=1
933 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
934 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
935 ,l,impf,jmpf,nframe,z0)
936 if(debugprint)print*,'sample ',varname,' = ',z0(im/2,(jsta+jend)/2)
937
938 varname='ths'
939 vcoordname='sfc'
940 l=1
941 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
942 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
943 ,l,impf,jmpf,nframe,ths)
944 if(debugprint)print*,'sample ',varname,' = ',ths(im/2,(jsta+jend)/2)
945
946 varname='qsh'
947 vcoordname='sfc'
948 l=1
949 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
950 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
951 ,l,impf,jmpf,nframe,qs)
952 if(debugprint)print*,'sample ',varname,' = ',qs(im/2,(jsta+jend)/2)
953
954 varname='twbs'
955 vcoordname='sfc'
956 l=1
957 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
958 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
959 ,l,impf,jmpf,nframe,twbs)
960 if(debugprint)print*,'sample ',varname,' = ',twbs(im/2,(jsta+jend)/2)
961
962 varname='qwbs'
963 vcoordname='sfc'
964 l=1
965 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
966 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
967 ,l,impf,jmpf,nframe,qwbs)
968 if(debugprint)print*,'sample ',varname,' = ',qwbs(im/2,(jsta+jend)/2)
969
970 varname='prec' ! instantaneous precip rate?
971 vcoordname='sfc'
972 l=1
973 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
974 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
975 ,l,impf,jmpf,nframe,prec)
976 if(debugprint)print*,'sample ',varname,' = ',prec(im/2,(jsta+jend)/2)
977
978 varname='acprec' ! accum total precip
979 vcoordname='sfc'
980 l=1
981 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
982 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
983 ,l,impf,jmpf,nframe,acprec)
984 if(debugprint)print*,'sample ',varname,' = ',acprec(im/2,(jsta+jend)/2)
985
986 varname='cuprec' ! accum cumulus precip
987 vcoordname='sfc'
988 l=1
989 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
990 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
991 ,l,impf,jmpf,nframe,cuprec)
992 if(debugprint)print*,'sample ',varname,' = ',cuprec(im/2,(jsta+jend)/2)
993
994! compute grid scale precip
995 do j=jsta,jend
996 do i=1,im
997 ancprc(i,j)=acprec(i,j)-cuprec(i,j)
998 end do
999 end do
1000
1001 varname='lspa'
1002 vcoordname='sfc'
1003 l=1
1004 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1005 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1006 ,l,impf,jmpf,nframe,lspa)
1007 if(debugprint)print*,'sample ',varname,' = ',lspa(im/2,(jsta+jend)/2)
1008
1009 varname='sno'
1010 vcoordname='sfc'
1011 l=1
1012 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1013 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1014 ,l,impf,jmpf,nframe,sno)
1015 if(debugprint)print*,'sample ',varname,' = ',sno(im/2,(jsta+jend)/2)
1016
1017 varname='snoavg'
1018 vcoordname='sfc'
1019 l=1
1020 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1021 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1022 ,l,impf,jmpf,nframe,snoavg)
1023 if(debugprint)print*,'sample ',varname,' = ',snoavg(im/2,(jsta+jend)/2)
1024
1025 varname='psfcavg'
1026 vcoordname='sfc'
1027 l=1
1028 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1029 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1030 ,l,impf,jmpf,nframe,psfcavg)
1031 if(debugprint)print*,'sample ',varname,' = ',psfcavg(im/2,(jsta+jend)/2)
1032
1033 varname='t10avg'
1034 vcoordname='10 m above gnd'
1035 l=1
1036 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1037 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1038 ,l,impf,jmpf,nframe,t10avg)
1039 if(debugprint)print*,'sample ',varname,' = ',t10avg(im/2,(jsta+jend)/2)
1040
1041 varname='t10'
1042 vcoordname='10 m above gnd'
1043 l=1
1044 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1045 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1046 ,l,impf,jmpf,nframe,t10m)
1047 if(debugprint)print*,'sample ',varname,' = ',t10m(im/2,(jsta+jend)/2)
1048
1049 varname='akhsavg'
1050 vcoordname='sfc'
1051 l=1
1052 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1053 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1054 ,l,impf,jmpf,nframe,akhsavg)
1055 if(debugprint)print*,'sample ',varname,' = ',akhsavg(im/2,(jsta+jend)/2)
1056
1057 varname='akmsavg'
1058 vcoordname='sfc'
1059 l=1
1060 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1061 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1062 ,l,impf,jmpf,nframe,akmsavg)
1063 if(debugprint)print*,'sample ',varname,' = ',akmsavg(im/2,(jsta+jend)/2)
1064
1065 varname='refdmax'
1066 vcoordname='sfc'
1067 l=1
1068 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1069 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1070 ,l,impf,jmpf,nframe,refd_max)
1071 if(debugprint)print*,'sample ',varname,' = ',refd_max(im/2,(jsta+jend)/2)
1072
1073 varname='upvvelmax'
1074 vcoordname='sfc' ! wrong
1075 l=1
1076 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1077 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1078 ,l,impf,jmpf,nframe,w_up_max)
1079 if(debugprint)print*,'sample ',varname,' = ',w_up_max(im/2,(jsta+jend)/2)
1080
1081 varname='dnvvelmax'
1082 vcoordname='sfc' ! wrong
1083 l=1
1084 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1085 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1086 ,l,impf,jmpf,nframe,w_dn_max)
1087 if(debugprint)print*,'sample ',varname,' = ',w_dn_max(im/2,(jsta+jend)/2)
1088
1089 varname='uphlmax'
1090 vcoordname='sfc' ! wrong
1091 l=1
1092 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1093 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1094 ,l,impf,jmpf,nframe,up_heli_max)
1095 if(debugprint)print*,'sample ',varname,' = ',up_heli_max(im/2,(jsta+jend)/2)
1096
1097 varname='pratemax' ! max hourly precip rate
1098 vcoordname='sfc' ! wrong
1099 l=1
1100 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1101 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1102 ,l,impf,jmpf,nframe,prate_max)
1103 if(debugprint)print*,'sample ',varname,' = ',prate_max(im/2,(jsta+jend)/2)
1104
1105 varname='fpratemax' ! max hourly frozen precip rate
1106 vcoordname='sfc' ! wrong
1107 l=1
1108 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1109 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1110 ,l,impf,jmpf,nframe,fprate_max)
1111 if(debugprint)print*,'sample ',varname,' = ',fprate_max(im/2,(jsta+jend)/2)
1112
1113 varname='si'
1114 vcoordname='sfc'
1115 l=1
1116 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1117 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1118 ,l,impf,jmpf,nframe,si)
1119 if(debugprint)print*,'sample ',varname,' = ',si(im/2,(jsta+jend)/2)
1120
1121 varname='cldefi'
1122 vcoordname='sfc'
1123 l=1
1124 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1125 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1126 ,l,impf,jmpf,nframe,cldefi)
1127 if(debugprint)print*,'sample ',varname,' = ',cldefi(im/2,(jsta+jend)/2)
1128
1129 varname='th10'
1130 vcoordname='10 m above gnd'
1131 l=1
1132 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1133 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1134 ,l,impf,jmpf,nframe,th10)
1135 if(debugprint)print*,'sample ',varname,' = ',th10(im/2,(jsta+jend)/2)
1136
1137 varname='q10'
1138 vcoordname='10 m above gnd'
1139 l=1
1140 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1141 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1142 ,l,impf,jmpf,nframe,q10)
1143 if(debugprint)print*,'sample ',varname,' = ',q10(im/2,(jsta+jend)/2)
1144
1145 varname='pshltr'
1146 vcoordname='sfc'
1147 l=1
1148 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1149 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1150 ,l,impf,jmpf,nframe,pshltr)
1151 if(debugprint)print*,'sample ',varname,' = ',pshltr(im/2,(jsta+jend)/2), &
1152 'max=',maxval(pshltr(1:im,jsta:jend)),minval(pshltr(1:im,jsta:jend))
1153
1154 varname='tshltr'
1155 vcoordname='sfc'
1156 l=1
1157 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1158 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1159 ,l,impf,jmpf,nframe,tshltr)
1160 if(debugprint)print*,'sample ',varname,' = ',tshltr(im/2,(jsta+jend)/2)
1161
1162 varname='qshltr'
1163 vcoordname='sfc'
1164 l=1
1165 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1166 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1167 ,l,impf,jmpf,nframe,qshltr)
1168 if(debugprint)print*,'sample ',varname,' = ',qshltr(im/2,(jsta+jend)/2)
1169
1170 tmaxmin=1.
1171! call mpi_bcast(TMAXMIN,1,MPI_REAL,0,mpi_comm_comp,iret)
1172
1173 varname='t02max'
1174 vcoordname='sfc'
1175 l=1
1176 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1177 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1178 ,l,impf,jmpf,nframe,maxtshltr)
1179 if(debugprint)print*,'sample ',varname,' = ',maxtshltr(im/2,(jsta+jend)/2)
1180
1181 varname='t02min'
1182 vcoordname='sfc'
1183 l=1
1184 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1185 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1186 ,l,impf,jmpf,nframe,mintshltr)
1187 if(debugprint)print*,'sample ',varname,' = ',mintshltr(im/2,(jsta+jend)/2)
1188
1189 varname='rh02max'
1190 vcoordname='sfc'
1191 l=1
1192 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1193 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1194 ,l,impf,jmpf,nframe,maxrhshltr)
1195 if(debugprint)print*,'sample ',varname,' = ',maxrhshltr(im/2,(jsta+jend)/2)
1196
1197 varname='rh02min'
1198 vcoordname='sfc'
1199 l=1
1200 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1201 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1202 ,l,impf,jmpf,nframe,minrhshltr)
1203 if(debugprint)print*,'sample ',varname,' = ',minrhshltr(im/2,(jsta+jend)/2)
1204
1205! model level q2
1206 varname='q2'
1207 vcoordname='mid layer'
1208 do l=1,lm
1209! ll=lm-l+1
1210 ll=l
1211 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1212 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1213 ,l,impf,jmpf,nframe,q2(1,jsta_2l,ll))
1214 if(debugprint)print*,'sample l ',varname,' = ',ll,q2(im/2,(jsta+jend)/2,ll)
1215 end do ! do loop for l
1216
1217 varname='akhs_out'
1218 vcoordname='sfc'
1219 l=1
1220 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1221 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1222 ,l,impf,jmpf,nframe,akhs)
1223 if(debugprint)print*,'sample ',varname,' = ',akhs(im/2,(jsta+jend)/2)
1224
1225 varname='akms_out'
1226 vcoordname='sfc'
1227 l=1
1228 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1229 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1230 ,l,impf,jmpf,nframe,akms)
1231 if(debugprint)print*,'sample ',varname,' = ',akms(im/2,(jsta+jend)/2)
1232
1233 varname='albase'
1234 vcoordname='sfc'
1235 l=1
1236 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1237 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1238 ,l,impf,jmpf,nframe,albase)
1239 if(debugprint)print*,'sample ',varname,' = ',albase(im/2,(jsta+jend)/2)
1240
1241 varname='albedo'
1242 vcoordname='sfc'
1243 l=1
1244 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1245 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1246 ,l,impf,jmpf,nframe,albedo)
1247 if(debugprint)print*,'sample ',varname,' = ',albedo(im/2,(jsta+jend)/2)
1248
1249 varname='czen'
1250 vcoordname='sfc'
1251 l=1
1252 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1253 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1254 ,l,impf,jmpf,nframe,czen)
1255 if(debugprint)print*,'sample ',varname,' = ',czen(im/2,(jsta+jend)/2)
1256
1257 varname='czmean'
1258 vcoordname='sfc'
1259 l=1
1260 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1261 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1262 ,l,impf,jmpf,nframe,czmean)
1263 if(debugprint)print*,'sample ',varname,' = ',czmean(im/2,(jsta+jend)/2)
1264
1265 varname='mxsnal'
1266 vcoordname='sfc'
1267 l=1
1268 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1269 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1270 ,l,impf,jmpf,nframe,mxsnal)
1271 if(debugprint)print*,'sample ',varname,' = ',mxsnal(im/2,(jsta+jend)/2)
1272
1273 varname='radot'
1274 vcoordname='sfc'
1275 l=1
1276 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1277 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1278 ,l,impf,jmpf,nframe,radot)
1279 if(debugprint)print*,'sample ',varname,' = ',radot(im/2,(jsta+jend)/2)
1280
1281 varname='sigt4'
1282 vcoordname='sfc'
1283 l=1
1284 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1285 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1286 ,l,impf,jmpf,nframe,sigt4)
1287 if(debugprint)print*,'sample ',varname,' = ',sigt4(im/2,(jsta+jend)/2)
1288
1289 varname='tg'
1290 vcoordname='sfc'
1291 l=1
1292 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1293 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1294 ,l,impf,jmpf,nframe,tg)
1295 if(debugprint)print*,'sample ',varname,' = ',tg(im/2,(jsta+jend)/2)
1296
1297! Chuang: mods to read variables from other mp schemes
1298 if(imp_physics==5 .or. imp_physics==85 .or. imp_physics==95 .or. imp_physics==99)then
1299! model level cwm
1300 varname='clwmr'
1301 vcoordname='mid layer'
1302 do l=1,lm
1303! ll=lm-l+1
1304 ll=l
1305 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1306 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1307 ,l,impf,jmpf,nframe,cwm(1,jsta_2l,ll))
1308 if(debugprint)print*,'sample l ',varname,' = ',ll,cwm(im/2,(jsta+jend)/2,ll)
1309 end do ! do loop for l
1310
1311 varname='f_ice'
1312 vcoordname='mid layer'
1313 do l=1,lm
1314! ll=lm-l+1
1315 ll=l
1316 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1317 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1318 ,l,impf,jmpf,nframe,f_ice(1,jsta_2l,ll))
1319 if(debugprint)print*,'sample l ',varname,' = ',ll,f_ice(im/2,(jsta+jend)/2,ll)
1320 if(debugprint)print*,'max min ',varname,' = ',ll,maxval(f_ice(:,:,ll)),minval(f_ice(:,:,ll))
1321 end do ! do loop for l
1322
1323 varname='f_rain'
1324 vcoordname='mid layer'
1325 do l=1,lm
1326! ll=lm-l+1
1327 ll=l
1328 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1329 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1330 ,l,impf,jmpf,nframe,f_rain(1,jsta_2l,ll))
1331 if(debugprint)print*,'sample l ',varname,' = ',ll,f_rain(im/2,(jsta+jend)/2,ll)
1332 if(debugprint)print*,'max min ',varname,' = ',ll,maxval(f_rain(:,:,ll)),minval(f_rain(:,:,ll))
1333 end do ! do loop for l
1334
1335 if(imp_physics /= 99)then
1336 varname='f_rimef'
1337 vcoordname='mid layer'
1338 do l=1,lm
1339! ll=lm-l+1
1340 ll=l
1341 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1342 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1343 ,l,impf,jmpf,nframe,f_rimef(1,jsta_2l,ll))
1344 if(debugprint)print*,'sample l ',varname,' = ',ll,f_rimef(im/2,(jsta+jend)/2,ll)
1345 if(debugprint)print*,'max min ',varname,' = ',ll,maxval(f_rimef(:,:,ll)),minval(f_rimef(:,:,ll))
1346 end do ! do loop for l
1347 endif
1348
1349 if(imp_physics==5)then
1350 varname='refl_10cm'
1351 vcoordname='mid layer'
1352 do l=1,lm
1353 ll=l
1354 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1355 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1356 ,l,impf,jmpf,nframe,ref_10cm(1,jsta_2l,ll))
1357 if(debugprint)print*,'sample l ',varname,' = ' &
1358 ,ll,ref_10cm(im/2,(jsta+jend)/2,ll)
1359 if(debugprint)print*,'max min ',varname,' = ' &
1360 ,ll,maxval(ref_10cm(:,:,ll)),minval(ref_10cm(:,:,ll))
1361 end do ! do loop for l
1362 endif
1363
1364 else ! retrieve hydrometeo fields directly for non-Ferrier
1365 if(imp_physics==8 .or. imp_physics==9)then
1366 varname='ni'
1367 vcoordname='mid layer'
1368 do l=1,lm
1369 ll=l
1370 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1371 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1372 ,l,impf,jmpf,nframe,qqni(1,jsta_2l,ll))
1373 if(debugprint)print*,'sample l ',varname,' = ' &
1374 ,ll,qqni(im/2,(jsta+jend)/2,ll)
1375 if(debugprint)print*,'max min ',varname,' = ' &
1376 ,ll,maxval(qqni(:,:,ll)),minval(qqni(:,:,ll))
1377 end do ! do loop for l
1378
1379 varname='nr'
1380 vcoordname='mid layer'
1381 do l=1,lm
1382 ll=l
1383 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1384 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1385 ,l,impf,jmpf,nframe,qqnr(1,jsta_2l,ll))
1386 if(debugprint)print*,'sample l ',varname,' = ' &
1387 ,ll,qqnr(im/2,(jsta+jend)/2,ll)
1388 if(debugprint)print*,'max min ',varname,' = ' &
1389 ,ll,maxval(qqnr(:,:,ll)),minval(qqnr(:,:,ll))
1390 end do ! do loop for l
1391 end if !imp_physics==8 or 9
1392
1393! water vapor, will be converted to specific humidity
1394 varname='qv'
1395 vcoordname='mid layer'
1396 do l=1,lm
1397 ll=l
1398 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1399 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1400 ,l,impf,jmpf,nframe,q(1,jsta_2l,ll))
1401!$omp parallel do private(i,j)
1402 do j=jsta_2l,jend_2u
1403 do i=1,im
1404 q(i,j,ll) = max(10e-12,q(i,j,ll))
1405 q(i,j,ll) = q(i,j,ll) / (1.0+q(i,j,ll))
1406 end do
1407 end do
1408 if(debugprint)print*,'sample l ',varname,' = ' &
1409 ,ll,q(im/2,(jsta+jend)/2,ll)
1410 if(debugprint)print*,'max min ',varname,' = ' &
1411 ,ll,maxval(q(:,:,ll)),minval(q(:,:,ll))
1412 end do ! ,j,l) do loop for l
1413
1414!$omp parallel do private(i,j,l)
1415 do l=1,lm
1416 do j=jsta_2l,jend_2u
1417 do i=1,im
1418 qqw(i,j,l) = spval
1419 qqr(i,j,l) = spval
1420 qqs(i,j,l) = spval
1421 qqi(i,j,l) = spval
1422 qqg(i,j,l) = spval
1423 cwm(i,j,l) = spval
1424 enddo
1425 enddo
1426 enddo
1427! cloud water
1428 if(imp_physics/=0)then
1429 varname='qc'
1430 vcoordname='mid layer'
1431 do l=1,lm
1432 ll=l
1433 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1434 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1435 ,l,impf,jmpf,nframe,qqw(1,jsta_2l,ll))
1436 if(imp_physics==3)then !WSM3
1437!$omp parallel do private(i,j)
1438 do j=jsta_2l,jend_2u
1439 do i=1,im
1440 if(t(i,j,ll)<tfrz) qqi(i,j,ll) = qqw(i,j,ll)
1441 end do
1442 end do
1443 end if
1444!$omp parallel do private(i,j)
1445 do j=jsta_2l,jend_2u
1446 do i=1,im
1447 cwm(i,j,ll) = qqw(i,j,ll)
1448 end do
1449 end do
1450 if(debugprint)print*,'sample l ',varname,' = ' &
1451 ,ll,qqw(im/2,(jsta+jend)/2,ll)
1452 if(debugprint)print*,'max min ',varname,' = ' &
1453 ,ll,maxval(qqw(:,:,ll)),minval(qqw(:,:,ll))
1454 end do ! do loop for l
1455 end if ! for non-zero mp_physics
1456
1457! cloud ice
1458 if(imp_physics/=0 .and. imp_physics/=1 .and. imp_physics/=3)then
1459 varname='qi'
1460 vcoordname='mid layer'
1461 do l=1,lm
1462 ll=l
1463 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1464 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1465 ,l,impf,jmpf,nframe,qqi(1,jsta_2l,ll))
1466!$omp parallel do private(i,j)
1467 do j=jsta_2l,jend_2u
1468 do i=1,im
1469 cwm(i,j,ll) = cwm(i,j,ll) + qqi(i,j,ll)
1470 end do
1471 end do
1472 if(debugprint)print*,'sample l ',varname,' = ' &
1473 ,ll,qqi(im/2,(jsta+jend)/2,ll)
1474 if(debugprint)print*,'max min ',varname,' = ' &
1475 ,ll,maxval(qqi(:,:,ll)),minval(qqi(:,:,ll))
1476 end do ! do loop for l
1477 end if
1478! rain
1479 if(imp_physics/=0)then
1480 varname='qr'
1481 vcoordname='mid layer'
1482 do l=1,lm
1483 ll=l
1484 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1485 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1486 ,l,impf,jmpf,nframe,qqr(1,jsta_2l,ll))
1487 if(imp_physics==3)then !WSM3
1488!$omp parallel do private(i,j)
1489 do j=jsta_2l,jend_2u
1490 do i=1,im
1491 if(t(i,j,ll)<tfrz) qqs(i,j,ll) = qqr(i,j,ll)
1492 end do
1493 end do
1494 end if
1495!$omp parallel do private(i,j)
1496 do j=jsta_2l,jend_2u
1497 do i=1,im
1498 cwm(i,j,ll) = cwm(i,j,ll) + qqr(i,j,ll)
1499 end do
1500 end do
1501 if(debugprint)print*,'sample l ',varname,' = ' &
1502 ,ll,qqr(im/2,(jsta+jend)/2,ll)
1503 if(debugprint)print*,'max min ',varname,' = ' &
1504 ,ll,maxval(qqr(:,:,ll)),minval(qqr(:,:,ll))
1505 end do ! do loop for l
1506 end if ! for non-zero mp_physics
1507!snow
1508 if(imp_physics/=0 .and. imp_physics/=1 .and. imp_physics/=3)then
1509 varname='qs'
1510 vcoordname='mid layer'
1511 do l=1,lm
1512 ll=l
1513 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1514 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1515 ,l,impf,jmpf,nframe,qqs(1,jsta_2l,ll))
1516!$omp parallel do private(i,j)
1517 do j=jsta_2l,jend_2u
1518 do i=1,im
1519 cwm(i,j,ll) = cwm(i,j,ll) + qqs(i,j,ll)
1520 end do
1521 end do
1522 if(debugprint)print*,'sample l ',varname,' = ' &
1523 ,ll,qqs(im/2,(jsta+jend)/2,ll)
1524 if(debugprint)print*,'max min ',varname,' = ' &
1525 ,ll,maxval(qqs(:,:,ll)),minval(qqs(:,:,ll))
1526 end do ! do loop for l
1527 end if
1528!graupel
1529 if(imp_physics==2 .or. imp_physics==6 .or. imp_physics==8)then
1530 varname='qg'
1531 vcoordname='mid layer'
1532 do l=1,lm
1533 ll=l
1534 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1535 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1536 ,l,impf,jmpf,nframe,qqg(1,jsta_2l,ll))
1537!$omp parallel do private(i,j)
1538 do j=jsta_2l,jend_2u
1539 do i=1,im
1540 cwm(i,j,ll) = cwm(i,j,ll) + qqg(i,j,ll)
1541 end do
1542 end do
1543 if(debugprint)print*,'sample l ',varname,' = ' &
1544 ,ll,qqg(im/2,(jsta+jend)/2,ll)
1545 if(debugprint)print*,'max min ',varname,' = ' &
1546 ,ll,maxval(qqg(:,:,ll)),minval(qqg(:,:,ll))
1547 end do ! do loop for l
1548 end if
1549! radar ref and effective radius for Thompson
1550 if(imp_physics==8)then
1551 varname='refl_10cm'
1552 vcoordname='mid layer'
1553 do l=1,lm
1554 ll=l
1555 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1556 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1557 ,l,impf,jmpf,nframe,ref_10cm(1,jsta_2l,ll))
1558 if(debugprint)print*,'sample l ',varname,' = ' &
1559 ,ll,ref_10cm(im/2,(jsta+jend)/2,ll)
1560 if(debugprint)print*,'max min ',varname,' = ' &
1561 ,ll,maxval(ref_10cm(:,:,ll)),minval(ref_10cm(:,:,ll))
1562 end do ! do loop for l
1563! effective radius
1564 varname='re_cloud'
1565 vcoordname='mid layer'
1566 do l=1,lm
1567 ll=l
1568 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1569 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1570 ,l,impf,jmpf,nframe,radius_cloud(1,jsta_2l,ll))
1571 if(debugprint)print*,'sample l ',varname,' = ' &
1572 ,ll,radius_cloud(im/2,(jsta+jend)/2,ll)
1573 if(debugprint)print*,'max min ',varname,' = ' &
1574 ,ll,maxval(radius_cloud(:,:,ll)),minval(radius_cloud(:,:,ll))
1575 end do ! do loop for l
1576
1577 varname='re_ice'
1578 vcoordname='mid layer'
1579 do l=1,lm
1580 ll=l
1581 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1582 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1583 ,l,impf,jmpf,nframe,radius_ice(1,jsta_2l,ll))
1584 if(debugprint)print*,'sample l ',varname,' = ' &
1585 ,ll,radius_ice(im/2,(jsta+jend)/2,ll)
1586 if(debugprint)print*,'max min ',varname,' = ' &
1587 ,ll,maxval(radius_ice(:,:,ll)),minval(radius_ice(:,:,ll))
1588 end do ! do loop for l
1589
1590 varname='re_snow'
1591 vcoordname='mid layer'
1592 do l=1,lm
1593 ll=l
1594 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1595 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1596 ,l,impf,jmpf,nframe,radius_snow(1,jsta_2l,ll))
1597 if(debugprint)print*,'sample l ',varname,' = ' &
1598 ,ll,radius_snow(im/2,(jsta+jend)/2,ll)
1599 if(debugprint)print*,'max min ',varname,' = ' &
1600 ,ll,maxval(radius_snow(:,:,ll)),minval(radius_snow(:,:,ll))
1601 end do ! do loop for l
1602
1603 end if
1604 ! end of retrieving hydrometeo for non Ferrier's MP options
1605 end if ! end of retrieving hydrometeo for all MP options
1606
1607
1608 varname='cldfra'
1609 vcoordname='mid layer'
1610 do l=1,lm
1611! ll=lm-l+1
1612 ll=l
1613 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1614 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1615 ,l,impf,jmpf,nframe,cfr(1,jsta_2l,ll))
1616 if(debugprint)print*,'sample l ',varname,' = ',ll,cfr(im/2,(jsta+jend)/2,ll)
1617 end do ! do loop for l
1618
1619 varname='sr'
1620 vcoordname='sfc'
1621 l=1
1622 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1623 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1624 ,l,impf,jmpf,nframe,sr)
1625 if(debugprint)print*,'sample ',varname,' = ',sr(im/2,(jsta+jend)/2)
1626
1627 varname='cfrach'
1628 vcoordname='sfc'
1629 l=1
1630 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1631 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1632 ,l,impf,jmpf,nframe,cfrach)
1633 if(debugprint)print*,'sample ',varname,' = ',cfrach(im/2,(jsta+jend)/2)
1634
1635 varname='cfracl'
1636 vcoordname='sfc'
1637 l=1
1638 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1639 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1640 ,l,impf,jmpf,nframe,cfracl)
1641 if(debugprint)print*,'sample ',varname,' = ',cfracl(im/2,(jsta+jend)/2)
1642
1643 varname='cfracm'
1644 vcoordname='sfc'
1645 l=1
1646 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1647 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1648 ,l,impf,jmpf,nframe,cfracm)
1649 if(debugprint)print*,'sample ',varname,' = ',cfracm(im/2,(jsta+jend)/2)
1650
1651 varname='islope' !???
1652 vcoordname='sfc'
1653 l=1
1654 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1655 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1656 ,l,impf,jmpf,nframe,islope)
1657 if(debugprint)print*,'sample ',varname,' = ',islope(im/2,(jsta+jend)/2)
1658
1659! either assign SLDPTH to be the same as eta (which is original
1660! setup in WRF LSM) or extract thickness of soil layers from wrf
1661! output
1662
1663! or get SLDPTH from wrf output
1664 varname='sldpth'
1665 if(me==0)then
1666 call nemsio_getheadvar(nfile,trim(varname),sldpth,iret)
1667! if (iret /= 0) then
1668! print*,VarName," not found in file-Assigned NOAH default values"
1669! SLDPTH(1)=0.10
1670! SLDPTH(2)=0.3
1671! SLDPTH(3)=0.6
1672! SLDPTH(4)=1.0
1673! end if
1674 end if
1675 call mpi_bcast(sldpth,4,mpi_real,0,mpi_comm_comp,iret)
1676 if(me==0)print*,'SLDPTH= ',(sldpth(n),n=1,nsoil)
1677
1678 varname='cmc'
1679 vcoordname='sfc'
1680 l=1
1681 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1682 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1683 ,l,impf,jmpf,nframe,cmc)
1684 if(debugprint)print*,'sample ',varname,' = ',cmc(im/2,(jsta+jend)/2)
1685
1686 varname='grnflx'
1687 vcoordname='sfc'
1688 l=1
1689 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1690 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1691 ,l,impf,jmpf,nframe,grnflx)
1692 if(debugprint)print*,'sample ',varname,' = ',grnflx(im/2,(jsta+jend)/2)
1693
1694 varname='pctsno'
1695 vcoordname='sfc'
1696 l=1
1697 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1698 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1699 ,l,impf,jmpf,nframe,pctsno)
1700 if(debugprint)print*,'sample ',varname,' = ',pctsno(im/2,(jsta+jend)/2)
1701
1702 varname='soiltb'
1703 vcoordname='sfc'
1704 l=1
1705 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1706 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1707 ,l,impf,jmpf,nframe,soiltb)
1708 if(debugprint)print*,'sample ',varname,' = ',soiltb(im/2,(jsta+jend)/2)
1709 if(debugprint)then
1710 do j=jsta,jend
1711 do i=1,im
1712 if(soiltb(i,j)>350.)print*,'large soiltb='
1713 end do
1714 end do
1715 end if
1716
1717 varname='vegfrc'
1718 vcoordname='sfc'
1719 l=1
1720 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1721 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1722 ,l,impf,jmpf,nframe,vegfrc)
1723 if(debugprint)print*,'sample ',varname,' = ',vegfrc(im/2,(jsta+jend)/2)
1724
1725 varname='sh2o'
1726 vcoordname='soil layer'
1727 l=1
1728 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1729 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1730 ,l,impf,jmpf,nframe,sh2o(1,jsta_2l,l))
1731 if(debugprint)print*,'sample l ',varname,' = ',l,sh2o(im/2,(jsta+jend)/2,l)
1732
1733 varname='sh2o'
1734 vcoordname='soil layer'
1735 l=2
1736 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1737 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1738 ,l,impf,jmpf,nframe,sh2o(1,jsta_2l,l))
1739 if(debugprint)print*,'sample l ',varname,' = ',l,sh2o(im/2,(jsta+jend)/2,l)
1740
1741 varname='sh2o'
1742 vcoordname='soil layer'
1743 l=3
1744 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1745 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1746 ,l,impf,jmpf,nframe,sh2o(1,jsta_2l,l))
1747 if(debugprint)print*,'sample l ',varname,' = ',l,sh2o(im/2,(jsta+jend)/2,l)
1748
1749 varname='sh2o'
1750 vcoordname='soil layer'
1751 l=4
1752 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1753 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1754 ,l,impf,jmpf,nframe,sh2o(1,jsta_2l,l))
1755 if(debugprint)print*,'sample l ',varname,' = ',l,sh2o(im/2,(jsta+jend)/2,l)
1756
1757 varname='smc'
1758 vcoordname='soil layer'
1759 l=1
1760 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1761 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1762 ,l,impf,jmpf,nframe,smc(1,jsta_2l,l))
1763 if(debugprint)print*,'sample l ',varname,' = ',l,smc(im/2,(jsta+jend)/2,l)
1764
1765 varname='smc'
1766 vcoordname='soil layer'
1767 l=2
1768 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1769 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1770 ,l,impf,jmpf,nframe,smc(1,jsta_2l,l))
1771 if(debugprint)print*,'sample l ',varname,' = ',l,smc(im/2,(jsta+jend)/2,l)
1772
1773 varname='smc'
1774 vcoordname='soil layer'
1775 l=3
1776 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1777 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1778 ,l,impf,jmpf,nframe,smc(1,jsta_2l,l))
1779 if(debugprint)print*,'sample l ',varname,' = ',l,smc(im/2,(jsta+jend)/2,l)
1780
1781 varname='smc'
1782 vcoordname='soil layer'
1783 l=4
1784 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1785 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1786 ,l,impf,jmpf,nframe,smc(1,jsta_2l,l))
1787 if(debugprint)print*,'sample l ',varname,' = ',l,smc(im/2,(jsta+jend)/2,l)
1788
1789 varname='stc'
1790 vcoordname='soil layer'
1791 l=1
1792 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1793 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1794 ,l,impf,jmpf,nframe,stc(1,jsta_2l,l))
1795 if(debugprint)print*,'sample l ',varname,' = ',l,stc(im/2,(jsta+jend)/2,l)
1796
1797 varname='stc'
1798 vcoordname='soil layer'
1799 l=2
1800 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1801 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1802 ,l,impf,jmpf,nframe,stc(1,jsta_2l,l))
1803 if(debugprint)print*,'sample l ',varname,' = ',l,stc(im/2,(jsta+jend)/2,l)
1804
1805 varname='stc'
1806 vcoordname='soil layer'
1807 l=3
1808 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1809 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1810 ,l,impf,jmpf,nframe,stc(1,jsta_2l,l))
1811 if(debugprint)print*,'sample l ',varname,' = ',l,stc(im/2,(jsta+jend)/2,l)
1812
1813 varname='stc'
1814 vcoordname='soil layer'
1815 l=4
1816 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1817 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1818 ,l,impf,jmpf,nframe,stc(1,jsta_2l,l))
1819 if(debugprint)print*,'sample l ',varname,' = ',l,stc(im/2,(jsta+jend)/2,l)
1820
1821 varname='pres'
1822 vcoordname='layer'
1823 do l=1,lp1
1824! ll=lp1-l+1
1825 ll=l
1826 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1827 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1828 ,l,impf,jmpf,nframe,pint(1,jsta_2l,ll))
1829 if(debugprint)print*,'sample l ',varname,' = ',ll,pint(im/2,(jsta+jend)/2,ll)
1830 if(l /= 1)then ! assuming post counts from top down
1831 do j=jsta,jend
1832 do i=1,im
1833 alpint(i,j,ll)=alog(pint(i,j,ll))
1834 end do
1835 end do
1836 end if
1837 end do ! do loop for l
1838
1839! do l = 1, lp1
1840 l=1
1841 do j = jsta, jend
1842 do i = 1, im
1843 if(pint(i,j,l) /= 0.0)then
1844 alpint(i,j,l)=alog(pint(i,j,l))
1845 else
1846 alpint(i,j,l)=spval
1847 end if
1848 end do
1849 end do
1850! end do
1851
1852 do l = 2, lp1
1853 do j = jsta_2l, jend_2u
1854 do i = 1, im
1855 pmid(i,j,l-1 ) = (pint(i,j,l-1)+ &
1856 pint(i,j,l))*0.5 ! representative of what model does
1857 end do
1858 end do
1859 if(debugprint)print*,'sample l, PMID = ',l-1,pmid(im/2,(jsta+jend)/2,l-1)
1860 end do
1861
1862 if(gridtype=='E')then
1863 do l = 1, lm
1864 call exch(pmid(1:im,jsta_2l:jend_2u,l))
1865 do j = jsta, jend
1866 do i = 1, im-mod(j,2)
1867 IF(j == 1 .AND. i < im)THEN !SOUTHERN BC
1868 pmidv(i,j,l)=0.5*(pmid(i,j,l)+pmid(i+1,j,l))
1869 ELSE IF(j==jm .AND. i<im)THEN !NORTHERN BC
1870 pmidv(i,j,l)=0.5*(pmid(i,j,l)+pmid(i+1,j,l))
1871 ELSE IF(i == 1 .AND. mod(j,2) == 0) THEN !WESTERN EVEN BC
1872 pmidv(i,j,l)=0.5*(pmid(i,j-1,l)+pmid(i,j+1,l))
1873 ELSE IF(i == im .AND. mod(j,2) == 0 &
1874 .AND. j < jm) THEN !EASTERN EVEN BC
1875 pmidv(i,j,l)=0.5*(pmid(i,j-1,l)+pmid(i,j+1,l))
1876 ELSE IF (mod(j,2) < 1) THEN
1877 pmidv(i,j,l)=0.25*(pmid(i,j,l)+pmid(i-1,j,l) &
1878 +pmid(i,j+1,l)+pmid(i,j-1,l))
1879 ELSE
1880 pmidv(i,j,l)=0.25*(pmid(i,j,l)+pmid(i+1,j,l) &
1881 +pmid(i,j+1,l)+pmid(i,j-1,l))
1882 END IF
1883 end do
1884 end do
1885 end do
1886 else if(gridtype=='B')then
1887 do l = 1, lm
1888 call exch(pmid(1:im,jsta_2l:jend_2u,l))
1889 do j = jsta, jend_m
1890 do i = 1, im-1
1891 pmidv(i,j,l)=0.25*(pmid(i,j,l) + pmid(i+1,j,l) &
1892 + pmid(i,j+1,l) + pmid(i+1,j+1,l))
1893 end do
1894 end do
1895 end do
1896 end if
1897
1898!!!!! COMPUTE Z
1899 allocate(fi(im,jsta_2l:jend_2u,2))
1900 do j = jsta, jend
1901 do i = 1, im
1902 zint(i,j,lm+1) = fis(i,j)/g
1903 if (debugprint .and. i == im/2 .and. j ==(jsta+jend)/2 ) then
1904 write(6,*) 'G,ZINT: ', g,zint(i,j,lm+1)
1905 endif
1906 fi(i,j,1) = fis(i,j)
1907 end do
1908 end do
1909
1910! SECOND, INTEGRATE HEIGHT HYDROSTATICLY
1911 ii = im/2
1912 jj = (jsta+jend)/2
1913 DO l=lm,1,-1
1914!$omp parallel do private(i,j)
1915 do j = jsta, jend
1916 do i = 1, im
1917 fi(i,j,2) = htm(i,j,l)*t(i,j,l)*(q(i,j,l)*d608+1.0)*rd* &
1918 (alpint(i,j,l+1)-alpint(i,j,l))+fi(i,j,1)
1919 zint(i,j,l) = fi(i,j,2)/g
1920 fi(i,j,1) = fi(i,j,2)
1921 ENDDO
1922 ENDDO
1923 if(debugprint)print*,'L,sample HTM,T,Q,ALPINT(L+1),ALPINT(l)', &
1924 ',ZINT= ',l,htm(ii,jj,l),t(ii,jj,l), &
1925 q(ii,jj,l),alpint(ii,jj,l+1), &
1926 alpint(ii,jj,l),zint(ii,jj,l)
1927 END DO
1928 if (me==0) then
1929 print*,'finish deriving geopotential in nmm'
1930 write(*,*)' after ZINT lm=',lm,' js=',js,' je=',je,' im=',im
1931 write(*,*)' zmid lbounds=',lbound(zmid),' ubounds=',ubound(zmid)
1932 write(*,*)' zint lbounds=',lbound(zint),' ubounds=',ubound(zint)
1933 write(*,*)' pmid lbounds=',lbound(pmid),' ubounds=',ubound(pmid)
1934 write(*,*)' pint lbounds=',lbound(pint),' ubounds=',ubound(pint)
1935 endif
1936 deallocate(fi)
1937!
1938 DO l=1,lm
1939! write(*,*)' zmid l=',l
1940!$omp parallel do private(i,j,fact)
1941 DO j=jsta,jend
1942! write(*,*)' zmid j=',j
1943 DO i=1,im
1944! write(*,*)' zmid i=',i
1945! ZMID(I,J,L)=(ZINT(I,J,L+1)+ZINT(I,J,L))*0.5 ! ave of z
1946! write(*,*)' pmid=',pmid(i,j,l)
1947! write(*,*)' pint=',pint(i,j,l),pint(i,j,l+1)
1948! write(*,*)' zint=',zint(i,j,l),zint(i,j,l+1)
1949 fact = (log(pmid(i,j,l))-log(pint(i,j,l))) &
1950 / (log(pint(i,j,l+1))-log(pint(i,j,l)))
1951 zmid(i,j,l) = zint(i,j,l) + (zint(i,j,l+1)-zint(i,j,l))*fact
1952 ENDDO
1953 ENDDO
1954 ENDDO
1955
1956 varname='vvel'
1957 vcoordname='mid layer'
1958 do l=1,lm
1959! ll=lm-l+1
1960 ll=l
1961 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1962 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1963 ,l,impf,jmpf,nframe,wh(1,jsta_2l,ll))
1964 if(debugprint)print*,'sample l ',varname,' = ',ll,wh(im/2,(jsta+jend)/2,ll)
1965 end do ! do loop for l
1966
1967 varname='acfrcv'
1968 vcoordname='sfc'
1969 l=1
1970 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1971 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1972 ,l,impf,jmpf,nframe,acfrcv)
1973 if(debugprint)print*,'sample ',varname,' = ',acfrcv(im/2,(jsta+jend)/2)
1974
1975 varname='acfrst'
1976 vcoordname='sfc'
1977 l=1
1978 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1979 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1980 ,l,impf,jmpf,nframe,acfrst)
1981 if(debugprint)print*,'sample ',varname,' = ',acfrst(im/2,(jsta+jend)/2)
1982
1983!insert-mp
1984 varname='ssroff'
1985 vcoordname='sfc'
1986 l=1
1987 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1988 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1989 ,l,impf,jmpf,nframe,ssroff)
1990 if(debugprint)print*,'sample ',varname,' = ',ssroff(im/2,(jsta+jend)/2)
1991
1992! reading UNDERGROUND RUNOFF
1993 varname='bgroff'
1994 vcoordname='sfc'
1995 l=1
1996 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
1997 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
1998 ,l,impf,jmpf,nframe,bgroff)
1999 if(debugprint)print*,'sample ',varname,' = ',bgroff(im/2,(jsta+jend)/2)
2000
2001 varname='rlwin'
2002 vcoordname='sfc'
2003 l=1
2004 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2005 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2006 ,l,impf,jmpf,nframe,rlwin)
2007 if(debugprint)print*,'sample ',varname,' = ',rlwin(im/2,(jsta+jend)/2)
2008
2009 varname='rlwtoa'
2010 vcoordname='sfc'
2011 l=1
2012 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2013 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2014 ,l,impf,jmpf,nframe,rlwtoa)
2015 if(debugprint)print*,'sample ',varname,' = ',rlwtoa(im/2,(jsta+jend)/2)
2016
2017 varname='alwin'
2018 vcoordname='sfc'
2019 l=1
2020 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2021 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2022 ,l,impf,jmpf,nframe,alwin)
2023 if(debugprint)print*,'sample ',varname,' = ',alwin(im/2,(jsta+jend)/2)
2024
2025 varname='alwout'
2026 vcoordname='sfc'
2027 l=1
2028 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2029 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2030 ,l,impf,jmpf,nframe,alwout)
2031 if(debugprint)print*,'sample ',varname,' = ',alwout(im/2,(jsta+jend)/2)
2032
2033 varname='alwtoa'
2034 vcoordname='sfc'
2035 l=1
2036 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2037 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2038 ,l,impf,jmpf,nframe,alwtoa)
2039 if(debugprint)print*,'sample ',varname,' = ',alwtoa(im/2,(jsta+jend)/2)
2040
2041 varname='rswin'
2042 vcoordname='sfc'
2043 l=1
2044 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2045 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2046 ,l,impf,jmpf,nframe,rswin)
2047 if(debugprint)print*,'sample ',varname,' = ',rswin(im/2,(jsta+jend)/2)
2048
2049 varname='rswinc'
2050 vcoordname='sfc'
2051 l=1
2052 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2053 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2054 ,l,impf,jmpf,nframe,rswinc)
2055 if(debugprint)print*,'sample ',varname,' = ',rswinc(im/2,(jsta+jend)/2)
2056
2057 varname='rswout'
2058 vcoordname='sfc'
2059 l=1
2060 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2061 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2062 ,l,impf,jmpf,nframe,rswout)
2063 if(debugprint)print*,'sample ',varname,' = ',rswout(im/2,(jsta+jend)/2)
2064
2065 varname='aswin'
2066 vcoordname='sfc'
2067 l=1
2068 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2069 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2070 ,l,impf,jmpf,nframe,aswin)
2071 if(debugprint)print*,'sample ',varname,' = ',aswin(im/2,(jsta+jend)/2)
2072
2073 varname='aswout'
2074 vcoordname='sfc'
2075 l=1
2076 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2077 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2078 ,l,impf,jmpf,nframe,aswout)
2079 if(debugprint)print*,'sample ',varname,' = ',aswout(im/2,(jsta+jend)/2)
2080
2081 varname='aswtoa'
2082 vcoordname='sfc'
2083 l=1
2084 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2085 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2086 ,l,impf,jmpf,nframe,aswtoa)
2087 if(debugprint)print*,'sample ',varname,' = ',aswtoa(im/2,(jsta+jend)/2)
2088
2089 varname='sfcshx'
2090 vcoordname='sfc'
2091 l=1
2092 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2093 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2094 ,l,impf,jmpf,nframe,sfcshx)
2095 if(debugprint)print*,'sample ',varname,' = ',sfcshx(im/2,(jsta+jend)/2)
2096
2097 varname='sfclhx'
2098 vcoordname='sfc'
2099 l=1
2100 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2101 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2102 ,l,impf,jmpf,nframe,sfclhx)
2103 if(debugprint)print*,'sample ',varname,' = ',sfclhx(im/2,(jsta+jend)/2)
2104
2105 varname='subshx'
2106 vcoordname='sfc'
2107 l=1
2108 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2109 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2110 ,l,impf,jmpf,nframe,subshx)
2111 if(debugprint)print*,'sample ',varname,' = ',subshx(im/2,(jsta+jend)/2)
2112
2113 varname='snopcx'
2114 vcoordname='sfc'
2115 l=1
2116 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2117 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2118 ,l,impf,jmpf,nframe,snopcx)
2119 if(debugprint)print*,'sample ',varname,' = ',snopcx(im/2,(jsta+jend)/2)
2120
2121 varname='sfcuvx'
2122 vcoordname='sfc'
2123 l=1
2124 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2125 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2126 ,l,impf,jmpf,nframe,sfcuvx)
2127 if(debugprint)print*,'sample ',varname,' = ',sfcuvx(im/2,(jsta+jend)/2)
2128
2129 varname='potevp'
2130 vcoordname='sfc'
2131 l=1
2132 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2133 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2134 ,l,impf,jmpf,nframe,potevp)
2135 if(debugprint)print*,'sample ',varname,' = ',potevp(im/2,(jsta+jend)/2)
2136
2137 varname='rlwtt'
2138 vcoordname='mid layer'
2139 do l=1,lm
2140! ll=lm-l+1
2141 ll=l
2142 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2143 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2144 ,l,impf,jmpf,nframe,rlwtt(1,jsta_2l,ll))
2145 if(debugprint)print*,'sample l ',varname,' = ',ll,rlwtt(im/2,(jsta+jend)/2,ll)
2146 end do ! do loop for l
2147
2148 varname='rswtt'
2149 vcoordname='mid layer'
2150 do l=1,lm
2151! ll=lm-l+1
2152 ll=l
2153 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2154 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2155 ,l,impf,jmpf,nframe,rswtt(1,jsta_2l,ll))
2156 if(debugprint)print*,'sample l ',varname,' = ',ll,rswtt(im/2,(jsta+jend)/2,ll)
2157 end do ! do loop for l
2158 where(rlwtt/=spval .and. rswtt/=spval)ttnd=rswtt+rlwtt
2159
2160 varname='tcucn'
2161 vcoordname='mid layer'
2162 do l=1,lm
2163! ll=lm-l+1
2164 ll=l
2165 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2166 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2167 ,l,impf,jmpf,nframe,tcucn(1,jsta_2l,ll))
2168 if(debugprint)print*,'sample l ',varname,' = ',ll,tcucn(im/2,(jsta+jend)/2,ll)
2169 end do ! do loop for l
2170
2171 varname='train'
2172 vcoordname='mid layer'
2173 do l=1,lm
2174! ll=lm-l+1
2175 ll=l
2176 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2177 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2178 ,l,impf,jmpf,nframe,train(1,jsta_2l,ll))
2179 if(debugprint)print*,'sample l ',varname,' = ',ll,train(im/2,(jsta+jend)/2,ll)
2180 end do ! do loop for l
2181
2182 varname='CFRCV'
2183 vcoordname='sfc'
2184 l=1
2185 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2186 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2187 ,l,impf,jmpf,nframe,ncfrcv)
2188 if(debugprint)print*,'sample ',varname,' = ',ncfrcv(im/2,(jsta+jend)/2)
2189
2190 varname='CFRST'
2191 vcoordname='sfc'
2192 l=1
2193 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2194 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2195 ,l,impf,jmpf,nframe,ncfrst)
2196 if(debugprint)print*,'sample ',varname,' = ',ncfrst(im/2,(jsta+jend)/2)
2197
2198! set default to not empty buket
2199 nsrfc=0
2200 nrdlw=0
2201 nrdsw=0
2202 nheat=0
2203 nclod=0
2204 nprec=0
2205 nphs=0
2206
2207 varname='nprec'
2208 if(me==0)then
2209 call nemsio_getheadvar(nfile,trim(varname),nprec,iret)
2210 if (iret /= 0) then
2211 print*,varname," not found in file-Assigned zero"
2212 end if
2213 end if
2214 call mpi_bcast(nprec,1,mpi_integer,0,mpi_comm_comp,iret)
2215 if(debugprint)print*,'sample ',varname,' = ',nprec
2216
2217 varname='nphs'
2218 if(me==0)then
2219 call nemsio_getheadvar(nfile,trim(varname),nphs,iret)
2220 if (iret /= 0) then
2221 print*,varname," not found in file-Assigned zero"
2222 end if
2223 end if
2224 call mpi_bcast(nphs,1,mpi_integer,0,mpi_comm_comp,iret)
2225 if(debugprint)print*,'sample ',varname,' = ',nphs
2226
2227 varname='nclod'
2228 if(me==0)then
2229 call nemsio_getheadvar(nfile,trim(varname),nclod,iret)
2230 if (iret /= 0) then
2231 print*,varname," not found in file-Assigned zero"
2232 end if
2233 end if
2234 call mpi_bcast(nclod,1,mpi_integer,0,mpi_comm_comp,iret)
2235 if(debugprint)print*,'sample ',varname,' = ',nclod
2236
2237 varname='nheat'
2238 if(me==0)then
2239 call nemsio_getheadvar(nfile,trim(varname),nheat,iret)
2240 if (iret /= 0) then
2241 print*,varname," not found in file-Assigned zero"
2242 end if
2243 end if
2244 call mpi_bcast(nheat,1,mpi_integer,0,mpi_comm_comp,iret)
2245 if(debugprint)print*,'sample ',varname,' = ',nheat
2246
2247 varname='nrdlw'
2248 if(me==0)then
2249 call nemsio_getheadvar(nfile,trim(varname),nrdlw,iret)
2250 if (iret /= 0) then
2251 print*,varname," not found in file-Assigned zero"
2252 end if
2253 end if
2254 call mpi_bcast(nrdlw,1,mpi_integer,0,mpi_comm_comp,iret)
2255 if(debugprint)print*,'sample ',varname,' = ',nrdlw
2256
2257 varname='nrdsw'
2258 if(me==0)then
2259 call nemsio_getheadvar(nfile,trim(varname),nrdsw,iret)
2260 if (iret /= 0) then
2261 print*,varname," not found in file-Assigned zero"
2262 end if
2263 end if
2264 call mpi_bcast(nrdsw,1,mpi_integer,0,mpi_comm_comp,iret)
2265 if(debugprint)print*,'sample ',varname,' = ',nrdsw
2266
2267 varname='nsrfc'
2268 if(me==0)then
2269 call nemsio_getheadvar(nfile,trim(varname),nsrfc,iret)
2270 if (iret /= 0) then
2271 print*,varname," not found in file-Assigned zero"
2272 end if
2273 end if
2274 call mpi_bcast(nsrfc,1,mpi_integer,0,mpi_comm_comp,iret)
2275 if(debugprint)print*,'sample ',varname,' = ',nsrfc
2276
2277! VarName='avrain'
2278! if(me==0)then
2279! call nemsio_getheadvar(nfile,trim(varname),avrain,iret)
2280! if (iret /= 0) then
2281! print*,VarName," not found in file-Assigned zero"
2282! end if
2283! end if
2284! call mpi_bcast(avrain,1,MPI_REAL,0,mpi_comm_comp,iret)
2285! if(debugprint)print*,'sample ',VarName,' = ',avrain
2286
2287! VarName='avcnvc'
2288! if(me==0)then
2289! call nemsio_getheadvar(nfile,trim(varname),avcnvc,iret)
2290! if (iret /= 0) then
2291! print*,VarName," not found in file-Assigned zero"
2292! end if
2293! end if
2294! call mpi_bcast(avcnvc,1,MPI_REAL,0,mpi_comm_comp,iret)
2295! if(debugprint)print*,'sample ',VarName,' = ',avcnvc
2296
2297! VarName='ardlw'
2298! if(me==0)then
2299! call nemsio_getheadvar(nfile,trim(varname),ardlw,iret)
2300! if (iret /= 0) then
2301! print*,VarName," not found in file-Assigned zero"
2302! end if
2303! end if
2304! call mpi_bcast(ardlw,1,MPI_REAL,0,mpi_comm_comp,iret)
2305! if(debugprint)print*,'sample ',VarName,' = ',ardlw
2306
2307! VarName='ardsw'
2308! if(me==0)then
2309! call nemsio_getheadvar(nfile,trim(varname),ardsw,iret)
2310! if (iret /= 0) then
2311! print*,VarName," not found in file-Assigned zero"
2312! end if
2313! end if
2314! call mpi_bcast(ardsw,1,MPI_REAL,0,mpi_comm_comp,iret)
2315! if(debugprint)print*,'sample ',VarName,' = ',ardsw
2316
2317! VarName='asrfc'
2318! if(me==0)then
2319! call nemsio_getheadvar(nfile,trim(varname),asrfc,iret)
2320! if (iret /= 0) then
2321! print*,VarName," not found in file-Assigned zero"
2322! end if
2323! end if
2324! call mpi_bcast(asrfc,1,MPI_REAL,0,mpi_comm_comp,iret)
2325! if(debugprint)print*,'sample ',VarName,' = ',asrfc
2326
2327! ! setting all counters to 1 because Ratko has divided all fluxes by counters within NEMS model
2328! avrain=1.0
2329! avcnvc=1.0
2330! ardlw=1.0
2331! ardsw=1.0
2332! asrfc=1.0
2333! if(debugprint)print*,'sample avrain, avcnvc, ardlw, ardsw, asrfc = ', &
2334! avrain, avcnvc, ardlw, ardsw, asrfc
2335
2336!-- Changes to NMMB to allow for counters to vary during the forecast by making them
2337! 2D arrays in order to get around an ESMF limitation (Ferrier 13 Aug 2009)
2338
2339 varname='AVRAIN'
2340 vcoordname='sfc'
2341 l=1
2342 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2343 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2344 ,l,impf,jmpf,nframe,buf)
2345 avrain=buf(im/2,(jsta+jend)/2)
2346 if(debugprint)print*,'sample ',varname,' = ',avrain
2347
2348 varname='AVCNVC'
2349 vcoordname='sfc'
2350 l=1
2351 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2352 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2353 ,l,impf,jmpf,nframe,buf)
2354 avcnvc=buf(im/2,(jsta+jend)/2)
2355 if(debugprint)print*,'sample ',varname,' = ',avcnvc
2356
2357 varname='ARDLW'
2358 vcoordname='sfc'
2359 l=1
2360 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2361 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2362 ,l,impf,jmpf,nframe,buf)
2363 ardlw=buf(im/2,(jsta+jend)/2)
2364 if(debugprint)print*,'sample ',varname,' = ',ardlw
2365
2366 varname='ARDSW'
2367 vcoordname='sfc'
2368 l=1
2369 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2370 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2371 ,l,impf,jmpf,nframe,buf)
2372 ardsw=buf(im/2,(jsta+jend)/2)
2373 if(debugprint)print*,'sample ',varname,' = ',ardsw
2374
2375 varname='ASRFC'
2376 vcoordname='sfc'
2377 l=1
2378 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2379 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2380 ,l,impf,jmpf,nframe,buf)
2381 asrfc=buf(im/2,(jsta+jend)/2)
2382 if(debugprint)print*,'sample ',varname,' = ',asrfc
2383
2384! reading 10 m wind
2385 varname='u10'
2386 vcoordname='10 m above gnd'
2387 l=1
2388 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2389 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2390 ,l,impf,jmpf,nframe,u10h)
2391! Chuang Aug 2012: 10 m winds are computed on mass points in the model
2392! post interpolates them onto V points because copygb interpolates
2393! wind points differently and 10 m winds are identified as 33/34
2394 call h2u(u10h,u10)
2395 if(debugprint)print*,'sample ',varname,' = ',u10(im/2,(jsta+jend)/2)
2396
2397 varname='v10'
2398 vcoordname='10 m above gnd'
2399 l=1
2400 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2401 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2402 ,l,impf,jmpf,nframe,v10h)
2403 call h2u(v10h,v10)
2404 if(debugprint)print*,'sample ',varname,' = ',v10(im/2,(jsta+jend)/2)
2405
2406 varname='u10max'
2407 vcoordname='10 m above gnd'
2408 l=1
2409 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2410 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2411 ,l,impf,jmpf,nframe,u10max)
2412 if(debugprint)print*,'sample ',varname,' = ',u10max(im/2,(jsta+jend)/2)
2413
2414 varname='v10max'
2415 vcoordname='10 m above gnd'
2416 l=1
2417 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2418 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2419 ,l,impf,jmpf,nframe,v10max)
2420 if(debugprint)print*,'sample ',varname,' = ',v10max(im/2,(jsta+jend)/2)
2421
2422! reading SMSTAV
2423 varname='smstav'
2424 vcoordname='sfc'
2425 l=1
2426 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2427 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2428 ,l,impf,jmpf,nframe,smstav)
2429 if(debugprint)print*,'sample ',varname,' = ',smstav(im/2,(jsta+jend)/2)
2430 if(debugprint)print*,'MAX/MIN ',varname,' = ' &
2431 ,maxval(smstav),minval(smstav)
2432
2433 varname='smstot'
2434 vcoordname='sfc'
2435 l=1
2436 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2437 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2438 ,l,impf,jmpf,nframe,smstot)
2439 if(debugprint)print*,'sample ',varname,' = ',smstot(im/2,(jsta+jend)/2)
2440
2441! reading VEGETATION TYPE
2442 varname='VGTYP'
2443 vcoordname='sfc'
2444 l=1
2445 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2446 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2447 ,l,impf,jmpf,nframe,sfcevp) ! temporary use sfcevp because it's real in nemsio
2448
2449! do j=jsta,jend
2450! do i=1,im
2451! if(sfcevp(i,j)> 27.0 .or. sfcevp(i,j)<1.0)print*, &
2452! 'bad vegtype=',i,j,sfcevp(i,j)
2453! end do
2454! end do
2455
2456 where(sfcevp /= spval)ivgtyp=nint(sfcevp)
2457 if(debugprint)print*,'sample ',varname,' = ',ivgtyp(im/2,(jsta+jend)/2)
2458
2459 varname='SLTYP'
2460 vcoordname='sfc'
2461 l=1
2462 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2463 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2464 ,l,impf,jmpf,nframe,sfcevp)
2465 where(sfcevp /= spval)isltyp=nint(sfcevp)
2466 if(debugprint)print*,'sample ',varname,' = ',isltyp(im/2,(jsta+jend)/2)
2467
2468! call retrieve_index(index,VarName,varname_all,nrecs,iret)
2469! if (iret /= 0) then
2470! print*,VarName," not found in file-Assigned missing values"
2471! ISLTYP=NINT(SPVAL)
2472! else
2473! this_offset=file_offset(index+1)+(jsta_2l-1)*4*im
2474! this_length=im*(jend_2u-jsta_2l+1)
2475! call mpi_file_read_at(iunit,this_offset &
2476! ,isltyp,this_length,mpi_integer4, mpi_status_ignore, ierr)
2477! if (ierr /= 0) then
2478! print*,"Error reading ", VarName,"Assigned missing values"
2479! ISLTYP=NINT(SPVAL)
2480! end if
2481! end if
2482
2483 varname='sfcevp'
2484 vcoordname='sfc'
2485 l=1
2486 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2487 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2488 ,l,impf,jmpf,nframe,sfcevp)
2489 if(debugprint)print*,'sample ',varname,' = ',sfcevp(im/2,(jsta+jend)/2)
2490
2491 varname='sfcexc'
2492 vcoordname='sfc'
2493 l=1
2494 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2495 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2496 ,l,impf,jmpf,nframe,sfcexc)
2497 if(debugprint)print*,'sample ',varname,' = ',sfcexc(im/2,(jsta+jend)/2)
2498 if(debugprint)print*,'MAX/MIN ',varname,' = ' &
2499 ,maxval(sfcexc),minval(sfcexc)
2500
2501 varname='acsnow'
2502 vcoordname='sfc'
2503 l=1
2504 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2505 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2506 ,l,impf,jmpf,nframe,acsnow)
2507 if(debugprint)print*,'sample ',varname,' = ',acsnow(im/2,(jsta+jend)/2)
2508
2509 varname='acsnom'
2510 vcoordname='sfc'
2511 l=1
2512 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2513 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2514 ,l,impf,jmpf,nframe,acsnom)
2515 if(debugprint)print*,'sample ',varname,' = ',acsnom(im/2,(jsta+jend)/2)
2516
2517 varname='tsea'
2518 vcoordname='sfc'
2519 l=1
2520 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2521 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2522 ,l,impf,jmpf,nframe,sst)
2523 if(debugprint)print*,'sample ',varname,' = ',sst(im/2,(jsta+jend)/2)
2524
2525! VarName='EL_PBL' ! not in nems io yet
2526 varname='xlen_mix'
2527 vcoordname='mid layer'
2528 do l=1,lm
2529! ll=lm-l+1
2530 ll=l
2531 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2532 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2533 ,l,impf,jmpf,nframe,el_pbl(1,jsta_2l,ll))
2534 if(debugprint)print*,'sample l ',varname,' = ',ll,el_pbl(im/2,(jsta+jend)/2,ll)
2535 end do ! do loop for l
2536
2537 varname='exch_h'
2538 vcoordname='mid layer'
2539 do l=1,lm
2540! ll=lm-l+1
2541 ll=l
2542 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2543 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2544 ,l,impf,jmpf,nframe,exch_h(1,jsta_2l,ll))
2545 if(debugprint)print*,'sample l ',varname,' = ',ll,exch_h(im/2,(jsta+jend)/2,ll)
2546 end do ! do loop for l
2547
2548 varname='thz0'
2549 vcoordname='sfc'
2550 l=1
2551 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2552 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2553 ,l,impf,jmpf,nframe,thz0)
2554 if(debugprint)print*,'sample ',varname,' = ',thz0(im/2,(jsta+jend)/2)
2555
2556 varname='qz0'
2557 vcoordname='sfc'
2558 l=1
2559 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2560 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2561 ,l,impf,jmpf,nframe,qz0)
2562 if(debugprint)print*,'sample ',varname,' = ',qz0(im/2,(jsta+jend)/2)
2563
2564 varname='uz0'
2565 vcoordname='sfc'
2566 l=1
2567 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2568 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2569 ,l,impf,jmpf,nframe,uz0)
2570 if(debugprint)print*,'sample ',varname,' = ',uz0(im/2,(jsta+jend)/2)
2571
2572 varname='vz0'
2573 vcoordname='sfc'
2574 l=1
2575 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2576 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2577 ,l,impf,jmpf,nframe,vz0)
2578 if(debugprint)print*,'sample ',varname,' = ',vz0(im/2,(jsta+jend)/2)
2579
2580!
2581! Very confusing story ...
2582!
2583! Retrieve htop and hbot => They are named CNVTOP, CNVBOT in the model and
2584! with HBOTS,HTOPS (shallow conv) and HBOTD,HTOPD (deep conv) represent
2585! the 3 sets of convective cloud base/top arrays tied to the frequency
2586! that history files are written.
2587!
2588! IN THE *MODEL*, arrays HBOT,HTOP are similar to CNVTOP,CNVBOT but are
2589! used in radiation and are tied to the frequency of radiation updates.
2590!
2591! For historical reasons model arrays CNVTOP,CNVBOT are renamed HBOT,HTOP
2592! and manipulated throughout the post.
2593
2594! retrieve htop and hbot
2595! VarName='HTOP'
2596 varname='cnvtop'
2597 vcoordname='sfc'
2598 l=1
2599 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2600 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2601 ,l,impf,jmpf,nframe,htop)
2602 where(htop /= spval)htop=float(lm)-htop+1.0
2603! where(htop /= spval .and. htop > lm)htop=lm*1.0
2604 if(debugprint)print*,'sample ',varname,' = ',htop(im/2,(jsta+jend)/2)
2605
2606! VarName='HBOT'
2607 varname='cnvbot'
2608 vcoordname='sfc'
2609 l=1
2610 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2611 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2612 ,l,impf,jmpf,nframe,hbot)
2613 where(hbot /= spval)hbot=float(lm)-hbot+1.0
2614! where(hbot /= spval .and. hbot > lm)hbot=lm*1.0
2615 if(debugprint)print*,'sample ',varname,' = ',hbot(im/2,(jsta+jend)/2)
2616
2617 varname='htopd'
2618 vcoordname='sfc'
2619 l=1
2620 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2621 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2622 ,l,impf,jmpf,nframe,htopd)
2623 where(htopd /= spval)htopd=float(lm)-htopd+1.0
2624! where(htopd /= spval .and. htopd > lm)htopd=lm*1.0
2625 if(debugprint)print*,'sample ',varname,' = ',htopd(im/2,(jsta+jend)/2)
2626
2627 varname='hbotd'
2628 vcoordname='sfc'
2629 l=1
2630 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2631 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2632 ,l,impf,jmpf,nframe,hbotd)
2633 where(hbotd /= spval)hbotd=float(lm)-hbotd+1.0
2634! where(hbotd /= spval .and. hbotd > lm)hbotd=lm*1.0
2635 if(debugprint)print*,'sample ',varname,' = ',hbotd(im/2,(jsta+jend)/2)
2636
2637 varname='htops'
2638 vcoordname='sfc'
2639 l=1
2640 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2641 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2642 ,l,impf,jmpf,nframe,htops)
2643 where(htops /= spval)htops=float(lm)-htops+1.0
2644! where(htops /= spval .and. htops > lm)htops=lm*1.0
2645 if(debugprint)print*,'sample ',varname,' = ',htops(im/2,(jsta+jend)/2)
2646
2647 varname='hbots'
2648 vcoordname='sfc'
2649 l=1
2650 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2651 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2652 ,l,impf,jmpf,nframe,hbots)
2653 where(hbots /= spval)hbots=float(lm)-hbots+1.0
2654! where(hbots /= spval .and. hbots > lm)hbots=lm*1.0
2655 if(debugprint)print*,'sample ',varname,' = ',hbots(im/2,(jsta+jend)/2)
2656
2657 varname='cuppt'
2658 vcoordname='sfc'
2659 l=1
2660 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2661 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2662 ,l,impf,jmpf,nframe,cuppt)
2663 if(debugprint)print*,'sample ',varname,' = ',cuppt(im/2,(jsta+jend)/2)
2664
2665 varname='cprate'
2666 vcoordname='sfc'
2667 l=1
2668 call getnemsandscatter(me,nfile,im,jm,jsta,jsta_2l &
2669 ,jend_2u,mpi_comm_comp,icnt,idsp,spval,varname,vcoordname &
2670 ,l,impf,jmpf,nframe,cprate)
2671 if(debugprint)print*,'sample ',varname,' = ',cprate(im/2,(jsta+jend)/2)
2672
2673!!!! DONE GETTING
2674
2675!$omp parallel do private(i,j,l)
2676 do l = 1, lm
2677 do j = jsta, jend
2678 do i = 1, im
2679 IF(abs(t(i,j,l))>1.0e-3 .and. (wh(i,j,1) < spval)) &
2680 omga(i,j,l) = -wh(i,j,l)*pmid(i,j,l)*g &
2681 / (rd*t(i,j,l)*(1.+d608*q(i,j,l)))
2682 end do
2683 end do
2684 end do
2685
2686 thl=210.
2687 plq=70000.
2688
2689 CALL table(ptbl,ttbl,pt, &
2690 rdq,rdth,rdp,rdthe,pl,thl,qs0,sqs,sthe,the0)
2691
2692 CALL tableq(ttblq,rdpq,rdtheq,plq,thl,stheq,the0q)
2693
2694!
2695!
2696 IF(me==0)THEN
2697 WRITE(6,*)' SPL (POSTED PRESSURE LEVELS) BELOW: '
2698 WRITE(6,51) (spl(l),l=1,lsm)
2699 50 FORMAT(14(f4.1,1x))
2700 51 FORMAT(8(f8.1,1x))
2701 ENDIF
2702!
2703! COMPUTE DERIVED TIME STEPPING CONSTANTS.
2704!
2705!MEB need to get DT
2706! DT = 120. !MEB need to get DT
2707! NPHS = 4 !MEB need to get physics DT
2708 dtq2 = dt * nphs !MEB need to get physics DT
2709 tsph = 3600./dt !MEB need to get DT
2710
2711 IF (pthresh > 0.) THEN
2712 pthresh = 0.01*dtq2/3.6e6 !-- Precip rate >= 0.01 mm/h
2713! PTHRESH = 0.01*DTQ2/(3600.*39.37) !-- Precip rate >= 0.01 inches/h
2714 ENDIF
2715
2716 tsrfc=float(nsrfc)/tsph
2717 if(me==0)write(6,*)'tsfrc ',tsrfc,nsrfc,tsph
2718 IF(nsrfc==0)tsrfc=float(ifhr) !in case buket does not get emptied
2719 trdlw=float(nrdlw)/tsph
2720 IF(nrdlw==0)trdlw=float(ifhr) !in case buket does not get emptied
2721 trdsw=float(nrdsw)/tsph
2722 IF(nrdsw==0)trdsw=float(ifhr) !in case buket does not get emptied
2723 theat=float(nheat)/tsph
2724 IF(nheat==0)theat=float(ifhr) !in case buket does not get emptied
2725 tclod=float(nclod)/tsph
2726 IF(nclod==0)tclod=float(ifhr) !in case buket does not get emptied
2727 tprec=float(nprec)/tsph
2728 IF(nprec==0)tprec=float(ifhr) !in case buket does not get emptied
2729! TPREC=float(ifhr)
2730 if(me==0)print*,'TSRFC TRDLW TRDSW THEAT TCLOD TPREC= ' &
2731 ,tsrfc, trdlw, trdsw, theat, tclod, tprec
2732!MEB need to get DT
2733
2734!how am i going to get this information?
2735! NPREC = INT(TPREC *TSPH+D50)
2736! NHEAT = INT(THEAT *TSPH+D50)
2737! NCLOD = INT(TCLOD *TSPH+D50)
2738! NRDSW = INT(TRDSW *TSPH+D50)
2739! NRDLW = INT(TRDLW *TSPH+D50)
2740! NSRFC = INT(TSRFC *TSPH+D50)
2741!how am i going to get this information?
2742!
2743! IF(ME==0)THEN
2744! WRITE(6,*)' '
2745! WRITE(6,*)'DERIVED TIME STEPPING CONSTANTS'
2746! WRITE(6,*)' NPREC,NHEAT,NSRFC : ',NPREC,NHEAT,NSRFC
2747! WRITE(6,*)' NCLOD,NRDSW,NRDLW : ',NCLOD,NRDSW,NRDLW
2748! ENDIF
2749!
2750! COMPUTE DERIVED MAP OUTPUT CONSTANTS.
2751 DO l = 1,lsm
2752 alsl(l) = alog(spl(l))
2753 END DO
2754!
2755!HC WRITE IGDS OUT FOR WEIGHTMAKER TO READ IN AS KGDSIN
2756 if(me==0)then
2757 print*,'writing out igds'
2758 igdout=110
2759! open(igdout,file='griddef.out',form='unformatted'
2760! + ,status='unknown')
2761 IF(maptype==203)THEN !A STAGGERED E-GRID
2762 WRITE(igdout)203
2763 WRITE(igdout)im
2764 WRITE(igdout)jm
2765 WRITE(igdout)latstart
2766 WRITE(igdout)lonstart
2767 WRITE(igdout)136
2768 WRITE(igdout)cenlat
2769 WRITE(igdout)cenlon
2770 WRITE(igdout)dxval
2771 WRITE(igdout)dyval
2772 WRITE(igdout)64
2773 WRITE(igdout)0
2774 WRITE(igdout)0
2775 WRITE(igdout)0
2776 WRITE(igdout)latlast
2777 WRITE(igdout)lonlast
2778 ELSE IF(maptype==205)THEN !A STAGGERED B-GRID
2779 WRITE(igdout)205
2780 WRITE(igdout)im
2781 WRITE(igdout)jm
2782 WRITE(igdout)latstart
2783 WRITE(igdout)lonstart
2784 WRITE(igdout)136
2785 WRITE(igdout)cenlat
2786 WRITE(igdout)cenlon
2787 WRITE(igdout)dxval
2788 WRITE(igdout)dyval
2789 WRITE(igdout)64
2790 WRITE(igdout)latlast
2791 WRITE(igdout)lonlast
2792 WRITE(igdout)0
2793 WRITE(igdout)0
2794 WRITE(igdout)0
2795 END IF
2796 open(111,file='copygb_gridnav.txt',form='formatted' &
2797 ,status='unknown')
2798 IF(maptype==203)THEN !A STAGGERED E-GRID
2799 write(111,1000) 2*im-1,jm,latstart,lonstart,cenlon, &
2800 nint(dxval*107.),nint(dyval*110.),cenlat,cenlat
2801 ELSE IF(maptype==205)THEN !A STAGGERED B-GRID
2802 if(grib=="grib2") then
2803 write(111,1000) im,jm,latstart/1000,lonstart/1000,cenlon/1000, &
2804 nint(dxval*107.)/1000,nint(dyval*110.)/1000, &
2805 cenlat/1000,cenlat/1000, &
2806 latlast/1000,lonlast/1000
2807 endif
2808 END IF
28091000 format('255 3 ',2(i4,x),i6,x,i7,x,'8 ',i7,x,2(i6,x),'0 64', &
2810 3(x,i6),x,i7)
2811 close(111)
2812!
2813 IF (maptype==205)THEN !A STAGGERED B-GRID
2814 open(112,file='latlons_corners.txt',form='formatted' &
2815 ,status='unknown')
2816 if(grib=="grib2") then
2817 write(112,1001)latstart/1000,(lonstart/1000)-360000, &
2818 latse/1000, &
2819 lonse/1000,latnw/1000,lonnw/1000,latlast/1000, &
2820 (lonlast/1000)-360000
2821 endif
28221001 format(4(i6,x,i7,x))
2823 close(112)
2824 ENDIF
2825
2826 end if
2827
2828! close all files
2829 call nemsio_close(nfile,iret=status)
2830!
2831 if(me==0)write(*,*)'end of INIT_NEMS'
2832
2833 RETURN
2834 END
initpost_nems
subroutine initpost_nems(nrec, nfile)
This routine initializes constants and variables at the start of an NEMS model or post processor run.
Definition INITPOST_NEMS.f:21