!= NetCDF のラッパープログラム
!
!= NetCDF wrapper
!
! Authors::   Yoshiyuki O. Takahashi
! Version::   $Id: netcdf_wrapper.F90,v 1.4 2015/01/29 12:02:28 yot Exp $
! Tag Name::  $Name:  $
! Copyright:: Copyright (C) GFD Dennou Club, 2008. All rights reserved.
! License::   See COPYRIGHT[link:../../../COPYRIGHT]
!

module netcdf_wrapper
  !
  != NetCDF のラッパープログラム
  !
  != NetCDF wrapper
  !
  ! <b>Note that Japanese and English are described in parallel.</b>
  !
  ! gtool で扱わない部分をカバーする NetCDF のラッパープログラム. 
  !
  ! The wrapper routines in this module treat NetCDF input/output that are not 
  ! covered by gtool.
  !
  !== Procedures List
  !
  ! NWInqDimLen :: 軸の長さを問い合わせる
  ! ----------- :: ------------
  ! NWInqDimLen :: Inquire length of a dimension
  !
  !== NAMELIST
  !
  ! NAMELIST#radiation_DennouAGCM_nml
  !

  ! モジュール引用 ; USE statements
  !

  ! 種別型パラメタ
  ! Kind type parameter
  !
  use dc_types, only: DP, &      ! 倍精度実数型. Double precision.
    &                 STRING     ! 文字列.       Strings.

  ! メッセージ出力
  ! Message output
  !
  use dc_message, only: MessageNotify

  ! NetCDF
  !
  use netcdf

  implicit none

  private


  public :: NWChkDim
  public :: NWInqDimLen
  public :: NWGetAtt
  public :: NWPresentAVarInFile


  interface NWGetAtt
    module procedure NWGetAttChar, NWGetAttInteger
  end interface


  character(*), parameter:: module_name = 'netcdf_wrapper'
                              ! モジュールの名称. 
                              ! Module name
  character(*), parameter:: version = &
    & '$Name:  $' // &
    & '$Id: netcdf_wrapper.F90,v 1.4 2015/01/29 12:02:28 yot Exp $'
                              ! モジュールのバージョン
                              ! Module version

  !--------------------------------------------------------------------------

