SCALE-RM
Functions/Subroutines
mod_atmos_phy_rd_driver Module Reference

module ATMOSPHERE / Physics Radiation More...

Functions/Subroutines

subroutine, public atmos_phy_rd_driver_setup
 Setup. More...
 
subroutine, public atmos_phy_rd_driver_resume
 Resume. More...
 
subroutine, public atmos_phy_rd_driver (update_flag)
 Driver. More...
 

Detailed Description

module ATMOSPHERE / Physics Radiation

Description
Atmospheric radiation transfer process driver
Author
Team SCALE
History
  • 2013-12-06 (S.Nishizawa) [new]
NAMELIST
  • No namelist group
History Output
namedescriptionunitvariable
COSZ cos(solar zenith angle) 1 cosSZA
OLR TOA net longwave radiation flux W/m2 flux_net_toa
OSR TOA net shortwave radiation flux W/m2 flux_net_toa
RADFLUX_LW net longwave radiation flux W/m2 flux_net
RADFLUX_LWDN downward longwave radiation flux W/m2 flux_dn
RADFLUX_LWUP upward longwave radiation flux W/m2 flux_up
RADFLUX_SW net shortwave radiation flux W/m2 flux_net
RADFLUX_SWDN downward shortwave radiation flux W/m2 flux_dn
RADFLUX_SWUP upward shortwave radiation flux W/m2 flux_up
RFLX_LW_dn downward longwave radiation flux (cell face) W/m2 flux_rad
RFLX_LW_up upward longwave radiation flux (cell face) W/m2 flux_rad
RFLX_SW_dn downward shortwave radiation flux (cell face) W/m2 flux_rad
RFLX_SW_up upward shortwave radiation flux (cell face) W/m2 flux_rad
RHOT_t_RD tendency of RHOT in rd K.kg/m3/s RHOT_t_RD
SFLX_LW_dn SFC downward longwave radiation flux W/m2 SFCFLX_LW_dn
SFLX_LW_dn_c SFC downward longwave radiation flux (clr) W/m2 SFCFLX_LW_dn_c
SFLX_LW_dn_dif SFC downward longwave flux (diffuse) W/m2 SFLX_rad_dn
SFLX_LW_dn_dir SFC downward longwave flux (direct ) W/m2 SFLX_rad_dn
SFLX_LW_up SFC upward longwave radiation flux W/m2 SFCFLX_LW_up
SFLX_LW_up_c SFC upward longwave radiation flux (clr) W/m2 SFCFLX_LW_up_c
SFLX_SW_dn SFC downward shortwave radiation flux W/m2 SFCFLX_SW_dn
SFLX_SW_dn_c SFC downward shortwave radiation flux (clr) W/m2 SFCFLX_SW_dn_c
SFLX_SW_dn_dif SFC downward shortwave flux (diffuse) W/m2 SFLX_rad_dn
SFLX_SW_dn_dir SFC downward shortwave flux (direct ) W/m2 SFLX_rad_dn
SFLX_SW_up SFC upward shortwave radiation flux W/m2 SFCFLX_SW_up
SFLX_SW_up_c SFC upward shortwave radiation flux (clr) W/m2 SFCFLX_SW_up_c
SLR SFC net longwave radiation flux W/m2 flux_net_sfc
SOLINS solar insolation W/m2 solins
SSR SFC net shortwave radiation flux W/m2 flux_net_sfc
TEMP_t_rd tendency of temp in rd K/day TEMP_t
TEMP_t_rd_LW tendency of temp in rd(LW) K/day TEMP_t
TEMP_t_rd_SW tendency of temp in rd(SW) K/day TEMP_t
TOAFLX_LW_dn TOA downward longwave radiation flux W/m2 TOAFLX_LW_dn
TOAFLX_LW_dn_c TOA downward longwave radiation flux (clr) W/m2 TOAFLX_LW_dn_c
TOAFLX_LW_up TOA upward longwave radiation flux W/m2 TOAFLX_LW_up
TOAFLX_LW_up_c TOA upward longwave radiation flux (clr) W/m2 TOAFLX_LW_up_c
TOAFLX_SW_dn TOA downward shortwave radiation flux W/m2 TOAFLX_SW_dn
TOAFLX_SW_dn_c TOA downward shortwave radiation flux (clr) W/m2 TOAFLX_SW_dn_c
TOAFLX_SW_up TOA upward shortwave radiation flux W/m2 TOAFLX_SW_up
TOAFLX_SW_up_c TOA upward shortwave radiation flux (clr) W/m2 TOAFLX_SW_up_c

Function/Subroutine Documentation

◆ atmos_phy_rd_driver_setup()

subroutine, public mod_atmos_phy_rd_driver::atmos_phy_rd_driver_setup ( )

Setup.

Definition at line 52 of file mod_atmos_phy_rd_driver.f90.

