CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_cloud_diagnostics.F90
1
3
5 use machine, only: kind_phys
7
8 ! Module parameters (imported directly from radiation_cloud.f)
9 integer, parameter :: &
10 nf_clds = 9, & ! Number of fields in cloud array
11 nk_clds = 3 ! Number of cloud vertical domains
12 real(kind_phys), parameter :: &
13 climit = 0.001, & ! Lowest allowable cloud-fraction
14 ovcst = 1.0 - 1.0e-8 ! Overcast cloud-fraction 0.999999999
15 real(kind_phys), parameter, dimension(NK_CLDS+1,2) :: &
16 ptopc = reshape(source=(/ 1050., 650., 400., 0.0, 1050., 750., 500., 0.0 /), &
17 shape=(/nk_clds+1,2/))
18
19 ! Version tag and last revision date
20 character(40), parameter :: vtagcld='UFS-cloud-diagnostics vX.x May 2020 '
21
23
24contains
25
34 subroutine gfs_cloud_diagnostics_run(nCol, nLev, iovr, iovr_rand, iovr_maxrand, &
35 iovr_max, iovr_dcorr, iovr_exp, iovr_exprand, lsswr, lslwr, lat, de_lgth, p_lay, &
36 cld_frac, p_lev, deltaZ, cloud_overlap_param, precip_overlap_param, con_pi, &
37 top_at_1, si, mtopa, mbota, cldsa, errmsg, errflg)
38 implicit none
39
40 ! Inputs
41 integer, intent(in) :: &
42 ncol, & ! Number of horizontal grid-points
43 nlev ! Number of vertical-layers
44 integer, intent(in) :: &
45 iovr, & ! Choice of cloud-overlap method
46 iovr_rand, & ! Flag for random cloud overlap method
47 iovr_maxrand, & ! Flag for maximum-random cloud overlap method
48 iovr_max, & ! Flag for maximum cloud overlap method
49 iovr_dcorr, & ! Flag for decorrelation-length cloud overlap method
50 iovr_exp, & ! Flag for exponential cloud overlap method
51 iovr_exprand ! Flag for exponential-random cloud overlap method
52 logical, intent(in) :: &
53 lsswr, & ! Call SW radiation?
54 lslwr, & ! Call LW radiation?
55 top_at_1 ! Vertical ordering flag
56 real(kind_phys), intent(in) :: &
57 con_pi ! Physical constant: pi
58 real(kind_phys), dimension(:), intent(in) :: &
59 lat, & ! Latitude
60 de_lgth, & ! Decorrelation length
61 si ! Vertical sigma coordinate
62 real(kind_phys), dimension(:,:), intent(in), optional :: &
63 p_lay ! Pressure at model-layer
64 real(kind_phys), dimension(:,:), intent(in) :: &
65 cld_frac ! Total cloud fraction
66 real(kind_phys), dimension(:,:), intent(in), optional :: &
67 p_lev ! Pressure at model interfaces
68 real(kind_phys), dimension(:,:), intent(in), optional :: &
69 deltaz, & ! Layer thickness (m)
70 cloud_overlap_param, & ! Cloud-overlap parameter
71 precip_overlap_param ! Precipitation overlap parameter
72
73 ! Outputs
74 character(len=*), intent(out) :: &
75 errmsg ! Error message
76 integer, intent(out) :: &
77 errflg ! Error flag
78 integer,dimension(:,:),intent(out) :: &
79 mbota, & ! Vertical indices for cloud tops
80 mtopa ! Vertical indices for cloud bases
81 real(kind_phys),dimension(:,:), intent(out) :: &
82 cldsa ! Fraction of clouds for low, middle, high, total and BL
83
84 ! Local variables
85 integer i,id,icol,ilay,icld
86 real(kind_phys) :: tem1
87 real(kind_phys),dimension(nCol,NK_CLDS+1) :: ptop1
88 real(kind_phys),dimension(nCol) :: rlat
89 real(kind_phys),dimension(nCol,nLev) :: cldcnv
90
91 if (.not. (lsswr .or. lslwr)) return
92
93 ! Initialize CCPP error handling variables
94 errmsg = ''
95 errflg = 0
96
97 ! This is set to zero in all of the progcld() routines and passed to gethml().
98 cldcnv(:,:) = 0._kind_phys
99
100 do icld = 1, nk_clds+1
101 tem1 = ptopc(icld,2) - ptopc(icld,1)
102 do i=1,ncol
103 rlat(i) = abs(lat(i) / con_pi )
104 ptop1(i,icld) = ptopc(icld,1) + tem1*max( 0.0, 4.0*rlat(i)-1.0 )
105 enddo
106 enddo
107
108 ! Compute low, mid, high, total, and boundary layer cloud fractions and clouds top/bottom
109 ! layer indices for low, mid, and high clouds. The three cloud domain boundaries are
110 ! defined by ptopc. The cloud overlapping method is defined by control flag 'iovr', which may
111 ! be different for lw and sw radiation programs.
112 call gethml(p_lay*0.01, ptop1, cld_frac, cldcnv, deltaz, de_lgth, cloud_overlap_param,&
113 ncol, nlev, iovr, iovr_rand, iovr_maxrand, iovr_max, iovr_dcorr, iovr_exp, &
114 iovr_exprand, top_at_1, si, cldsa, mtopa, mbota)
115
116 end subroutine gfs_cloud_diagnostics_run
118end module gfs_cloud_diagnostics
subroutine, public gfs_cloud_diagnostics_run(ncol, nlev, iovr, iovr_rand, iovr_maxrand, iovr_max, iovr_dcorr, iovr_exp, iovr_exprand, lsswr, lslwr, lat, de_lgth, p_lay, cld_frac, p_lev, deltaz, cloud_overlap_param, precip_overlap_param, con_pi, top_at_1, si, mtopa, mbota, cldsa, errmsg, errflg)
real(kind=kind_phys), dimension(nk_clds+1, 2), save ptopc
pressure limits of cloud domain interfaces (low, mid, high) in mb (0.1kPa)
integer, parameter, public nk_clds
number of cloud vertical domains
real(kind=kind_phys), parameter ovcst
subroutine, public gethml(plyr, ptop1, cldtot, cldcnv, dz, de_lgth, alpha, ix, nlay, iovr, iovr_rand, iovr_maxrand, iovr_max, iovr_dcorr, iovr_exp, iovr_exprand, top_at_1, si, clds, mtop, mbot)
This subroutine computes high, mid, low, total, and boundary cloud fractions and cloud top/bottom lay...
integer, parameter, public nf_clds
number of fields in cloud array
real(kind=kind_phys), parameter climit
This module computes cloud related quantities for radiation computations.