!! Fortran version of different diagnostics
! L. Fita. LMD May 2016
! gfortran module_generic.o -c module_ForDiagnostics.F90
!
! f2py -m module_ForDiagnostics --f90exec=/usr/bin/gfortran-4.7 -c module_generic.F90 module_ForDiagnostics.F90

MODULE module_ForDiagnosticsVars

  USE module_generic

  IMPLICIT NONE

!  INTEGER, PARAMETER                                     :: r_k = KIND(1.d0)
  REAL(r_k), PARAMETER                                   :: ZEPSEC=1.0D-12
! Low limit pressure for medium clouds [Pa]
  REAL(r_k), PARAMETER                                   :: prmhc = 44000.d0
! Low limit pressure for High clouds [Pa]
  REAL(r_k), PARAMETER                                   :: prmlc = 68000.d0

  REAL(r_k), PARAMETER                                   :: zero=0.d0
  REAL(r_k), PARAMETER                                   :: one=1.d0
  REAL(r_k), PARAMETER                                   :: two=2.d0

  CONTAINS

!!!!!!! Variables
! var_cllmh: low, medium, high-cloud [0,1]
! var_clt: total cloudiness [0,1]

!!!!!!! Calculations
! compute_clt: Computation of total cloudiness

!!!
! Variables
!!!

  FUNCTION var_cllmh(clfra, p, dz)
! Function to compute cllmh on a 1D column 1: low-cloud; 2: medium-cloud; 3: high-cloud [1]

    IMPLICIT NONE

!    INTEGER, PARAMETER                                   :: r_k = KIND(1.d0)
    INTEGER, INTENT(in)                                  :: dz
    REAL(r_k), DIMENSION(dz), INTENT(in)                 :: clfra, p
    REAL(r_k), DIMENSION(3)                              :: var_cllmh

! Local
    INTEGER                                              :: iz
    CHARACTER(LEN=50)                                    :: fname
    REAL(r_k)                                            :: zclearl, zcloudl, zclearm, zcloudm,       &
      zclearh, zcloudh

!!!!!!! Variables
! clfra: cloudfraction as 1D verical-column [1]
! p: pressure values of the column
    fname = 'var_cllmh'

    zclearl = one
    zcloudl = zero
    zclearm = one
    zcloudm = zero
    zclearh = one
    zcloudh = zero

    var_cllmh = one

    DO iz=1, dz
      IF (p(iz) < prmhc) THEN
        var_cllmh(3) = var_cllmh(3)*(one-MAX(clfra(iz),zcloudh))/(one-MIN(zcloudh,one-ZEPSEC))
        zcloudh = clfra(iz)
      ELSE IF ( (p(iz) >= prmhc) .AND. (p(iz) < prmlc) ) THEN
        var_cllmh(2) = var_cllmh(2)*(one-MAX(clfra(iz),zcloudm))/(one-MIN(zcloudm,one-ZEPSEC))
        zcloudm = clfra(iz)
      ELSE IF (p(iz) >= prmlc) THEN
        var_cllmh(1) = var_cllmh(1)*(one-MAX(clfra(iz),zcloudl))/(one-MIN(zcloudl,one-ZEPSEC))
        zcloudl = clfra(iz)
      ELSE
        PRINT *,'  ' // TRIM(fname) // ': This is weird, pressure:', p(iz), ' Pa fails out!!'
        PRINT *,'    from high, low cloud pressures:', prmhc, ' ,', prmlc,' Pa at z-level:', iz
        PRINT *,'    p_high > p:', prmhc,'> ',p(iz),' Pa'
        PRINT *,'    p_low > p >= p_high:', prmlc,'> ',p(iz),' >=', prmhc,' Pa'
        PRINT *,'    p_low >= p:', prmlc,'>= ',p(iz),' Pa'
        STOP
      END IF
    END DO

    var_cllmh = one - var_cllmh

    RETURN 

  END FUNCTION var_cllmh

  REAL(r_k) FUNCTION var_clt(clfra, dz)
! Function to compute the total cloud following 'newmicro.F90' from LMDZ using 1D vertical 
!   column values

    IMPLICIT NONE

!    INTEGER, PARAMETER                                   :: r_k = KIND(1.d0)
    INTEGER, INTENT(in)                                  :: dz
    REAL(r_k), DIMENSION(dz), INTENT(in)                 :: clfra
! Local
    INTEGER                                              :: iz
    REAL(r_k)                                            :: zclear, zcloud
    CHARACTER(LEN=50)                                    :: fname
!!!!!!! Variables
! cfra: 1-column cloud fraction values

    fname = 'var_clt'

    zclear = one
    zcloud = zero

    DO iz=1,dz
      zclear = zclear*(one-MAX(clfra(iz),zcloud))/(one-MIN(zcloud,1.-ZEPSEC))
      var_clt = one - zclear
      zcloud = clfra(iz)
    END DO

    RETURN

  END FUNCTION var_clt

END MODULE module_ForDiagnosticsVars