mpp_io_read.inc

Go to the documentation of this file.

! -*-f90-*-

!***********************************************************************
!*                   GNU Lesser General Public License
!*
!* This file is part of the GFDL Flexible Modeling System (FMS).
!*
!* FMS is free software: you can redistribute it and/or modify it under
!* the terms of the GNU Lesser General Public License as published by
!* the Free Software Foundation, either version 3 of the License, or (at
!* your option) any later version.
!*
!* FMS is distributed in the hope that it will be useful, but WITHOUT
!* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
!* FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
!* for more details.
!*
!* You should have received a copy of the GNU Lesser General Public
!* License along with FMS.  If not, see <http://www.gnu.org/licenses/>.
!***********************************************************************

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!                                                                      !
!                               MPP_READ                               !
!                                                                      !
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

#undef MPP_READ_2DDECOMP_2D_
#undef READ_RECORD_CORE_
#define READ_RECORD_CORE_ read_record_core
#undef READ_RECORD_
#define READ_RECORD_ read_record
#define MPP_READ_2DDECOMP_2D_ mpp_read_2ddecomp_r2d
#undef MPP_READ_2DDECOMP_3D_
#define MPP_READ_2DDECOMP_3D_ mpp_read_2ddecomp_r3d
#undef MPP_READ_2DDECOMP_4D_
#define MPP_READ_2DDECOMP_4D_ mpp_read_2ddecomp_r4d
#undef MPP_TYPE_
#define MPP_TYPE_ real
#include <mpp_read_2Ddecomp.h>

#ifdef OVERLOAD_R8
#undef READ_RECORD_CORE_
#define READ_RECORD_CORE_ read_record_core_r8
#undef READ_RECORD_
#define READ_RECORD_ read_record_r8
#undef MPP_READ_2DDECOMP_2D_
#define MPP_READ_2DDECOMP_2D_ mpp_read_2ddecomp_r2d_r8
#undef MPP_READ_2DDECOMP_3D_
#define MPP_READ_2DDECOMP_3D_ mpp_read_2ddecomp_r3d_r8
#undef MPP_READ_2DDECOMP_4D_
#define MPP_READ_2DDECOMP_4D_ mpp_read_2ddecomp_r4d_r8
#undef MPP_TYPE_
#define MPP_TYPE_ real(DOUBLE_KIND)
#include <mpp_read_2Ddecomp.h>
#endif

#undef MPP_READ_COMPRESSED_1D_
#define MPP_READ_COMPRESSED_1D_ mpp_read_compressed_r1d
#undef MPP_READ_COMPRESSED_2D_
#define MPP_READ_COMPRESSED_2D_ mpp_read_compressed_r2d
#undef MPP_READ_COMPRESSED_3D_
#define MPP_READ_COMPRESSED_3D_ mpp_read_compressed_r3d
#undef MPP_TYPE_
#define MPP_TYPE_ real
#include <mpp_read_compressed.h>

#include <mpp_read_distributed_ascii.inc>

! <SUBROUTINE NAME="mpp_read_r4D" INTERFACE="mpp_read">
!   <IN NAME="unit" TYPE="integer"></IN>
!   <IN NAME="field" TYPE="type(fieldtype)"></IN>
!   <INOUT NAME="data" TYPE="real" DIM="(:,:,:,:)"></INOUT>
!   <IN NAME="tindex" TYPE="integer"></IN>
! </SUBROUTINE>
    subroutine mpp_read_r4D( unit, field, data, tindex)
      integer, intent(in) :: unit
      type(fieldtype), intent(in) :: field
      real, intent(inout) :: data(:,:,:,:)
      integer, intent(in), optional :: tindex

      call read_record( unit, field, size(data(:,:,:,:)), data, tindex )
    end subroutine mpp_read_r4D


! <SUBROUTINE NAME="mpp_read_r3D" INTERFACE="mpp_read">
!   <IN NAME="unit" TYPE="integer"></IN>
!   <IN NAME="field" TYPE="type(fieldtype)"></IN>
!   <INOUT NAME="data" TYPE="real" DIM="(:,:,:)"></INOUT>
!   <IN NAME="tindex" TYPE="integer"></IN>
! </SUBROUTINE>
    subroutine mpp_read_r3D( unit, field, data, tindex)
      integer, intent(in) :: unit
      type(fieldtype), intent(in) :: field
      real, intent(inout) :: data(:,:,:)
      integer, intent(in), optional :: tindex

      call read_record( unit, field, size(data(:,:,:)), data, tindex )
    end subroutine mpp_read_r3D

    subroutine mpp_read_r2D( unit, field, data, tindex )
      integer, intent(in) :: unit
      type(fieldtype), intent(in) :: field
      real, intent(inout) :: data(:,:)
      integer, intent(in), optional :: tindex

      call read_record( unit, field, size(data(:,:)), data, tindex )
    end subroutine mpp_read_r2D

    subroutine mpp_read_r1D( unit, field, data, tindex )
      integer, intent(in) :: unit
      type(fieldtype), intent(in) :: field
      real, intent(inout) :: data(:)
      integer, intent(in), optional :: tindex

      call read_record( unit, field, size(data(:)), data, tindex )
    end subroutine mpp_read_r1D

    subroutine mpp_read_r0D( unit, field, data, tindex )
      integer, intent(in) :: unit
      type(fieldtype), intent(in) :: field
      real, intent(inout) :: data
      integer, intent(in), optional :: tindex
      real, dimension(1) :: data_tmp

      data_tmp(1)=data
      call read_record( unit, field, 1, data_tmp, tindex )
      data=data_tmp(1)
    end subroutine mpp_read_r0D

    subroutine mpp_read_region_r2D(unit, field, data, start, nread)
      integer,         intent(in) :: unit
      type(fieldtype), intent(in) :: field
      real,         intent(inout) :: data(:,:)
      integer,         intent(in) :: start(:), nread(:)

      if(size(start(:)) .NE. 4 .OR. size(nread(:)) .NE. 4) call mpp_error(FATAL, &
          "mpp_io_read.inc(mpp_read_region_r2D): size of start and nread must be 4")

      if(size(data,1) .NE. nread(1) .OR. size(data,2) .NE. nread(2)) then
         call mpp_error( FATAL, 'mpp_io_read.inc(mpp_read_block_r2D): size mismatch between data and nread')
      endif
      if(nread(3) .NE. 1 .OR. nread(4) .NE. 1) call mpp_error(FATAL, &
          "mpp_io_read.inc(mpp_read_region_r2D): nread(3) and nread(4) must be 1")
      call read_record_core(unit, field, nread(1)*nread(2), data, start, nread)

      return


    end subroutine mpp_read_region_r2D

    subroutine mpp_read_region_r3D(unit, field, data, start, nread)
      integer,         intent(in) :: unit
      type(fieldtype), intent(in) :: field
      real,         intent(inout) :: data(:,:,:)
      integer,         intent(in) :: start(:), nread(:)

      if(size(start(:)) .NE. 4 .OR. size(nread(:)) .NE. 4) call mpp_error(FATAL, &
          "mpp_io_read.inc(mpp_read_region_r3D): size of start and nread must be 4")

      if(size(data,1) .NE. nread(1) .OR. size(data,2) .NE. nread(2) .OR. size(data,3) .NE. nread(3) ) then
         call mpp_error( FATAL, 'mpp_io_read.inc(mpp_read_block_r3D): size mismatch between data and nread')
      endif
      if(nread(4) .NE. 1) call mpp_error(FATAL, &
          "mpp_io_read.inc(mpp_read_region_r3D): nread(4) must be 1")
      call read_record_core(unit, field, nread(1)*nread(2)*nread(3), data, start, nread)

      return
    end subroutine mpp_read_region_r3D