References mod_atmos_phy_rd_vars::atmos_phy_rd_cossza, scale_atmos_phy_rd::atmos_phy_rd_setup(), mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_downall, mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_lw_dn, mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_lw_up, mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_sw_dn, mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_sw_up, mod_atmos_phy_rd_vars::atmos_phy_rd_solins, mod_atmos_phy_rd_vars::atmos_phy_rd_toaflx_lw_dn, mod_atmos_phy_rd_vars::atmos_phy_rd_toaflx_lw_up, mod_atmos_phy_rd_vars::atmos_phy_rd_toaflx_sw_dn, mod_atmos_phy_rd_vars::atmos_phy_rd_toaflx_sw_up, mod_atmos_admin::atmos_phy_rd_type, mod_atmos_admin::atmos_sw_phy_rd, scale_stdio::io_fid_log, and scale_stdio::io_l.

Referenced by mod_atmos_driver::atmos_driver_setup().

52  use scale_atmos_phy_rd, only: &
54  use mod_atmos_admin, only: &
57  use mod_atmos_phy_rd_vars, only: &
58  sfcflx_lw_up => atmos_phy_rd_sflx_lw_up, &
59  sfcflx_lw_dn => atmos_phy_rd_sflx_lw_dn, &
60  sfcflx_sw_up => atmos_phy_rd_sflx_sw_up, &
61  sfcflx_sw_dn => atmos_phy_rd_sflx_sw_dn, &
62  toaflx_lw_up => atmos_phy_rd_toaflx_lw_up, &
63  toaflx_lw_dn => atmos_phy_rd_toaflx_lw_dn, &
64  toaflx_sw_up => atmos_phy_rd_toaflx_sw_up, &
65  toaflx_sw_dn => atmos_phy_rd_toaflx_sw_dn, &
66  sflx_rad_dn => atmos_phy_rd_sflx_downall, &
67  solins => atmos_phy_rd_solins, &
68  cossza => atmos_phy_rd_cossza
69  implicit none
70  !---------------------------------------------------------------------------
71 
72  if( io_l ) write(io_fid_log,*)
73  if( io_l ) write(io_fid_log,*) '++++++ Module[DRIVER] / Categ[ATMOS PHY_RD] / Origin[SCALE-RM]'
74 
75  if ( atmos_sw_phy_rd ) then
76 
77  ! setup library component
79 
80  else
81 
82  if( io_l ) write(io_fid_log,*) '*** this component is never called.'
83  if( io_l ) write(io_fid_log,*) '*** radiation fluxes are set to zero.'
84  sfcflx_lw_up(:,:) = 0.0_rp
85  sfcflx_lw_dn(:,:) = 0.0_rp
86  sfcflx_sw_up(:,:) = 0.0_rp
87  sfcflx_sw_dn(:,:) = 0.0_rp
88  toaflx_lw_up(:,:) = 0.0_rp
89  toaflx_lw_dn(:,:) = 0.0_rp
90  toaflx_sw_up(:,:) = 0.0_rp
91  toaflx_sw_dn(:,:) = 0.0_rp
92  sflx_rad_dn(:,:,:,:) = 0.0_rp
93  solins(:,:) = 0.0_rp
94  cossza(:,:) = 0.0_rp
95 
96  endif
97 
98  return
module ATMOS admin
logical, public atmos_sw_phy_rd
character(len=h_short), public atmos_phy_rd_type
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_sflx_lw_up
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_toaflx_sw_up
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_sflx_lw_dn
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_cossza
module Atmosphere / Physics Radiation
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_sflx_sw_up
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_sflx_sw_dn
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_toaflx_lw_dn
module ATMOSPHERE / Physics Radiation
real(rp), dimension(:,:,:,:), allocatable, public atmos_phy_rd_sflx_downall
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_toaflx_sw_dn
subroutine, public atmos_phy_rd_setup(RD_TYPE)
Setup.
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_solins
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_toaflx_lw_up
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_phy_rd_driver_resume()

subroutine, public mod_atmos_phy_rd_driver::atmos_phy_rd_driver_resume ( )

Resume.

Definition at line 104 of file mod_atmos_phy_rd_driver.f90.

References atmos_phy_rd_driver(), mod_atmos_admin::atmos_sw_phy_rd, scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

Referenced by mod_atmos_driver::atmos_driver_resume1().

104  use mod_atmos_admin, only: &
106  implicit none
107 
108  if ( atmos_sw_phy_rd ) then
109 
110  ! run once (only for the diagnostic value)
111  call prof_rapstart('ATM_Radiation', 1)
112  call atmos_phy_rd_driver( update_flag = .true. )
113  call prof_rapend ('ATM_Radiation', 1)
114 
115  end if
116 
117  return
module ATMOS admin
logical, public atmos_sw_phy_rd
subroutine, public atmos_phy_rd_driver(update_flag)
Driver.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_phy_rd_driver()

subroutine, public mod_atmos_phy_rd_driver::atmos_phy_rd_driver ( logical, intent(in)  update_flag)

Driver.

Definition at line 123 of file mod_atmos_phy_rd_driver.f90.

