!------------------------------------------------------------------------------ ! Global Modeling and Assimilation Office (GMAO) ! ! Goddard Earth Observing System (GEOS) ! ! MAPL Component ! !------------------------------------------------------------------------------ #include "unused_dummy.H" ! !> !### MODULE: `ESMF_CFIOUtilMod` ! ! Author: GMAO SI-Team ! ! `ESMF_CFIOUtilMod` - utility subroutines for CFIO ! ! The module `ESMF_CFIOUtilMod` provides all the necessary utility subroutines for CFIO. ! !#### History !- Sep2004 Baoyu Yin Modified from ESMF_CFIOMod.F90. !- Mar2005 Baoyu Yin Integrated all routines of GFIO into CFIO. !- Apr2005 Baoyu Yin Added Get_Missing routine !- Apr2006 da Silva Incorporated m_StrTemplate to eliminate mpeu dependency. !- Jun2006 Baoyu yin Added cyclic option to GetVar. !- Jul2006 da Silva Eliminated read(str,fmt) to parse time; replaced ! with more robust mod() calculations. ! module ESMF_CFIOUtilMod ! ! !USES: use ESMF_CFIOBaseMod use netcdf use, intrinsic :: iso_fortran_env, only: INT16, INT32, INT64 use, intrinsic :: iso_fortran_env, only: REAL32, REAL64 implicit none #if defined(HDFEOS) || defined(HDFSD) include "hdf.f90" #endif !------------------------------------------------------------------------------ ! !EOP !============================================================================== integer, parameter :: GCTP_GEO=0 integer, parameter :: HDFE_NOMERGE=0 integer, parameter :: HDFE_GD_LL=2 integer, parameter :: NDIMS_MAX = 4 integer, parameter :: MAX_VAR_DIMS = 32 character(len=7), parameter :: GRID_NAME='EOSGRID' integer, parameter :: MAXCHR = 256 integer, parameter :: PACK_BITS = 32766 integer, parameter :: PACK_FILL = 32767 integer, parameter :: MLEN = 1024 ! Max. length of an attribute integer, parameter :: MVARLEN = 256 ! Max. length of a variable name ! Define a new data type "List" -- private data type for variable and ! global attributes type iNode integer :: count integer, pointer :: intData(:) character(len=MVARLEN) :: name character(len=MVARLEN) :: vName type (iNode), pointer :: next => null() end type iNode type rNode integer :: count real, pointer :: realData(:) character(len=MVARLEN) :: name character(len=MVARLEN) :: vName type (rNode), pointer :: next => null() end type rNode type cNode integer :: count character(len=MLEN) :: charData character(len=MVARLEN) :: name character(len=MVARLEN) :: vName type (cNode), pointer :: next => null() end type cNode ! StrTemplate parameters (adapted from MPEU) private :: die, perr integer, parameter :: stderr = 6 character(len=*),parameter :: myname='m_StrTemplate' character(len=3),parameter,dimension(12) :: mon_lc =(/& 'jan','feb','mar','apr','may','jun',& 'jul','aug','sep','oct','nov','dec'/) character(len=3),parameter,dimension(12) :: mon_wd =(/& 'Jan','Feb','Mar','Apr','May','Jun',& 'Jul','Aug','Sep','Oct','Nov','Dec'/) character(len=3),parameter,dimension(12) :: mon_uc =(/& 'JAN','FEB','MAR','APR','MAY','JUN',& 'JUL','AUG','SEP','OCT','NOV','DEC'/) interface GetDate module procedure GetDateInt8 module procedure GetDateInt4 end interface !------------------------------------------------------------------------------ contains !------------------------------------------------------------------------------ !> ! `addList` -- put user defined attribute into a list. ! subroutine addList(name, count, vName, attInt, attChar, attReal, & iList, rList, cList) ! ! !INPUT PARAMETERS: ! character(len=*), intent(in) :: name integer, intent(in) :: count character(len=*), intent(in), OPTIONAL :: vName integer, intent(in), OPTIONAL :: attInt(*) real, intent(in), OPTIONAL :: attReal(*) character(len=MLEN), intent(in), OPTIONAL :: attChar ! ! !OUTPUT PARAMETERS: ! type(iNode), pointer, OPTIONAL :: iList !! int list type(rNode), pointer, OPTIONAL :: rList !! real list type(cNode), pointer, OPTIONAL :: cList !! char list ! !------------------------------------------------------------------------------ type(iNode), pointer :: p, q type(rNode), pointer :: rp, rq type(cNode), pointer :: cp, cq integer :: i ! put int attribute into an integer list if ( present(attInt) .and. present(iList)) then if ( .not. associated(iList) ) then allocate(iList) iList%count = count iList%name = name if (present(vName)) iList%vName = vName allocate(iList%intData(count)) do i =1, count iList%intData(i) = attInt(i) end do iList%next => null() else q => iList p => iList%next do while ( associated(p) ) q => p p => p%next end do allocate(p) p%count = count p%name = name if (present(vName)) p%vName = vName allocate(p%intData(count)) do i =1, p%count p%intData(i) = attInt(i) end do q%next => p end if end if ! put real attribute into a real list if ( present(attReal) .and. present(rList)) then if ( .not. associated(rList) ) then allocate(rList) rList%count = count rList%name = name if (present(vName)) rList%vName = vName allocate(rList%realData(count)) do i =1, count rList%realData(i) = attReal(i) end do rList%next => null() else rq => rList rp => rList%next do while ( associated(rp) ) rq => rp rp => rp%next end do allocate(rp) rp%count = count rp%name = name if (present(vName)) rp%vName = vName allocate(rp%RealData(count)) do i =1, rp%count rp%realData(i) = attReal(i) end do rq%next => rp end if end if ! put char attribute into a char list if ( present(attChar) .and. present(cList)) then if ( .not. associated(cList) ) then allocate(cList) cList%count = count cList%name = name if (present(vName)) cList%vName = vName cList%charData = attChar(1:count) cList%next => null() else cq => cList cp => cList%next do while ( associated(cp) ) cq => cp cp => cp%next end do allocate(cp) cp%count = count cp%name = name if (present(vName)) cp%vName = vName cp%charData = attChar(1:count) cq%next => cp end if end if end subroutine addList !------------------------------------------------------------------------------ !> ! `getList` -- retrieve defined attributes from a list ! subroutine getList(iList, nIntAtt, intAttNames, intAttCnts, intAtts, & rList, nRealAtt, realAttNames, realAttCnts, realAtts, & cList, nCharAtt, charAttNames, charAttCnts, charAtts, vNames) ! ! !INPUT PARAMETERS: ! type(iNode), pointer, OPTIONAL :: iList !! int list type(cNode), pointer, OPTIONAL :: cList !! real list type(rNode), pointer, OPTIONAL :: rList !! char list ! ! !OUTPUT PARAMETERS: ! integer, OPTIONAL :: nIntAtt !! num of int att integer, OPTIONAL :: nRealAtt !! num of real att integer, OPTIONAL :: nCharAtt !! num of char att character(len=MLEN), pointer, OPTIONAL :: intAttNames(:) character(len=MLEN), pointer, OPTIONAL :: realAttNames(:) character(len=MLEN), pointer, OPTIONAL :: charAttNames(:) integer, OPTIONAL, pointer :: intAttCnts(:) !!data count in int att integer, OPTIONAL, pointer :: realAttCnts(:)!!data count in real att integer, OPTIONAL, pointer :: charAttCnts(:)!!data count in char att integer, OPTIONAL, pointer :: intAtts(:,:) !!int attribute real, OPTIONAL, pointer :: realAtts(:,:) !!real attribute character(len=MLEN), pointer, OPTIONAL :: charAtts(:) !! char att character(len=MLEN), pointer, OPTIONAL :: vNames(:) ! !------------------------------------------------------------------------------ type(iNode), pointer :: p ! pointer for integer list type(rNode), pointer :: rp ! pointer for real list type(cNode), pointer :: cp ! pointer for char list integer :: maxLen ! length of a list integer :: cnt ! max number of data in nodes integer :: i maxLen = 0 cnt = 0 ! get attributes from an integer list if ( present(iList) ) then allocate(p) p = iList call getMaxLenCnt(maxLen, cnt, iList=iList) if ( present(nIntAtt) ) nIntAtt = cnt allocate(intAttCnts(cnt), intAttNames(cnt), intAtts(cnt,maxLen)) if (present(vNames)) allocate(vNames(cnt)) i = 1 do while ( associated(p) ) intAttNames(i) = trim(p%name) if (present(vNames)) vNames(i) = trim(p%vName) intAttCnts(i) = p%count intAtts(i,:) = 0 intAtts(i,1:size(p%intData)) = p%intData p => p%next i = i + 1 end do end if ! get attributes from a real list if ( present(rList) ) then allocate(rp) rp = rList call getMaxLenCnt(maxLen, cnt, rList=rList) if (present(nRealAtt)) nRealAtt = cnt allocate(realAttCnts(cnt),realAttNames(cnt),realAtts(cnt,maxLen)) if (present(vNames)) allocate(vNames(cnt)) i = 1 do while ( associated(rp) ) realAttNames(i) = trim(rp%name) realAttCnts(i) = rp%count if (present(vNames)) vNames(i) = trim(rp%vName) realAtts(i,1:size(rp%realData)) = rp%realData rp => rp%next i = i + 1 end do end if ! get attributes from a char list if ( present(cList) ) then allocate(cp) cp = cList call getMaxLenCnt(maxLen, cnt, cList=cList) if ( present(nCharAtt) ) nCharAtt = cnt allocate(charAttCnts(cnt), charAttNames(cnt), charAtts(cnt)) if (present(vNames)) allocate(vNames(cnt)) i = 1 do while ( associated(cp) ) charAttNames(i) = trim(cp%name) if (present(vNames)) vNames(i) = trim(cp%vName) charAttCnts(i) = cp%count charAtts(i) = cp%charData cp => cp%next i = i + 1 end do end if end subroutine getList !------------------------------------------------------------------------------ !> ! `getMaxLenCnt` -- get length of a list and max number of data in the nodes ! ! Get length of a list and max number of data in the nodes so that ! maxLen/count can been used to allocate array. ! subroutine getMaxLenCnt(maxLen, count, iList, rList, cList) ! ! !INPUT PARAMETERS: ! type(iNode), pointer, OPTIONAL :: iList type(cNode), pointer, OPTIONAL :: cList type(rNode), pointer, OPTIONAL :: rList ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: maxLen integer, intent(out) :: count ! !------------------------------------------------------------------------------ type(iNode), pointer :: p ! int pointer for int list type(rNode), pointer :: rp ! real pointer for real list type(cNode), pointer :: cp ! char pointer for char list count = 0 maxLen = 0 ! get maxLen/count from a integer list if ( present(iList) ) then allocate(p) p = iList do while ( associated(p) ) if (p%count .gt. maxLen) maxLen = p%count count = count + 1 p => p%next end do end if ! get maxLen/count from a real list if ( present(rList) ) then allocate(rp) rp = rList do while ( associated(rp) ) if (rp%count .gt. maxLen) maxLen = rp%count count = count + 1 rp => rp%next end do end if ! get maxLen/count from a character list if ( present(cList) ) then allocate(cp) cp = cList do while ( associated(cp) ) if (cp%count .gt. maxLen) maxLen = cp%count count = count + 1 cp => cp%next end do end if end subroutine getMaxLenCnt !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! ! `CFIO_DimInquire` -- Gets dimension information from a CFIO file. ! ! This routine is used to get dimension information from ! an existing CFIO file. This dimension information can ! subsequently be used to allocate arrays for reading data ! from the file. For more complete information about the ! contents of a file, Cfio\_Inquire should be used. ! !#### History !- 1998.07.02 Lucchesi Initial interface design. !- 1998.08.05 Lucchesi Added "ngatts" !- 1998.09.24 Lucchesi Revamped error codes !- 1998.12.22 Lucchesi Added IdentifyDim and associated code !- 1999.01.04 Lucchesi Changed variable initialization !- 2008.03.14 Kokron Initialize stationFile to false ! subroutine CFIO_DimInquire (fid,im,jm,km,lm,nvars,ngatts,vdir,rc) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid ! File handle ! ! !OUTPUT PARAMETERS: ! integer im !! Size of longitudinal dimension integer jm !! Size of latitudinal dimension integer km !! Size of vertical dimension !! km=0 if surface-only file integer lm !! Number of times integer nvars !! Number of variables integer ngatts !! Number of global attributes integer, optional :: vdir !! Positive vertical direction. !! If `-1`, level 1 in the file is TOA. !! If `1`, level 1 in the file is the surface. !! If `0`, the file has no vertical co-ordinate (default). integer, optional :: rc !! Error return code: !! rc = 0 all is well !! rc = -19 unable to identify coordinate variable !! !! NetCDF Errors !! ------------- !! rc = -40 error from NF90_INQ_VARID !! rc = -41 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (lat or lon) !! rc = -42 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (lev) !! rc = -43 error from NF90_INQ_VARID (time variable) !! rc = -47 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (time) !! rc = -48 error from NF90_INQUIRE !! rc = -53 error from NF90_GET_ATT ! !------------------------------------------------------------------------- integer dimId, i character(len=MAXCHR) dimName character(len=MAXCHR) stdName character(len=MAXCHR) dimUnits character(len=MAXCHR) posStr character(len=MAXCHR) vname integer dimSize integer nDims logical surfaceOnly logical stationFile integer myIndex integer varType, nvDims, vDims(MAXVDIMS), nvAtts integer tmpNvar logical :: cs_found integer :: vdir_ integer :: found_xc, found_yc,vid ! Initialize variables surfaceOnly = .FALSE. stationFile = .false. vdir_ = -1 ! Assume 3D and same orientation found_xc = 0 found_yc = 0 ! Check FID here. ! Check to make sure max string lengths are large enough. NetCDF defines ! MAXNCNAM, but it can't be used in a character(len=MAXNCNAM) statement. if (MAXCHR .LT. MAXNCNAM) then print *, 'CFIO_DimInquire warning: MAXNCNAM is larger than ', & 'dimName array size.' endif ! Get basic information from file. rc = NF90_INQUIRE(fid, nDims, nvars, ngatts, dimId) if (err("DimInqure: NF90_INQUIRE failed",rc,-48) .NE. 0)return ! Subtract dimension variables from the variable count. cs_found = .false. tmpNvar = nvars do i=1,nvars rc = NF90_INQUIRE_VARIABLE (fid,i,vname,varType,nvDims,vDims,nvAtts) if (err("DimInquire: variable inquire error",rc,-52) .NE. 0) & return if (nvDims .EQ. 1 .or. trim(vname) .eq. 'time_bnds') then tmpNvar = tmpNvar - 1 endif if (vname == 'nf') then cs_found = .true. end if enddo if (cs_found) then tmpNvar = tmpNvar - 4 found_xc = NF90_INQ_VARID(fid,"corner_lons",vid) if (found_xc ==0) tmpNvar = tmpNvar - 1 found_yc = NF90_INQ_VARID(fid,"corner_lats",vid) if (found_yc ==0) tmpNvar = tmpNvar - 1 end if nvars = tmpNvar ! Extract dimension information do i=1,nDims rc = NF90_INQUIRE_DIMENSION (fid, i, dimName, dimSize) if (err("DimInqure: can't get dim info",rc,-41) .NE. 0) return if (index(dimName,'station') .gt. 0) then stationFile = .true. im = dimSize jm = dimSize cycle end if if (trim(dimName) .eq. 'nv') cycle if (trim(dimName) .eq. 'nf') cycle if (trim(dimName) .eq. 'ncontact') cycle if (trim(dimName) .eq. 'XCdim') cycle if (trim(dimName) .eq. 'YCdim') cycle if (trim(dimName) .eq. 'orientationStrLen') cycle rc = NF90_INQ_VARID (fid, dimName, dimId) if (err("DimInqure: NF90_INQ_VARID failed",rc,-40) .NE. 0) return ! If it has the standard_name attribute, use that instead rc = NF90_GET_ATT(fid,dimId,'standard_name',stdName) if (rc .ne. 0) stdName = Trim(dimName) dimunits = '' rc = NF90_GET_ATT(fid,dimId,'units',dimUnits) if (.not. cs_found) then !ALT: the new cubed-sphere format might be missing some units if (err("DimInqure: could not get units for dimension",rc,-53)& .NE. 0) return end if myIndex = IdentifyDim (dimName, dimUnits) if ( myIndex .EQ. 0 ) then im = dimSize else if ( myIndex .EQ. 1 ) then jm = dimSize if(cs_found) jm = jm*6 else if ( myIndex .EQ. 2 ) then km = dimSize if (km.eq.1) then ! 2D vdir_ = 0 else rc = NF90_GET_ATT(fid,dimId,'positive',posStr) if (rc.eq.0) then if ((Trim(posStr)=="up").or.(Trim(posStr)=="Up")) then ! Level 1 = surface vdir_ = 1 elseif ((Trim(posStr)=="down").or.(Trim(posStr)=="Down")) then ! Level 1 = TOA vdir_ = -1 endif endif endif else if ( myIndex .EQ. 3 ) then lm = dimSize !print *, "dimUnits: ", trim(dimUnits) !print *, "dimName: ", trim(dimName) ! else ! print *, 'CFIO_DimInquire: Coordinate variable ', & ! TRIM(dimName),' with units of ',TRIM(dimUnits), & ! ' is not understood.' ! rc = -19 ! return endif enddo if (cs_found .and. nDims == 6) surfaceOnly = .TRUE. if (cs_found .and. nDims == 8) surfaceOnly = .TRUE. if (nDims .EQ. 3 .and. .NOT. stationFile) then surfaceOnly = .TRUE. endif if (surfaceOnly) then km = 0 vdir_ = 0 endif if (present(vdir)) vdir = vdir_ rc=0 return end subroutine CFIO_DimInquire !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `GetDateTimeVec` - Get date/time of file samples ! ! This routine returns the date/times on file. ! !#### History !- 1999.11.01 da Silva Initial code. !- 1999.11.08 da Silva Generic time coordinate variable (no name assumed). !- 2000.10.18 Lucchesi Added ParseTimeUnits subroutine to handle a wider variety of Units formats. !- 2016.12.06 Eastham Adapted to return the vector of offsets instead ! subroutine GetDateTimeVec ( fid, begDate, begTime, incVec, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! file ID ! ! !OUTPUT PARAMETERS: ! integer :: begDate !! Beginning date integer :: begTime !! Beginning time integer(Kind=INT64) :: incVec(:) !! Vector of offsets (seconds) integer :: rc !! error return code ! !------------------------------------------------------------------------- integer i, timeId, hour, min, sec, corner(1) !integer incSecs integer year, month, day character(len=MAXCHR) timeUnits, dimUnits !character(len=MAXCHR) strTmp character(len=MAXCHR) varName, dimName, stdName integer type, nvDims, vdims(MAXVDIMS), nvAtts, dimSize integer nDims, nvars, ngatts, dimId ! Time conversion local variables real(kind=REAL32) rtime, rtime_array(1) real(kind=REAL64) dtime, dtime_array(1) integer(kind=INT16) itime, itime_array(1) integer(kind=REAL32) ltime, ltime_array(1) !integer t1 integer newDate, newTime ! We now have the possibility of a very large interval integer(Kind=INT64) :: t1Long, t2Long, tMultLong, incSecsLong integer(Kind=INT64),allocatable :: incVecLong(:) ! Vector of offsets (seconds) ! Get the starting date and time ! --------------------------------------------------------- ! Start by determing the ID of the time coordinate variable ! --------------------------------------------------------- timeId = -1 rc = NF90_INQUIRE (fid, nDims, nvars, ngatts, dimId) if (err("GetDateTimeVec: NF90_INQUIRE failed",rc,-48) .NE. 0)return do i=1,nDims rc = NF90_INQUIRE_DIMENSION (fid, i, dimName, dimSize) if (err("GetDateTimeVec: can't get dim info",rc,-41) .NE. 0) return if (index(dimName,'station') .gt. 0) cycle if (trim(dimName) .eq. 'nv') cycle if ( index(dimName,'edges') .gt. 0 ) cycle if (trim(dimName) .eq. 'nv') cycle if (trim(dimName) .eq. 'nf') cycle if (trim(dimName) .eq. 'ncontact') cycle if (trim(dimName) .eq. 'XCdim') cycle if (trim(dimName) .eq. 'YCdim') cycle if (trim(dimName) .eq. 'orientationStrLen') cycle rc = NF90_INQ_VARID (fid, dimName, dimId) if (err("GetDateTimeVec: NF90_INQ_VARID failed",rc,-40) .NE. 0) return ! If it has the standard_name attribute, use that instead rc = NF90_GET_ATT(fid,dimId,'standard_name',stdName) if (rc .ne. 0) stdName = Trim(dimName) rc = NF90_GET_ATT(fid,dimId,'units',dimUnits) if (err("GetDateTimeVec: could not get units for dimension",rc,-53)& .NE. 0) return if ( IdentifyDim (stdName, dimUnits) .eq. 3 ) then timeId = dimId timeUnits = dimUnits exit end if end do if ( timeId .lt. 0 ) then rc = -43 print *, "GetDateTimeVec: could not find time coord" return end if ! Attempt to parse the COARDS compliant time units ! -------------------------------------------------- !ams rc = NF90_GET_ATT(fid,timeId,'units',timeUnits) !ams if (err("GetDateTimeVec: missing time.units",rc,-44) .NE. 0) return i = index(timeUnits,'since') if ( i .le. 0 ) then if (err("GetDateTimeVec: invalid time units",1,-44) .NE. 0) return endif ! Call to ParseTimeUnits replaces an internal read, that made assumptions ! about the format of the Time Units string that were not always true. ! (RL: 10/2000) call ParseTimeUnits ( timeUnits, year, month, day, hour, min, sec, rc ) begDate = year*10000 + month*100 + day begTime = hour*10000 + min*100 + sec ! Retreive time vector. ! ------------------------- rc = NF90_INQUIRE_VARIABLE (fid, timeID, varName, type, nvDims, vDims, & nvAtts) if (err("GetDateTimeVec: error in time variable inquire",& rc,-52) .NE. 0) return allocate(incVecLong(dimSize)) incVecLong(:) = 0 do i=1,dimsize if ( type .eq. NF90_FLOAT ) then corner(1) = i rc = NF90_GET_VAR(fid,timeID,rtime_array,corner,(/1/)) rtime = rtime_array(1) incVecLong(i) = int(rtime,INT64) else if ( type .eq. NF90_DOUBLE ) then corner(1) = i rc = NF90_GET_VAR(fid,timeID,dtime_array,corner,(/1/)) dtime = dtime_array(1) incVecLong(i) = int(dtime,INT64) else if ( type .eq. NF90_SHORT ) then corner(1) = i rc = NF90_GET_VAR(fid,timeID,itime_array,corner,(/1/)) itime = itime_array(1) incVecLong(i) = int(itime,INT64) else if ( type .eq. NF90_INT ) then corner(1) = i rc = NF90_GET_VAR(fid,timeID,ltime_array,corner,(/1/)) ltime = ltime_array(1) incVecLong(i) = int(ltime,INT64) else if (err("GetDateTimeVec: invalid time data type",& 1,-44) .NE. 0) return endif end do ! Convert time increment to seconds if necessary ! ---------------------------------------------- if ( timeUnits(1:6) .eq. 'minute' ) then tMultLong = 60 else if ( timeUnits(1:4) .eq. 'hour' ) then tMultLong = 60 * 60 else if ( timeUnits(1:3) .eq. 'day' ) then tMultLong = 60 * 60 * 24 else if (err("GetDateTimeVec: invalid time unit name",& 1,-44) .NE. 0) return endif ! Combine the first time offset with the reference time to get the beginning ! date and time ! ----------------------------------------------------------------------------- t1Long = incVecLong(1) Call GetDate(begDate,begTime,t1Long*tMultLong,newDate,newTime,rc) begDate=newDate begTime=newTime ! Convert all the offsets to reference the first sample ! ----------------------------------------------------------------------------- incVec(:) = 0 do i=1,dimsize t2Long = incVecLong(i) incSecsLong = (t2Long - t1Long)*tMultLong ! If (incSecsLong.gt.huge(t1)) Then ! print *, 'Time interval too large' ! rc = -10 ! return ! End If ! Sometimes we actually do need to return this as a longlong.. !incSecs = int(incSecsLong,4) !incVec(i) = incSecs incVec(i) = incSecsLong end do deallocate(incVecLong) !ams print *, 'begdate, begtime, incsecs: ',begdate, begtime, incsecs rc = 0 ! all done return end subroutine GetDateTimeVec !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `GetBegDateTime` - Get begin date/time of file ! ! This routine returns the begin date/begin time on file. ! For a native CFIO file, it simply returns the value of ! global attributes begin_date/begin_time. If these do no exist then ! it attempts to parse the COARDS compliant time unit. ! !#### History !- 1999.11.01 da Silva Initial code. !- 1999.11.08 da Silva Generic time coordinate variable (no name assumed). !- 2000.10.18 Lucchesi Added ParseTimeUnits subroutine to handle a wider variety of Units formats. ! subroutine GetBegDateTime ( fid, begDate, begTime, incSecs, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! file ID ! ! !OUTPUT PARAMETERS: ! integer begDate !! beginning date integer begTime !! beginning time integer incSecs !! time increment in secs integer rc !! error return code ! !------------------------------------------------------------------------- integer i, timeId, hour, min, sec, corner(1), timInc integer year, month, day character(len=MAXCHR) timeUnits, dimUnits, stdName character(len=MAXCHR) varName, dimName integer type, nvDims, vdims(MAXVDIMS), nvAtts, dimSize integer nDims, nvars, ngatts, dimId ! Time conversion local variables real(kind=REAL32) :: rtime, rtime_array(1) real(kind=REAL64) :: dtime, dtime_array(1) integer(kind=INT16) :: itime, itime_array(1) integer(kind=INT32) :: ltime, ltime_array(1) integer t1, t2, tMult, newDate, newTime ! Start by determing the ID of the time coordinate variable ! --------------------------------------------------------- timeId = -1 rc = NF90_INQUIRE (fid, nDims, nvars, ngatts, dimId) if (err("GetBegDateTime: NF90_INQUIRE failed",rc,-48) .NE. 0)return do i=1,nDims rc = NF90_INQUIRE_DIMENSION (fid, i, dimName, dimSize) if (err("GetBegDateTime: can't get dim info",rc,-41) .NE. 0) return if (index(dimName,'station') .gt. 0) cycle if (trim(dimName) .eq. 'nv') cycle if ( index(dimName,'edges') .gt. 0 ) cycle if (dimName=='nf') cycle if (dimName=='orientationStrLen') cycle if (dimName=='ncontact') cycle if (trim(dimName) .eq. 'XCdim') cycle if (trim(dimName) .eq. 'YCdim') cycle rc = NF90_INQ_VARID (fid, dimName, dimId) if (err("GetBegDateTime: NF90_INQ_VARID failed",rc,-40) .NE. 0) return ! If it has the standard_name attribute, use that instead rc = NF90_GET_ATT(fid,dimId,'standard_name',stdName) if (rc .ne. 0) stdName = Trim(dimName) rc = NF90_GET_ATT(fid,dimId,'units',dimUnits) if (err("GetBegDateTime: could not get units for dimension",rc,-53)& .NE. 0) return if ( IdentifyDim (stdName, dimUnits) .eq. 3 ) then timeId = dimId timeUnits = dimUnits exit end if end do if ( timeId .lt. 0 ) then rc = -43 print *, "GetBegDateTime: could not find time coord" return end if ! Try assuming this file has been written with CFIO ! ------------------------------------------------- rc = NF90_GET_ATT(fid,timeId,'begin_date',begDate) if ( rc .eq. 0 ) then rc = NF90_GET_ATT(fid,timeId,'begin_time',begTime) end if ! Well, it must be a native CFIO file ! ----------------------------------- if ( rc .eq. 0 ) then rc = NF90_GET_ATT(fid,timeId,'time_increment',timInc) if (err("GetBegDateTime: missing time increment",rc,-44) .NE. 0) & return !ams write (strTmp,'(i6)') timinc !ams read (strTmp,'(3I2)') hour, min, sec call CFIO_parseIntTime ( timinc, hour, min, sec ) incSecs = hour*3600 + min*60 + sec !ams print *, 'begdate, begtime, incsecs: ',begdate, begtime, incsecs return ! all done. end if ! If could not find begin_date/begin_time attributes ! then this is not a native CFIO file. In this case ! attempt to parse the COARDS compliant time units ! -------------------------------------------------- !ams rc = NF90_GET_ATT(fid,timeId,'units',timeUnits) !ams if (err("GetBegDateTime: missing time.units",rc,-44) .NE. 0) return i = index(timeUnits,'since') if ( i .le. 0 ) then if (err("GetBegDateTime: invalid time units",1,-44) .NE. 0) return endif ! Call to ParseTimeUnits replaces an internal read, that made assumptions ! about the format of the Time Units string that were not always true. ! (RL: 10/2000) call ParseTimeUnits ( timeUnits, year, month, day, hour, min, sec, rc ) begDate = year*10000 + month*100 + day begTime = hour*10000 + min*100 + sec ! Determine time increment. ! ------------------------- rc = NF90_INQUIRE_VARIABLE (fid, timeID, varName, type, nvDims, vDims, & nvAtts) if (err("GetBegDateTime: error in time variable inquire",& rc,-52) .NE. 0) return if ( type .eq. NF90_FLOAT ) then corner(1) = 1 rc = NF90_GET_VAR(fid,timeID,rtime_array,corner,(/1/)) rtime = rtime_array(1) t1 = int(rtime) corner(1) = 2 rc = NF90_GET_VAR(fid,timeID,rtime_array,corner,(/1/)) rtime = rtime_array(1) t2 = int(rtime) else if ( type .eq. NF90_DOUBLE ) then corner(1) = 1 rc = NF90_GET_VAR(fid,timeID,dtime_array,corner,(/1/)) dtime = dtime_array(1) t1 = int(dtime) !ams print *, t1, dtime, rc corner(1) = 2 rc = NF90_GET_VAR(fid,timeID,dtime_array,corner,(/1/)) dtime = dtime_array(1) t2 = int(dtime) !ams print *, t2, dtime, rc else if ( type .eq. NF90_SHORT ) then corner(1) = 1 rc = NF90_GET_VAR(fid,timeID,itime_array,corner,(/1/)) itime = itime_array(1) t1 = itime corner(1) = 2 rc = NF90_GET_VAR(fid,timeID,itime_array,corner,(/1/)) itime = itime_array(1) t2 = itime else if ( type .eq. NF90_INT ) then corner(1) = 1 rc = NF90_GET_VAR(fid,timeID,ltime_array,corner,(/1/)) ltime = ltime_array(1) t1 = ltime corner(1) = 2 rc = NF90_GET_VAR(fid,timeID,ltime_array,corner,(/1/)) ltime = ltime_array(1) t2 = ltime else if (err("GetBegDateTime: invalid time data type",& 1,-44) .NE. 0) return endif ! Convert time increment to seconds if necessary ! ---------------------------------------------- incSecs = t2 - t1 if ( timeUnits(1:6) .eq. 'minute' ) then tMult = 60 else if ( timeUnits(1:4) .eq. 'hour' ) then tMult = 60 * 60 else if ( timeUnits(1:3) .eq. 'day' ) then tMult = 60 * 60 * 24 else if (err("GetBegDateTime: invalid time unit name",& 1,-44) .NE. 0) return endif incSecs = incSecs*tMult ! Combine the first time offset with the reference time to get the beginning ! date and time ! ----------------------------------------------------------------------------- Call GetDate(begDate,begTime,t1*tMult,newDate,newTime,rc) begDate=newDate begTime=newTime !ams print *, 'begdate, begtime, incsecs: ',begdate, begtime, incsecs incSecs = max ( 1, incSecs ) rc = 0 ! all done return end subroutine GetBegDateTime !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `IdentifyDim` - Identify a cooridate variable ! ! This function attempts to identify a coordiante variable ! from the name or units of the variable. It does so by ! attempting to match the units specified in the COARDS ! conventions or by checking the name against commonly used ! names. ! ! RETURN VALUES: !* 0 = X dimension (longitude) !* 1 = Y dimension (latitude) !* 2 = Z dimension (level) !* 3 = Time !* -1 = Unable to determine dimension ! !#### History !- 1998.12.22 Lucchesi Initial coding. !- 1999.11.02 da Silva Made LATS4D compatible, !- 2001.01.02 da Silva Cimcurventing PGI bugs. ! integer function IdentifyDim (dimName, dimUnits) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! character(len=*) dimName !! Name of the coordinate variable character(len=*) dimUnits !! Units of the coordinate variable ! !------------------------------------------------------------------------- if (TRIM(dimUnits) .EQ. "degrees_north" ) then IdentifyDim = 1 return end if if (TRIM(dimUnits) .EQ. "hPa" ) then IdentifyDim = 2 return end if if ( trim(dimName) .eq. "time" ) then IdentifyDim = 3 return end if if (TRIM(dimUnits) .EQ. "degrees_east" .OR. & trim(dimName) .eq. "longitude" .OR. & trim(dimName) .eq. "lon" ) then IdentifyDim = 0 else if (TRIM(dimUnits) .EQ. "degrees_north" ) then IdentifyDim = 1 else if ( trim(dimName) .eq. "latitude" .OR. & trim(dimName) .eq. "lat" ) then IdentifyDim = 1 else if (INDEX(dimName,"lev") .NE. 0 .OR. & INDEX(dimName,"Height") .NE. 0) then IdentifyDim = 2 else if (TRIM(dimUnits) .EQ. "mb" .OR. & TRIM(dimUnits) .EQ. "millibar" .OR. & TRIM(dimUnits) .EQ. "sigma_level" .OR. & TRIM(dimUnits) .EQ. "hPa") then IdentifyDim = 2 else if (trim(dimName) .eq. "TIME" .OR. & trim(dimName) .eq. "TIME:EOSGRID" .OR. & trim(dimName) .eq. "time" .OR. & trim(dimName) .eq. "Time") then IdentifyDim = 3 else IdentifyDim = -1 endif if (TRIM(dimUnits) .EQ. "degrees_east" .OR. & INDEX(dimName,"XDim") .NE. 0 .OR. & INDEX(dimName,"lon") .NE. 0) then IdentifyDim = 0 else if (TRIM(dimUnits) .EQ. "degrees_north" .OR. & INDEX(dimName,"YDim") .NE. 0 .OR. & INDEX(dimName,"lat") .NE. 0) then IdentifyDim = 1 else if (INDEX(dimName,"lev") .NE. 0 .OR. & INDEX(dimName,"Height") .NE. 0) then IdentifyDim = 2 else if (TRIM(dimUnits) .EQ. "mb" .OR. & TRIM(dimUnits) .EQ. "millibar" .OR. & TRIM(dimUnits) .EQ. "sigma_level" .OR. & TRIM(dimUnits) .EQ. "hPa") then IdentifyDim = 2 else if (INDEX(dimName,"TIME") .NE. 0 .OR. & INDEX(dimName,"time") .NE. 0 .OR. & INDEX(dimName,"Time") .NE. 0) then IdentifyDim = 3 else IdentifyDim = -1 endif end function IdentifyDim subroutine CFIO_parseIntTime ( hhmmss, hour, min, sec ) integer, intent(in) :: hhmmss integer, intent(out) :: hour, min, sec hour = hhmmss / 10000 min = mod(hhmmss,10000)/100 sec = mod(hhmmss,100) end subroutine CFIO_parseIntTime !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `strToInt` - convert timeString to integer date and time ! ! This function attempts to identify a coordiante variable ! from the name or units of the variable. It does so by ! attempting to match the units specified in the COARDS ! conventions or by checking the name against commonly used ! names. ! !#### History !- 2004.10.04 B. Yin first version. ! subroutine strToInt(timeString, date, begTime) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! character(len=*), intent(in) :: timeString !! string expression of data and time ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: date integer, intent(out) :: begTime ! !------------------------------------------------------------------------- integer :: iCnt, jCnt, dLen character(len=MVARLEN) :: sDate, sTime character(len=MVARLEN) :: strDate, strTime character :: char dLen = index(timeString, 'T' ) sDate = timeString(1:dLen-1) sTime = timeString(dLen+1:len(trim(timeString))) jCnt = 1 strDate = '' do iCnt = 1, len(sDate) char = sDate(iCnt:iCnt) if (char .ne. ':' .and. char .ne. '-') then strDate(jCnt:jCnt) = char jCnt = jCnt + 1 end if end do jCnt = 1 strTime = '' do iCnt = 1, len(sTime) char = sTime(iCnt:iCnt) if (char .ne. ':' .and. char .ne. '-') then strTime(jCnt:jCnt) = char jCnt = jCnt + 1 end if end do read(strDate,*) date read(strTime,*) begTime return end subroutine strToInt !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `CFIO_Open` -- Opens an existing DAO gridded file ! ! This routine opens an existing DAO gridded file. The file ! mode will be read/write. If the application already knows ! the contents of the file, it may begin interaction with the ! file using the returned file handle. Otherwise, the file ! handle can be used with the "inquire" routines to gather ! information about the contents. A negative return code ! indicates there were problems opening the file. ! !#### History !- 1998.07.02 Lucchesi Initial interface design. !- 1998.07.07 Lucchesi Initial coding. !- 1998.12.09 Lucchesi Corrected for ncopn bug. ! subroutine CFIO_Open ( fname, fmode, fid, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! character(len=*) fname !! File name integer fmode !! File mode: !! 0 for READ-WRITE !! non-zero for READ-ONLY ! ! !OUTPUT PARAMETERS: ! integer fid !! File handle integer rc !! Error return code: !! rc = 0 All is well !! rc = -39 error from ncopn (file open) ! !------------------------------------------------------------------------- if ( fmode .EQ. 0) then fid = ncopn (fname, NCWRITE, rc) else fid = ncopn (fname, NCNOWRIT, rc) endif if (fid .LT. 0) then ! ncopn has a bug. error codes should rc = fid ! be returned in rc, but in reality they endif ! are returned in fid. if (err("Open: error opening file",rc,-39) .NE. 0) return rc = 0 return end subroutine CFIO_Open !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `CFIO_Close` -- Closes file ! ! This routine is used to close an open CFIO stream. ! !#### History !- 1997.10.13 da Silva/Lucchesi Initial interface design. !- 1998.03.30 Lucchesi Documentation expanded. Clean-up of code. Added rc. ! subroutine CFIO_Close ( fid, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle ! ! !OUTPUT PARAMETERS: ! integer rc !! Error return code: !! rc = 0 all is well !! !! NetCDF Errors !! ------------- !! rc = -54 error from ncclos (file close) ! !------------------------------------------------------------------------- call ncclos (fid, rc) if (err("Close: error closing file",rc,-54) .NE. 0) return rc = 0 return end subroutine CFIO_Close !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `CFIO_PutIntAtt` -- Write a user-defined integer attribute ! ! This routine allows the user to define an integer ! attribute in an open CFIO file. ! !#### History !- 1998.07.30 Lucchesi Initial interface design. !- 1998.07.30 Lucchesi Initial coding. !- 1998.09.24 Lucchesi Changed error handling. !- 1998.09.28 Lucchesi Added support for multiple precisions ! subroutine CFIO_PutIntAtt ( fid, name, count, buf, prec, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) name !! Name of attribute integer count !! Number of integers to write integer buf(count) !! Buffer with integer values integer prec !! Desired precision of attribute value: !! 0 = 32 bit !! 1 = 64 bit ! ! !OUTPUT PARAMETERS: ! integer rc !! Error return code: !! rc = 0 all is well !! rc = -12 error determining default precision !! !! NetCDF Errors !! ------------- !! rc = -36 error from NF90_PUT_ATT (global attribute) !! rc = -55 error from NF90_REDEF (enter define mode) !! rc = -56 error from NF90_ENDDEF (exit define mode) ! !------------------------------------------------------------------------- integer(kind=INT32) dummy32 integer(kind=INT64) dummy64 integer i integer(kind=INT32), allocatable :: buf32(:) integer(kind=INT64), allocatable :: buf64(:) rc = NF90_REDEF ( fid ) if (err("PutIntAtt: could not enter define mode",rc,-55) .NE. 0) & return if ( HUGE(dummy32) .EQ. HUGE(i) .AND. prec .EQ. 0 ) then ! -i4 rc = NF90_PUT_ATT ( fid, NF90_GLOBAL, name, buf) ! 32-bit out else if ( HUGE(dummy32) .EQ. HUGE(i) .AND. prec .EQ. 1 ) then ! -i4 allocate ( buf64(count) ) ! 64-bit out do i=1,count buf64(i) = buf(i) enddo rc = NF90_PUT_ATT ( fid, NF90_GLOBAL, name, buf64 ) deallocate (buf64) else if (HUGE(dummy64) .EQ. HUGE(i) .AND. prec .EQ. 0 ) then ! -i8 allocate ( buf32(count) ) ! 32-bit out do i=1,count buf32(i) = buf(i) enddo rc = NF90_PUT_ATT ( fid, NF90_GLOBAL, name, buf32 ) deallocate (buf32) else if (HUGE(dummy64) .EQ. HUGE(i) .AND. prec .EQ. 1 ) then ! -i8 rc = NF90_PUT_ATT ( fid, NF90_GLOBAL, name, buf ) ! 64-bit out else rc = -12 return endif if (err("PutIntAtt: error writing attribute",rc,-36) .NE. 0) & return rc = NF90_ENDDEF(fid) if (err("PutIntAtt: could not exit define mode",rc,-56) .NE. 0) & return rc = 0 return end subroutine CFIO_PutIntAtt !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `CFIO_PutRealAtt` -- Write a user-defined real attribute ! ! This routine allows the user to define a real ! attribute in an open CFIO file. ! !#### History !- 1998.07.30 Lucchesi Initial interface design. !- 1998.07.30 Lucchesi Initial coding. !- 1998.09.24 Lucchesi Changed error handling. !- 1998.09.28 Lucchesi Added support for multiple precisions ! subroutine CFIO_PutRealAtt ( fid, name, count, buf, prec, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) name !! Name of attribute integer count !! Number of integers to write real buf(count) !! Buffer with real values integer prec !! Desired precision of attribute value: !! 0 = 32 bit !! 1 = 64 bit ! ! !OUTPUT PARAMETERS: ! integer rc !! Error return code: !! rc = 0 all is well !! rc = -12 error determining default precision !! !! NetCDF Errors !! ------------- !! rc = -36 error from NF90_PUT_ATT (global attribute) !! rc = -55 error from NF90_REDEF (enter define mode) !! rc = -56 error from NF90_ENDDEF (exit define mode) ! !------------------------------------------------------------------------- real(kind=REAL32) dummy32 real(kind=REAL64) dummy64 real r integer i real(kind=REAL32), allocatable :: buf32(:) real(kind=REAL64), allocatable :: buf64(:) rc = NF90_REDEF ( fid ) if (err("PutRealAtt: could not enter define mode",rc,-55) .NE. 0) & return if (HUGE(dummy32) .EQ. HUGE(r) .AND. prec .EQ. 0) then ! -r4 rc = NF90_PUT_ATT ( fid, NF90_GLOBAL, name, buf ) ! 32-bit out else if (HUGE(dummy32) .EQ. HUGE(r) .AND. prec .EQ. 1) then ! -r4 allocate (buf64(count)) ! 64-bit out do i=1,count buf64(i) = buf(i) enddo rc = NF90_PUT_ATT ( fid, NF90_GLOBAL, name, buf64 ) deallocate (buf64) else if (HUGE(dummy64) .EQ. huge(r) .AND. prec .EQ. 0) then ! -r8 allocate (buf32(count)) ! 32-bit out do i=1,count buf32(i) = buf(i) enddo rc = NF90_PUT_ATT ( fid, NF90_GLOBAL, name, buf32 ) deallocate (buf32) else if (HUGE(dummy64) .EQ. huge(r) .AND. prec .EQ. 1) then ! -r8 rc = NF90_PUT_ATT ( fid, NF90_GLOBAL, name, buf ) ! 64-bit out else rc = -12 return endif if (err("PutRealAtt: error writing attribute",rc,-36) .NE. 0) & return rc = NF90_ENDDEF(fid) if (err("PutRealAtt: could not exit define mode",rc,-56) .NE. 0) & return rc = 0 return end subroutine CFIO_PutRealAtt !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `CFIO_PutCharAtt` -- Write a user-defined character attribute ! ! This routine allows the user to define a character (string) ! attribute in an open CFIO file. ! !#### History !- 1998.07.30 Lucchesi Initial interface design. !- 1998.07.30 Lucchesi Initial coding. !- 1998.09.24 Lucchesi Changed error handling. ! subroutine CFIO_PutCharAtt ( fid, name, count, buf, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) name !! Name of attribute integer count !! Number of characters to write character(len=MLEN) :: buf !! Buffer containing string ! ! !OUTPUT PARAMETERS: ! integer rc !! Error return code: !! rc = 0 all is well !! !! NetCDF Errors !! ------------- !! rc = -36 error from NF90_PUT_ATT (global attribute) !! rc = -55 error from NF90_REDEF (enter define mode) !! rc = -56 error from NF90_ENDDEF (exit define mode) ! !------------------------------------------------------------------------- _UNUSED_DUMMY(count) rc = NF90_REDEF ( fid ) if (err("PutCharAtt: could not enter define mode",rc,-55) .NE. 0) & return rc = NF90_PUT_ATT ( fid, NF90_GLOBAL, name, buf ) if (err("PutCharAtt: error writing attribute",rc,-36) .NE. 0) & return rc = NF90_ENDDEF(fid) if (err("PutCharAtt: could not exit define mode",rc,-56) .NE. 0) & return rc = 0 return end subroutine CFIO_PutCharAtt !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `CFIO_GetAttNames` -- Get global attribute names ! ! This routine allows the user to get the names of ! global attributes. ! !#### History !- 1998.08.05 Lucchesi Initial interface design. !- 1998.08.05 Lucchesi Initial coding. !- 1998.09.24 Lucchesi Changed error handling. ! subroutine CFIO_GetAttNames ( fid, ngatts, aname, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle ! ! !INPUT/OUTPUT PARAMETERS: ! integer ngatts !! Expected number of attributes (input) !! Actual number of attributes (output if rc=-2) ! ! !OUTPUT PARAMETERS: ! character(len=*) aname(ngatts) !! Array of attribute names integer rc !! Error return code: !! rc = 0 all is well !! rc = -10 ngatts is incompatible with file !! rc = -11 character string not long enough !! !! NetCDF Errors !! ------------- !! rc = -48 error from NF90_INQUIRE !! rc = -57 error from NF90_INQ_ATTNAME ! !------------------------------------------------------------------------- integer ngattsFile, i integer nDims,dimSize,recDim ! Check number of attributes against file rc = NF90_INQUIRE (fid,nDims,dimSize,ngattsFile,recdim) if (err("GetAttNames: NF90_INQUIRE failed",rc,-48) .NE. 0) return if (ngattsFile .NE. ngatts) then rc = -10 ngatts = ngattsFile return endif ! Check length of aname string if (LEN(aname(1)) .lt. MAXNCNAM) then print *,'CFIO_GetAttNames: length of aname array must be at ', & 'least ',MAXNCNAM,' bytes.' rc = -11 return endif ! Read global attribute names do i=1,ngatts rc = NF90_INQ_ATTNAME(fid,NF90_GLOBAL,i,aname(i)) if (err("GetAttNames: error reading attribute name",rc,-57) & .NE. 0) return enddo rc = 0 return end subroutine CFIO_GetAttNames !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! ! `CFIO_AttInquire` -- Get information about an attribute ! ! This routine allows the user to get information about ! a global attribute of an open CFIO file. This is most ! useful for determining the number of values stored in an ! attribute. ! ! @note ! The returned integer "type" for 64-bit integer is not supported ! in the current implementation of netCDF/HDF. When a user writes a ! 64-bit integer attribute using PutIntAtt, it is actually saved as ! a 64-bit real by the HDF library. Thus, upon reading the attribute, ! there is no way for HDF/CFIO to distinguish it from a REAL number. ! The user must realize this variable is really an integer and call ! GetIntAtt to read it. Even for a 64-bit integer, type=4 will never ! be returned unless there are changed to HDF/netCDF. !@endnote ! !#### History !- 1998.07.30 Lucchesi Initial interface design. !- 1998.07.30 Lucchesi Initial coding. !- 1998.09.24 Lucchesi Changed error codes, added type assignment. ! subroutine CFIO_AttInquire ( fid, name, type, count, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) name !! Name of attribute ! ! !OUTPUT PARAMETERS: ! integer type !! Code for attribute type !! 0 = integer !! 1 = real !! 2 = character !! 3 = 64-bit real !! 4 = 64-bit integer !! -1 = other integer count !! Number of items (length of array) integer rc !! Error return code: !! rc = 0 all is well !! !! NetCDF Errors !! ------------- !! rc = -58 error from NF90_INQUIRE_ATTRIBUTE ! !------------------------------------------------------------------------- integer nctype rc = NF90_INQUIRE_ATTRIBUTE (fid, NF90_GLOBAL, name, nctype, count) if (err("AttInquire: error reading attribute info",rc,-58) & .NE. 0) return if (nctype .EQ. NF90_INT) then type = 0 elseif (nctype .EQ. NF90_FLOAT) then type = 1 elseif (nctype .EQ. NCCHAR) then type = 2 elseif (nctype .EQ. NF90_DOUBLE) then type = 3 else type = -1 endif rc = 0 return end subroutine CFIO_AttInquire !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `CFIO_GetIntAtt` -- Read a user-defined integer attribute ! ! This routine allows the user to read an integer ! attribute from an open CFIO file. ! !#### History !- 1998.07.30 Lucchesi Initial interface design. !- 1998.07.30 Lucchesi Initial coding. !- 1998.09.29 Lucchesi Changed error handling. Added 64-bit support. ! subroutine CFIO_GetIntAtt ( fid, name, count, buf, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) name !! Name of attribute ! ! !INPUT/OUTPUT PARAMETERS: ! integer count !! On input: Number of items in attribute !! On output: If rc = -1, count will contain !! the correct count of this attribute ! ! !OUTPUT PARAMETERS: ! integer buf(count) !! Buffer with integer values integer rc !! Error return code: !! rc = 0 all is well !! rc = -1 invalid count !! rc = -2 type mismatch !! rc = -12 error determining default precision !! !! NetCDF Errors !! ------------- !! rc = -36 error from NF90_PUT_ATT (global attribute) !! rc = -51 error from NF90_GET_ATT (global attribute) ! !------------------------------------------------------------------------- integer length, type integer i integer(kind=INT32) dummy32 integer(kind=INT64) dummy64 integer(kind=INT32), allocatable :: buf32(:) integer(kind=INT64), allocatable :: buf64(:) rc = NF90_INQUIRE_ATTRIBUTE (fid, NF90_GLOBAL, name, type, length) if (err("GetIntAtt: error reading attribute info",rc,-58) & .NE. 0) return if ( count .NE. length ) then rc = -1 count = length return endif if ( type .NE. NF90_INT .AND. type .NE. NF90_DOUBLE) then rc = -2 return endif if ( HUGE(dummy32) .EQ. HUGE(i)) then if ( type .EQ. NF90_INT ) then ! -i4 32bit rc = NF90_GET_ATT(fid,NF90_GLOBAL,name,buf) else ! type .EQ. NF90_DOUBLE allocate (buf64(count)) ! -i4 64bit rc = NF90_GET_ATT(fid,NF90_GLOBAL,name,buf64) do i=1,count buf(i) = buf64(i) enddo deallocate (buf64) endif else if (HUGE(dummy64) .EQ. HUGE(i)) then if ( type .EQ. NF90_INT ) then allocate (buf32(count)) ! -i8 32bit rc = NF90_GET_ATT(fid,NF90_GLOBAL,name,buf32) do i=1,count buf(i) = buf32(i) enddo deallocate (buf32) else ! type .EQ. NF90_DOUBLE rc = NF90_GET_ATT(fid,NF90_GLOBAL,name,buf) endif else rc = -12 return endif if (err("GetIntAtt: error reading attribute value",rc,-51) & .NE. 0) return rc = 0 return end subroutine CFIO_GetIntAtt !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `CFIO_GetRealAtt` -- Read a user-defined real attribute ! ! This routine allows the user to read a real ! attribute from an open CFIO file. ! !#### History !- 1998.07.30 Lucchesi Initial interface design. !- 1998.07.30 Lucchesi Initial coding. !- 1998.09.29 Lucchesi Changed error handling. Added 64-bit support. !- 1999.08.23 Lucchesi Changed .OR. to .AND. ! subroutine CFIO_GetRealAtt ( fid, name, count, buf, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) name !! Name of attribute ! ! !INPUT/OUTPUT PARAMETERS: ! integer count !! On input: Number of items in attribute !! On output: If rc = -1, count will contain !! the correct number of attributes ! ! !OUTPUT PARAMETERS: ! real buf(count) !! Buffer with real values integer rc !! Error return code: !! rc = 0 all is well !! rc = -1 invalid count !! rc = -2 type mismatch !! rc = -12 error determining default precision !! !! NetCDF Errors !! ------------- !! rc = -36 error from NF90_PUT_ATT (global attribute) !! rc = -51 error from NF90_GET_ATT (global attribute) ! !------------------------------------------------------------------------- integer length, type real r integer i real(kind=REAL32) dummy32 real(kind=REAL64) dummy64 real(kind=REAL32), allocatable :: buf32(:) real(kind=REAL64), allocatable :: buf64(:) rc = NF90_INQUIRE_ATTRIBUTE (fid, NF90_GLOBAL, name, type, length) if (err("GetRealAtt: error reading attribute info",rc,-58) & .NE. 0) return if ( count .NE. length ) then rc = -1 count = length return endif if ( type .NE. NF90_FLOAT .AND. type .NE. NF90_DOUBLE) then rc = -2 return endif if ( HUGE(dummy32) .EQ. HUGE(r)) then if ( type .EQ. NF90_FLOAT ) then ! -r4 32bit rc = NF90_GET_ATT(fid,NF90_GLOBAL,name,buf) else ! type .EQ. NF90_DOUBLE allocate (buf64(count)) ! -r4 64bit rc = NF90_GET_ATT(fid,NF90_GLOBAL,name,buf64) do i=1,count buf(i) = buf64(i) enddo deallocate (buf64) endif else if (HUGE(dummy64) .EQ. HUGE(r)) then if ( type .EQ. NF90_FLOAT ) then allocate (buf32(count)) ! -r8 32bit rc = NF90_GET_ATT(fid,NF90_GLOBAL,name,buf32) do i=1,count buf(i) = buf32(i) enddo deallocate (buf32) else ! type .EQ. NF90_DOUBLE rc = NF90_GET_ATT(fid,NF90_GLOBAL,name,buf) endif else rc = -12 return endif if (err("GetRealAtt: error reading attribute value",rc,-51) & .NE. 0) return rc = 0 return end subroutine CFIO_GetRealAtt !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `CFIO_GetCharAtt` -- Read a user-defined character attribute ! ! This routine allows the user to read a character ! attribute from an open CFIO file. ! !#### History !- 1998.07.30 Lucchesi Initial interface design. !- 1998.07.30 Lucchesi Initial coding. !- 1998.09.29 Lucchesi Changed error handling. ! subroutine CFIO_GetCharAtt ( fid, name, count, buf, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) name !! Name of attribute ! ! !INPUT/OUTPUT PARAMETERS: ! integer count !! On input: Number of items in attribute !! On output: If rc = -1, count will contain !! the correct number of attributes ! ! !OUTPUT PARAMETERS: ! character :: buf(count) !! Buffer with character values ! character(len=MLEN) :: buf !! Buffer with character values integer rc !! Error return code: !! rc = 0 all is well !! rc = -1 invalid count !! rc = -2 type mismatch !! !! NetCDF Errors !! ------------- !! rc = -36 error from NF90_PUT_ATT (global attribute) !! rc = -51 error from NF90_GET_ATT (global attribute) ! !------------------------------------------------------------------------- integer length, type character(len=count) :: chartmp rc = NF90_INQUIRE_ATTRIBUTE (fid, NF90_GLOBAL, name, type, length) if (err("GetCharAtt: error reading attribute info",rc,-58) & .NE. 0) return if ( count .NE. length ) then rc = -1 count = length return endif if ( type .NE. NCCHAR) then rc = -2 return endif rc = NF90_GET_ATT(fid,NF90_GLOBAL,name,chartmp) if (err("GetCharAtt: error reading attribute value",rc,-51) & .NE. 0) return buf = chartmp rc = 0 return end subroutine CFIO_GetCharAtt ! The following function was taken from the book "Numerical Recipes in ! FORTRAN, the art of scientific computing (2nd Ed.), by William H. Press, ! Saul A. Teukolsky, William T. Vetterling, and Brian P. Flannery (Cambridge ! University Press, 1992). INTEGER FUNCTION julday(mm,id,iyyy) INTEGER id,iyyy,mm,IGREG PARAMETER (IGREG=15+31*(10+12*1582)) INTEGER ja,jm,jy jy=iyyy if (jy.eq.0) then print *, 'julday: there is no year zero' return endif if (jy.lt.0) jy=jy+1 if (mm.gt.2) then jm=mm+1 else jy=jy-1 jm=mm+13 endif julday=int(365.25*jy)+int(30.6001*jm)+id+1720995 if (id+31*(mm+12*iyyy).ge.IGREG) then ja=int(0.01*jy) julday=julday+2-ja+int(0.25*ja) endif return END function julday !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !BOI ! ! !TITLE: Calculate difference in seconds between two times. ! ! !AUTHORS: Rob Lucchesi ! ! !AFFILIATION: Data Assimilation Office, NASA/GSFC, Greenbelt, MD 20771 ! ! !DATE: October 17, 1997 ! !EOI !------------------------------------------------------------------------- !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `DiffDate` --- Calculates the number of seconds between two times. ! ! This function returns the number of seconds between two ! times. Each time is specified with two integers, one ! representing a date in the format YYYYMMDD and one ! representing a time in the format HHMMSS. This function ! determines the Julian day of each date using the "julday" ! function from the book "Numerical Recipes in FORTRAN, the ! art of scientific computing (2nd Ed.), by William H. Press, ! Saul A. Teukolsky, William T. Vetterling, and Brian P. ! Flannery (Cambridge University Press, 1992). The difference ! between the two times is then calculated and returned. The ! times need not be in chronological order as the function returns ! the abs value. -1 is returned in the event of an error. ! ! RETURNED VALUE: ! ! Integer function returns number of seconds between the ! the times given as input. -1 is returned in the event ! of an error. ! !#### History !- 17Oct97 Lucchesi Initial version. !- 2010.05.11 Lucchesi Integer for julian seconds changed to 64-bit. StartDate ! constant no longer needed. ! integer(kind=INT64) function DiffDate (yyyymmhh_1,hhmmss_1,yyyymmhh_2,hhmmss_2) implicit none ! ! !INPUT PARAMETERS: ! integer yyyymmhh_1 ! First date in YYYYYMMDD format integer hhmmss_1 ! First time in HHMMSS format integer yyyymmhh_2 ! Second date in YYYYMMDD format integer hhmmss_2 ! Second time in HHMMSS format ! !------------------------------------------------------------------------- integer year1,mon1,day1,hour1,min1,sec1 integer year2,mon2,day2,hour2,min2,sec2 integer(INT64) julian1, julian2, julsec1, julsec2 ! Error checking. !rl if (yyyymmhh_1 .lt. 19000000 .or. yyyymmhh_1 .gt. 21000000 ) then !rl DiffDate=-1 !rl return !rl endif !rl if (yyyymmhh_2 .lt. 19000000 .or. yyyymmhh_2 .gt. 21000000 ) then !rl DiffDate=-1 !rl return !rl endif if (hhmmss_1 .lt. 0 .or. hhmmss_1 .ge. 240000 ) then DiffDate=-1 return endif if (hhmmss_2 .lt. 0 .or. hhmmss_2 .ge. 240000 ) then DiffDate=-1 return endif ! Convert Date/Time strings to integer variables. !ams write (dateString, 200) yyyymmhh_1 !ams 200 format (I8) !ams read (dateString, 201) !ams 201 format (I4,2I2) !ams write (dateString, 200) yyyymmhh_2 !ams read (dateString, 201) year2, mon2, day2 call CFIO_parseIntTime ( yyyymmhh_1, year1, mon1, day1 ) call CFIO_parseIntTime ( yyyymmhh_2, year2, mon2, day2 ) !ams write (dateString, 202) hhmmss_1 !ams 202 format (I6) !ams read (dateString, 203) hour1, min1, sec1 !ams 203 format (3I2) !ams write (dateString, 202) hhmmss_2 !ams read (dateString, 203) hour2, min2, sec2 call CFIO_parseIntTime ( hhmmss_1, hour1, min1, sec1 ) call CFIO_parseIntTime ( hhmmss_2, hour2, min2, sec2 ) ! Get Julian Days and subtract off a constant (Julian days since 7/14/66) julian1 = julday (mon1, day1, year1) julian2 = julday (mon2, day2, year2) ! Calculcate Julian seconds julsec1 = (julian1-1)*86400 + hour1*3600 + min1*60 + sec1 julsec2 = (julian2-1)*86400 + hour2*3600 + min2*60 + sec2 !!! DiffDate = iabs (julsec2 - julsec1) DiffDate = julsec2 - julsec1 return end function DiffDate SUBROUTINE CALDAT (JULIAN,MM,ID,IYYY) ! !C ROUTINE CONVERTS JULIAN DAY TO MONTH, DAY, & YEAR. !C THIS CODE IS LIFTED FROM THE BOOK: !C W.H. PRESS ET AL., NUMERICAL RECIPES, CAMBRIDGE UNIV. PRESS, 1986. !C THE ONLY MODIFICATION IS THAT REAL ARITHMETIC IS DONE IN R*8. !C !C THE ROUTINE OUTPUTS THE MONTH, DAY, AND YEAR ON WHICH THE !C SPECIFIED JULIAN DAY STARTED AT NOON. !C !C TO CONVERT MODIFIED JULIAN DAY, CALL THIS ROUTINE WITH !C JULIAN = MJD + 2400001 !C integer(INT64) JULIAN integer IGREG integer JALPHA, JA, JB, JC, JD, JE, ID, MM, IYYY PARAMETER (IGREG=2299161) IF (JULIAN.GE.IGREG) THEN JALPHA=INT((DBLE(JULIAN-1867216)-0.25D0)/36524.25D0) JA=JULIAN+1+JALPHA-INT(0.25D0*DBLE(JALPHA)) ELSE JA=JULIAN ENDIF JB=JA+1524 JC=INT(6680.D0+(DBLE(JB-2439870)-122.1D0)/365.25D0) JD=365*JC+INT(0.25D0*DBLE(JC)) JE=INT(DBLE(JB-JD)/30.6001D0) ID=JB-JD-INT(30.6001D0*DBLE(JE)) MM=JE-1 IF (MM.GT.12) MM=MM-12 IYYY=JC-4715 IF (MM.GT.2) IYYY=IYYY-1 IF (IYYY.LE.0) IYYY=IYYY-1 RETURN END subroutine CALDAT integer function err ( outstring, iret, ec ) character(len=*) outstring integer ec, iret if (iret .EQ. 0) then err = iret else print *,"CFIO_",outstring iret = ec err = ec endif return end function err !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `ParseTimeUnits` -- Parse the COARDS time units string ! ! This subroutine takes as input the COARDS metadata time ! units string and parses out the date included in the string. ! This date typically represents the first time in a COARDS ! HDF files. ! !#### History !- 2000.10.18 Lucchesi Initial coding. !- 2001.08.13 Lucchesi Modified to better parse time:units string (needed for lats4d support) ! subroutine ParseTimeUnits ( TimeUnits, year, month, day, hour, min, sec, rc ) ! ! !USES: ! implicit none ! ! !INPUT PARAMETERS: ! character(len=MAXCHR) TimeUnits !! Units metadata string from the Time coord var ! ! !OUTPUT PARAMETERS: ! integer year !! 4-digit year integer month !! month integer day !! day integer hour !! hour integer min !! minute integer sec !! second integer rc !! return code !! 0 = no error !! -1 = problem parsing string ! !------------------------------------------------------------------------- ! Local variables integer ypos(2), mpos(2), dpos(2), hpos(2), spos(2) integer strlen integer firstdash, lastdash integer firstcolon, lastcolon integer lastspace strlen = LEN_TRIM (TimeUnits) firstdash = index(TimeUnits, '-') lastdash = index(TimeUnits, '-', BACK=.TRUE.) if (firstdash .LE. 0 .OR. lastdash .LE. 0) then rc = -1 return endif ypos(2) = firstdash - 1 mpos(1) = firstdash + 1 ypos(1) = ypos(2) - 4 mpos(2) = lastdash - 1 dpos(1) = lastdash + 1 dpos(2) = dpos(1) + 2 read ( TimeUnits(ypos(1):ypos(2)), * ) year read ( TimeUnits(mpos(1):mpos(2)), * ) month read ( TimeUnits(dpos(1):dpos(2)), * ) day firstcolon = index(TimeUnits, ':') if (firstcolon .LE. 0) then ! If no colons, check for hour. ! Logic below assumes a null character or something else is after the hour ! if we do not find a null character add one so that it correctly parses time if (TimeUnits(strlen:strlen) /= char(0)) then TimeUnits = trim(TimeUnits)//char(0) strlen=len_trim(TimeUnits) endif lastspace = index(TRIM(TimeUnits), ' ', BACK=.TRUE.) if ((strlen-lastspace).eq.2 .or. (strlen-lastspace).eq.3) then hpos(1) = lastspace+1 hpos(2) = strlen-1 read (TimeUnits(hpos(1):hpos(2)), * ) hour min = 0 sec = 0 else print *, 'ParseTimeUnits: Assuming a starting time of 00z' hour = 0 min = 0 sec = 0 endif else hpos(1) = firstcolon - 2 hpos(2) = firstcolon - 1 lastcolon = index(TimeUnits, ':', BACK=.TRUE.) if ( lastcolon .EQ. firstcolon ) then mpos(1) = firstcolon + 1 mpos(2) = firstcolon + 2 read (TimeUnits(hpos(1):hpos(2)), * ) hour read (TimeUnits(mpos(1):mpos(2)), * ) min sec = 0 else mpos(1) = firstcolon + 1 mpos(2) = lastcolon - 1 spos(1) = lastcolon + 1 spos(2) = lastcolon + 2 read (TimeUnits(hpos(1):hpos(2)), * ) hour read (TimeUnits(mpos(1):mpos(2)), * ) min read (TimeUnits(spos(1):spos(2)), * ) sec endif endif rc = 0 return end subroutine ParseTimeUnits !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `CFIO_SPutVar` -- Write a variable to the file ! ! This routine is used to write a variable to an open CFIO ! stream. Multiple vertical levels can be written at one ! time provided they are contiguous in memory. Date and time ! must be consistent with the time increment and the starting ! date/time as defined in CFIO\_Create. Times must fall on ! minute boundaries to allow GrADS to work. Error checking is ! done for dimensions that are out of bounds. ! !#### History !- 1997.10.13 da Silva/Lucchesi Initial interface design. !- 1998.02.10 Lucchesi Added support for applications running with 64-bit reals. !- 1998.03.30 Lucchesi Documentation expanded. Clean-up of code. !- 1998.07.02 Lucchesi Replaced vid with vname in argument list & made related mods to code. !- 1998.09.24 Lucchesi Changed err codes, removed DIM_CHECK if-def !- 1998.10.27 Lucchesi Added support for packing and range checks !- 1998.12.15 Lucchesi Added support for skipping times (allTimes) !- 1999.01.04 Lucchesi Fixed bug in skipping times (allTimes)/also changed variable initialization. !- 1999.07.13 Lucchesi Changes for REAL or INT time dimension ! subroutine CFIO_SPutVar ( fid, vname, yyyymmdd, hhmmss, & im, jm, kbeg, kount, grid, & rc ) ! ! !USES: Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) vname !! Variable name integer yyyymmdd !! Year-month-day, e.g., 19971003 integer hhmmss !! Hour-minute-second, e.g., 120000 integer im !! size of longitudinal dimension integer jm !! size of latitudinal dimension integer kbeg !! first level to write; if 2-D grid !! use kbeg = 0. integer kount !! number of levels to write real grid(im,kount) !! Gridded data to write at this time ! !OUTPUT PARAMETERS: integer rc !! Error return code: !! rc = 0 all is well !! rc = -2 time is inconsistent with increment !! rc = -3 number of levels is incompatible with file !! rc = -4 im is incompatible with file !! rc = -5 jm is incompatible with file !! rc = -6 time must fall on a minute boundary !! rc = -7 error in diffdate !! rc = -12 error determining default precision !! rc = -13 error determining variable type !! rc = -15 data outside of valid range !! rc = -16 data outside of packing range !! rc = -17 data outside of pack and valid range !! !! NetCDF Errors !! ------------- !! rc = -38 error from NF90_PUT_VAR (dimension variable) !! rc = -40 error from NF90_INQ_VARID !! rc = -41 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (lat or lon) !! rc = -42 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (lev) !! rc = -43 error from NF90_INQ_VARID (time variable) !! rc = -44 error from NF90_GET_ATT (time attribute) !! rc = -45 error from NF90_PUT_VAR !! rc = -46 error from NF90_GET_VAR !! rc = -52 error from NF90_INQUIRE_VARIABLE !! rc = -53 error from NF90_GET_ATT ! !------------------------------------------------------------------------- integer timeid, timeDimId, dimSize, timeType character(len=MAXCHR) dimName integer corner(3), edges(3) integer vid integer(INT64) seconds integer timeIndex integer minutes ! added as a work-around integer i, k integer begDate, begTime, timInc integer hour,min,sec,incSecs integer, allocatable :: allTimes(:) integer fillTime ! Variables for dealing with precision real(kind=REAL32), allocatable :: grid_32(:,:) real(kind=REAL64), allocatable :: grid_64(:,:) real(kind=REAL32) dummy32 real(kind=REAL64) dummy64 real dummy ! Variables for NF90_INQUIRE_VARIABLE character(len=MAXCHR) varName integer type, nvDims, vdims(MAXVDIMS), nvAtts ! Variables for packing and range checking integer(kind=INT16), allocatable :: grid_16(:,:) real(kind=REAL32), allocatable :: fminutes_32(:) real(kind=REAL32) high_32, low_32, amiss_32 real(kind=REAL32) scale_32, offset_32 logical outRange logical outPRange _UNUSED_DUMMY(jm) ! Variable initialization outRange = .FALSE. outPRange = .FALSE. ! Check to make sure max string lengths are large enough. NetCDF defines ! MAXNCNAM, but it can't be used in a character(len=MAXNCNAM) statement. if (MAXCHR .LT. MAXNCNAM) then print *, 'CFIO_PutVar warning: MAXNCNAM is larger than ', & 'dimName array size.' endif ! Determine NetCDF variable ID. rc = NF90_INQ_VARID (fid, vname, vid) if (err("PutVar: variable not defined",rc,-40) .NE. 0) return ! Determine number of seconds since starting date/time. rc = NF90_INQ_VARID (fid, 'time', timeId) if (err("PutVar: time not defined",rc,-43) .NE. 0) return rc = NF90_GET_ATT(fid,timeId,'begin_date',begDate) if (err("PutVar: missing begin_date",rc,-44) .NE. 0) return rc = NF90_GET_ATT(fid,timeId,'begin_time',begTime) if (err("PutVar: missing begin_time",rc,-44) .NE. 0) return seconds = DiffDate (begDate, begTime, yyyymmdd, hhmmss) if (seconds .lt. 0) then print *, 'CFIO_PutVar: Error code from diffdate. Problem with', & ' date/time.' rc = -7 return endif if ( MOD (int(seconds),60) .eq. 0 ) then minutes = seconds / 60 else print *, 'CFIO_PutVar: Currently, times must fall on minute ',& 'boundaries.' rc = -6 return endif ! Confirm that this time is consistent with the starting time coupled with ! the time increment. rc = NF90_GET_ATT(fid,timeId,'time_increment',timInc) if (err("PutVar: missing time increment",rc,-44) .NE. 0) return ! Convert time increment to seconds. !ams write (strBuf,203) timinc !ams 203 format (I6) !ams read (strBuf,204) hour, min, sec !ams 204 format (3I2) call CFIO_parseIntTime ( timinc, hour, min, sec ) incSecs = hour*3600 + min*60 + sec if ( MOD (int(seconds), incSecs) .ne. 0 ) then print *, 'CFIO_putvar: Absolute time of ',seconds,' not ',& 'possible with an interval of ',incSecs rc = -2 return else timeIndex = seconds/incSecs + 1 endif ! Load starting indicies. if ( kbeg .eq. 0 ) then corner(1)=1 corner(2)=timeIndex edges(1)=im edges(2)=1 else corner(1)=1 corner(2)=kbeg corner(3)=timeIndex edges(1)=im edges(2)=kount edges(3)=1 endif ! Check variable against valid range. rc = NF90_GET_ATT(fid,vid,'vmin',low_32) if (err("PutVar: can't get vmin",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'vmax',high_32) if (err("PutVar: can't get vmax",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'fmissing_value',amiss_32) if (err("PutVar: can't get fmissing_value",rc,-53) .NE. 0) return if (abs(low_32) .NE. amiss_32 .OR. high_32 .NE. amiss_32) then do k=1,kount do i=1,im if (grid(i,k) .GT. high_32 .OR. grid(i,k) .LT. & low_32) then outRange = .TRUE. goto 100 endif enddo enddo 100 continue endif ! Determine if we are writing single- or double-precision. rc = NF90_INQUIRE_VARIABLE (fid, vid, varName, type, nvDims, vDims, nvAtts) if (err("PutVar: error in variable inquire",rc,-52) .NE. 0) return ! Write variable in the appropriate precision. if (HUGE(dummy) .EQ. HUGE(dummy32)) then ! -r4 if (type .EQ. NF90_FLOAT) then ! 32-bit rc = NF90_PUT_VAR(fid,vid,grid,corner,edges) else if (type .EQ. NF90_DOUBLE) then ! 64-bit allocate (grid_64(im,kount)) do k=1,kount do i=1,im grid_64(i,k) = grid(i,k) enddo enddo rc = NF90_PUT_VAR(fid,vid,grid_64,corner,edges) deallocate (grid_64) else if (type .EQ. NF90_SHORT) then rc = NF90_GET_ATT(fid,vid,'packmax',high_32) if (err("PutVar: error getting packmax",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'packmin',low_32) if (err("PutVar: error getting packmin",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'scale_factor',scale_32) if (err("PutVar: error getting scale",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'add_offset',offset_32) if (err("PutVar: error getting offset",rc,-53) .NE. 0) return allocate (grid_16(im,kount)) do k=1,kount do i=1,im if ( grid(i,k) .LT. low_32 .OR. grid(i,k) .GT. & high_32) then grid_16(i,k) = PACK_FILL outPRange = .TRUE. else grid_16(i,k) = (grid(i,k) - offset_32)/scale_32 endif enddo enddo rc = NF90_PUT_VAR(fid,vid,grid_16,corner,edges) deallocate (grid_16) else rc = -13 return endif else if (HUGE(dummy) .EQ. HUGE(dummy64)) then ! -r8 if (type .EQ. NF90_FLOAT) then ! 32-bit allocate (grid_32(im,kount)) do k=1,kount do i=1,im grid_32(i,k) = grid(i,k) enddo enddo rc = NF90_PUT_VAR(fid, vid, grid_32, corner, edges) deallocate (grid_32) else if (type .EQ. NF90_DOUBLE) then ! 64-bit rc = NF90_PUT_VAR(fid,vid,grid,corner,edges) else if (type .EQ. NF90_SHORT) then rc = NF90_GET_ATT(fid,vid,'packmax',high_32) if (err("PutVar: error getting packmax",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'packmin',low_32) if (err("PutVar: error getting packmin",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'scale_factor',scale_32) if (err("PutVar: error getting scale",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'add_offset',offset_32) if (err("PutVar: error getting offset",rc,-53) .NE. 0) return allocate (grid_16(im,kount)) do k=1,kount do i=1,im if ( grid(i,k) .LT. low_32 .OR. grid(i,k) .GT. & high_32) then grid_16(i,k) = PACK_FILL outPRange = .TRUE. else grid_16(i,k) = (grid(i,k) - offset_32)/scale_32 endif enddo enddo rc = NF90_PUT_VAR(fid,vid,grid_16,corner,edges) deallocate (grid_16) else rc = -13 return endif else rc = -12 return endif if (err("PutVar: error writing variable",rc,-45) .NE. 0) return ! Read time dimension scale and fill all values up to the current time. ! This will insure missing times are defined with the proper time value. rc = NF90_INQ_DIMID(fid, 'time', timeDimId) rc = NF90_INQUIRE_DIMENSION (fid, timeDimId, dimName, dimSize) dimSize = dimSize - 1 ! We've already written the ! the new time. allocate ( allTimes (MAX(timeIndex,dimSize)) ) allocate ( fminutes_32 (MAX(timeIndex,dimSize)) ) if (dimSize .GT. 0) then ! Depending on the version of CFIO used to write the file, the Time ! dimension variable can either be floating point or integer. corner(1)=1 edges(1)=dimSize rc = NF90_INQUIRE_VARIABLE (fid,timeId,dimName,timeType,nvDims,vDims,nvAtts) if (timeType .EQ. NF90_FLOAT) then rc = NF90_GET_VAR(fid,timeId,fminutes_32,corner,edges) do i=1,dimSize allTimes(i) = INT(fminutes_32(i)) enddo else if (timeType .EQ. NF90_INT) then rc = NF90_GET_VAR(fid,timeId,allTimes,corner,edges) endif if (err("SPutVar: error reading times from file",rc,-46) .NE. 0) & return endif ! This loop fills the time dimension scale based on the time increment ! specified in CFIO_Create. If CFIO ever changes to support variable ! time increments, this code MUST be changed. do i=1,timeIndex-1 fillTime = (i-1) * incSecs/60 allTimes(i) = fillTime enddo allTimes(timeIndex) = minutes ! Write filled time array to file. corner(1)=1 edges(1)=timeIndex if (timeType .EQ. NF90_FLOAT) then do i=1,timeIndex fminutes_32(i) = INT(allTimes(i)) enddo rc = NF90_PUT_VAR(fid,timeId,fminutes_32,corner,edges) else if (timeType .EQ. NF90_INT) then rc = NF90_PUT_VAR(fid,timeId,allTimes,corner,edges) endif if (err("PutVar: error writing time",rc,-38) .NE. 0) return if (outRange .AND. outPRange) then rc = -17 else if (outPRange) then rc = -16 else if (outRange) then rc = -15 else rc = 0 endif deallocate ( allTimes ) deallocate ( fminutes_32 ) return end subroutine CFIO_SPutVar !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `CFIO_SGetVar` -- Read a variable from the file ! ! This routine will read one or more levels of "vname" ! into the buffer passed in as "grid." "fid" is the file ! handle returned by Cfio_open. ! !#### History !- 1997.10.13 da Silva/Lucchesi Initial interface design. !- 1998.07.07 Lucchesi Combined two GetVar routines into this one. !- 1998.09.24 Lucchesi Updated error codes. !- 1999.06.21 Lucchesi Bug fixed. Unable to read HDF-EOS files because was still looking for "lon" and "lat" !- 1999.06.21 Lucchesi Added a check for time too early. !- 1999.11.02 da Silva Made LATS4D compatible. !- 2006.06.08 byin Added cyclic option ! subroutine CFIO_SGetVar ( fid, vname, yyyymmdd, hhmmss,& im, jm, kbeg, kount, lm, grid, cyclic, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) vname !! Variable name integer yyyymmdd !! Year-month-day, e.g., 19971003 integer hhmmss !! Hour-minute-second, e.g., 120000 integer im !! size of longitudinal dimension integer jm !! size of latitudinal dimension integer kbeg !! first level to read; if 2-D grid !! set kbeg = 0. integer kount !! number of levels to read integer lm !! number of time steps logical cyclic !! input file is cyclic or not ! ! !OUTPUT PARAMETERS: ! real grid(im,kount) !! Gridded data read for this time integer rc !! Error return code: !! rc = 0 all is well !! rc = -2 time is inconsistent with increment !! rc = -3 number of levels is incompatible with file !! rc = -4 im is incompatible with file !! rc = -5 jm is incompatible with file !! rc = -6 time must fall on a minute boundary !! rc = -7 error in diffdate !! rc = -12 error determining default precision !! rc = -13 error determining variable type !! rc = -19 unable to identify coordinate variable !! !! NetCDF Errors !! ------------- !! rc = -38 error from NF90_PUT_VAR (dimension variable) !! rc = -40 error from NF90_INQ_VARID !! rc = -41 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (lat or lon) !! rc = -42 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (lev) !! rc = -43 error from NF90_INQ_VARID (time variable) !! rc = -44 error from NF90_GET_ATT (time attribute) !! rc = -46 error from NF90_GET_VAR !! rc = -48 error from NF90_INQUIRE !! rc = -52 error from NF90_INQUIRE_VARIABLE ! !------------------------------------------------------------------------- integer begDate, begTime, seconds integer corner(3), edges(3), timeIndex integer vid integer i,j,k !integer incSecs integer(INT64), allocatable :: incVec(:) ! Variables for working with dimensions character(len=MAXCHR) dimName character(len=MAXCHR) stdName character(len=MAXCHR) dimUnits character(len=MAXCHR) varName integer dimSize, dimId integer nDims,nvars,ngatts integer varType, myIndex !integer timeShift ! Variables for dealing with precision real(kind=REAL32), allocatable :: grid_32(:,:) real(kind=REAL64), allocatable :: grid_64(:,:) real(kind=REAL32) dummy32 real(kind=REAL64) dummy64 real dummy ! Variables for NF90_INQUIRE_VARIABLE integer type, nvDims, vdims(MAXVDIMS), nvAtts ! Variables for packing integer(kind=INT16), allocatable :: grid_16(:,:) integer(kind=INT16) amiss_16 real(kind=REAL32) amiss_32 real(kind=REAL32) scale_32, offset_32 ! Check to make sure max string lengths are large enough. NetCDF defines ! MAXNCNAM, but it can't be used in a character(len=MAXNCNAM) statement. if (MAXCHR .LT. MAXNCNAM) then print *, 'CFIO_GetVar warning: MAXNCNAM is larger than ', & 'dimName array size.' endif ! Get basic information from file. rc = NF90_INQUIRE (fid, nDims, nvars, ngatts, dimId) if (err("DimInqure: NF90_INQUIRE failed",rc,-48) .NE. 0)return ! Subtract dimension variables from the variable count. do i=1,nvars rc = NF90_INQUIRE_VARIABLE (fid,i,varName,varType,nvDims,vDims,nvAtts) if (err("DimInquire: variable inquire error",rc,-52) .NE. 0) & return if (nvDims .EQ. 1) then nvars = nvars - 1 endif enddo ! Extract dimension information do i=1,nDims rc = NF90_INQUIRE_DIMENSION (fid, i, dimName, dimSize) if (err("DimInqure: can't get dim info",rc,-41) .NE. 0) return if (index(dimName,'station') .gt. 0) then im = dimSize jm = dimSize cycle end if rc = NF90_INQ_VARID (fid, dimName, dimId) if (err("DimInqure: NF90_INQ_VARID failed",rc,-40) .NE. 0) return ! If it has the standard_name attribute, use that instead rc = NF90_GET_ATT(fid,dimId,'standard_name',stdName) if (rc .ne. 0) stdName = Trim(dimName) rc = NF90_GET_ATT(fid,dimId,'units',dimUnits) if (err("DimInqure: could not get units for dimension",rc,-53) & .NE. 0) return myIndex = IdentifyDim (stdName, dimUnits) if ( myIndex .EQ. 0 ) then if (dimSize .ne. im) then rc = -4 im = dimSize ! return endif else if ( myIndex .EQ. 1 ) then if (dimSize .ne. jm) then rc = -5 jm = dimSize ! return endif else if ( myIndex .EQ. 2 ) then if (kount .gt. dimSize) then rc = -3 ! return endif else if ( myIndex .EQ. 3 ) then else print *, 'CFIO_GetVar: Coordinate variable ', & TRIM(dimName),' with units of ',TRIM(dimUnits),& ' is not understood.' rc = -19 ! return endif enddo ! Determine NetCDF variable ID. rc = NF90_INQ_VARID (fid, vname, vid) if (err("GetVar: variable not defined",rc,-40) .NE. 0) return ! Get beginning time & date. Calculate offset seconds from start. !ams rc = NF90_GET_ATT(fid,timeId,'begin_date',begDate) !ams if (err("GetVar: missing begin_date",rc,-44) .NE. 0) return !ams rc = NF90_GET_ATT(fid,timeId,'begin_time',begTime) !ams if (err("GetVar: missing begin_time",rc,-44) .NE. 0) return allocate(incVec(lm)) call GetDateTimeVec( fID, begDate, begTime, incVec, rc ) if (err("GetVar: could not determine begin_date/begin_time",rc,-44)& .NE. 0) return seconds = DiffDate (begDate, begTime, yyyymmdd, hhmmss) ! Make sure input time are valid, if time is not periodic !ams print *, '+++ incSecs, begDate, begTime: ', incsecs, begDate, begTime !ams print *, '+++ seconds, yyyymmdd, hhmmss: ', seconds, yyyymmdd, hhmmss if ( .not. cyclic ) then if (seconds .LT. 0) then print *, 'CFIO_SGetVar: Error code from diffdate. Problem with', & ' date/time.' rc = -7 return endif if (yyyymmdd .LT. begDate .OR. (begDate .EQ. yyyymmdd .AND. & hhmmss .LT. begTime) ) then print *, 'CFIO_GetVar: Requested time earlier than first time.' rc = -7 return endif end if if ( MOD (seconds,60) .eq. 0 ) then else print *, 'CFIO_GetVar: Currently, times must fall on minute ',& 'boundaries.' rc = -6 return endif ! Determine the time index from the time vector. timeIndex=-1 j=0 Do While (timeIndex.lt.1) j=j+1 do i=1,lm if ((seconds.le.incVec(i)).and.(timeIndex.lt.1)) timeIndex=i end do if (timeIndex.lt.1) then if (cyclic) then seconds = seconds - (incVec(2)-incVec(1))*lm else j=j+1000 end if end if if (j.gt.1000) then Write(*,'(a,L1,a,I0.6)') 'CFIO_Sgetvar: failed to find a valid sample (C', cyclic, & '). Total iterations: ', j rc = -20 return end if End Do deallocate(incVec) !ams print *, '+++ Time Index, timeShift: ', timeIndex, timeShift ! Load starting indices. if ( kbeg .eq. 0 ) then corner(1)=1 corner(2)=timeIndex edges(1)=im edges(2)=1 else corner(1)=1 corner(2)=kbeg corner(3)=timeIndex edges(1)=im edges(2)=kount edges(3)=1 endif ! Determine data type. rc = NF90_INQUIRE_VARIABLE (fid, vid, varName, type, nvDims, vDims, nvAtts) if (err("GetVar: error in variable inquire",rc,-52) .NE. 0) return ! Read variable in the appropriate precision. if (HUGE(dummy) .EQ. HUGE(dummy32)) then ! -r4 if (type .EQ. NF90_FLOAT) then ! 32-bit rc = NF90_GET_VAR(fid,vid,grid,corner,edges) else if (type .EQ. NF90_DOUBLE) then ! 64-bit allocate (grid_64(im,kount)) rc = NF90_GET_VAR(fid,vid,grid_64,corner,edges) do k=1,kount do i=1,im grid(i,k) = grid_64(i,k) enddo enddo deallocate (grid_64) else if (type .EQ. NF90_SHORT) then rc = NF90_GET_ATT(fid,vid,'scale_factor',scale_32) if (err("GetVar: error getting scale",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'add_offset',offset_32) if (err("GetVar: error getting offset",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'missing_value',amiss_16) if (err("GetVar: error getting missing",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'fmissing_value',amiss_32) if (err("GetVar: error getting fmissing",rc,-53) .NE. 0) return allocate (grid_16(im,kount)) rc = NF90_GET_VAR(fid,vid,grid_16,corner,edges) do k=1,kount do i=1,im if ( grid_16(i,k) .EQ. amiss_16 ) then grid(i,k) = amiss_32 else grid(i,k) = scale_32*grid_16(i,k) + offset_32 endif enddo enddo deallocate (grid_16) else rc = -13 return endif else if (HUGE(dummy) .EQ. HUGE(dummy64)) then ! -r8 if (type .EQ. NF90_FLOAT) then ! 32-bit allocate (grid_32(im,kount)) ! print *, "im,kount, varName,rc: ",im,kount,trim(varName), rc ! print *, "corner: ",corner ! print *, "edges: ", edges rc = NF90_GET_VAR(fid,vid,grid_32,corner,edges) ! print *, "ts: ",grid_32 do k=1,kount do i=1,im grid(i,k) = grid_32(i,k) enddo enddo deallocate (grid_32) elseif (type .EQ. NF90_DOUBLE) then ! 64-bit rc= NF90_GET_VAR(fid,vid,grid,corner,edges) else if (type .EQ. NF90_SHORT) then rc = NF90_GET_ATT(fid,vid,'scale_factor',scale_32) if (err("GetVar: error getting scale",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'add_offset',offset_32) if (err("GetVar: error getting offset",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'missing_value',amiss_16) if (err("GetVar: error getting missing",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'fmissing_value',amiss_32) if (err("GetVar: error getting fmissing",rc,-53) .NE. 0) return allocate (grid_16(im,kount)) rc = NF90_GET_VAR(fid,vid,grid_16,corner,edges) do k=1,kount do i=1,im if ( grid_16(i,k) .EQ. amiss_16 ) then grid(i,k) = amiss_32 else grid(i,k) = scale_32*grid_16(i,k) + offset_32 endif enddo enddo deallocate (grid_16) else rc = -13 return endif else rc = -12 return endif if (err("GetVar: error reading variable",rc,-46) .NE. 0) return rc = 0 return end subroutine CFIO_SGetVar !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> !`CFIO_GetVar` -- Read a variable from the file ! ! This routine will read one or more levels of "vname" ! into the buffer passed in as "grid." "fid" is the file ! handle returned by Cfio\_open. ! !#### History !- 1997.10.13 da Silva/Lucchesi Initial interface design. !- 1998.07.07 Lucchesi Combined two GetVar routines into this one. !- 1998.09.24 Lucchesi Updated error codes. !- 1999.06.21 Lucchesi Bug fixed. Unable to read HDF-EOS files because was still looking for "lon" and "lat" !- 1999.06.21 Lucchesi Added a check for time too early. !- 1999.11.02 da Silva Made LATS4D compatible. !- 2006.06.08 byin Added cyclic option ! subroutine CFIO_GetVar ( fid, vname, yyyymmdd, hhmmss,& im, jm, kbeg, kount, lm, grid, cyclic, & useFaceDim, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) vname !! Variable name integer yyyymmdd !! Year-month-day, e.g., 19971003 integer hhmmss !! Hour-minute-second, e.g., 120000 integer im !! size of longitudinal dimension integer jm !! size of latitudinal dimension integer kbeg !! first level to read; if 2-D grid !! set kbeg = 0. integer kount !! number of levels to read integer lm !! number of time steps logical cyclic !! input file is cyclic or not logical :: useFaceDim ! ! !OUTPUT PARAMETERS: ! real grid(im,jm,kount) !! Gridded data read for this time integer rc !! Error return code: !! rc = 0 all is well !! rc = -2 time is inconsistent with increment !! rc = -3 number of levels is incompatible with file !! rc = -4 im is incompatible with file !! rc = -5 jm is incompatible with file !! rc = -6 time must fall on a minute boundary !! rc = -7 error in diffdate !! rc = -12 error determining default precision !! rc = -13 error determining variable type !! rc = -19 unable to identify coordinate variable !! !! NetCDF Errors !! ------------- !! rc = -38 error from NF90_PUT_VAR (dimension variable) !! rc = -40 error from NF90_INQ_VARID !! rc = -41 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (lat or lon) !! rc = -42 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (lev) !! rc = -43 error from NF90_INQ_VARID (time variable) !! rc = -44 error from NF90_GET_ATT (time attribute) !! rc = -46 error from NF90_GET_VAR !! rc = -48 error from NF90_INQUIRE !! rc = -52 error from NF90_INQUIRE_VARIABLE ! !------------------------------------------------------------------------- integer begDate, begTime, seconds !integer timeShift integer corner(5), edges(5), timeIndex integer vid integer i,j,k !integer incSecs integer(INT64), allocatable :: incVec(:) ! Variables for working with dimensions character(len=MAXCHR) dimName character(len=MAXCHR) dimUnits character(len=MAXCHR) varName integer dimSize, dimId integer nDims,nvars,ngatts integer varType ! Variables for dealing with precision real(kind=REAL32), allocatable :: grid_32(:,:,:) real(kind=REAL64), allocatable :: grid_64(:,:,:) real(kind=REAL32) dummy32 real(kind=REAL64) dummy64 real dummy ! Variables for NF90_INQUIRE_VARIABLE integer type, nvDims, vdims(MAXVDIMS), nvAtts ! Variables for packing integer(kind=INT16), allocatable :: grid_16(:,:,:) integer(kind=INT16) amiss_16 real(kind=REAL32) amiss_32 real(kind=REAL32) scale_32, offset_32 ! Initialize these, just in case corner = 1 edges = 1 ! Check to make sure max string lengths are large enough. NetCDF defines ! MAXNCNAM, but it can't be used in a character(len=MAXNCNAM) statement. if (MAXCHR .LT. MAXNCNAM) then print *, 'CFIO_GetVar warning: MAXNCNAM is larger than ',& 'dimName array size.' endif ! Get basic information from file. rc = NF90_INQUIRE (fid, nDims, nvars, ngatts, dimId) if (err("DimInqure: NF90_INQUIRE failed",rc,-48) .NE. 0)return ! Subtract dimension variables from the variable count. do i=1,nvars rc = NF90_INQUIRE_VARIABLE (fid,i,varName,varType,nvDims,vDims,nvAtts) if (err("DimInquire: variable inquire error",rc,-52) .NE. 0)& return if (nvDims .EQ. 1 .or. trim(vname) .eq. 'time_bnds') then nvars = nvars - 1 endif enddo ! Extract dimension information do i=1,nDims rc = NF90_INQUIRE_DIMENSION (fid, i, dimName, dimSize) if (err("DimInqure: can't get dim info",rc,-41) .NE. 0) return if (trim(dimName) .eq. 'nv' ) cycle if (index(dimName,'edges') .gt. 0 ) cycle if (index(dimName,'station') .gt. 0 ) cycle if (dimName=='nf') cycle if (dimName=='orientationStrLen') cycle if (dimName=='ncontact') cycle if (trim(dimName) .eq. 'XCdim') cycle if (trim(dimName) .eq. 'YCdim') cycle rc = NF90_INQ_VARID (fid, dimName, dimId) if (err("DimInqure: NF90_INQ_VARID failed",rc,-40) .NE. 0) return rc = NF90_GET_ATT(fid,dimId,'units',dimUnits) if (err("DimInqure: could not get units for dimension",rc,-53)& .NE. 0) return ! myIndex = IdentifyDim (dimName, dimUnits) ! if ( myIndex .EQ. 0 ) then ! if (dimSize .ne. im) then ! rc = -4 ! im = dimSize ! return ! endif ! else if ( myIndex .EQ. 1 ) then ! if (dimSize .ne. jm) then ! rc = -5 ! jm = dimSize ! return ! endif ! else if ( myIndex .EQ. 2 ) then ! if (kount .gt. dimSize) then ! rc = -3 ! return ! endif ! else if ( myIndex .EQ. 3 ) then ! timeId = dimId ! else ! print *, 'CFIO_GetVar: Coordinate variable ',& ! TRIM(dimName),' with units of ',TRIM(dimUnits),& ! ' is not understood.' ! rc = -19 ! return ! endif enddo ! Determine NetCDF variable ID. rc = NF90_INQ_VARID (fid, vname, vid) if (err("GetVar: variable not defined",rc,-40) .NE. 0) return ! Get beginning time & date. Calculate offset seconds from start. !ams rc = NF90_GET_ATT(fid,timeId,'begin_date',begDate) !ams if (err("GetVar: missing begin_date",rc,-44) .NE. 0) return !ams rc = NF90_GET_ATT(fid,timeId,'begin_time',begTime) !ams if (err("GetVar: missing begin_time",rc,-44) .NE. 0) return allocate(incVec(lm)) call GetDateTimeVec ( fID, begDate, begTime, incVec, rc ) if (err("GetVar: could not determine begin_date/begin_time",rc,-44) & .NE. 0) return seconds = DiffDate (begDate, begTime, yyyymmdd, hhmmss) ! Make sure input time are valid, if time is not periodic !ams print *, '+++ incSecs, begDate, begTime: ', incsecs, begDate, begTime !ams print *, '+++ seconds, yyyymmdd, hhmmss: ', seconds, yyyymmdd, hhmmss if ( .not. cyclic ) then if (seconds .LT. 0) then print *, 'CFIO_GetVar: Error code from diffdate. Problem with', & ' date/time.' rc = -7 return endif if (yyyymmdd .LT. begDate .OR. (begDate .EQ. yyyymmdd .AND. & hhmmss .LT. begTime) ) then print *, 'CFIO_GetVar: Requested time earlier than first time.' rc = -7 return endif end if if ( MOD (seconds,60) .eq. 0 ) then else print *, 'CFIO_GetVar: Currently, times must fall on minute ',& 'boundaries.' rc = -6 return endif ! Determine the time index from the time vector. timeIndex=-1 j=0 Do While (timeIndex.lt.1) j=j+1 do i=1,lm if ((seconds.le.incVec(i)).and.(timeIndex.lt.1)) timeIndex=i end do if (timeIndex.lt.1) then if (cyclic) then seconds = seconds - (incVec(2)-incVec(1))*lm else j=j+1000 end if end if if (j.gt.1000) then Write(*,'(a,L1,a,I0.6)') 'CFIO_getvar: failed to find a valid sample (C', cyclic, & '). Total iterations: ', j rc = -20 return end if End Do deallocate(incVec) !ams print *, '+++ Time Index, timeShift: ', timeIndex, timeShift ! Load starting indicies. if (useFaceDim) then if ( kbeg .eq. 0 ) then corner(4)=timeIndex edges(1)=im edges(2)=jm/6 edges(3)=6 else corner(4)=kbeg corner(5)=timeIndex edges(1)=im edges(2)=jm/6 edges(3)=6 edges(4)=kount end if else if ( kbeg .eq. 0 ) then corner(1)=1 corner(2)=1 corner(3)=timeIndex edges(1)=im edges(2)=jm edges(3)=1 else corner(1)=1 corner(2)=1 corner(3)=kbeg corner(4)=timeIndex edges(1)=im edges(2)=jm edges(3)=kount edges(4)=1 end if endif ! Determine data type. rc = NF90_INQUIRE_VARIABLE (fid, vid, varName, type, nvDims, vDims, nvAtts) if (err("GetVar: error in variable inquire",rc,-52) .NE. 0) return ! Read variable in the appropriate precision. if (HUGE(dummy) .EQ. HUGE(dummy32)) then ! -r4 if (type .EQ. NF90_FLOAT) then ! 32-bit rc = NF90_GET_VAR(fid, vid, grid, corner, edges) if(rc /=0) then print*,'Error reading variable using NF90_GET_VAR',rc print*, NF90_STRERROR(rc) return endif else if (type .EQ. NF90_DOUBLE) then ! 64-bit allocate (grid_64(im,jm,kount)) rc = NF90_GET_VAR(fid, vid, grid_64, corner, edges) do k=1,kount do j=1,jm do i=1,im grid(i,j,k) = grid_64(i,j,k) enddo enddo enddo deallocate (grid_64) else if (type .EQ. NF90_SHORT) then rc = NF90_GET_ATT(fid,vid,'scale_factor',scale_32) if (err("GetVar: error getting scale",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'add_offset',offset_32) if (err("GetVar: error getting offset",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'missing_value',amiss_16) if (err("GetVar: error getting missing",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'fmissing_value',amiss_32) if (err("GetVar: error getting fmissing",rc,-53) .NE. 0) return allocate (grid_16(im,jm,kount)) rc = NF90_GET_VAR(fid, vid, grid_16, corner, edges) do k=1,kount do j=1,jm do i=1,im if ( grid_16(i,j,k) .EQ. amiss_16 ) then grid(i,j,k) = amiss_32 else grid(i,j,k) = scale_32*grid_16(i,j,k) + offset_32 endif enddo enddo enddo deallocate (grid_16) else rc = -13 return endif else if (HUGE(dummy) .EQ. HUGE(dummy64)) then ! -r8 if (type .EQ. NF90_FLOAT) then ! 32-bit allocate (grid_32(im,jm,kount)) rc = NF90_GET_VAR(fid, vid, grid_32, corner, edges) do k=1,kount do j=1,jm do i=1,im grid(i,j,k) = grid_32(i,j,k) enddo enddo enddo deallocate (grid_32) elseif (type .EQ. NF90_DOUBLE) then ! 64-bit rc = NF90_GET_VAR(fid, vid, grid, corner, edges) else if (type .EQ. NF90_SHORT) then rc = NF90_GET_ATT(fid,vid,'scale_factor',scale_32) if (err("GetVar: error getting scale",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'add_offset',offset_32) if (err("GetVar: error getting offset",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'missing_value',amiss_16) if (err("GetVar: error getting missing",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'fmissing_value',amiss_32) if (err("GetVar: error getting fmissing",rc,-53) .NE. 0) return allocate (grid_16(im,jm,kount)) rc = NF90_GET_VAR(fid, vid, grid_16, corner, edges) do k=1,kount do j=1,jm do i=1,im if ( grid_16(i,j,k) .EQ. amiss_16 ) then grid(i,j,k) = amiss_32 else grid(i,j,k) = scale_32*grid_16(i,j,k) + offset_32 endif enddo enddo enddo deallocate (grid_16) else rc = -13 return endif else rc = -12 return endif if (err("GetVar: error reading variable",rc,-46) .NE. 0) return rc = 0 return end subroutine CFIO_GetVar !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `CFIO_PutVar` -- Write a variable to the file ! ! This routine is used to write a variable to an open CFIO ! stream. Multiple vertical levels can be written at one ! time provided they are contiguous in memory. Date and time ! must be consistent with the time increment and the starting ! date/time as defined in CFIO\_Create. Times must fall on ! minute boundaries to allow GrADS to work. Error checking is ! done for dimensions that are out of bounds. ! !#### HIstory !- 1997.10.13 da Silva/Lucchesi Initial interface design. !- 1998.02.10 Lucchesi Added support for applications running with 64-bit reals. !- 1998.03.30 Lucchesi Documentation expanded. Clean-up of code. !- 1998.07.02 Lucchesi Replaced vid with vname in argument list & made related mods to code. !- 1998.09.24 Lucchesi Changed error codes, removed DIM_CHECK if-def !- 1998.10.27 Lucchesi Added support for packing and range checks !- 1998.12.15 Lucchesi Added support for skipping times (allTimes) !- 1999.01.04 Lucchesi Fixed bug in skipping times (allTimes)/also changed variable initialization. !- 1999.07.13 Lucchesi Changes for REAL or INT time dimension ! subroutine CFIO_PutVar ( fid, vname, yyyymmdd, hhmmss, & im, jm, kbeg, kount, grid, & useFaceDim, & rc ) ! ! !USES: Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) vname !! Variable name integer yyyymmdd !! Year-month-day, e.g., 19971003 integer hhmmss !! Hour-minute-second, e.g., 120000 integer im !! size of longitudinal dimension integer jm !! size of latitudinal dimension integer kbeg !! first level to write; if 2-D grid !! use kbeg = 0. integer kount !! number of levels to write real grid(im,jm,kount) !! Gridded data to write at this time logical :: useFaceDim ! !OUTPUT PARAMETERS: integer rc !! Error return code: !! rc = 0 all is well !! rc = -2 time is inconsistent with increment !! rc = -3 number of levels is incompatible with file !! rc = -4 im is incompatible with file !! rc = -5 jm is incompatible with file !! rc = -6 time must fall on a minute boundary !! rc = -7 error in diffdate !! rc = -12 error determining default precision !! rc = -13 error determining variable type !! rc = -15 data outside of valid range !! rc = -16 data outside of packing range !! rc = -17 data outside of pack and valid range !! !! NetCDF Errors !! ------------- !! rc = -38 error from NF90_PUT_VAR (dimension variable) !! rc = -40 error from NF90_INQ_VARID !! rc = -41 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (lat or lon) !! rc = -42 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (lev) !! rc = -43 error from NF90_INQ_VARID (time variable) !! rc = -44 error from NF90_GET_ATT (time attribute) !! rc = -45 error from NF90_PUT_VAR !! rc = -46 error from NF90_GET_VAR !! rc = -52 error from NF90_INQUIRE_VARIABLE !! rc = -53 error from NF90_GET_ATT ! !------------------------------------------------------------------------- integer timeid, timeDimId, dimSize, timeType character(len=MAXCHR) dimName integer corner(5), edges(5) integer vid integer(INT64) seconds integer timeIndex integer minutes ! added as a work-around integer i, j, k integer begDate, begTime, timInc integer hour,min,sec,incSecs integer, allocatable :: allTimes(:) integer fillTime ! Variables for dealing with precision real(kind=REAL32), allocatable :: grid_32(:,:,:) real(kind=REAL64), allocatable :: grid_64(:,:,:) real(kind=REAL32) dummy32 real(kind=REAL64) dummy64 real dummy ! Variables for NF90_INQUIRE_VARIABLE character(len=MAXCHR) varName integer type, nvDims, vdims(MAXVDIMS), nvAtts ! Variables for packing and range checking integer(kind=INT16), allocatable :: grid_16(:,:,:) real(kind=REAL32), allocatable :: fminutes_32(:) real(kind=REAL32) high_32, low_32, amiss_32 real(kind=REAL32) scale_32, offset_32 logical outRange logical outPRange ! Variable initialization outRange = .FALSE. outPRange = .FALSE. ! Check to make sure max string lengths are large enough. NetCDF defines ! MAXNCNAM, but it can't be used in a character(len=MAXNCNAM) statement. if (MAXCHR .LT. MAXNCNAM) then print *, 'CFIO_PutVar warning: MAXNCNAM is larger than ',& 'dimName array size.' endif ! Determine NetCDF variable ID. rc = NF90_INQ_VARID (fid, vname, vid) if (err("PutVar: variable not defined",rc,-40) .NE. 0) return ! Determine number of seconds since starting date/time. rc = NF90_INQ_VARID (fid, 'time', timeId) if (err("PutVar: time not defined",rc,-43) .NE. 0) return rc = NF90_INQ_DIMID(fid, 'time', timeDimId) rc = NF90_INQUIRE_DIMENSION(fid, timeDimId, dimName, dimSize) rc = NF90_GET_ATT(fid,timeId,'begin_date',begDate) if (err("PutVar: missing begin_date",rc,-44) .NE. 0) return rc = NF90_GET_ATT(fid,timeId,'begin_time',begTime) if (err("PutVar: missing begin_time",rc,-44) .NE. 0) return seconds = DiffDate (begDate, begTime, yyyymmdd, hhmmss) if (seconds .lt. 0) then print *, 'CFIO_PutVar: Error code from diffdate. Problem with',& ' date/time.' rc = -7 return endif if ( MOD (int(seconds),60) .eq. 0 ) then minutes = seconds / 60 else print *, 'CFIO_PutVar: Currently, times must fall on minute ',& 'boundaries.' rc = -6 return endif ! Confirm that this time is consistent with the starting time coupled with ! the time increment. rc = NF90_GET_ATT(fid,timeId,'time_increment',timInc) if (err("PutVar: missing time increment",rc,-44) .NE. 0) return ! Convert time increment to seconds. !ams write (strBuf,203) timinc !ams 203 format (I6) !ams read (strBuf,204) hour, min, sec !ams 204 format (3I2) call CFIO_parseIntTime ( timinc, hour, min, sec ) incSecs = hour*3600 + min*60 + sec if ( MOD (int(seconds), incSecs) .ne. 0 ) then print *, 'CFIO_putvar: Absolute time of ',seconds,' not ',& 'possible with an interval of ',incSecs rc = -2 return else timeIndex = seconds/incSecs + 1 endif ! Load starting indicies. corner = 1 edges = 1 if ( kbeg .eq. 0 ) then edges(1)=im if (useFaceDim) then edges(2)=jm/6 edges(3)=6 corner(4)=timeIndex else edges(2)=jm corner(3)=timeIndex end if else edges(1)=im if (useFaceDim) then corner(4)=kbeg corner(5)=timeIndex edges(2)=jm/6 edges(3)=6 edges(4)=kount else corner(3)=kbeg corner(4)=timeIndex edges(2)=jm edges(3)=kount end if endif ! Check variable against valid range. rc = NF90_GET_ATT(fid,vid,'vmin',low_32) if (err("PutVar: can't get vmin",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'vmax',high_32) if (err("PutVar: can't get vmax",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'fmissing_value',amiss_32) if (err("PutVar: can't get fmissing_value",rc,-53) .NE. 0) return if (abs(low_32) .NE. amiss_32 .OR. high_32 .NE. amiss_32) then do k=1,kount do j=1,jm do i=1,im if (grid(i,j,k) .GT. high_32 .OR. grid(i,j,k) .LT. & low_32) then outRange = .TRUE. goto 100 endif enddo enddo enddo 100 continue endif ! Determine if we are writing single- or double-precision. rc = NF90_INQUIRE_VARIABLE (fid, vid, varName, type, nvDims, vDims, nvAtts) if (err("PutVar: error in variable inquire",rc,-52) .NE. 0) return ! Write variable in the appropriate precision. if (HUGE(dummy) .EQ. HUGE(dummy32)) then ! -r4 if (type .EQ. NF90_FLOAT) then ! 32-bit rc = NF90_PUT_VAR(fid, vid, grid, corner, edges) else if (type .EQ. NF90_DOUBLE) then ! 64-bit allocate (grid_64(im,jm,kount)) do k=1,kount do j=1,jm do i=1,im grid_64(i,j,k) = grid(i,j,k) enddo enddo enddo rc = NF90_PUT_VAR(fid, vid, grid_64, corner, edges) deallocate (grid_64) else if (type .EQ. NF90_SHORT) then rc = NF90_GET_ATT(fid,vid,'packmax',high_32) if (err("PutVar: error getting packmax",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'packmin',low_32) if (err("PutVar: error getting packmin",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'scale_factor',scale_32) if (err("PutVar: error getting scale",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'add_offset',offset_32) if (err("PutVar: error getting offset",rc,-53) .NE. 0) return allocate (grid_16(im,jm,kount)) do k=1,kount do j=1,jm do i=1,im if ( grid(i,j,k) .LT. low_32 .OR. grid(i,j,k) .GT. & high_32) then grid_16(i,j,k) = PACK_FILL outPRange = .TRUE. else grid_16(i,j,k) = (grid(i,j,k) - offset_32)/scale_32 endif enddo enddo enddo rc = NF90_PUT_VAR(fid, vid, grid_16, corner, edges) deallocate (grid_16) else rc = -13 return endif else if (HUGE(dummy) .EQ. HUGE(dummy64)) then ! -r8 if (type .EQ. NF90_FLOAT) then ! 32-bit allocate (grid_32(im,jm,kount)) do k=1,kount do j=1,jm do i=1,im grid_32(i,j,k) = grid(i,j,k) enddo enddo enddo rc= NF90_PUT_VAR(fid, vid, grid_32, corner, edges) deallocate (grid_32) else if (type .EQ. NF90_DOUBLE) then ! 64-bit rc = NF90_PUT_VAR(fid, vid, grid, corner, edges) else if (type .EQ. NF90_SHORT) then rc = NF90_GET_ATT(fid,vid,'packmax',high_32) if (err("PutVar: error getting packmax",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'packmin',low_32) if (err("PutVar: error getting packmin",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'scale_factor',scale_32) if (err("PutVar: error getting scale",rc,-53) .NE. 0) return rc = NF90_GET_ATT(fid,vid,'add_offset',offset_32) if (err("PutVar: error getting offset",rc,-53) .NE. 0) return allocate (grid_16(im,jm,kount)) do k=1,kount do j=1,jm do i=1,im if ( grid(i,j,k) .LT. low_32 .OR. grid(i,j,k) .GT. & high_32) then grid_16(i,j,k) = PACK_FILL outPRange = .TRUE. else grid_16(i,j,k) = (grid(i,j,k) - offset_32)/scale_32 endif enddo enddo enddo rc = NF90_PUT_VAR(fid, vid, grid_16, corner, edges) deallocate (grid_16) else rc = -13 return endif else rc = -12 return endif if (err("PutVar: error writing variable",rc,-45) .NE. 0) return ! Read time dimension scale and fill all values up to the current time. ! This will insure missing times are defined with the proper time value. rc = NF90_INQUIRE_DIMENSION (fid, timeDimId, dimName, dimSize) dimSize = dimSize - 1 ! We've already written the ! the new time. allocate ( allTimes (MAX(timeIndex,dimSize)) ) allocate ( fminutes_32 (MAX(timeIndex,dimSize)) ) rc = NF90_INQUIRE_VARIABLE (fid,timeId,dimName,timeType,nvDims,vDims,nvAtts) if (dimSize .GT. 0) then ! Depending on the version of CFIO used to write the file, the Time ! dimension variable can either be floating point or integer. corner(1)=1 edges(1)=dimSize if (timeType .EQ. NF90_FLOAT) then rc = NF90_GET_VAR(fid, timeId, fminutes_32, corner, edges) do i=1,dimSize allTimes(i) = INT(fminutes_32(i)) enddo else if (timeType .EQ. NF90_INT) then rc = NF90_GET_VAR(fid, timeId, allTimes, corner, edges) endif if (err("PutVar: error reading times from file",rc,-46) .NE. 0)& return endif ! This loop fills the time dimension scale based on the time increment ! specified in CFIO_Create. If CFIO ever changes to support variable ! time increments, this code MUST be changed. do i=1,timeIndex-1 fillTime = (i-1) * incSecs/60 allTimes(i) = fillTime enddo allTimes(timeIndex) = minutes ! Write filled time array to file. corner(1)=1 edges(1)=timeIndex if (timeType .EQ. NF90_FLOAT) then do i=1,timeIndex fminutes_32(i) = INT(allTimes(i)) enddo rc = NF90_PUT_VAR(fid,timeId,fminutes_32,corner,edges) else if (timeType .EQ. NF90_INT) then rc = NF90_PUT_VAR(fid,timeId,allTimes,corner,edges) endif if (err("PutVar: error writing time",rc,-38) .NE. 0) return if (outRange .AND. outPRange) then rc = -17 else if (outPRange) then rc = -16 else if (outRange) then rc = -15 else rc = 0 endif deallocate ( allTimes ) deallocate ( fminutes_32 ) return end subroutine CFIO_PutVar !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !BOI ! ! !TITLE: Returns a new date/time from an initial date/time and offset ! ! !AUTHORS: Rob Lucchesi ! ! !AFFILIATION: Data Assimilation Office, NASA/GSFC, Greenbelt, MD 20771 ! ! !DATE: July 20, 1998 ! !EOI !------------------------------------------------------------------------- !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `GetDateInt8` --- Returns a new date/time from an initial date/time and offset ! ! This subroutine returns a new date and time in yyyymmdd ! and hhmmss format given and initial date, time, and ! offset in seconds. The routine converts the input date ! and time to julian seconds, adds the offset, and converts ! back to yyyymmdd and hhmmss format. This routine has been ! tested for Y2K compiance. ! !#### History !- 1998.07.20 Lucchesi Initial version. !- 2010.05.11 Lucchesi Integer for julian seconds changed to 64-bit. StartDate ! constant no longer needed. ! subroutine GetDateInt8 (yyyymmdd_1,hhmmss_1,offset, & yyyymmdd_2,hhmmss_2,rc) implicit none ! ! !INPUT PARAMETERS: ! integer yyyymmdd_1 !! Initial date in YYYYYMMDD format integer hhmmss_1 !! Initial time in HHMMSS format integer(kind=INT64) offset !! Offset to add (in seconds) ! ! !OUTPUT PARAMETERS: ! integer yyyymmdd_2 !! New date in YYYYMMDD format integer hhmmss_2 !! New time in HHMMSS format integer rc !! Return code. (<0 = error) ! !------------------------------------------------------------------------- integer year1,mon1,day1,hour1,min1,sec1 integer year2,mon2,day2,hour2,min2,sec2 integer(kind=INT64) julian1 integer(kind=INT64) julsec, remainder !character(len=8) dateString ! Error checking. if (yyyymmdd_1 .lt. 19000000 .or. yyyymmdd_1 .gt. 21000000 ) then rc=-1 return endif if (hhmmss_1 .lt. 0 .or. hhmmss_1 .ge. 240000 ) then rc=-1 return endif ! Convert Date/Time strings to integer variables. !ams write (dateString, 200) yyyymmdd_1 !ams 200 format (I8) !ams read (dateString, 201) year1, mon1, day1 !ams 201 format (I4,2I2) !ams write (dateString, 202) hhmmss_1 !ams 202 format (I6) !ams read (dateString, 203) hour1, min1, sec1 !ams 203 format (3I2) call CFIO_parseIntTime ( yyyymmdd_1, year1, mon1, day1 ) call CFIO_parseIntTime ( hhmmss_1, hour1, min1, sec1 ) ! Get Julian Day and subtract off a constant (Julian days since 7/14/66) julian1 = julday (mon1, day1, year1) ! Calculcate Julian seconds julsec = (julian1-1)*86400 + hour1*3600 + min1*60 + sec1 ! Add offset and calculate new julian day. julsec = julsec + offset julian1 = INT(julsec/86400) + 1 remainder = MOD(julsec,86400_INT64) ! Convert julian day to YYYYMMDD. call caldat (julian1, mon2, day2, year2) ! Calculate HHMMSS from the remainder. hour2 = INT(remainder/3600) remainder = MOD(remainder,3600_INT64) min2 = INT(remainder/60) sec2 = MOD(remainder,60_INT64) ! Build YYYYMMDD and HHMMSS variables. yyyymmdd_2 = year2*10000 + mon2*100 + day2 hhmmss_2 = hour2*10000 + min2*100 + sec2 rc = 0 return end subroutine GetDateInt8 !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !> ! `GetDateInt4` --- Returns a new date/time from an initial date/time and offset ! ! This subroutine returns a new date and time in yyyymmdd ! and hhmmss format given and initial date, time, and ! offset in seconds. The routine converts the input date ! and time to julian seconds, adds the offset, and converts ! back to yyyymmdd and hhmmss format. This routine has been ! tested for Y2K compiance. ! !#### History !- 1998.07.20 Lucchesi Initial version. !- 2010.05.11 Lucchesi Integer for julian seconds changed to 64-bit. StartDate ! constant no longer needed. ! subroutine GetDateInt4 (yyyymmdd_1,hhmmss_1,offset, & yyyymmdd_2,hhmmss_2,rc) implicit none ! ! !INPUT PARAMETERS: ! integer yyyymmdd_1 !! Initial date in YYYYYMMDD format integer hhmmss_1 !! Initial time in HHMMSS format integer offset !! Offset to add (in seconds) ! ! !OUTPUT PARAMETERS: ! integer yyyymmdd_2 !! New date in YYYYMMDD format integer hhmmss_2 !! New time in HHMMSS format integer rc !! Return code. (<0 = error) ! !------------------------------------------------------------------------- integer(Kind=INT64) offsetLong offsetLong = INT(offset,Kind=INT64) call GetDateInt8(yyyymmdd_1,hhmmss_1,offsetLong,yyyymmdd_2,hhmmss_2,rc) return end subroutine GetDateInt4 !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !------------------------------------------------------------------------- !BOP ! ! !ROUTINE: CFIO_GetMissing --- Return missing value ! ! !DESCRIPTION: Returns missing value on file ! ! !INTERFACE: real function CFIO_GetMissing ( fid, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: integer fid ! file id integer rc ! Error code ! ! !REVISION HISTORY: ! ! 1999.11.01 da Silva Initial code. ! !EOP !------------------------------------------------------------------------- integer nDims, recdim, ngatts integer varType, nvDims, vDims(MAXVDIMS), nvAtts character(len=MAXCHR) vnameTemp integer i integer attType, attLen integer allVars ! all variables - includes dimension vars real(kind=REAL32) amiss_32 ! Get basic information from the file rc = NF90_INQUIRE (fid,nDims,allVars,ngatts,recdim) if (err("Inqure: NF90_INQUIRE failed",rc,-48) .NE. 0) return do i= 1, allVars rc = NF90_INQUIRE_VARIABLE (fid,i,vnameTemp,varType,nvDims,vDims,nvAtts) if (err("CFIO_GetMissing: variable inquire error",rc,-52) .NE. 0) return if (nvDims .EQ. 1) then ! coord variable cycle else ! noon-coord variable rc = NF90_GET_ATT(fid,i,'fmissing_value',amiss_32) if (rc .NE. 0) then rc = NF90_INQUIRE_ATTRIBUTE (fid, i, 'missing_value', attType, attLen) if (rc.eq.0 .and. attType .EQ. NF90_FLOAT) then rc = NF90_GET_ATT(fid,allVars,'missing_value',amiss_32) if (err("CFIO_GetMissing: error getting missing value",rc,-53) & .NE. 0) return else print *, & 'CFIO_GetMissing: Cannot find missing value, assuming 1E+15' amiss_32 = 1.0E+15 end if endif exit ! just check first non-ccordinate variable endif end do CFIO_GetMissing = amiss_32 rc = 0 end function CFIO_GetMissing ! ! The next routine adapted from mpeu... ! !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !BOP ------------------------------------------------------------------- ! ! !MODULE: m_StrTemplate - A template formatting a string with variables ! ! !DESCRIPTION: ! ! A template resolver formatting a string with a string variable ! and time variables. The format descriptors are similar to those ! used in the GrADS. ! ! "%y4" substitute with a 4 digit year ! "%y2" a 2 digit year ! "%m1" a 1 or 2 digit month ! "%m2" a 2 digit month ! "%mc" a 3 letter month in lower cases ! "%Mc" a 3 letter month with a leading letter in upper case ! "%MC" a 3 letter month in upper cases ! "%d1" a 1 or 2 digit day ! "%d2" a 2 digit day ! "%h1" a 1 or 2 digit hour ! "%h2" a 2 digit hour ! "%h3" a 3 digit hour (?) ! "%n2" a 2 digit minute ! "%s" a string variable ! "%%" a "%" ! ! !INTERFACE: ! !REVISION HISTORY: ! 19Dec06 - Jing Guo <guo@gmao.gsfc.nasa.gov> ! - Merged changes between 1.1.2.6 and 1.1.2.9 to 1.2, ! including a fix at bug nymd==0 and environment ! variable ($env or ${env}) support if getenv() is ! available from the system. ! 01Jun99 - Jing Guo <guo@dao.gsfc.nasa.gov> ! - initial prototype/prolog/code !EOP ___________________________________________________________________ !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !BOP ------------------------------------------------------------------- !> ! `strTemplate_` - expanding a format template to a string ! !#### History !- 03Jun99 - Jing Guo <guo@dao.gsfc.nasa.gov> - initial prototype/prolog/code !- 08Jan01 - da Silva: moved uppercase() to outside select() to avoid coredump on Linux/PGI. ! subroutine strTemplate_(str,tmpl,class,xid,nymd,nhms,stat) implicit none character(len=*),intent(out) :: str ! the output character(len=*),intent(in ) :: tmpl ! a "format" character(len=*),intent(in ),optional :: class ! choose a UNIX or a GrADS(defulat) type format character(len=*),intent(in ),optional :: xid ! a string substituting a "%s". Trailing ! spaces will be ignored integer,intent(in ),optional :: nymd ! yyyymmdd, substituting "%y4", "%y2", "%m1", ! "%m2", "%mc", "%Mc', and "%MC" integer,intent(in ),optional :: nhms ! hhmmss, substituting "%h1", "%h2", "%h3", ! and "%n2" integer,intent(out),optional :: stat ! error code character(len=*),parameter :: myname_=myname//'::strTemplate_' character(len=16) :: tmpl_class,uc_class tmpl_class="GX" if(present(class)) tmpl_class=class !!! uc_class=uppercase(tmpl_class) uc_class=trim(tmpl_class) ! removed this dependency on mpeu select case(uc_class) case("GX","GRADS") call GX_(str,tmpl,xid,nymd,nhms,stat) !case("UX","UNIX") ! yet to be implemented ! call UX_(str,tmpl,xid,nymd,nhms,stat) case default write(stderr,'(4a)') myname_,': unknown class, "',trim(tmpl_class),'"' if(.not.present(stat)) call die(myname_) stat=-1 return end select end subroutine strTemplate_ !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! !BOP ------------------------------------------------------------------- !> ! `GX_` - evaluate a GrADS style string template ! !#### History !- 01Jun99 - Jing Guo <guo@dao.gsfc.nasa.gov> - initial prototype/prolog/code ! subroutine GX_(str,tmpl,xid,nymd,nhms,stat) implicit none character(len=*),intent(out) :: str character(len=*),intent(in ) :: tmpl character(len=*),optional,intent(in) :: xid integer,optional,intent(in) :: nymd integer,optional,intent(in) :: nhms integer,optional,intent(out) :: stat character(len=*),parameter :: myname_=myname//'::GX_' integer :: iy4,iy2,imo,idy integer :: ihr,imn integer :: i,i1,i2,m,k integer :: ln_tmpl,ln_str integer :: istp,kstp integer :: ier character(len=1) :: c0,c1,c2 character(len=4) :: sbuf !________________________________________ ! Determine iyr, imo, and idy iy4=-1 iy2=-1 imo=-1 idy=-1 if(present(nymd)) then if(nymd <= 0) then call perr(myname_,'nymd <= 0',nymd) if(.not.present(stat)) call die(myname_) stat=1 return endif i=nymd iy4=i/10000 iy2=mod(iy4,100) i=mod(i,10000) imo=i/100 i=mod(i,100) idy=i endif !________________________________________ ! Determine ihr and imn ihr=-1 imn=-1 if(present(nhms)) then if(nhms < 0) then call perr(myname_,'nhms < 0',nhms) if(.not.present(stat)) call die(myname_) stat=1 return endif i=nhms ihr=i/10000 i=mod(i,10000) imn=i/100 endif !________________________________________ ln_tmpl=len_trim(tmpl) ! size of the format template ln_str =len(str) ! size of the output string !________________________________________ if(present(stat)) stat=0 str="" i=0; istp=1 k=1; kstp=1 do while( i+istp <= ln_tmpl ) ! A loop over all tokens in (tmpl) if(k>ln_Str) exit ! truncate the output here. i=i+istp c0=tmpl(i:i) select case(c0) case ("$") call genv_(tmpl,ln_tmpl,i,istp,str,ln_str,k,ier) if(ier/=0) then call perr(myname_,'genv_("'//tmpl(i:ln_tmpl)//'"',ier) if(.not.present(stat)) call die(myname_) stat=1 return endif case ("%") !________________________________________ c1="" i1=i+1 if(i1 <= ln_Tmpl) c1=tmpl(i1:i1) !________________________________________ select case(c1) case("s") if(.not.present(xid)) then write(stderr,'(2a)') myname_, & ': optional argument expected, "xid="' if(.not.present(stat)) call die(myname_) stat=1 return endif istp=2 m=min(k+len_trim(xid)-1,ln_str) str(k:m)=xid k=m+1 cycle case("%","$") istp=2 str(k:k)=c1 k=k+1 ! kstp=1 cycle case default c2="" i2=i+2 if(i2 <= ln_Tmpl) c2=tmpl(i2:i2) !________________________________________ select case(c1//c2) case("y4","y2","m1","m2","mc","Mc","MC","d1","d2") if(.not.present(nymd)) then write(stderr,'(2a)') myname_, & ': optional argument expected, "nymd="' if(.not.present(stat)) call die(myname_) stat=1 return endif istp=3 case("h1","h2","h3","n2") if(.not.present(nhms)) then write(stderr,'(2a)') myname_, & ': optional argument expected, "nhms="' if(.not.present(stat)) call die(myname_) stat=1 return endif istp=3 case default write(stderr,'(4a)') myname_, & ': invalid template entry, "',trim(tmpl(i:)),'"' if(.not.present(stat)) call die(myname_) stat=2 return end select ! case(c1//c2) end select ! case(c1) !________________________________________ select case(c1) case("y") select case(c2) case("2") write(sbuf,'(i2.2)') iy2 kstp=2 case("4") write(sbuf,'(i4.4)') iy4 kstp=4 case default write(stderr,'(4a)') myname_, & ': invalid template entry, "',trim(tmpl(i:)),'"' if(.not.present(stat)) call die(myname_) stat=2 return end select case("m") select case(c2) case("1") if(imo < 10) then write(sbuf,'(i1)') imo kstp=1 else write(sbuf,'(i2)') imo kstp=2 endif case("2") write(sbuf,'(i2.2)') imo kstp=2 case("c") sbuf=mon_lc(imo) kstp=3 case default write(stderr,'(4a)') myname_, & ': invalid template entry, "',trim(tmpl(i:)),'"' if(.not.present(stat)) call die(myname_) stat=2 return end select case("M") select case(c2) case("c") sbuf=mon_wd(imo) kstp=3 case("C") sbuf=mon_uc(imo) kstp=3 case default write(stderr,'(4a)') myname_, & ': invalid template entry, "',trim(tmpl(i:)),'"' if(.not.present(stat)) call die(myname_) stat=2 return end select case("d") select case(c2) case("1") if(idy < 10) then write(sbuf,'(i1)') idy kstp=1 else write(sbuf,'(i2)') idy kstp=2 endif case("2") write(sbuf,'(i2.2)') idy kstp=2 case default write(stderr,'(4a)') myname_, & ': invalid template entry, "',trim(tmpl(i:)),'"' if(.not.present(stat)) call die(myname_) stat=2 return end select case("h") select case(c2) case("1") if(ihr < 10) then write(sbuf,'(i1)') ihr kstp=1 else write(sbuf,'(i2)') ihr kstp=2 endif case("2") write(sbuf,'(i2.2)') ihr kstp=2 case("3") write(sbuf,'(i3.3)') ihr kstp=3 case default write(stderr,'(4a)') myname_, & ': invalid template entry, "',trim(tmpl(i:)),'"' if(.not.present(stat)) call die(myname_) stat=2 return end select case("n") select case(c2) case("2") write(sbuf,'(i2.2)') imn kstp=2 case default write(stderr,'(4a)') myname_, & ': invalid template entry, "',trim(tmpl(i:)),'"' if(.not.present(stat)) call die(myname_) stat=2 return end select case default write(stderr,'(4a)') myname_, & ': invalid template entry, "',trim(tmpl(i:)),'"' if(.not.present(stat)) call die(myname_) stat=2 return end select ! case(c1) m=min(k+kstp-1,ln_Str) str(k:m)=sbuf k=m+1 case default istp=1 str(k:k)=tmpl(i:i) k=k+1 end select ! case(c0) end do contains subroutine genv_(tmpl,lnt,i,istp,str,lns,k,ier) implicit none character(len=*),intent(in) :: tmpl integer,intent(in) :: lnt integer,intent(in) :: i integer,intent(out) :: istp character(len=*),intent(inout) :: str integer ,intent(in) :: lns integer ,intent(inout) :: k integer,intent(out) :: ier integer :: j,jb,je integer :: l,m logical :: bracket,more character(len=256) :: env j=i+1 ! skip "$" ier=0 if(j>lnt) then ier=1 return endif bracket = tmpl(j:j)=='{' if(bracket) j=j+1 ! There is at least one a letter (including "_") to start a ! variable name select case(tmpl(j:j)) case ("A":"Z","a":"z","_") case default ier=2 return end select jb=j je=j if(bracket) then more=.true. do while(more) select case(tmpl(j:j)) case ("A":"Z","a":"z","_","0":"9") je=j j=j+1 case ("}") ! End if "}" or eos j=j+1 exit case default ier=3 return end select more=j<=lnt enddo else more=.true. do while(more) select case(tmpl(j:j)) case ("A":"Z","a":"z","_","0":"9") je=j j=j+1 case default exit end select more=j<=lnt enddo endif istp=j-i call get_environment_variable(tmpl(jb:je),env) l=len_trim(env) m=min(k+l-1,lns) str(k:m)=env k=m+1 end subroutine genv_ end subroutine GX_ ! ! MPEU Stubs ! subroutine perr ( name, msg, i ) character(len=*) :: name, msg integer :: i print *, trim(name)//':'//trim(msg), i return end subroutine perr subroutine die ( name ) character(len=*) :: name print *, trim(name)//': fatal error in CFIO, but not aborting...' return end subroutine die #if defined(HDFEOS) || defined(HDFSD) !------------------------------------------------------------------------------- !> ! The routine `EOS_Close` ... ! !#### History !- 1997.10.13 da Silva/Lucchesi Initial interface design. !- 1998.03.30 Lucchesi Documentation expanded. Clean-up of code. Added rc. ! subroutine EOS_Close ( fid, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle ! ! !OUTPUT PARAMETERS: ! integer rc !! Error return code: !! rc = 0 all is well !! !! NetCDF Errors !! ------------- !! rc = -54 error from ncclos (file close) ! !------------------------------------------------------------------------- integer i integer sfend integer GDclose ! Make NetCDF errors non-fatal, but issue warning messages. #if defined(HDFSD) rc = sfend (fid) #endif #if defined(HDFEOS) rc = GDclose (fid) #endif rc = 0 return end subroutine EOS_Close !------------------------------------------------------------------------- !> ! The routine `EOS_PutIntAtt` ... ! !#### History !- 1998.07.30 Lucchesi Initial interface design. !- 1998.07.30 Lucchesi Initial coding. !- 1998.09.24 Lucchesi Changed error handling. !- 1998.09.28 Lucchesi Added support for multiple precisions !- 1999.01.29 Lucchesi Converted API to SD for HDFEOS ! subroutine EOS_PutIntAtt ( fid, name, count, buf, prec, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) name !! Name of attribute integer count !! Number of integers to write integer buf(count) !! Buffer with integer values integer prec !! Desired precision of attribute value: !! 0 = 32 bit !! 1 = 64 bit ! ! !OUTPUT PARAMETERS: ! integer rc !! Error return code: !! rc = 0 all is well !! rc = -12 error determining default precision !! !! NetCDF Errors !! ------------- !! rc = -36 error from NF90_PUT_ATT (global attribute) !! rc = -55 error from NF90_REDEF (enter define mode) !! rc = -56 error from NF90_ENDDEF (exit define mode) ! !------------------------------------------------------------------------- integer(kind=INT32) dummy32 integer(kind=INT64) dummy64 integer i integer(kind=INT32), allocatable :: buf32(:) integer(kind=INT64), allocatable :: buf64(:) integer rct #if defined(HDFEOS) ! Variables for HDF-EOS integer hdfFid, hdfeos_fid integer gridId ! Functions for HDF-EOS integer GDattach, EHidinfo, GDdetach #endif ! Functions for SD interface integer sfsnatt rc = 0 #if defined (HDFEOS) ! Attach EOS-HDF grid and obtain SD itentifier. Place this ! in 'fid' after saving the HDF-EOS fid in 'hdfeos_fid.' hdfeos_fid = fid gridId = GDattach (hdfeos_fid, GRID_NAME) rct = EHidinfo (hdfeos_fid, hdfFid, fid) #endif if ( HUGE(dummy32) .EQ. HUGE(i) .AND. prec .EQ. 0 ) then ! -i4 rct = sfsnatt( fid, name, DFNT_INT32, count, buf) else if ( HUGE(dummy32) .EQ. HUGE(i) .AND. prec .EQ. 1 ) then ! -i4 allocate ( buf64(count) ) ! 64-bit out do i=1,count buf64(i) = buf(i) enddo rct = sfsnatt( fid, name, DFNT_FLOAT64, count, buf) deallocate (buf64) else if (HUGE(dummy64) .EQ. HUGE(i) .AND. prec .EQ. 0 ) then ! -i8 allocate ( buf32(count) ) ! 32-bit out do i=1,count buf32(i) = buf(i) enddo rct = sfsnatt( fid, name, DFNT_INT32, count, buf) deallocate (buf32) else if (HUGE(dummy64) .EQ. HUGE(i) .AND. prec .EQ. 1 ) then ! -i8 rct = sfsnatt( fid, name, DFNT_FLOAT64, count, buf) else rc = -12 goto 999 endif if (err("PutIntAtt: error writing attribute",rc,-36) .NE. 0) return 999 continue ! Return processing. I hate GOTOs, but we need to detach ! the grid if in HDF-EOS mode before returning and this seemed ! like the best way without having messy if-defs all through ! the code. #if defined (HDFEOS) fid = hdfeos_fid rct = GDdetach (gridId) #endif return end subroutine EOS_PutIntAtt !---------------------------------------------------------------------------- !> ! The routine `EOS_PutRealAtt` ... ! !#### History !- 1998.07.30 Lucchesi Initial interface design. !- 1998.07.30 Lucchesi Initial coding. !- 1998.09.24 Lucchesi Changed error handling. !- 1998.09.28 Lucchesi Added support for multiple precisions !- 1999.01.29 Lucchesi Converted API to SD for HDFEOS ! subroutine EOS_PutRealAtt ( fid, name, count, buf, prec, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) name !! Name of attribute integer count !! Number of integers to write real buf(count) !! Buffer with real values integer prec !! Desired precision of attribute value: !! 0 = 32 bit !! 1 = 64 bit ! ! !OUTPUT PARAMETERS: ! integer rc !! Error return code: !! rc = 0 all is well !! rc = -12 error determining default precision !! !! NetCDF Errors !! ------------- !! rc = -36 error from NF90_PUT_ATT (global attribute) !! rc = -55 error from NF90_REDEF (enter define mode) !! rc = -56 error from NF90_ENDDEF (exit define mode) ! !------------------------------------------------------------------------- real(kind=REAL32) dummy32 real(kind=REAL64) dummy64 real r integer i real(kind=REAL32), allocatable :: buf32(:) real(kind=REAL64), allocatable :: buf64(:) integer rct #if defined(HDFEOS) ! Variables for HDF-EOS integer hdfFid, hdfeos_fid integer gridId ! Functions for HDF-EOS integer GDattach, EHidinfo, GDdetach #endif ! Functions for SD interface integer sfsnatt rc = 0 #if defined (HDFEOS) ! Attach EOS-HDF grid and obtain SD itentifier. Place this ! in 'fid' after saving the HDF-EOS fid in 'hdfeos_fid.' hdfeos_fid = fid gridId = GDattach (hdfeos_fid,GRID_NAME) rct = EHidinfo (hdfeos_fid, hdfFid, fid) #endif if (HUGE(dummy32) .EQ. HUGE(r) .AND. prec .EQ. 0) then ! -r4 rct = sfsnatt( fid, name, DFNT_FLOAT32, count, buf) else if (HUGE(dummy32) .EQ. HUGE(r) .AND. prec .EQ. 1) then ! -r4 allocate (buf64(count)) ! 64-bit out do i=1,count buf64(i) = buf(i) enddo rct = sfsnatt( fid, name, DFNT_FLOAT64, count, buf64) deallocate (buf64) else if (HUGE(dummy64) .EQ. huge(r) .AND. prec .EQ. 0) then ! -r8 allocate (buf32(count)) ! 32-bit out do i=1,count buf32(i) = buf(i) enddo rct = sfsnatt( fid, name, DFNT_FLOAT32, count, buf32) deallocate (buf32) else if (HUGE(dummy64) .EQ. huge(r) .AND. prec .EQ. 1) then ! -r8 rct = sfsnatt( fid, name, DFNT_FLOAT64, count, buf) else rc = -12 goto 999 endif if (err("PutRealAtt: error writing attribute",rc,-36) .NE. 0) & goto 999 999 continue ! Return processing. I hate GOTOs, but we need to detach ! the grid if in HDF-EOS mode before returning and this seemed ! like the best way without having messy if-defs all through ! the code. #if defined (HDFEOS) fid = hdfeos_fid rct = GDdetach (gridId) #endif return end subroutine EOS_PutRealAtt !------------------------------------------------------------------------------------ !> ! The routine `EOS_PutCharAtt` .... ! !#### History !- 1998.07.30 Lucchesi Initial interface design. !- 1998.07.30 Lucchesi Initial coding. !- 1998.09.24 Lucchesi Changed error handling. !- 1999.01.29 Lucchesi Converted API to SD for HDFEOS ! subroutine EOS_PutCharAtt ( fid, name, count, buf, rc ) ! ! !USES: ! Implicit NONE ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) name !! Name of attribute integer count !! Number of characters to write character(len=MLEN) buf !! Buffer containing string ! ! !OUTPUT PARAMETERS: ! integer rc !! Error return code: !! rc = 0 all is well !! !! NetCDF Errors !! ------------- !! rc = -36 error from NF90_PUT_ATT (global attribute) !! rc = -55 error from NF90_REDEF (enter define mode) !! rc = -56 error from NF90_ENDDEF (exit define mode) ! !------------------------------------------------------------------------- integer rct #if defined(HDFEOS) ! Variables for HDF-EOS integer hdfFid, hdfeos_fid integer gridId ! Functions for HDF-EOS integer GDattach, EHidinfo, GDdetach #endif ! Functions for SD interface integer sfscatt rc = 0 #if defined (HDFEOS) ! Attach EOS-HDF grid and obtain SD itentifier. Place this ! in 'fid' after saving the HDF-EOS fid in 'hdfeos_fid.' hdfeos_fid = fid gridId = GDattach (hdfeos_fid, GRID_NAME) rct = EHidinfo (hdfeos_fid, hdfFid, fid) #endif rct = sfscatt (fid, name, DFNT_CHAR8, count, buf) 999 continue ! Return processing. I hate GOTOs, but we need to detach ! the grid if in HDF-EOS mode before returning and this seemed ! like the best way without having messy if-defs all through ! the code. #if defined (HDFEOS) fid = hdfeos_fid rct = GDdetach (gridId) #endif return end subroutine EOS_PutCharAtt !------------------------------------------------------------------------- !> ! The routine `EOS_PutVar` ... ! !#### History ! 1997.10.13 da Silva/Lucchesi Initial interface design. ! 1998.02.10 Lucchesi Added support for applications running with ! 64-bit reals. ! 1998.03.30 Lucchesi Documentation expanded. Clean-up of code. ! 1998.07.02 Lucchesi Replaced vid with vname in argument list & ! made related mods to code. ! 1998.09.24 Lucchesi Changed error codes, removed DIM_CHECK if-def ! 1998.10.27 Lucchesi Added support for packing and range checks ! 1999.01.29 Lucchesi Converted API to SD for HDFEOS ! 1999.05.27 Lucchesi Updated error codes. ! subroutine EOS_PutVar ( fid, vname, yyyymmdd, hhmmss, & im, jm, kbeg, kount, grid, & do_comp, do_chunk, rc ) ! ! !USES: Implicit NONE !#include "hlimits.h" ! ! !INPUT PARAMETERS: ! integer fid !! File handle character(len=*) vname !! Variable name integer yyyymmdd !! Year-month-day, e.g., 19971003 integer hhmmss !! Hour-minute-second, e.g., 120000 integer im !! size of longitudinal dimension integer jm !! size of latitudinal dimension integer kbeg !! first level to write; if 2-D grid !! use kbeg = 0. integer kount !! number of levels to write real grid(im,jm,kount) !! Gridded data to write at this time logical do_comp logical do_chunk integer comp_num !! 1 -- COMP_CODE_RLE; 2 -- COMP_CODE_NBIT !! 3 --COMP_CODE_SKPHUFF; 4 -- COMP_CODE_DEFLATE !! 5 --COMP_CODE_SZIP ! !OUTPUT PARAMETERS: integer rc !! Error return code: !! rc = 0 all is well !! rc = -2 time is inconsistent with increment !! rc = -3 number of levels is incompatible with file !! rc = -4 im is incompatible with file !! rc = -5 jm is incompatible with file !! rc = -6 time must fall on a minute boundary !! rc = -7 error in diffdate !! rc = -12 error determining default precision !! rc = -13 error determining variable type !! rc = -15 data outside of valid range !! rc = -16 data outside of packing range !! rc = -17 data outside of pack and valid range !! !! NetCDF Errors !! ------------- !! rc = -32 error detaching from grid !! rc = -37 error attaching to grid (HDFEOS) !! rc = -38 error from NF90_PUT_VAR (dimension variable) NOTUSED !! rc = -40 variable not defined !! rc = -41 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (lat or lon) NOTUSED !! rc = -42 error from NF90_INQ_DIMID or NF90_INQUIRE_DIMENSION (lev) NOTUSED !! rc = -43 error from NF90_INQ_VARID (time variable) !! rc = -44 error reading time information !! rc = -45 error writing data !! rc = -52 error from NF90_INQUIRE_VARIABLE NOTUSED !! rc = -53 error getting variable attributes ! !------------------------------------------------------------------------- integer timeid, dimSize, dimId character(len=MAXCHR) dimName character(len=MAXCHR) attrName character(len=MAXCHR) dimUnits integer attrType, attrCount integer corner(4), edges(4), stride(4) integer dim_chunk(4), origin(4) integer vid integer seconds, timeIndex integer minutes ! added as a work-around integer idx, i, j, k integer begDate, begTime, timInc character(len=8) strBuf integer hour,min,sec,incSecs integer rct ! Variables for dealing with precision real(kind=REAL32), allocatable :: grid_32(:,:,:) real(kind=REAL64), allocatable :: grid_64(:,:,:) real(kind=REAL32) dummy32 real(kind=REAL64) dummy64 real dummy ! Variables for NF90_INQUIRE_VARIABLE character(len=MAXCHR) varName integer type, nvDims, dimSizes(MAX_VAR_DIMS), nvAtts ! Variables for packing and range checking integer(kind=INT16), allocatable :: grid_16(:,:,:) real(kind=REAL32) high_32, low_32, amiss_32 real(kind=REAL32) scale_32, offset_32 logical outRange, outPRange ! Variables for SD interface integer sdsId, sdsIndex, attrIdx, dsdsId integer dataType, numAttr real(kind=REAL64), allocatable :: timeScale(:) real(kind=REAL64), allocatable :: eosTimeScale(:) integer numTimes, fillTime ! Functions for SD interface integer sfgainfo integer sfgdinfo, sffattr, sfrnatt, sfginfo, sfwdata, sfdimid integer sfgdscale integer sfsdscale integer sfn2index, sfselect integer sfrcatt integer sfscompress ! Variables for HDF-EOS integer hdfFid, hdfeos_fid integer gridId integer eosTimeId integer EHidinfo integer GDattach integer GDdetach integer GDwrfld integer comp_type, comp_arg(2) integer COMP_CODE_SZIP, SZ_NN_OPTION_MASK parameter (COMP_CODE_SZIP=5, SZ_NN_OPTION_MASK=32) integer COMP_CODE_DEFLATE parameter (COMP_CODE_DEFLATE = 4) integer COMP_CODE_RLE, COMP_CODE_NBIT, COMP_CODE_SKPHUFF parameter (COMP_CODE_RLE=1, COMP_CODE_NBIT=2, COMP_CODE_SKPHUFF=3) integer DEFLATE_LEVEL parameter (DEFLATE_LEVEL = 6) integer pixels_per_block parameter (pixels_per_block = 12) integer flag, maxcache integer sfwchnk, sfscchnk, sfschnk ! Internal CFIO functions rc = 0 outRange = .FALSE. outPRange = .FALSE. ! Check to make sure max string lengths are large enough. NetCDF defines ! MAXNCNAM, but it can't be used in a character(len=MAXNCNAM) statement. if (MAXCHR .LT. MAXNCNAM) then print *, 'EOS_PutVar warning: MAXNCNAM is larger than ', & 'dimName array size.' endif #if defined (HDFEOS) ! Attach EOS-HDF grid and obtain SD itentifier. Place this ! in 'fid' after saving the HDF-EOS fid in 'hdfeos_fid.' hdfeos_fid = fid gridId = GDattach (hdfeos_fid, GRID_NAME) if (err("PutVar: error in GDattach",rc,-37) .NE. 0) & goto 999 rct = EHidinfo (hdfeos_fid, hdfFid, fid) if (err("PutVar: error in EHidinfo",rc,-37) .NE. 0) goto 999 #endif ! Get SDS id. sdsIndex = sfn2index (fid, vname) if (err("PutVar: can't find variable",rc,-40).NE.0) then print *, 'Error details: no index for ',TRIM(vname) goto 999 endif sdsId = sfselect (fid, sdsIndex) if (err("PutVar: can't find variable",rc,-40).NE.0) then print *, 'Error details: no id for ',TRIM(vname) goto 999 endif ! Dimension error checking do i=0,NDIMS_MAX-1 ! get the dimension ID and infomation about the dimension dimid = sfdimid (sdsId, i) rct = sfgdinfo (dimid, dimName, dimSize, type, numAttr) if (rct .NE. 0) then ! No dimensions left. exit endif ! look for a "units" attribute. if none, use a dummy string for units. dimUnits = REPEAT(' ',MAXCHR) ! zero out dimUnits: sfrcatt won't pad ! with blanks. attrIdx = sffattr (dimid, 'units') if (attrIdx .GE. 0) then rct = sfgainfo (dimid, attrIdx, attrName, attrType, attrCount) if (attrCount .LE. MAXCHR) then rct = sfrcatt (dimid, attrIdx, dimUnits) endif else dimUnits='dummy' endif ! try to identify the dimension based on the name and the units ! then do appropriate error checking idx = IdentifyDim (dimName, dimUnits) if (idx .EQ. 0) then if (dimSize .ne. im) then rc = -4 goto 999 endif else if (idx .EQ. 1) then if (dimSize .ne. jm) then rc = -5 goto 999 endif else if (idx .EQ. 2) then if (kbeg-1 + kount .gt. dimSize) then rc = -3 goto 999 endif else if (idx .EQ. 3) then ! for the time dimension, extract the CFIO-specific time information ! stored as attributes of the time coordinate variable timeid = dimid attrIdx = sffattr (timeid, 'begin_date') rct = sfrnatt (timeid, attrIdx, begDate) if (err("PutVar: error getting begin_date",rc,-44) & .NE. 0) goto 999 attrIdx = sffattr (timeid, 'begin_time') rct = sfrnatt (timeid, attrIdx, begTime) if (err("PutVar: error getting begin_time",rc,-44) & .NE. 0) goto 999 seconds = DiffDate (begDate, begTime, yyyymmdd, hhmmss) if (seconds .lt. 0) then print *, 'EOS_PutVar: Error code from diffdate. Problem ', & 'with date/time.' rc = -7 goto 999 endif if ( MOD (seconds,60) .eq. 0 ) then minutes = seconds / 60 else print *, 'EOS_PutVar: Currently, times must fall on ', & 'minute boundaries.' rc = -6 goto 999 endif attrIdx = sffattr (timeid, 'time_increment') rct = sfrnatt (timeid, attrIdx, timInc) if (err("PutVar: error getting time_increment",rc,-44) & .NE.0) goto 999 ! Convert time increment to seconds. write (strBuf,203) timinc 203 format (I6) read (strBuf,204) hour, min, sec 204 format (3I2) incSecs = hour*3600 + min*60 + sec if ( MOD (seconds, incSecs) .ne. 0 ) then print *, 'CFIO_putvar: Absolute time of ',seconds,' not ', & 'possible with an interval of ',timInc rc = -2 goto 999 else timeIndex = seconds/incSecs + 1 endif else print *, 'EOS_Putvar: WARNING. Dimension ',TRIM(dimName), & ' is unknown.' endif enddo ! Load starting indicies. if ( kbeg .eq. 0 ) then corner(1)=0 corner(2)=0 corner(3)=timeIndex-1 stride(1)=1 stride(2)=1 stride(3)=1 edges(1)=im edges(2)=jm edges(3)=1 dim_chunk(1) = im dim_chunk(2) = jm dim_chunk(3) = 1 else corner(1)=0 corner(2)=0 corner(3)=kbeg-1 corner(4)=timeIndex-1 stride(1)=1 stride(2)=1 stride(3)=1 stride(4)=1 edges(1)=im edges(2)=jm edges(3)=kount edges(4)=1 dim_chunk(1) = im dim_chunk(2) = jm dim_chunk(3) = 1 dim_chunk(4) = 1 endif ! Check variable against valid range. attrIdx = sffattr (sdsId, 'vmin') rct = sfrnatt (sdsId, attrIdx, low_32) if (err("PutVar: error getting vmin",rc,-53) .NE. 0) goto 999 attrIdx = sffattr (sdsId, 'vmax') rct = sfrnatt (sdsId, attrIdx, high_32) if (err("PutVar: error getting vmax",rc,-53) .NE. 0) goto 999 attrIdx = sffattr (sdsId, 'fmissing_value') rct = sfrnatt (sdsId, attrIdx, amiss_32) if (err("PutVar: error getting FILL",rc,-53) .NE. 0) goto 999 if (abs(low_32) .NE. amiss_32 .OR. high_32 .NE. amiss_32) then do k=1,kount do j=1,jm do i=1,im if (grid(i,j,k) .GT. high_32 .OR. grid(i,j,k) .LT. & low_32) then outRange = .TRUE. goto 100 endif enddo enddo enddo 100 continue endif ! Determine if we are writing single- or double-precision from the "type" ! flag. Also get the size of the unlimited dimension (time) for later use. rct = sfginfo (sdsId, varName, nvDims, dimSizes, type, nvAtts) if (err("PutVar: can't get # of times",rc,-44) .NE. 0) & goto 999 numTimes = dimSizes(nvDims) ! DO szip compression comp_num = 4 if (comp_num .eq. 5) then comp_type = COMP_CODE_SZIP comp_arg(1) = SZ_NN_OPTION_MASK comp_arg(2) = pixels_per_block end if if (comp_num .eq. 4) then comp_type = COMP_CODE_DEFLATE comp_arg(1) = DEFLATE_LEVEL end if if ( do_comp ) then if ( do_chunk ) then if ( timeIndex .eq. 1 ) then rct = sfschnk(sdsId, dim_chunk, comp_type, comp_arg) if( rct .ne. 0 ) then print *, 'sfschnk failed for chunk' endif end if ! flag = 0 ! maxcache = 2 ! rct = sfscchnk(sdsId, maxcache, flag) ! if( rct .ne. 0 ) then ! print *, 'sfscchnk failed for chunk' ! endif else rct = sfscompress(sdsId, comp_type, comp_arg) if( rct .ne. 0 ) then print *, 'sfscompress failed ' endif end if end if ! Write variable in the appropriate precision. if (HUGE(dummy) .EQ. HUGE(dummy32)) then ! -r4 if (type .EQ. DFNT_FLOAT32) then ! 32-bit if ( do_chunk ) then if ( kbeg .eq. 0 ) then origin(1) = 1 origin(2) = 1 origin(3) = timeIndex rct = sfwchnk(sdsId, origin, grid(:,:,1)) else do k=1,kount origin(1) = 1 origin(2) = 1 origin(3) = k origin(4) = timeIndex rct = sfwchnk(sdsId, origin, grid(:,:,k)) end do end if if( rct .ne. 0 ) then print *, 'sfwchnk failed for chunk' endif else rct = sfwdata (sdsId, corner, stride, edges, grid) end if if( rct .ne. 0 ) then print *, "return code from sfwdata in CFIO: ", rct end if else if (type .EQ. DFNT_FLOAT64) then ! 64-bit allocate (grid_64(im,jm,kount)) do k=1,kount do j=1,jm do i=1,im grid_64(i,j,k) = grid(i,j,k) enddo enddo enddo if ( do_chunk ) then if ( kbeg .eq. 0 ) then origin(1) = 1 origin(2) = 1 origin(3) = timeIndex rct = sfwchnk(sdsId, origin, grid_64(:,:,1)) else do k=1,kount origin(1) = 1 origin(2) = 1 origin(3) = k origin(4) = timeIndex rct = sfwchnk(sdsId, origin, grid_64(:,:,k)) end do end if if( rct .ne. 0 ) then print *, 'sfwchnk failed for chunk' endif else rct = sfwdata (sdsId, corner, stride, edges, grid_64) end if if( rct .ne. 0 ) then print *, "return code from sfwdata in CFIO: ", rct end if deallocate (grid_64) else if (type .EQ. DFNT_INT16) then attrIdx = sffattr (sdsId, 'packmax') rct = sfrnatt (sdsId, attrIdx, high_32) if (err("PutVar: error getting packmax",rc,-53) .NE. 0) & goto 999 attrIdx = sffattr (sdsId, 'packmin') rct = sfrnatt (sdsId, attrIdx, low_32) if (err("PutVar: error getting packmin",rc,-53) .NE. 0) & goto 999 attrIdx = sffattr (sdsId, 'scale_factor') rct = sfrnatt (sdsId, attrIdx, scale_32) if (err("PutVar: error getting scale_factor",rc,-53) & .NE. 0) goto 999 attrIdx = sffattr (sdsId, 'add_offset') rct = sfrnatt (sdsId, attrIdx, offset_32) if (err("PutVar: error getting add_offset",rc,-53) & .NE. 0) goto 999 allocate (grid_16(im,jm,kount)) do k=1,kount do j=1,jm do i=1,im if ( grid(i,j,k) .LT. low_32 .OR. grid(i,j,k) .GT. & high_32) then grid_16(i,j,k) = PACK_FILL outPRange = .TRUE. else grid_16(i,j,k) = (grid(i,j,k) - offset_32)/scale_32 endif enddo enddo enddo if ( do_chunk ) then if ( kbeg .eq. 0 ) then origin(1) = 1 origin(2) = 1 origin(3) = timeIndex rct = sfwchnk(sdsId, origin, grid_16(:,:,1)) else do k=1,kount origin(1) = 1 origin(2) = 1 origin(3) = k origin(4) = timeIndex rct = sfwchnk(sdsId, origin, grid_16(:,:,k)) end do end if if( rct .ne. 0 ) then print *, 'sfwchnk failed for chunk' endif else rct = sfwdata (sdsId, corner, stride, edges, grid_16) end if deallocate (grid_16) else rc = -13 goto 999 endif else if (HUGE(dummy) .EQ. HUGE(dummy64)) then ! -r8 if (type .EQ. DFNT_FLOAT32) then ! 32-bit allocate (grid_32(im,jm,kount)) do k=1,kount do j=1,jm do i=1,im grid_32(i,j,k) = grid(i,j,k) enddo enddo enddo if ( do_chunk ) then if ( kbeg .eq. 0 ) then origin(1) = 1 origin(2) = 1 origin(3) = timeIndex rct = sfwchnk(sdsId, origin, grid_32(:,:,1)) else do k=1,kount origin(1) = 1 origin(2) = 1 origin(3) = k origin(4) = timeIndex rct = sfwchnk(sdsId, origin, grid_32(:,:,k)) end do end if if( rct .ne. 0 ) then print *, 'sfwchnk failed for chunk' endif else rct = sfwdata (sdsId, corner, stride, edges, grid_32) end if if( rct .ne. 0 ) then print *, "return code from sfwdata in CFIO: ", rct end if deallocate (grid_32) else if (type .EQ. DFNT_FLOAT64) then ! 64-bit if ( do_chunk ) then if ( kbeg .eq. 0 ) then origin(1) = 1 origin(2) = 1 origin(3) = timeIndex rct = sfwchnk(sdsId, origin, grid(:,:,1)) else do k=1,kount origin(1) = 1 origin(2) = 1 origin(3) = k origin(4) = timeIndex rct = sfwchnk(sdsId, origin, grid(:,:,k)) end do end if if( rct .ne. 0 ) then print *, 'sfwchnk failed for chunk' endif else rct = sfwdata (sdsId, corner, stride, edges, grid) end if if( rct .ne. 0 ) then print *, "return code from sfwdata in CFIO: ", rct end if else if (type .EQ. DFNT_INT16) then attrIdx = sffattr (sdsId, 'packmax') rct = sfrnatt (sdsId, attrIdx, high_32) if (err("PutVar: error getting packmax",rc,-53) .NE. 0) & goto 999 attrIdx = sffattr (sdsId, 'packmin') rct = sfrnatt (sdsId, attrIdx, low_32) if (err("PutVar: error getting packmin",rc,-53) .NE. 0) & goto 999 attrIdx = sffattr (sdsId, 'scale_factor') rct = sfrnatt (sdsId, attrIdx, scale_32) if (err("PutVar: error getting scale_factor",rc,-53) & .NE. 0) goto 999 attrIdx = sffattr (sdsId, 'add_offset') rct = sfrnatt (sdsId, attrIdx, offset_32) if (err("PutVar: error getting add_offset",rc,-53) & .NE. 0) goto 999 allocate (grid_16(im,jm,kount)) do k=1,kount do j=1,jm do i=1,im if ( grid(i,j,k) .LT. low_32 .OR. grid(i,j,k) .GT. & high_32) then grid_16(i,j,k) = PACK_FILL else grid_16(i,j,k) = (grid(i,j,k) - offset_32)/scale_32 endif enddo enddo enddo if ( do_chunk ) then if ( kbeg .eq. 0 ) then origin(1) = 1 origin(2) = 1 origin(3) = timeIndex rct = sfwchnk(sdsId, origin, grid_16(:,:,1)) else do k=1,kount origin(1) = 1 origin(2) = 1 origin(3) = k origin(4) = timeIndex rct = sfwchnk(sdsId, origin, grid_16(:,:,k)) end do end if if( rct .ne. 0 ) then print *, 'sfwchnk failed for chunk' endif else rct = sfwdata (sdsId, corner, stride, edges, grid_16) end if deallocate (grid_16) else rc = -13 goto 999 endif else rc = -12 goto 999 endif if (err("PutVar: error writing variable",rc,-45) .NE. 0) & goto 999 ! Write time to file. if (numTimes .GT. timeIndex) then allocate ( timeScale(numTimes) ) allocate ( eosTimeScale(numTimes) ) else allocate ( timeScale(timeIndex) ) allocate ( eosTimeScale(timeIndex) ) endif rct = sfgdscale (timeid, timeScale) if (err("PutVar: can't read times",rc,-44) .NE. 0) goto 999 do i=1,timeIndex-1 fillTime = (i-1) * incSecs/60 if ( timeScale(i) .NE. fillTime ) then timeScale(i) = fillTime endif eosTimeScale(i) = timeScale(i) enddo timeScale(timeIndex) = minutes eosTimeScale(timeIndex) = minutes rct = sfsdscale (timeid, timeIndex, DFNT_FLOAT64, timeScale) if (err("PutVar: can't write times",rc,-44) .NE. 0) goto 999 #if defined (HDFEOS) ! This code will write out placeholders for the EOS time. EOS ! requires seconds from a date in 1993 and also requires leap seconds. ! A post-processing job will use the EOS Toolkit to add this information. sdsIndex = sfn2index (fid, "Time") eosTimeId = sfselect (fid, sdsIndex) if (err("PutVar: can't find EOS time coordinate", & rc,-44) .NE. 0) goto 999 corner(1)=0 stride(1)=1 edges(1)=timeIndex rct = sfwdata (eosTimeId, corner, stride, edges, eosTimeScale) if (err("PutVar: can't write EOS times",rc,-44) .NE. 0) & goto 999 #endif deallocate (timeScale) deallocate (eosTimeScale) if (outRange .AND. outRange) then rc = -17 else if (outPRange) then rc = -16 else if (outRange) then rc = -15 else rc = 0 endif 999 continue #if defined (HDFEOS) fid = hdfeos_fid rct = GDdetach (gridId) #endif return end subroutine EOS_PutVar #if defined (HDFEOS) INTEGER FUNCTION GetSDSid (fid, varName) IMPLICIT NONE integer fid character(len=*) varName integer sdid, rc, idx integer HDFfid, sd_id integer EHidinfo, sfn2index, sfselect rc = EHidinfo (fid, HDFfid, sd_id) idx = sfn2index (sd_id, varName) GetSDSid = sfselect (sd_id, idx) return end FUNCTION GetSDSid #endif #endif end module ESMF_CFIOUtilMod