CFIO_SGetVar Subroutine

public subroutine CFIO_SGetVar(fid, vname, yyyymmdd, hhmmss, im, jm, kbeg, kount, lm, grid, cyclic, rc)

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

Arguments

Type IntentOptional Attributes Name
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

real :: grid(im,kount)

Gridded data read for this time

logical :: cyclic

input file is cyclic or not

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


Calls

proc~~cfio_sgetvar~~CallsGraph proc~cfio_sgetvar CFIO_SGetVar nf90_get_att nf90_get_att proc~cfio_sgetvar->nf90_get_att nf90_get_var nf90_get_var proc~cfio_sgetvar->nf90_get_var nf90_inq_varid nf90_inq_varid proc~cfio_sgetvar->nf90_inq_varid nf90_inquire nf90_inquire proc~cfio_sgetvar->nf90_inquire nf90_inquire_dimension nf90_inquire_dimension proc~cfio_sgetvar->nf90_inquire_dimension nf90_inquire_variable nf90_inquire_variable proc~cfio_sgetvar->nf90_inquire_variable proc~diffdate DiffDate proc~cfio_sgetvar->proc~diffdate proc~err err proc~cfio_sgetvar->proc~err proc~getdatetimevec GetDateTimeVec proc~cfio_sgetvar->proc~getdatetimevec proc~identifydim IdentifyDim proc~cfio_sgetvar->proc~identifydim proc~cfio_parseinttime CFIO_parseIntTime proc~diffdate->proc~cfio_parseinttime proc~julday julday proc~diffdate->proc~julday proc~getdatetimevec->nf90_get_att proc~getdatetimevec->nf90_get_var proc~getdatetimevec->nf90_inq_varid proc~getdatetimevec->nf90_inquire proc~getdatetimevec->nf90_inquire_dimension proc~getdatetimevec->nf90_inquire_variable proc~getdatetimevec->proc~err proc~getdatetimevec->proc~identifydim interface~getdate GetDate proc~getdatetimevec->interface~getdate proc~parsetimeunits ParseTimeUnits proc~getdatetimevec->proc~parsetimeunits proc~getdateint4 GetDateInt4 interface~getdate->proc~getdateint4 proc~getdateint8 GetDateInt8 interface~getdate->proc~getdateint8 proc~getdateint4->proc~getdateint8 proc~getdateint8->proc~cfio_parseinttime proc~getdateint8->proc~julday proc~caldat CALDAT proc~getdateint8->proc~caldat

Source Code

      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