#ifdef OVERLOAD_R8
    subroutine mpp_read_region_r2D_r8(unit, field, data, start, nread)
      integer,                intent(in) :: unit
      type(fieldtype),        intent(in) :: field
      real(kind=DOUBLE_KIND), intent(inout) :: data(:,:)
      integer,                intent(in) :: start(:), nread(:)

      if(size(start(:)) .NE. 4 .OR. size(nread(:)) .NE. 4) call mpp_error(FATAL, &
          "mpp_io_read.inc(mpp_read_region_r2D_r8): size of start and nread must be 4")

      if(size(data,1).NE.nread(1) .OR. size(data,2).NE.nread(2)) then
         call mpp_error( FATAL, 'mpp_io_read.inc(mpp_read_block_r2D_r8): size mismatch between data and nread')
      endif
      if(nread(3) .NE. 1 .OR. nread(4) .NE. 1) call mpp_error(FATAL, &
          "mpp_io_read.inc(mpp_read_region_r2D_r8): nread(3) and nread(4) must be 1")
      call read_record_core_r8(unit, field, nread(1)*nread(2), data, start, nread)

      return
    end subroutine mpp_read_region_r2D_r8

    subroutine mpp_read_region_r3D_r8(unit, field, data, start, nread)
      integer,                intent(in) :: unit
      type(fieldtype),        intent(in) :: field
      real(kind=DOUBLE_KIND), intent(inout) :: data(:,:,:)
      integer,                intent(in) :: start(:), nread(:)

      if(size(start(:)) .NE. 4 .OR. size(nread(:)) .NE. 4) call mpp_error(FATAL, &
          "mpp_io_read.inc(mpp_read_region_r3D_r8): size of start and nread must be 4")

      if(size(data,1).NE.nread(1) .OR. size(data,2).NE.nread(2) .OR. size(data,3).NE.nread(3) ) then
         call mpp_error( FATAL, 'mpp_io_read.inc(mpp_read_block_r3D_r8): size mismatch between data and nread')
      endif
      if(nread(4) .NE. 1) call mpp_error(FATAL, &
          "mpp_io_read.inc(mpp_read_region_r3D_r8): nread(4) must be 1")
      call read_record_core_r8(unit, field, nread(1)*nread(2)*nread(3), data, start, nread)

      return
    end subroutine mpp_read_region_r3D_r8
#endif


    !--- Assume the text field is at most two-dimensional
    !--- the level is always for the first dimension
    subroutine mpp_read_text( unit, field, data, level )
      integer,             intent(in) :: unit
      type(fieldtype),     intent(in) :: field
      character(len=*), intent(inout) :: data
      integer, intent(in), optional   :: level
      integer                         :: lev, n
      character(len=256)              :: error_msg
      integer, dimension(size(field%axes(:))) :: start, axsiz
      character(len=len(data))        :: text

