CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_rrtmgp_post.F90
1
5 use machine, only: kind_phys
8 use mo_heating_rates, only: compute_heating_rate
9 use radiation_tools, only: check_error_msg
10 implicit none
11
12 public gfs_rrtmgp_post_run
13
14contains
15
25 subroutine gfs_rrtmgp_post_run (nCol, nLev, nDay, iSFC, iTOA, idxday, doLWrad, doSWrad, &
26 do_lw_clrsky_hr, do_sw_clrsky_hr, save_diag, fhlwr, fhswr, sfc_alb_nir_dir, &
27 sfc_alb_nir_dif, sfc_alb_uvvis_dir, sfc_alb_uvvis_dif, p_lev, tsfa, coszen, coszdg, &
28 fluxlwDOWN_clrsky, fluxlwUP_allsky, fluxlwDOWN_allsky, fluxlwUP_clrsky, &
29 fluxswDOWN_clrsky, fluxswUP_allsky, fluxswDOWN_allsky, fluxswUP_clrsky, &
30 raddt, aerodp, cldsa, mtopa, mbota, cld_frac, cldtaulw, cldtausw, scmpsw, fluxr, &
31 sfcdlw, sfculw, sfcflw, tsflw, htrlw, htrlwu, topflw, nirbmdi, nirdfdi, visbmdi, &
32 visdfdi, nirbmui, nirdfui, visbmui, visdfui, sfcnsw, sfcdsw, htrsw, sfcfsw, topfsw, &
33 htrswc, htrlwc, errmsg, errflg)
34
35 ! Inputs
36 integer, intent(in) :: &
37 ncol, & !< Horizontal loop extent
38 nlev, & !< Number of vertical layers
39 nday, & !< Number of daylit columns
40 isfc, & !< Vertical index for surface level
41 itoa
42 integer, intent(in), dimension(:) :: &
43 idxday
44 integer, intent(in), dimension(:,:) :: &
45 mbota, & !< Vertical indices for low, middle and high cloud tops
46 mtopa
47 logical, intent(in) :: &
48 dolwrad, & !< Logical flags for lw radiation calls
49 doswrad, & !< Logical flags for sw radiation calls
50 do_lw_clrsky_hr, & !< Output clear-sky LW heating-rate?
51 do_sw_clrsky_hr, & !< Output clear-sky SW heating-rate?
52 save_diag
53 real(kind_phys), intent(in) :: &
54 fhlwr, & !< Frequency for LW radiation calls
55 fhswr
56 real(kind_phys), dimension(:), intent(in) :: &
57 tsfa, & !< Lowest model layer air temperature for radiation (K)
58 coszen, & !< Cosine(SZA)
59 coszdg, & !< Cosine(SZA), daytime
60 sfc_alb_nir_dir, & !< Surface albedo (direct)
61 sfc_alb_nir_dif, & !< Surface albedo (diffuse)
62 sfc_alb_uvvis_dir, & !< Surface albedo (direct)
63 sfc_alb_uvvis_dif
64 real(kind_phys), dimension(:,:), intent(in), optional :: &
65 p_lev, & !< Pressure @ model layer-interfaces (Pa)
66 fluxlwup_allsky, & !< RRTMGP longwave all-sky flux (W/m2)
67 fluxlwdown_allsky, & !< RRTMGP longwave all-sky flux (W/m2)
68 fluxlwup_clrsky, & !< RRTMGP longwave clear-sky flux (W/m2)
69 fluxlwdown_clrsky, & !< RRTMGP longwave clear-sky flux (W/m2)
70 fluxswup_allsky, & !< RRTMGP shortwave all-sky flux (W/m2)
71 fluxswdown_allsky, & !< RRTMGP shortwave all-sky flux (W/m2)
72 fluxswup_clrsky, & !< RRTMGP shortwave clear-sky flux (W/m2)
73 fluxswdown_clrsky
74 real(kind_phys), intent(in) :: &
75 raddt
76 real(kind_phys), dimension(:,:), intent(in) :: &
77 aerodp, & !< Vertical integrated optical depth for various aerosol species
78 cldsa, & !< Fraction of clouds for low, middle, high, total and BL
79 cld_frac, & !< Total cloud fraction in each layer
80 cldtaulw, & !< approx 10.mu band layer cloud optical depth
81 cldtausw
82 type(cmpfsw_type), dimension(:), intent(in) :: &
83 scmpsw
90
91
92 real(kind=kind_phys), dimension(:,:), intent(inout) :: fluxr
93
94 ! Outputs (mandatory)
95 real(kind_phys), dimension(:), intent(inout) :: &
96 tsflw, & !< LW sfc air temp during calculation (K)
97 sfcdlw, & !< LW sfc all-sky downward flux (W/m2)
98 sfculw, & !< LW sfc all-sky upward flux (W/m2)
99 nirbmdi, & !< SW sfc nir beam downward flux (W/m2)
100 nirdfdi, & !< SW sfc nir diff downward flux (W/m2)
101 visbmdi, & !< SW sfc uv+vis beam downward flux (W/m2)
102 visdfdi, & !< SW sfc uv+vis diff downward flux (W/m2)
103 nirbmui, & !< SW sfc nir beam upward flux (W/m2)
104 nirdfui, & !< SW sfc nir diff upward flux (W/m2)
105 visbmui, & !< SW sfc uv+vis beam upward flux (W/m2)
106 visdfui, & !< SW sfc uv+vis diff upward flux (W/m2)
107 sfcnsw, & !< SW sfc all-sky net flux (W/m2) flux into ground
108 sfcdsw
109 real(kind_phys), dimension(:,:), intent(inout) :: &
110 htrlw, & !< LW all-sky heating rate (K/s)
111 htrsw
112 real(kind_phys), dimension(:,:), intent(inout), optional :: &
113 htrlwu
114 type(sfcflw_type), dimension(:), intent(inout) :: &
115 sfcflw
116 type(sfcfsw_type), dimension(:), intent(inout) :: &
117 sfcfsw
118 type(topfsw_type), dimension(:), intent(inout) :: &
119 topfsw
120 type(topflw_type), dimension(:), intent(inout) :: &
121 topflw
122 character(len=*), intent(out) :: &
123 errmsg
124 integer, intent(out) :: &
125 errflg
126
127 ! Outputs (optional)
128 real(kind_phys),dimension(:,:),intent(inout),optional :: &
129 htrlwc, & !< LW clear-sky heating-rate (K/s)
130 htrswc
131
132 ! Local variables
133 integer :: i, j, k, itop, ibtc
134 real(kind_phys) :: tem0d, tem1, tem2
135 real(kind_phys), dimension(nDay, nLev) :: thetatendclrsky, thetatendallsky
136
137 ! Initialize CCPP error handling variables
138 errmsg = ''
139 errflg = 0
140
141 if (.not. (dolwrad .or. doswrad)) return
142
143 if (dolwrad) then
144 ! #######################################################################################
145 ! Compute LW heating-rates.
146 ! #######################################################################################
147
148 ! Clear-sky heating-rate (optional)
149 if (do_lw_clrsky_hr) then
150 call check_error_msg('GFS_rrtmgp_post',compute_heating_rate( &
151 fluxlwup_clrsky, & ! IN - RRTMGP upward longwave clear-sky flux profiles (W/m2)
152 fluxlwdown_clrsky, & ! IN - RRTMGP downward longwave clear-sky flux profiles (W/m2)
153 p_lev, & ! IN - Pressure @ layer-interfaces (Pa)
154 htrlwc)) ! OUT - Longwave clear-sky heating rate (K/sec)
155 endif
156
157 ! All-sky heating-rate (mandatory)
158 call check_error_msg('GFS_rrtmgp_post',compute_heating_rate( &
159 fluxlwup_allsky, & ! IN - RRTMGP upward longwave all-sky flux profiles (W/m2)
160 fluxlwdown_allsky, & ! IN - RRTMGP downward longwave all-sky flux profiles (W/m2)
161 p_lev, & ! IN - Pressure @ layer-interfaces (Pa)
162 htrlw)) ! OUT - Longwave all-sky heating rate (K/sec)
163
164 ! #######################################################################################
165 ! Save LW outputs.
166 ! (Copy fluxes from RRTMGP types into model radiation types.)
167 ! #######################################################################################
168 ! TOA fluxes
169
170 topflw(:)%upfxc = fluxlwup_allsky(:,itoa)
171 topflw(:)%upfx0 = fluxlwup_clrsky(:,itoa)
172
173 ! Surface fluxes
174 sfcflw(:)%upfxc = fluxlwup_allsky(:,isfc)
175 sfcflw(:)%upfx0 = fluxlwup_clrsky(:,isfc)
176 sfcflw(:)%dnfxc = fluxlwdown_allsky(:,isfc)
177 sfcflw(:)%dnfx0 = fluxlwdown_clrsky(:,isfc)
178
179 ! Save surface air temp for diurnal adjustment at model t-steps
180 tsflw(:) = tsfa(:)
181
182 ! Radiation fluxes for other physics processes
183 sfcdlw(:) = sfcflw(:)%dnfxc
184 sfculw(:) = sfcflw(:)%upfxc
185
186 ! Heating-rate at radiation timestep, used for adjustment between radiation calls.
187 htrlwu = htrlw
188
189 ! #######################################################################################
190 ! Save LW diagnostics
191 ! - For time averaged output quantities (including total-sky and clear-sky SW and LW
192 ! fluxes at TOA and surface; conventional 3-domain cloud amount, cloud top and base
193 ! pressure, and cloud top temperature; aerosols AOD, etc.), store computed results in
194 ! corresponding slots of array fluxr with appropriate time weights.
195 ! - Collect the fluxr data for wrtsfc
196 ! #######################################################################################
197 if (save_diag) then
198 do i=1,ncol
199 ! LW all-sky fluxes
200 fluxr(i,1 ) = fluxr(i,1 ) + fhlwr * fluxlwup_allsky( i,itoa) ! total sky top lw up
201 fluxr(i,19) = fluxr(i,19) + fhlwr * fluxlwdown_allsky(i,isfc) ! total sky sfc lw dn
202 fluxr(i,20) = fluxr(i,20) + fhlwr * fluxlwup_allsky( i,isfc) ! total sky sfc lw up
203 ! LW clear-sky fluxes
204 fluxr(i,28) = fluxr(i,28) + fhlwr * fluxlwup_clrsky( i,itoa) ! clear sky top lw up
205 fluxr(i,30) = fluxr(i,30) + fhlwr * fluxlwdown_clrsky(i,isfc) ! clear sky sfc lw dn
206 fluxr(i,33) = fluxr(i,33) + fhlwr * fluxlwup_clrsky( i,isfc) ! clear sky sfc lw up
207 enddo
208
209 ! Save cld frac,toplyr,botlyr and top temp, note that the order of h,m,l cloud is reversed for
210 ! the fluxr output. save interface pressure (pa) of top/bot
211 do j = 1, 3
212 do i = 1, ncol
213 tem0d = raddt * cldsa(i,j)
214 itop = mtopa(i,j)
215 ibtc = mbota(i,j)
216
217 ! Add optical depth and emissivity output
218 tem2 = 0.
219 do k=ibtc,itop
220 tem2 = tem2 + cldtaulw(i,k) ! approx 10. mu channel
221 enddo
222 fluxr(i,46-j) = fluxr(i,46-j) + tem0d * (1.0-exp(-tem2))
223 enddo
224 enddo
225 endif
226 endif
227 ! #######################################################################################
228 ! #######################################################################################
229 ! #######################################################################################
230 ! #######################################################################################
231 ! #######################################################################################
232 ! #######################################################################################
233 if (doswrad) then
234 if (nday .gt. 0) then
235 ! #################################################################################
236 ! Compute SW heating-rates
237 ! #################################################################################
238
239 ! Clear-sky heating-rate (optional)
240 if (do_sw_clrsky_hr) then
241 htrswc(:,:) = 0._kind_phys
242 call check_error_msg('GFS_rrtmgp_post',compute_heating_rate( &
243 fluxswup_clrsky(idxday(1:nday),:), & ! IN - Shortwave upward clear-sky flux profiles (W/m2)
244 fluxswdown_clrsky(idxday(1:nday),:), & ! IN - Shortwave downward clear-sky flux profiles (W/m2)
245 p_lev(idxday(1:nday),:), & ! IN - Pressure at model-interface (Pa)
246 thetatendclrsky)) ! OUT - Clear-sky heating-rate (K/sec)
247 htrswc(idxday(1:nday),:)=thetatendclrsky !**NOTE** GP doesn't use radiation levels, it uses the model fields. Not sure if this is necessary
248 endif
249
250 ! All-sky heating-rate (mandatory)
251 htrsw(:,:) = 0._kind_phys
252 call check_error_msg('GFS_rrtmgp_post',compute_heating_rate( &
253 fluxswup_allsky(idxday(1:nday),:), & ! IN - Shortwave upward all-sky flux profiles (W/m2)
254 fluxswdown_allsky(idxday(1:nday),:), & ! IN - Shortwave downward all-sky flux profiles (W/m2)
255 p_lev(idxday(1:nday),:), & ! IN - Pressure at model-interface (Pa)
256 thetatendallsky)) ! OUT - All-sky heating-rate (K/sec)
257 htrsw(idxday(1:nday),:) = thetatendallsky
258
259 ! #################################################################################
260 ! Save SW outputs
261 ! (Copy fluxes from RRTMGP types into model radiation types.)
262 ! #################################################################################
263
264 ! TOA fluxes
265 topfsw(:)%upfxc = fluxswup_allsky(:,itoa)
266 topfsw(:)%upfx0 = fluxswup_clrsky(:,itoa)
267 topfsw(:)%dnfxc = fluxswdown_allsky(:,itoa)
268
269 ! Surface fluxes
270 sfcfsw(:)%upfxc = fluxswup_allsky(:,isfc)
271 sfcfsw(:)%upfx0 = fluxswup_clrsky(:,isfc)
272 sfcfsw(:)%dnfxc = fluxswdown_allsky(:,isfc)
273 sfcfsw(:)%dnfx0 = fluxswdown_clrsky(:,isfc)
274
275 ! Surface down and up spectral component fluxes
276 ! - Save two spectral bands' surface downward and upward fluxes for output.
277 do i=1,ncol
278 nirbmdi(i) = scmpsw(i)%nirbm
279 nirdfdi(i) = scmpsw(i)%nirdf
280 visbmdi(i) = scmpsw(i)%visbm
281 visdfdi(i) = scmpsw(i)%visdf
282 nirbmui(i) = scmpsw(i)%nirbm * sfc_alb_nir_dir(i)
283 nirdfui(i) = scmpsw(i)%nirdf * sfc_alb_nir_dif(i)
284 visbmui(i) = scmpsw(i)%visbm * sfc_alb_uvvis_dir(i)
285 visdfui(i) = scmpsw(i)%visdf * sfc_alb_uvvis_dif(i)
286 enddo
287 else ! if_nday_block
288 ! #################################################################################
289 ! Dark everywhere
290 ! #################################################################################
291 htrsw(:,:) = 0.0
292 sfcfsw = sfcfsw_type( 0.0, 0.0, 0.0, 0.0 )
293 topfsw = topfsw_type( 0.0, 0.0, 0.0 )
294 do i=1,ncol
295 nirbmdi(i) = 0.0
296 nirdfdi(i) = 0.0
297 visbmdi(i) = 0.0
298 visdfdi(i) = 0.0
299 nirbmui(i) = 0.0
300 nirdfui(i) = 0.0
301 visbmui(i) = 0.0
302 visdfui(i) = 0.0
303 enddo
304
305 if (do_sw_clrsky_hr) then
306 htrswc(:,:) = 0
307 endif
308 endif ! end_if_nday
309
310 ! Radiation fluxes for other physics processes
311 do i=1,ncol
312 sfcnsw(i) = sfcfsw(i)%dnfxc - sfcfsw(i)%upfxc
313 sfcdsw(i) = sfcfsw(i)%dnfxc
314 enddo
315
316 ! #################################################################################
317 ! Save SW diagnostics
318 ! - For time averaged output quantities (including total-sky and clear-sky SW and LW
319 ! fluxes at TOA and surface; conventional 3-domain cloud amount, cloud top and base
320 ! pressure, and cloud top temperature; aerosols AOD, etc.), store computed results in
321 ! corresponding slots of array fluxr with appropriate time weights.
322 ! - Collect the fluxr data for wrtsfc
323 ! #################################################################################
324 if (save_diag) then
325 do i=1,ncol
326 fluxr(i,34) = aerodp(i,1) ! total aod at 550nm
327 fluxr(i,35) = aerodp(i,2) ! DU aod at 550nm
328 fluxr(i,36) = aerodp(i,3) ! BC aod at 550nm
329 fluxr(i,37) = aerodp(i,4) ! OC aod at 550nm
330 fluxr(i,38) = aerodp(i,5) ! SU aod at 550nm
331 fluxr(i,39) = aerodp(i,6) ! SS aod at 550nm
332 if (coszen(i) > 0.) then
333 ! SW all-sky fluxes
334 tem0d = fhswr * coszdg(i) / coszen(i)
335 fluxr(i,2 ) = fluxr(i,2) + topfsw(i)%upfxc * tem0d ! total sky top sw up
336 fluxr(i,3 ) = fluxr(i,3) + sfcfsw(i)%upfxc * tem0d
337 fluxr(i,4 ) = fluxr(i,4) + sfcfsw(i)%dnfxc * tem0d ! total sky sfc sw dn
338 ! SW uv-b fluxes
339 fluxr(i,21) = fluxr(i,21) + scmpsw(i)%uvbfc * tem0d ! total sky uv-b sw dn
340 fluxr(i,22) = fluxr(i,22) + scmpsw(i)%uvbf0 * tem0d ! clear sky uv-b sw dn
341 ! SW TOA incoming fluxes
342 fluxr(i,23) = fluxr(i,23) + topfsw(i)%dnfxc * tem0d ! top sw dn
343 ! SW SFC flux components
344 fluxr(i,24) = fluxr(i,24) + visbmdi(i) * tem0d ! uv/vis beam sw dn
345 fluxr(i,25) = fluxr(i,25) + visdfdi(i) * tem0d ! uv/vis diff sw dn
346 fluxr(i,26) = fluxr(i,26) + nirbmdi(i) * tem0d ! nir beam sw dn
347 fluxr(i,27) = fluxr(i,27) + nirdfdi(i) * tem0d ! nir diff sw dn
348 ! SW clear-sky fluxes
349 fluxr(i,29) = fluxr(i,29) + topfsw(i)%upfx0 * tem0d
350 fluxr(i,31) = fluxr(i,31) + sfcfsw(i)%upfx0 * tem0d
351 fluxr(i,32) = fluxr(i,32) + sfcfsw(i)%dnfx0 * tem0d
352 endif
353 enddo
354
355 ! Save total and boundary-layer clouds
356 do i=1,ncol
357 fluxr(i,17) = fluxr(i,17) + raddt * cldsa(i,4)
358 fluxr(i,18) = fluxr(i,18) + raddt * cldsa(i,5)
359 enddo
360
361 ! Save cld frac,toplyr,botlyr and top temp, note that the order of h,m,l cloud
362 ! is reversed for the fluxr output. save interface pressure (pa) of top/bot
363 do j = 1, 3
364 do i = 1, ncol
365 tem0d = raddt * cldsa(i,j)
366 itop = mtopa(i,j)
367 ibtc = mbota(i,j)
368 fluxr(i, 8-j) = fluxr(i, 8-j) + tem0d
369 fluxr(i,11-j) = fluxr(i,11-j) + tem0d * p_lev(i,itop)
370 fluxr(i,14-j) = fluxr(i,14-j) + tem0d * p_lev(i,ibtc)
371 fluxr(i,17-j) = fluxr(i,17-j) + tem0d * p_lev(i,itop)
372
373 ! Add optical depth and emissivity output
374 tem1 = 0.
375 do k=ibtc,itop
376 tem1 = tem1 + cldtausw(i,k) ! approx .55 mu channel
377 enddo
378 fluxr(i,43-j) = fluxr(i,43-j) + tem0d * tem1
379 enddo
380 enddo
381 endif
382 endif
383
384 end subroutine gfs_rrtmgp_post_run
385end module gfs_rrtmgp_post
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
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 special components of surface SW fluxes
derived type for SW fluxes at surface
Definition radsw_param.f:89
derived type for SW fluxes at TOA
Definition radsw_param.f:79