CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
radiation_cloud_overlap.F90
1
3
7 use machine, only : kind_phys
8 implicit none
9 public :: cmp_dcorr_lgth
10
12 module procedure cmp_dcorr_lgth_hogan
13 module procedure cmp_dcorr_lgth_oreopoulos
14 end interface
15
16contains
17
18 ! ######################################################################################
19 ! Hogan et al. (2010)
20 ! "Effect of improving representation of horizontal and vertical cloud structure on the
21 ! Earth's global radiation budget. Part I: Review and parametrization"
22 ! https://rmets.onlinelibrary.wiley.com/doi/full/10.1002/qj.647
23 ! ######################################################################################
25 subroutine cmp_dcorr_lgth_hogan(nCol, lat, con_pi, dcorr_lgth)
26 ! Inputs
27 integer, intent(in) :: &
28 nCol ! Number of horizontal grid-points
29 real(kind_phys), intent(in) :: &
30 con_pi ! Physical constant: Pi
31 real(kind_phys), dimension(:), intent(in) :: &
32 lat ! Latitude
33 ! Outputs
34 real(kind_phys), dimension(:),intent(out) :: &
35 dcorr_lgth ! Decorrelation length
36
37 ! Local variables
38 integer :: iCol
39
40 ! Parameters
41 real(kind_phys),parameter :: min_dcorr = 0.6 ! (see section 2.3)
42
43 do icol =1,ncol
44 dcorr_lgth(icol) = max(min_dcorr, 2.78-4.6*abs(lat(icol)/con_pi)) ! (eq. 13)
45 enddo
46
47 end subroutine cmp_dcorr_lgth_hogan
48 ! ######################################################################################
49 ! Oreopoulos et al. (2012)
50 ! "Radiative impacts of cloud heterogeneity and overlap in an
51 ! atmospheric General Circulation Model"
52 ! 10.5194/acp-12-9097-2012
53 ! ######################################################################################
55 subroutine cmp_dcorr_lgth_oreopoulos(nCol, lat, juldat, yearlength, dcorr_lgth)
56 ! Inputs
57 integer, intent(in) :: &
58 nCol, & ! Number of horizontal grid-points
59 yearlength ! Number of days in year
60
61 real(kind_phys), intent(in) :: &
62 juldat ! Julian date
63 real(kind_phys), dimension(:), intent(in) :: &
64 lat ! Latitude
65
66 ! Outputs
67 real(kind_phys), dimension(:),intent(out) :: &
68 dcorr_lgth ! Decorrelation length (km)
69
70 ! Parameters for the Gaussian fits per Eqs. (10) and (11) (See Table 1)
71 real(kind_phys), parameter :: & !
72 am1 = 1.4315_kind_phys, & !
73 am2 = 2.1219_kind_phys, & !
74 am4 = -25.584_kind_phys, & !
75 amr = 7.0_kind_phys !
76
77 ! Local variables
78 integer :: iCol
79 real(kind_phys) :: am3
80
81 do icol = 1, ncol
82 if (juldat .gt. 181._kind_phys) then
83 am3 = -4._kind_phys * amr * (juldat - 272._kind_phys) / yearlength ! (eq. 11a)
84 else
85 am3 = 4._kind_phys * amr * (juldat - 91._kind_phys) / yearlength ! (eq. 11b)
86 endif
87 dcorr_lgth(icol) = am1 + am2 * exp( -(lat(icol) - am3)**2 / am4**2) ! (eq. 10)
88 enddo
89
90 end subroutine cmp_dcorr_lgth_oreopoulos
91
93 subroutine get_alpha_exper(nCol, nLay, iovr, iovr_exprand, dzlay, &
94 dcorr_lgth, cld_frac, alpha)
95
96 ! Inputs
97 integer, intent(in) :: &
98 nCol, & ! Number of horizontal grid points
99 nLay ! Number of vertical grid points
100 integer, intent(in) :: &
101 iovr, &
102 iovr_exprand
103 real(kind_phys), dimension(:), intent(in) :: &
104 dcorr_lgth ! Decorrelation length (km)
105 real(kind_phys), dimension(:,:), intent(in) :: &
106 dzlay !
107 real(kind_phys), dimension(:,:), intent(in) :: &
108 cld_frac
109
110 ! Outputs
111 real(kind_phys), dimension(:,:) :: &
112 alpha ! Cloud overlap parameter
113
114 ! Local variables
115 integer :: iCol,iLay
116
117 do icol = 1, ncol
118 alpha(icol,1) = 0.0d0
119 do ilay = 2, nlay
120 alpha(icol,ilay) = exp( -(dzlay(icol,ilay)) / dcorr_lgth(icol))
121 enddo
122 enddo
123
124 ! Revise alpha for exponential-random cloud overlap
125 ! Decorrelate layers when a clear layer follows a cloudy layer to enforce
126 ! random correlation between non-adjacent blocks of cloudy layers
127 if (iovr == iovr_exprand) then
128 do ilay = 2, nlay
129 do icol = 1, ncol
130 if (cld_frac(icol,ilay) == 0.0 .and. cld_frac(icol,ilay-1) > 0.0) then
131 alpha(icol,ilay) = 0.0
132 endif
133 enddo
134 enddo
135 endif
136
137 return
138
139 end subroutine get_alpha_exper