contains

  !--------------------------------------------------------------------------

  function NWChkDim( &
    & ncfn,    & ! (in )
    & dimname  & ! (in )
    & ) result( FlagExist )
    !
    ! 軸の長さを問い合わせます. 
    !
    ! Inquire length of a dimension. 
    !

    character(*), intent(in ) :: ncfn
                                  ! NetCDF filename
    character(*), intent(in ) :: dimname
                                  ! Dimension name
    logical :: FlagExist

    !
    ! Local variables
    !
    character(STRING) :: mode
    integer           :: ncid
    integer           :: dimid
    integer           :: status
    character(STRING) :: err_mes


    err_mes = "In NWChkDim: " // trim( dimname ) // " in " // trim( ncfn )

    mode = "read"
    call NWOpen( NWMkNCFN(ncfn), mode, ncid )

    status = NF90_inq_dimid( ncid, dimname, dimid )

    if ( status == NF90_EBADDIM ) then
      FlagExist = .false.
    else
      FlagExist = .true.
    end if

    call NWClose( ncid )


  end function NWChkDim

  !--------------------------------------------------------------------------

  subroutine NWInqDimLen( &
    & ncfn,    & ! (in )
    & dimname, & ! (in )
    & len      & ! (out)
    & )
    !
    ! 軸の長さを問い合わせます. 
    !
    ! Inquire length of a dimension. 
    !

    character(*), intent(in ) :: ncfn
                                  ! NetCDF filename
    character(*), intent(in ) :: dimname
                                  ! Dimension name
    integer     , intent(out) :: len
                                  ! Length of a dimension


    !
    ! Local variables
    !
    character(STRING) :: mode
    integer           :: ncid
    integer           :: dimid
    integer           :: status
    character(STRING) :: err_mes


    err_mes = "In NWInqDimLen: " // trim( dimname ) // " in " // trim( ncfn )

    mode = "read"
    call NWOpen( NWMkNCFN(ncfn), mode, ncid )

    status = NF90_inq_dimid( ncid, dimname, dimid )
    call NWHandleErr( status, err_mes )

    status = NF90_Inquire_Dimension( ncid, dimid, len = len )
    call NWHandleErr( status, err_mes )

    call NWClose( ncid )


  end subroutine NWInqDimLen

  !--------------------------------------------------------------------------

  subroutine NWGetAttChar( ncfn, varname, attname, att )

    character(*), intent(in ) :: ncfn
                                  ! NetCDF filename
    character(*), intent(in ) :: varname
                                  ! Variable name
    character(*), intent(in ) :: attname
                                  ! Attribute name
    character(*), intent(out) :: att
                                  ! Attribute


    ! Local variables
    !
    character(STRING) :: mode
    integer           :: ncid
    integer           :: varid
    integer           :: status
    character(STRING) :: err_mes


    err_mes = "In NWGetAttChar"
    err_mes = trim( err_mes ) &
      // ', Var: ' // trim( varname ) // ', Att: ' // trim( attname )

    mode = "read"
    call NWOpen( NWMkNCFN(ncfn), mode, ncid )

    if ( ( varname == 'global' ) .or. ( varname == 'GLOBAL' ) ) then
      varid = NF90_GLOBAL
    else
      status = nf90_inq_varid( ncid, varname, varid )
      call NWHandleErr( status, err_mes )
    end if

    status = nf90_get_att( ncid, varid, attname, att )
    call NWHandleErr( status, err_mes )

    call NWClose( ncid )

  end subroutine NWGetAttChar

  !--------------------------------------------------------------------------

  subroutine NWGetAttInteger( ncfn, varname, attname, att )

    character(*), intent(in ) :: ncfn
                                  ! NetCDF filename
    character(*), intent(in ) :: varname
                                  ! Variable name
    character(*), intent(in ) :: attname
                                  ! Attribute name
    integer     , intent(out) :: att
                                  ! Attribute


    ! Local variables
    !
    character(STRING) :: mode
    integer           :: ncid
    integer           :: varid
    integer           :: status
    character(STRING) :: err_mes


    err_mes = "In NWGetAttInteger"
    err_mes = trim( err_mes ) &
      // ', Var: ' // trim( varname ) // ', Att: ' // trim( attname )

    mode = "read"
    call NWOpen( NWMkNCFN(ncfn), mode, ncid )

    if ( ( varname == 'global' ) .or. ( varname == 'GLOBAL' ) ) then
      varid = NF90_GLOBAL
    else
      status = nf90_inq_varid( ncid, varname, varid )
      call NWHandleErr( status, err_mes )
    end if

    status = nf90_get_att( ncid, varid, attname, att )
    call NWHandleErr( status, err_mes )

    call NWClose( ncid )

  end subroutine NWGetAttInteger

  !--------------------------------------------------------------------------

  function NWPresentAVarInFile( ncfn, varname ) result( Flag )

    character(*), intent(in ) :: ncfn
                                  ! NetCDF filename
    character(*), intent(in ) :: varname
                                  ! Variable name

    logical                   :: Flag

    ! Local variables
    !
    character(STRING) :: mode
    integer           :: ncid
    integer           :: varid
    integer           :: status
    character(STRING) :: err_mes


    err_mes = "In NWChkAVarInFile"
    err_mes = trim( err_mes ) // ', Var: ' // trim( varname )

    mode = "read"
    call NWOpen( NWMkNCFN(ncfn), mode, ncid )

    status = nf90_inq_varid( ncid, varname, varid )
    if( status == nf90_noerr ) then
      Flag = .true.
    else
      Flag = .false.
      call MessageNotify( 'M', module_name, &
        & 'Variable, %c, cannot be found in file, %c.', &
        & c1 = trim(varname), c2 = trim(NWMkNCFN(ncfn)) &
        & )
    end if

    call NWClose( ncid )

  end function NWPresentAVarInFile

  !--------------------------------------------------------------------------
  ! Routines below are internal ones.
  !--------------------------------------------------------------------------

  character(STRING) function NWMkNCFN( ncfn_in ) result( ncfn_out )

    ! MPI wrapper
    !
