84 subroutine rrtmgp_lw_main_run(doLWrad, doLWclrsky, top_at_1, doGP_lwscat, &
85 use_LW_jacobian, doGP_sgs_cnv, doGP_sgs_pbl, nCol, nLay, nGases,rrtmgp_phys_blksz,&
86 nGauss_angles, icseed_lw, iovr, iovr_convcld, iovr_max, iovr_maxrand, iovr_rand, &
87 iovr_dcorr, iovr_exp, iovr_exprand, isubc_lw, semis, tsfg, p_lay, p_lev, t_lay, &
88 t_lev, vmr_o2, vmr_h2o, vmr_o3, vmr_ch4, vmr_n2o, vmr_co2, &
89 cld_frac, cld_lwp, cld_reliq, cld_iwp, cld_reice, cld_swp, cld_resnow, &
90 cld_rwp, cld_rerain, precip_frac, cld_cnv_lwp, cld_cnv_reliq, cld_cnv_iwp, &
91 cld_cnv_reice, cld_pbl_lwp, cld_pbl_reliq, cld_pbl_iwp, cld_pbl_reice, &
92 cloud_overlap_param, active_gases_array, aerlw_tau, aerlw_ssa, aerlw_g, &
93 fluxlwUP_allsky, fluxlwDOWN_allsky, fluxlwUP_clrsky, fluxlwDOWN_clrsky, &
94 fluxlwUP_jac, fluxlwUP_radtime, fluxlwDOWN_radtime, errmsg, errflg)
97 logical,
intent(in) :: &
105 integer,
intent(in) :: &
120 integer,
intent(in),
dimension(:),
optional :: &
122 real(kind_phys),
dimension(:),
intent(in) :: &
125 real(kind_phys),
dimension(:,:),
intent(in),
optional :: &
136 real(kind_phys),
dimension(:,:),
intent(in) :: &
146 real(kind_phys),
dimension(:,:),
intent(in),
optional :: &
157 real(kind_phys),
dimension(:,:,:),
intent(in) :: &
161 character(len=*),
dimension(:),
intent(in),
optional :: &
165 real(kind_phys),
dimension(:,:),
intent(inout),
optional :: &
173 character(len=*),
intent(out) :: &
175 integer,
intent(out) :: &
179 type(ty_fluxes_byband) :: flux_allsky, flux_clrsky
180 integer :: icol, ilay, igas, iband, icol2, ix, iblck
181 integer,
dimension(rrtmgp_phys_blksz) :: ipseed_lw
183 real(kind_phys),
dimension(rrtmgp_phys_blksz) :: zcf0, zcf1
184 logical,
dimension(rrtmgp_phys_blksz,nLay,lw_gas_props%get_ngpt()) :: maskmcica
185 real(kind_phys),
dimension(rrtmgp_phys_blksz) :: tau_rain, tau_snow
186 real(kind_dbl_prec),
dimension(lw_gas_props%get_ngpt()) :: rng1d
187 real(kind_dbl_prec),
dimension(lw_gas_props%get_ngpt(),nLay,rrtmgp_phys_blksz) :: rng3d,rng3d2
188 real(kind_dbl_prec),
dimension(lw_gas_props%get_ngpt()*nLay) :: rng2d
189 real(kind_phys),
dimension(rrtmgp_phys_blksz,nLay+1,lw_gas_props%get_nband()),
target :: &
190 fluxlw_up_allsky, fluxlw_up_clrsky, fluxlw_dn_allsky, fluxlw_dn_clrsky
191 real(kind_phys),
dimension(rrtmgp_phys_blksz,lw_gas_props%get_ngpt()) :: lw_ds
192 real(kind_phys),
dimension(lw_gas_props%get_nband(),rrtmgp_phys_blksz) :: sfc_emiss_byband
195 type(ty_gas_concs) :: gas_concs
196 type(ty_optical_props_1scl) :: lw_optical_props_clrsky, lw_optical_props_aerosol_local
197 type(ty_optical_props_2str) :: lw_optical_props_clouds, lw_optical_props_cloudsbyband, &
198 lw_optical_props_cnvcloudsbyband, lw_optical_props_pblcloudsbyband, &
199 lw_optical_props_precipbyband
200 type(ty_source_func_lw) :: sources
206 if (.not. dolwrad)
return
213 call check_error_msg(
'rrtmgp_lw_main_gas_concs_run',gas_concs%init(active_gases_array))
216 call check_error_msg(
'rrtmgp_lw_main_gas_optics_run',&
217 lw_optical_props_clrsky%alloc_1scl(rrtmgp_phys_blksz, nlay, lw_gas_props))
218 call check_error_msg(
'rrtmgp_lw_main_sources_run',&
219 sources%alloc(rrtmgp_phys_blksz, nlay, lw_gas_props))
220 call check_error_msg(
'rrtmgp_lw_main_cloud_optics_run',&
221 lw_optical_props_cloudsbyband%alloc_2str(rrtmgp_phys_blksz, nlay, lw_gas_props%get_band_lims_wavenumber()))
222 call check_error_msg(
'rrtmgp_lw_main_precip_optics_run',&
223 lw_optical_props_precipbyband%alloc_2str(rrtmgp_phys_blksz, nlay, lw_gas_props%get_band_lims_wavenumber()))
224 call check_error_msg(
'rrtmgp_lw_mian_cloud_sampling_run', &
225 lw_optical_props_clouds%alloc_2str(rrtmgp_phys_blksz, nlay, lw_gas_props))
226 call check_error_msg(
'rrtmgp_lw_main_aerosol_optics_run',&
227 lw_optical_props_aerosol_local%alloc_1scl(rrtmgp_phys_blksz, nlay, lw_gas_props%get_band_lims_wavenumber()))
228 if (dogp_sgs_cnv)
then
229 call check_error_msg(
'rrtmgp_lw_main_cnv_cloud_optics_run',&
230 lw_optical_props_cnvcloudsbyband%alloc_2str(rrtmgp_phys_blksz, nlay, lw_gas_props%get_band_lims_wavenumber()))
232 if (dogp_sgs_pbl)
then
233 call check_error_msg(
'rrtmgp_lw_main_pbl_cloud_optics_run',&
234 lw_optical_props_pblcloudsbyband%alloc_2str(rrtmgp_phys_blksz, nlay, lw_gas_props%get_band_lims_wavenumber()))
242 do icol=1,ncol,rrtmgp_phys_blksz
243 icol2 = icol + rrtmgp_phys_blksz - 1
248 lw_optical_props_clrsky%tau = 0._kind_phys
249 lw_optical_props_precipbyband%tau = 0._kind_phys
250 lw_optical_props_precipbyband%ssa = 0._kind_phys
251 lw_optical_props_precipbyband%g = 0._kind_phys
252 lw_optical_props_cloudsbyband%tau = 0._kind_phys
253 lw_optical_props_cloudsbyband%ssa = 0._kind_phys
254 lw_optical_props_cloudsbyband%g = 0._kind_phys
255 lw_optical_props_clouds%tau = 0._kind_phys
256 lw_optical_props_clouds%ssa = 0._kind_phys
257 lw_optical_props_clouds%g = 0._kind_phys
258 sources%sfc_source = 0._kind_phys
259 sources%lay_source = 0._kind_phys
260 sources%lev_source_inc = 0._kind_phys
261 sources%lev_source_dec = 0._kind_phys
262 sources%sfc_source_Jac = 0._kind_phys
263 fluxlw_up_allsky = 0._kind_phys
264 fluxlw_dn_allsky = 0._kind_phys
265 fluxlw_up_clrsky = 0._kind_phys
266 fluxlw_dn_clrsky = 0._kind_phys
267 if (dogp_sgs_cnv) lw_optical_props_cnvcloudsbyband%tau = 0._kind_phys
268 if (dogp_sgs_pbl) lw_optical_props_pblcloudsbyband%tau = 0._kind_phys
271 fluxlw_up_allsky = 0._kind_phys
272 fluxlw_dn_allsky = 0._kind_phys
273 fluxlw_up_clrsky = 0._kind_phys
274 fluxlw_dn_clrsky = 0._kind_phys
275 flux_allsky%bnd_flux_up => fluxlw_up_allsky
276 flux_allsky%bnd_flux_dn => fluxlw_dn_allsky
277 flux_clrsky%bnd_flux_up => fluxlw_up_clrsky
278 flux_clrsky%bnd_flux_dn => fluxlw_dn_clrsky
285 call check_error_msg(
'rrtmgp_lw_main_set_vmr_o2', &
286 gas_concs%set_vmr(trim(active_gases_array(istr_o2)), vmr_o2(icol:icol2,:)))
287 call check_error_msg(
'rrtmgp_lw_main_set_vmr_co2', &
288 gas_concs%set_vmr(trim(active_gases_array(istr_co2)),vmr_co2(icol:icol2,:)))
289 call check_error_msg(
'rrtmgp_lw_main_set_vmr_ch4', &
290 gas_concs%set_vmr(trim(active_gases_array(istr_ch4)),vmr_ch4(icol:icol2,:)))
291 call check_error_msg(
'rrtmgp_lw_main_set_vmr_n2o', &
292 gas_concs%set_vmr(trim(active_gases_array(istr_n2o)),vmr_n2o(icol:icol2,:)))
293 call check_error_msg(
'rrtmgp_lw_main_set_vmr_h2o', &
294 gas_concs%set_vmr(trim(active_gases_array(istr_h2o)),vmr_h2o(icol:icol2,:)))
295 call check_error_msg(
'rrtmgp_lw_main_set_vmr_o3', &
296 gas_concs%set_vmr(trim(active_gases_array(istr_o3)), vmr_o3(icol:icol2,:)))
304 do iblck=1,rrtmgp_phys_blksz
305 if (semis(icol+iblck-1) > eps .and. semis(icol+iblck-1) <= 1._kind_phys)
then
306 do iband=1,lw_gas_props%get_nband()
307 sfc_emiss_byband(iband,iblck) = semis(icol+iblck-1)
310 sfc_emiss_byband(1:lw_gas_props%get_nband(),iblck) = 1.0
319 call check_error_msg(
'rrtmgp_lw_main_gas_optics',lw_gas_props%gas_optics(&
320 p_lay(icol:icol2,:), &
321 p_lev(icol:icol2,:), &
322 t_lay(icol:icol2,:), &
325 lw_optical_props_clrsky, &
327 tlev=t_lev(icol:icol2,:)))
335 zcf0(:) = 1._kind_phys
336 zcf1(:) = 1._kind_phys
337 do iblck = 1, rrtmgp_phys_blksz
339 zcf0(iblck) = min(zcf0(iblck), 1._kind_phys - cld_frac(icol+iblck-1,ilay))
341 if (zcf0(iblck) <= ftiny) zcf0(iblck) = 0._kind_phys
342 if (zcf0(iblck) > oneminus) zcf0(iblck) = 1._kind_phys
343 zcf1(iblck) = 1._kind_phys - zcf0(iblck)
346 if (any(zcf1 .gt. eps))
then
348 call check_error_msg(
'rrtmgp_lw_main_cloud_optics',lw_cloud_props%cloud_optics(&
349 cld_lwp(icol:icol2,:), &
350 cld_iwp(icol:icol2,:), &
351 cld_reliq(icol:icol2,:), &
352 cld_reice(icol:icol2,:), &
353 lw_optical_props_cloudsbyband))
356 if (dogp_sgs_cnv)
then
358 call check_error_msg(
'rrtmgp_lw_main_cnv_cloud_optics',lw_cloud_props%cloud_optics(&
359 cld_cnv_lwp(icol:icol2,:), &
360 cld_cnv_iwp(icol:icol2,:), &
361 cld_cnv_reliq(icol:icol2,:), &
362 cld_cnv_reice(icol:icol2,:), &
363 lw_optical_props_cnvcloudsbyband))
366 call check_error_msg(
'rrtmgp_lw_main_increment_cnvclouds_to_clouds',&
367 lw_optical_props_cnvcloudsbyband%increment(lw_optical_props_cloudsbyband))
371 if (dogp_sgs_pbl)
then
373 call check_error_msg(
'rrtmgp_lw_main_pbl_cloud_optics',lw_cloud_props%cloud_optics(&
374 cld_pbl_lwp(icol:icol2,:), &
375 cld_pbl_iwp(icol:icol2,:), &
376 cld_pbl_reliq(icol:icol2,:), &
377 cld_pbl_reice(icol:icol2,:), &
378 lw_optical_props_pblcloudsbyband))
381 call check_error_msg(
'rrtmgp_lw_main_increment_pblclouds_to_clouds',&
382 lw_optical_props_pblcloudsbyband%increment(lw_optical_props_cloudsbyband))
391 tau_rain(:) = 0._kind_phys
392 tau_snow(:) = 0._kind_phys
393 do ix=1,rrtmgp_phys_blksz
395 if (cld_frac(icol+ix-1,ilay) .gt. eps)
then
397 tau_rain(ix) = absrain*cld_rwp(icol+ix-1,ilay)
400 if (cld_swp(icol+ix-1,ilay) .gt. 0. .and. cld_resnow(icol+ix-1,ilay) .gt. 10._kind_phys)
then
401 tau_snow(ix) = abssnow0*1.05756*cld_swp(icol+ix-1,ilay)/cld_resnow(icol+ix-1,ilay)
405 do iband=1,lw_gas_props%get_nband()
406 lw_optical_props_precipbyband%tau(ix,ilay,iband) = tau_rain(ix) + tau_snow(ix)
412 call check_error_msg(
'rrtmgp_lw_main_increment_precip_to_clouds',&
413 lw_optical_props_precipbyband%increment(lw_optical_props_cloudsbyband))
421 if (any(zcf1 .gt. eps))
then
423 if(isubc_lw == 1)
then
424 do ix=1,rrtmgp_phys_blksz
425 ipseed_lw(ix) = lw_gas_props%get_ngpt() + icol + ix - 1
427 elseif (isubc_lw == 2)
then
428 do ix=1,rrtmgp_phys_blksz
429 ipseed_lw(ix) = icseed_lw(icol+ix-1)
434 do ix=1,rrtmgp_phys_blksz
437 if (iovr == iovr_max)
then
440 rng3d(:,ilay,ix) = rng1d
445 rng3d(:,ilay,ix) = rng1d
452 if (iovr == iovr_maxrand .or. iovr == iovr_rand .or. iovr == iovr_max)
then
453 call sampled_mask(real(rng3d,kind=kind_phys), cld_frac(icol:icol2,:), maskmcica)
456 if (iovr == iovr_dcorr)
then
457 do ix=1,rrtmgp_phys_blksz
461 rng3d2(:,:,ix) = reshape(source = rng2d,shape=[lw_gas_props%get_ngpt(),nlay])
464 call sampled_mask(real(rng3d,kind=kind_phys), cld_frac(icol:icol2,:), maskmcica, &
465 overlap_param = cloud_overlap_param(icol:icol2,1:nlay-1), randoms2 = real(rng3d2, kind=kind_phys))
468 if (iovr == iovr_exp .or. iovr == iovr_exprand)
then
469 call sampled_mask(real(rng3d,kind=kind_phys), cld_frac(icol:icol2,:), maskmcica, &
470 overlap_param = cloud_overlap_param(icol:icol2,1:nlay-1))
473 call check_error_msg(
'rrtmgp_lw_main_cloud_sampling',&
474 draw_samples(maskmcica, .true., &
475 lw_optical_props_cloudsbyband, lw_optical_props_clouds))
484 lw_optical_props_aerosol_local%tau = aerlw_tau(icol:icol2,:,:)
485 call check_error_msg(
'rrtmgp_lw_main_increment_aerosol_to_clrsky',&
486 lw_optical_props_aerosol_local%increment(lw_optical_props_clrsky))
490 call check_error_msg(
'rrtmgp_lw_main_opt_angle',&
491 lw_gas_props%compute_optimal_angles(lw_optical_props_clrsky,lw_ds))
492 if (ngauss_angles .gt. 1)
then
493 call check_error_msg(
'rrtmgp_lw_main_lw_rte_clrsky',rte_lw( &
494 lw_optical_props_clrsky, &
499 n_gauss_angles = ngauss_angles))
501 call check_error_msg(
'rrtmgp_lw_main_lw_rte_clrsky',rte_lw( &
502 lw_optical_props_clrsky, &
511 fluxlwup_clrsky(icol:icol2,:) = sum(flux_clrsky%bnd_flux_up, dim=3)
512 fluxlwdown_clrsky(icol:icol2,:) = sum(flux_clrsky%bnd_flux_dn, dim=3)
514 fluxlwup_clrsky(icol:icol2,:) = 0.0
515 fluxlwdown_clrsky(icol:icol2,:) = 0.0
534 if (dogp_lwscat)
then
536 call check_error_msg(
'rrtmgp_lw_main_increment_clrsky_to_clouds',&
537 lw_optical_props_clrsky%increment(lw_optical_props_clouds))
539 if (use_lw_jacobian)
then
541 call check_error_msg(
'rrtmgp_lw_main_lw_rte_allsky',rte_lw( &
542 lw_optical_props_clouds, &
547 n_gauss_angles = ngauss_angles, &
548 flux_up_jac = fluxlwup_jac))
550 call check_error_msg(
'rrtmgp_lw_main_lw_rte_allsky',rte_lw( &
551 lw_optical_props_clouds, &
556 n_gauss_angles = ngauss_angles))
561 call check_error_msg(
'rrtmgp_lw_main_increment_clouds_to_clrsky', &
562 lw_optical_props_clouds%increment(lw_optical_props_clrsky))
564 if (use_lw_jacobian)
then
566 call check_error_msg(
'rrtmgp_lw_rte_run',rte_lw( &
567 lw_optical_props_clrsky, &
572 n_gauss_angles = ngauss_angles, &
573 flux_up_jac = fluxlwup_jac))
575 call check_error_msg(
'rrtmgp_lw_rte_run',rte_lw( &
576 lw_optical_props_clrsky, &
581 n_gauss_angles = ngauss_angles))
586 fluxlwup_allsky(icol:icol2,:) = sum(flux_allsky%bnd_flux_up, dim=3)
587 fluxlwdown_allsky(icol:icol2,:) = sum(flux_allsky%bnd_flux_dn, dim=3)
590 fluxlwup_radtime(icol:icol2,:) = fluxlwup_allsky(icol:icol2,:)
591 fluxlwdown_radtime(icol:icol2,:) = fluxlwdown_allsky(icol:icol2,:)