#include <misc.h>
#include <params.h>
subroutine cldnrh(lchnk ,ncol , & 1,5
pmid ,t ,q ,omga , &
cnt ,cnb ,cldn ,clc ,pdel , &
cmfmc ,cmfmc2 ,landfrac,snowh ,concld ,cldst , &
ts ,sst ,ps ,zdu ,ocnfrac , &
rhdfda ,rhu00 ,phis)
!-----------------------------------------------------------------------
!
! Purpose:
! Interface to compute cloud, and its response to rh perturbation
!
!
! Author: W. Lin
!
!-----------------------------------------------------------------------
use shr_kind_mod
, only: r8 => shr_kind_r8
use ppgrid
use cloud_fraction
, only: cldfrc
implicit none
!------------------------------Parameters-------------------------------
real(r8) pnot ! reference pressure
parameter (pnot = 1.e5)
!------------------------------Arguments--------------------------------
!
! Input arguments
!
integer, intent(in) :: lchnk ! chunk identifier
integer, intent(in) :: ncol ! number of atmospheric columns
real(r8), intent(in) :: pmid(pcols,pver) ! midpoint pressures
real(r8), intent(in) :: t(pcols,pver) ! temperature
real(r8), intent(in) :: q(pcols,pver) ! specific humidity
real(r8), intent(in) :: omga(pcols,pver) ! vertical pressure velocity
real(r8), intent(in) :: cnt(pcols) ! top level of convection
real(r8), intent(in) :: cnb(pcols) ! bottom level of convection
real(r8), intent(in) :: phis(pcols) ! surface geopotential
real(r8), intent(in) :: cmfmc(pcols,pverp) ! convective mass flux--m sub c
real(r8), intent(in) :: cmfmc2(pcols,pverp) ! shallow convective mass flux--m sub c
real(r8), intent(in) :: snowh(pcols) ! snow depth (liquid water equivalent)
real(r8), intent(in) :: pdel(pcols,pver) ! pressure depth of layer
real(r8), intent(in) :: landfrac(pcols) ! Land fraction
real(r8), intent(in) :: ocnfrac(pcols) ! Ocean fraction
real(r8), intent(in) :: ts(pcols) ! surface temperature
real(r8), intent(in) :: sst(pcols) ! sea surface temperature
real(r8), intent(in) :: ps(pcols) ! surface pressure
real(r8), intent(in) :: zdu(pcols,pver) ! detrainment rate from deep convection
!
! Output arguments
!
real(r8), intent(out) :: cldn(pcols,pver) ! cloud fraction
real(r8), intent(out) :: clc(pcols) ! column convective cloud amount
real(r8), intent(out) :: cldst(pcols,pver) ! cloud fraction
real(r8), intent(out) :: rhdfda(pcols,pver) ! d_RH/d_cloud_fraction ====wlin
real(r8), intent(out) :: rhu00(pcols,pver) ! RH limit, U00 ====wlin
real(r8), intent(out) :: concld(pcols,pver) ! convective cloud cover
!
! local work space
!
real(r8) relhum(pcols,pver) ! RH, output to determine drh/da
real(r8) rhu002(pcols,pver) ! same as rhu00 but for perturbed rh
real(r8) cldn2(pcols,pver) ! same as cldn but for perturbed rh
real(r8) concld2(pcols,pver) ! same as concld but for perturbed rh
real(r8) cldst2(pcols,pver) ! same as cldst but for perturbed rh
real(r8) relhum2(pcols,pver) ! RH after perturbation
integer i,k ! longitude, level indices
!
call t_startf("cldfrc")
call cldfrc
(lchnk, ncol, &
pmid, t, q, omga, phis, &
cnt, cnb, cldn, clc, pdel, &
cmfmc, cmfmc2, landfrac,snowh, concld, cldst, &
ts, sst, ps, zdu, ocnfrac, rhu00, &
relhum, 0 )
! re-calculate cloud with perturbed rh add call cldfrc
call cldfrc
(lchnk, ncol, &
pmid, t, q, omga, phis, &
cnt, cnb, cldn2, clc, pdel, &
cmfmc, cmfmc2 ,landfrac,snowh, concld2, cldst2, &
ts, sst, ps, zdu, ocnfrac, rhu002, &
relhum2, 1 )
call t_stopf("cldfrc")
! cldfrc does not define layer cloud for model layer at k=1
! so set rhu00(k=1)=2.0 to not calculate cme for this layer
rhu00(:ncol,1) = 2.0
! Add following to estimate rhdfda
do k=1,pver
do i=1,ncol
if(relhum(i,k) < rhu00(i,k) ) then
rhdfda(i,k)=0.0
else if (relhum(i,k) >= 1.0 ) then
rhdfda(i,k)=0.0
else
!under certain circumstances, rh+ cause cld not to changed
!when at an upper limit, or w/ strong subsidence
!need to further check whether this if-block is necessary
if((cldn2(i,k) - cldn(i,k) ) < 1.e-4 ) then
rhdfda(i,k) = 0.01*relhum(i,k)*1.e+4 !instead of 0.0
else
rhdfda(i,k)=0.01*relhum(i,k)/(cldn2(i,k)-cldn(i,k))
endif
endif
enddo
enddo
return
end subroutine cldnrh