#ifdef LIB_MPI
    use mpi_wrapper, only : myrank
#endif

    character(*), intent(in ) :: ncfn_in
                                  ! NetCDF filename


#ifdef LIB_MPI
    if ( ncfn_in(len_trim(ncfn_in)-2:len_trim(ncfn_in)) /= '.nc' ) then
      write( 6, * ) trim(ncfn_in), ' is inappropriate.'
      stop
    end if
    write( ncfn_out, '(a,a,i6.6,a)' ) ncfn_in(1:len_trim(ncfn_in)-3), '_rank', myrank, '.nc'
#else
    ncfn_out = ncfn_in
#endif

  end function NWMkNCFN

  !--------------------------------------------------------------------------

  subroutine NWHandleErr( &
    & status, &
    & err_mes &
    & )

    integer     , intent(in)           :: status
    character(*), intent(in), optional :: err_mes


    if( status .ne. nf90_noerr ) then
      if( present( err_mes ) ) print *, trim( err_mes )
      print *, trim( nf90_strerror( status ) )
      stop "STOP"
    end if


  end subroutine NWHandleErr

  !--------------------------------------------------------------------------

  subroutine NWOpen( &
    & path, &
    & mode, &
    & ncid &
    & )

    character(*), intent(in ) :: path
    character(*), intent(in ) :: mode
    integer     , intent(out) :: ncid


    ! Local variables
    !
    integer           :: cmode, omode
    integer           :: status
    character(STRING) :: err_mes


    err_mes = "In NWOpen"
    err_mes = trim( err_mes ) // ', path: ' // trim( path ) &
      // ', mode: ' // trim( mode )


    if( ( mode .eq. "new" ) .or. ( mode .eq. "NEW" ) ) then

      cmode = NF90_CLOBBER
      status = nf90_create( path, cmode, ncid )
      call NWHandleErr( status, err_mes )

    else

      if( ( mode .eq. "read" ) .or. ( mode .eq. "READ" ) ) then
        omode = NF90_NOWRITE
      else if( ( mode .eq. "write" ) .or. ( mode .eq. "WRITE" ) ) then
        omode = NF90_WRITE
      else
        write( 6, * ) "Inproper argument, mode"
        write( 6, * ) "Argument mode should be either of words below."
        write( 6, * ) " 'new', 'NEW', 'read', 'READ', 'write', or 'WRITE'"
        stop
      end if

      status = nf90_open( path, omode, ncid )
      call NWHandleErr( status, err_mes )

    end if


  end subroutine NWOpen

  !--------------------------------------------------------------------------

  subroutine NWClose( &
    & ncid &
    & )

    integer, intent(in) :: ncid


    ! Local variables
    !
    integer           :: status
    character(STRING) :: err_mes


    err_mes = "In NWClose"


    status = nf90_close( ncid )
    call NWHandleErr( status, err_mes )

  end subroutine NWClose

  !--------------------------------------------------------------------------

  subroutine NWReDef( &
    & ncid &
    & )

    integer, intent(in) :: ncid


    ! Local variables
    !
    integer           :: status
    character(STRING) :: err_mes


    err_mes = "In NWReDef"


    status = nf90_redef( ncid )
    if( ( status .ne. nf90_noerr ) .and. ( status .ne. nf90_eindefine ) ) &
      call NWHandleErr( status, err_mes )


  end subroutine NWReDef

  !--------------------------------------------------------------------------

  subroutine NWEndDef( &
    & ncid &
    & )

    integer, intent(in) :: ncid


    ! Local variables
    !
    integer           :: status
    character(STRING) :: err_mes


    err_mes = "In NWEndDef"


    status = nf90_enddef( ncid )
    if( ( status .ne. nf90_noerr ) .and. ( status .ne. nf90_enotindefine ) ) &
      call NWHandleErr( status, err_mes )


  end subroutine NWEndDef

  !--------------------------------------------------------------------------

end module netcdf_wrapper
