module opt_prop

  use dc_types
  use gtool_history

  implicit none

  private

  public :: OptPropInit
  public :: OptPropFinalize
  public :: OptPropGetNWN
  public :: OptPropGetNPress
  public :: OptPropGetNMol
  public :: OptPropGetNBand
  public :: OptPropGetDelWN
  public :: OptPropGetWN
  public :: OptPropGetBandWNBnds
  public :: OptPropGetBandBinIndexBnds
  public :: OptPropGetAbsCoefProf
  public :: OptPropGetPFInted
  public :: OptPropGetPFIntedRatio


  integer, save               :: IDMethod
  integer, parameter          :: IDMethodLBL        =  1
  integer, parameter          :: IDMethodLBLConstPF =  2
  integer, parameter          :: IDMethodKD         = 10


  character(128)       , save :: OptPropNcFNSave
  integer              , save :: NPressSave
  integer              , save :: NMolSave

  ! LBL
  integer              , save :: NWaveNumSave
  real(DP), allocatable, save :: w_WaveNumSave(:)

  ! KD
  integer              , save :: NBinSave
  integer , allocatable, save :: a_BinNumSave(:)
  integer , parameter         :: NBndNumSave = 2
  integer              , save :: NBandSave
  real(DP), allocatable, save :: aa_BandWaveNumBndsSave     (:,:)
!!$  integer , allocatable, save :: aa_BandWaveNumIndexBndsSave(:,:)
  integer , allocatable, save :: aa_BandBinIndexBndsSave    (:,:)

  character(128), save :: WaveNumAxisName
  character(128), save :: AbsCoefVarName
  character(128), save :: PFRatioVarName


  logical              , save :: FlagHaveAllOptProp = .true.
  real(DP), allocatable, save :: za_PFRatio  (:,:)
  real(DP), allocatable, save :: zma_AbsCoef (:,:,:)
  real(DP), allocatable, save :: a_DelWaveNum(:)


  character(128), save :: ModuleName = 'opt_prop'

  logical, save :: FlagInited = .false.

