UPP (develop)
Loading...
Searching...
No Matches
WRFPOST.F
Go to the documentation of this file.
1
43!---------------------------------------------------------------------
45!---------------------------------------------------------------------
46 PROGRAM wrfpost
47
48!
49!
50!============================================================================================================
51!
52! This is an MPI code. All array indexing is with respect to the global indices. Loop indices
53! look as follows for N MPI tasks.
54!
55!
56!
57! Original New
58! Index Index
59!
60! JM ----------------------------------------------- JEND
61! JM-1 - - JEND_M
62! JM-2 - MPI TASK N-1 - JEND_M2
63! - -
64! - -
65! ----------------------------------------------- JSTA, JSTA_M, JSTA_M2
66! ----------------------------------------------- JEND, JEND_M, JEND_M2
67! - -
68! - MPI TASK N-2 -
69! - -
70! - -
71! ----------------------------------------------- JSTA, JSTA_M, JSTA_M2
72!
73! .
74! .
75! .
76!
77! ----------------------------------------------- JEND, JEND_M, JEND_M2
78! - -
79! - MPI TASK 1 -
80! - -
81! - -
82! ----------------------------------------------- JSTA, JSTA_M, JSTA_M2
83! ----------------------------------------------- JEND, JEND_M, JEND_M2
84! - -
85! - MPI TASK 0 -
86! 3 - - JSTA_M2
87! 2 - - JSTA_M
88! 1 ----------------------------------------------- JSTA
89!
90! 1 IM
91!
92!
93! Jim Tuccillo
94! Jan 2000
95!
96! README - Jim Tuccillo Feb 2001
97!
98! Many common blocks have been replaced by modules to support Fortran
99! "allocate" commands. Many of the 3-D arrays are now allocated to be the
100! exact size required based on the number of MPI tasks. The dimensioning will be
101! x ( im,jsta_2l:jend_2u,lm)
102! Most 2-D arrays continue to be dimensioned (im,jm). This is fine but please be aware
103! that the EXCH routine for arrays dimensioned (im,jm) is different than arrays dimensioned
104! (im,jsta_2l:jend_2u). Also, be careful about passing any arrays dimensioned
105! (im,jst_2l:jend_2u,lm). See examples in the code as to the correct calling sequence and
106! EXCH routine to use.
107!
108!
109! ASYNCHRONOUS I/O HAS BEEN ADDED. THE LAST MPI TASK DOES THE I/O. IF THERE IS
110! ONLY ONE MPI TASK THN TASK ) DOES THE I/O.
111! THE CODE HAS GOTTEN A LITTLE KLUDGY. BASICLY, IM, IMX and IMOUT MUST BE EQUAL
112! AND REPRESENT THE VALUE USED IN THE MODEL. THE SAME HOLDS FOR JM, JMX and JMOUT.
113!
114! Jim Tuccillo June 2001
115!
116!
117!===========================================================================================
118!
119 use netcdf
120#if defined(BUILD_WITH_NEMSIO)
121 use nemsio_module, only: nemsio_getheadvar, nemsio_gfile, nemsio_init, nemsio_open, &
122 nemsio_getfilehead,nemsio_close
123#endif
124 use ctlblk_mod, only: filenameaer, me, num_procs, num_servers, mpi_comm_comp, datestr, &
125 mpi_comm_inter, filename, ioform, grib, idat, filenameflux, filenamed3d, gdsdegr, &
126 spldef, modelname, ihrst, lsmdef,vtimeunits, tprec, pthresh, datahandle, im, jm, lm, &
127 lp1, lm1, im_jm, isf_surface_physics, nsoil, spl, lsmp1, global, imp_physics, &
128 ista, iend, ista_m, iend_m, ista_2l, iend_2u, &
129 jsta, jend, jsta_m, jend_m, jsta_2l, jend_2u, novegtype, icount_calmict, npset, datapd,&
132 fixed_tim, time_output, imin, surfce2_tim, komax, ivegsrc, d3d_on, gocart_on,rdaod, &
133 readxml_tim, spval, fullmodelname, submodelname, hyb_sigp, filenameflat, aqf_on,numx, &
134 run_ifi_tim, slrutah_on, d2d_chem, gtg_on, method_blsn
135 use grib2_module, only: gribit2,num_pset,nrecout,first_grbtbl,grib_info_finalize
136 use upp_ifi_mod, only: write_ifi_debug_files
137!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
138 implicit none
139!
140#if defined(BUILD_WITH_NEMSIO)
141 type(nemsio_gfile) :: nfile,ffile,rfile
142#endif
143 include "mpif.h"
144!
145! DECLARE VARIABLES.
146!
147! SET HEADER WRITER FLAGS TO TRUE.
148!
149!temporary vars
150!
151 real(kind=8) :: time_initpost=0.,initpost_tim=0.,btim,bbtim
152 real rinc(5), untcnvt
153 integer :: status=0,iostatusd3d=0,iostatusflux=0
154 integer i,j,iii,l,k,ierr,nrec,ist,lusig,idrt,ncid3d,ncid2d,varid
155 integer :: prntsec,iim,jjm,llm,ioutcount,itmp,iret,iunit, &
156 iunitd3d,iyear,imn,iday,lcntrl,ieof
157 integer :: iostatusaer
158 logical :: popascal
159!
160 integer :: kpo,kth,kpv
161 real,dimension(komax) :: po,th,pv
162 namelist/nampgb/kpo,po,kth,th,kpv,pv,filenameaer,d3d_on,gocart_on,gccpp_on, nasa_on,gtg_on,method_blsn,popascal &
163 ,hyb_sigp,rdaod,d2d_chem, aqf_on,slrutah_on, vtimeunits,numx,write_ifi_debug_files
164 integer :: itag_ierr
165 namelist/model_inputs/filename,ioform,grib,datestr,modelname,submodelname &
166 ,filenameflux,filenameflat
167
168 character startdate*19,sysdepinfo*80,iowrfname*3,post_fname*255
169 character cgar*1,cdum*4,line*10
170!
171!------------------------------------------------------------------------------
172! START HERE
173!
174 call start()
175!
176! INITIALIZE MPI
177
178 CALL setup_servers(me, &
179 & num_procs, &
180 & num_servers, &
181 & mpi_comm_comp, &
182 & mpi_comm_inter)
183!
184! ME IS THE RANK
185! NUM_PROCS IS THE NUMBER OF TASKS DOING POSTING
186! NUM_SERVERS IS ONE IF THERE ARE MORE THAN ONE TOTAL MPI TASKS, OTHERWISE ZERO
187! MPI_COMM_COMP IS THE INTRACOMMUNICATOR
188! MPI_COMM_INTER IS THE INTERCOMMUNICATOR FOR COMMUNCATION BETWEEN TASK 0 OF THE
189! TASKS DOING THE POSTING AND THE I/O SERVER
190!
191!
192! IF WE HAVE MORE THAN 1 MPI TASK THEN WE WILL FIRE UP THE IO SERVER
193! THE LAST TASK ( IN THE CONTEXT OF MPI_COMM_WORLD ) IS THE I/O SERVER
194!
195 if (me == 0) CALL w3tagb('nems ',0000,0000,0000,'np23 ')
196
197 if ( me >= num_procs ) then
198!
199 call server
200!
201 else
202 spval = 9.9e10
203!
204!**************************************************************************
205!KaYee: Read itag in Fortran Namelist format
206!Set default
207 submodelname='NONE'
208!Set control file name
209 filenameflat='postxconfig-NT.txt'
210 numx=1
211!open namelist
212 open(5,file='itag')
213 read(5,nml=model_inputs,iostat=itag_ierr,err=888)
214888 if (itag_ierr /= 0) then
215 print*,'Incorrect namelist variable(s) found in the itag file,stopping.'
216 stop
217 endif
218 if (me == 0) write(6, model_inputs)
219
220! if(MODELNAME == 'NMM')then
221! read(5,1114) VTIMEUNITS
222! 1114 format(a4)
223! if (me==0) print*,'VALID TIME UNITS = ', VTIMEUNITS
224! endif
225!
226 303 format('MODELNAME="',a,'" SUBMODELNAME="',a,'"')
227
228 if(me==0) write(*,*)'MODELNAME: ', modelname, submodelname
229
230! assume for now that the first date in the stdin file is the start date
231 read(datestr,300) iyear,imn,iday,ihrst,imin
232 if (me==0) write(*,*) 'in WRFPOST iyear,imn,iday,ihrst,imin', &
233 iyear,imn,iday,ihrst,imin
234 300 format(i4,1x,i2,1x,i2,1x,i2,1x,i2)
235
236 idat(1) = imn
237 idat(2) = iday
238 idat(3) = iyear
239 idat(4) = ihrst
240 idat(5) = imin
241
242 111 format(a256)
243 112 format(a19)
244 113 format(a20)
245 114 format(a8)
246 120 format(a5)
247 121 format(a4)
248
249!KaYee: Read in GFS/FV3 runs in Fortran Namelist Format.
250 if(grib=='grib2') then
251 gdsdegr = 1.d6
252 endif
253!
254! set default for kpo, kth, th, kpv, pv
255 kpo = 0
256 po = 0
257 kth = 6
258 th = (/310.,320.,350.,450.,550.,650.,(0.,k=kth+1,komax)/) ! isentropic level to output
259 kpv = 8
260 pv = (/0.5,-0.5,1.0,-1.0,1.5,-1.5,2.0,-2.0,(0.,k=kpv+1,komax)/)
261
262 hyb_sigp = .true.
263 d3d_on = .false.
264 gocart_on = .false.
265 gccpp_on = .false.
266 nasa_on = .false.
267 aqf_on = .false.
268 slrutah_on = .false.
269 gtg_on = .false.
270 method_blsn = .true.
271 popascal = .false.
272 filenameaer = ''
273 rdaod = .false.
274 d2d_chem = .false.
275 vtimeunits = ''
276
277 read(5,nampgb,iostat=iret,end=119)
278 119 continue
279 if (me == 0) write(6, nampgb)
280 if(mod(num_procs,numx)/=0) then
281 if (me==0) then
282 print*,'total proces, num_procs=', num_procs
283 print*,'number of subdomain in x direction, numx=', numx
284 print*,'remainder of num_procs/numx = ', mod(num_procs,numx)
285 print*,'Warning!!! the remainder of num_procs/numx is not 0, reset numx=1 &
286 & in this run or you adjust numx in the itag file to restart'
287 endif
288! stop 9999
289 numx=1
290 if(me == 0) print*,'Warning!!! Reset numx as 1, numx=',numx
291 endif
292 if(numx>num_procs/2) then
293 if (me==0) then
294 print*,'total proces, num_procs=', num_procs
295 print*,'number of subdomain in x direction, numx=', numx
296 print*,'Warning!!! numx cannot exceed num_procs/2, reset numx=1 in this run'
297 print*,'or you adjust numx in the itag file to restart'
298 endif
299 numx=1
300 if(me == 0) print*,'Warning!!! Reset numx as 1, numx=',numx
301 endif
302 if(me == 0) then
303 print*,'komax,iret for nampgb= ',komax,iret
304 print*,'komax,kpo,kth,th,kpv,pv,fileNameAER,nasa_on,popascal= ',komax,kpo &
305 & ,kth,th(1:kth),kpv,pv(1:kpv),trim(filenameaer),nasa_on,popascal
306 print*,'NUM_PROCS=',num_procs
307 print*,'numx= ',numx
308 endif
309
310 IF(trim(ioform) /= 'netcdfpara' .AND. trim(ioform) /= 'netcdf' ) THEN
311 numx=1
312 if(me == 0) print*,'2D decomposition only supports netcdfpara IO.'
313 if(me == 0) print*,'Reset numx= ',numx
314 ENDIF
315
316 IF(modelname /= 'FV3R' .AND. modelname /= 'GFS') THEN
317 numx=1
318 if(me == 0) print*,'2D decomposition only supports GFS and FV3R.'
319 if(me == 0) print*,'Reset numx= ',numx
320 ENDIF
321
322! set up pressure level from POSTGPVARS or DEFAULT
323 if(kpo == 0) then
324! use default pressure levels
325 if(me == 0) then
326 print*,'using default pressure levels,spldef=',(spldef(l),l=1,lsmdef)
327 endif
328 lsm = lsmdef
329 do l=1,lsm
330 spl(l) = spldef(l)
331 end do
332 else
333! use POSTGPVARS
334 if(me == 0) then
335 print*,'using pressure levels from POSTGPVARS'
336 endif
337 lsm = kpo
338 if( .not. popascal ) then
339 untcnvt = 100.
340 else
341 untcnvt = 1.
342 endif
343 if(po(lsm) < po(1))then ! post logic assumes asscending
344 do l=1,lsm
345 spl(l) = po(lsm-l+1)*untcnvt
346 end do
347 else
348 do l=1,lsm
349 spl(l) = po(l)*untcnvt
350 end do
351 end if
352 end if
353 lsmp1 = lsm+1
354
355 116 continue
356
357! set PTHRESH for different models
358 if(modelname == 'NMM')then
359 pthresh = 0.000004
360 else
361 pthresh = 0.000001
362 end if
363!Chuang: add dynamical allocation
364 if(trim(ioform) == 'netcdf' .OR. trim(ioform) == 'netcdfpara') THEN
365 IF(modelname == 'NCAR' .OR. modelname == 'RAPR' .OR. modelname == 'NMM') THEN
366 call ext_ncd_ioinit(sysdepinfo,status)
367 call ext_ncd_open_for_read( trim(filename), 0, 0, " ", &
368 datahandle, status)
369 if ( status /= 0 ) then
370 print*,'error opening ',filename, ' Status = ', status ; stop
371 endif
372 call ext_ncd_get_dom_ti_integer(datahandle &
373 ,'WEST-EAST_GRID_DIMENSION',iim,1,ioutcount, status )
374 im = iim-1
375 call ext_ncd_get_dom_ti_integer(datahandle &
376 ,'SOUTH-NORTH_GRID_DIMENSION',jjm,1,ioutcount, status )
377 jm = jjm-1
378 call ext_ncd_get_dom_ti_integer(datahandle &
379 ,'BOTTOM-TOP_GRID_DIMENSION',llm,1,ioutcount, status )
380 lm = llm-1
381 lp1 = lm+1
382 lm1 = lm-1
383 im_jm = im*jm
384
385! Read and set global value for surface physics scheme
386 call ext_ncd_get_dom_ti_integer(datahandle &
387 ,'SF_SURFACE_PHYSICS',itmp,1,ioutcount, status )
388
389 ! J. Kenyon / 23 Aug 2024: GSL MPAS output erroneously
390 ! indicates SF_SURFACE_PHYSICS = 0 (should be 3), so overwrite
391 ! here:
392 if (modelname == 'RAPR' .and. submodelname == 'MPAS') itmp = 3
393
394 isf_surface_physics = itmp
395! set NSOIL to 4 as default for NOAH but change if using other
396! SFC scheme
397 nsoil = 4
398 IF(itmp == 1) then !thermal diffusion scheme
399 nsoil = 5
400 ELSE IF(itmp == 3) then ! RUC LSM
401 nsoil = 9
402 ELSE IF(itmp == 7) then ! Pleim Xu
403 nsoil = 2
404 END IF
405
406 call ext_ncd_ioclose ( datahandle, status )
407 ELSE
408! use parallel netcdf lib directly to read FV3 output in netCDF
409 spval = 9.99e20
410 status = nf90_open(trim(filename),ior(nf90_nowrite,nf90_mpiio), &
411 ncid3d,comm=mpi_comm_world,info=mpi_info_null)
412 if ( status /= 0 ) then
413 print*,'error opening ',filename, ' Status = ', status
414 stop
415 endif
416 status = nf90_open(trim(filenameflux),ior(nf90_nowrite,nf90_mpiio), &
417 ncid2d,comm=mpi_comm_world,info=mpi_info_null)
418 if ( status /= 0 ) then
419 print*,'error opening ',filenameflux, ' Status = ', status
420 stop
421 endif
422! read in LSM index and nsoil here
423 status=nf90_get_att(ncid2d,nf90_global,'landsfcmdl', isf_surface_physics)
424 if(status/=0)then
425 print*,'landsfcmdl not found; assigning to 2'
426 isf_surface_physics=2 !set LSM physics to 2 for NOAH
427 endif
428 if(isf_surface_physics<2)then
429 isf_surface_physics=2 !set LSM physics to 2 for NOAH
430 endif
431 status=nf90_get_att(ncid2d,nf90_global,'nsoil', nsoil)
432 if(status/=0)then
433 print*,'nsoil not found; assigning to 4'
434 nsoil=4 !set nsoil to 4 for NOAH
435 endif
436! read imp_physics
437 status=nf90_get_att(ncid2d,nf90_global,'imp_physics',imp_physics)
438 if(status/=0)then
439 print*,'imp_physics not found; assigning to GFDL 11'
440 imp_physics=11
441 endif
442! get dimesions
443 status = nf90_inq_dimid(ncid3d,'grid_xt',varid)
444 if ( status /= 0 ) then
445 print*,status,varid
446 stop 1
447 end if
448 status = nf90_inquire_dimension(ncid3d,varid,len=im)
449 if ( status /= 0 ) then
450 print*,status
451 stop 1
452 end if
453 status = nf90_inq_dimid(ncid3d,'grid_yt',varid)
454 if ( status /= 0 ) then
455 print*,status,varid
456 stop 1
457 end if
458 status = nf90_inquire_dimension(ncid3d,varid,len=jm)
459 if ( status /= 0 ) then
460 print*,status
461 stop 1
462 end if
463 status = nf90_inq_dimid(ncid3d,'pfull',varid)
464 if ( status /= 0 ) then
465 print*,status,varid
466 stop 1
467 end if
468 status = nf90_inquire_dimension(ncid3d,varid,len=lm)
469 if ( status /= 0 ) then
470 print*,status
471 stop 1
472 end if
473 lp1 = lm+1
474 lm1 = lm-1
475 im_jm = im*jm
476! set NSOIL to 4 as default for NOAH but change if using other
477! SFC scheme
478! NSOIL = 4
479 END IF
480
481 ELSE IF(trim(ioform) == 'binary' .OR. &
482 trim(ioform) == 'binarympiio' ) THEN
483 print*,'WRF Binary format is no longer supported'
484 stop 9996
485! NEMSIO format
486#if defined(BUILD_WITH_NEMSIO)
487 ELSE IF(trim(ioform) == 'binarynemsio' .or. &
488 trim(ioform) == 'binarynemsiompiio' )THEN
489
490 spval = 9.99e20
491 IF(me == 0)THEN
492 call nemsio_init(iret=status)
493 call nemsio_open(nfile,trim(filename),'read',iret=status)
494 if ( status /= 0 ) then
495 print*,'error opening ',filename, ' Status = ', status ; stop
496 endif
497!---
498 call nemsio_getfilehead(nfile,iret=status,nrec=nrec &
499 ,dimx=im,dimy=jm,dimz=lm,nsoil=nsoil)
500 if ( status /= 0 ) then
501 print*,'error finding model dimensions '; stop
502 endif
503 call nemsio_getheadvar(nfile,'global',global,iret)
504 if (iret /= 0)then
505 print*,"global not found in file-Assigned false"
506 global = .false.
507 end if
508 IF(modelname == 'GFS') global = .true.
509! global NMMB has i=1 overlap with i=im so post will cut off i=im
510 if(global .and. modelname == 'NMM') im = im-1
511
512 END IF
513
514 CALL mpi_bcast(im, 1,mpi_integer,0, mpi_comm_comp,status)
515 call mpi_bcast(jm, 1,mpi_integer,0, mpi_comm_comp,status)
516 call mpi_bcast(lm, 1,mpi_integer,0, mpi_comm_comp,status)
517 call mpi_bcast(nsoil,1,mpi_integer,0, mpi_comm_comp,status)
518
519 call mpi_bcast(global,1,mpi_logical,0,mpi_comm_comp,status)
520 lp1 = lm+1
521 lm1 = lm-1
522 im_jm = im*jm
523
524! opening GFS flux file
525 IF(modelname == 'GFS') THEN
526! iunit=33
527 call nemsio_open(ffile,trim(filenameflux),'read',iret=iostatusflux)
528 if ( iostatusflux /= 0 ) then
529 print*,'error opening ',filenameflux, ' Status = ', iostatusflux
530 endif
531 iostatusd3d = -1
532 iunitd3d = -1
533!
534! opening GFS aer file
535 call nemsio_open(rfile,trim(filenameaer),'read',iret=iostatusaer)
536 if ( iostatusaer /= 0 .and. me == 0) then
537 print*,'error opening AER ',filenameaer, ' Status = ', iostatusaer
538 endif
539!
540! print*,'iostatusD3D in WRFPOST= ',iostatusD3D
541
542 END IF
543#endif
544 ELSE
545 print*,'UNKNOWN MODEL OUTPUT FORMAT, STOPPING'
546 stop 9999
547 END IF
548
549
550 CALL mpi_first()
551 CALL allocate_all()
552
553!
554! INITIALIZE POST COMMON BLOCKS
555!
556 lcntrl = 14
557 rewind(lcntrl)
558
559! EXP. initialize netcdf here instead
560 bbtim = mpi_wtime()
561 btim = mpi_wtime()
562! set default novegtype
563 if(modelname == 'GFS')THEN
564 novegtype = 13
565 ivegsrc = 2
566 else if(modelname=='NMM' .and. trim(ioform)=='binarynemsio')then
567 novegtype = 20
568 ivegsrc = 1
569 else if(modelname=='RAPR')then
570 novegtype = 20
571 ivegsrc = 1
572 else ! USGS
573 novegtype = 24
574 ivegsrc = 0
575 end if
576
577! Reading model output for different models and IO format
578
579 IF(trim(ioform) == 'netcdf' .OR. trim(ioform) == 'netcdfpara') THEN
580 IF(modelname == 'NCAR' .OR. modelname == 'RAPR') THEN
581 IF (submodelname == 'MPAS') THEN
582 CALL initpost_mpas
583 ELSE
584 CALL initpost
585 ENDIF
586 ELSE IF (modelname == 'FV3R' .OR. modelname == 'GFS') THEN
587! use parallel netcdf library to read output directly
588 CALL initpost_netcdf(ncid2d,ncid3d)
589 ELSE
590 print*,'POST does not have netcdf option for model,',modelname,' STOPPING,'
591 stop 9998
592 END IF
593 ELSE IF(trim(ioform) == 'binarympiio') THEN
594 IF(modelname == 'NCAR' .OR. modelname == 'RAPR' .OR. modelname == 'NMM') THEN
595 print*,'WRF BINARY IO FORMAT IS NO LONGER SUPPORTED, STOPPING'
596 stop 9996
597 ELSE IF(modelname == 'RSM') THEN
598 print*,'MPI BINARY IO IS NOT YET INSTALLED FOR RSM, STOPPING'
599 stop 9997
600 ELSE
601 print*,'POST does not have mpiio option for this model, STOPPING'
602 stop 9998
603 END IF
604#if defined(BUILD_WITH_NEMSIO)
605 ELSE IF(trim(ioform) == 'binarynemsio') THEN
606 IF(modelname == 'NMM') THEN
607 CALL initpost_nems(nrec,nfile)
608 ELSE
609 print*,'POST does not have nemsio option for model,',modelname,' STOPPING,'
610 stop 9998
611
612 END IF
613
614 ELSE IF(trim(ioform) == 'binarynemsiompiio')THEN
615 IF(modelname == 'GFS') THEN
616! close nemsio file for serial read
617 call nemsio_close(nfile,iret=status)
618 call nemsio_close(ffile,iret=status)
619 call nemsio_close(rfile,iret=status)
620 CALL initpost_gfs_nems_mpiio(iostatusaer)
621 ELSE
622 print*,'POST does not have nemsio mpi option for model,',modelname, &
623 'STOPPING,'
624 stop 9999
625
626 END IF
627#endif
628 ELSE
629 print*,'UNKNOWN MODEL OUTPUT FORMAT, STOPPING'
630 stop 9999
631 END IF
632 initpost_tim = initpost_tim +(mpi_wtime() - btim)
633 IF(me == 0)THEN
634 WRITE(6,*)'WRFPOST: INITIALIZED POST COMMON BLOCKS'
635 ENDIF
636!
637! IF GRIB2 read out post aviable fields xml file and post control file
638!
639 if(grib == "grib2") then
640 btim=mpi_wtime()
641 call read_xml()
642 readxml_tim = readxml_tim + (mpi_wtime() - btim)
643 endif
644!
645! LOOP OVER THE OUTPUT GRID(S). FIELD(S) AND OUTPUT GRID(S) ARE SPECIFIED
646! IN THE CONTROL FILE. WE PROCESS ONE GRID AND ITS FIELDS AT A TIME.
647! THAT'S WHAT THIS LOOP DOES.
648!
649 icount_calmict = 0
650 first_grbtbl = .true.
651 npset = 0
652!10 CONTINUE
653!
654! READ CONTROL FILE DIRECTING WHICH FIELDS ON WHICH
655! LEVELS AND TO WHICH GRID TO INTERPOLATE DATA TO.
656! VARIABLE IEOF/=0 WHEN THERE ARE NO MORE GRIDS TO PROCESS.
657!
658! -------- grib1 processing ---------------
659! ------------------
660! if (grib == "grib1") then !DO NOT REVERT TO GRIB1. GRIB1 NOT SUPPORTED ANYMORE
661! IEOF = 0
662! do while (ieof == 0)
663! CALL READCNTRL(kth,IEOF)
664! IF(ME == 0)THEN
665! WRITE(6,*)'POST: RETURN FROM READCNTRL. ', 'IEOF=',IEOF
666! ENDIF
667!
668! PROCESS SELECTED FIELDS. FOR EACH SELECTED FIELD/LEVEL
669! WE GO THROUGH THE FOLLOWING STEPS:
670! (1) COMPUTE FIELD IF NEED BE
671! (2) WRITE FIELD TO OUTPUT FILE IN GRIB.
672!
673! if (ieof == 0) then
674! CALL PROCESS(kth,kpv,th(1:kth),pv(1:kpv),iostatusD3D)
675! IF(ME == 0)THEN
676! WRITE(6,*)' '
677! WRITE(6,*)'WRFPOST: PREPARE TO PROCESS NEXT GRID'
678! ENDIF
679! endif
680!
681! PROCESS NEXT GRID.
682!
683! enddo
684! -------- grib2 processing ---------------
685! ------------------
686! elseif (grib == "grib2") then
687 if (me==0) write(*,*) ' in WRFPOST OUTFORM= ',grib
688 if (me==0) write(*,*) ' GRIB1 IS NOT SUPPORTED ANYMORE'
689 if (grib == "grib2") then
690 do while (npset < num_pset)
691 npset = npset+1
692 if (me==0) write(*,*)' in WRFPOST npset=',npset,' num_pset=',num_pset
693 CALL set_outflds(kth,th,kpv,pv)
694 if (me==0) write(*,*)' in WRFPOST size datapd',size(datapd)
695 if(allocated(datapd)) deallocate(datapd)
696!Jesse x-decomposition
697! allocate(datapd(im,1:jend-jsta+1,nrecout+100))
698 allocate(datapd(1:iend-ista+1,1:jend-jsta+1,nrecout+100))
699!$omp parallel do private(i,j,k)
700 do k=1,nrecout+100
701 do j=1,jend+1-jsta
702!Jesse x-decomposition
703! do i=1,im
704 do i =1,iend+1-ista
705 datapd(i,j,k) = 0.
706 enddo
707 enddo
708 enddo
709 call get_postfilename(post_fname)
710 if (me==0) write(*,*)'post_fname=',trim(post_fname)
711 if (me==0) write(*,*)'get_postfilename,post_fname=',trim(post_fname), &
712 'npset=',npset, 'num_pset=',num_pset, &
713 'iSF_SURFACE_PHYSICS=',isf_surface_physics
714!
715! PROCESS SELECTED FIELDS. FOR EACH SELECTED FIELD/LEVEL
716! WE GO THROUGH THE FOLLOWING STEPS:
717! (1) COMPUTE FIELD IF NEED BE
718! (2) WRITE FIELD TO OUTPUT FILE IN GRIB.
719!
720 CALL process(kth,kpv,th(1:kth),pv(1:kpv),iostatusd3d)
721 IF(me == 0) WRITE(6,*)'WRFPOST: PREPARE TO PROCESS NEXT GRID'
722!
723! write(*,*)'enter gribit2 before mpi_barrier'
724 call mpi_barrier(mpi_comm_comp,ierr)
725
726! if(me==0)call w3tage('bf grb2 ')
727! write(*,*)'enter gribit2 after mpi barrier'
728 call gribit2(post_fname)
729 deallocate(datapd)
730 deallocate(fld_info)
731!
732! PROCESS NEXT GRID.
733!
734 enddo
735
736 endif
737!
738!-------
739 call grib_info_finalize()
740!
741 IF(me == 0) THEN
742 WRITE(6,*)' '
743 WRITE(6,*)'ALL GRIDS PROCESSED.'
744 WRITE(6,*)' '
745 ENDIF
746!
747 call de_allocate
748
749! GO TO 98
750 1000 CONTINUE
751!exp call ext_ncd_ioclose ( DataHandle, Status )
752!
753 IF(me == 0) THEN
754 print*, 'INITPOST_tim = ', initpost_tim
755 print*, 'MDLFLD_tim = ', etafld2_tim
756 print*, 'MDL2P_tim = ',eta2p_tim
757 print*, 'MDL2SIGMA_tim = ',mdl2sigma_tim
758 print*, 'MDL2AGL_tim = ',mdl2agl_tim
759 print*, 'SURFCE_tim = ',surfce2_tim
760 print*, 'CLDRAD_tim = ',cldrad_tim
761 print*, 'MISCLN_tim = ',miscln_tim
762 print*, 'MDL2STD_tim = ',mdl2std_tim
763 print*, 'FIXED_tim = ',fixed_tim
764 print*, 'MDL2THANDPV_tim = ',mdl2thandpv_tim
765 print*, 'CALRAD_WCLOUD_tim = ',calrad_wcloud_tim
766 print*, 'RUN_IFI_tim = ',run_ifi_tim
767 print*, 'Total time = ',(mpi_wtime() - bbtim)
768 print*, 'Time for OUTPUT = ',time_output
769 print*, 'Time for READxml = ',readxml_tim
770 endif
771!
772! END OF PROGRAM.
773!
774!
775! MPI_LAST WILL SHUTDOWN THE IO SERVER, IF IT EXISTS
776!
777 CALL mpi_last
778!
779!
780 end if
781!
782!
783!
784! call summary()
785 if (me == 0) CALL w3tage('UNIFIED_POST')
786 CALL mpi_finalize(ierr)
787
788
789 stop 0
790
791 END
792
subroutine allocate_all()
SET UP MESSGAE PASSING INFO.
subroutine de_allocate
2023-08-16 | Yali Mao | Add CIT to GTG fields.
Definition DEALLOCATE.f:19
subroutine initpost
This routine initializes constants and variables at the start of an ETA model or post processor run.
Definition INITPOST.F:26
subroutine initpost_gfs_nems_mpiio(iostatusaer)
initializes constants and variables at the start of GFS model or post processor run.
subroutine initpost_mpas
This routine initializes constants and variables when using UPP on MPAS model output.
subroutine initpost_nems(nrec, nfile)
INITPOST_NEMS This routine initializes constants and variables at the start of an NEMS model or post ...
subroutine initpost_netcdf(ncid2d, ncid3d)
2023-04-17 | Eric James | Read in unified ext550 extinction (and remove aodtot) for RRFS 2023-04-21 |...
subroutine mpi_first()
MPI_FIRST() Sets up message passing info (MPI).
Definition MPI_FIRST.f:46
subroutine mpi_last
SUBPROGRAM: MPI_LAST SHUTS DOWN THE IO SERVER PRGRMMR: TUCCILLO ORG: IBM.
Definition MPI_LAST.f:32
subroutine process(kth, kpv, th, pv, iostatusd3d)
process() is a driver for major post routines.
Definition PROCESS.f:40
subroutine read_xml()
Definition READ_xml.f:14
subroutine server
SUBPROGRAM: SERVER PERFORMS IO TO DISK PRGRMMR: TUCCILLO ORG: IBM.
Definition SERVER.f:37
subroutine setup_servers(mype, npes, inumq, mpi_comm_comp, mpi_comm_inter)
This subroutine is to setup I/O servers.
subroutine set_outflds(kth, th, kpv, pv)
Reads post XML control file.
Definition SET_OUTFLDS.f:24
program wrfpost
Definition WRFPOST.F:46
real(kind=8) etafld2_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209
real(kind=8) mdl2agl_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209
real(kind=8) surfce2_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209
real(kind=8) eta2p_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209
real(kind=8) calrad_wcloud_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209
real(kind=8) mdl2std_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209
real(kind=8) run_ifi_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209
real(kind=8) time_output
Initialized as 0, but never used.
Definition CTLBLK.f:219
real(kind=8) fixed_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209
real(kind=8) readxml_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209
real(kind=8) miscln_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209
real(kind=8) mdl2thandpv_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209
real(kind=8) cldrad_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209
real(kind=8) mdl2sigma_tim
Time to execute named routine; note that ETAFLD2 and ETA2P refer to MDLFLD and MDL2P routines respect...
Definition CTLBLK.f:209