UPP (develop)
Loading...
Searching...
No Matches
CLDFRAC_ZHAO.f
1!-----------------------------------
2 subroutine progcld1 &
3!...................................
4
5! --- inputs:
6 ( plyr,tlyr,qlyr,qstl,clw, &
7! xlat,xlon, &
8 ix, nlay, iflip, &
9! --- outputs:
10 cldtot &
11 )
12
13! ================= subprogram documentation block ================ !
14! !
15! subprogram: progcld1 computes cloud related quantities using !
16! zhao/moorthi's prognostic cloud microphysics scheme. !
17! !
18! abstract: this program computes cloud fractions from cloud !
19! condensates, calculates liquid/ice cloud droplet effective radius, !
20! and computes the low, mid, high, total and boundary layer cloud !
21! fractions and the vertical indices of low, mid, and high cloud !
22! top and base. the three vertical cloud domains are set up in the !
23! initial subroutine "cldinit". !
24! !
25! program history log: !
26! 11-xx-1992 y.h., k.a.c, a.k. - cloud parameterization !
27! 'cldjms' patterned after slingo and slingo's work (jgr, !
28! 1992), stratiform clouds are allowed in any layer except !
29! the surface and upper stratosphere. the relative humidity !
30! criterion may cery in different model layers. !
31! 10-25-1995 kenneth campana - tuned cloud rh curves !
32! rh-cld relation from tables created using mitchell-hahn !
33! tuning technique on airforce rtneph observations. !
34! 11-02-1995 kenneth campana - the bl relationships used !
35! below llyr, except in marine stratus regions. !
36! 04-11-1996 kenneth campana - save bl cld amt in cld(,5) !
37! 12-29-1998 s. moorthi - prognostic cloud method !
38! 04-15-2003 yu-tai hou - rewritten in frotran 90 !
39! modulized form, seperate prognostic and diagnostic methods !
40! into two packages. !
41! !
42! usage: call progcld1 !
43! !
44! subprograms called: gethml !
45! !
46! attributes: !
47! language: fortran 90 !
48! machine: ibm-sp, sgi !
49! !
50! !
51! ==================== defination of variables ==================== !
52! !
53! input variables: !
54! plyr (IX,NLAY) : model layer mean pressure in mb (100Pa) !
55! tlyr (IX,NLAY) : model layer mean temperature in k !
56! qlyr (IX,NLAY) : layer specific humidity in gm/gm !
57! qstl (IX,NLAY) : layer saturate humidity in gm/gm !
58! rhly (IX,NLAY) : layer relative humidity (=qlyr/qstl) !
59! clw (IX,NLAY) : layer cloud condensate amount !
60! xlat (IX) : grid latitude in radians !
61! xlon (IX) : grid longitude in radians (not used) !
62! slmsk (IX) : sea/land mask array (sea:0,land:1,sea-ice:2) !
63! IX : horizontal dimention !
64! NLAY : vertical layer/level dimensions !
65! iflip : control flag for in/out vertical indexing !
66! =0: index from toa to surface !
67! =1: index from surface to toa !
68! iovr : control flag for cloud overlap !
69! =0 random overlapping clouds !
70! =1 max/ran overlapping clouds !
71! !
72! output variables: !
73! clouds(IX,NLAY,NF_CLDS) : cloud profiles !
74! clouds(:,:,1) - layer total cloud fraction !
75! clouds(:,:,2) - layer cloud liq water path (g/m**2) !
76! clouds(:,:,3) - mean eff radius for liq cloud (micron) !
77! clouds(:,:,4) - layer cloud ice water path (g/m**2) !
78! clouds(:,:,5) - mean eff radius for ice cloud (micron) !
79! clouds(:,:,6) - layer rain drop water path not assigned !
80! clouds(:,:,7) - mean eff radius for rain drop (micron) !
81! *** clouds(:,:,8) - layer snow flake water path not assigned !
82! clouds(:,:,9) - mean eff radius for snow flake (micron) !
83! *** fu's scheme need to be normalized by snow density (g/m**3/1.0e6) !
84! clds (IX,5) : fraction of clouds for low, mid, hi, tot, bl !
85! mtop (IX,3) : vertical indices for low, mid, hi cloud tops !
86! mbot (IX,3) : vertical indices for low, mid, hi cloud bases !
87! !
88! ==================== end of description ===================== !
89!
90 use kinds, only: r_kind
91 implicit none
92 real, parameter:: climit=0.001
93! --- inputs
94 integer, intent(in) :: IX, NLAY, iflip
95
96 real (kind=r_kind), dimension(ix,nlay), intent(in) :: plyr, &
97 tlyr, qlyr, qstl, clw
98
99! real (kind=r_kind), dimension(ix), intent(in) :: xlat, xlon
100
101! --- outputs
102! real (kind=kind_phys), dimension(:,:,:), intent(out) :: clouds
103
104 real (kind=r_kind), dimension(ix,nlay), intent(out) :: cldtot
105
106! integer, dimension(:,:), intent(out) :: mtop,mbot
107
108! --- local variables:
109! real (kind=kind_phys), dimension(IX,NLAY) :: cldtot, cldcnv, &
110! & cwp, cip, crp, csp, rew, rei, res, rer, delp, tem2d
111
112! real (kind=kind_phys) :: ptop1(IX,4)
113 real (kind=r_kind) :: rhly(ix,nlay)
114 real (kind=r_kind) :: clwmin, clwm, clwt, onemrh, value, &
115 tem1, tem2, tem3
116
117 integer, dimension(IX) :: kinver
118
119 integer :: i, k, id, id1
120
121 logical :: inversn(IX)
122
123!
124!===> ... begin here
125!
126! clouds(:,:,:) = 0.0
127
128!$omp parallel do private(i,k)
129 do k = 1, nlay
130 do i = 1, ix
131 cldtot(i,k) = 0.0
132 rhly(i,k) = qlyr(i,k)/qstl(i,k) ! Chuang: add RH computation here
133 rhly(i,k) = max(0.0,min(1.0,rhly(i,k))) ! moorthi
134 enddo
135 enddo
136
137! --- layer cloud fraction
138
139 if (iflip == 0) then ! input data from toa to sfc
140
141 do i = 1, ix
142 inversn(i) = .false.
143 kinver(i) = 1
144 enddo
145
146 do k = nlay-1, 2, -1
147!$omp parallel do private(i,tem1)
148 do i = 1, ix
149 if (plyr(i,k) > 600.0 .and. (.not.inversn(i))) then
150 tem1 = tlyr(i,k-1) - tlyr(i,k)
151
152 if (tem1 > 0.1 .and. tlyr(i,k) > 278.0) then
153 inversn(i) = .true.
154 kinver(i) = k
155 endif
156 endif
157 enddo
158 enddo
159
160 clwmin = 0.0
161 do k = nlay, 1, -1
162!$omp parallel do private(i,tem1,tem2,clwt,clwm,onemrh,value)
163 do i = 1, ix
164 clwt = 1.0e-6 * (plyr(i,k)*0.001)
165! clwt = 2.0e-6 * (plyr(i,k)*0.001)
166
167 if (clw(i,k) > clwt .or. &
168 (inversn(i) .and. k >= kinver(i)) ) then
169
170 onemrh= max( 1.e-10, 1.0-rhly(i,k) )
171 clwm = clwmin / max( 0.01, plyr(i,k)*0.001 )
172
173 tem1 = min(max(sqrt(sqrt(onemrh*qstl(i,k))),0.0001),1.0)
174 tem1 = 2000.0 / tem1
175! tem1 = 1000.0 / tem1
176! if (inversn(i) .and. k >= kinver(i)) tem1 = tem1 * 5.0
177
178 value = max( min( tem1*(clw(i,k)-clwm), 50.0 ), 0.0 )
179 tem2 = sqrt( sqrt(rhly(i,k)) )
180
181 cldtot(i,k) = max( tem2*(1.0-exp(-value)), 0.0 )
182 endif
183 enddo
184 enddo
185
186 else ! input data from sfc to toa
187
188 do i = 1, ix
189 inversn(i) = .false.
190 kinver(i) = nlay
191 enddo
192
193 do k = 2, nlay
194!$omp parallel do private(i,tem1)
195 do i = 1, ix
196 if (plyr(i,k) > 600.0 .and. (.not.inversn(i))) then
197 tem1 = tlyr(i,k+1) - tlyr(i,k)
198
199 if (tem1 > 0.1 .and. tlyr(i,k) > 278.0) then
200 inversn(i) = .true.
201 kinver(i) = k
202 endif
203 endif
204 enddo
205 enddo
206
207 clwmin = 0.0
208 do k = 1, nlay
209!$omp parallel do private(i,tem1,tem2,clwt,clwm,onemrh,value)
210 do i = 1, ix
211 clwt = 1.0e-6 * (plyr(i,k)*0.001)
212! clwt = 2.0e-6 * (plyr(i,k)*0.001)
213
214 if (clw(i,k) > clwt .or. &
215 (inversn(i) .and. k <= kinver(i)) ) then
216
217 onemrh= max( 1.e-10, 1.0-rhly(i,k) )
218 clwm = clwmin / max( 0.01, plyr(i,k)*0.001 )
219
220 tem1 = min(max(sqrt(sqrt(onemrh*qstl(i,k))),0.0001),1.0)
221 tem1 = 2000.0 / tem1
222! tem1 = 1000.0 / tem1
223! if (inversn(i) .and. k <= kinver(i)) tem1 = tem1 * 5.0
224
225 value = max( min( tem1*(clw(i,k)-clwm), 50.0 ), 0.0 )
226 tem2 = sqrt( sqrt(rhly(i,k)) )
227
228 cldtot(i,k) = max( tem2*(1.0-exp(-value)), 0.0 )
229 endif
230 enddo
231 enddo
232
233 endif ! end_if_flip
234
235 where (cldtot < climit)
236 cldtot = 0.0
237 endwhere
238!
239 return
240!...................................
241 end subroutine progcld1
242!-----------------------------------
243
244