#ifdef use_netCDF
      if( .NOT.module_is_initialized )call mpp_error( FATAL, 'READ_RECORD: must first call mpp_io_init.' )
      if( .NOT.mpp_file(unit)%opened )call mpp_error( FATAL, 'READ_RECORD: invalid unit number.' )
      if( mpp_file(unit)%threading.EQ.MPP_SINGLE .AND. pe.NE.mpp_root_pe() )return

      if( .NOT.mpp_file(unit)%initialized ) call mpp_error( FATAL, 'MPP_READ: must first call mpp_read_meta.' )
      lev = 1
      if(present(level)) lev = level

      if( verbose )print '(a,2i6,2i5)', 'MPP_READ: PE, unit, %id, level =', pe, unit, mpp_file(unit)%id, lev

      if( mpp_file(unit)%format.EQ.MPP_NETCDF )then
         start = 1
         axsiz(:) = field%size(:)
         if(len(data) < field%size(1) ) call mpp_error(FATAL, &
                'mpp_io(mpp_read_text): the first dimension size is greater than data length')
         select case( field%ndim)
         case(1)
            if(lev .NE. 1) call mpp_error(FATAL,'mpp_io(mpp_read_text): level should be 1 when ndim is 1')
         case(2)
            if(lev<1 .OR. lev > field%size(2)) then
               write(error_msg,'(I5,"/",I5)') lev, field%size(2)
               call mpp_error(FATAL,'mpp_io(mpp_read_text): level out of range, level/max_level='//trim(error_msg))
            end if
            start(2) = lev
            axsiz(2) = 1
         case default
            call mpp_error( FATAL, 'MPP_READ: ndim of text field should be at most 2')
         end select

         if( verbose )print '(a,2i6,i6,12i4)', 'mpp_read_text: PE, unit, nwords, start, axsiz=', pe, unit, len(data), start, axsiz

          select case (field%type)
             case(NF_CHAR)
                if(field%ndim==1) then
                   error = NF_GET_VAR_TEXT(mpp_file(unit)%ncid, field%id, text)
                else
                   error = NF_GET_VARA_TEXT(mpp_file(unit)%ncid, field%id, start, axsiz, text)
                end if
                call netcdf_err( error, mpp_file(unit), field=field )
                do n = 1, len_trim(text)
                   if(text(n:n) == CHAR(0) ) exit
                end do
                n = n-1
                data = text(1:n)
             case default
                call mpp_error( FATAL, 'mpp_read_text: the field type should be NF_CHAR' )
          end select
      else                      !non-netCDF
          call mpp_error( FATAL, 'Currently dont support non-NetCDF mpp read' )

      end if
#else
      call mpp_error( FATAL, 'mpp_read_text currently requires use_netCDF option' )
#endif
      return
    end subroutine mpp_read_text

! <SUBROUTINE NAME="mpp_read_meta">

!   <OVERVIEW>
!     Read metadata.
!   </OVERVIEW>
!   <DESCRIPTION>
!     This routine is used to read the <LINK SRC="#metadata">metadata</LINK>
!     describing the contents of a file. Each file can contain any number of
!     fields, which are functions of 0-3 space axes and 0-1 time axes. (Only
!     one time axis can be defined per file). The basic metadata defined <LINK
!     SRC="#metadata">above</LINK> for <TT>axistype</TT> and
!     <TT>fieldtype</TT> are stored in <TT>mpp_io_mod</TT> and
!     can be accessed outside of <TT>mpp_io_mod</TT> using calls to
!     <TT>mpp_get_info</TT>, <TT>mpp_get_atts</TT>,
!     <TT>mpp_get_vars</TT> and
!     <TT>mpp_get_times</TT>.
!   </DESCRIPTION>
!   <TEMPLATE>
!     call mpp_read_meta(unit)
!   </TEMPLATE>
!   <IN NAME="unit" TYPE="integer"> </IN>
!   <NOTE>
!     <TT>mpp_read_meta</TT> must be called prior to <TT>mpp_read</TT>.
!   </NOTE>
! </SUBROUTINE>
    subroutine mpp_read_meta(unit, read_time)
!
! read file attributes including dimension and variable attributes
! and store in filetype structure.  All of the file information
! with the exception of the (variable) data is stored.  Attributes
! are supplied to the user by get_info,get_atts,get_axes and get_fields
!
! every PE is eligible to call mpp_read_meta
!
      integer, intent(in) :: unit
      logical, intent(in), optional :: read_time ! read_time is set to false when file action is appending or writing.
                                                 ! This is temporary fix for MOM6 reopen_file issue.
      integer         :: ncid,ndim,nvar_total,natt,recdim,nv,nvar,len
      integer :: error, i, j, istat, check_exist
      integer         :: type, nvdims, nvatts, dimid
      integer, allocatable, dimension(:) :: dimids
      character(len=128) :: name, attname, unlimname, attval, bounds_name
      logical :: isdim, found_bounds, get_time_info
      integer(LONG_KIND) :: checksumf
      character(len=64)  :: checksum_char
      integer :: num_checksumf, last, is, k

      integer(SHORT_KIND), allocatable :: i2vals(:)
      integer(INT_KIND), allocatable :: ivals(:)
      real(FLOAT_KIND), allocatable  :: rvals(:)
      real(DOUBLE_KIND), allocatable :: r8vals(:)

      get_time_info = .TRUE.
      if(present(read_time)) get_time_info = read_time

#ifdef use_netCDF

      if( mpp_file(unit)%format.EQ.MPP_NETCDF )then
        ncid = mpp_file(unit)%ncid
        error = NF_INQ(ncid,ndim, nvar_total,&
                      natt, recdim);call netcdf_err( error, mpp_file(unit) )


        mpp_file(unit)%ndim = ndim
        mpp_file(unit)%natt = natt
        mpp_file(unit)%recdimid = recdim
!
! if no recdim exists, recdimid = -1
! variable id of unlimdim and length
!
        if( recdim.NE.-1 )then
           error = NF_INQ_DIM( ncid, recdim, unlimname, mpp_file(unit)%time_level )
           call netcdf_err( error, mpp_file(unit) )
           error = NF_INQ_VARID( ncid, unlimname, mpp_file(unit)%id )
           call netcdf_err( error, mpp_file(unit), string='Field='//unlimname )
        else
           mpp_file(unit)%time_level = -1 ! set to zero so mpp_get_info returns ntime=0 if no time axis present
        endif

        allocate(mpp_file(unit)%Att(natt))
        allocate(dimids(ndim))
        allocate(mpp_file(unit)%Axis(ndim))

!
! initialize fieldtype and axis type
!


        do i=1,ndim
           mpp_file(unit)%Axis(i) = default_axis
        enddo

        do i=1,natt
           mpp_file(unit)%Att(i) = default_att
        enddo

!
! assign global attributes
!
        do i=1,natt
           error=NF_INQ_ATTNAME(ncid,NF_GLOBAL,i,name);call netcdf_err( error, mpp_file(unit), string=' Global attribute error.' )
           error=NF_INQ_ATT(ncid,NF_GLOBAL,trim(name),type,len);call netcdf_err( error, mpp_file(unit), string=' Attribute='//name )
           mpp_file(unit)%Att(i)%name = name
           mpp_file(unit)%Att(i)%len = len
           mpp_file(unit)%Att(i)%type = type
!
!  allocate space for att data and assign
!
           select case (type)
              case (NF_CHAR)
                 if (len.gt.MAX_ATT_LENGTH) then
                    call mpp_error(NOTE,'GLOBAL ATT too long - not reading this metadata')
                    len=7
                    mpp_file(unit)%Att(i)%len=len
                    mpp_file(unit)%Att(i)%catt = 'unknown'
                 else
                     error=NF_GET_ATT_TEXT(ncid,NF_GLOBAL,name,mpp_file(unit)%Att(i)%catt)
                     call netcdf_err( error, mpp_file(unit), attr=mpp_file(unit)%att(i) )
                     if (verbose.and.pe == 0) print *, 'GLOBAL ATT ',trim(name),' ',mpp_file(unit)%Att(i)%catt(1:len)
                 endif
!
! store integers in float arrays
!
              case (NF_SHORT)
                 allocate(mpp_file(unit)%Att(i)%fatt(len), STAT=istat)
                 if ( istat .ne. 0 ) then
                    write(text,'(A)') istat
                    call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Att%fatt, NF_SHORT case.  "//&
                         & "STAT = "//trim(text))
                 end if
                 allocate(i2vals(len), STAT=istat)
                 if ( istat .ne. 0 ) then
                    write(text,'(A)') istat
                    call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array i2vals.  STAT = "&
                         & //trim(text))
                 end if
                 error=NF_GET_ATT_INT2(ncid,NF_GLOBAL,name,i2vals)
                 call netcdf_err( error, mpp_file(unit), attr=mpp_file(unit)%att(i) )
                 if( verbose .and. pe == 0 )print *, 'GLOBAL ATT ',trim(name),' ',i2vals(1:len)
                 mpp_file(unit)%Att(i)%fatt(1:len)=i2vals(1:len)
                 deallocate(i2vals)
              case (NF_INT)
                 allocate(mpp_file(unit)%Att(i)%fatt(len), STAT=istat)
                 if ( istat .ne. 0 ) then
                    write(text,'(A)') istat
                    call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Att%fatt, NF_INT case.  "//&
                         & "STAT = "//trim(text))
                 end if
                 allocate(ivals(len), STAT=istat)
                 if ( istat .ne. 0 ) then
                    write(text,'(A)') istat
                    call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array ivals.  STAT = "&
                         & //trim(text))
                 end if
                 error=NF_GET_ATT_INT(ncid,NF_GLOBAL,name,ivals)
                 call netcdf_err( error, mpp_file(unit), attr=mpp_file(unit)%att(i) )
                 if( verbose .and. pe == 0 )print *, 'GLOBAL ATT ',trim(name),' ',ivals(1:len)
                 mpp_file(unit)%Att(i)%fatt(1:len)=ivals(1:len)
                 if(lowercase(trim(name)) == 'time_axis' .and. ivals(1)==0) &
                        mpp_file(unit)%time_level = -1 ! This file is an unlimited axis restart
                 deallocate(ivals)
              case (NF_FLOAT)
                 allocate(mpp_file(unit)%Att(i)%fatt(len), STAT=istat)
                 if ( istat .ne. 0 ) then
                    write(text,'(A)') istat
                    call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Att%fatt, NF_FLOAT case.  "//&
                         & "STAT = "//trim(text))
                 end if
                 allocate(rvals(len), STAT=istat)
                 if ( istat .ne. 0 ) then
                    write(text,'(A)') istat
                    call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array rvals.  STAT = "&
                         & //trim(text))
                 end if
                 error=NF_GET_ATT_REAL(ncid,NF_GLOBAL,name,rvals)
                 call netcdf_err( error, mpp_file(unit), attr=mpp_file(unit)%att(i) )
                 mpp_file(unit)%Att(i)%fatt(1:len)=rvals(1:len)
                 if( verbose .and. pe == 0)print *, 'GLOBAL ATT ',trim(name),' ',mpp_file(unit)%Att(i)%fatt(1:len)
                 deallocate(rvals)
              case (NF_DOUBLE)
                 allocate(mpp_file(unit)%Att(i)%fatt(len), STAT=istat)
                 if ( istat .ne. 0 ) then
                    write(text,'(A)') istat
                    call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Att%fatt, NF_DOUBLE case.  "//&
                         & "STAT = "//trim(text))
                 end if
                 allocate(r8vals(len), STAT=istat)
                 if ( istat .ne. 0 ) then
                    write(text,'(A)') istat
                    call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array r8vals.  STAT = "&
                         & //trim(text))
                 end if
                 error=NF_GET_ATT_DOUBLE(ncid,NF_GLOBAL,name,r8vals)
                 call netcdf_err( error, mpp_file(unit), attr=mpp_file(unit)%att(i) )
                 mpp_file(unit)%Att(i)%fatt(1:len)=r8vals(1:len)
                 if( verbose .and. pe == 0)print *, 'GLOBAL ATT ',trim(name),' ',mpp_file(unit)%Att(i)%fatt(1:len)
                 deallocate(r8vals)
           end select

        enddo
!
! assign dimension name and length
!
        do i=1,ndim
           error = NF_INQ_DIM(ncid,i,name,len);call netcdf_err( error, mpp_file(unit) )
           mpp_file(unit)%Axis(i)%name = name
           mpp_file(unit)%Axis(i)%len = len
        enddo

        nvar=0
        do i=1, nvar_total
           error=NF_INQ_VAR(ncid,i,name,type,nvdims,dimids,nvatts);call netcdf_err( error, mpp_file(unit) )
           isdim=.false.
           do j=1,ndim
              if( trim(lowercase(name)).EQ.trim(lowercase(mpp_file(unit)%Axis(j)%name)) )isdim=.true.
           enddo
           if (.not.isdim) nvar=nvar+1
        enddo
        mpp_file(unit)%nvar = nvar
        allocate(mpp_file(unit)%Var(nvar))

        do i=1,nvar
           mpp_file(unit)%Var(i) = default_field
        enddo

!
! assign dimension info
!
        do i=1, nvar_total
           error=NF_INQ_VAR(ncid,i,name,type,nvdims,dimids,nvatts);call netcdf_err( error, mpp_file(unit) )
           isdim=.false.
           do j=1,ndim
              if( trim(lowercase(name)).EQ.trim(lowercase(mpp_file(unit)%Axis(j)%name)) )isdim=.true.
           enddo

           if( isdim )then
              error=NF_INQ_DIMID(ncid,name,dimid);call netcdf_err( error, mpp_file(unit), string=' Axis='//name )
              mpp_file(unit)%Axis(dimid)%type = type
              mpp_file(unit)%Axis(dimid)%did = dimid
              mpp_file(unit)%Axis(dimid)%id = i
              mpp_file(unit)%Axis(dimid)%natt = nvatts
              ! get axis values
              if( i.NE.mpp_file(unit)%id )then   ! non-record dims
                 select case (type)
                 case (NF_INT)
                    len=mpp_file(unit)%Axis(dimid)%len
                    allocate(mpp_file(unit)%Axis(dimid)%data(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Axis%data, NF_INT case.  "//&
                            & "STAT = "//trim(text))
                    end if
                    allocate(ivals(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array ivals.  STAT = "&
                            & //trim(text))
                    end if
                    error = NF_GET_VAR_INT(ncid,i,ivals);call netcdf_err( error, mpp_file(unit), mpp_file(unit)%Axis(dimid) )
                    mpp_file(unit)%Axis(dimid)%data(1:len)=ivals(1:len)
                    deallocate(ivals)
                 case (NF_FLOAT)
                    len=mpp_file(unit)%Axis(dimid)%len
                    allocate(mpp_file(unit)%Axis(dimid)%data(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Axis%data, "//&
                            & "NF_FLOAT case.  STAT = "//trim(text))
                    end if
                    allocate(rvals(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array rvals.  STAT = "&
                            & //trim(text))
                    end if
                    error = NF_GET_VAR_REAL(ncid,i,rvals);call netcdf_err( error, mpp_file(unit), mpp_file(unit)%Axis(dimid) )
                    mpp_file(unit)%Axis(dimid)%data(1:len)=rvals(1:len)
                    deallocate(rvals)
                 case (NF_DOUBLE)
                    len=mpp_file(unit)%Axis(dimid)%len
                    allocate(mpp_file(unit)%Axis(dimid)%data(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Axis%data, "//&
                            & "NF_DOUBLE case.  STAT = "//trim(text))
                    end if
                    allocate(r8vals(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array r8vals.  STAT = "&
                            & //trim(text))
                    end if
                    error = NF_GET_VAR_DOUBLE(ncid,i,r8vals);call netcdf_err( error, mpp_file(unit), mpp_file(unit)%Axis(dimid) )
                    mpp_file(unit)%Axis(dimid)%data(1:len) = r8vals(1:len)
                    deallocate(r8vals)
                 case default
                    call mpp_error( FATAL, 'Invalid data type for dimension' )
                 end select
             else if(get_time_info) then
                 len = mpp_file(unit)%time_level
                 if( len > 0 ) then
                    allocate(mpp_file(unit)%time_values(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%time_valuse.  STAT = "&
                            & //trim(text))
                    end if
                    select case (type)
                    case (NF_FLOAT)
                       allocate(rvals(len), STAT=istat)
                       if ( istat .ne. 0 ) then
                          write(text,'(A)') istat
                          call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array rvals.  STAT = "&
                               & //trim(text))
                       end if
                       !z1l read from root pe and broadcast to other processor.
                       !In the future we will modify the code if there is performance issue for very high MPI ranks.
                       if(mpp_pe()==mpp_root_pe()) then
                          error = NF_GET_VAR_REAL(ncid,i,rvals)
                          call netcdf_err( error, mpp_file(unit), mpp_file(unit)%Axis(dimid) )
                       endif
                       call mpp_broadcast(rvals, len, mpp_root_pe())
                       mpp_file(unit)%time_values(1:len) = rvals(1:len)
                       deallocate(rvals)
                    case (NF_DOUBLE)
                       allocate(r8vals(len), STAT=istat)
                       if ( istat .ne. 0 ) then
                          write(text,'(A)') istat
                          call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array r8vals.  STAT = "&
                            & //trim(text))
                       end if
                       !z1l read from root pe and broadcast to other processor.
                       !In the future we will modify the code if there is performance issue for very high MPI ranks.
                       if(mpp_pe()==mpp_root_pe()) then
                         error = NF_GET_VAR_DOUBLE(ncid,i,r8vals)
                         call netcdf_err( error, mpp_file(unit), mpp_file(unit)%Axis(dimid) )
                       endif
                       call mpp_broadcast(r8vals, len, mpp_root_pe())
                       mpp_file(unit)%time_values(1:len) = r8vals(1:len)
                       deallocate(r8vals)
                    case default
                       call mpp_error( FATAL, 'Invalid data type for dimension' )
                    end select
                 endif
              endif
              ! assign dimension atts
              if( nvatts.GT.0 )allocate(mpp_file(unit)%Axis(dimid)%Att(nvatts))

              do j=1,nvatts
                 mpp_file(unit)%Axis(dimid)%Att(j) = default_att
              enddo

              do j=1,nvatts
                 error=NF_INQ_ATTNAME(ncid,i,j,attname);call netcdf_err( error, mpp_file(unit) )
                 error=NF_INQ_ATT(ncid,i,trim(attname),type,len)
                 call netcdf_err( error, mpp_file(unit), string=' Attribute='//attname )

                 mpp_file(unit)%Axis(dimid)%Att(j)%name = trim(attname)
                 mpp_file(unit)%Axis(dimid)%Att(j)%type = type
                 mpp_file(unit)%Axis(dimid)%Att(j)%len = len

                 select case (type)
                 case (NF_CHAR)
                    if (len.gt.MAX_ATT_LENGTH) call mpp_error(FATAL,'DIM ATT too long')
                    error=NF_GET_ATT_TEXT(ncid,i,trim(attname),mpp_file(unit)%Axis(dimid)%Att(j)%catt);
                    call netcdf_err( error, mpp_file(unit), attr=mpp_file(unit)%Axis(dimid)%att(j) )
                    if( verbose .and. pe == 0 ) &
                         print *, 'AXIS ',trim(mpp_file(unit)%Axis(dimid)%name),' ATT ',trim(attname),' ',mpp_file(unit)%Axis(dimid)%Att(j)%catt(1:len)
                    ! store integers in float arrays
                    ! assume dimension data not packed
                 case (NF_SHORT)
                    allocate(mpp_file(unit)%Axis(dimid)%Att(j)%fatt(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Axis%Att%fatt, "//&
                            & "NF_SHORT CASE.  STAT = "//trim(text))
                    end if
                    allocate(i2vals(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array i2vals.  STAT = "&
                            & //trim(text))
                    end if
                    error=NF_GET_ATT_INT2(ncid,i,trim(attname),i2vals);
                    call netcdf_err( error, mpp_file(unit), attr=mpp_file(unit)%Axis(dimid)%att(j) )
                    mpp_file(unit)%Axis(dimid)%Att(j)%fatt(1:len)=i2vals(1:len)
                    if( verbose .and. pe == 0  ) &
                         print *, 'AXIS ',trim(mpp_file(unit)%Axis(dimid)%name),' ATT ',trim(attname),' ',mpp_file(unit)&
                         &%Axis(dimid)%Att(j)%fatt
                    deallocate(i2vals)
                 case (NF_INT)
                    allocate(mpp_file(unit)%Axis(dimid)%Att(j)%fatt(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Axis%Att%fatt, "//&
                            & "NF_INT CASE.  STAT = "//trim(text))
                    end if
                    allocate(ivals(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array ivals.  STAT = "&
                            & //trim(text))
                    end if
                    error=NF_GET_ATT_INT(ncid,i,trim(attname),ivals);
                    call netcdf_err( error, mpp_file(unit), attr=mpp_file(unit)%Axis(dimid)%att(j) )
                    mpp_file(unit)%Axis(dimid)%Att(j)%fatt(1:len)=ivals(1:len)
                    if( verbose .and. pe == 0  ) &
                         print *, 'AXIS ',trim(mpp_file(unit)%Axis(dimid)%name),' ATT ',trim(attname),' ',&
                         & mpp_file(unit)%Axis(dimid)%Att(j)%fatt
                    deallocate(ivals)
                 case (NF_FLOAT)
                    allocate(mpp_file(unit)%Axis(dimid)%Att(j)%fatt(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Axis%Att%fatt, "//&
                            & "NF_FLOAT CASE.  STAT = "//trim(text))
                    end if
                    allocate(rvals(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array rvals.  STAT = "&
                            & //trim(text))
                    end if
                    error=NF_GET_ATT_REAL(ncid,i,trim(attname),rvals);
                    call netcdf_err( error, mpp_file(unit), attr=mpp_file(unit)%Axis(dimid)%att(j) )
                    mpp_file(unit)%Axis(dimid)%Att(j)%fatt(1:len)=rvals(1:len)
                    if( verbose  .and. pe == 0 ) &
                         print *, 'AXIS ',trim(mpp_file(unit)%Axis(dimid)%name),' ATT ',trim(attname),' ',&
                         & mpp_file(unit)%Axis(dimid)%Att(j)%fatt
                    deallocate(rvals)
                 case (NF_DOUBLE)
                    allocate(mpp_file(unit)%Axis(dimid)%Att(j)%fatt(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Axis%Att%fatt, "//&
                            & "NF_DOUBLE CASE.  STAT = "//trim(text))
                    end if
                    allocate(r8vals(len), STAT=istat)
                    if ( istat .ne. 0 ) then
                       write(text,'(A)') istat
                       call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array r8vals.  STAT = "&
                            & //trim(text))
                    end if
                    error=NF_GET_ATT_DOUBLE(ncid,i,trim(attname),r8vals);
                    call netcdf_err( error, mpp_file(unit), attr=mpp_file(unit)%Axis(dimid)%att(j) )
                    mpp_file(unit)%Axis(dimid)%Att(j)%fatt(1:len)=r8vals(1:len)
                    if( verbose  .and. pe == 0 ) &
                         print *, 'AXIS ',trim(mpp_file(unit)%Axis(dimid)%name),' ATT ',trim(attname),' ',&
                         & mpp_file(unit)%Axis(dimid)%Att(j)%fatt
                    deallocate(r8vals)
                 case default
                    call mpp_error( FATAL, 'Invalid data type for dimension at' )
                 end select
                 ! assign pre-defined axis attributes
                 select case(trim(attname))
                 case('long_name')
                    mpp_file(unit)%Axis(dimid)%longname=mpp_file(unit)%Axis(dimid)%Att(j)%catt(1:len)
                 case('units')
                    mpp_file(unit)%Axis(dimid)%units=mpp_file(unit)%Axis(dimid)%Att(j)%catt(1:len)
                 case('cartesian_axis')
                    mpp_file(unit)%Axis(dimid)%cartesian=mpp_file(unit)%Axis(dimid)%Att(j)%catt(1:len)
                 case('calendar')
                    mpp_file(unit)%Axis(dimid)%calendar=mpp_file(unit)%Axis(dimid)%Att(j)%catt(1:len)
                    mpp_file(unit)%Axis(dimid)%calendar = lowercase(cut0(mpp_file(unit)%Axis(dimid)%calendar))
                    if (trim(mpp_file(unit)%Axis(dimid)%calendar) == 'none') &
                         mpp_file(unit)%Axis(dimid)%calendar = 'no_calendar'
                    if (trim(mpp_file(unit)%Axis(dimid)%calendar) == 'no_leap') &
                         mpp_file(unit)%Axis(dimid)%calendar = 'noleap'
                    if (trim(mpp_file(unit)%Axis(dimid)%calendar) == '365_days') &
                         mpp_file(unit)%Axis(dimid)%calendar = '365_day'
                    if (trim(mpp_file(unit)%Axis(dimid)%calendar) == '360_days') &
                         mpp_file(unit)%Axis(dimid)%calendar = '360_day'
                 case('calendar_type')
                    mpp_file(unit)%Axis(dimid)%calendar=mpp_file(unit)%Axis(dimid)%Att(j)%catt(1:len)
                    mpp_file(unit)%Axis(dimid)%calendar = lowercase(cut0(mpp_file(unit)%Axis(dimid)%calendar))
                    if (trim(mpp_file(unit)%Axis(dimid)%calendar) == 'none') &
                         mpp_file(unit)%Axis(dimid)%calendar = 'no_calendar'
                    if (trim(mpp_file(unit)%Axis(dimid)%calendar) == 'no_leap') &
                         mpp_file(unit)%Axis(dimid)%calendar = 'noleap'
                    if (trim(mpp_file(unit)%Axis(dimid)%calendar) == '365_days') &
                         mpp_file(unit)%Axis(dimid)%calendar = '365_day'
                    if (trim(mpp_file(unit)%Axis(dimid)%calendar) == '360_days') &
                         mpp_file(unit)%Axis(dimid)%calendar = '360_day'
                 case('compress')
                    mpp_file(unit)%Axis(dimid)%compressed=mpp_file(unit)%Axis(dimid)%Att(j)%catt(1:len)
                 case('positive')
                    attval = mpp_file(unit)%Axis(dimid)%Att(j)%catt(1:len)
                    if( attval.eq.'down' )then
                       mpp_file(unit)%Axis(dimid)%sense=-1
                    else if( attval.eq.'up' )then
                       mpp_file(unit)%Axis(dimid)%sense=1
                    endif
                 end select

              enddo
           endif
        enddo

        ! assign axis bounds
        do j = 1, mpp_file(unit)%ndim
           if(.not. associated(mpp_file(unit)%Axis(j)%data)) cycle
           len = size(mpp_file(unit)%Axis(j)%data(:))
           allocate(mpp_file(unit)%Axis(j)%data_bounds(len+1))
           mpp_file(unit)%Axis(j)%name_bounds = 'none'
           bounds_name = 'none'
           found_bounds = .false.
           do i = 1, mpp_file(unit)%Axis(j)%natt
              if(trim(mpp_file(unit)%Axis(j)%Att(i)%name) == 'bounds' .OR. &
                   trim(mpp_file(unit)%Axis(j)%Att(i)%name) == 'edges' ) then
                 bounds_name = mpp_file(unit)%Axis(j)%Att(i)%catt
                 found_bounds = .true.
                 exit
              endif
           enddo
           !-- loop through all the fields to locate bounds_name
           if( found_bounds ) then
              found_bounds = .false.
              do i = 1, mpp_file(unit)%ndim
                 if(.not. associated(mpp_file(unit)%Axis(i)%data)) cycle
                 if(trim(mpp_file(unit)%Axis(i)%name) == trim(bounds_name)) then
                    found_bounds = .true.
                    if(size(mpp_file(unit)%Axis(i)%data(:)) .NE. len+1) &
                         call mpp_error(FATAL, "mpp_read_meta: improperly size bounds for field "// &
                            trim(bounds_name)//" in file "// trim(mpp_file(unit)%name) )
                    mpp_file(unit)%Axis(j)%data_bounds(:) = mpp_file(unit)%Axis(i)%data(:)
                    exit
                 endif
              enddo
              if( .not. found_bounds ) then
                 do i=1, nvar_total
                    error=NF_INQ_VAR(ncid,i,name,type,nvdims,dimids,nvatts);call netcdf_err( error, mpp_file(unit) )
                    if(trim(name) == trim(bounds_name)) then
                       found_bounds = .true.
                       if(nvdims .NE. 2) &
                            call mpp_error(FATAL, "mpp_read_meta: field "//trim(bounds_name)//" in file "//&
                              trim(mpp_file(unit)%name)//" must be 2-D field")
                       if(mpp_file(unit)%Axis(dimids(1))%len .NE. 2) &
                            call mpp_error(FATAL, "mpp_read_meta: first dimension size of field "// &
                            trim(mpp_file(unit)%Var(i)%name)//" from file "//trim(mpp_file(unit)%name)// &
                            " must be 2")
                       if(mpp_file(unit)%Axis(dimids(2))%len .NE. len) &
                            call mpp_error(FATAL, "mpp_read_meta: second dimension size of field "// &
                            trim(mpp_file(unit)%Var(i)%name)//" from file "//trim(mpp_file(unit)%name)// &
                            " is not correct")
                       select case (type)
                       case (NF_INT)
                          allocate(ivals(2*len), STAT=istat)
                          if ( istat .ne. 0 ) then
                             write(text,'(A)') istat
                             call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate array ivals."//&
                                  "  STAT = "//trim(text))
                          end if
                          error = NF_GET_VAR_INT(ncid,i,ivals)
                          call netcdf_err( error, mpp_file(unit), string=" Field="//trim(bounds_name) )
                          mpp_file(unit)%Axis(j)%data_bounds(1:len) =ivals(1:(2*len-1):2)
                          mpp_file(unit)%Axis(j)%data_bounds(len+1) = ivals(2*len)
                          deallocate(ivals)
                       case (NF_FLOAT)
                          allocate(rvals(2*len), STAT=istat)
                          if ( istat .ne. 0 ) then
                             write(text,'(A)') istat
                             call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate array rvals. "// &
                                  "  STAT = "//trim(text))
                          end if
                          error = NF_GET_VAR_REAL(ncid,i,rvals)
                          call netcdf_err( error, mpp_file(unit), string=" Field="//trim(bounds_name) )
                          mpp_file(unit)%Axis(j)%data_bounds(1:len) =rvals(1:(2*len-1):2)
                          mpp_file(unit)%Axis(j)%data_bounds(len+1) = rvals(2*len)
                          deallocate(rvals)
                       case (NF_DOUBLE)
                          allocate(r8vals(2*len), STAT=istat)
                          if ( istat .ne. 0 ) then
                             write(text,'(A)') istat
                             call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate array r8vals. "//&
                                  "  STAT = "//trim(text))
                          end if
                          error = NF_GET_VAR_DOUBLE(ncid,i,r8vals)
                          call netcdf_err( error, mpp_file(unit), string=" Field="//trim(bounds_name) )
                          mpp_file(unit)%Axis(j)%data_bounds(1:len) =r8vals(1:(2*len-1):2)
                          mpp_file(unit)%Axis(j)%data_bounds(len+1) = r8vals(2*len)
                          deallocate(r8vals)
                       case default
                          call mpp_error( FATAL, 'mpp_io_mod(mpp_read_meta): Invalid data type for dimension' )
                       end select
                       exit
                    endif
                 enddo
              endif
           endif
           if (found_bounds) then
              mpp_file(unit)%Axis(j)%name_bounds = trim(bounds_name)
           else
              deallocate(mpp_file(unit)%Axis(j)%data_bounds)
              mpp_file(unit)%Axis(j)%data_bounds =>NULL()
           endif
        enddo

! assign variable info
        nv = 0
        do i=1, nvar_total
           error=NF_INQ_VAR(ncid,i,name,type,nvdims,dimids,nvatts);call netcdf_err( error, mpp_file(unit) )
!
! is this a dimension variable?
!
           isdim=.false.
           do j=1,ndim
              if( trim(lowercase(name)).EQ.trim(lowercase(mpp_file(unit)%Axis(j)%name)) )isdim=.true.
           enddo

           if( .not.isdim )then
! for non-dimension variables
              nv=nv+1; if( nv.GT.mpp_file(unit)%nvar )call mpp_error( FATAL, 'variable index exceeds number of defined variables' )
              mpp_file(unit)%Var(nv)%type = type
              mpp_file(unit)%Var(nv)%id = i
              mpp_file(unit)%Var(nv)%name = name
              mpp_file(unit)%Var(nv)%natt = nvatts
! determine packing attribute based on NetCDF variable type
             select case (type)
             case(NF_SHORT)
                 mpp_file(unit)%Var(nv)%pack = 4
             case(NF_FLOAT)
                 mpp_file(unit)%Var(nv)%pack = 2
             case(NF_DOUBLE)
                 mpp_file(unit)%Var(nv)%pack = 1
             case (NF_INT)
                 mpp_file(unit)%Var(nv)%pack = 2
             case (NF_CHAR)
                 mpp_file(unit)%Var(nv)%pack = 1
             case default
                   call mpp_error( FATAL, 'Invalid variable type in NetCDF file' )
             end select
! assign dimension ids
              mpp_file(unit)%Var(nv)%ndim = nvdims
              allocate(mpp_file(unit)%Var(nv)%axes(nvdims))
              do j=1,nvdims
                 mpp_file(unit)%Var(nv)%axes(j) = mpp_file(unit)%Axis(dimids(j))
              enddo
              allocate(mpp_file(unit)%Var(nv)%size(nvdims))

              do j=1,nvdims
                 if(dimids(j).eq.mpp_file(unit)%recdimid  .and. mpp_file(unit)%time_level/=-1)then
                    mpp_file(unit)%Var(nv)%time_axis_index = j   !dimids(j). z1l: Should be j. 
                                                                 !This will cause problem when appending to existed file.
                    mpp_file(unit)%Var(nv)%size(j)=1    ! dimid length set to 1 here for consistency w/ mpp_write
                 else
                    mpp_file(unit)%Var(nv)%size(j)=mpp_file(unit)%Axis(dimids(j))%len
                 endif
              enddo
! assign variable atts
              if( nvatts.GT.0 )allocate(mpp_file(unit)%Var(nv)%Att(nvatts))

              do j=1,nvatts
                 mpp_file(unit)%Var(nv)%Att(j) = default_att
              enddo

              do j=1,nvatts
                 error=NF_INQ_ATTNAME(ncid,i,j,attname);call netcdf_err( error, mpp_file(unit), field=mpp_file(unit)%Var(nv) )
                 error=NF_INQ_ATT(ncid,i,attname,type,len)
                 call netcdf_err( error, mpp_file(unit),field= mpp_file(unit)%Var(nv), string=' Attribute='//attname )
                 mpp_file(unit)%Var(nv)%Att(j)%name = trim(attname)
                 mpp_file(unit)%Var(nv)%Att(j)%type = type
                 mpp_file(unit)%Var(nv)%Att(j)%len = len

                 select case (type)
                   case (NF_CHAR)
                     if (len.gt.512) call mpp_error(FATAL,'VAR ATT too long')
                     error=NF_GET_ATT_TEXT(ncid,i,trim(attname),mpp_file(unit)%Var(nv)%Att(j)%catt(1:len))
                     call netcdf_err( error, mpp_file(unit), field=mpp_file(unit)%var(nv), attr=mpp_file(unit)%var(nv)%att(j) )
                     if (verbose .and. pe == 0 )&
                           print *, 'Var ',nv,' ATT ',trim(attname),' ',mpp_file(unit)%Var(nv)%Att(j)%catt(1:len)
! store integers as float internally
                   case (NF_SHORT)
                     allocate(mpp_file(unit)%Var(nv)%Att(j)%fatt(len), STAT=istat)
                     if ( istat .ne. 0 ) then
                        write(text,'(A)') istat
                        call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Var%Att%fatt, "//&
                             & "NF_SHORT CASE.  STAT = "//trim(text))
                     end if
                     allocate(i2vals(len), STAT=istat)
                     if ( istat .ne. 0 ) then
                        write(text,'(A)') istat
                        call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array i2vals.  STAT = "&
                             & //trim(text))
                     end if
                     error=NF_GET_ATT_INT2(ncid,i,trim(attname),i2vals)
                     call netcdf_err( error, mpp_file(unit), field=mpp_file(unit)%var(nv), attr=mpp_file(unit)%var(nv)%att(j) )
                     mpp_file(unit)%Var(nv)%Att(j)%fatt(1:len)= i2vals(1:len)
                     if( verbose  .and. pe == 0 )&
                          print *, 'Var ',nv,' ATT ',trim(attname),' ',mpp_file(unit)%Var(nv)%Att(j)%fatt
                     deallocate(i2vals)
                   case (NF_INT)
                     allocate(mpp_file(unit)%Var(nv)%Att(j)%fatt(len), STAT=istat)
                     if ( istat .ne. 0 ) then
                        write(text,'(A)') istat
                        call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Var%Att%fatt, "//&
                             & "NF_INT CASE.  STAT = "//trim(text))
                     end if
                     allocate(ivals(len), STAT=istat)
                     if ( istat .ne. 0 ) then
                        write(text,'(A)') istat
                        call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array ivals.  STAT = "&
                             & //trim(text))
                     end if
                     error=NF_GET_ATT_INT(ncid,i,trim(attname),ivals)
                     call netcdf_err( error, mpp_file(unit), field=mpp_file(unit)%var(nv), attr=mpp_file(unit)%var(nv)%att(j) )
                     mpp_file(unit)%Var(nv)%Att(j)%fatt(1:len)=ivals(1:len)
                     if( verbose .and. pe == 0  )&
                          print *, 'Var ',nv,' ATT ',trim(attname),' ',mpp_file(unit)%Var(nv)%Att(j)%fatt
                     deallocate(ivals)
                   case (NF_FLOAT)
                     allocate(mpp_file(unit)%Var(nv)%Att(j)%fatt(len), STAT=istat)
                     if ( istat .ne. 0 ) then
                        write(text,'(A)') istat
                        call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Var%Att%fatt, "//&
                             & "NF_FLOAT CASE.  STAT = "//trim(text))
                     end if
                     allocate(rvals(len), STAT=istat)
                     if ( istat .ne. 0 ) then
                        write(text,'(A)') istat
                        call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array rvals.  STAT = "&
                             & //trim(text))
                     end if
                     error=NF_GET_ATT_REAL(ncid,i,trim(attname),rvals)
                     call netcdf_err( error, mpp_file(unit), field=mpp_file(unit)%var(nv), attr=mpp_file(unit)%var(nv)%att(j) )
                     mpp_file(unit)%Var(nv)%Att(j)%fatt(1:len)=rvals(1:len)
                     if( verbose  .and. pe == 0 )&
                          print *, 'Var ',nv,' ATT ',trim(attname),' ',mpp_file(unit)%Var(nv)%Att(j)%fatt
                     deallocate(rvals)
                   case (NF_DOUBLE)
                     allocate(mpp_file(unit)%Var(nv)%Att(j)%fatt(len), STAT=istat)
                     if ( istat .ne. 0 ) then
                        write(text,'(A)') istat
                        call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate mpp_file%Var%Att%fatt, "//&
                             & "NF_DOUBLE CASE.  STAT = "//trim(text))
                     end if
                     allocate(r8vals(len), STAT=istat)
                     if ( istat .ne. 0 ) then
                        write(text,'(A)') istat
                        call mpp_error(FATAL, "mpp_io_mod(mpp_read_meta): Unable to allocate temporary array r8vals.  STAT = "&
                             & //trim(text))
                     end if
                      error=NF_GET_ATT_DOUBLE(ncid,i,trim(attname),r8vals)
                     call netcdf_err( error, mpp_file(unit), field=mpp_file(unit)%var(nv), attr=mpp_file(unit)%var(nv)%att(j) )
                     mpp_file(unit)%Var(nv)%Att(j)%fatt(1:len)=r8vals(1:len)
                     if( verbose .and. pe == 0  ) &
                          print *, 'Var ',nv,' ATT ',trim(attname),' ',mpp_file(unit)%Var(nv)%Att(j)%fatt
                     deallocate(r8vals)
                   case default
                        call mpp_error( FATAL, 'Invalid data type for variable att' )
                 end select
! assign pre-defined field attributes
                 select case (trim(attname))
                    case ('long_name')
                      mpp_file(unit)%Var(nv)%longname=mpp_file(unit)%Var(nv)%Att(j)%catt(1:len)
                    case('units')
                      mpp_file(unit)%Var(nv)%units=mpp_file(unit)%Var(nv)%Att(j)%catt(1:len)
                    case('scale_factor')
                       mpp_file(unit)%Var(nv)%scale=mpp_file(unit)%Var(nv)%Att(j)%fatt(1)
                    case('missing')
                       mpp_file(unit)%Var(nv)%missing=mpp_file(unit)%Var(nv)%Att(j)%fatt(1)
                    case('missing_value')
                       mpp_file(unit)%Var(nv)%missing=mpp_file(unit)%Var(nv)%Att(j)%fatt(1)
                    case('_FillValue')
                       mpp_file(unit)%Var(nv)%fill=mpp_file(unit)%Var(nv)%Att(j)%fatt(1)
                    case('add_offset')
                       mpp_file(unit)%Var(nv)%add=mpp_file(unit)%Var(nv)%Att(j)%fatt(1)
                    case('packing')
                       mpp_file(unit)%Var(nv)%pack=mpp_file(unit)%Var(nv)%Att(j)%fatt(1)
                    case('valid_range')
                       mpp_file(unit)%Var(nv)%min=mpp_file(unit)%Var(nv)%Att(j)%fatt(1)
                       mpp_file(unit)%Var(nv)%max=mpp_file(unit)%Var(nv)%Att(j)%fatt(2)
                    case('checksum')
                         checksum_char = mpp_file(unit)%Var(nv)%Att(j)%catt
!                        Scan  checksum attribute for , delimiter. If found implies multiple time levels.
                         checksumf = 0
                         num_checksumf = 1

                         last = len_trim(checksum_char)
                         is = index (trim(checksum_char),",") ! A value of 0 implies only 1 checksum value
                         do while ((is > 0) .and. (is < (last-15)))
                           is = is + scan(checksum_char(is:last), "," ) ! move starting pointer after ","
                           num_checksumf = num_checksumf + 1
                         enddo
                         is =1
                         do k = 1, num_checksumf
                           read (checksum_char(is:is+15),'(Z16)') checksumf
                           mpp_file(unit)%Var(nv)%checksum(k) = checksumf
                           is = is + 17 ! Move index past the ,
                         enddo
                 end select
              enddo
           endif
        enddo   ! end variable loop
      else
        call mpp_error( FATAL,  'MPP READ CURRENTLY DOES NOT SUPPORT NON-NETCDF' )
      endif

      mpp_file(unit)%initialized = .TRUE.
#else
      call mpp_error( FATAL, 'MPP_READ currently requires use_netCDF option' )
#endif
      return
    end subroutine mpp_read_meta


    function cut0(string)
      character(len=256) :: cut0
      character(len=*), intent(in) :: string
      integer :: i

      cut0 = string
      i = index(string,achar(0))
      if(i > 0) cut0(i:i) = ' '

      return
    end function cut0


    subroutine mpp_get_tavg_info(unit, field, fields, tstamp, tstart, tend, tavg)
      implicit none
      integer, intent(in) :: unit
      type(fieldtype), intent(in) :: field
      type(fieldtype), intent(in), dimension(:) :: fields
      real, intent(inout), dimension(:) :: tstamp, tstart, tend, tavg
!balaji: added because mpp_read can only read default reals
!      when running with -r4 this will read a default real and then cast double
      real :: t_default_real


      integer :: n, m
      logical :: tavg_info_exists

      tavg = -1.0


      if (size(tstamp,1) /= size(tstart,1)) call mpp_error(FATAL,&
            'size mismatch in mpp_get_tavg_info')

      if ((size(tstart,1) /= size(tend,1)) .OR. (size(tstart,1) /= size(tavg,1))) then
          call mpp_error(FATAL,'size mismatch in mpp_get_tavg_info')
      endif

      tstart = tstamp
      tend = tstamp

      tavg_info_exists = .false.

#ifdef use_netCDF
      do n= 1, field%natt
         if (field%Att(n)%type .EQ. NF_CHAR) then
             if (field%Att(n)%name(1:13) == 'time_avg_info') then
                 tavg_info_exists = .true.
                 exit
             endif
         endif
      enddo
#endif
      if (tavg_info_exists) then
          do n = 1, size(fields(:))
             if (trim(fields(n)%name) == 'average_T1') then
                 do m = 1, size(tstart(:))
                    call mpp_read(unit, fields(n),t_default_real, m)
                    tstart(m) = t_default_real
                 enddo
             endif
             if (trim(fields(n)%name) == 'average_T2') then
                 do m = 1, size(tend(:))
                    call mpp_read(unit, fields(n),t_default_real, m)
                    tend(m) = t_default_real
                 enddo
             endif
             if (trim(fields(n)%name) == 'average_DT') then
                 do m = 1, size(tavg(:))
                    call mpp_read(unit, fields(n),t_default_real, m)
                    tavg(m) = t_default_real
                 enddo
             endif
          enddo

      end if
      return
    end subroutine mpp_get_tavg_info

!#######################################################################