#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 ! ! !DESCRIPTION: Write and read restart files. !\\ !\\ ! !INTERFACE: ! module Restart ! ! !USES: ! use GO , only : gol, goPr, goErr use dims , only : nregions implicit none private ! ! !PUBLIC MEMBER FUNCTIONS: ! public :: Restart_Init public :: Restart_Done public :: Restart_Save public :: Restart_Write public :: Restart_Read public :: rs_write ! ! !PRIVATE DATA MEMBERS: ! character(len=*), parameter :: mname = 'Restart' character(len=256) :: rs_write_dir logical :: rs_write logical :: rs_write_extra integer :: rs_write_extra_dhour, rs_write_extra_hour integer :: fid ! file id for IF_NOTOK_MDF macro ! ! !REVISION HISTORY: ! 8 Apr 2011 - P. Le Sager - Close MDF file if error occurs. This is ! needed for mpi_abort not to hang. See TM5_MPI_Abort in ! partools, and remarks below. Made IF_NOTOK_MDF macro for ! that purpose. ! 28 Apr 2011 - P. Le Sager - Read method : handle restart file with extra ! tracers. ! ! !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_Init( status ) use GO , only : TrcFile, Init, Done, ReadRc use global_data, only : rcfile use global_data, only : outdir use meteodata , only : lli ! --- in/out ------------------------------- integer, intent(out) :: status ! --- const -------------------------------- character(len=*), parameter :: rname = 'Restart_Init' ! --- local -------------------------------- type(TrcFile) :: rcF ! --- begin -------------------------------- ! ! read settings from rcfile ! ! open: call Init( rcF, rcfile, status ) IF_NOTOK_RETURN(status=1) ! write restart files at all ? call ReadRc( rcF, 'restart.write', rs_write, status, default=.false. ) IF_ERROR_RETURN(status=1) ! further settings ... if ( rs_write ) then ! output directory: call ReadRc( rcF, 'restart.write.dir', rs_write_dir, status, default=outdir ) IF_ERROR_RETURN(status=1) ! extra restart files ? call ReadRc( rcF, 'restart.write.extra', rs_write_extra, status, default=.false. ) IF_ERROR_RETURN(status=1) if ( rs_write_extra ) then call ReadRc( rcF, 'restart.write.extra.hour', rs_write_extra_hour, status, default=0 ) IF_ERROR_RETURN(status=1) call ReadRc( rcF, 'restart.write.extra.dhour', rs_write_extra_dhour, status, default=24 ) IF_ERROR_RETURN(status=1) end if end if ! write restart files ! close: call Done( rcF, status ) IF_NOTOK_RETURN(status=1) ! ! done ! ! ok status = 0 end subroutine Restart_Init ! *** subroutine Restart_Done( status ) ! --- in/out ------------------------------- integer, intent(out) :: status ! --- const -------------------------------- character(len=*), parameter :: rname = 'Restart_Done' ! --- begin -------------------------------- ! nothing to be done ... ! ok status = 0 end subroutine Restart_Done ! *** subroutine Restart_Save( status, extra, isfirst ) use dims, only : idate ! --- in/out ------------------------------- integer, intent(out) :: status logical, intent(in), optional :: extra logical, intent(in), optional :: isfirst ! --- const -------------------------------- character(len=*), parameter :: rname = 'Restart_Save' ! --- local -------------------------------- logical :: is_extra ! --- begin -------------------------------- ! options ... is_extra = .false. if ( present(extra) ) is_extra = extra ! write restart files at all ? if ( rs_write ) then ! end or extra ? if ( is_extra ) then ! save extra restart files ? if ( rs_write_extra ) then ! every hour+n*dhour only : if ( modulo( idate(4) - rs_write_extra_hour, rs_write_extra_dhour ) == 0 .and. & all( idate(5:6) == 0 ) ) then ! write restart file for this time: call Restart_Write( status, isfirst=isfirst ) IF_NOTOK_RETURN(status=1) end if ! for this hour end if ! extra restart files ? else ! write restart file : call Restart_Write( status, isfirst=isfirst ) IF_NOTOK_RETURN(status=1) end if ! not extra end if ! write at all ! ok status = 0 end subroutine Restart_Save ! *** 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 !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Restart_Write ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! subroutine Restart_Write( status, key, region, isfirst ) use GO , only : Get use dims , only : nregions, at, bt ! use dims , only : nlon360, nlat180, iglbsfc use chem_param , only : ntracet, ntrace_chem, ntrace, names use partools , only : localComm use partools , only : myid, root use partools , only : ntracetloc, offsetn use partools , only : lmloc, offsetl use partools , only : MPI_INFO_NULL use partools , only : previous_par use tracer_data , only : mass_dat, chem_dat #ifdef with_tendencies use tm5_tendency, only : plc_ntr, plc_trname use tm5_tendency, only : plc_npr, plc_prname use tracer_data , only : plc_dat #endif use meteodata , only : lli, levi use meteodata , only : sp_dat, phlb_dat, m_dat !use meteodata , only : sp1_dat, sp2_dat ! debug !use meteodata , only : mfu_dat, mfv_dat, mfw_dat ! debug !use meteodata , only : pu_dat, pv_dat, pw_dat ! debug !use meteodata , only : sshf_dat, slhf_dat use MDF , only : MDF_Create, MDF_EndDef, MDF_Close use MDF , only : MDF_Var_Par_Access, MDF_INDEPENDENT use MDF , only : MDF_Def_Dim, MDF_Def_Var use MDF , only : MDF_Put_Att use MDF , only : MDF_Put_Var use MDF , only : MDF_REPLACE, MDF_NETCDF, MDF_NETCDF4 use MDF , only : MDF_FLOAT, MDF_DOUBLE, MDF_CHAR ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !INPUT PARAMETERS: ! character(len=*), intent(in), optional :: key integer, intent(in), optional :: region logical, intent(in), optional :: isfirst ! ! !REVISION HISTORY: ! 8 Apr 2011 - P. Le Sager - use IF_NOTOK_MDF macro ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC ! --- local -------------------------------- character(len=*), parameter :: rname = 'Restart_Write' integer :: imr, jmr, lmr integer :: iregion character(len=256) :: fname integer :: ftype integer :: ncid integer :: dimid_lon, dimid_lat, dimid_lev, dimid_hlev integer :: dimid_lon_sfc, dimid_lat_sfc integer :: dimid_trace, dimid_trace_transp, dimid_trace_chem integer :: dimid_name integer :: varid integer :: varid_sp, varid_ph, varid_m, varid_at, varid_bt !integer :: varid_sp1, varid_sp2 ! debug !integer :: varid_mfu, varid_mfv, varid_mfw ! debug !integer :: varid_pu, varid_pv, varid_pw ! debug !integer :: varid_sshf, varid_slhf integer :: varid_names integer :: varid_rm #ifdef slopes integer :: varid_rxm, varid_rym, varid_rzm #endif integer :: varid_rmc #ifdef with_tendencies integer :: varid_plc(plc_ntr,plc_npr) integer :: itr, ipr integer :: time6(6) #endif integer :: rtype ! --- begin -------------------------------- write (gol,'("write restart file(s) ...")'); call goPr ! loop over regions: do iregion = 1, nregions ! only selected region ? if ( present(region) ) then if ( iregion /= region ) cycle end if ! grid size imr = lli(iregion)%nlon jmr = lli(iregion)%nlat lmr = levi%nlev ! name of restart file call Restart_FileName( iregion, fname, status, key=key, dir=rs_write_dir, isfirst=isfirst ) IF_NOTOK_RETURN(status=1) write (gol,'(" destination : ",a)') trim(fname); call goPr ! o open netcdf file #ifdef MPI ! overwrite existing files (clobber), provide MPI stuff: call MDF_Create( trim(fname), MDF_NETCDF4, MDF_REPLACE, ncid, status, & mpi_comm=localComm, mpi_info=MPI_INFO_NULL ) if (status/=0) then write (gol,'("from creating NetCDF4 file for writing in parallel;")'); call goErr write (gol,'("MDF module not compiled with netcdf4 support ?")'); call goErr TRACEBACK; status=1; return end if #else ! overwrite existing files (clobber) call MDF_Create( trim(fname), MDF_NETCDF4, MDF_REPLACE, ncid, status ) IF_NOTOK_RETURN(status=1) #endif ! o define dimensions call MDF_Def_Dim( ncid, 'lon', lli(iregion)%nlon, dimid_lon, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Def_Dim( ncid, 'lat', lli(iregion)%nlat, dimid_lat, status ) IF_NOTOK_MDF(fid=ncid) ! These Dims are not used, so commented (pls, 30-03-2011) ! call MDF_Def_Dim( ncid, 'lon_sfc', lli(iglbsfc)%nlon, dimid_lon_sfc, status ) ! IF_NOTOK_MDF(fid=ncid) ! call MDF_Def_Dim( ncid, 'lat_sfc', lli(iglbsfc)%nlat, dimid_lat_sfc, status ) ! IF_NOTOK_MDF(fid=ncid) call MDF_Def_Dim( ncid, 'lev', levi%nlev, dimid_lev, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Def_Dim( ncid, 'hlev', levi%nlev+1, dimid_hlev, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Def_Dim( ncid, 'trace_transp', ntracet, dimid_trace_transp, status ) IF_NOTOK_MDF(fid=ncid) if ( ntrace_chem > 0 ) then call MDF_Def_Dim( ncid, 'trace_chem', ntrace_chem, dimid_trace_chem, status ) IF_NOTOK_MDF(fid=ncid) else dimid_trace_chem = -1 end if call MDF_Def_Dim( ncid, 'trace', ntrace, dimid_trace, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Def_Dim( ncid, 'name', len(names(1)), dimid_name, status ) IF_NOTOK_MDF(fid=ncid) ! o define variables select case ( kind(m_dat(iregion)%data) ) case ( 4 ) ; rtype = MDF_FLOAT case ( 8 ) ; rtype = MDF_DOUBLE case default write (gol,'("unsupported real kind : ",i6)') kind(m_dat(iregion)%data) TRACEBACK; status=1; return end select ! surface pressure call MDF_Def_Var( ncid, 'sp', rtype, & (/dimid_lon,dimid_lat/), varid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'long_name', 'surface pressure', status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'unit', 'Pa', status ) IF_NOTOK_MDF(fid=ncid) varid_sp = varid ! at, bt coefficients for hybrid grid call MDF_Def_Var( ncid, 'at', rtype, (/dimid_hlev/), varid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'long_name', 'hybrid grid a_t coefficient', status ) IF_NOTOK_MDF(fid=ncid) varid_at = varid call MDF_Def_Var( ncid, 'bt', rtype, (/dimid_hlev/), varid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'long_name', 'hybrid grid b_t coefficient', status ) IF_NOTOK_MDF(fid=ncid) varid_bt = varid !! debug ... !MDF_Def_Var( ncid, 'sp1', rtype, (/dimid_lon,dimid_lat/), varid_sp1, status ) !IF_NOTOK_MDF(fid=ncid) !MDF_Def_Var( ncid, 'sp2', rtype, (/dimid_lon,dimid_lat/), varid_sp2, status ) !IF_NOTOK_MDF(fid=ncid) ! half level pressure call MDF_Def_Var( ncid, 'ph', rtype, & (/dimid_lon,dimid_lat,dimid_hlev/), varid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'long_name', 'half level pressure', status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'unit', 'Pa', status ) IF_NOTOK_MDF(fid=ncid) varid_ph = varid ! air mass call MDF_Def_Var( ncid, 'm', rtype, & (/dimid_lon,dimid_lat,dimid_lev/), varid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'long_name', 'air mass', status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'unit', 'kg', status ) IF_NOTOK_MDF(fid=ncid) varid_m = varid !! debug !call MDF_Def_Var( ncid, 'mfu', rtype, (/dimid_lon,dimid_lat,dimid_lev/), varid_mfu, status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Def_Var( ncid, 'mfv', rtype, (/dimid_lon,dimid_lat,dimid_lev/), varid_mfv, status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Def_Var( ncid, 'mfw', rtype, (/dimid_lon,dimid_lat,dimid_lev/), varid_mfw, status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Def_Var( ncid, 'pu', rtype, (/dimid_lon,dimid_lat,dimid_lev/), varid_pu, status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Def_Var( ncid, 'pv', rtype, (/dimid_lon,dimid_lat,dimid_lev/), varid_pv, status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Def_Var( ncid, 'pw', rtype, (/dimid_lon,dimid_lat,dimid_lev/), varid_pw, status ) !IF_NOTOK_MDF(fid=ncid) !! accumulated surface fluxes !! !call MDF_Def_Var( ncid, 'slhf', rtype, (/dimid_lon_sfc,dimid_lat_sfc/), varid, status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Put_Att( ncid, varid, 'long_name', 'surface latent heat flux', status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Put_Att( ncid, varid, 'unit', 'W/m2', status ) !IF_NOTOK_MDF(fid=ncid) !varid_slhf = varid !! !call MDF_Def_Var( ncid, 'sshf', rtype, (/dimid_lon_sfc,dimid_lat_sfc/), varid, status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Put_Att( ncid, varid, 'long_name', 'surface sensible heat flux', status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Put_Att( ncid, varid, 'unit', 'W/m2', status ) !IF_NOTOK_MDF(fid=ncid) !varid_sshf = varid ! tracer names call MDF_Def_Var( ncid, 'names', MDF_CHAR, (/dimid_name,dimid_trace/), varid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'long_name', 'tracer names', status ) IF_NOTOK_MDF(fid=ncid) varid_names = varid ! tracer mass call MDF_Def_Var( ncid, 'rm', rtype, & (/dimid_lon,dimid_lat,dimid_lev,dimid_trace_transp/), varid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'long_name', 'transported tracer mass', status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'unit', 'kg', status ) IF_NOTOK_MDF(fid=ncid) varid_rm = varid ! tracer mass slopes: #ifdef slopes call MDF_Def_Var( ncid, 'rxm', rtype, & (/dimid_lon,dimid_lat,dimid_lev,dimid_trace_transp/), varid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'long_name', 'tracer mass slope in x direction', status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'unit', 'kg/(half cell)', status ) IF_NOTOK_MDF(fid=ncid) varid_rxm = varid call MDF_Def_Var( ncid, 'rym', rtype, & (/dimid_lon,dimid_lat,dimid_lev,dimid_trace_transp/), varid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'long_name', 'tracer mass slope in y direction', status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'unit', 'kg/(half cell)', status ) IF_NOTOK_MDF(fid=ncid) varid_rym = varid call MDF_Def_Var( ncid, 'rzm', rtype, & (/dimid_lon,dimid_lat,dimid_lev,dimid_trace_transp/), varid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'long_name', 'tracer mass slope in z direction', status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'unit', 'kg/(half cell)', status ) IF_NOTOK_MDF(fid=ncid) varid_rzm = varid #endif ! non-transported tracers: if ( ntrace_chem > 0 ) then call MDF_Def_Var( ncid, 'rmc', rtype, & (/dimid_lon,dimid_lat,dimid_lev,dimid_trace_chem/), varid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'long_name', 'non-transported tracer mass', status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'unit', 'kg', status ) IF_NOTOK_MDF(fid=ncid) varid_rmc = varid end if #ifdef with_tendencies ! production, loss, and concentration: do itr = 1, plc_ntr do ipr = 1, plc_npr ! define netcdf variable: call MDF_Def_Var( ncid, trim(plc_trname(itr))//'_'//trim(plc_prname(ipr)), rtype, & (/dimid_lon,dimid_lat,dimid_lev/), varid, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'long_name', 'chemical tendency', status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Att( ncid, varid, 'unit', trim(plc_dat(region,itr,ipr)%unit), status ) IF_NOTOK_MDF(fid=ncid) ! extract time as 6 integers: call Get( plc_dat(region,itr,ipr)%t, time6=time6 ) ! add time attribute: call MDF_Put_Att( ncid, varid, 'time', time6, status ) IF_NOTOK_MDF(fid=ncid) ! store variable id: varid_plc(itr,ipr) = varid end do end do #endif ! o end definition mode call MDF_EndDef( ncid, status ) IF_NOTOK_MDF(fid=ncid) ! o write variables ! set independent data mode (not all processes are writing): call MDF_Var_Par_Access( ncid, varid_sp , MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Var_Par_Access( ncid, varid_ph , MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Var_Par_Access( ncid, varid_m , MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Var_Par_Access( ncid, varid_at , MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Var_Par_Access( ncid, varid_bt , MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) !call MDF_Var_Par_Access( ncid, varid_slhf, MDF_INDEPENDENT, status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Var_Par_Access( ncid, varid_sshf, MDF_INDEPENDENT, status ) !IF_NOTOK_MDF(fid=ncid) ! write meteo etc from root: if ( myid == root ) then ! surface pressure call MDF_Put_Var( ncid, varid_sp, sp_dat(iregion)%data(1:imr,1:jmr,1), status ) IF_NOTOK_MDF(fid=ncid) !! debug ... !call MDF_Put_Var( ncid, varid_sp1, sp1_dat(iregion)%data(1:imr,1:jmr,1), status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Put_Var( ncid, varid_sp2, sp2_dat(iregion)%data(1:imr,1:jmr,1), status ) !IF_NOTOK_MDF(fid=ncid) ! half level pressure call MDF_Put_Var( ncid, varid_ph, phlb_dat(iregion)%data(1:imr,1:jmr,1:lmr+1), status ) IF_NOTOK_MDF(fid=ncid) ! at, bt coefficients call MDF_Put_Var( ncid, varid_at, at(1:lmr+1), status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Var( ncid, varid_bt, bt(1:lmr+1), status ) IF_NOTOK_MDF(fid=ncid) ! air mass call MDF_Put_Var( ncid, varid_m, m_dat(iregion)%data(1:imr,1:jmr,1:lmr), status ) IF_NOTOK_MDF(fid=ncid) !! debug ... !call MDF_Put_Var( ncid, varid_mfu, mfu_dat(iregion)%data(1:imr,1:jmr,1:lmr), status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Put_Var( ncid, varid_mfv, mfv_dat(iregion)%data(1:imr,1:jmr,1:lmr), status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Put_Var( ncid, varid_mfw, mfw_dat(iregion)%data(1:imr,1:jmr,1:lmr), status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Put_Var( ncid, varid_pu, pu_dat(iregion)%data(1:imr,1:jmr,1:lmr), status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Put_Var( ncid, varid_pv, pv_dat(iregion)%data(1:imr,1:jmr,1:lmr), status ) !IF_NOTOK_MDF(fid=ncid) !call MDF_Put_Var( ncid, varid_pw, pw_dat(iregion)%data(1:imr,1:jmr,1:lmr), status ) !IF_NOTOK_MDF(fid=ncid) !! surface latent heat flux; global surface field ! !call MDF_Put_Var( ncid, varid_slhf, slhf_dat(iglbsfc)%data(1:nlon360,1:nlat180,1), status ) !IF_NOTOK_MDF(fid=ncid) ! !! surface sensible heat flux; global surface field ! !call MDF_Put_Var( ncid, varid_sshf, sshf_dat(iglbsfc)%data(1:nlon360,1:nlat180,1), status ) !IF_NOTOK_MDF(fid=ncid) ! tracer names: call MDF_Put_Var( ncid, varid_names, names, status ) IF_NOTOK_MDF(fid=ncid) end if ! root ! set independent data mode (not all processes are writing): call MDF_Var_Par_Access( ncid, varid_rm , MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) #ifdef slopes call MDF_Var_Par_Access( ncid, varid_rxm, MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Var_Par_Access( ncid, varid_rym, MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Var_Par_Access( ncid, varid_rzm, MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) #endif ! write transported tracers, distributed over tracers or levels: select case ( previous_par(iregion) ) case ( 'tracer' ) ! distributed over tracers; might be empty on this pe if ( ntracetloc > 0 ) then call MDF_Put_Var( ncid, varid_rm, & mass_dat(iregion)%rm_t(1:imr,1:jmr,1:lmr,1:ntracetloc), status, & start=(/1,1,1,offsetn+1/), & count=(/imr,jmr,lmr,ntracetloc/) ) IF_NOTOK_MDF(fid=ncid) #ifdef slopes call MDF_Put_Var( ncid, varid_rxm, & mass_dat(iregion)%rxm_t(1:imr,1:jmr,1:lmr,1:ntracetloc), status, & start=(/1,1,1,offsetn+1/), & count=(/imr,jmr,lmr,ntracetloc/) ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Var( ncid, varid_rym, & mass_dat(iregion)%rym_t(1:imr,1:jmr,1:lmr,1:ntracetloc), status, & start=(/1,1,1,offsetn+1/), & count=(/imr,jmr,lmr,ntracetloc/) ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Var( ncid, varid_rzm, & mass_dat(iregion)%rzm_t(1:imr,1:jmr,1:lmr,1:ntracetloc), status, & start=(/1,1,1,offsetn+1/), & count=(/imr,jmr,lmr,ntracetloc/) ) IF_NOTOK_MDF(fid=ncid) #endif end if ! ntraceloc > 0 case ( 'levels' ) ! distributed over levels; might be empty on this pe if ( lmloc > 0 ) then call MDF_Put_Var( ncid, varid_rm, & mass_dat(iregion)%rm_k(1:imr,1:jmr,1:lmloc,1:ntracet), status, & start=(/1,1,offsetl+1,1/), & count=(/imr,jmr,lmloc,ntracet/) ) IF_NOTOK_MDF(fid=ncid) #ifdef slopes call MDF_Put_Var( ncid, varid_rxm, & mass_dat(iregion)%rxm_k(1:imr,1:jmr,1:lmloc,1:ntracet), status, & start=(/1,1,offsetl+1,1/), & count=(/imr,jmr,lmloc,ntracet/) ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Var( ncid, varid_rym, & mass_dat(iregion)%rym_k(1:imr,1:jmr,1:lmloc,1:ntracet), status, & start=(/1,1,offsetl+1,1/), & count=(/imr,jmr,lmloc,ntracet/) ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Var( ncid, varid_rzm, & mass_dat(iregion)%rzm_k(1:imr,1:jmr,1:lmloc,1:ntracet), status, & start=(/1,1,offsetl+1,1/), & count=(/imr,jmr,lmloc,ntracet/) ) IF_NOTOK_MDF(fid=ncid) #endif end if ! ntraceloc > 0 case default write (gol,'("unsupported distribution:")'); call goPr write (gol,'(" region : ",i2)') iregion ;call goPr write (gol,'(" previous_par : ",a)') previous_par(iregion) ;call goPr TRACEBACK; status=1; return end select ! write non-transported tracers levels on this pe: if ( (ntrace_chem > 0) .and. (lmloc > 0) ) then ! set independent data mode (not all processes are writing): call MDF_Var_Par_Access( ncid, varid_rmc, MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Put_Var( ncid, varid_rmc, & chem_dat(iregion)%rm_k(1:imr,1:jmr,1:lmloc,ntracet+1:ntracet+ntrace_chem), status, & start=(/1,1,offsetl+1,1/), & count=(/imr,jmr,lmloc,ntrace_chem/) ) IF_NOTOK_MDF(fid=ncid) end if ! lmloc > 0 #ifdef with_tendencies ! set independent data mode (not all processes are writing): do itr = 1, plc_ntr do ipr = 1, plc_npr call MDF_Var_Par_Access( ncid, varid_plc(itr,ipr) , MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) end do end do ! write production/loss/concentration levels on this pe: if ( lmloc > 0 ) then do itr = 1, plc_ntr do ipr = 1, plc_npr call MDF_Put_Var( ncid, varid_plc(itr,ipr), & plc_dat(region,itr,ipr)%rm_k(1:imr,1:jmr,1:lmloc), status, & start=(/1,1,offsetl+1/), count=(/imr,jmr,lmloc/) ) IF_NOTOK_MDF(fid=ncid) end do end do end if #endif ! o close file call MDF_Close( ncid, status ) IF_NOTOK_RETURN(status=1) end do ! regions !write (gol,'(" ok")'); call goPr ! ok status = 0 end subroutine Restart_Write !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: Restart_Read ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! subroutine Restart_Read( status, & region, & surface_pressure, pressure, air_mass, surface_fluxes, & tracer_mass, tendencies, megan_history, nox_pulsing, & override_istart) ! ! !USES: ! use GO , only : TrcFile, Init, Done, ReadRc use GO , only : Set use GO , only : goMatchValue use dims , only : nregions 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 #ifdef with_tendencies use tm5_tendency, only : plc_ntr, plc_trname use tm5_tendency, only : plc_npr, plc_prname use tracer_data , only : plc_dat #endif 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 ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !INPUT PARAMETERS: ! integer, intent(in), optional :: region logical, intent(in), optional :: surface_pressure, pressure, air_mass, surface_fluxes logical, intent(in), optional :: tracer_mass, tendencies, megan_history, nox_pulsing logical, intent(in), optional :: override_istart ! ! !REVISION HISTORY: ! 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 = mname//'/Restart_Read' ! --- local -------------------------------- type(TrcFile) :: rcF character(len=256) :: rs_read_dir logical :: exist logical :: do_sp, do_ph, do_m, do_sflux, do_rm, do_plc, do_megan, do_pulse integer :: imr, jmr, lmr integer :: n, region1, region2 character(len=256) :: fname integer :: ncid integer :: varid_sp integer :: varid_ph integer :: varid_m !integer :: varid_slhf, varid_sshf 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) #ifdef slopes integer :: varid_rxm, varid_rym, varid_rzm #endif integer :: varid_rmc #ifdef with_tendencies integer :: varid_plc(plc_ntr,plc_npr) integer :: itr, ipr integer :: time6(6) #endif ! --- begin -------------------------------- write (gol,'("read restart file ...")'); call goPr ! correct istart ? if ( istart /= 33 .and. .not. present(override_istart) ) then write (gol,'(" skip; istart not 33 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 ! 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) ! region range: if ( present(region) ) then region1 = region region2 = region else region1 = 1 region2 = nregions end if ! data sets: do_sp = .false. ; if ( present(surface_pressure ) ) do_sp = surface_pressure do_ph = .false. ; if ( present(pressure ) ) do_ph = pressure do_m = .false. ; if ( present(air_mass ) ) do_m = air_mass do_sflux = .false. ; if ( present(surface_fluxes ) ) do_sflux = surface_fluxes do_rm = .false. ; if ( present(tracer_mass ) ) do_rm = tracer_mass do_plc = .false. ; if ( present(tendencies ) ) do_plc = tendencies do_megan = .false. ; if ( present(megan_history ) ) do_megan = megan_history do_pulse = .false. ; if ( present(nox_pulsing ) ) do_pulse = nox_pulsing ! sorry .. if ( do_sflux ) then write (gol,'("no surface fluxes in restart files until somebody")'); call goErr write (gol,'("has a good idea on what should be storred:")'); call goErr write (gol,'(" o global surface field (1x1 ?)")'); call goErr write (gol,'(" o zoom regions")'); call goErr write (gol,'(" o both")'); call goErr TRACEBACK; status=1; return 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,'(" restore 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 #if defined MPI && defined with_netcdf4_par ! 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 ! o define variables ! surface pressure if ( do_sp ) then call MDF_Inq_VarID( ncid, 'sp', varid_sp, status ) IF_NOTOK_MDF(fid=ncid) end if ! half level pressure if ( do_ph ) then call MDF_Inq_VarID( ncid, 'ph', varid_ph, status ) IF_NOTOK_MDF(fid=ncid) end if ! air mass if ( do_m ) then call MDF_Inq_VarID( ncid, 'm', varid_m, status ) IF_NOTOK_MDF(fid=ncid) end if !! surface fluxes; only once is ok !if ( do_sflux .and. (n==region1) ) then ! call MDF_Inq_VarID( ncid, 'slhf', varid_slhf, status ) ! IF_NOTOK_MDF(fid=ncid) ! call MDF_Inq_VarID( ncid, 'sshf', varid_sshf, status ) ! IF_NOTOK_MDF(fid=ncid) !end if ! tracer mass if ( do_rm ) then 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) #ifdef slopes call MDF_Inq_VarID( ncid, 'rxm', varid_rxm, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Inq_VarID( ncid, 'rym', varid_rym, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Inq_VarID( ncid, 'rzm', varid_rzm, status ) IF_NOTOK_MDF(fid=ncid) #endif ! read non-transported tracers levels on this pe: if ( ntrace_chem > 0 ) then call MDF_Inq_VarID( ncid, 'rmc', varid_rmc, status ) IF_NOTOK_MDF(fid=ncid) end if end if #ifdef with_tendencies ! production/loss/concentration if ( do_plc ) then do itr = 1, plc_ntr do ipr = 1, plc_npr call MDF_Inq_VarID( ncid, trim(plc_trname(itr))//'_'//trim(plc_prname(ipr)), varid_plc(itr,ipr), status ) IF_NOTOK_MDF(fid=ncid) end do end do end if #endif ! o read variables ! read meteo on each processor: ! surface pressure if ( do_sp ) then write (gol,'(" restore surface pressure ...")'); call goPr call MDF_Var_Par_Access( ncid, varid_sp, MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Get_Var( ncid, varid_sp, sp_dat(n)%data(1:imr,1:jmr,1), status ) IF_NOTOK_MDF(fid=ncid) end if ! half level pressure if ( do_ph ) then write (gol,'(" restore half level pressure ...")'); call goPr call MDF_Var_Par_Access( ncid, varid_ph, MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Get_Var( ncid, varid_ph, phlb_dat(n)%data(1:imr,1:jmr,1:lmr+1), status ) IF_NOTOK_MDF(fid=ncid) end if ! air mass if ( do_m ) then write (gol,'(" restore air mass ...")'); call goPr call MDF_Var_Par_Access( ncid, varid_m, MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Get_Var( ncid, varid_m, m_dat(n)%data(1:imr,1:jmr,1:lmr), status ) IF_NOTOK_MDF(fid=ncid) end if !! surface fluxes !if ( do_sflux ) then ! write (gol,'(" restore surface fluxes ...")'); call goPr ! if ( slhf_dat(iglbsfc)%used ) then ! call MDF_Var_Par_Access( ncid, varid_slhf, MDF_INDEPENDENT, status ) ! IF_NOTOK_MDF(fid=ncid) ! call MDF_Get_Var( ncid, varid_slhf, slhf_dat(iglbsfc)%data(1:nlon360,1:nlat180,1), status ) ! IF_NOTOK_MDF(fid=ncid) ! end if ! if ( sshf_dat(iglbsfc)%used ) then ! call MDF_Var_Par_Access( ncid, varid_sshf, MDF_INDEPENDENT, status ) ! IF_NOTOK_MDF(fid=ncid) ! call MDF_Get_Var( ncid, varid_sshf, sshf_dat(iglbsfc)%data(1:nlon360,1:nlat180,1), status ) ! IF_NOTOK_MDF(fid=ncid) ! end if !end if ! tracer mass if ( do_rm ) then write (gol,'(" restore tracer mass ...")'); call goPr ! set independent data mode: call MDF_Var_Par_Access( ncid, varid_rm , MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) #ifdef slopes call MDF_Var_Par_Access( ncid, varid_rxm, MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Var_Par_Access( ncid, varid_rym, MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) call MDF_Var_Par_Access( ncid, varid_rzm, MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) #endif if ( ntrace_chem > 0 ) then call MDF_Var_Par_Access( ncid, varid_rmc, MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) end if ! 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) ! loop over all tracers: do itr = 1, ntrace !WP! We skip this step as the tracer names in CTDAS are not unique !WP! Instead, simply assume the same number for each tracer in file and simulation ! search in list: !WP! call goMatchValue( names(itr), values_names, itr_file, status ) itr_file=itr status = 0 if ( status < 0 ) then !WP! kind of obsolete now write(gol,'("Requested tracer `",a,"` not found in restart file!")') trim(names(itr)); call goErr IF_NOTOK_MDF(fid=ncid) !! TODO : !! insert (optional) default startup values here ?? !! User can always add those in the restart file first, !! outside TM5... probably easier, faster, and more flexible. endif ! check ... if ( trim(names(itr)) /= trim(values_names(itr_file)) ) then write (gol,'("tracer `",a,"` does not have the same name as its counterpart `",a"` in the restart file")') trim(names(itr)), trim(values_names(itr_file)); call goErr status=1 IF_NOTOK_MDF(fid=ncid) end if ! 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) ! read tracer mass: call MDF_Get_Var( ncid, varid_rm, & mass_dat(n)%rm_t(1:imr,1:jmr,1:lmr,itr_loc:itr_loc), & status, & start=(/1,1,1,itr_file/), & count=(/imr,jmr,lmr,1/) ) IF_NOTOK_MDF(fid=ncid) #ifdef slopes ! read slopes: call MDF_Get_Var( ncid, varid_rxm, & mass_dat(n)%rxm_t(1:imr,1:jmr,1:lmr,itr_loc:itr_loc), & status, & start=(/1,1,1,itr_file/), & count=(/imr,jmr,lmr,1/) ) IF_NOTOK_MDF(fid=ncid) call MDF_Get_Var( ncid, varid_rym, & mass_dat(n)%rym_t(1:imr,1:jmr,1:lmr,itr_loc:itr_loc), & status, & start=(/1,1,1,itr_file/), & count=(/imr,jmr,lmr,1/) ) IF_NOTOK_MDF(fid=ncid) call MDF_Get_Var( ncid, varid_rzm, & mass_dat(n)%rzm_t(1:imr,1:jmr,1:lmr,itr_loc:itr_loc), & status, & start=(/1,1,1,itr_file/), & count=(/imr,jmr,lmr,1/) ) IF_NOTOK_MDF(fid=ncid) #endif end if ! tracer active case ( 'levels' ) ! distributed over layers; might be empty on this pe if ( lmloc > 0 ) then call MDF_Get_Var( ncid, varid_rm, & mass_dat(n)%rm_k(1:imr,1:jmr,1:lmloc,itr:itr), & status, & start=(/1,1,offsetl+1,itr_file/), & count=(/imr,jmr,lmloc,1/) ) IF_NOTOK_MDF(fid=ncid) #ifdef slopes call MDF_Get_Var( ncid, varid_rxm, & mass_dat(n)%rxm_k(1:imr,1:jmr,1:lmloc,itr:itr), & status, & start=(/1,1,offsetl+1,itr_file/), & count=(/imr,jmr,lmloc,1/) ) IF_NOTOK_MDF(fid=ncid) call MDF_Get_Var( ncid, varid_rym, & mass_dat(n)%rym_k(1:imr,1:jmr,1:lmloc,itr:itr), & status, & start=(/1,1,offsetl+1,itr_file/), & count=(/imr,jmr,lmloc,1/) ) IF_NOTOK_MDF(fid=ncid) call MDF_Get_Var( ncid, varid_rzm, & mass_dat(n)%rzm_k(1:imr,1:jmr,1:lmloc,itr:itr), & status, & start=(/1,1,offsetl+1,itr_file/), & count=(/imr,jmr,lmloc,1/) ) IF_NOTOK_MDF(fid=ncid) #endif end if ! local layers > 0 case default write (gol,'("unsupported distribution:")'); call goPr write (gol,'(" n : ",i2)') n ;call goPr write (gol,'(" previous_par : ",a)') previous_par(n) ;call goPr write (gol,'("in ",a)') rname; call goErr; status=1; return end select else ! short lived tracer: ! check ... if ( itr_file <= ntracet_restart ) then write (gol,'("tracer `",a,"` is not-transported but seems to be transported in restart file")') trim(names(itr)); call goErr status=1 IF_NOTOK_MDF(fid=ncid) end if ! non-transported tracer; ! distributed over layers, any layers on this pe ? if ( lmloc > 0 ) then ! read non-transported tracer into: ! chem_dat%rm_k(:,:,1:lmloc,ntracet+1:ntracet+ntrace_chem) call MDF_Get_Var( ncid, varid_rmc, & chem_dat(n)%rm_k(1:imr,1:jmr,1:lmloc,itr:itr), & status, & start=(/1,1,offsetl+1,itr_file-ntracet_restart/), & count=(/imr,jmr,lmloc,1/) ) IF_NOTOK_MDF(fid=ncid) end if ! lmloc > 0 end if ! transported or short-lived end do ! tracers ! free mem for next region DeAllocate( values_names) end if ! do_rm #ifdef with_tendencies ! production/loss/concentration if ( do_plc ) then write (gol,'(" restore tendency data ...")'); call goPr ! set independent data mode (not all processes are reading): do itr = 1, plc_ntr do ipr = 1, plc_npr call MDF_Var_Par_Access( ncid, varid_plc(itr,ipr), MDF_INDEPENDENT, status ) IF_NOTOK_MDF(fid=ncid) end do end do ! distributed over levels: if ( lmloc > 0 ) then do itr = 1, plc_ntr do ipr = 1, plc_npr ! read plc slab: call MDF_Get_Var( ncid, varid_plc(itr,ipr), & plc_dat(n,itr,ipr)%rm_k(1:imr,1:jmr,1:lmloc), status, & start=(/1,1,offsetl+1/), count=(/imr,jmr,lmloc/) ) IF_NOTOK_MDF(fid=ncid) ! get time attribute: call MDF_Get_Att( ncid, varid_plc(itr,ipr), 'time', time6, status ) IF_NOTOK_MDF(fid=ncid) ! store time: call Set( plc_dat(n,itr,ipr)%t, time6=time6 ) end do end do end if end if ! do_plc #endif ! o close file call MDF_Close( ncid, status ) IF_NOTOK_RETURN(status=1) end do ! regions ! ok status = 0 end subroutine Restart_Read end module Restart