CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_rrtmg_pre.F90
1
3
5
7
8 contains
9
13
14 ! Attention - the output arguments lm, im, lmk, lmp must not be set
15 ! in the CCPP version - they are defined in the interstitial_create routine
20 subroutine gfs_rrtmg_pre_run (im, levs, lm, lmk, lmp, n_var_lndp, lextop,&
21 ltp, imfdeepcnv, imfdeepcnv_gf, imfdeepcnv_c3, me, ncnd, ntrac, &
22 num_p3d, npdf3d, xr_cnvcld, &
23 ncnvcld3d,ntqv, ntcw,ntiw, ntlnc, ntinc, ntrnc, ntsnc, ntccn, top_at_1,&
24 ntrw, ntsw, ntgl, nthl, ntwa, ntoz, ntsmoke, ntdust, ntcoarsepm, &
25 ntclamt, nleffr, nieffr, nseffr, lndp_type, kdt, &
26 ntdu1, ntdu2, ntdu3, ntdu4, ntdu5, ntss1, ntss2, &
27 ntss3, ntss4, ntss5, ntsu, ntbcb, ntbcl, ntocb, ntocl, ntchm, &
28 imp_physics,imp_physics_nssl, nssl_ccn_on, nssl_invertccn, &
29 imp_physics_thompson, imp_physics_gfdl, imp_physics_zhao_carr, &
30 imp_physics_zhao_carr_pdf, imp_physics_mg, imp_physics_wsm6, &
31 imp_physics_fer_hires, iovr, iovr_rand, iovr_maxrand, iovr_max, &
32 iovr_dcorr, iovr_exp, iovr_exprand, idcor, idcor_con, idcor_hogan, &
33 idcor_oreopoulos, dcorr_con, julian, yearlen, lndp_var_list, lsswr, &
34 lslwr, ltaerosol, mraerosol, lgfdlmprad, uni_cld, effr_in, do_mynnedmf,&
35 lmfshal, lcnorm, lmfdeep2, lcrick, fhswr, fhlwr, solhr, sup, con_eps, &
36 epsm1, fvirt, rog, rocp, con_rd, xlat_d, xlat, xlon, coslat, sinlat, &
37 tsfc, slmsk, prsi, prsl, prslk, tgrs, sfc_wts, mg_cld, effrr_in, &
38 pert_clds, sppt_wts, sppt_amp, cnvw_in, cnvc_in, qgrs, aer_nm, dx, &
39 icloud, iaermdl, iaerflg, con_pi, con_g, con_ttp, con_thgni, si, & !inputs from here and above
40 coszen, coszdg, effrl_inout, effri_inout, effrs_inout, &
41 clouds1, clouds2, clouds3, clouds4, clouds5, qci_conv, & !in/out from here and above
42 kd, kt, kb, mtopa, mbota, raddt, tsfg, tsfa, de_lgth, alb1d, delp, dz, & !output from here and below
43 plvl, plyr, tlvl, tlyr, qlyr, olyr, gasvmr_co2, gasvmr_n2o, gasvmr_ch4,&
44 gasvmr_o2, gasvmr_co, gasvmr_cfc11, gasvmr_cfc12, gasvmr_cfc22, &
45 gasvmr_ccl4, gasvmr_cfc113, aerodp,ext550, clouds6, clouds7, clouds8, &
46 clouds9, cldsa, cldfra, cldfra2d, lwp_ex,iwp_ex, lwp_fc,iwp_fc, &
47 faersw1, faersw2, faersw3, faerlw1, faerlw2, faerlw3, alpha, rrfs_sd, &
48 aero_dir_fdb, fdb_coef, spp_wts_rad, spp_rad, ico2, ozphys, &
49 errmsg, errflg)
50
51 use machine, only: kind_phys
52
53 use radcons, only: itsfc, qmin, qme5, qme6, epsq, prsmin
54 use funcphys, only: fpvs
55
56 use module_radiation_astronomy,only: coszmn ! sol_init, sol_update
57 use module_radiation_gases, only: nf_vgas, getgases ! gas_init, gas_update,
58 use module_radiation_aerosols, only: nf_aesw, nf_aelw, setaer, & ! aer_init, aer_update,
59 & nspc1
60 use module_radiation_clouds, only: nf_clds, & ! cld_init
62 & cal_cldfra3, &
67
69 & profsw_type, nbdsw
71 & proflw_type, nbdlw
73
74 ! For Thompson MP
76 nt_c_l, nt_c_o, &
77 re_qc_min, re_qc_max, &
78 re_qi_min, re_qi_max, &
79 re_qs_min, re_qs_max
84 ! For NRL Ozone
85 use module_ozphys, only: ty_ozphys
86 implicit none
87
88 integer, intent(in) :: im, levs, lm, lmk, lmp, ltp, &
89 n_var_lndp, imfdeepcnv, &
90 imfdeepcnv_gf, imfdeepcnv_c3, &
91 me, ncnd, ntrac, &
92 num_p3d, npdf3d, ncnvcld3d, ntqv, &
93 ntcw, ntiw, ntlnc, ntinc, &
94 ntrnc, ntsnc,ntccn, &
95 ntrw, ntsw, ntgl, nthl, ntwa, ntoz, &
96 ntsmoke, ntdust, ntcoarsepm, &
97 ntclamt, nleffr, nieffr, nseffr, &
98 lndp_type, &
99 kdt, imp_physics, &
100 imp_physics_thompson, &
101 imp_physics_gfdl, &
102 imp_physics_zhao_carr, &
103 imp_physics_zhao_carr_pdf, &
104 imp_physics_mg, imp_physics_wsm6, &
105 imp_physics_nssl, &
106 imp_physics_fer_hires, &
107 yearlen, icloud, iaermdl, iaerflg
108
109 integer, intent(in) :: &
110 iovr, & ! choice of cloud-overlap method
111 iovr_rand, & ! Flag for random cloud overlap method
112 iovr_maxrand, & ! Flag for maximum-random cloud overlap method
113 iovr_max, & ! Flag for maximum cloud overlap method
114 iovr_dcorr, & ! Flag for decorrelation-length cloud overlap method
115 iovr_exp, & ! Flag for exponential cloud overlap method
116 iovr_exprand, & ! Flag for exponential-random cloud overlap method
117 idcor_con, &
118 idcor, &
119 idcor_hogan, &
120 idcor_oreopoulos, &
121 ico2 ! Flag for co2 source used in radiation
122
123 integer, intent(in) :: ntdu1, ntdu2, ntdu3, ntdu4, ntdu5, ntss1, ntss2, ntss3, &
124 ntss4, ntss5, ntsu, ntbcb, ntbcl, ntocb, ntocl, ntchm
125
126 character(len=3), dimension(:), intent(in), optional :: lndp_var_list
127
128 logical, intent(in) :: lsswr, lslwr, ltaerosol, lgfdlmprad, &
129 uni_cld, effr_in, do_mynnedmf, &
130 lmfshal, lmfdeep2, pert_clds, lcrick,&
131 lcnorm, top_at_1, lextop, mraerosol
132 logical, intent(in) :: rrfs_sd, aero_dir_fdb, xr_cnvcld
133
134 logical, intent(in) :: nssl_ccn_on, nssl_invertccn
135 integer, intent(in) :: spp_rad
136 real(kind_phys), intent(in), optional :: spp_wts_rad(:,:)
137
138 real(kind=kind_phys), intent(in) :: fhswr, fhlwr, solhr, sup, julian, sppt_amp, dcorr_con
139 real(kind=kind_phys), intent(in) :: con_eps, epsm1, fvirt, rog, rocp, con_rd, con_pi, con_g, con_ttp, con_thgni
140
141 real(kind=kind_phys), dimension(:), intent(in) :: xlat_d, xlat, xlon, &
142 coslat, sinlat, tsfc, &
143 slmsk, dx, si
144
145 real(kind=kind_phys), dimension(:,:), intent(in) :: prsi, prsl, prslk, &
146 tgrs
147 real(kind=kind_phys), dimension(:,:), intent(in), optional :: &
148 mg_cld, effrr_in, &
149 cnvw_in, cnvc_in, &
150 sppt_wts, sfc_wts
151
152 real(kind=kind_phys), dimension(:,:,:), intent(in) :: qgrs
153 real(kind=kind_phys), dimension(:,:,:), intent(inout) :: aer_nm
154
155 real(kind=kind_phys), dimension(:), intent(inout) :: coszen, coszdg
156
157 real(kind=kind_phys), dimension(:,:), intent(inout), optional :: effrl_inout, &
158 effri_inout, &
159 effrs_inout
160 real(kind=kind_phys), dimension(:,:), intent(inout) :: clouds1, &
161 clouds2, clouds3, &
162 clouds4, clouds5
163 real(kind=kind_phys), dimension(:,:), intent(in), optional :: qci_conv
164 real(kind=kind_phys), dimension(:), intent(in) :: fdb_coef
165 real(kind=kind_phys), dimension(:), intent(out) :: lwp_ex,iwp_ex, &
166 lwp_fc,iwp_fc
167
168 integer, intent(out) :: kd, kt, kb
169
170 integer, dimension(:,:), intent(out) :: mbota, mtopa
171
172 real(kind=kind_phys), intent(out) :: raddt
173
174 real(kind=kind_phys), dimension(:), intent(out) :: tsfg, tsfa
175 real(kind=kind_phys), dimension(:), intent(out) :: de_lgth, &
176 alb1d
177
178 real(kind=kind_phys), dimension(:,:), intent(out) :: delp, dz, &
179 plyr, tlyr, &
180 qlyr, olyr
181
182 real(kind=kind_phys), dimension(:,:), intent(out) :: plvl, tlvl
183
184 real(kind=kind_phys), dimension(:,:), intent(out) :: gasvmr_co2, &
185 gasvmr_n2o, &
186 gasvmr_ch4, &
187 gasvmr_o2, &
188 gasvmr_co, &
189 gasvmr_cfc11,&
190 gasvmr_cfc12,&
191 gasvmr_cfc22,&
192 gasvmr_ccl4,&
193 gasvmr_cfc113
194 real(kind=kind_phys), dimension(:,:), intent(out) :: aerodp
195 real(kind=kind_phys), dimension(:,:), intent(out) :: ext550
196 real(kind=kind_phys), dimension(:,:), intent(out) :: clouds6, &
197 clouds7, &
198 clouds8, &
199 clouds9, &
200 cldfra
201 real(kind=kind_phys), dimension(:), intent(out) :: cldfra2d
202 real(kind=kind_phys), dimension(:,:), intent(out) :: cldsa
203
204 real(kind=kind_phys), dimension(:,:,:), intent(out) :: faersw1,&
205 faersw2,&
206 faersw3
207
208 real(kind=kind_phys), dimension(:,:,:), intent(out) :: faerlw1,&
209 faerlw2,&
210 faerlw3
211 real(kind=kind_phys), dimension(:,:), intent(out) :: alpha
212 character(len=*), intent(out) :: errmsg
213 integer, intent(out) :: errflg
214
215 ! Local variables
216 integer :: ncndl
217
218 integer :: i, j, k, k1, k2, lsk, lv, n, itop, ibtc, lp1, lla, llb, lya,lyb
219
220 real(kind=kind_phys) :: es, qs, delt, tem0d, pfac
221 real(kind=kind_phys), dimension(im) :: gridkm
222
223 real(kind=kind_phys), dimension(im) :: cvt1, cvb1, tem1d, tskn, xland
224
225 real(kind=kind_phys), dimension(im,lm+LTP) :: &
226 htswc, htlwc, gcice, grain, grime, htsw0, htlw0, &
227 rhly, tvly,qstl, vvel, clw, ciw, prslk1, tem2da, &
228 dzb, hzb, cldcov, deltaq, cnvc, cnvw, &
229 effrl, effri, effrr, effrs, rho, orho, plyrpa
230
231 ! for Thompson MP
232 real(kind=kind_phys), dimension(im,lm+LTP) :: &
233 qv_mp, qc_mp, qi_mp, qs_mp, &
234 nc_mp, ni_mp, nwfa
235 real (kind=kind_phys), dimension(lm) :: cldfra1d, qv1d, &
236 & qc1d, qi1d, qs1d, dz1d, p1d, t1d
237
238 ! for F-A MP
239 real(kind=kind_phys), dimension(im,lm+LTP+1) :: tem2db, hz
240
241 real(kind=kind_phys), dimension(im,lm+LTP,min(4,ncnd)) :: ccnd
242 real(kind=kind_phys), dimension(im,lm+LTP,2:ntrac) :: tracer1
243 real(kind=kind_phys), dimension(im,lm+LTP) :: &
244 & cld_frac, cld_lwp, cld_reliq, cld_iwp, cld_reice, &
245 & cld_rwp, cld_rerain, cld_swp, cld_resnow
246 real(kind=kind_phys), dimension(im,lm+LTP,NF_VGAS) :: gasvmr
247 real(kind=kind_phys), dimension(im,lm+LTP,NBDSW,NF_AESW) :: faersw
248 real(kind=kind_phys), dimension(im,lm+LTP,NBDLW,NF_AELW) :: faerlw
249
250 ! for stochastic cloud perturbations
251 real(kind=kind_phys), dimension(im) :: cldp1d
252 real (kind=kind_phys) :: alpha0,beta0,m,s,cldtmp,tmp_wt,cdfz
253 real (kind=kind_phys) :: max_relh
254 integer :: iflag
255 integer :: islmsk
256
257 ! For NRL Ozone
258 type(ty_ozphys),intent(in) :: ozphys
259
260 integer :: ids, ide, jds, jde, kds, kde, &
261 ims, ime, jms, jme, kms, kme, &
262 its, ite, jts, jte, kts, kte
263
264 real(kind=kind_phys) :: qvs
265!
266!===> ... begin here
267!
268 ! Initialize CCPP error handling variables
269 errmsg = ''
270 errflg = 0
271
272 if (.not. (lsswr .or. lslwr)) return
273
274 !--- set commonly used integers
275 ncndl = min(ncnd,4)
276
277 lp1 = lm + 1 ! num of in/out levels
278
279 if (imp_physics == imp_physics_thompson) then
280 max_relh = 1.5
281 else
282 max_relh = 1.1
283 endif
284
285 do i = 1, im
286 gridkm(i) = dx(i)*0.001
287 lwp_ex(i) = 0.0
288 iwp_ex(i) = 0.0
289 lwp_fc(i) = 0.0
290 iwp_fc(i) = 0.0
291 enddo
292
293! --- ... set local /level/layer indexes corresponding to in/out
294! variables
295
296 if ( lextop ) then
297 if (.not. top_at_1) then ! vertical from sfc upward
298 kd = 0 ! index diff between in/out and local
299 kt = 1 ! index diff between lyr and upper bound
300 kb = 0 ! index diff between lyr and lower bound
301 lla = lmk ! local index at the 2nd level from top
302 llb = lmp ! local index at toa level
303 lya = lm ! local index for the 2nd layer from top
304 lyb = lp1 ! local index for the top layer
305 else ! vertical from toa downward
306 kd = 1 ! index diff between in/out and local
307 kt = 0 ! index diff between lyr and upper bound
308 kb = 1 ! index diff between lyr and lower bound
309 lla = 2 ! local index at the 2nd level from top
310 llb = 1 ! local index at toa level
311 lya = 2 ! local index for the 2nd layer from top
312 lyb = 1 ! local index for the top layer
313 endif ! end if_top_at_1_block
314 else
315 kd = 0
316 if (.not. top_at_1) then ! vertical from sfc upward
317 kt = 1 ! index diff between lyr and upper bound
318 kb = 0 ! index diff between lyr and lower bound
319 else ! vertical from toa downward
320 kt = 0 ! index diff between lyr and upper bound
321 kb = 1 ! index diff between lyr and lower bound
322 endif ! end if_top_at_1_block
323 endif ! end if_lextop_block
324
325 raddt = min(fhswr, fhlwr)
326! print *,' in grrad : raddt=',raddt
327
328
331
332 if ( itsfc == 0 ) then ! use same sfc skin-air/ground temp
333 do i = 1, im
334 tskn(i) = tsfc(i)
335 tsfg(i) = tsfc(i)
336 enddo
337 else ! use diff sfc skin-air/ground temp
338 do i = 1, im
339 tskn(i) = tsfc(i)
340 tsfg(i) = tsfc(i)
341 enddo
342 endif
343
344
346!
347
348 lsk = 0
349 if (top_at_1 .and. lm < levs) lsk = levs - lm
350
351! convert pressure unit from pa to mb
352 do k = 1, lm
353 k1 = k + kd
354 k2 = k + lsk
355 do i = 1, im
356 plvl(i,k1+kb) = prsi(i,k2+kb) * 0.01 ! pa to mb (hpa)
357 plyr(i,k1) = prsl(i,k2) * 0.01 ! pa to mb (hpa)
358 tlyr(i,k1) = tgrs(i,k2)
359 prslk1(i,k1) = prslk(i,k2)
360 rho(i,k1) = prsl(i,k2)/(con_rd*tlyr(i,k1))
361 orho(i,k1) = 1.0/rho(i,k1)
362
364 es = min( prsl(i,k2), fpvs( tgrs(i,k2) ) ) ! fpvs and prsl in pa
365 qs = max( qmin, con_eps * es / (prsl(i,k2) + epsm1*es) )
366 rhly(i,k1) = max( 0.0, min( 1.0, max(qmin, qgrs(i,k2,ntqv))/qs ) )
367 qstl(i,k1) = qs
368 enddo
369 enddo
370
372 do j = 2, ntrac
373 do k = 1, lm
374 k1 = k + kd
375 k2 = k + lsk
376 tracer1(:,k1,j) = max(0.0, qgrs(:,k2,j))
377 enddo
378 enddo
379!
380 if (top_at_1) then ! input data from toa to sfc
381 if (lsk > 0) then
382 k1 = 1 + kd
383 k2 = k1 + kb
384 do i = 1, im
385 plvl(i,k2) = 0.01 * prsi(i,1+kb) ! pa to mb (hpa)
386 plyr(i,k1) = 0.5 * (plvl(i,k2+1) + plvl(i,k2))
387 prslk1(i,k1) = (plyr(i,k1)*0.001) ** rocp
388 enddo
389 else
390 k1 = 1 + kd
391 do i = 1, im
392 plvl(i,k1) = prsi(i,1) * 0.01 ! pa to mb (hpa)
393 enddo
394 endif
395 else ! input data from sfc to top
396 if (levs > lm) then
397 k1 = lm + kd
398 do i = 1, im
399 plvl(i,k1+1) = 0.01 * prsi(i,levs+1) ! pa to mb (hpa)
400 plyr(i,k1) = 0.5 * (plvl(i,k1+1) + plvl(i,k1))
401 prslk1(i,k1) = (plyr(i,k1)*0.001) ** rocp
402 enddo
403 else
404 k1 = lp1 + kd
405 do i = 1, im
406 plvl(i,k1) = prsi(i,lp1) * 0.01 ! pa to mb (hpa)
407 enddo
408 endif
409 endif
410!
411 if ( lextop ) then ! values for extra top layer
412 do i = 1, im
413 plvl(i,llb) = prsmin
414 if ( plvl(i,lla) <= prsmin ) plvl(i,lla) = 2.0*prsmin
415 plyr(i,lyb) = 0.5 * plvl(i,lla)
416 tlyr(i,lyb) = tlyr(i,lya)
417 prslk1(i,lyb) = (plyr(i,lyb)*0.001) ** rocp ! plyr in hPa
418 rho(i,lyb) = plyr(i,lyb) *100.0/(con_rd*tlyr(i,lyb))
419 orho(i,lyb) = 1.0/rho(i,lyb)
420 rhly(i,lyb) = rhly(i,lya)
421 qstl(i,lyb) = qstl(i,lya)
422 enddo
423
424! --- note: may need to take care the top layer amount
425 tracer1(:,lyb,:) = tracer1(:,lya,:)
426 endif
427
428
430
431 if (ntoz > 0) then ! interactive ozone generation
432 do k=1,lmk
433 do i=1,im
434 olyr(i,k) = max( qmin, tracer1(i,k,ntoz) )
435 enddo
436 enddo
437 else ! climatological ozone
438 call ozphys%run_o3clim(xlat, prslk1, con_pi, olyr)
439 endif ! end_if_ntoz
440
442 if (lsswr) then
443 call coszmn (xlon,sinlat,coslat,solhr,im,me, & ! --- inputs
444 coszen, coszdg) ! --- outputs
445 endif
446
449! - gasvmr(:,:,1) - co2 volume mixing ratio
450! - gasvmr(:,:,2) - n2o volume mixing ratio
451! - gasvmr(:,:,3) - ch4 volume mixing ratio
452! - gasvmr(:,:,4) - o2 volume mixing ratio
453! - gasvmr(:,:,5) - co volume mixing ratio
454! - gasvmr(:,:,6) - cf11 volume mixing ratio
455! - gasvmr(:,:,7) - cf12 volume mixing ratio
456! - gasvmr(:,:,8) - cf22 volume mixing ratio
457! - gasvmr(:,:,9) - ccl4 volume mixing ratio
458! - gasvmr(:,:,10) - cfc113 volumne mixing ratio
459
460! --- ... set up non-prognostic gas volume mixing ratioes
461
462 call getgases (plvl, xlon, xlat, im, lmk, ico2, top_at_1,& ! --- inputs
463 con_pi, gasvmr) ! --- outputs
464
465!CCPP: re-assign gasvmr(:,:,NF_VGAS) to gasvmr_X(:,:)
466 do k = 1, lmk
467 do i = 1, im
468 gasvmr_co2(i,k) = gasvmr(i,k,1)
469 gasvmr_n2o(i,k) = gasvmr(i,k,2)
470 gasvmr_ch4(i,k) = gasvmr(i,k,3)
471 gasvmr_o2(i,k) = gasvmr(i,k,4)
472 gasvmr_co(i,k) = gasvmr(i,k,5)
473 gasvmr_cfc11(i,k) = gasvmr(i,k,6)
474 gasvmr_cfc12(i,k) = gasvmr(i,k,7)
475 gasvmr_cfc22(i,k) = gasvmr(i,k,8)
476 gasvmr_ccl4(i,k) = gasvmr(i,k,9)
477 gasvmr_cfc113(i,k) = gasvmr(i,k,10)
478 enddo
479 enddo
480
482 do k = 2, lmk
483 do i = 1, im
484 tem2da(i,k) = log( plyr(i,k) )
485 tem2db(i,k) = log( plvl(i,k) )
486 enddo
487 enddo
488
489 if (top_at_1) then ! input data from toa to sfc
490 do i = 1, im
491 tem1d(i) = qme6
492 tem2da(i,1) = log( plyr(i,1) )
493 tem2db(i,1) = log( max(prsmin, plvl(i,1)) )
494 tem2db(i,lmp) = log( plvl(i,lmp) )
495 tsfa(i) = tlyr(i,lmk) ! sfc layer air temp
496 tlvl(i,1) = tlyr(i,1)
497 tlvl(i,lmp) = tskn(i)
498 enddo
499
500 do k = 1, lm
501 k1 = k + kd
502 do i = 1, im
503 qlyr(i,k1) = max( tem1d(i), qgrs(i,k,ntqv) )
504 tem1d(i) = min( qme5, qlyr(i,k1) )
505 tvly(i,k1) = tgrs(i,k) * (1.0 + fvirt*qlyr(i,k1)) ! virtual T (K)
506 delp(i,k1) = plvl(i,k1+1) - plvl(i,k1)
507 enddo
508 enddo
509
510 if ( lextop ) then
511 do i = 1, im
512 qlyr(i,lyb) = qlyr(i,lya)
513 tvly(i,lyb) = tvly(i,lya)
514 delp(i,lyb) = plvl(i,lla) - plvl(i,llb)
515 enddo
516 endif
517
518 do k = 2, lmk
519 do i = 1, im
520 tlvl(i,k) = tlyr(i,k) + (tlyr(i,k-1) - tlyr(i,k)) &
521 & * (tem2db(i,k) - tem2da(i,k)) &
522 & / (tem2da(i,k-1) - tem2da(i,k))
523 enddo
524 enddo
525
526! --- ... level height and layer thickness (km)
527! dz: Layer thickness between layer boundaries
528! dzb: Layer thickness between layer centers (lowest is from surface to lowest layer center)
529! hz: Height of each level (i.e. layer boundary)
530! hzb: Height of each layer center
531
532 tem0d = 0.001 * rog
533 do i = 1, im
534 do k = 1, lmk
535 dz(i,k) = tem0d * (tem2db(i,k+1) - tem2db(i,k)) * tvly(i,k)
536 enddo
537
538 hz(i,lmp) = 0.0
539 do k = lmk, 1, -1
540 hz(i,k) = hz(i,k+1) + dz(i,k)
541 enddo
542
543 do k = lmk, 1, -1
544 pfac = (tem2db(i,k+1) - tem2da(i,k)) / (tem2db(i,k+1) - tem2db(i,k))
545 hzb(i,k) = hz(i,k+1) + pfac * (hz(i,k) - hz(i,k+1))
546 enddo
547
548 do k = lmk-1, 1, -1
549 dzb(i,k) = hzb(i,k) - hzb(i,k+1)
550 enddo
551 dzb(i,lmk) = hzb(i,lmk) - hz(i,lmp)
552 enddo
553
554 else ! input data from sfc to toa
555
556 do i = 1, im
557 tem1d(i) = qme6
558 tem2da(i,1) = log( plyr(i,1) )
559 tem2db(i,1) = log( plvl(i,1) )
560 tem2db(i,lmp) = log( max(prsmin, plvl(i,lmp)) )
561 tsfa(i) = tlyr(i,1) ! sfc layer air temp
562 tlvl(i,1) = tskn(i)
563 tlvl(i,lmp) = tlyr(i,lmk)
564 enddo
565
566 do k = lm, 1, -1
567 do i = 1, im
568 qlyr(i,k) = max( tem1d(i), qgrs(i,k,ntqv) )
569 tem1d(i) = min( qme5, qlyr(i,k) )
570 tvly(i,k) = tgrs(i,k) * (1.0 + fvirt*qlyr(i,k)) ! virtual T (K)
571 delp(i,k) = plvl(i,k) - plvl(i,k+1)
572 enddo
573 enddo
574
575 if ( lextop ) then
576 do i = 1, im
577 qlyr(i,lyb) = qlyr(i,lya)
578 tvly(i,lyb) = tvly(i,lya)
579 delp(i,lyb) = plvl(i,lla) - plvl(i,llb)
580 enddo
581 endif
582
583 do k = 1, lmk-1
584 do i = 1, im
585 tlvl(i,k+1) = tlyr(i,k) + (tlyr(i,k+1) - tlyr(i,k)) &
586 & * (tem2db(i,k+1) - tem2da(i,k)) &
587 & / (tem2da(i,k+1) - tem2da(i,k))
588 enddo
589 enddo
590
591! --- ... level height and layer thickness (km)
592! dz: Layer thickness between layer boundaries
593! dzb: Layer thickness between layer centers (lowest is from surface to lowest layer center)
594! hz: Height of each level (i.e. layer boundary)
595! hzb: Height of each layer center
596
597 tem0d = 0.001 * rog
598 do i = 1, im
599 do k = lmk, 1, -1
600 dz(i,k) = tem0d * (tem2db(i,k) - tem2db(i,k+1)) * tvly(i,k)
601 enddo
602
603 hz(i,1) = 0.0
604 do k = 1, lmk
605 hz(i,k+1) = hz(i,k) + dz(i,k)
606 enddo
607
608 do k = 1, lmk
609 pfac = (tem2db(i,k) - tem2da(i,k)) / (tem2db(i,k) - tem2db(i,k+1))
610 hzb(i,k) = hz(i,k) + pfac * (hz(i,k+1) - hz(i,k))
611 enddo
612
613 do k = 2, lmk
614 dzb(i,k) = hzb(i,k) - hzb(i,k-1)
615 enddo
616 dzb(i,1) = hzb(i,1) - hz(i,1)
617 enddo
618
619 endif ! end_if_top_at_1
620
621
622!check print *,' in grrad : calling setaer '
623
625 if (ntchm>0 .and. iaermdl==2) then
626 do k=1,levs
627 do i=1,im
628 aer_nm(i,k,1) = qgrs(i,k,ntdu1)*1.e-9_kind_phys
629 aer_nm(i,k,2) = qgrs(i,k,ntdu2)*1.e-9_kind_phys
630 aer_nm(i,k,3) = qgrs(i,k,ntdu3)*1.e-9_kind_phys
631 aer_nm(i,k,4) = qgrs(i,k,ntdu4)*1.e-9_kind_phys
632 aer_nm(i,k,5) = qgrs(i,k,ntdu5)*1.e-9_kind_phys
633 aer_nm(i,k,6) = qgrs(i,k,ntss1)*1.e-9_kind_phys
634 aer_nm(i,k,7) = qgrs(i,k,ntss2)*1.e-9_kind_phys
635 aer_nm(i,k,8) = qgrs(i,k,ntss3)*1.e-9_kind_phys
636 aer_nm(i,k,9) = qgrs(i,k,ntss4)*1.e-9_kind_phys
637 aer_nm(i,k,10) = qgrs(i,k,ntss5)*1.e-9_kind_phys
638 aer_nm(i,k,11) = qgrs(i,k,ntsu)*1.e-9_kind_phys
639 aer_nm(i,k,12) = qgrs(i,k,ntbcb)*1.e-9_kind_phys
640 aer_nm(i,k,13) = qgrs(i,k,ntbcl)*1.e-9_kind_phys
641 aer_nm(i,k,14) = qgrs(i,k,ntocb)*1.e-9_kind_phys
642 aer_nm(i,k,15) = qgrs(i,k,ntocl)*1.e-9_kind_phys
643 enddo
644 enddo
645 endif
646
648 if (rrfs_sd .and. aero_dir_fdb) then
649 do k=1,lmk
650 do i=1,im
651 aer_nm(i,k,1 )=aer_nm(i,k,1 )+ qgrs(i,k,ntdust)*fdb_coef(1)*1.e-9 ! dust bin1
652 aer_nm(i,k,2 )=aer_nm(i,k,2 )+(qgrs(i,k,ntdust)*fdb_coef(2) &
653 +qgrs(i,k,ntcoarsepm)*fdb_coef(3))*1.e-9 ! dust bin2
654 aer_nm(i,k,3 )=aer_nm(i,k,3 )+qgrs(i,k,ntcoarsepm)*fdb_coef(4)*1.e-9 ! dust bin3
655 aer_nm(i,k,4 )=aer_nm(i,k,4 )+qgrs(i,k,ntcoarsepm)*fdb_coef(5)*1.e-9 ! dust bin4
656 aer_nm(i,k,12)=aer_nm(i,k,12)+qgrs(i,k,ntsmoke)*fdb_coef(6)*1.e-9 ! Smoke BC
657 aer_nm(i,k,14)=aer_nm(i,k,14)+qgrs(i,k,ntsmoke)*fdb_coef(7)*1.e-9 ! Smoke OA
658 enddo
659 enddo
660 endif
661
662
665 call setaer (plvl, plyr, prslk1, tvly, rhly, slmsk, & ! --- inputs
666 tracer1, aer_nm, xlon, xlat, im, lmk, lmp,&
667 lsswr, lslwr, iaermdl, iaerflg, top_at_1, con_pi, &
668 con_rd, con_g, faersw, faerlw, aerodp, ext550, errflg, errmsg) ! --- outputs
669
670! CCPP
671 do j = 1,nbdsw
672 do k = 1, lmk
673 do i = 1, im
674 ! NF_AESW = 3
675 faersw1(i,k,j) = faersw(i,k,j,1)
676 faersw2(i,k,j) = faersw(i,k,j,2)
677 faersw3(i,k,j) = faersw(i,k,j,3)
678 enddo
679 enddo
680 enddo
681
682 do j = 1,nbdlw
683 do k = 1, lmk
684 do i = 1, im
685 ! NF_AELW = 3
686 faerlw1(i,k,j) = faerlw(i,k,j,1)
687 faerlw2(i,k,j) = faerlw(i,k,j,2)
688 faerlw3(i,k,j) = faerlw(i,k,j,3)
689 enddo
690 enddo
691 enddo
692
695
696! --- ... obtain cloud information for radiation calculations
697
698! if (ntcw > 0) then ! prognostic cloud schemes
699 ccnd = 0.0_kind_phys
700 if (ncnd == 1) then ! Zhao_Carr_Sundqvist
701 do k=1,lmk
702 do i=1,im
703 ccnd(i,k,1) = tracer1(i,k,ntcw) ! liquid water/ice
704 enddo
705 enddo
706 elseif (ncnd == 2) then ! MG
707 do k=1,lmk
708 do i=1,im
709 ccnd(i,k,1) = tracer1(i,k,ntcw) ! liquid water
710 ccnd(i,k,2) = tracer1(i,k,ntiw) ! ice water
711 enddo
712 enddo
713 elseif (ncnd == 4) then ! MG2
714 do k=1,lmk
715 do i=1,im
716 ccnd(i,k,1) = tracer1(i,k,ntcw) ! liquid water
717 ccnd(i,k,2) = tracer1(i,k,ntiw) ! ice water
718 ccnd(i,k,3) = tracer1(i,k,ntrw) ! rain water
719 ccnd(i,k,4) = tracer1(i,k,ntsw) ! snow water
720 enddo
721 enddo
722 elseif (ncnd == 5 .or. ncnd == 6) then ! GFDL MP, Thompson, MG3, NSSL
723 do k=1,lmk
724 do i=1,im
725 ccnd(i,k,1) = tracer1(i,k,ntcw) ! liquid water
726 ccnd(i,k,2) = tracer1(i,k,ntiw) ! ice water
727 ccnd(i,k,3) = tracer1(i,k,ntrw) ! rain water
728 if (imp_physics == imp_physics_fer_hires ) then
729 ccnd(i,k,4) = 0.0
730 else
731 IF ( ncnd == 5 ) THEN
732 ccnd(i,k,4) = tracer1(i,k,ntsw) + tracer1(i,k,ntgl) ! snow + graupel
733 ELSEIF ( ncnd == 6 ) THEN
734 ccnd(i,k,4) = tracer1(i,k,ntsw) + tracer1(i,k,ntgl) + tracer1(i,k,nthl) ! snow + graupel + hail
735 ENDIF
736 endif
737 enddo
738 enddo
739 ! for Thompson MP - prepare variables for calc_effr
740 if_thompson: if (imp_physics == imp_physics_thompson .and. (ltaerosol .or. mraerosol)) then
741 do k=1,lmk
742 do i=1,im
743 qvs = qlyr(i,k)
744 qv_mp(i,k) = qvs/(1.-qvs)
745 rho(i,k) = con_eps*plyr(i,k)*100./(con_rd*tlyr(i,k)*(qv_mp(i,k)+con_eps))
746 orho(i,k) = 1.0/rho(i,k)
747 qc_mp(i,k) = tracer1(i,k,ntcw)/(1.-qvs)
748 qi_mp(i,k) = tracer1(i,k,ntiw)/(1.-qvs)
749 qs_mp(i,k) = tracer1(i,k,ntsw)/(1.-qvs)
750 nc_mp(i,k) = tracer1(i,k,ntlnc)/(1.-qvs)
751 ni_mp(i,k) = tracer1(i,k,ntinc)/(1.-qvs)
752 nwfa(i,k) = tracer1(i,k,ntwa)
753 enddo
754 enddo
755 elseif (imp_physics == imp_physics_thompson) then
756 do k=1,lmk
757 do i=1,im
758 qvs = qlyr(i,k)
759 qv_mp(i,k) = qvs/(1.-qvs)
760 rho(i,k) = con_eps*plyr(i,k)*100./(con_rd*tlyr(i,k)*(qv_mp(i,k)+con_eps))
761 orho(i,k) = 1.0/rho(i,k)
762 qc_mp(i,k) = tracer1(i,k,ntcw)/(1.-qvs)
763 qi_mp(i,k) = tracer1(i,k,ntiw)/(1.-qvs)
764 qs_mp(i,k) = tracer1(i,k,ntsw)/(1.-qvs)
765 if(nint(slmsk(i)) == 1) then
766 nc_mp(i,k) = nt_c_l*orho(i,k)
767 else
768 nc_mp(i,k) = nt_c_o*orho(i,k)
769 endif
770 ni_mp(i,k) = tracer1(i,k,ntinc)/(1.-qvs)
771 enddo
772 enddo
773 endif if_thompson
774 endif
775 do n=1,ncndl
776 do k=1,lmk
777 do i=1,im
778 if (ccnd(i,k,n) < epsq) ccnd(i,k,n) = 0.0
779 enddo
780 enddo
781 enddo
782 if (imp_physics == imp_physics_gfdl ) then
783 if (.not. lgfdlmprad) then
784
785
786! rsun the summation methods and order make the difference in calculation
787
788! clw(:,:) = clw(:,:) + tracer1(:,1:LMK,ntcw) &
789! + tracer1(:,1:LMK,ntiw) &
790! + tracer1(:,1:LMK,ntrw) &
791! + tracer1(:,1:LMK,ntsw) &
792! + tracer1(:,1:LMK,ntgl)
793 ccnd(:,:,1) = tracer1(:,1:lmk,ntcw)
794 ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:lmk,ntrw)
795 ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:lmk,ntiw)
796 ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:lmk,ntsw)
797 ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:lmk,ntgl)
798 endif
799 do k=1,lmk
800 do i=1,im
801 if (ccnd(i,k,1) < epsq ) ccnd(i,k,1) = 0.0
802 enddo
803 enddo
804 endif
805!
806 if (uni_cld) then
807 if (effr_in) then
808 do k=1,lm
809 k1 = k + kd
810 do i=1,im
811 cldcov(i,k1) = mg_cld(i,k)
812 effrl(i,k1) = effrl_inout(i,k)
813 effri(i,k1) = effri_inout(i,k)
814 effrr(i,k1) = effrr_in(i,k)
815 effrs(i,k1) = effrs_inout(i,k)
816 enddo
817 enddo
818 else
819 do k=1,lm
820 k1 = k + kd
821 do i=1,im
822 cldcov(i,k1) = mg_cld(i,k)
823 enddo
824 enddo
825 endif
826 elseif (imp_physics == imp_physics_gfdl) then ! GFDL MP
827 if ((imfdeepcnv==imfdeepcnv_gf .or. imfdeepcnv==imfdeepcnv_c3) .and. kdt>1) then
828 do k=1,lm
829 k1 = k + kd
830 do i=1,im
831 if (qci_conv(i,k)>0.) then
832 ! GF sub-grid cloud fraction
833 cldcov(i,k1) = clouds1(i,k1)
834 else
835 cldcov(i,k1) = tracer1(i,k1,ntclamt)
836 endif
837 enddo
838 enddo
839 else
840 ! GFDL cloud fraction
841 cldcov(1:im,1+kd:lm+kd) = tracer1(1:im,1:lm,ntclamt)
842 endif
843 if(effr_in) then
844 do k=1,lm
845 k1 = k + kd
846 do i=1,im
847 effrl(i,k1) = effrl_inout(i,k)
848 effri(i,k1) = effri_inout(i,k)
849 effrr(i,k1) = effrr_in(i,k)
850 effrs(i,k1) = effrs_inout(i,k)
851! if(me==0) then
852! if(effrl(i,k1)> 5.0) then
853! write(6,*) 'rad driver:cloud radii:',kdt, i,k1, &
854! effrl(i,k1)
855! endif
856! if(effrs(i,k1)==0.0) then
857! write(6,*) 'rad driver:snow mixing ratio:',kdt, i,k1, &
858! tracer1(i,k,ntsw)
859! endif
860! endif
861 enddo
862 enddo
863 endif
864
865 elseif (imp_physics == imp_physics_nssl ) then ! NSSL MP
866 cldcov = 0.0
867 if(effr_in) then
868 do k=1,lm
869 k1 = k + kd
870 do i=1,im
871 effrl(i,k1) = effrl_inout(i,k)! re_cloud (i,k)
872 effri(i,k1) = effri_inout(i,k)! re_ice (i,k)
873 effrr(i,k1) = effrr_in(i,k)
874 effrs(i,k1) = effrs_inout(i,k) ! re_snow(i,k)
875 enddo
876 enddo
877 else
878 ! not used yet -- effr_in should always be true for now
879 endif
880
881 elseif (imp_physics == imp_physics_thompson) then ! Thompson MP
882 !
883 ! Compute effective radii for QC, QI, QS with (GF, MYNN) or without (all others) sub-grid clouds
884 !
885 ! Update number concentration, consistent with sub-grid clouds (GF, MYNN) or without (all others)
886 do k=1,lm
887 do i=1,im
888 if ((ltaerosol .or. mraerosol) .and. qc_mp(i,k)>1.e-12 .and. nc_mp(i,k)<100.) then
889 nc_mp(i,k) = make_dropletnumber(qc_mp(i,k)*rho(i,k), nwfa(i,k)*rho(i,k)) * orho(i,k)
890 endif
891 if (qi_mp(i,k)>1.e-12 .and. ni_mp(i,k)<100.) then
892 ni_mp(i,k) = make_icenumber(qi_mp(i,k)*rho(i,k), tlyr(i,k)) * orho(i,k)
893 endif
894 end do
895 end do
897 do i=1,im
898 islmsk = nint(slmsk(i))
899 ! Effective radii [m] are now intent(out), bounds applied in calc_effectRad
900 !tgs: progclduni has different limits for ice radii (10.0-150.0) than
901 ! calc_effectRad (4.99-125.0 for WRFv3.8.1; 2.49-125.0 for WRFv4+)
902 ! it will raise the low limit from 5 to 10, but the high limit will remain 125.
903 call calc_effectrad (tlyr(i,:), plyr(i,:)*100., qv_mp(i,:), qc_mp(i,:), &
904 nc_mp(i,:), qi_mp(i,:), ni_mp(i,:), qs_mp(i,:), &
905 effrl(i,:), effri(i,:), effrs(i,:), islmsk, 1, lm )
906 ! Scale Thompson's effective radii from meter to micron
907 do k=1,lm
908 effrl(i,k) = max(re_qc_min, min(effrl(i,k), re_qc_max))*1.e6
909 effri(i,k) = max(re_qi_min, min(effri(i,k), re_qi_max))*1.e6
910 effrs(i,k) = max(re_qs_min, min(effrs(i,k), re_qs_max))*1.e6
911 end do
912 effrl(i,lmk) = re_qc_min*1.e6
913 effri(i,lmk) = re_qi_min*1.e6
914 effrs(i,lmk) = re_qs_min*1.e6
915 end do
916 effrr(:,:) = 1000. ! rrain_def=1000.
917 ! Update global arrays
918 do k=1,lm
919 k1 = k + kd
920 do i=1,im
921 effrl_inout(i,k) = effrl(i,k1)
922 effri_inout(i,k) = effri(i,k1)
923 effrs_inout(i,k) = effrs(i,k1)
924 enddo
925 enddo
926 else ! all other cases
927 cldcov = 0.0
928 endif
929
930!
931! --- add suspended convective cloud water to grid-scale cloud water
932! only for cloud fraction & radiation computation
933! it is to enhance cloudiness due to suspended convec cloud water
934! for zhao/moorthi's (imp_phys=99) &
935! ferrier's (imp_phys=5) microphysics schemes
936
937 if ((num_p3d == 4) .and. (npdf3d == 3)) then ! same as imp_physics = imp_physics_zhao_carr_pdf
938 do k=1,lm
939 k1 = k + kd
940 do i=1,im
941 !GJF: this is not consistent with GFS_typedefs,
942 ! but it looks like the Zhao-Carr-PDF scheme is not in the CCPP
943 deltaq(i,k1) = 0.0!Tbd%phy_f3d(i,k,5) !GJF: this variable is not in phy_f3d anymore
944 cnvw(i,k1) = cnvw_in(i,k)
945 cnvc(i,k1) = cnvc_in(i,k)
946 enddo
947 enddo
948 elseif ((npdf3d == 0) .and. (ncnvcld3d == 1)) then ! all other microphysics with pdfcld = .false. and cnvcld = .true.
949 do k=1,lm
950 k1 = k + kd
951 do i=1,im
952 deltaq(i,k1) = 0.0
953 cnvw(i,k1) = cnvw_in(i,k)
954 cnvc(i,k1) = 0.0
955 enddo
956 enddo
957 else ! all the rest
958 do k=1,lmk
959 do i=1,im
960 deltaq(i,k) = 0.0
961 cnvw(i,k) = 0.0
962 cnvc(i,k) = 0.0
963 enddo
964 enddo
965 endif
966
967 if (imp_physics == imp_physics_zhao_carr) then
968 ccnd(1:im,1:lmk,1) = ccnd(1:im,1:lmk,1) + cnvw(1:im,1:lmk)
969 endif
970
973 & ( plyr, plvl, tlyr, tvly, qlyr, qstl, rhly, & ! --- inputs:
974 & ccnd, ncndl, cnvw, cnvc, tracer1, &
975 & xlat, xlon, slmsk, dz, delp, im, lm, lmk, lmp, &
976 & deltaq, sup, dcorr_con, me, icloud, kdt, &
977 & ntrac, ntcw, ntiw, ntrw, ntsw, ntgl, ntclamt, &
978 & imp_physics, imp_physics_nssl, imp_physics_fer_hires, &
979 & imp_physics_gfdl, imp_physics_thompson, imp_physics_wsm6, &
980 & imp_physics_zhao_carr, imp_physics_zhao_carr_pdf, &
981 & imp_physics_mg, iovr, iovr_rand, iovr_maxrand, iovr_max, &
982 & iovr_dcorr, iovr_exp, iovr_exprand, idcor, idcor_con, &
983 & idcor_hogan, idcor_oreopoulos, lcrick, lcnorm, &
984 & imfdeepcnv, imfdeepcnv_gf, imfdeepcnv_c3, do_mynnedmf, &
985 & lgfdlmprad, xr_cnvcld, &
986 & uni_cld, lmfshal, lmfdeep2, cldcov, clouds1, &
987 & effrl, effri, effrr, effrs, effr_in, &
988 & effrl_inout, effri_inout, effrs_inout, &
989 & lwp_ex, iwp_ex, lwp_fc, iwp_fc, &
990 & dzb, xlat_d, julian, yearlen, gridkm, top_at_1, si, &
991 & con_ttp, con_pi, con_g, con_rd, con_thgni, &
992 & cld_frac, cld_lwp, cld_reliq, cld_iwp, cld_reice, & ! --- outputs:
993 & cld_rwp, cld_rerain, cld_swp, cld_resnow, & ! --- outputs:
994 & cldsa, mtopa, mbota, de_lgth, alpha & ! --- outputs:
995 & )
996
997! endif ! end_if_ntcw
998
1000 if (pert_clds) then
1001 do i=1,im
1002 tmp_wt= -1*log( ( 2.0 / ( sppt_wts(i,38) ) ) - 1 )
1003 call cdfnor(tmp_wt,cdfz)
1004 cldp1d(i) = cdfz
1005 enddo
1006 do k = 1, lmk
1007 do i = 1, im
1008 ! compute beta distribution parameters
1009 m = cld_frac(i,k)
1010 if (m<0.99 .AND. m > 0.01) then
1011 s = sppt_amp*m*(1.-m)
1012 alpha0 = m*m*(1.-m)/(s*s)-m
1013 beta0 = alpha0*(1.-m)/m
1014 ! compute beta distribution value corresponding
1015 ! to the given percentile albPpert to use as new albedo
1016 call ppfbet(cldp1d(i),alpha0,beta0,iflag,cldtmp)
1017 cld_frac(i,k) = cldtmp
1018 else
1019 cld_frac(i,k) = m
1020 endif
1021 enddo ! end_do_i_loop
1022 enddo ! end_do_k_loop
1023 endif
1024 do k = 1, lm
1025 do i = 1, im
1026 clouds1(i,k) = cld_frac(i,k)
1027 clouds2(i,k) = cld_lwp(i,k)
1028 clouds3(i,k) = cld_reliq(i,k)
1029 clouds4(i,k) = cld_iwp(i,k)
1030 clouds5(i,k) = cld_reice(i,k)
1031 clouds6(i,k) = cld_rwp(i,k)
1032 clouds7(i,k) = cld_rerain(i,k)
1033 clouds8(i,k) = cld_swp(i,k)
1034 clouds9(i,k) = cld_resnow(i,k)
1035 cldfra(i,k) = cld_frac(i,k)
1036 enddo
1037 enddo
1038 do i = 1, im
1039 cldfra2d(i) = 0.0
1040 do k = 1, lm-1
1041 cldfra2d(i) = max(cldfra2d(i), cldfra(i,k))
1042 enddo
1043 enddo
1044
1045 if ( spp_rad == 1 ) then
1046 do k=1,lm
1047 if (k < levs) then
1048 do i=1,im
1049 clouds3(i,k) = clouds3(i,k) - spp_wts_rad(i,k) * clouds3(i,k)
1050 clouds5(i,k) = clouds5(i,k) - spp_wts_rad(i,k) * clouds5(i,k)
1051 clouds9(i,k) = clouds9(i,k) - spp_wts_rad(i,k) * clouds9(i,k)
1052 enddo
1053 else
1054 do i=1,im
1055 clouds3(i,k) = clouds3(i,k) - spp_wts_rad(i,levs) * clouds3(i,k)
1056 clouds5(i,k) = clouds5(i,k) - spp_wts_rad(i,levs) * clouds5(i,k)
1057 clouds9(i,k) = clouds9(i,k) - spp_wts_rad(i,levs) * clouds9(i,k)
1058 enddo
1059 endif
1060 enddo
1061 endif
1062
1063! mg, sfc-perts
1064! --- scale random patterns for surface perturbations with
1065! perturbation size
1066! --- turn vegetation fraction pattern into percentile pattern
1068 alb1d(:) = 0.
1069 if (lndp_type==1) then
1070 do k =1,n_var_lndp
1071 if (lndp_var_list(k) == 'alb') then
1072 do i=1,im
1073 call cdfnor(sfc_wts(i,k),alb1d(i))
1074 !lndp_alb = lndp_prt_list(k)
1075 enddo
1076 endif
1077 enddo
1078 endif
1079! mg, sfc-perts
1080
1081 end subroutine gfs_rrtmg_pre_run
1083 end module gfs_rrtmg_pre
subroutine, public gfs_rrtmg_pre_run(im, levs, lm, lmk, lmp, n_var_lndp, lextop, ltp, imfdeepcnv, imfdeepcnv_gf, imfdeepcnv_c3, me, ncnd, ntrac, num_p3d, npdf3d, xr_cnvcld, ncnvcld3d, ntqv, ntcw, ntiw, ntlnc, ntinc, ntrnc, ntsnc, ntccn, top_at_1, ntrw, ntsw, ntgl, nthl, ntwa, ntoz, ntsmoke, ntdust, ntcoarsepm, ntclamt, nleffr, nieffr, nseffr, lndp_type, kdt, ntdu1, ntdu2, ntdu3, ntdu4, ntdu5, ntss1, ntss2, ntss3, ntss4, ntss5, ntsu, ntbcb, ntbcl, ntocb, ntocl, ntchm, imp_physics, imp_physics_nssl, nssl_ccn_on, nssl_invertccn, imp_physics_thompson, imp_physics_gfdl, imp_physics_zhao_carr, imp_physics_zhao_carr_pdf, imp_physics_mg, imp_physics_wsm6, imp_physics_fer_hires, iovr, iovr_rand, iovr_maxrand, iovr_max, iovr_dcorr, iovr_exp, iovr_exprand, idcor, idcor_con, idcor_hogan, idcor_oreopoulos, dcorr_con, julian, yearlen, lndp_var_list, lsswr, lslwr, ltaerosol, mraerosol, lgfdlmprad, uni_cld, effr_in, do_mynnedmf, lmfshal, lcnorm, lmfdeep2, lcrick, fhswr, fhlwr, solhr, sup, con_eps, epsm1, fvirt, rog, rocp, con_rd, xlat_d, xlat, xlon, coslat, sinlat, tsfc, slmsk, prsi, prsl, prslk, tgrs, sfc_wts, mg_cld, effrr_in, pert_clds, sppt_wts, sppt_amp, cnvw_in, cnvc_in, qgrs, aer_nm, dx, icloud, iaermdl, iaerflg, con_pi, con_g, con_ttp, con_thgni, si, coszen, coszdg, effrl_inout, effri_inout, effrs_inout, clouds1, clouds2, clouds3, clouds4, clouds5, qci_conv, kd, kt, kb, mtopa, mbota, raddt, tsfg, tsfa, de_lgth, alb1d, delp, dz, plvl, plyr, tlvl, tlyr, qlyr, olyr, gasvmr_co2, gasvmr_n2o, gasvmr_ch4, gasvmr_o2, gasvmr_co, gasvmr_cfc11, gasvmr_cfc12, gasvmr_cfc22, gasvmr_ccl4, gasvmr_cfc113, aerodp, ext550, clouds6, clouds7, clouds8, clouds9, cldsa, cldfra, cldfra2d, lwp_ex, iwp_ex, lwp_fc, iwp_fc, faersw1, faersw2, faersw3, faerlw1, faerlw2, faerlw3, alpha, rrfs_sd, aero_dir_fdb, fdb_coef, spp_wts_rad, spp_rad, ico2, ozphys, errmsg, errflg)
subroutine calc_effectrad(t1d, p1d, qv1d, qc1d, nc1d, qi1d, ni1d, qs1d, re_qc1d, re_qi1d, re_qs1d, lsml, kts, kte)
Compute radiation effective radii of cloud water, ice, and snow. These are entirely consistent with m...
elemental real function, public make_icenumber(q_ice, temp)
Table of lookup values of radiative effective radius of ice crystals as a function of Temperature fro...
subroutine, public ppfbet(pr, p, q, iflag, x)
This subroutine computes the beta distribution value that matches the percentile from the random patt...
subroutine, public cdfnor(z, cdfz)
This subrtouine calculates the CDF of the standard normal distribution evaluated at z.
subroutine, public coszmn(xlon, sinlat, coslat, solhr, im, me, coszen, coszdg)
This subroutine computes mean cos solar zenith angle over SW calling interval.
subroutine, public radiation_clouds_prop(plyr, plvl, tlyr, tvly, qlyr, qstl, rhly, ccnd, ncndl, cnvw, cnvc, tracer1, xlat, xlon, slmsk, dz, delp, ix, lm, nlay, nlp1, deltaq, sup, dcorr_con, me, icloud, kdt, ntrac, ntcw, ntiw, ntrw, ntsw, ntgl, ntclamt, imp_physics, imp_physics_nssl, imp_physics_fer_hires, imp_physics_gfdl, imp_physics_thompson, imp_physics_wsm6, imp_physics_zhao_carr, imp_physics_zhao_carr_pdf, imp_physics_mg, iovr, iovr_rand, iovr_maxrand, iovr_max, iovr_dcorr, iovr_exp, iovr_exprand, idcor, idcor_con, idcor_hogan, idcor_oreopoulos, lcrick, lcnorm, imfdeepcnv, imfdeepcnv_gf, imfdeepcnv_c3, do_mynnedmf, lgfdlmprad, xr_cnvcld, uni_cld, lmfshal, lmfdeep2, cldcov, clouds1, effrl, effri, effrr, effrs, effr_in, effrl_inout, effri_inout, effrs_inout, lwp_ex, iwp_ex, lwp_fc, iwp_fc, dzlay, latdeg, julian, yearlen, gridkm, top_at_1, si, con_ttp, con_pi, con_g, con_rd, con_thgni, cld_frac, cld_lwp, cld_reliq, cld_iwp, cld_reice, cld_rwp, cld_rerain, cld_swp, cld_resnow, clds, mtop, mbot, de_lgth, alpha)
Subroutine radiation_clouds_prop computes cloud related quantities for different cloud microphysics s...
subroutine, public adjust_cloudice(cfr, qi, qs, qvs, t, dz, entr, k1, k2, kts, kte)
integer, parameter, public nf_clds
number of fields in cloud array
subroutine, public adjust_cloudh2o(cfr, qc, qvs, t, dz, entr, k1, k2, kts, kte)
subroutine, public adjust_cloudfinal(cfr, qc, qi, rho, dz, kts, kte)
Do not alter any grid-explicitly resolved hydrometeors, rather only the supposed amounts due to the c...
subroutine, public cal_cldfra3(cldfra, qv, qc, qi, qs, dz, p, t, xland, gridkm, modify_qvapor, max_relh, kts, kte, debug_flag)
Cloud fraction scheme by G. Thompson (NCAR-RAL), not intended for combining with any cumulus or shall...
subroutine, public find_cloudlayers(qvs1d, cfr1d, t1d, p1d, dz1d, entrmnt, debugfl, qc1d, qi1d, qs1d, kts, kte)
From cloud fraction array, find clouds of multi-level depth and compute a reasonable value of LWP or ...
subroutine, public getgases(plvl, xlon, xlat, imax, lmax, ico2flg, top_at_1, con_pi, gasdat)
This subroutine sets up global distribution of radiation absorbing gases in volume mixing ratio....
integer, parameter, public nf_vgas
This module provides an Application Program Interface (API) for computing basic thermodynamic physics...
Definition funcphys.f90:26
The operational GFS currently parameterizes ozone production and destruction based on monthly mean co...
This module sets up astronomy quantities for solar radiation calculations.
This module computes cloud related quantities for radiation computations.
This module sets up constant gas rofiles, such as co2, ch4, n2o, o2, and those of cfc gases.
This module contains LW band parameters set up.
Definition radlw_param.f:61
This module is for specifying the band structures and program parameters used by the RRTMG-SW scheme.
Definition radsw_param.f:62
This module contains some of the most frequently used math and physics constants for RRTMG.
Definition radcons.f90:7
Derived type containing data and procedures needed by ozone photochemistry parameterization Note All ...
define type construct for optional radiation flux profiles
Definition radlw_param.f:94
derived type for LW fluxes at surface
Definition radlw_param.f:87
derived type for LW fluxes at top of atmosphere
Definition radlw_param.f:78
derived type for SW fluxes' column profiles (at layer interfaces)
derived type for SW fluxes at surface
Definition radsw_param.f:89
derived type for SW fluxes at TOA
Definition radsw_param.f:79