subroutine read_files(mype) !$$$ subprogram documentation block ! . . . . ! subprogram: read_files get info about atm & sfc guess files ! prgmmr: derber org: np23 date: 2002-11-14 ! ! abstract: This routine determines how many global atmospheric and ! surface guess files are present. The valid time for each ! guess file is determine. The time are then sorted in ! ascending order. This information is broadcast to all ! mpi tasks. ! ! program history log: ! 2002-11-14 derber ! 2004-06-16 treadon - update documentation ! 2004-08-02 treadon - add only to module use, add intent in/out ! 2004-12-02 treadon - replace mpe_ibcast (IBM extension) with ! standard mpi_bcast ! 2005-01-27 treadon - make use of sfcio module ! 2005-02-18 todling - no need to read entire sfc file; only head needed ! 2005-03-30 treadon - clean up formatting of write statements ! 2006-01-09 treadon - use sigio to read gfs spectral coefficient file header ! 2007-03-01 tremolet - measure time from beginning of assimilation window ! 2007-04-17 todling - getting nhr_assimilation from gsi_4dvar ! 2008-05-27 safford - rm unused vars ! 2009-01-07 todling - considerable revamp (no pre-assigned dims) ! 2009-08-26 li - add variables to handle nst guess files ! 2010-04-20 jing - set hrdifsig_all and hrdifsfc_all for non-ESMF cases. ! 2010-12-06 Huang - make use of nemsio_module to check whether atm and sfc files ! are in NEMSIO format and get header informaion 'lpl' ! 2010-12-16 li - (1) set nfldnst, ntguesnst, ifilenst, hrdifnst, hrdifnst_all for nst files. ! (2) add zero initialization of nfldsfc ! 2011-04-04 Huang - change looping over 0,99 to find existed sigf and sfcf ! twice. Use fcst_hr_sig and fcst_hr_sfc to store info of ! files found in first loop of 0,99. Use nfldsig and nfldsfc ! to access needed sigf and sfcf w/ fcst_hr_sig and *_sfc. ! 2015-02-23 Rancic/Thomas - add l4densvar to time window logical ! 2017-09-08 li - add sfcnst_comb to get nfldnst and control when sfc & nst combined ! 2019-03-21 Wei/Martin - add capability to read in aerosol guess from NEMS ! 2019-09-24 martin - add support for use_gfs_ncio ! 2020-07-08 Wei - fix the capability to count external aerosol files ! ! input argument list: ! mype - mpi task id ! ! output argument list: ! ! comments: ! The difference of time Info between operational GFS IO (gfshead%, sfc_head%), ! analysis time (iadate), and NEMSIO (idate=) ! ! gfshead & sfc_head NEMSIO Header Analysis time (obsmod) ! =================== ============================ ========================== ! %idate(1) Hour idate(1) Year iadate(1) Year ! %idate(2) Month idate(2) Month iadate(2) Month ! %idate(3) Day idate(3) Day iadate(3) Day ! %idate(4) Year idate(4) Hour iadate(4) Hour ! idate(5) Minute iadate(5) Minute ! idate(6) Scaled seconds ! idate(7) Seconds multiplier ! ! The difference of header forecasting hour Info bewteen operational GFS IO ! (gfshead%, sfc_head%) and NEMSIO ! ! gfshead & sfc_head NEMSIO Header ! ========================== ============================ ! %fhour FCST Hour (r_kind) nfhour FCST Hour (i_kind) ! nfminute FCST Mins (i_kind) ! nfsecondn FCST Secs (i_kind) numerator ! nfsecondd FCST Secs (i_kind) denominator ! ! %fhour = float(nfhour) + float(nfminute)/r60 + float(nfsecondn)/float(nfsecondd)/r3600 ! attributes: ! language: f90 ! machine: ibm RS/6000 SP ! !$$$ use kinds, only: r_kind,r_single,i_kind,i_llong use mpimod, only: mpi_rtype,mpi_comm_world,ierror,npe,mpi_itype use guess_grids, only: nfldsig,nfldsfc,nfldnst,ntguessig,ntguessfc,ntguesnst,& ifilesig,ifilesfc,ifilenst,hrdifsig,hrdifsfc,hrdifnst,create_gesfinfo use guess_grids, only: hrdifsig_all,hrdifsfc_all,hrdifnst_all use guess_grids, only: nfldaer, ntguesaer, ifileaer, hrdifaer, hrdifaer_all !for aerosol use gsi_4dvar, only: l4dvar,l4densvar,iwinbgn,winlen,nhr_assimilation,nhr_obsbin use hybrid_ensemble_parameters, only: ntlevs_ens use gridmod, only: nlat_sfc,nlon_sfc,lpl_gfs,dx_gfs,use_gfs_nemsio,sfcnst_comb,use_gfs_ncio use constants, only: zero,r60inv,r60,r3600,i_missing use obsmod, only: iadate use gsi_nstcouplermod, only: nst_gsi use sfcio_module, only: sfcio_head,sfcio_sropen,& sfcio_sclose,sfcio_srhead use nstio_module, only: nstio_head,nstio_sropen,& nstio_srclose,nstio_srhead use sigio_module, only: sigio_head,sigio_sropen,& sigio_sclose,sigio_srhead use nemsio_module, only: nemsio_init,nemsio_open,nemsio_close use nemsio_module, only: nemsio_gfile,nemsio_getfilehead,nemsio_getheadvar use read_obsmod, only: gsi_inquire use gsi_io, only: verbose use module_ncio, only: Dataset, Dimension, open_dataset, get_dim, & read_vardata, get_idate_from_time_units, & close_dataset use chemmod, only: lread_ext_aerosol implicit none ! Declare passed variables integer(i_kind),intent(in ) :: mype ! Declare local parameters integer(i_kind),parameter:: lunsfc=11 integer(i_kind),parameter:: lunatm=12 integer(i_kind),parameter:: lunnst=13 integer(i_kind),parameter:: num_lpl=2000 integer(i_kind),parameter:: max_file = 100 real(r_kind),parameter:: r0_001=0.001_r_kind ! Declare local variables logical(4) fexist logical:: present character(6) filename integer(i_kind) i,j,iwan,npem1,iret integer(i_kind) nhr_half,ihr integer(i_kind) iamana(4) ! changed to 4 from 3 for aer files integer(i_kind) nminanl,nmings,nming2,ndiff integer(i_kind),dimension(4):: idateg integer(i_kind),dimension(2):: i_ges integer(i_kind),allocatable,dimension(:):: nst_ges integer(i_kind),dimension(5):: idate5 integer(i_kind),dimension(6):: idate6 integer(i_kind),dimension(num_lpl):: lpl_dum integer(i_kind),dimension(7):: idate integer(i_kind) :: nfhour, nfminute, nfsecondn, nfsecondd integer(i_kind),dimension(:,:),allocatable:: irec integer(i_llong) :: lenbytes real(r_single) hourg4 real(r_kind), allocatable, dimension(:) :: fhour real(r_kind) hourg,t4dv real(r_kind),allocatable,dimension(:,:):: time_atm real(r_kind),allocatable,dimension(:,:):: time_sfc real(r_kind),allocatable,dimension(:,:):: time_nst real(r_kind),allocatable,dimension(:,:):: time_aer type(sfcio_head):: sfc_head type(sigio_head):: sigatm_head type(nstio_head):: nst_head type(nemsio_gfile) :: gfile_atm,gfile_sfc,gfile_nst,gfile_aer logical :: print_verbose type(Dataset) :: atmges, sfcges, nstges type(Dimension) :: ncdim print_verbose=.false. if(verbose)print_verbose=.true. !----------------------------------------------------------------------------- ! Initialize variables nhr_half=nhr_assimilation/2 if(nhr_half*2 < nhr_assimilation) nhr_half=nhr_half+1 npem1=npe-1 fexist=.true. nfldsig=0 nfldsfc=0 nfldnst=0 nfldaer=0 iamana=0 ! Check for non-zero length atm, sfc, aer, and nst files on single task if(mype==npem1)then allocate( irec(max_file,4) ) irec=i_missing ! Check for atm files with non-zero length do i=0,max_file-1 write(filename,'(''sigf'',i2.2)')i call gsi_inquire(lenbytes,fexist,filename,mype) if(fexist .and. lenbytes>0) then nfldsig=nfldsig+1 irec(nfldsig,1) = i end if enddo if(nfldsig==0) then write(6,*)'READ_FILES: ***ERROR*** NO atm fields; aborting' call stop2(169) end if ! Check for sfc files with non-zero length do i=0,max_file-1 write(filename,'(''sfcf'',i2.2)')i call gsi_inquire(lenbytes,fexist,filename,mype) if(fexist .and. lenbytes>0) then nfldsfc=nfldsfc+1 irec(nfldsfc,2) = i end if enddo if(nfldsfc==0) then write(6,*)'READ_FILES: ***ERROR* NO sfc fields; aborting' call stop2(170) end if allocate(time_atm(nfldsig,2),time_sfc(nfldsfc,2)) if(nst_gsi > 0 ) then ! nst application is an option ! Check for nsf files with non-zero length if ( sfcnst_comb ) then nfldnst = nfldsfc else do i=0,max_file-1 write(filename,'(''nstf'',i2.2)')i call gsi_inquire(lenbytes,fexist,filename,mype) if(fexist .and. lenbytes>0) then nfldnst=nfldnst+1 irec(nfldnst,3) = i end if enddo if(nfldnst==0) then write(6,*)'READ_FILES: ***ERROR*** NO nst fields; aborting' call stop2(170) end if end if allocate(time_nst(nfldnst,2)) end if if(lread_ext_aerosol) then ! Check for aer files with non-zero length do i=0,max_file-1 write(filename,'(''aerf'',i2.2)')i call gsi_inquire(lenbytes,fexist,filename,mype) if(fexist .and. lenbytes>0) then nfldaer=nfldaer+1 irec(nfldaer,4) = i end if enddo if(nfldaer==0) then write(6,*)'READ_FILES: ***ERROR*** NO aer fields; aborting' call stop2(170) end if allocate(time_aer(nfldaer,2)) end if ! Let a single task query the guess files. ! Convert analysis time to minutes relative to fixed date call w3fs21(iadate,nminanl) write(6,*)'READ_FILES: analysis date,minutes ',iadate,nminanl ! Check for consistency of times from atmospheric guess files. iwan=0 do i=1,nfldsig write(filename,'(''sigf'',i2.2)')irec(i,1) if(print_verbose)write(6,*)'READ_FILES: process ',trim(filename) if ( (.not. use_gfs_nemsio) .and. (.not. use_gfs_ncio) ) then call sigio_sropen(lunatm,filename,iret) call sigio_srhead(lunatm,sigatm_head,iret) hourg4=sigatm_head%fhour idateg=sigatm_head%idate call sigio_sclose(lunatm,iret) else if (use_gfs_ncio) then atmges = open_dataset(filename,errcode=iret) if (iret /=0 .and. mype==0) & write(6,*)'READ_FILES: ***WARNING*** problem reading atm file ',trim(filename),iret idate6 = get_idate_from_time_units(atmges) call read_vardata(atmges, 'time', fhour) hourg4 = float(nint(fhour(1))) ! going to make this nearest integer for now idateg(1) = idate6(4) idateg(2) = idate6(2) idateg(3) = idate6(3) idateg(4) = idate6(1) call close_dataset(atmges) else call nemsio_init(iret=iret) call nemsio_open(gfile_atm,filename,'READ',iret=iret) idate = i_missing nfhour = i_missing; nfminute = i_missing nfsecondn = i_missing; nfsecondd = i_missing call nemsio_getfilehead(gfile_atm, nfhour=nfhour, nfminute=nfminute, & nfsecondn=nfsecondn, nfsecondd=nfsecondd, idate=idate, iret=iret) call nemsio_close(gfile_atm,iret=iret) if ( nfhour == i_missing .or. nfminute == i_missing .or. & nfsecondn == i_missing .or. nfsecondd == i_missing ) then write(6,*)'READ_FILES: ***ERROR*** some forecast hour info ', & 'are not defined in ', trim(filename) write(6,*)'READ_FILES: nfhour, nfminute, nfsecondn, and nfsecondd = ', & nfhour, nfminute, nfsecondn, nfsecondd call stop2(80) endif hourg4 = float(nfhour) + float(nfminute)/r60 + float(nfsecondn)/float(nfsecondd)/r3600 idateg(1) = idate(4) !hour idateg(2) = idate(2) !month idateg(3) = idate(3) !day idateg(4) = idate(1) !year endif hourg = hourg4 idate5(1)=idateg(4); idate5(2)=idateg(2) idate5(3)=idateg(3); idate5(4)=idateg(1); idate5(5)=0 call w3fs21(idate5,nmings) nming2=nmings+60*hourg write(6,*)'READ_FILES: atm guess file ',filename,hourg,idateg,nming2 t4dv=real((nming2-iwinbgn),r_kind)*r60inv if (l4dvar.or.l4densvar) then if (t4dvwinlen) cycle else ndiff=nming2-nminanl if(abs(ndiff) > 60*nhr_half ) cycle endif iwan=iwan+1 if(nminanl==nming2) iamana(1)=iwan time_atm(iwan,1) = t4dv time_atm(iwan,2) = irec(i,1)+r0_001 end do ! Check for consistency of times from surface guess files. iwan=0 do i=1,nfldsfc write(filename,'(''sfcf'',i2.2)')irec(i,2) if(print_verbose)write(6,*)'READ_FILES: process ',trim(filename) if ( (.not. use_gfs_nemsio) .and. (.not. use_gfs_ncio) ) then call sfcio_sropen(lunsfc,filename,iret) call sfcio_srhead(lunsfc,sfc_head,iret) hourg4=sfc_head%fhour idateg=sfc_head%idate i_ges(1)=sfc_head%lonb i_ges(2)=sfc_head%latb+2 if(sfc_head%latb/2>num_lpl)then write(6,*)'READ_FILES: increase dimension of variable lpl_dum' call stop2(80) endif lpl_dum=0 lpl_dum(1:sfc_head%latb/2)=sfc_head%lpl call sfcio_sclose(lunsfc,iret) if(i == 1 .and. print_verbose)write(6,*)' READ_FILES: in sfcio sfc_head%lpl = ', sfc_head%lpl else if (use_gfs_ncio) then sfcges = open_dataset(filename,errcode=iret) if (iret /=0 .and. mype==0) & write(6,*)'READ_FILES: ***WARNING*** problem reading sfc file ',trim(filename),iret ncdim = get_dim(sfcges, 'grid_xt'); sfc_head%lonb = ncdim%len ncdim = get_dim(sfcges, 'grid_yt'); sfc_head%latb = ncdim%len idate6 = get_idate_from_time_units(sfcges) call read_vardata(sfcges, 'time', fhour) hourg4 = float(nint(fhour(1))) ! going to make this nearest integer for now idateg(1) = idate6(4) idateg(2) = idate6(2) idateg(3) = idate6(3) idateg(4) = idate6(1) i_ges(1)=sfc_head%lonb i_ges(2)=sfc_head%latb+2 if((sfc_head%latb+1)/2>num_lpl)then write(6,*)'READ_FILES: increase dimension of variable lpl_dum' call stop2(80) endif if ( (sfc_head%latb+1)/2 /= sfc_head%latb/2 ) then write(6,*) 'READ_FILES: ****WARNING**** (sfc_head%latb+1)/2 = ', & (sfc_head%latb+1)/2, 'sfc_head%latb/2 = ', sfc_head%latb/2 end if if (allocated(sfc_head%lpl)) deallocate(sfc_head%lpl) allocate(sfc_head%lpl((sfc_head%latb+1)/2)) sfc_head%lpl=sfc_head%lonb call close_dataset(sfcges) lpl_dum=0 lpl_dum(1:sfc_head%latb/2)=sfc_head%lpl deallocate(sfc_head%lpl) else call nemsio_init(iret=iret) call nemsio_open(gfile_sfc,filename,'READ',iret=iret) idate = i_missing nfhour = i_missing; nfminute = i_missing nfsecondn = i_missing; nfsecondd = i_missing call nemsio_getfilehead(gfile_sfc, nfhour=nfhour, nfminute=nfminute, & nfsecondn=nfsecondn, nfsecondd=nfsecondd, idate=idate, & dimx=sfc_head%lonb, dimy=sfc_head%latb, iret=iret) if ( nfhour == i_missing .or. nfminute == i_missing .or. & nfsecondn == i_missing .or. nfsecondd == i_missing ) then write(6,*)'READ_FILES: ***ERROR*** some forecast hour info ', & 'are not defined in ', trim(filename) write(6,*)'READ_FILES: nfhour, nfminute, nfsecondn, and nfsecondd = ', & nfhour, nfminute, nfsecondn, nfsecondd call stop2(80) endif hourg4 = float(nfhour) + float(nfminute)/r60 + float(nfsecondn)/float(nfsecondd)/r3600 idateg(1) = idate(4) !hour idateg(2) = idate(2) !month idateg(3) = idate(3) !day idateg(4) = idate(1) !year i_ges(1)=sfc_head%lonb i_ges(2)=sfc_head%latb+2 if((sfc_head%latb+1)/2>num_lpl)then write(6,*)'READ_FILES: increase dimension of variable lpl_dum' call stop2(80) endif if ( (sfc_head%latb+1)/2 /= sfc_head%latb/2 ) then write(6,*) 'READ_FILES: ****WARNING**** (sfc_head%latb+1)/2 = ', & (sfc_head%latb+1)/2, 'sfc_head%latb/2 = ', sfc_head%latb/2 end if if (allocated(sfc_head%lpl)) deallocate(sfc_head%lpl) allocate(sfc_head%lpl((sfc_head%latb+1)/2)) sfc_head%lpl=i_missing call nemsio_getheadvar(gfile_sfc,'lpl',sfc_head%lpl,iret=iret) if ( iret /= 0 ) then write(6,*)' READ_FILES: ****ERROR**** reading sfc_head%lpl from ', & trim(filename),', iret = ', iret write(6,*)' sfc_head%lpl = ', sfc_head%lpl(1:20) call nemsio_close(gfile_sfc,iret=iret) call stop2(80) end if call nemsio_close(gfile_sfc,iret=iret) lpl_dum=0 lpl_dum(1:sfc_head%latb/2)=sfc_head%lpl deallocate(sfc_head%lpl) endif hourg = hourg4 idate5(1)=idateg(4); idate5(2)=idateg(2) idate5(3)=idateg(3); idate5(4)=idateg(1); idate5(5)=0 call w3fs21(idate5,nmings) nming2=nmings+60*hourg write(6,*)'READ_FILES: sfc guess file ',filename,hourg,idateg,nming2 t4dv=real((nming2-iwinbgn),r_kind)*r60inv if (l4dvar.or.l4densvar) then if (t4dvwinlen) cycle else ndiff=nming2-nminanl if(abs(ndiff) > 60*nhr_half ) cycle endif iwan=iwan+1 if(nminanl==nming2) iamana(2)=iwan time_sfc(iwan,1) = t4dv time_sfc(iwan,2) = irec(i,2)+r0_001 end do ! Check for consistency of times from nst guess files. if ( nst_gsi > 0 ) then allocate(nst_ges(2)) if ( sfcnst_comb ) then time_nst = time_sfc iamana(3) = iamana(2) else iwan=0 do i=1,nfldnst write(filename,'(''nstf'',i2.2)')irec(i,3) write(6,*)'READ_FILES: process ',trim(filename) if ( (.not. use_gfs_nemsio) .and. (.not. use_gfs_ncio) ) then call nstio_sropen(lunnst,filename,iret) call nstio_srhead(lunnst,nst_head,iret) hourg4=nst_head%fhour idateg=nst_head%idate nst_ges(1)=nst_head%lonb nst_ges(2)=nst_head%latb+2 else if (use_gfs_ncio) then nstges = open_dataset(filename) ncdim = get_dim(nstges, 'grid_xt'); nst_head%lonb = ncdim%len ncdim = get_dim(nstges, 'grid_yt'); nst_head%latb = ncdim%len idate6 = get_idate_from_time_units(nstges) call read_vardata(nstges, 'time', fhour) hourg4 = fhour(1) idateg(1) = idate6(4) idateg(2) = idate6(2) idateg(3) = idate6(3) idateg(4) = idate6(1) else call nemsio_init(iret=iret) call nemsio_open(gfile_nst,filename,'READ',iret=iret) idate = i_missing nfhour = i_missing; nfminute = i_missing nfsecondn = i_missing; nfsecondd = i_missing call nemsio_getfilehead(gfile_nst, nfhour=nfhour, nfminute=nfminute, & nfsecondn=nfsecondn, nfsecondd=nfsecondd, idate=idate, & dimx=nst_head%lonb, dimy=nst_head%latb, iret=iret) call nemsio_close(gfile_nst,iret=iret) if ( nfhour == i_missing .or. nfminute == i_missing .or. & nfsecondn == i_missing .or. nfsecondd == i_missing ) then write(6,*)'READ_FILES: ***ERROR*** some forecast hour info ', & 'are not defined in ', trim(filename) write(6,*)'READ_FILES: nfhour, nfminute, nfsecondn, and nfsecondd = ', & nfhour, nfminute, nfsecondn, nfsecondd call stop2(80) endif hourg4 = float(nfhour) + float(nfminute)/r60 + float(nfsecondn)/float(nfsecondd)/r3600 idateg(1) = idate(4) !hour idateg(2) = idate(2) !month idateg(3) = idate(3) !day idateg(4) = idate(1) !year nst_ges(1)=nst_head%lonb nst_ges(2)=nst_head%latb+2 endif if ( i_ges(1) /= nst_ges(1) .or. i_ges(2) /= nst_ges(2) ) then write(6,'(''READ_FILES: sfc file lat,lon '',2i5,'' do not match with nst file lat,lon '',2i5)') & i_ges(2)-2,i_ges(1),nst_ges(2)-2,nst_ges(1) call stop2(80) end if hourg = hourg4 idate5(1)=idateg(4); idate5(2)=idateg(2) idate5(3)=idateg(3); idate5(4)=idateg(1); idate5(5)=0 call w3fs21(idate5,nmings) nming2=nmings+60*hourg write(6,*)'READ_FILES: nst guess file',filename,hourg,idateg,nming2 t4dv=real((nming2-iwinbgn),r_kind)*r60inv if (l4dvar.or.l4densvar) then if (t4dvwinlen) cycle else ndiff=nming2-nminanl if(abs(ndiff) > 60*nhr_half ) cycle endif iwan=iwan+1 if(nminanl==nming2) iamana(3)=iwan time_nst(iwan,1) = t4dv time_nst(iwan,2) = irec(i,3)+r0_001 end do deallocate(nst_ges) endif ! if ( sfcnst_comb ) then endif ! if ( nst_gsi > 0 ) then ! for external aerosol files only ! Check for consistency of times from aer guess files. if ( lread_ext_aerosol ) then write(6,*) 'READ_FILES: nfldaer ', nfldaer iwan=0 do i=1,nfldaer write(filename,'(''aerf'',i2.2)')irec(i,4) write(6,*)'READ_FILES: process ',trim(filename) if ( .not. use_gfs_nemsio ) then write(6,*)'READ_FILES: ***ERROR*** aerosol files only work with nemsio' else call nemsio_init(iret=iret) call nemsio_open(gfile_aer,filename,'READ',iret=iret) idate = i_missing nfhour = i_missing; nfminute = i_missing nfsecondn = i_missing; nfsecondd = i_missing call nemsio_getfilehead(gfile_aer, nfhour=nfhour, nfminute=nfminute, & nfsecondn=nfsecondn, nfsecondd=nfsecondd, idate=idate ) call nemsio_close(gfile_aer,iret=iret) if ( nfhour == i_missing .or. nfminute == i_missing .or. & nfsecondn == i_missing .or. nfsecondd == i_missing ) then write(6,*)'READ_FILES: ***ERROR*** some forecast hour info ', & 'are not defined in ', trim(filename) write(6,*)'READ_FILES: nfhour, nfminute, nfsecondn, and nfsecondd = ', & nfhour, nfminute, nfsecondn, nfsecondd call stop2(80) endif hourg4 = float(nfhour) + float(nfminute)/r60 + float(nfsecondn)/float(nfsecondd)/r3600 idateg(1) = idate(4) !hour idateg(2) = idate(2) !month idateg(3) = idate(3) !day idateg(4) = idate(1) !year endif hourg = hourg4 idate5(1)=idateg(4); idate5(2)=idateg(2) idate5(3)=idateg(3); idate5(4)=idateg(1); idate5(5)=0 call w3fs21(idate5,nmings) nming2=nmings+60*hourg write(6,*)'READ_FILES: aer guess file',filename,hourg,idateg,nming2 t4dv=real((nming2-iwinbgn),r_kind)*r60inv if (l4dvar.or.l4densvar) then if (t4dvwinlen) cycle else ndiff=nming2-nminanl if(abs(ndiff) > 60*nhr_half ) cycle endif iwan=iwan+1 if(nminanl==nming2) iamana(4)=iwan time_aer(iwan,1) = t4dv time_aer(iwan,2) = irec(i,4)+r0_001 end do endif ! if ( lread_ext_aerosol ) then deallocate( irec ) end if ! Broadcast guess file information to all tasks call mpi_bcast(nfldsig,1,mpi_itype,npem1,mpi_comm_world,ierror) call mpi_bcast(nfldsfc,1,mpi_itype,npem1,mpi_comm_world,ierror) print_verbose=.false. if(verbose)print_verbose=.true. if (nst_gsi > 0) call mpi_bcast(nfldnst,1,mpi_itype,npem1,mpi_comm_world,ierror) if (lread_ext_aerosol) call mpi_bcast(nfldaer,1,mpi_itype,npem1,mpi_comm_world,ierror)! for external aerosol files if(.not.allocated(time_atm)) allocate(time_atm(nfldsig,2)) if(.not.allocated(time_sfc)) allocate(time_sfc(nfldsfc,2)) call mpi_bcast(time_atm,2*nfldsig,mpi_rtype,npem1,mpi_comm_world,ierror) call mpi_bcast(time_sfc,2*nfldsfc,mpi_rtype,npem1,mpi_comm_world,ierror) if(.not.allocated(time_nst)) allocate(time_nst(nfldnst,2)) if (nst_gsi > 0 ) call mpi_bcast(time_nst,2*nfldnst,mpi_rtype,npem1,mpi_comm_world,ierror) ! for external aerosol files if(.not.allocated(time_aer)) allocate(time_aer(nfldaer,2)) if (lread_ext_aerosol) call mpi_bcast(time_aer,2*nfldaer,mpi_rtype,npem1,mpi_comm_world,ierror) call mpi_bcast(iamana,3,mpi_rtype,npem1,mpi_comm_world,ierror) call mpi_bcast(i_ges,2,mpi_itype,npem1,mpi_comm_world,ierror) nlon_sfc=i_ges(1) nlat_sfc=i_ges(2) call mpi_bcast(lpl_dum,num_lpl,mpi_itype,npem1,mpi_comm_world,ierror) allocate(lpl_gfs(nlat_sfc/2)) allocate(dx_gfs(nlat_sfc/2)) lpl_gfs(1)=1 ! singularity at pole dx_gfs(1) = 360._r_kind / lpl_gfs(1) do j=2,nlat_sfc/2 lpl_gfs(j)=lpl_dum(j-1) dx_gfs(j) = 360._r_kind / lpl_gfs(j) enddo ! Allocate space for guess information files call create_gesfinfo ! Load time information for atm guess field sinfo into output arrays ntguessig = iamana(1) do i=1,nfldsig hrdifsig(i) = time_atm(i,1) ifilesig(i) = nint(time_atm(i,2)) hrdifsig_all(i) = hrdifsig(i) end do if(mype == 0) write(6,*)'READ_FILES: atm fcst files used in analysis : ',& (ifilesig(i),i=1,nfldsig),(hrdifsig(i),i=1,nfldsig),ntguessig if (ntguessig==0) then write(6,*)'READ_FILES: ***ERROR*** center atm fcst NOT AVAILABLE: PROGRAM STOPS' call stop2(99) endif if (l4densvar .and. nfldsig/=ntlevs_ens) then if (mype==0) then write(6,*)'READ_FILES: ***FATAL ERROR*** insufficient atm fcst for 4densvar: PROGRAM STOPS' do i=1,ntlevs_ens ihr=nhr_obsbin*(i-1)+nhr_half present=.false. do j=1,nfldsig if (ihr == ifilesig(j)) present=.true. end do if (.not.present) then write(filename,'(''sigf'',i2.2)')ihr write(6,*)'READ_FILES: ***FATAL ERROR*** file ',trim(filename),' missing: PROGRAM STOPS' endif end do endif call mpi_barrier(mpi_comm_world,ierror) call stop2(99) endif ! Load time information for surface guess field info into output arrays ntguessfc = iamana(2) do i=1,nfldsfc hrdifsfc(i) = time_sfc(i,1) ifilesfc(i) = nint(time_sfc(i,2)) hrdifsfc_all(i) = hrdifsfc(i) end do if(mype == 0) write(6,*)'READ_FILES: sfc fcst files used in analysis: ',& (ifilesfc(i),i=1,nfldsfc),(hrdifsfc(i),i=1,nfldsfc),ntguessfc if (ntguessfc==0) then write(6,*)'READ_FILES: ***ERROR*** center sfc fcst NOT AVAILABLE: PROGRAM STOPS' call stop2(99) endif if (l4densvar .and. nfldsfc/=ntlevs_ens) then if (mype==0) then write(6,*)'READ_FILES: ***FATAL ERROR*** insufficient sfc fcst for 4densvar: PROGRAM STOPS' do i=1,ntlevs_ens ihr=nhr_obsbin*(i-1)+nhr_half present=.false. do j=1,nfldsfc if (ihr == ifilesfc(j)) present=.true. end do if (.not.present) then write(filename,'(''sfcf'',i2.2)')ihr write(6,*)'READ_FILES: ***FATAL ERROR*** file ',trim(filename),' missing: PROGRAM STOPS' endif end do endif call mpi_barrier(mpi_comm_world,ierror) call stop2(99) endif deallocate(time_atm,time_sfc) ! Load time information for nst guess field info into output arrays ntguesnst = iamana(3) if ( nst_gsi > 0 ) then if ( sfcnst_comb ) then hrdifnst = hrdifsfc hrdifnst_all = hrdifnst if(mype == 0) write(6,*)'READ_FILES: hrdifnst = hrdifsfc' else do i=1,nfldnst hrdifnst(i) = time_nst(i,1) ifilenst(i) = nint(time_nst(i,2)) hrdifnst_all(i) = hrdifnst(i) end do if(mype == 0) write(6,*)'READ_FILES: nst fcst files used in analysis: ',& (ifilenst(i),i=1,nfldnst),(hrdifnst(i),i=1,nfldnst),ntguesnst if (ntguesnst==0) then write(6,*)'READ_FILES: ***ERROR*** center nst fcst NOT AVAILABLE: PROGRAM STOPS' call stop2(99) endif if (l4densvar .and. nfldnst/=ntlevs_ens) then write(6,*)'READ_FILES: ***ERROR*** insufficient nst fcst for 4densvar: PROGRAM STOPS' call stop2(99) endif endif deallocate(time_nst) endif ! for external aerosol files ! Load time information for aer guess field info into output arrays ntguesaer = iamana(4) if ( lread_ext_aerosol ) then do i=1,nfldaer hrdifaer(i) = time_aer(i,1) ifileaer(i) = nint(time_aer(i,2)) hrdifaer_all(i) = hrdifaer(i) end do if(mype == 0) write(6,*)'READ_FILES: aer fcst files used in analysis: ',& (ifileaer(i),i=1,nfldaer),(hrdifaer(i),i=1,nfldaer),ntguesaer if (ntguesaer==0) then write(6,*)'READ_FILES: ***ERROR*** center aer fcst NOT AVAILABLE: PROGRAM STOPS' call stop2(99) endif if (l4densvar .and. nfldaer/=ntlevs_ens) then write(6,*)'READ_FILES: ***ERROR*** insufficient aer fcst for 4densvar: PROGRAM STOPS' call stop2(99) endif deallocate(time_aer) endif ! End of routine return end subroutine read_files