contains

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

  subroutine OptPropInit( &
    & Type, &
    & OptPropNcFN &
    & )

    character(*), intent(in) :: Type
    character(*), intent(in) :: OptPropNcFN

    if ( FlagInited ) return


    select case ( Type )
    case ( 'LBL' )
      IDMethod = IDMethodLBL
    case ( 'LBLConstPF' )
      IDMethod = IDMethodLBLConstPF
    case ( 'KD' )
      IDMethod = IDMethodKD
    case default
      write( 6, * ) 'In module, ', trim( ModuleName )
      write( 6, * ) 'Unexpected type: ', trim( Type )
      stop
    end select


    call OptPropInitCore( &
      & OptPropNcFN &
      & )


    FlagInited = .true.

  end subroutine OptPropInit

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

  subroutine OptPropFinalize


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    call OptPropFinalizeCore


    FlagInited = .false.

  end subroutine OptPropFinalize

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

  function OptPropGetNWN() result( NWN )


    integer :: NWN

    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    NWN = OptPropGetNWNCore()


  end function OptPropGetNWN

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

  function OptPropGetNPress() result( NPress )


    integer :: NPress

    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    NPress = OptPropGetNPressCore()


  end function OptPropGetNPress

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

  function OptPropGetNMol() result( NMol )


    integer :: NMol

    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    NMol = OptPropGetNMolCore()


  end function OptPropGetNMol

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

  function OptPropGetNBand() result( NBand )


    integer :: NBand

    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    NBand = OptPropGetNBandCore()


  end function OptPropGetNBand

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

  function OptPropGetDelWN( iWaveNum ) result( DelWaveNum )


    integer, intent(in) :: iWaveNum


    real(DP) :: DelWaveNum

    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    DelWaveNum = OptPropGetDelWNCore( iWaveNum )


  end function OptPropGetDelWN

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

  function OptPropGetWN( iWaveNum ) result( WaveNum )


    integer, intent(in) :: iWaveNum

    real(DP) :: WaveNum

    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    WaveNum = OptPropGetWNCore( iWaveNum )


  end function OptPropGetWN

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

  function OptPropGetBandWNBnds() result( aa_BandWNBnds )


    real(DP) :: aa_BandWNBnds(2,NBandSave)


    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    aa_BandWNBnds = OptPropGetBandWNBndsCore()


  end function OptPropGetBandWNBnds

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

  function OptPropGetBandBinIndexBnds() result( aa_BandBinIndexBnds )


    integer :: aa_BandBinIndexBnds(2,NBandSave)


    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    aa_BandBinIndexBnds = OptPropGetBandBinIndexBndsCore()


  end function OptPropGetBandBinIndexBnds

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

  subroutine OptPropGetAbsCoefProf( &
    & iWaveNum, NPress, NMol, &
    & zm_AbsCoef &
    & )


    integer     , intent(in ) :: iWaveNum
    integer     , intent(in ) :: NPress
    integer     , intent(in ) :: NMol
    real(DP)    , intent(out) :: zm_AbsCoef(1:NPress,1:NMol)


    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    call OptPropGetAbsCoefProfCore( &
      & iWaveNum, NPress, NMol, &
      & zm_AbsCoef &
      & )


  end subroutine OptPropGetAbsCoefProf

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

  subroutine OptPropGetPFInted( &
    & iWaveNum, &
    & kmax, &
    & r_Temp, SurfTemp, &
    & r_PFInted, SurfPFInted, SurfDPFDTInted &
    & )

    integer , intent(in ) :: iWaveNum
    integer , intent(in ) :: kmax
    real(DP), intent(in ) :: r_Temp(0:kmax)
    real(DP), intent(in ) :: SurfTemp
    real(DP), intent(out) :: r_PFInted  (0:kmax)
    real(DP), intent(out) :: SurfPFInted
    real(DP), intent(out) :: SurfDPFDTInted


    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    call OptPropGetPFIntedCore( &
      & iWaveNum, &
      & kmax, &
      & r_Temp, SurfTemp, &
      & r_PFInted, SurfPFInted, SurfDPFDTInted &
      & )


  end subroutine OptPropGetPFInted

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

  subroutine OptPropGetPFIntedRatio( &
    & iWaveNum, NPress, &
    & z_PFRatio, SurfPFRatio &
    & )

    integer     , intent(in ) :: iWaveNum
    integer     , intent(in ) :: NPress
    real(DP)    , intent(out) :: z_PFRatio(1:NPress)
    real(DP)    , intent(out) :: SurfPFRatio


    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    select case ( IDMethod )
    case ( IDMethodLBL )
      write( 6, * ) 'In ', trim( ModuleName )
      write( 6, * ) '  function OptPropGetPFIntedRatio cannot be used when IDMethod is IDMethodLBL'
      stop
    end select


    call OptPropGetPFIntedRatioCore( &
      & iWaveNum, NPress, &
      & z_PFRatio, SurfPFRatio &
      & )


  end subroutine OptPropGetPFIntedRatio

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

  subroutine OptPropInitCore( &
    & OptPropNcFN &
    & )

    use rad_utils, only : InqDimLength


    character(*), intent(in) :: OptPropNcFN


    ! local variables
    !
    integer :: NWaveNumBin


!!$    if ( .not. FlagInited ) then
!!$      write( 6, * ) trim( ModuleName ), " is not initialized."
!!$      stop
!!$    end if


    OptPropNcFNSave = OptPropNcFN

    NPressSave = InqDimLength( OptPropNcFNSave, 'r_Press'   )
    NMolSave   = InqDimLength( OptPropNcFNSave, 'm_MolNum'  )

    select case ( IDMethod )
    case ( IDMethodLBL, IDMethodLBLConstPF )
      WaveNumAxisName = 'w_WaveNum'
      AbsCoefVarName  = 'rmw_TAbsorpK'
      PFRatioVarName  = 'rw_PFRatio'
      !
      NWaveNumSave = InqDimLength( OptPropNcFNSave, WaveNumAxisName )
      allocate( w_WaveNumSave( 1:NWaveNumSave ) )
      call HistoryGet( OptPropNcFNSave, WaveNumAxisName, w_WaveNumSave, quiet = .true. )
      !
      NWaveNumBin = NWaveNumSave
    case ( IDMethodKD )
      WaveNumAxisName = 'BinNum'
      AbsCoefVarName  = 'rma_TAbsorpK'
      PFRatioVarName  = 'ra_PFRatio'
      !
      NBinSave    = InqDimLength( OptPropNcFNSave, WaveNumAxisName )
      NBandSave   = InqDimLength( OptPropNcFNSave, 'BandNum' )
      allocate( a_BinNumSave               (NBinSave) )
      allocate( aa_BandWaveNumBndsSave     (NBndNumSave,NBandSave) )
!!$      allocate( aa_BandWaveNumIndexBndsSave(NBndNumSave,NBandSave) )
      allocate( aa_BandBinIndexBndsSave    (NBndNumSave,NBandSave) )

      call HistoryGet( OptPropNcFNSave, WaveNumAxisName       , a_BinNumSave, quiet = .true. )
      call HistoryGet( OptPropNcFNSave, 'BandWaveNumBnds'     , aa_BandWaveNumBndsSave, quiet = .true. )
!!$      call HistoryGet( OptPropNcFNSave, 'BandWaveNumIndexBnds', aa_BandWaveNumIndexBndsSave, quiet = .true. )
      call HistoryGet( OptPropNcFNSave, 'BandBinIndexBnds'    , aa_BandBinIndexBndsSave, quiet = .true. )
      !
      NWaveNumBin = NBinSave
    end select

    if ( FlagHaveAllOptProp ) then
      select case ( IDMethod )
      case ( IDMethodKD )
        allocate( a_DelWaveNum(NWaveNumBin)                     )
        call HistoryGet( OptPropNcFNSave, "BinWeight", a_DelWaveNum, &
          & quiet = .true. )
      end select
      !
      allocate( zma_AbsCoef (NPressSave,NMolSave,NWaveNumBin) )
      call HistoryGet( OptPropNcFNSave, AbsCoefVarName, zma_AbsCoef, &
        & quiet = .true. )
      !
      select case ( IDMethod )
      case ( IDMethodLBL )
      case default
        allocate( za_PFRatio  (NPressSave,NWaveNumBin)          )
        call HistoryGet( OptPropNcFNSave, PFRatioVarName, za_PFRatio, &
          & quiet = .true. )
      end select
    end if


!!$    allocate( pmw_LBLAbsCoef( 1:NPressSave, 1:NMolSave, 1:NWaveNumSave ) )
!!$    call HistoryGet( OptPropNcFNSave, AbsCoefVarName, pmw_LBLAbsCoef, quiet = .true. )

  end subroutine OptPropInitCore

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

  subroutine OptPropFinalizeCore


    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if

    NPressSave = 0
    NMolSave   = 0

    select case ( IDMethod )
    case ( IDMethodLBL, IDMethodLBLConstPF )
      NWaveNumSave    = 0
      deallocate( w_WaveNumSave )
    case ( IDMethodKD )
      NBinSave    = 0
      deallocate( a_BinNumSave                )
      deallocate( aa_BandWaveNumBndsSave      )