References scale_atmos_phy_rd::atmos_phy_rd, mod_atmos_phy_rd_vars::atmos_phy_rd_cossza, scale_atmos_phy_rd_common::atmos_phy_rd_heating(), mod_atmos_phy_rd_vars::atmos_phy_rd_rhot_t, mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_downall, mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_lw_dn, mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_lw_up, mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_sw_dn, mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_sw_up, mod_atmos_phy_rd_vars::atmos_phy_rd_solins, mod_atmos_phy_rd_vars::atmos_phy_rd_toaflx_lw_dn, mod_atmos_phy_rd_vars::atmos_phy_rd_toaflx_lw_up, mod_atmos_phy_rd_vars::atmos_phy_rd_toaflx_sw_dn, mod_atmos_phy_rd_vars::atmos_phy_rd_toaflx_sw_up, mod_atmos_admin::atmos_phy_rd_type, mod_atmos_phy_sf_vars::atmos_phy_sf_sfc_albedo, mod_atmos_phy_sf_vars::atmos_phy_sf_sfc_temp, scale_const::const_cpdry, scale_const::const_pre00, scale_const::const_rdry, mod_atmos_vars::dens, mod_atmos_vars::dens_av, scale_atmos_phy_rd_common::i_diffuse, scale_atmos_phy_rd_common::i_direct, scale_atmos_phy_rd_common::i_dn, scale_atmos_phy_rd_common::i_lw, scale_atmos_hydrometeor::i_qc, scale_atmos_hydrometeor::i_qg, scale_atmos_hydrometeor::i_qi, scale_atmos_hydrometeor::i_qr, scale_atmos_hydrometeor::i_qs, scale_atmos_hydrometeor::i_qv, scale_atmos_phy_rd_common::i_sw, scale_atmos_phy_rd_common::i_up, scale_grid_index::ia, scale_grid_index::ie, scale_grid_index::is, scale_grid_index::ja, scale_grid_index::je, scale_grid_index::js, scale_grid_index::ka, scale_grid_index::ke, scale_grid_index::ks, scale_landuse::landuse_fact_land, scale_landuse::landuse_fact_ocean, scale_landuse::landuse_fact_urban, scale_process::prc_mpistop(), mod_atmos_vars::pres, mod_atmos_vars::qtrc, mod_atmos_vars::qtrc_av, scale_grid_real::real_cz, scale_grid_real::real_fz, scale_grid_real::real_lat, scale_grid_real::real_lon, mod_atmos_vars::rhot, mod_atmos_vars::rhot_av, mod_atmos_vars::rhot_tp, scale_rm_statistics::statistics_checktotal, scale_atmos_phy_rd_mm5sw::swrad(), mod_atmos_vars::temp, scale_time::time_dtsec_atmos_phy_rd, scale_time::time_nowdate, and scale_time::time_offset_year.

Referenced by mod_atmos_driver::atmos_driver(), and atmos_phy_rd_driver_resume().

