#define TRACEBACK write (gol,'("in ",a," (",a,", line",i5,")")') rname, __FILE__, __LINE__; call goErr #define IF_NOTOK_RETURN(action) if (status/=0) then; TRACEBACK; action; return; end if #define IF_NOTOK_MDF(action) if (status/=0) then; TRACEBACK; action; call MDF_CLose(fid,status); status=1; return; end if #define IF_ERROR_RETURN(action) if (status> 0) then; TRACEBACK; action; return; end if #include "tm5.inc" !---------------------------------------------------------------------------- ! TM5 ! !---------------------------------------------------------------------------- !BOP ! ! !MODULE: Restart_remap ! ! !DESCRIPTION: Read 'rm' from restart files at one vertical resolution, ! use fill3d to expand to higher (current) vertical resolution. !\\ !\\ ! !INTERFACE: ! module Restart_remap ! ! !USES: ! use GO , only : gol, goPr, goErr use dims , only : nregions implicit none private ! ! !PUBLIC MEMBER FUNCTIONS: ! public :: Restart_remap_rm ! ! !PRIVATE DATA MEMBERS: ! character(len=*), parameter :: mname = 'Restart_remap' integer :: fid ! file id for IF_NOTOK_MDF macro ! ! !REVISION HISTORY: ! 3 Jan 2012 - Andy Jacobson, borrowed wholesale from Le Sager restart.F90 ! most comments and code are his. ! ! !REMARKS: ! (1) when an error occurs when accessing MDF files, you should first close ! the file before returning. The IF_NOTOK_MDF macro takes care of that. ! The only thing you need is to call it like that : ! ! IF_NOTOK_MDF(fid=xxxx) ! ! where you replace xxxx with the integer id (file handler) of the file ! you are accessing. ! !EOP !------------------------------------------------------------------------ contains ! ================================================================ subroutine Restart_remap_rm( status, region ) use GO , only : TrcFile, Init, Done, ReadRc use global_data, only : rcfile use global_data, only : outdir use meteodata , only : lli use GO , only : Set use GO , only : goMatchValue use dims , only : nregions, at, bt use dims , only : istart, idate, idatei ! use dims , only : nlon360, nlat180, iglbsfc use chem_param , only : ntracet, ntrace_chem, ntrace use chem_param , only : names, tracer_name_len use global_data , only : rcfile use partools , only : localComm use partools , only : myid, root use partools , only : ntracetloc, offsetn, tracer_active, tracer_loc use partools , only : lmloc, offsetl use partools , only : MPI_INFO_NULL use partools , only : previous_par use tracer_data , only : mass_dat, chem_dat use meteodata , only : lli, levi use meteodata , only : sp_dat, phlb_dat, m_dat !use meteodata , only : slhf_dat, sshf_dat use MDF , only : MDF_Open, MDF_Close, MDF_Inquire_Dimension use MDF , only : MDF_Inq_VarID, MDF_Inquire_Variable, MDF_Inq_DimID use MDf , only : MDF_Var_Par_Access, MDF_INDEPENDENT use MDF , only : MDF_Get_Att use MDF , only : MDF_Get_Var use MDF , only : MDF_READ, MDF_NETCDF4 ! for restart_remap (ARJ Jan 2012) ! use meteo, only : pressure_to_mass use Grid , only : TllGridInfo, TLevelInfo, Init, Done use Grid , only : HPressure use Grid , only : Fill3D ! ! !REVISION HISTORY: ! 3 Jan 2012 - Andy Jacobson - remap rm ! 8 Apr 2011 - P. Le Sager - use IF_NOTOK_MDF macro ! 28 Apr 2011 - P. Le Sager - Check on tracer availability in restart file. ! - Allows for more tracers in restart file than needed ! 10 May 2011 - P. Le Sager - Added deallocate statement to work with zoom regions ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC ! --- const -------------------------------- character(len=*), parameter :: rname = 'Restart_remap_rm' ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !INPUT PARAMETERS: ! integer, intent(in), optional :: region ! --- local -------------------------------- type(TrcFile) :: rcF character(len=256) :: rs_read_dir logical :: exist integer :: imr, jmr, lmr integer :: n, region1, region2 character(len=256) :: fname integer :: ncid integer :: varid_names integer :: iname character(len=tracer_name_len), allocatable :: values_names(:) integer :: itr, itr_loc, itr_file integer :: varid_rm integer :: ntracet_restart, dimid integer :: shp(2) integer :: lmr_restart, varid_at, varid_bt real, allocatable :: at_restart(:),bt_restart(:) type(TllGridInfo) :: lli_restart type(TLevelInfo) :: levi_restart real, allocatable :: inbuf(:,:,:) ! --- begin -------------------------------- write (gol,'("[restart_remap_rm] ")'); call goPr ! ! read settings from rcfile ! ! input directory: call Init( rcF, rcfile, status ) IF_NOTOK_RETURN(status=1) call ReadRc( rcF, 'restart.read.dir', rs_read_dir, status ) IF_NOTOK_RETURN(status=1) call Done( rcF, status ) IF_NOTOK_RETURN(status=1) ! correct istart ? if ( istart /= 61 ) then write (gol,'(" - skip; istart not 61 but ",i3)') istart; call goPr status=0; return end if ! init time ? if ( any( idate /= idatei ) ) then write (gol,'(" - skip; idate not idatei but ",i4,5i2.2)') idate; call goPr status=0; return end if ! region range: if ( present(region) ) then region1 = region region2 = region else region1 = 1 region2 = nregions end if ! loop over ns do n = region1, region2 ! grid size imr = lli(n)%nlon jmr = lli(n)%nlat lmr = levi%nlev ! name of restart file call Restart_FileName( n, fname, status, dir=trim(rs_read_dir) ) IF_NOTOK_RETURN(status=1) write (gol,'(" - reading tracer mass from ",a)') trim(fname); call goPr ! test ... inquire( file=fname, exist=exist ) if ( .not. exist ) then write (gol,'(" - restart file not found : ",a)') trim(fname); call goErr TRACEBACK; status=1; return end if ! o open netcdf file #ifdef MPI ! open existing file: call MDF_Open( trim(fname), MDF_NETCDF4, MDF_READ, ncid, status, & mpi_comm=localComm, mpi_info=MPI_INFO_NULL ) IF_NOTOK_RETURN(status=1) #else ! open existing file: call MDF_Open( trim(fname), MDF_NETCDF4, MDF_READ, ncid, status ) IF_NOTOK_RETURN(status=1) #endif ! check for vertical dimension compatibility call MDF_Inq_DimID( ncid, 'lev', dimid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Inquire_Dimension( ncid, dimid, status, length=lmr_restart ) IF_NOTOK_MDF(fid=ncid) write (gol,'(" - [restart_remap] Model compiled for ",i3," levels; restart file has ",i3" levels.")') lmr, lmr_restart; call goPr ! init levi_restart allocate(at_restart(lmr_restart+1)) allocate(bt_restart(lmr_restart+1)) ! call MDF_Inq_VarID( ncid, 'at', varid_at, status ) IF_NOTOK_MDF(fid=ncid) ! call MDF_Get_Var( ncid, varid_at, at_restart(1:(lmr_restart+1)), status ) IF_NOTOK_MDF(fid=ncid) ! call MDF_Inq_VarID( ncid, 'bt', varid_bt, status ) IF_NOTOK_MDF(fid=ncid) ! call MDF_Get_Var( ncid, varid_bt, bt_restart(1:(lmr_restart+1)), status ) IF_NOTOK_MDF(fid=ncid) ! call Init( levi_restart, lmr_restart, at_restart, bt_restart, status ) IF_NOTOK_RETURN(status=1) ! init lli_restart call Init( lli_restart, lli(n)%lon_deg(1), lli(n)%dlon_deg, lli(n)%im, & lli(n)%lat_deg(1), lli(n)%dlat_deg, lli(n)%jm, status ) IF_NOTOK_RETURN(status=1) call MDF_Inq_VarID( ncid, 'names', varid_names, status ) if ( status /= 0 ) then write (gol,'("could not find variable `names` in restart file;")'); call goErr write (gol,'(" using an old restart file to initialize the model ?")'); call goErr status=1 end if IF_NOTOK_MDF(fid=ncid) ! get dimension of "names" call MDF_Inquire_Variable( ncid, varid_names, status, shp=shp ) IF_NOTOK_MDF(fid=ncid) ! get number of transported tracer in restart file call MDF_Inq_DimID( ncid, 'trace_transp', dimid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Inquire_Dimension( ncid, dimid, status, length=ntracet_restart ) IF_NOTOK_MDF(fid=ncid) ! various id call MDF_Inq_VarID( ncid, 'rm', varid_rm, status ) IF_NOTOK_MDF(fid=ncid) ! set independent data mode: call MDF_Var_Par_Access( ncid, varid_rm , MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) ! read list with tracer names in rcfile: allocate( values_names(shp(2)) ) call MDF_Get_Var( ncid, varid_names, values_names, status ) IF_NOTOK_MDF(fid=ncid) allocate(inbuf(imr,jmr,lmr_restart)) ! loop over all tracers: do itr = 1, ntrace ! search in list: call goMatchValue( names(itr), values_names, itr_file, status ) if ( status < 0 ) then write(gol,'("Requested tracer `",a,"` not found in restart file!")') trim(names(itr)); call goErr IF_NOTOK_MDF(fid=ncid) endif ! transported or short lived tracer ? if ( itr <= ntracet ) then ! check ... if ( itr_file > ntracet_restart ) then write (gol,'("tracer `",a,"` is transported but seems to be not-transported in restart file")') trim(names(itr)); call goErr status=1 IF_NOTOK_MDF(fid=ncid) end if ! transported tracers distributed over tracers or levels: select case ( previous_par(n) ) case ( 'tracer' ) ! on this pe ? if ( tracer_active(itr) ) then ! local tracer index: itr_loc = tracer_loc(itr) call MDF_Get_Var( ncid, varid_rm, inbuf, & start=(/1,1,1,itr_file/), & count=(/imr,jmr,lmr_restart,1/), status=status ) IF_NOTOK_MDF(fid=ncid) call Fill3D( lli(n), levi, 'n', sp_dat(n)%data(1:imr,1:jmr,1), & mass_dat(n)%rm_t(1:imr,1:jmr,1:lmr,itr_loc), & lli_restart, levi_restart, inbuf(1:imr,1:jmr,1:lmr_restart), 'sum', status ) IF_NOTOK_RETURN(status=1) #ifdef slopes ! Ensure that slopes are initialized to "unset" values of 0.0. ! Wouter says that we could remap levels for rxm et al., but ! 0.0 will also work. The noise induced from remapping the ! rm array is almost certainly bigger than any issues from ! having this "default=0.0" slopes information. -ARJ 1 Jan 12 mass_dat(n)%rxm_t(1:imr,1:jmr,1:lmr,itr_loc:itr_loc) = 0.0 mass_dat(n)%rym_t(1:imr,1:jmr,1:lmr,itr_loc:itr_loc) = 0.0 mass_dat(n)%rzm_t(1:imr,1:jmr,1:lmr,itr_loc:itr_loc) = 0.0 #endif end if ! tracer active case default write (gol,'("remapping only supported for distribution of tracer by PE (not levels or otherwise).")'); call goPr write (gol,'("in ",a)') rname; call goErr; status=1; return end select end if end do ! tracers ! free mem for next region DeAllocate( values_names) deallocate(inbuf) deallocate(at_restart) deallocate(bt_restart) call MDF_Close( ncid, status ) IF_NOTOK_RETURN(status=1) call Done(levi_restart, status) IF_NOTOK_RETURN(status=1) call Done(lli_restart, status) IF_NOTOK_RETURN(status=1) end do ! regions ! ok status = 0 end subroutine Restart_remap_rm ! *** ! direct copy of subroutine from restart module. - ARJ subroutine Restart_FileName( region, fname, status, key, dir, isfirst ) use dims , only : idate use global_data, only : outdir use meteodata , only : lli ! --- in/out ------------------------------- integer, intent(in) :: region character(len=*), intent(out) :: fname integer, intent(out) :: status character(len=*), intent(in), optional :: dir character(len=*), intent(in), optional :: key logical, intent(in), optional :: isfirst ! --- local ------------------------------- character(len=*), parameter :: rname = 'Restart_FileName' character(len=256) :: adir character(len=32) :: akey ! --- begin -------------------------------- ! destination directory: adir = trim(outdir) if ( present(dir) ) adir = trim(dir) ! extra key, for example '_x' to denote that ! a restart file was dumped after process 'x': akey = '' if ( present(key) ) akey = trim(key) ! if this is the initial time, add an extra key to avoid ! that the restart file for this hour from the previous ! run is overwritten: if ( present(isfirst) ) then if ( isfirst ) akey = trim(akey)//'_initial' end if ! write filename: write (fname,'(a,"/TM5_restart_",i4.4,2i2.2,"_",2i2.2,"_",a,a,".nc")') & trim(adir), idate(1:5), trim(lli(region)%name), trim(akey) ! ok status = 0 end subroutine Restart_FileName end module Restart_remap