MODULE jules_netcdf ! Contains procedures used to read/write from netCDF files. USE inout, ONLY : & ! imported scalar parameters formatNc & ! imported scalars ,undefOut USE netcdf, ONLY : & ! imported procedures nf90_def_dim,nf90_def_var,nf90_enddef & ,nf90_float,nf90_get_var,nf90_global,nf90_inq_dimid & ,nf90_inq_varid,nf90_put_att,nf90_put_var & ! imported scalar parameters ,nf90_noErr,nf90_noWrite,nf90_unlimited USE rwErr_mod, ONLY : & ! imported procedures rwErr IMPLICIT NONE !------------------------------------------------------------------------------- ! Global scalar parameters. CHARACTER(len=7), PARAMETER :: ncTypeGswp2 = 'gswp2' ! value of ncType (qv) ! indicating GSWP2 format CHARACTER(len=7), PARAMETER :: ncTypeJules = 'jules' ! value of ncType (qv) ! indicating JULES format CHARACTER(len=7), PARAMETER :: ncTypePilps2e = 'pilps2e' ! value of ncType ! (qv) indicating PILPS2e format CHARACTER(len=7), PARAMETER :: ncTypePrincet = 'princet' ! value of ncType ! (qv) indicating Princeton format CHARACTER(len=7), PARAMETER :: ncTypeTseries = 'tseries' ! value of ncType ! (qv) indicating time series format. CHARACTER(len=7), PARAMETER :: ncTypeWatch = 'watch' ! value of ncType (qv) ! indicating WATCH format ! Global scalars. CHARACTER(len=7) :: ncType ! indicates what "type" of netCDF files are to be read ! Does not apply to driving (meteorological) data. ! Only used if netCDF files are to be read. ! This flag is used so that netCDF files from different experiments can be ! read by the same not-really-good netCDF-reading code. CHARACTER(len=7) :: ncTypeDrive ! as ncType, but for driving (meteorological) data CONTAINS !############################################################################### !############################################################################### !############################################################################### SUBROUTINE checkNcType( inType,errMess ) ! Simply checks that the declared "type" (i.e. format) of netCDf files is ! recognised, and a does a few other checks. USE inout, ONLY : & ! imported scalars with intent(in) nxIn,nyIn IMPLICIT NONE ! Scalar arguments with intent(in) CHARACTER(len=*), INTENT(in) :: inType ! indicated "type" of netCDf files CHARACTER(len=*), INTENT(in) :: errMess ! message to be printed on error LOGICAL :: errFound ! error flag !------------------------------------------------------------------------------- errFound = .FALSE. SELECT CASE ( inType ) CASE ( ncTypePilps2e,ncTypePrincet ) ! Recognised, and no restrictions. CASE ( ncTypeGswp2,ncTypeWatch ) ! Can only be used with vector input. IF ( nyIn /= 1 ) THEN errFound = .TRUE. WRITE(*,*)'ERROR: GSWP2- or WATCH-style netCDF files must be used with nyIn=1.' ENDIF CASE ( ncTypeTseries ) ! Can only be used with single point input. IF ( nxIn*nyIn /= 1 ) THEN errFound = .TRUE. WRITE(*,*)'ERROR: "Simple time series netCDF files must be used with nxIn=nyIn=1.' ENDIF CASE default errFound = .TRUE. WRITE(*,*)'ERROR: do not recognise ncType=',TRIM(inType) END SELECT IF ( errFound ) THEN WRITE(*,*)'ERROR in checkNcType.' WRITE(*,*)'ErrMess=',TRIM(errMess) WRITE(*,*)'Stopping in subroutine checkNcType.' STOP ENDIF END SUBROUTINE checkNcType !############################################################################### !############################################################################### ! subroutine createNcFile ! Internal procedure in module jules_netcdf. ! Set up a new, but already existing, netCDF file (create dimensions etc). SUBROUTINE createNcFile( iout,ncID,gridNx,gridNy,pointsOut,nvar & ,gridDx,gridDy,gridX1,gridY1,forGrads & ,haveSCpool,havePFT,haveSnow & ,haveSoil,haveTile,haveType,haveScalar,haveNOlevs & ,haveNfarray,haveSeed & ,outCompress,callType,compressGridFile & ,dateTimeString,fileName,globalAtt,compressMap & ,varName,varNcAtt,varType,varID ) USE ancil_info, ONLY : & ! imported scalars with intent(in) dim_cs1,nsmax,nsoil=>sm_levels,ntiles USE inout, ONLY : & ! imported arrays with intent(out) outTimeID USE nstypes, ONLY : & ! imported scalars with intent(in) npft,ntype USE time_loc, ONLY : & ! imported scalars with intent(in) timestep USE imogen_constants, ONLY : n_olevs,nfarray IMPLICIT NONE !------------------------------------------------------------------------------- ! Scalar arguments with intent(in) !------------------------------------------------------------------------------- INTEGER, INTENT(in) :: iout ! output profile number (only ! used for callType='diag') INTEGER, INTENT(in) :: ncID ! netCDF ID of file INTEGER, INTENT(in) :: gridNx ! x size of full output grid INTEGER, INTENT(in) :: gridNy ! y size of full output grid INTEGER, INTENT(in) :: pointsOut ! number of output points INTEGER, INTENT(in) :: nvar ! number of variables to be written REAL, INTENT(in) :: gridDx ! gridbox size in x (deg lon) REAL, INTENT(in) :: gridDy ! gridbox size in y (deg lat) REAL, INTENT(in) :: gridX1 ! x coord of point (1,1) in ! output grid (deg lon) REAL, INTENT(in) :: gridY1 ! y coord of point (1,1) in ! output grid (deg lat) LOGICAL, INTENT(in) :: forGrads ! T means netCDF files should be readable ! by GrADS ! ! F no requirement (butmight be readable) LOGICAL, INTENT(in) :: haveSCpool, havePFT,haveSnow,haveSoil,haveTile, & haveType,haveScalar,haveNOlevs,haveNfarray,haveSeed ! Flags indicating what "types" of variables are to be output. LOGICAL, INTENT(in) :: outCompress ! T if output is a compressed grid CHARACTER(len=*), INTENT(in) :: callType ! type of call ! This should be 'diag' or 'dump'. CHARACTER(len=*), INTENT(in) :: compressGridFile ! name of compression ! mapping data file (if used) (for information only) CHARACTER(len=*), INTENT(in) :: dateTimeString ! date and time as ! yyyy-mm-dd hh:mm:ss CHARACTER(len=*), INTENT(in) :: fileName ! name of file (for diagnostic ! purposes only) CHARACTER(len=*), INTENT(in) :: globalAtt ! a line to use as a global ! attribute (for a "note") !------------------------------------------------------------------------------- ! Array arguments with intent(in) !------------------------------------------------------------------------------- INTEGER, INTENT(in) :: compressMap(pointsOut) ! index in grid of output points CHARACTER(len=*), INTENT(in) :: varName(:) ! name of each output variable CHARACTER(len=*), INTENT(in) :: varNcAtt(:,:) ! attributes for each variable CHARACTER(len=*), INTENT(in) :: varType(:) ! type of each output variable !------------------------------------------------------------------------------- ! Array arguments with intent(out) !------------------------------------------------------------------------------- INTEGER, INTENT(out) :: varID(:) ! netCDF ID for each created variable !------------------------------------------------------------------------------- ! Local scalar parameters. !------------------------------------------------------------------------------- INTEGER, PARAMETER :: ndimMax = 10 ! number of possible "vertical" ! dimensions. These are pft,csPool,snow,soil,tile,type. INTEGER, PARAMETER :: nvarDimMax = 5 ! max possible number of netCDF ! for any one variable. These are x,y,z,pseudo and t. !------------------------------------------------------------------------------- ! Local scalar variables. !------------------------------------------------------------------------------- INTEGER :: idim,ierrLast ! loop counters/work INTEGER :: ierrSum ! work INTEGER :: ierr ! work INTEGER :: indexDimID ! netCDF ID for "index" dimension INTEGER :: indexID ! netCDF ID for "index" variable INTEGER :: idimPft ! index for PFT dimension INTEGER :: idimSCpool ! index for soil carbon pool dimension INTEGER :: idimSnow ! index for snow dimension INTEGER :: idimSoil ! index for soil dimension INTEGER :: idimTile ! index for tile dimension INTEGER :: idimType ! index for type dimension INTEGER :: idimScalar ! index for 'scalar' dimension INTEGER :: idimNOlevs ! index for n_olevs (see IMOGEN) dimension INTEGER :: idimNfarray ! index for nfarray (see IMOGEN) dimension INTEGER :: idimSeed ! index for seed dimension (see IMOGEN) INTEGER :: ivar,ix,iy,iz ! loop counters/work INTEGER :: kvar ! work INTEGER :: latID ! netCDF ID for latitude variable INTEGER :: lonID ! netCDF ID for longitude variable INTEGER :: ndim ! # of netCDF dimensions for a variable INTEGER :: nvarCreate ! number of variables created INTEGER :: timeDimID ! netCDF ID for time dimension INTEGER :: timeID ! netCDF ID for a time variable INTEGER :: vID ! netCDF ID for a variable INTEGER :: xDimID ! netCDF ID for x dimension INTEGER :: yDimID ! netCDF ID for y dimension CHARACTER(len=100) :: useName ! work CHARACTER(len=100) :: cwork ! workspace - hopefully large enough! !------------------------------------------------------------------------------- ! Local array variables. !------------------------------------------------------------------------------- INTEGER :: dimID(ndimMax) ! netCDF ID of each possible "levels" dimension INTEGER :: dimSize(ndimMax) ! size of each possible "levels" dimension INTEGER :: dimVarID(ndimMax) ! ID for variables related to "levels" dimensions INTEGER :: varDims(nvarDimMax) ! netCDF IDs for dimensions for variables REAL :: lat(gridNx,gridNy) ! latitudes of points in full output grid REAL :: lon(gridNx,gridNy) ! longitudes of points in full output grid LOGICAL :: useDim(ndimMax) ! flag indicating if a dimension used CHARACTER(len=10) :: dimName(ndimMax) ! name of each possible "levels" dimension CHARACTER(len=30) :: dimLongName(ndimMax)! name of each possible "levels" dimension !------------------------------------------------------------------------------- ! Of course we should probably try to follow the CF conventions...but I'm in a ! hurry. Instead I'm basing this on GSWP2 and WATCH Forcing Data files. !------------------------------------------------------------------------------- ! Set up all possible "levels" dimensions. idim = 0 idim=idim+1; idimPFT=idim; dimName(idim)='pft' dimLongName(idim) = 'plant functional types' dimSize(idim) = npft idim=idim+1; idimSCpool=idim; dimName(idim)='SCpool' dimLongName(idim) = 'soil carbon pools' dimSize(idim) = dim_cs1 idim=idim+1; idimSnow=idim; dimName(idim)='snow' dimLongName(idim) = 'snow layers' dimSize(idim) = nsmax idim=idim+1; idimSoil=idim; dimName(idim)='soil' dimLongName(idim) = 'soil layers' dimSize(idim) = nsoil idim=idim+1; idimTile=idim; dimName(idim)='tile' dimLongName(idim) = 'surface tiles' dimSize(idim) = ntiles idim=idim+1; idimType=idim; dimName(idim)='type' dimLongName(idim) = 'surface types' dimSize(idim) = ntype idim=idim+1; idimScalar=idim; dimName(idim)='scalar' dimLongName(idim) = 'fake dimension to allow scalars in dump code' dimSize(idim) = 1 idim=idim+1; idimNOlevs=idim; dimName(idim)='n_olevs' dimLongName(idim) = 'number of IMOGEN ocean levels' dimSize(idim) = n_olevs idim=idim+1; idimNfarray=idim; dimName(idim)='nfarray' dimLongName(idim) = 'nfarray' dimSize(idim) = nfarray idim=idim+1; idimSeed=idim; dimName(idim)='seed' dimLongName(idim) = 'number of components in IMOGEN seed variables' dimSize(idim) = 4 !------------------------------------------------------------------------------- ! Calculate lat and lon for every gridpoint (in the full, uncompressed output grid). !------------------------------------------------------------------------------- DO iy=1,gridNy DO ix=1,gridNx lon(ix,iy) = gridX1 + REAL(ix-1) * gridDx lat(ix,iy) = gridY1 + REAL(iy-1) * gridDy ENDDO ENDDO !------------------------------------------------------------------------------- ! Create dimensions. ! Allow a few errors to accumulate before stopping! !------------------------------------------------------------------------------- ierrSum = 0 ierrLast = 0 IF ( callType == 'diag' ) THEN ! Unlimited dimension. ierr = nf90_def_dim( ncID,'tstep',nf90_unlimited,timeDimID ) ELSE ierr = nf90_def_dim( ncID,'tstep',1,timeDimID ) ENDIF ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! x dimension. ierr = nf90_def_dim( ncID,'x',gridNx,xDimID ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! y dimension. ierr = nf90_def_dim( ncID,'y',gridNy,yDimID ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Index dimension (e.g. land points) - only used if grid is compressed. IF ( outCompress ) THEN ierr = nf90_def_dim( ncID,'index',pointsOut,indexDimID ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ENDIF !------------------------------------------------------------------------------- ! Create vertical/levels/pseudo dimensions. ! Loop over all possible vertical dimensions. ! This isn't great but works. Keeping everything in arrays reduces repetition, ! but still requires an "IF haveXX" test. !------------------------------------------------------------------------------- useDim(:) = .FALSE. DO idim=1,ndimMax IF ( idim==idimPFT .AND. havePFT ) THEN useDim(idim) = .TRUE. ELSEIF ( idim==idimSCpool .AND. haveSCpool ) THEN useDim(idim) = .TRUE. ELSEIF ( idim==idimSnow .AND. haveSnow ) THEN useDim(idim) = .TRUE. ELSEIF ( idim==idimSoil .AND. haveSoil ) THEN useDim(idim) = .TRUE. ELSEIF ( idim==idimTile .AND. haveTile ) THEN useDim(idim) = .TRUE. ELSEIF ( idim==idimType .AND. haveType ) THEN useDim(idim) = .TRUE. ELSEIF ( idim==idimScalar .AND. haveScalar ) THEN useDim(idim) = .TRUE. ELSEIF ( idim==idimNOlevs .AND. haveNOlevs ) THEN useDim(idim) = .TRUE. ELSEIF ( idim==idimNfarray .AND. haveNfarray ) THEN useDim(idim) = .TRUE. ELSEIF ( idim==idimSeed .AND. haveSeed ) THEN useDim(idim) = .TRUE. ENDIF IF ( useDim(idim) ) THEN ! Create dimension. ierr = nf90_def_dim( ncID,TRIM(dimName(idim)),dimSize(idim),dimID(idim) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ENDIF ENDDO ! Snow layer variables also require a tile dimension, unless ! writing for GrADS. useDim=T means it has already been defined. IF ( haveSnow .AND. .NOT.forGrads .AND. .NOT.useDim(idimTile) ) THEN ! Create tile dimension. idim = idimTile ierr = nf90_def_dim( ncID,TRIM(dimName(idim)),dimSize(idim),dimID(idim) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ENDIF ! If error, call error routine, indicating a stop. IF ( ierrSum /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierrLast,fileName=fileName,ncID=ncid & ,errMess1='createNcFile nf90_def_dim - last error reported' ) !------------------------------------------------------------------------------- ! Create variables to describe the grid. !------------------------------------------------------------------------------- !------------------------------------------------------------------------------- ! Create time variable. !------------------------------------------------------------------------------- ierr = nf90_def_var( ncID,'time',nf90_float,(/timeDimID/),timeID ) ! Save ID for variable. IF ( callType == 'diag' ) outTimeID(iout,1) = timeID IF ( ierr /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,fileName=fileName,ncID=ncid,errMess1='createNcFile nf90_def_var time' ) ! Write attributes. ! Units. IF ( callType == 'diag' ) THEN cwork = 'seconds since ' // TRIM(dateTimeString) ELSE cwork = TRIM(dateTimeString) ENDIF ierr = nf90_put_att( ncID,timeID,'units',TRIM(cwork) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Title ierr = nf90_put_att( ncID,timeID,'title','Time' ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Long name. ierr = nf90_put_att( ncID,timeID,'long_name','Time axis' ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Time origin. ierr = nf90_put_att( ncID,timeID,'time_origin',TRIM(dateTimeString) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! If error, call error routine, indicating a stop. IF ( ierrSum /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierrLast,fileName=fileName,ncID=ncid & ,errMess1='createNcFile nf90_put_att, last error for var= time' ) !------------------------------------------------------------------------------- ! Create time step variable. !------------------------------------------------------------------------------- ! Get ID for variable. ierr = nf90_def_var( ncID,'timestp ',nf90_float,(/timeDimID/),timeID ) ! Save ID for variable. IF ( callType == 'diag' ) outTimeID(iout,2) = timeID IF ( ierr /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,fileName=fileName,ncID=ncid & ,errMess1='createNcFile nf90_def_var timestp' ) ! Write attributes. ! Units. cwork = 'timesteps starting ' // TRIM(dateTimeString) ierr = nf90_put_att( ncID,timeID,'units',TRIM(cwork) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Title ierr = nf90_put_att( ncID,timeID,'title','Time steps' ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Long name. ierr = nf90_put_att( ncID,timeID,'long_name','Time step axis' ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Time origin. ierr = nf90_put_att( ncID,timeID,'time_origin',TRIM(dateTimeString) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Timestep length. ierr = nf90_put_att( ncID,timeID,'tstep_sec',timestep ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! If error, call error routine, indicating a stop. IF ( ierrSum /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierrLast,fileName=fileName,ncID=ncid & ,errMess1='createNcFile nf90_put_att, last error for var= timestp' ) !------------------------------------------------------------------------------- ! Create longitude variable. ! This is longitude(x,y). !------------------------------------------------------------------------------- ! Create variable. ierr = nf90_def_var( ncID,'lon',nf90_float,(/xDimID,yDimID/),lonID ) IF ( ierr /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,fileName=fileName,ncID=ncid,errMess1='createNcFile nf90_def_var lon' ) ! Write attributes. ! Units ierr = nf90_put_att( ncID,lonID,'units','pixel centre degrees east' ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Min ierr = nf90_put_att( ncID,lonID,'valid_min',lon(1,1) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Max ierr = nf90_put_att( ncID,lonID,'valid_max',lon(gridNx,1) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Long name ierr = nf90_put_att( ncID,lonID,'long_name','Longitude' ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! If error, call error routine, indicating a stop. IF ( ierrSum /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierrLast,fileName=fileName,ncID=ncid & ,errMess1='createNcFile nf90_put_att, last error for var= lon' ) !------------------------------------------------------------------------------- ! Create latitude variable. ! This is lat(x,y). !------------------------------------------------------------------------------- ! Create variable. ierr = nf90_def_var( ncID,'lat',nf90_float,(/xDimID,yDimID/),latID ) IF ( ierr /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,fileName=fileName,ncID=ncid,errMess1='createNcFile nf90_def_var lat' ) ! Write attributes. ! Units ierr = nf90_put_att( ncID,latID,'units','pixel centre degrees north' ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Min ierr = nf90_put_att( ncID,latID,'valid_min',lat(1,1) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Max ierr = nf90_put_att( ncID,latID,'valid_max',lat(1,gridNy) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Long name ierr = nf90_put_att( ncID,latID,'long_name','Latitude' ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! If error, call error routine, indicating a stop. IF ( ierrSum /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierrLast,fileName=fileName,ncID=ncid & ,errMess1='createNcFile nf90_put_att, last error for var=latitude' ) IF ( outCompress ) THEN !------------------------------------------------------------------------------- ! Grid index (point number). Only needed for a compressed grid (and not used ! by GrADS at present, which uses a pdef file). ! Really integer, but leave as real for now. !------------------------------------------------------------------------------- ! Create index variable. ierr = nf90_def_var( ncID,'index',nf90_float,(/indexDimID/),indexID ) IF ( ierr /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,fileName=fileName,ncID=ncid,errMess1='createNcFile nf90_def_var index' ) ! Write attributes. ! Compress ierr = nf90_put_att( ncID,indexID,'compress','y x' ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Notes ierr = nf90_put_att( ncID,indexID,'note','index in grid' ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! If error, call error routine, indicating a stop. IF ( ierrSum /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierrLast,fileName=fileName,ncID=ncid & ,errMess1='createNcFile nf90_put_att, last error for var=index' ) ENDIF ! outCompress !------------------------------------------------------------------------------- ! Loop over all possible "layers" dimensions, creating a variable if required. ! Currently using an index - really integer, but leave as real for now. !------------------------------------------------------------------------------- DO idim=1,ndimMax IF ( useDim(idim) ) THEN ! Create variable. ierr = nf90_def_var( ncID,TRIM(dimName(idim)),nf90_float,(/dimID(idim)/),dimVarID(idim) ) IF ( ierr /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,fileName=fileName,ncID=ncid & ,errMess1='createNcFile nf90_def_var var='//TRIM(dimName(idim)) ) ! Write attributes for variable. ! Title ierr = nf90_put_att( ncID,dimVarID(idim),'title',TRIM(dimName(idim)) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Long name. ierr = nf90_put_att( ncID,dimVarID(idim),'long_name',TRIM(dimLongName(idim)) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! If error, call error routine, indicating a stop. IF ( ierrSum /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierrLast,fileName=fileName,ncID=ncid & ,errMess1='createNcFile nf90_put_att, last error for var='//TRIM(dimName(idim)) ) ENDIF ! useDim ENDDO ! idim !------------------------------------------------------------------------------- !------------------------------------------------------------------------------- ! Create variables and attributes for the "data" variables. !------------------------------------------------------------------------------- !------------------------------------------------------------------------------- DO ivar=1,nvar nvarCreate = 1 ! by default, each variable generates one netCDF var !------------------------------------------------------------------------------- ! Decide on dimensions for this variable. !------------------------------------------------------------------------------- SELECT CASE ( varType(ivar) ) CASE ( 'LA', 'RG', 'SI' ) IF ( outCompress ) THEN ndim = 2 ! index,t varDims(1:ndim) = (/ indexDimID,timeDimID /) ELSE ndim = 3 ! x,y,t varDims(1:ndim) = (/ xDimID,ydimID,timeDimID /) ENDIF CASE ( 'PF' ) IF ( outCompress ) THEN ndim = 3 ! index,pft,t varDims(1:ndim) = (/ indexDimID,dimID(idimPFT),timeDimID /) ELSE ndim = 4 ! x,y,pseudo,t varDims(1:ndim) = (/ xDimID,ydimID,dimID(idimPFT),timeDimID /) ENDIF CASE ( 'SC' ) IF ( outCompress ) THEN ndim = 3 ! index,pool,t varDims(1:ndim) = (/ indexDimID,dimID(idimSCpool),timeDimID /) ELSE ndim = 4 ! x,y,pool,t varDims(1:ndim) = (/ xDimID,ydimID,dimID(idimSCpool),timeDimID /) ENDIF CASE ( 'SO' ) IF ( outCompress ) THEN ndim = 3 ! index,z,t varDims(1:ndim) = (/ indexDimID,dimID(idimSoil),timeDimID /) ELSE ndim = 4 ! x,y,z,t varDims(1:ndim) = (/ xDimID,ydimID,dimID(idimSoil),timeDimID /) ENDIF CASE ( 'SN' ) ! Snow layer variables are awkward! IF ( .NOT. forGrads ) THEN ! Keep layer and pseudo dimensions. IF ( outCompress ) THEN ndim = 4 ! index,z,tile,t varDims(1:ndim) = (/ indexDimID,dimID(idimSnow) & ,dimID(idimTile),timeDimID /) ELSE ndim = 5 ! x,y,z,tile,t varDims(1:ndim) = (/ xDimID,ydimID,dimID(idimSnow) & ,dimID(idimTile),timeDimID /) ENDIF ELSE ! Each tile (pseudo dimension) is represented as a separate variable. nvarCreate = ntiles IF ( outCompress ) THEN ndim = 3 ! index,z,t varDims(1:ndim) = (/ indexDimID,dimID(idimSnow),timeDimID /) ELSE ndim = 4 ! x,y,z,t varDims(1:ndim) = (/ xDimID,ydimID,dimID(idimSnow),timeDimID /) ENDIF ENDIF ! forGrads CASE ( 'TI' ) IF ( outCompress ) THEN ndim = 3 ! index,tile,t varDims(1:ndim) = (/ indexDimID,dimID(idimTile),timeDimID /) ELSE ndim = 4 ! x,y,tile,t varDims(1:ndim) = (/ xDimID,ydimID,dimID(idimTile),timeDimID /) ENDIF CASE ( 'TY' ) IF ( outCompress ) THEN ndim = 3 ! index,type,t varDims(1:ndim) = (/ indexDimID,dimID(idimType),timeDimID /) ELSE ndim = 4 ! x,y,type,t varDims(1:ndim) = (/ xDimID,ydimID,dimID(idimType),timeDimID /) ENDIF !------------------------------------------------------------------------------- ! Compression is not enabled for IMOGEN variables since they are on a ! different grid !------------------------------------------------------------------------------- CASE ( 'SL' ) ndim = 2 ! scalar,t varDims(1:ndim) = (/ dimID(idimScalar),timeDimID /) CASE ( 'OC' ) ndim = 2 ! n_olevs,t varDims(1:ndim) = (/ dimID(idimNOlevs),timeDimID /) CASE ( 'FA' ) ndim = 2 ! nfarray,t varDims(1:ndim) = (/ dimID(idimNfarray),timeDimID /) CASE ( 'SD' ) ndim = 2 ! nfarray,t varDims(1:ndim) = (/ dimID(idimSeed),timeDimID /) CASE default WRITE(*,*)'ERROR: createNcFile: ndim: varType=',varType(ivar) WRITE(*,*)'No code for this varType.' STOP END SELECT !------------------------------------------------------------------------------- ! Create one or more netCDF variables. ! This is a painful way to get snow layer variables that GrADS can read! ! Otherwise we don't need the kvar loop. !------------------------------------------------------------------------------- DO kvar=1,nvarCreate IF ( nvarCreate>1 .AND. varType(ivar)/='SN' ) THEN WRITE(*,*)'ERROR: createNcFile: nvarCreate>1 .AND. varType(ivar)/=SN' STOP ENDIF !------------------------------------------------------------------------------- ! Get ID for each created variable. !------------------------------------------------------------------------------- useName = TRIM(varName(ivar)) IF ( nvarCreate>1 .AND. varType(ivar)=='SN' ) THEN ! For snow layer variables, add tile number to name. IF ( nvarCreate < 10 ) THEN WRITE(cwork,"(i1.1)") kvar ELSE WRITE(cwork,"(i2.2)") kvar ENDIF useName = TRIM(useName) // TRIM(cwork) ENDIF ierr = nf90_def_var( ncID,TRIM(useName),nf90_float & ,varDims(1:ndim),vID ) ! A hack for snow layers vars with gradsNC=T: only save the ! var ID for first level. IF ( kvar == 1 ) varID(ivar) = vID IF ( ierr /= nf90_noErr ) THEN WRITE(*,*)'Error for variable=',TRIM(useName) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,fileName=fileName,ncID=ncid & ,errMess1='createNcFile nf90_def_var' ) ENDIF !------------------------------------------------------------------------------- ! Write attributes for variable. !------------------------------------------------------------------------------- ierr = nf90_put_att( ncID,vID,'axis',TRIM(varNcAtt(ivar,1)) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ierr = nf90_put_att( ncID,vID,'units',TRIM(varNcAtt(ivar,2)) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ierr = nf90_put_att( ncID,vID,'missing_value',undefOut ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ierr = nf90_put_att( ncID,vID,'title',TRIM(useName) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! Here's another painful consequence of hacking snow variables for GrADS - ! varNcAtt has an entry for each var as shown in ctl file. ierr = nf90_put_att( ncID,vID,'long_name',TRIM(varNcAtt(ivar+kvar-1,4)) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! If error, call error routine, indicating a stop. IF ( ierrSum /= 0 ) THEN WRITE(*,*)'Error for variable=',TRIM(varName(ivar)) CALL rwErr( formatNc,.TRUE. & ,ierr=ierrLast,fileName=fileName,ncID=ncid & ,errMess1='createNcFile nf90_put_att, last error' ) ENDIF ENDDO ! kvar created variables ENDDO ! ivar !------------------------------------------------------------------------------- ! Global attributes. !------------------------------------------------------------------------------- ierr = nf90_put_att( ncID,nf90_global,'title','Output from the JULES model.' ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ierr = nf90_put_att( ncID,nf90_global,'note',TRIM(globalAtt) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr IF ( outCompress ) THEN cwork = 'Mapping information held in ' // TRIM(compressGridFile) ierr = nf90_put_att( ncID,nf90_global,'note1',TRIM(cwork) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ENDIF !------------------------------------------------------------------------------- ! Finished definitions, leave define mode. !------------------------------------------------------------------------------- ierr = nf90_enddef( ncID ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! If error, call error routine, indicating a stop. IF ( ierrSum /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierrLast,fileName=fileName,ncID=ncid & ,errMess1='createNcFile last error') !------------------------------------------------------------------------------- ! Write grid variables. !------------------------------------------------------------------------------- ! Write lat and lon variables. These are for the full (uncompressed) output grid. ierr = nf90_put_var( ncID,lonID,lon(:,:) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ierr = nf90_put_var( ncID,latID,lat(:,:) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! If error, call error routine, indicating a stop. IF ( ierrSum /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierrLast,fileName=fileName,ncID=ncid & ,errMess1='createNcFile last error from write lat/lon' ) ! Write index variable. IF ( outCompress ) THEN ierr = nf90_put_var( ncID,indexID,REAL( compressMap ) ) ierrSum = ierrSum + ierr IF ( ierr /= nf90_noErr ) ierrLast = ierr ! If error, call error routine, indicating a stop. IF ( ierrSum /= 0 ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierrLast,fileName=fileName,ncID=ncid & ,errMess1='createNcFile error from write indexVar' ) ENDIF !------------------------------------------------------------------------------- ! Write layer variables. ! Loop over all possible "layers" dimensions. ! Currently using an index - really integer, but leave as real for now. !------------------------------------------------------------------------------- DO idim=1,ndimMax IF ( useDim(idim) ) THEN ierr = nf90_put_var( ncID,dimVarID(idim),(/ (iz,iz=1,dimSize(idim)) /) ) ! If error, call error routine, indicating a stop. IF ( ierr /= nf90_noErr ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,fileName=fileName,ncID=ncid & ,errMess1='createNcFile error from write levels for var='//TRIM(dimName(idim)) ) ENDIF ENDDO END SUBROUTINE createNcFile !############################################################################### !############################################################################### !############################################################################### ! subroutine readVar2dReal_ncVector_gswp2 ! Internal procedure in module jules_netcdf. ! Read a pseudo-2D field from netCDF file - a temporary (hopefully) fix for ! GSWP2 and WATCH files. In fact, reads a vector (eg compressed to land points). ! I have combined these GSWP2 and WATCH cases because they are so similar, ! but it does make this routine slightly more complicated. SUBROUTINE readVar2dReal_ncVector_gswp2( ncCallType,ncTypeArg,ncID,varName & ,t,z,npoints,outval,start ) IMPLICIT NONE !------------------------------------------------------------------------------- ! Scalar arguments with intent(in). !------------------------------------------------------------------------------- INTEGER, INTENT(in) :: & ncID &! the ID (number) of the netCDF file ,npoints &! number of points to read ,t &! time number to read ,z ! level number to read CHARACTER(len=*), INTENT(in) :: & ncCallType &! type (location) of call ,ncTypeArg &! type (format) of netCDF files ,varName ! the name of the netCDF variable !------------------------------------------------------------------------------- ! Array arguments with intent(out) !------------------------------------------------------------------------------- REAL, INTENT(out) :: outval(npoints) ! output data !------------------------------------------------------------------------------- ! Optional scalar arguments with intent(in). !------------------------------------------------------------------------------- INTEGER, INTENT(in), OPTIONAL :: start ! starting point number !------------------------------------------------------------------------------- ! Local scalar variables. !------------------------------------------------------------------------------- INTEGER :: & ierr &! error flag ,varID &! the ID (number) of the netCDF variable ,xstart ! starting point number LOGICAL :: errFound ! error flag CHARACTER(len=len(ncCallType)+6) :: callType ! local version of ! ncCallType, used to select code where GSWP2 and WATCH differ. ! Length is sufficient for ncCallType + _gswp2 or _watch. !------------------------------------------------------------------------------- ! Look for optional arguments. xstart = 1 IF ( PRESENT(start) ) xstart=start ! Get the ID of the variable. ierr = nf90_inq_varid( ncID,varName,varID ) ! If error, call error routine, indicating a stop. IF( ierr /= nf90_noErr ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,ncID=ncID,varName=varName & ,errMess1='readVar2dReal_ncVector_gswp2 nf90_inq_varid' & ,errMess2=TRIM(ncCallType) ) !------------------------------------------------------------------------------- ! Create a new variable that can be used to select code for GSWP2 or WATCH ! "formats" when these need different code. !------------------------------------------------------------------------------- errFound = .FALSE. callType = ncCallType ! default SELECT CASE ( ncCallType ) ! Cases for which identical code can be used for GSWP2 and WATCH. ! Nothing more to do. CASE ( 'init_agric','init_grid','init_hgt','init_pdm','init_soil' & ,'init_top','veg_read' ) ! Cases that are only considered for GSWP2. CASE ( 'drive_read','init_ic' ) if ( ncTypeArg /= ncTypeGSWP2 ) errFound = .TRUE. ! Cases that require different code for GSWP2 and WATCH. CASE ( 'init_frac' ) if ( ncTypeArg == ncTypeGSWP2 ) then callType = trim(ncCallType) // '_gswp2' elseif ( ncTypeArg == ncTypeWatch ) then callType = trim(ncCallType) // '_watch' else errFound = .TRUE. endif CASE default errFound = .TRUE. END SELECT ! ncCallType !------------------------------------------------------------------------------- ! Stop if there is no code for this case. !------------------------------------------------------------------------------- if ( errFound ) then write(*,*)'ERROR: readVar2dReal_ncVector_gswp2: no code for these options.' write(*,*)'ncCallType=',trim(ncCallType) write(*,*)'ncTypeArg=',trim(ncTypeArg) stop endif !------------------------------------------------------------------------------- ! Read data. !------------------------------------------------------------------------------- SELECT CASE ( callType ) !------------------------------------------------------------------------------- ! First, cases other than 'init_ic'. !------------------------------------------------------------------------------- ! 1-D cases. CASE ( 'init_agric','init_grid','init_pdm','init_soil','init_top' ) ! Read a vector with single dimension ('Land' or 'land'). ierr = nf90_get_var( ncID,varID,outval,(/xstart/) ) ! 2-D cases. CASE ( 'drive_read' ) ! Read a vector from 1st dimension ('land'), ! for given value of 2nd dimension ('tstep'). ierr = nf90_get_var( ncID,varID,outval,(/xstart,t/),(/npoints,1/) ) CASE ( 'init_frac_watch','init_hgt' ) ! Read a vector from 1st dimension ('Land' or 'land'), ! for given value of 2nd dimension ('Psuedo' or 'pseudo'). ierr = nf90_get_var( ncID,varID,outval,(/xstart,z/),(/npoints,1/) ) ! 3-D cases. CASE ( 'init_frac_gswp2','veg_read' ) ! Read a vector from 1st dimension ('Land' or 'land'), ! for given value of 2nd dimension ('Psuedo' or 'pseudo') ! and 3rd dimension ('Time' or 'tstep'). ierr = nf90_get_var( ncID,varID,outval,(/xstart,z,t/),(/npoints,1,1/) ) !------------------------------------------------------------------------------- ! Now, init_ic cases, selected on variable name. !------------------------------------------------------------------------------- CASE ( 'init_ic' ) SELECT CASE ( varName ) CASE ( 'cs', 'gs' ) ! Read a vector with single dimension ('Land'). ierr = nf90_get_var( ncID,varID,outval,(/xstart/) ) CASE ( 'sthuf', 't_soil' ) ! Read a vector from 1st dimension ('Land'), ! for given value of 2nd dimension ('Soil'). ierr = nf90_get_var( ncID,varID,outval,(/xstart,z/),(/npoints,1/) ) CASE ( 'canopy', 'frac', 'rgrain', 'snow_tile', 'tstar_tile' ) ! Read a vector from 1st dimension ('Land'), ! for given value of 2nd dimension ('Psuedo'). ierr = nf90_get_var( ncID,varID,outval,(/xstart,z/),(/npoints,1/) ) CASE default WRITE(*,*)'ERROR: readVar2dReal_ncVector_gswp2: no code!' WRITE(*,*)'callType=',TRIM(callType),' varName=',TRIM(varName) STOP END SELECT !------------------------------------------------------------------------------- CASE default WRITE(*,*)'ERROR: readVar2dReal_ncVector_gswp2: no code for callType.' WRITE(*,*)'callType=',TRIM(callType) WRITE(*,*)'ncCallType=',TRIM(ncCallType) WRITE(*,*)'ncTypeArg=',TRIM(ncTypeArg) STOP END SELECT ! callType !------------------------------------------------------------------------------- ! If error, call error routine, indicating a stop. !------------------------------------------------------------------------------- IF ( ierr /= nf90_noErr ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,ncID=ncID,varName=varName,varID=varID & ,t=t,xstart=xstart,z=z & ,errMess1='readVar2dReal_ncVector_gswp2 nf90_get_var' & ,errMess2=TRIM(callType) ) END SUBROUTINE readVar2dReal_ncVector_gswp2 !############################################################################### !############################################################################### !############################################################################### ! subroutine readVar2dReal_nc_pilps2e ! Internal procedure in module jules_netcdf. ! Read a 2D field from netCDF file - a temporary (hopefully) fix for PILPS2e ! data. ! SUBROUTINE readVar2dReal_nc_pilps2e( ncCallType,ncID,varName & ,nx,xstart,ny,ystart,t,outval ) IMPLICIT NONE INTEGER, INTENT(in) :: &! in scalars ncID &! the ID (number) of the netCDF file ,nx &! size in x of data range to be read ,ny &! size in y of data range to be read ,t &! time number to read ,xstart &! first x location to read ,ystart ! first x location to read INTEGER :: &! local scalars ierr &! error flag ,varID ! the ID (number) of the netCDF variable REAL, INTENT(out) :: &! out arrays outval(nx*ny) CHARACTER(len=*), INTENT(in) :: &! in scalars ncCallType &! type of call ,varName ! the name of the netCDF variable !------------------------------------------------------------------------------- ! Get the ID of the variable. ierr = nf90_inq_varid( ncID,varName,varID ) ! If error, call error routine, indicating a stop. IF( ierr /= nf90_noErr ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,ncID=ncID,varName=varName & ,errMess1='readVar2dReal_nc_pilps2e nf90_inq_varid' & ,errMess2=TRIM(ncCallType) ) SELECT CASE ( ncCallType ) !------------------------------------------------------------------------------- ! Only driving (meteorological) data can be read. CASE ( 'drive_read' ) ! Read a grid of nx*ny points from dimensions #1 and 2 ('x' and 'y'), for given ! value of 3rd dimension ('tstep'). ierr = nf90_get_var( ncID,varID,outval,(/xstart,ystart,t/),(/nx,ny,1/) ) !-------------------------------------------------- CASE default WRITE(*,*)'ERROR: readVar2dReal_nc_pilps2e: no code for ncCallType.' WRITE(*,*)'ncCallType=',TRIM(ncCallType) STOP END SELECT ! If error, call error routine, indicating a stop. IF ( ierr /= nf90_noErr ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,ncID=ncID,varName=varName,varID=varID,t=t,xstart=xstart & ,ystart=ystart,errMess1='readVar2dReal_nc_pilps2e nf90_get_var' & ,errMess2=TRIM(ncCallType) ) END SUBROUTINE readVar2dReal_nc_pilps2e !############################################################################### !############################################################################### !############################################################################### ! subroutine readVar2dReal_nc_princet ! Internal procedure in module jules_netcdf. ! Read a 2D field from netCDF file - a temporary (hopefully) fix for Princeton ! data (Sheffield et al., 2006, J. Climate.) ! data. ! SUBROUTINE readVar2dReal_nc_princet( ncCallType,ncID,varName & ,nx,xstart,ny,ystart,t,outval ) IMPLICIT NONE ! Local scalar parameters. INTEGER, PARAMETER :: z = 1 ! index to use in "z" dimension INTEGER, INTENT(in) :: &! in scalars ncID &! the ID (number) of the netCDF file ,nx &! size in x of data range to be read ,ny &! size in y of data range to be read ,t &! time number to read ,xstart &! first x location to read ,ystart ! first x location to read INTEGER :: &! local scalars ierr &! error flag ,varID ! the ID (number) of the netCDF variable REAL, INTENT(out) :: &! out arrays outval(nx*ny) CHARACTER(len=*), INTENT(in) :: &! in scalars ncCallType &! type of call ,varName ! the name of the netCDF variable !------------------------------------------------------------------------------- ! Get the ID of the variable. ierr = nf90_inq_varid( ncID,varName,varID ) ! If error, call error routine, indicating a stop. IF( ierr /= nf90_noErr ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,ncID=ncID,varName=varName & ,errMess1='readVar2dReal_nc_princet nf90_inq_varid' & ,errMess2=TRIM(ncCallType) ) SELECT CASE ( ncCallType ) !------------------------------------------------------------------------------- ! Only driving (meteorological) data can be read. CASE ( 'drive_read' ) ! Read a grid of nx*ny points from dimensions #1 and 2 ('longitude' and ! 'latitude'), for given values of 3rd and 4th dimensions ('z','time'). ierr = nf90_get_var( ncID,varID,outval,(/xstart,ystart,z,t/),(/nx,ny,1,1/) ) !-------------------------------------------------- CASE default WRITE(*,*)'ERROR: readVar2dReal_nc_princet: no code for ncCallType.' WRITE(*,*)'ncCallType=',TRIM(ncCallType) STOP END SELECT ! If error, call error routine, indicating a stop. IF ( ierr /= nf90_noErr ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,ncID=ncID,varName=varName,varID=varID,t=t,xstart=xstart & ,ystart=ystart,z=z,errMess1='readVar2dReal_nc_princet nf90_get_var' & ,errMess2=TRIM(ncCallType) ) END SUBROUTINE readVar2dReal_nc_princet !############################################################################### !############################################################################### !############################################################################### ! subroutine readVar2dReal_nc_tseries ! Internal procedure in module jules_netcdf. ! Read a 2D (to JULES) field from a netCDF file that contains a single space point. ! SUBROUTINE readVar2dReal_nc_tseries( ncCallType,ncID,varName,t & ,outval ) IMPLICIT NONE INTEGER, INTENT(in) :: &! in scalars ncID &! the ID (number) of the netCDF file ,t ! time number to read INTEGER :: &! local scalars ierr &! error flag ,varID ! the ID (number) of the netCDF variable REAL, INTENT(out) :: &! out arrays outval(1) CHARACTER(len=*), INTENT(in) :: &! in scalars ncCallType &! type of call ,varName ! the name of the netCDF variable !------------------------------------------------------------------------------- ! Get the ID of the variable. ierr = nf90_inq_varid( ncID,varName,varID ) ! If error, call error routine, indicating a stop. IF( ierr /= nf90_noErr ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,ncID=ncID,varName=varName & ,errMess1='readVar2dReal_nc_tseries nf90_inq_varid' & ,errMess2=TRIM(ncCallType) ) SELECT CASE ( ncCallType ) !------------------------------------------------------------------------------- ! Only driving (meteorological) data can be read. CASE ( 'drive_read' ) ! Read a single value from dimension #1 ( time ). ierr = nf90_get_var( ncID,varID,outval,(/t/),(/1/) ) !-------------------------------------------------- CASE default WRITE(*,*)'ERROR: readVar2dReal_nc_tseries: no code for ncCallType.' WRITE(*,*)'ncCallType=',TRIM(ncCallType) STOP END SELECT ! If error, call error routine, indicating a stop. IF ( ierr /= nf90_noErr ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,ncID=ncID,varName=varName,varID=varID,t=t & ,errMess1='readVar2dReal_nc_tseries nf90_get_var' & ,errMess2=TRIM(ncCallType) ) END SUBROUTINE readVar2dReal_nc_tseries !############################################################################### !############################################################################### ! subroutine check_nc_dims ! Internal procedure in module jules_netcdf. ! Checks that the dimensions found in an existing netCDF file match those ! expected. ! Note that this is far from being a full check - it will establish if the ! dimensions exist at all in the file, but not if the required variables have ! these dimensions. This is bit of a faff, and could alternatively be dealt with ! by adding to error diagnostics if an error occurs on read. SUBROUTINE check_nc_dims( ncID,ncCallType,ncFileType ) IMPLICIT NONE ! Local scalar parameters. INTEGER, PARAMETER :: ndimMax = 4 ! maximum possible number of dimensions ! Scalar arguments with intent(in) INTEGER, INTENT(in) :: ncID ! ID (number) of an existing and ! open netCDf file CHARACTER(len=*), INTENT(in) :: ncCallType ! flag indicating what type of ! data to be read from this file ! This is used to identify what dimensions are expected. CHARACTER(len=*), INTENT(in) :: ncFileType ! indicates what "type" of netCDF ! files are to be read ! Local scalar variables. INTEGER :: & dimID &! dimension ID (number) ,i &! loop counter ,ierr &! error flag ,ndim ! number of dimensions for a particular case LOGICAL :: errFound ! flag indicating if an error has been found ! Local array variables. CHARACTER(len=10) :: dimName(ndimMax) ! names of the netCDF dimensions !------------------------------------------------------------------------------- errFound = .FALSE. ! Set up the number and names of dimensions in each case. ! Note that the Sun (Solaris) compiler currently requires "array constructor values" ! for character variables to all be of the same length - hence some of the ! trailing blanks in dimName below. SELECT CASE ( ncFileType ) !------------------------------------------------------------------------------- CASE ( ncTypeGswp2 ) ! NB The "pseudo" dimension in GSWP2 files was mistakenly called "Psuedo". ! To allow this code to be used with existing files, we have to indicate ! the misspelled name. SELECT CASE ( ncCallType ) CASE ( 'init_agric','init_grid','init_soil','init_top' ) ndim = 1 dimName(1:ndim) = (/ 'Land' /) CASE ( 'drive_read' ) ndim = 2 dimName(1:ndim) = (/ 'land ', 'tstep' /) CASE ( ' init_hgt' ) ndim = 2 dimName(1:ndim) = (/ 'Land ', 'Psuedo' /) CASE ( 'init_frac', 'veg_read' ) ndim = 3 dimName(1:ndim) = (/ 'Land ', 'Psuedo', 'Time ' /) CASE ( 'init_ic' ) ndim = 3 dimName(1:ndim) = (/ 'Land ', 'Psuedo', 'Soil ' /) CASE default errFound = .TRUE. END SELECT !------------------------------------------------------------------------------- CASE ( ncTypeJules ) SELECT CASE ( ncCallType ) CASE ( 'dump_io' ) ! In fact we do nothing here as this is handled in dump i/o code. ! Just leave the routine (rather untidy!) RETURN CASE default errFound = .TRUE. END SELECT !------------------------------------------------------------------------------- CASE ( ncTypePilps2e ) SELECT CASE ( ncCallType ) CASE ( 'drive_read' ) ndim = 3 dimName(1:ndim) = (/ 'x ', 'y ', 'tstep' /) CASE default errFound = .TRUE. END SELECT !------------------------------------------------------------------------------- CASE ( ncTypePrincet ) SELECT CASE ( ncCallType ) CASE ( 'drive_read' ) ndim = 4 dimName(1:ndim) = (/ 'longitude', 'latitude ', 'z ', 'time ' /) CASE default errFound = .TRUE. END SELECT !------------------------------------------------------------------------------- CASE ( ncTypeTseries ) SELECT CASE ( ncCallType ) CASE ( 'drive_read' ) ndim = 1 dimName(1:ndim) = (/ 'time' /) CASE default errFound = .TRUE. END SELECT !------------------------------------------------------------------------------- CASE ( ncTypeWatch ) ! Similar to ncTypeGSWP2 but with relatively subtle but important differences! SELECT CASE ( ncCallType ) CASE ( 'init_agric', 'init_grid', 'init_soil', 'init_top' ) ndim = 1 dimName(1:ndim) = (/ 'land' /) CASE ( 'init_frac',' init_hgt' ) ndim = 2 dimName(1:ndim) = (/ 'land ', 'pseudo' /) CASE ( 'veg_read' ) ndim = 3 dimName(1:ndim) = (/ 'land ', 'pseudo', 'tstep ' /) CASE default errFound = .TRUE. END SELECT !------------------------------------------------------------------------------- CASE default errFound = .TRUE. END SELECT ! Stop if error found. IF ( errFound ) THEN WRITE(*,*) 'ERROR: check_nc_dims: no code for this combination:' WRITE(*,*) 'ncFileType=',TRIM(ncFileType),' ncCallType=',TRIM(ncCallType) WRITE(*,*)'This code is checking that an existing netCDF file is of the' WRITE(*,*)'type claimed, with the expected DIMENSION names.' WRITE(*,*)'To introduce a new "type" of netCDF file, add more code!' WRITE(*,*)'Stopping in subroutine check_nc_dims' STOP ENDIF !------------------------------------------------------------------------------- ! Check that the dimensions exist in file. !------------------------------------------------------------------------------- DO i=1,ndim ierr = nf90_inq_dimid( ncID,TRIM(dimName(i)),dimID ) ! If error, call error routine. IF ( ierr /= nf90_noErr ) THEN errFound = .TRUE. WRITE(*,*)'ERROR: check_nc_dims: error checking for dimensions.' CALL rwErr( formatNc,.FALSE.,ierr=ierr,ncID=ncID,dimName=dimName(i) & ,errMess1='nf90_inq_dimid' ) ENDIF ENDDO ! Stop if error found. IF ( errFound ) THEN WRITE(*,*) WRITE(*,*)'The code expected to find the following dimensions in the file:' DO i=1,ndim WRITE(*,*) TRIM( dimName(i) ) ENDDO WRITE(*,*)'Stopping in subroutine check_nc_dims' STOP ENDIF END SUBROUTINE check_nc_dims !############################################################################### !############################################################################### SUBROUTINE readVarNcDump( ncID,varName,ncCount,outval ) ! Read a variable from a netCDF dump file. ! Module procedure in module jules_netcdf. IMPLICIT NONE !------------------------------------------------------------------------------- ! Scalar arguments with intent(in) !------------------------------------------------------------------------------- INTEGER, INTENT(in) :: ncID ! ID of netCDF dataset CHARACTER(len=*), INTENT(in) :: varName ! netCDF variable name !------------------------------------------------------------------------------- ! Array arguments with intent(in) !------------------------------------------------------------------------------- INTEGER, INTENT(in) :: ncCount(:) ! counts for variable !------------------------------------------------------------------------------- ! Array arguments with intent(out) !------------------------------------------------------------------------------- REAL, INTENT(out) :: outval(:,:) ! data read from file ! Note that a 3-D variable is read into this 2-D space, not because that ! is sensible, but because that means that existing 2-D code can be reused. !------------------------------------------------------------------------------- ! Local scalar variables. !------------------------------------------------------------------------------- INTEGER :: ierr ! return code INTEGER :: varID ! netCDF ID of a variable !------------------------------------------------------------------------------- ! Get ID of variable. !------------------------------------------------------------------------------- ierr = nf90_inq_varid( ncID,varName,varID ) ! If error, call error routine, indicating a stop. IF( ierr /= nf90_noErr ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,ncID=ncID,varName=varName & ,errMess1='Error from nf90_inq_varid' & ,errMess2='Called from readVarNcDump' ) !------------------------------------------------------------------------------- ! Read data. By default start=1. !------------------------------------------------------------------------------- ierr = nf90_get_var( ncID,varID,outval,count=ncCount ) ! If error, call error routine, indicating a stop. IF( ierr /= nf90_noErr ) CALL rwErr( formatNc,.TRUE. & ,ierr=ierr,ncID=ncID,varName=varName & ,errMess1='Error from nf90_get_var' & ,errMess2='Called from readVarNcDump' ) END SUBROUTINE readVarNcDump !############################################################################### !############################################################################### END MODULE jules_netcdf