123  use scale_grid_real, only: &
124  real_cz, &
125  real_fz, &
126  real_lon, &
127  real_lat
128  use scale_process, only: &
130  use scale_const, only: &
131  pre00 => const_pre00, &
132  rdry => const_rdry, &
133  cpdry => const_cpdry
134  use scale_landuse, only: &
138  use scale_time, only: &
139  dt_rd => time_dtsec_atmos_phy_rd, &
140  time_nowdate, &
142  use scale_rm_statistics, only: &
144  stat_total
145  use scale_history, only: &
146  hist_in
147  use mod_atmos_admin, only: &
149  use scale_atmos_solarins, only: &
150  solarins_insolation => atmos_solarins_insolation
151  use scale_atmos_phy_rd, only: &
153  use scale_atmos_phy_rd_mm5sw, only: &
154  swrad
155  use scale_atmos_phy_rd_common, only: &
156  rd_heating => atmos_phy_rd_heating, &
157  i_sw, &
158  i_lw, &
159  i_dn, &
160  i_up, &
161  i_direct, &
162  i_diffuse
163  use scale_atmos_hydrometeor, only: &
164  i_qv, &
165  i_qc, &
166  i_qr, &
167  i_qi, &
168  i_qs, &
169  i_qg
170  use mod_atmos_vars, only: &
171  temp, &
172  pres, &
173  dens => dens_av, &
174  rhot => rhot_av, &
175  qtrc => qtrc_av, &
176  rhot_t => rhot_tp
177  use mod_atmos_phy_sf_vars, only: &
178  sfc_temp => atmos_phy_sf_sfc_temp, &
179  sfc_albedo => atmos_phy_sf_sfc_albedo
180  use mod_atmos_phy_rd_vars, only: &
181  rhot_t_rd => atmos_phy_rd_rhot_t, &
182  sfcflx_lw_up => atmos_phy_rd_sflx_lw_up, &
183  sfcflx_lw_dn => atmos_phy_rd_sflx_lw_dn, &
184  sfcflx_sw_up => atmos_phy_rd_sflx_sw_up, &
185  sfcflx_sw_dn => atmos_phy_rd_sflx_sw_dn, &
186  toaflx_lw_up => atmos_phy_rd_toaflx_lw_up, &
187  toaflx_lw_dn => atmos_phy_rd_toaflx_lw_dn, &
188  toaflx_sw_up => atmos_phy_rd_toaflx_sw_up, &
189  toaflx_sw_dn => atmos_phy_rd_toaflx_sw_dn, &
190  sflx_rad_dn => atmos_phy_rd_sflx_downall, &
191  solins => atmos_phy_rd_solins, &
192  cossza => atmos_phy_rd_cossza
193  implicit none
194 
195  logical, intent(in) :: update_flag
196 
197  real(RP) :: TEMP_t (KA,IA,JA,3)
198  real(RP) :: flux_rad (KA,IA,JA,2,2,2)
199  real(RP) :: flux_rad_top( IA,JA,2,2,2)
200 
201  real(RP) :: flux_up (KA,IA,JA,2)
202  real(RP) :: flux_dn (KA,IA,JA,2)
203  real(RP) :: flux_net (KA,IA,JA,2)
204  real(RP) :: flux_net_toa( IA,JA,2)
205  real(RP) :: flux_net_sfc( IA,JA,2)
206 
207  real(RP) :: SFCFLX_LW_up_c(IA,JA)
208  real(RP) :: SFCFLX_LW_dn_c(IA,JA)
209  real(RP) :: SFCFLX_SW_up_c(IA,JA)
210  real(RP) :: SFCFLX_SW_dn_c(IA,JA)
211  real(RP) :: TOAFLX_LW_up_c(IA,JA)
212  real(RP) :: TOAFLX_LW_dn_c(IA,JA)
213  real(RP) :: TOAFLX_SW_up_c(IA,JA)
214  real(RP) :: TOAFLX_SW_dn_c(IA,JA)
215 
216  ! for WRF radiation scheme added by Adachi; array order is (i,k,j)
217  real(RP) :: RTHRATENSW(IA,KA,JA)
218  real(RP) :: SDOWN3D (IA,KA,JA) ! downward short wave flux (W/m2)
219  real(RP) :: GSW (IA,JA) ! net short wave flux at ground surface (W/m2)
220  real(RP) :: RHO3D (IA,KA,JA)
221  real(RP) :: T3D (IA,KA,JA)
222  real(RP) :: P3D (IA,KA,JA)
223  real(RP) :: PI3D (IA,KA,JA)
224  real(RP) :: DZ8W (IA,KA,JA)
225  real(RP) :: QV3D (IA,KA,JA)
226  real(RP) :: QC3D (IA,KA,JA)
227  real(RP) :: QR3D (IA,KA,JA)
228  real(RP) :: QI3D (IA,KA,JA)
229  real(RP) :: QS3D (IA,KA,JA)
230  real(RP) :: QG3D (IA,KA,JA)
231  real(RP) :: flux_rad_org(KA,IA,JA,2,2,2)
232 
233  real(RP) :: total ! dummy
234 
235  integer :: k, i, j
236  !---------------------------------------------------------------------------
237 
238  if ( update_flag ) then
239 
240  call solarins_insolation( solins(:,:), & ! [OUT]
241  cossza(:,:), & ! [OUT]
242  real_lon(:,:), & ! [IN]
243  real_lat(:,:), & ! [IN]
244  time_nowdate(:), & ! [IN]
245  time_offset_year ) ! [IN]
246 
247 
248  call atmos_phy_rd( dens, rhot, qtrc, & ! [IN]
249  real_cz, real_fz, & ! [IN]
250  landuse_fact_ocean, & ! [IN]
251  landuse_fact_land, & ! [IN]
252  landuse_fact_urban, & ! [IN]
253  sfc_temp, & ! [IN]
254  sfc_albedo, & ! [IN]
255  solins, cossza, & ! [IN]
256  flux_rad, & ! [OUT]
257  flux_rad_top, & ! [OUT]
258  sflx_rad_dn ) ! [OUT]
259 
260 
261  ! apply radiative flux convergence -> heating rate
262  call rd_heating( flux_rad(:,:,:,:,:,2), & ! [IN]
263  dens(:,:,:), & ! [IN]
264  rhot(:,:,:), & ! [IN]
265  qtrc(:,:,:,:), & ! [IN]
266  real_fz(:,:,:), & ! [IN]
267  dt_rd, & ! [IN]
268  temp_t(:,:,:,:), & ! [OUT]
269  rhot_t_rd(:,:,:) ) ! [OUT]
270 
271 
272 !OCL XFILL
273  do j = js, je
274  do i = is, ie
275  ! for clear-sky
276  sfcflx_lw_up_c(i,j) = flux_rad(ks-1,i,j,i_lw,i_up,1)
277  sfcflx_lw_dn_c(i,j) = flux_rad(ks-1,i,j,i_lw,i_dn,1)
278  sfcflx_sw_up_c(i,j) = flux_rad(ks-1,i,j,i_sw,i_up,1)
279  sfcflx_sw_dn_c(i,j) = flux_rad(ks-1,i,j,i_sw,i_dn,1)
280  ! for all-sky
281  sfcflx_lw_up(i,j) = flux_rad(ks-1,i,j,i_lw,i_up,2)
282  sfcflx_lw_dn(i,j) = flux_rad(ks-1,i,j,i_lw,i_dn,2)
283  sfcflx_sw_up(i,j) = flux_rad(ks-1,i,j,i_sw,i_up,2)
284  sfcflx_sw_dn(i,j) = flux_rad(ks-1,i,j,i_sw,i_dn,2)
285 
286  flux_net_sfc(i,j,i_lw) = sfcflx_lw_up(i,j) - sfcflx_lw_dn(i,j)
287  flux_net_sfc(i,j,i_sw) = sfcflx_sw_up(i,j) - sfcflx_sw_dn(i,j)
288  enddo
289  enddo
290 
291 !OCL XFILL
292  do j = js, je
293  do i = is, ie
294  ! for clear-sky
295  toaflx_lw_up_c(i,j) = flux_rad_top(i,j,i_lw,i_up,1)
296  toaflx_lw_dn_c(i,j) = flux_rad_top(i,j,i_lw,i_dn,1)
297  toaflx_sw_up_c(i,j) = flux_rad_top(i,j,i_sw,i_up,1)
298  toaflx_sw_dn_c(i,j) = flux_rad_top(i,j,i_sw,i_dn,1)
299  ! for all-sky
300  toaflx_lw_up(i,j) = flux_rad_top(i,j,i_lw,i_up,2)
301  toaflx_lw_dn(i,j) = flux_rad_top(i,j,i_lw,i_dn,2)
302  toaflx_sw_up(i,j) = flux_rad_top(i,j,i_sw,i_up,2)
303  toaflx_sw_dn(i,j) = flux_rad_top(i,j,i_sw,i_dn,2)
304 
305  flux_net_toa(i,j,i_lw) = toaflx_lw_up(i,j) - toaflx_lw_dn(i,j)
306  flux_net_toa(i,j,i_sw) = toaflx_sw_up(i,j) - toaflx_sw_dn(i,j)
307  enddo
308  enddo
309 
310 !OCL XFILL
311  do j = js, je
312  do i = is, ie
313  do k = ks, ke
314  flux_up(k,i,j,i_lw) = 0.5_rp * ( flux_rad(k-1,i,j,i_lw,i_up,2) + flux_rad(k,i,j,i_lw,i_up,2) )
315  flux_dn(k,i,j,i_lw) = 0.5_rp * ( flux_rad(k-1,i,j,i_lw,i_dn,2) + flux_rad(k,i,j,i_lw,i_dn,2) )
316  flux_up(k,i,j,i_sw) = 0.5_rp * ( flux_rad(k-1,i,j,i_sw,i_up,2) + flux_rad(k,i,j,i_sw,i_up,2) )
317  flux_dn(k,i,j,i_sw) = 0.5_rp * ( flux_rad(k-1,i,j,i_sw,i_dn,2) + flux_rad(k,i,j,i_sw,i_dn,2) )
318 
319  flux_net(k,i,j,i_lw) = flux_up(k,i,j,i_lw) - flux_dn(k,i,j,i_lw)
320  flux_net(k,i,j,i_sw) = flux_up(k,i,j,i_sw) - flux_dn(k,i,j,i_sw)
321  enddo
322  enddo
323  enddo
324 
325  if ( atmos_phy_rd_type == 'WRF' ) then
326 
327  flux_rad_org(:,:,:,:,:,:) = flux_rad(:,:,:,:,:,:)
328  rthratensw = 0.0_rp
329  sdown3d = 0.0_rp
330  gsw = 0.0_rp
331 
332  if ( i_qv > 0 ) then
333  qv3d(:,:,:) = qtrc(:,:,:,i_qv)
334  else
335  qv3d = 0.0_rp
336  end if
337  if ( i_qc > 0 ) then
338  qc3d(:,:,:) = qtrc(:,:,:,i_qc)
339  else
340  qc3d = 0.0_rp
341  end if
342  if ( i_qr > 0 ) then
343  qr3d(:,:,:) = qtrc(:,:,:,i_qr)
344  else
345  qr3d = 0.0_rp
346  end if
347  if ( i_qi > 0 ) then
348  qi3d(:,:,:) = qtrc(:,:,:,i_qi)
349  else
350  qi3d = 0.0_rp
351  end if
352  if ( i_qs > 0 ) then
353  qs3d(:,:,:) = qtrc(:,:,:,i_qs)
354  else
355  qs3d = 0.0_rp
356  end if
357  if ( i_qg > 0 ) then
358  qg3d(:,:,:) = qtrc(:,:,:,i_qg)
359  else
360  qg3d = 0.0_rp
361  end if
362 
363  do j = 1, ja
364  do i = 1, ia
365  do k = 1, ka
366  t3d(i,k,j) = temp(k,i,j) ! temperature (K)
367  rho3d(i,k,j) = dens(k,i,j) ! density (kg/m^3)
368  p3d(i,k,j) = pres(k,i,j) ! pressure (Pa)
369  pi3d(i,k,j) = (pres(k,i,j)/pre00)**(rdry/cpdry) ! exner function (dimensionless)
370  dz8w(i,k,j) = real_fz(k,i,j)-real_fz(k-1,i,j) ! dz between full levels(m)
371  enddo
372  enddo
373  enddo
374 
375  call swrad( dt_rd, & ! [IN]
376  rthratensw, & ! [INOUT]
377  sdown3d, & ! [INOUT]
378  gsw, & ! [INOUT]
379  real_lat, & ! [IN]
380  real_lon, & ! [IN]
381  sfc_albedo(:,:,i_sw), & ! [IN]
382  rho3d, & ! [IN]
383  t3d, & ! [IN]
384  p3d, & ! [IN]
385  pi3d, & ! [IN]
386  dz8w, & ! [IN]
387  solins(:,:), & ! [IN]
388  cossza(:,:), & ! [IN]
389  qv3d, & ! [IN]
390  qc3d, & ! [IN]
391  qr3d, & ! [IN]
392  qi3d, & ! [IN]
393  qs3d, & ! [IN]
394  qg3d, & ! [IN]
395  f_qv = .true., & ! [IN]
396  f_qc = .true., & ! [IN]
397  f_qr = .true., & ! [IN]
398  f_qi = .true., & ! [IN]
399  f_qs = .true., & ! [IN]
400  f_qg = .true., & ! [IN]
401  icloud = 1, & ! [IN]
402  warm_rain = .true. ) ! [IN]
403 
404  do j = js, je
405  do i = is, ie
406  flux_net_sfc(i,j,i_sw) = gsw(i,j)
407  do k = ks-1, ke
408  flux_rad(k,i,j,i_sw,i_up,2) = 0.0_rp
409  flux_rad(k,i,j,i_sw,i_dn,2) = sdown3d(i,k,j)
410  enddo
411 
412  do k = 1, ks-2
413  flux_rad(k,i,j,i_sw,i_dn,2) = sdown3d(i,ks-1,j)
414  enddo
415 
416  do k = ke+1, ka
417  flux_rad(k,i,j,i_sw,i_dn,2) = sdown3d(i,ke,j)
418  enddo
419  enddo
420  enddo
421 
422  do j = js, je
423  do i = is, ie
424  sfcflx_sw_up(i,j) = flux_rad(ks-1,i,j,i_sw,i_dn,2) * sfc_albedo(i,j,i_sw)
425  sfcflx_sw_dn(i,j) = flux_rad(ks-1,i,j,i_sw,i_dn,2)
426 
427  toaflx_sw_up(i,j) = 0.0_rp
428  toaflx_sw_dn(i,j) = flux_rad(ke,i,j,i_sw,i_dn,2) ! sometimes TOA altitude is very low
429  enddo
430  enddo
431 
432  do j = js, je
433  do i = is, ie
434  flux_net_toa(i,j,i_sw) = 0.0_rp
435  do k = ks, ke
436  flux_net(k,i,j,i_sw) = 0.0_rp
437  flux_up(k,i,j,i_sw) = 0.0_rp
438  flux_dn(k,i,j,i_sw) = 0.5_rp * ( flux_rad(k-1,i,j,i_sw,i_dn,2) + flux_rad(k,i,j,i_sw,i_dn,2) )
439  enddo
440  enddo
441  enddo
442 
443  do j = js, je
444  do i = is, ie
445  flux_net_sfc(i,j,i_sw) = flux_rad(ks-1,i,j,i_sw,i_dn,2)*sfc_albedo(i,j,i_sw)-flux_rad(ks-1,i,j,i_sw,i_dn,2)
446  enddo
447  enddo
448 
449  endif
450 
451  call hist_in( solins(:,:), 'SOLINS', 'solar insolation', 'W/m2', nohalo=.true. )
452  call hist_in( cossza(:,:), 'COSZ', 'cos(solar zenith angle)', '1', nohalo=.true. )
453 
454  call hist_in( sfcflx_lw_up_c(:,:), 'SFLX_LW_up_c', 'SFC upward longwave radiation flux (clr)', 'W/m2', nohalo=.true. )
455  call hist_in( sfcflx_lw_dn_c(:,:), 'SFLX_LW_dn_c', 'SFC downward longwave radiation flux (clr)', 'W/m2', nohalo=.true. )
456  call hist_in( sfcflx_sw_up_c(:,:), 'SFLX_SW_up_c', 'SFC upward shortwave radiation flux (clr)', 'W/m2', nohalo=.true. )
457  call hist_in( sfcflx_sw_dn_c(:,:), 'SFLX_SW_dn_c', 'SFC downward shortwave radiation flux (clr)', 'W/m2', nohalo=.true. )
458 
459  call hist_in( sfcflx_lw_up(:,:), 'SFLX_LW_up', 'SFC upward longwave radiation flux', 'W/m2', nohalo=.true. )
460  call hist_in( sfcflx_lw_dn(:,:), 'SFLX_LW_dn', 'SFC downward longwave radiation flux', 'W/m2', nohalo=.true. )
461  call hist_in( sfcflx_sw_up(:,:), 'SFLX_SW_up', 'SFC upward shortwave radiation flux', 'W/m2', nohalo=.true. )
462  call hist_in( sfcflx_sw_dn(:,:), 'SFLX_SW_dn', 'SFC downward shortwave radiation flux', 'W/m2', nohalo=.true. )
463 
464  call hist_in( toaflx_lw_up_c(:,:), 'TOAFLX_LW_up_c', 'TOA upward longwave radiation flux (clr)', 'W/m2', nohalo=.true. )
465  call hist_in( toaflx_lw_dn_c(:,:), 'TOAFLX_LW_dn_c', 'TOA downward longwave radiation flux (clr)', 'W/m2', nohalo=.true. )
466  call hist_in( toaflx_sw_up_c(:,:), 'TOAFLX_SW_up_c', 'TOA upward shortwave radiation flux (clr)', 'W/m2', nohalo=.true. )
467  call hist_in( toaflx_sw_dn_c(:,:), 'TOAFLX_SW_dn_c', 'TOA downward shortwave radiation flux (clr)', 'W/m2', nohalo=.true. )
468 
469  call hist_in( toaflx_lw_up(:,:), 'TOAFLX_LW_up', 'TOA upward longwave radiation flux', 'W/m2', nohalo=.true. )
470  call hist_in( toaflx_lw_dn(:,:), 'TOAFLX_LW_dn', 'TOA downward longwave radiation flux', 'W/m2', nohalo=.true. )
471  call hist_in( toaflx_sw_up(:,:), 'TOAFLX_SW_up', 'TOA upward shortwave radiation flux', 'W/m2', nohalo=.true. )
472  call hist_in( toaflx_sw_dn(:,:), 'TOAFLX_SW_dn', 'TOA downward shortwave radiation flux', 'W/m2', nohalo=.true. )
473 
474  call hist_in( flux_net_sfc(:,:,i_lw), 'SLR', 'SFC net longwave radiation flux', 'W/m2', nohalo=.true. )
475  call hist_in( flux_net_sfc(:,:,i_sw), 'SSR', 'SFC net shortwave radiation flux', 'W/m2', nohalo=.true. )
476  call hist_in( flux_net_toa(:,:,i_lw), 'OLR', 'TOA net longwave radiation flux', 'W/m2', nohalo=.true. )
477  call hist_in( flux_net_toa(:,:,i_sw), 'OSR', 'TOA net shortwave radiation flux', 'W/m2', nohalo=.true. )
478 
479  call hist_in( flux_up(:,:,:,i_lw), 'RADFLUX_LWUP', 'upward longwave radiation flux', 'W/m2', nohalo=.true. )
480  call hist_in( flux_dn(:,:,:,i_lw), 'RADFLUX_LWDN', 'downward longwave radiation flux', 'W/m2', nohalo=.true. )
481  call hist_in( flux_net(:,:,:,i_lw), 'RADFLUX_LW', 'net longwave radiation flux', 'W/m2', nohalo=.true. )
482  call hist_in( flux_up(:,:,:,i_sw), 'RADFLUX_SWUP', 'upward shortwave radiation flux', 'W/m2', nohalo=.true. )
483  call hist_in( flux_dn(:,:,:,i_sw), 'RADFLUX_SWDN', 'downward shortwave radiation flux', 'W/m2', nohalo=.true. )
484  call hist_in( flux_net(:,:,:,i_sw), 'RADFLUX_SW', 'net shortwave radiation flux', 'W/m2', nohalo=.true. )
485 
486  call hist_in( sflx_rad_dn(:,:,i_lw,i_direct ), 'SFLX_LW_dn_dir', 'SFC downward longwave flux (direct )', 'W/m2', nohalo=.true. )
487  call hist_in( sflx_rad_dn(:,:,i_lw,i_diffuse), 'SFLX_LW_dn_dif', 'SFC downward longwave flux (diffuse)', 'W/m2', nohalo=.true. )
488  call hist_in( sflx_rad_dn(:,:,i_sw,i_direct ), 'SFLX_SW_dn_dir', 'SFC downward shortwave flux (direct )', 'W/m2', nohalo=.true. )
489  call hist_in( sflx_rad_dn(:,:,i_sw,i_diffuse), 'SFLX_SW_dn_dif', 'SFC downward shortwave flux (diffuse)', 'W/m2', nohalo=.true. )
490 
491  call hist_in( temp_t(:,:,:,i_lw), 'TEMP_t_rd_LW', 'tendency of temp in rd(LW)', 'K/day', nohalo=.true. )
492  call hist_in( temp_t(:,:,:,i_sw), 'TEMP_t_rd_SW', 'tendency of temp in rd(SW)', 'K/day', nohalo=.true. )
493  call hist_in( temp_t(:,:,:,3 ), 'TEMP_t_rd', 'tendency of temp in rd', 'K/day', nohalo=.true. )
494  call hist_in( rhot_t_rd(:,:,:), 'RHOT_t_RD', 'tendency of RHOT in rd', 'K.kg/m3/s', nohalo=.true. )
495 
496  ! output of raw data, for offline output
497  call hist_in( flux_rad(:,:,:,i_lw,i_up,2), 'RFLX_LW_up', 'upward longwave radiation flux (cell face)', 'W/m2', nohalo=.true. )
498  call hist_in( flux_rad(:,:,:,i_lw,i_dn,2), 'RFLX_LW_dn', 'downward longwave radiation flux (cell face)', 'W/m2', nohalo=.true. )
499  call hist_in( flux_rad(:,:,:,i_sw,i_up,2), 'RFLX_SW_up', 'upward shortwave radiation flux (cell face)', 'W/m2', nohalo=.true. )
500  call hist_in( flux_rad(:,:,:,i_sw,i_dn,2), 'RFLX_SW_dn', 'downward shortwave radiation flux (cell face)', 'W/m2', nohalo=.true. )
501 
502  if ( atmos_phy_rd_type == 'WRF' ) then
503  ! revert all radiation flux from MM5 scheme to default
504  flux_rad(:,:,:,:,:,:) = flux_rad_org(:,:,:,:,:,:)
505 
506  do j = js, je
507  do i = is, ie
508  sfcflx_sw_up(i,j) = flux_rad(ks-1,i,j,i_sw,i_up,2)
509  sfcflx_sw_dn(i,j) = flux_rad(ks-1,i,j,i_sw,i_dn,2)
510 
511  toaflx_sw_up(i,j) = flux_rad_top(i,j,i_sw,i_up,2) ! mstrnx
512  toaflx_sw_dn(i,j) = flux_rad_top(i,j,i_sw,i_dn,2) ! mstrnx
513  enddo
514  enddo
515 
516  endif
517 
518  endif
519 
520  do j = js, je
521  do i = is, ie
522  do k = ks, ke
523  rhot_t(k,i,j) = rhot_t(k,i,j) + rhot_t_rd(k,i,j)
524  enddo
525  enddo
526  enddo
527 
528  if ( statistics_checktotal ) then
529  call stat_total( total, rhot_t_rd(:,:,:), 'RHOT_t_RD' )
530  endif
531 
532  return
module ATMOS admin
logical, public statistics_checktotal
calc&report variable totals to logfile?
subroutine, public atmos_phy_rd_heating(flux_rad, DENS, RHOT, QTRC, FZ, dt, TEMP_t, RHOT_t)
Calc heating rate.
real(rp), public const_cpdry
specific heat (dry air,constant pressure) [J/kg/K]
Definition: scale_const.F90:58
subroutine, public prc_mpistop
Abort MPI.
integer, parameter, public i_direct
real(rp), dimension(:,:), allocatable, public landuse_fact_urban
urban factor
real(dp), public time_dtsec_atmos_phy_rd
time interval of physics(radiation ) [sec]
Definition: scale_time.F90:42
real(rp), dimension(:,:,:), allocatable, target, public rhot
procedure(rd), pointer, public atmos_phy_rd
real(rp), dimension(:,:,:), allocatable, public pres
module ATMOSPHERIC Variables
real(rp), dimension(:,:,:,:), pointer, public qtrc_av
real(rp), dimension(:,:,:), allocatable, public rhot_tp
real(rp), dimension(:,:,:), allocatable, public atmos_phy_rd_rhot_t
real(rp), dimension(:,:,:), allocatable, public real_fz
geopotential height [m] (cell face )
real(rp), dimension(:,:), allocatable, public atmos_phy_sf_sfc_temp
real(rp), dimension(:,:,:), allocatable, public real_cz
geopotential height [m] (cell center)
real(rp), dimension(:,:,:), allocatable, target, public dens
integer, parameter, public i_diffuse
character(len=h_short), public atmos_phy_rd_type
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_sflx_lw_up
integer, parameter, public i_lw
real(rp), public const_rdry
specific gas constant (dry air) [J/kg/K]
Definition: scale_const.F90:57
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_toaflx_sw_up
module Statistics
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_sflx_lw_dn
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_cossza
module Atmosphere / Physics Radiation
module ATMOSPHERIC Surface Variables
integer, parameter, public i_sw
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_sflx_sw_up
real(rp), dimension(:,:,:), allocatable, public temp
module GRID (real space)
module LANDUSE
integer, parameter, public i_dn
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_sflx_sw_dn
real(rp), public const_pre00
pressure reference [Pa]
Definition: scale_const.F90:90
integer, public time_offset_year
time offset [year]
Definition: scale_time.F90:73
real(rp), dimension(:,:,:), allocatable, public atmos_phy_sf_sfc_albedo
real(rp), dimension(:,:), allocatable, public landuse_fact_ocean
ocean factor
module TIME
Definition: scale_time.F90:15
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_toaflx_lw_dn
module ATMOSPHERE / Physics Radiation
module PROCESS
real(rp), dimension(:,:,:), pointer, public dens_av
module CONSTANT
Definition: scale_const.F90:14
real(rp), dimension(:,:,:,:), allocatable, public atmos_phy_rd_sflx_downall
subroutine, public swrad(dt, RTHRATEN, SDOWN3D, GSW, XLAT, XLONG, ALBEDO, rho_phy, T3D, P3D, pi3D, dz8w, solins, cosSZA, QV3D, QC3D, QR3D, QI3D, QS3D, QG3D, F_QV, F_QC, F_QR, F_QI, F_QS, F_QG, icloud, warm_rain)
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_toaflx_sw_dn
real(rp), dimension(:,:), allocatable, public real_lon
longitude [rad,0-2pi]
module HISTORY
module ATMOSPHERE / Physics Radiation
integer, dimension(6), public time_nowdate
current time [YYYY MM DD HH MM SS]
Definition: scale_time.F90:65
real(rp), dimension(:,:,:), pointer, public rhot_av
real(rp), dimension(:,:), allocatable, public real_lat
latitude [rad,-pi,pi]
real(rp), dimension(:,:), allocatable, public landuse_fact_land
land factor
module ATMOSPHERE / Physics Radiation
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_solins
integer, parameter, public i_up
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_toaflx_lw_up
real(rp), dimension(:,:,:,:), allocatable, target, public qtrc
Here is the call graph for this function:
Here is the caller graph for this function: