module sfx_thermal_roughness
    !< @brief surface thermal roughness parameterizations

    ! modules used
    ! --------------------------------------------------------------------------------
    use sfx_phys_const
    use sfx_surface
    ! --------------------------------------------------------------------------------

    ! directives list
    ! --------------------------------------------------------------------------------
    implicit none
    ! --------------------------------------------------------------------------------

    ! public interface
    ! --------------------------------------------------------------------------------
    public :: get_thermal_roughness_kl
    public :: get_thermal_roughness_cz
    public :: get_thermal_roughness_zi
    public :: get_thermal_roughness_ca
    ! --------------------------------------------------------------------------------

    ! --------------------------------------------------------------------------------
    real, parameter, private :: kappa = 0.40         !< von Karman constant [n/d]
    ! --------------------------------------------------------------------------------

contains


    ! thermal roughness definition by (Kazakov, Lykosov)
    ! --------------------------------------------------------------------------------
    subroutine get_thermal_roughness_kl(z0_t, B, &
            z0_m, Re, surface_type)
        ! ----------------------------------------------------------------------------
        real, intent(out) :: z0_t               !< thermal roughness [m]
        real, intent(out) :: B                  !< = log(z0_m / z0_t) [n/d]

        real, intent(in) :: z0_m                !< aerodynamic roughness [m]
        real, intent(in) :: Re                  !< roughness Reynolds number [n/d]
        integer, intent(in) :: surface_type     !< = [ocean] || [land] || [lake]
        ! ----------------------------------------------------------------------------

        ! --- local variables
        ! ----------------------------------------------------------------------------

        !--- define B = log(z0_m / z0_t)
        if (Re <= Re_rough_min) then
            B = B1_rough * alog(B3_rough * Re) + B2_rough
        else
            ! *: B4 takes into account Re value at z' ~ O(10) z0
            B = B4_rough * (Re**B2_rough)
        end if

        ! --- apply max restriction based on surface type
        if (surface_type == surface_ocean) then
            B = min(B, B_max_ocean)
        else if (surface_type == surface_lake) then
            B = min(B, B_max_lake)
        else if (surface_type == surface_land) then
            B = min(B, B_max_land)
        end if

        ! --- define roughness [thermal]
        z0_t = z0_m / exp(B)

    end subroutine
    ! --------------------------------------------------------------------------------

    ! thermal roughness definition by (Chen, F., Zhang, Y., 2009)
    ! --------------------------------------------------------------------------------
    subroutine get_thermal_roughness_cz(z0_t, B, &
            z0_m, Re, surface_type)
        ! ----------------------------------------------------------------------------
        real, intent(out) :: z0_t               !< thermal roughness [m]
        real, intent(out) :: B                  !< = log(z0_m / z0_t) [n/d]

        real, intent(in) :: z0_m                !< aerodynamic roughness [m]
        real, intent(in) :: Re                  !< roughness Reynolds number [n/d]
        integer, intent(in) :: surface_type     !< = [ocean] || [land] || [lake]

        ! --- local variables
        ! ----------------------------------------------------------------------------

        !--- define B = log(z0_m / z0_t)
        B = (kappa * 10.0**(-0.4 * z0_m / 0.07)) * (Re**0.45)

        ! --- define roughness [thermal]
        z0_t = z0_m / exp(B)

    end subroutine
    ! --------------------------------------------------------------------------------
    
    ! thermal roughness definition by (Zilitinkevich, S., 1995)
    ! --------------------------------------------------------------------------------
    subroutine get_thermal_roughness_zi(z0_t, B, &
            z0_m, Re, surface_type)
        ! ----------------------------------------------------------------------------
        real, intent(out) :: z0_t               !< thermal roughness [m]
        real, intent(out) :: B                  !< = log(z0_m / z0_t) [n/d]

        real, intent(in) :: z0_m                !< aerodynamic roughness [m]
        real, intent(in) :: Re                  !< roughness Reynolds number [n/d]
        integer, intent(in) :: surface_type     !< = [ocean] || [land] || [lake]

        ! --- local variables
        ! ----------------------------------------------------------------------------

        !--- define B = log(z0_m / z0_t)
        B = 0.1 * kappa * (Re**0.5)

        ! --- define roughness [thermal]
        z0_t = z0_m / exp(B)

    end subroutine
    ! --------------------------------------------------------------------------------
    
    ! thermal roughness definition by (Cahill, A.T., Parlange, M.B., Albertson, J.D., 1997)
    ! --------------------------------------------------------------------------------
    subroutine get_thermal_roughness_ca(z0_t, B, &
            z0_m, Re, surface_type)
        ! ----------------------------------------------------------------------------
        real, intent(out) :: z0_t               !< thermal roughness [m]
        real, intent(out) :: B                  !< = log(z0_m / z0_t) [n/d]

        real, intent(in) :: z0_m                !< aerodynamic roughness [m]
        real, intent(in) :: Re                  !< roughness Reynolds number [n/d]
        integer, intent(in) :: surface_type     !< = [ocean] || [land] || [lake]

        ! --- local variables
        ! ----------------------------------------------------------------------------

        !--- define B = log(z0_m / z0_t)
        B = 2.46 * (Re**0.25) - 3.8

        ! --- define roughness [thermal]
        z0_t = z0_m / exp(B)

    end subroutine
    ! --------------------------------------------------------------------------------

end module sfx_thermal_roughness