!!$      deallocate( aa_BandWaveNumIndexBndsSave )
      deallocate( aa_BandBinIndexBndsSave     )
    end select
    !
    if ( FlagHaveAllOptProp ) then
      select case ( IDMethod )
      case ( IDMethodKD )
        deallocate( a_DelWaveNum )
      end select
      !
      deallocate( zma_AbsCoef  )
      !
      select case ( IDMethod )
      case ( IDMethodLBL )
      case default
        deallocate( za_PFRatio   )
      end select
    end if


  end subroutine OptPropFinalizeCore

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

  function OptPropGetNWNCore() result( NWN )


    integer :: NWN

    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    select case ( IDMethod )
    case ( IDMethodLBL, IDMethodLBLConstPF )
      NWN = NWaveNumSave
    case ( IDMethodKD )
      NWN = NBinSave
    case default
      write( 6, * ) 'In ', trim( ModuleName )
      write( 6, * ) '  Unexpected case default'
      stop
    end select

  end function OptPropGetNWNCore

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

  function OptPropGetNPressCore() result( NPress )


    integer :: NPress

    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    NPress = NPressSave


  end function OptPropGetNPressCore

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

  function OptPropGetNMolCore() result( NMol )


    integer :: NMol

    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    NMol = NMolSave


  end function OptPropGetNMolCore

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

  function OptPropGetNBandCore() result( NBand )


    integer :: NBand

    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    select case ( IDMethod )
    case ( IDMethodLBL, IDMethodLBLConstPF )
      write( 6, * ) 'In ', trim( ModuleName )
      write( 6, * ) '  function OptPropGetNBandCore cannot be used when IDMethod is IDMethodLBL or IDMethodLBLConstPF'
      stop
    case ( IDMethodKD )
      NBand = NBandSave
    case default
      write( 6, * ) 'In ', trim( ModuleName )
      write( 6, * ) '  Unexpected case default'
      stop
    end select


  end function OptPropGetNBandCore

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

  function OptPropGetDelWNCore( iWaveNum ) result( DelWaveNum )

    use rad_utils, only : MkRangeString

    integer, intent(in) :: iWaveNum

    real(DP) :: DelWaveNum

    ! local variables
    !
    character(128) :: String


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    select case ( IDMethod )
    case ( IDMethodLBL, IDMethodLBLConstPF )
      DelWaveNum = w_WaveNumSave( 2 ) - w_WaveNumSave( 1 )
    case ( IDMethodKD )
      if ( FlagHaveAllOptProp ) then
        DelWaveNum = a_DelWaveNum(iWaveNum)
      else
        String = MkRangeString( "BinNum", iWaveNum )
        call HistoryGet( OptPropNcFNSave, "BinWeight", DelWaveNum, &
          & range = string, quiet = .true. )
      end if
    end select


  end function OptPropGetDelWNCore

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

  function OptPropGetWNCore( iWaveNum ) result( WaveNum )


    integer, intent(in) :: iWaveNum

    real(DP) :: WaveNum

    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    select case ( IDMethod )
    case ( IDMethodLBL, IDMethodLBLConstPF )
      WaveNum = w_WaveNumSave( iWaveNum )
    case ( IDMethodKD )
      write( 6, * ) 'In ', trim( ModuleName )
      write( 6, * ) '  function OptPropGetWNCore cannot be used when IDMethod is IDMethodKD'
      stop
    case default
      write( 6, * ) 'In ', trim( ModuleName )
      write( 6, * ) '  Unexpected case default'
      stop
    end select


  end function OptPropGetWNCore

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

  function OptPropGetBandWNBndsCore( ) result( aa_BandWNBnds )


    real(DP) :: aa_BandWNBnds(2,NBandSave)

    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    select case ( IDMethod )
    case ( IDMethodLBL, IDMethodLBLConstPF )
      write( 6, * ) 'In ', trim( ModuleName )
      write( 6, * ) '  function OptPropGetBandWNBndsCore cannot be used when IDMethod is IDMethodLBL'
      stop
    case ( IDMethodKD )
      aa_BandWNBnds = aa_BandWaveNumBndsSave
    end select


  end function OptPropGetBandWNBndsCore

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

  function OptPropGetBandBinIndexBndsCore( ) result( aa_BandBinIndexBnds )


    integer :: aa_BandBinIndexBnds(2,NBandSave)


    ! local variables
    !


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    select case ( IDMethod )
    case ( IDMethodLBL, IDMethodLBLConstPF )
      write( 6, * ) 'In ', trim( ModuleName )
      write( 6, * ) '  function OptPropGetBandBinIndexBndsCore cannot be used when IDMethod is IDMethodLBL or IDMethodLBLConstPF'
      stop
    case ( IDMethodKD )
      aa_BandBinIndexBnds = aa_BandBinIndexBndsSave
    end select


  end function OptPropGetBandBinIndexBndsCore

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

  subroutine OptPropGetAbsCoefProfCore( &
    & iWaveNum, NPress, NMol, &
    & zm_AbsCoef &
    & )

    use rad_utils, only : MkRangeString


    integer     , intent(in ) :: iWaveNum
    integer     , intent(in ) :: NPress
    integer     , intent(in ) :: NMol
    real(DP)    , intent(out) :: zm_AbsCoef(1:NPress,1:NMol)


    ! local variables
    !
    character(128) :: String


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    if ( FlagHaveAllOptProp ) then
      zm_AbsCoef = zma_AbsCoef(:,:,iWaveNum)
    else
      String = MkRangeString( WaveNumAxisName, iWaveNum )
      call HistoryGet( OptPropNcFNSave, AbsCoefVarName, zm_AbsCoef, &
        & range = string, quiet = .true. )
    end if


  end subroutine OptPropGetAbsCoefProfCore

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

  subroutine OptPropGetPFIntedCore( &
    & iWaveNum, &
    & kmax, &
    & r_Temp, SurfTemp, &
    & r_PFInted, SurfPFInted, SurfDPFDTInted &
    & )

    use constants0, only : PI

    use planck_func_wrapper

    integer , intent(in ) :: iWaveNum
    integer , intent(in ) :: kmax
    real(DP), intent(in ) :: r_Temp(0:kmax)
    real(DP), intent(in ) :: SurfTemp
    real(DP), intent(out) :: r_PFInted  (0:kmax)
    real(DP), intent(out) :: SurfPFInted
    real(DP), intent(out) :: SurfDPFDTInted


    ! local variables
    !
    real(DP) :: WaveNum
    real(DP) :: DelWaveNum
    real(DP) :: WNs
    real(DP) :: WNe
    integer  :: GQNum


    select case ( IDMethod )
    case ( IDMethodLBL, IDMethodLBLConstPF )
      WaveNum    = OptPropGetWNCore( iWaveNum )
      DelWaveNum = OptPropGetDelWNCore( iWaveNum )
      WNs = WaveNum - DelWaveNum / 2.0d0
      WNe = WaveNum + DelWaveNum / 2.0d0
    case ( IDMethodKD )
      write( 6, * ) 'In ', trim( ModuleName )
      write( 6, * ) '  function OptPropGetPFIntedCore cannot be used when IDMethod is IDMethodKD'
      stop
    case default
      write( 6, * ) 'In ', trim( ModuleName )
      write( 6, * ) '  Unexpected case default'
      stop
    end select

!!$    GQNum = 1
    GQNum = int( DelWaveNum / 10.0d0 ) + 1

    call Integ_PF_GQ_Array1D( &
      & WNs, WNe, GQNum, &
      & 0, kmax, &
      & r_Temp, &
      & r_PFInted &
      & )
    r_PFInted = PI * r_PFInted
    call Integ_PF_GQ_Array0D( &
      & WNs, WNe, GQNum, &
      & SurfTemp, &
      & SurfPFInted &
      & )
    SurfPFInted = PI * SurfPFInted
    call Integ_DPFDT_GQ_Array0D( &
      & WNs, WNe, GQNum,                 & ! (in )
      & SurfTemp,                & ! (in )
      & SurfDPFDTInted                  & ! (out)
      & )
    SurfDPFDTInted = PI * SurfDPFDTInted


  end subroutine OptPropGetPFIntedCore

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

  subroutine OptPropGetPFIntedRatioCore( &
    & iWaveNum, NPress, &
    & z_PFRatio, SurfPFRatio &
    & )

    use rad_utils, only : MkRangeString


    integer     , intent(in ) :: iWaveNum
    integer     , intent(in ) :: NPress
    real(DP)    , intent(out) :: z_PFRatio(1:NPress)
    real(DP)    , intent(out) :: SurfPFRatio


    ! local variables
    !
    character(128) :: String

!!$    real(DP), allocatable :: za_PFRatio(:,:)


    if ( .not. FlagInited ) then
      write( 6, * ) trim( ModuleName ), " is not initialized."
      stop
    end if


    if ( FlagHaveAllOptProp ) then
      z_PFRatio = za_PFRatio(:,iWaveNum)
    else
      String = MkRangeString( WaveNumAxisName, iWaveNum )
      call HistoryGet( OptPropNcFNSave, PFRatioVarName, z_PFRatio, &
        & range = string, quiet = .true. )
    end if


!!$    allocate( za_PFRatio( NPressSave, NWaveNumSave ) )
!!$    call HistoryGet( OptPropNcFNSave, PFRatioVarName, za_PFRatio, &
!!$      & quiet = .true. )
!!$    z_PFRatio = za_PFRatio(:,iWaveNum)
!!$    deallocate( za_PFRatio )


    SurfPFRatio = z_PFRatio(1)


  end subroutine OptPropGetPFIntedRatioCore

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

end module opt_prop

