source: LMDZ6/trunk/libf/phylmd/ecrad/radiation_interface.F90 @ 4013

Last change on this file since 4013 was 3908, checked in by idelkadi, 4 years ago

Online implementation of the radiative transfer code ECRAD in the LMDZ model.

  • Inclusion of the ecrad directory containing the sources of the ECRAD code
    • interface routine : radiation_scheme.F90
  • Adaptation of compilation scripts :
    • compilation under CPP key CPP_ECRAD
    • compilation with option "-rad ecard" or "-ecard true"
    • The "-rad old/rtm/ecran" build option will need to replace the "-rrtm true" and "-ecrad true" options in the future.
  • Runing LMDZ simulations with ecrad, you need :
    • logical key iflag_rrtm = 2 in physiq.def
    • namelist_ecrad (DefLists?)
    • the directory "data" containing the configuration files is temporarily placed in ../libfphylmd/ecrad/
  • Compilation and execution are tested in the 1D case. The repository under svn would allow to continue the implementation work: tests, verification of the results, ...
File size: 26.6 KB
Line 
1! radiation_interface.F90 - Public interface to radiation scheme
2!
3! (C) Copyright 2014- ECMWF.
4!
5! This software is licensed under the terms of the Apache Licence Version 2.0
6! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
7!
8! In applying this licence, ECMWF does not waive the privileges and immunities
9! granted to it by virtue of its status as an intergovernmental organisation
10! nor does it submit to any jurisdiction.
11!
12! Author:  Robin Hogan
13! Email:   r.j.hogan@ecmwf.int
14!
15! Modifications
16!   2017-04-11  R. Hogan  Changes to enable generalized surface description
17!   2017-09-08  R. Hogan  Reverted some changes
18!
19! To use the radiation scheme, create a configuration_type object,
20! call "setup_radiation" on it once to load the look-up-tables and
21! data describing how gas and hydrometeor absorption/scattering are to
22! be represented, and call "radiation" multiple times on different
23! input profiles.
24
25module radiation_interface
26
27  implicit none
28
29  public  :: setup_radiation, set_gas_units, radiation
30  private :: radiation_reverse
31
32contains
33
34  !---------------------------------------------------------------------
35  ! Load the look-up-tables and data describing how gas and
36  ! hydrometeor absorption/scattering are to be represented
37  subroutine setup_radiation(config)
38
39    use parkind1,         only : jprb
40    use yomhook,          only : lhook, dr_hook
41    use radiation_config, only : config_type, ISolverMcICA, &
42         &   IGasModelMonochromatic, IGasModelIFSRRTMG
43
44    ! Currently there are two gas absorption models: RRTMG (default)
45    ! and monochromatic
46    use radiation_monochromatic,  only : &
47         &   setup_gas_optics_mono     => setup_gas_optics, &
48         &   setup_cloud_optics_mono   => setup_cloud_optics, &
49         &   setup_aerosol_optics_mono => setup_aerosol_optics
50    use radiation_ifs_rrtm,       only :  setup_gas_optics
51    use radiation_cloud_optics,   only :  setup_cloud_optics
52    use radiation_aerosol_optics, only :  setup_aerosol_optics
53
54    type(config_type), intent(inout) :: config
55
56    real(jprb) :: hook_handle
57
58    if (lhook) call dr_hook('radiation_interface:setup_radiation',0,hook_handle)
59
60    ! Consolidate configuration data, including setting data file
61    ! names
62    call config%consolidate()
63
64    ! Load the look-up tables from files in the specified directory
65    if (config%i_gas_model == IGasModelMonochromatic) then
66      call setup_gas_optics_mono(config, trim(config%directory_name))
67    else if (config%i_gas_model == IGasModelIFSRRTMG) then
68      call setup_gas_optics(config, trim(config%directory_name))
69    end if
70
71    ! Whether or not the "radiation" subroutine needs ssa_lw and g_lw
72    ! arrays depends on whether longwave scattering by aerosols is to
73    ! be included.  If not, one of the array dimensions will be set to
74    ! zero.
75    if (config%do_lw_aerosol_scattering) then
76      config%n_g_lw_if_scattering = config%n_g_lw
77    else
78      config%n_g_lw_if_scattering = 0
79    end if
80
81    ! Whether or not the "radiation" subroutine needs ssa_lw_cloud and
82    ! g_lw_cloud arrays depends on whether longwave scattering by
83    ! hydrometeors is to be included.  If not, one of the array
84    ! dimensions will be set to zero.
85    if (config%do_lw_cloud_scattering) then
86      config%n_bands_lw_if_scattering = config%n_bands_lw
87    else
88      config%n_bands_lw_if_scattering = 0
89    end if
90
91    ! If we have longwave scattering and McICA then even if there is
92    ! no aerosol, it is convenient if single scattering albedo and
93    ! g factor arrays are allocated before the call to
94    ! solver_lw as they will be needed.
95    if (config%do_lw_cloud_scattering &
96         &  .and. config%i_solver_lw == ISolverMcICA) then
97      config%n_g_lw_if_scattering = config%n_g_lw
98    end if
99
100    ! Consolidate the albedo/emissivity intervals with the shortwave
101    ! and longwave spectral bands
102    call config%consolidate_intervals(.true., &
103           &  config%do_nearest_spectral_sw_albedo, &
104           &  config%sw_albedo_wavelength_bound, config%i_sw_albedo_index, &
105           &  config%wavenumber1_sw, config%wavenumber2_sw, &
106           &  config%i_albedo_from_band_sw, config%sw_albedo_weights)
107    call config%consolidate_intervals(.false., &
108           &  config%do_nearest_spectral_lw_emiss, &
109           &  config%lw_emiss_wavelength_bound, config%i_lw_emiss_index, &
110           &  config%wavenumber1_lw, config%wavenumber2_lw, &
111           &  config%i_emiss_from_band_lw, config%lw_emiss_weights)
112
113    if (config%do_clouds) then
114      if (config%i_gas_model == IGasModelMonochromatic) then
115        !      call setup_cloud_optics_mono(config)
116      else
117        call setup_cloud_optics(config)
118      end if
119    end if
120
121    if (config%use_aerosols) then
122      if (config%i_gas_model == IGasModelMonochromatic) then
123!        call setup_aerosol_optics_mono(config)
124      else
125        call setup_aerosol_optics(config)
126      end if
127    end if
128
129    ! Load cloud water PDF look-up table for McICA
130    if (         config%i_solver_sw == ISolverMcICA &
131         &  .or. config%i_solver_lw == ISolverMcICA) then
132      call config%pdf_sampler%setup(config%cloud_pdf_file_name, &
133           &                        iverbose=config%iverbosesetup)
134    end if
135
136    if (lhook) call dr_hook('radiation_interface:setup_radiation',1,hook_handle)
137
138  end subroutine setup_radiation
139
140
141  !---------------------------------------------------------------------
142  ! Scale the gas mixing ratios so that they have the units (and
143  ! possibly scale factors) required by the specific gas absorption
144  ! model.  This subroutine simply passes the gas object on to the
145  ! module of the currently active gas model.
146  subroutine set_gas_units(config, gas)
147   
148    use radiation_config
149    use radiation_gas,           only : gas_type
150    use radiation_monochromatic, only : set_gas_units_mono  => set_gas_units
151    use radiation_ifs_rrtm,      only : set_gas_units_ifs   => set_gas_units
152
153    type(config_type), intent(in)    :: config
154    type(gas_type),    intent(inout) :: gas
155
156    if (config%i_gas_model == IGasModelMonochromatic) then
157      call set_gas_units_mono(gas)
158    else
159      call set_gas_units_ifs(gas)
160    end if
161
162  end subroutine set_gas_units
163
164
165  !---------------------------------------------------------------------
166  ! Run the radiation scheme according to the configuration in the
167  ! config object. There are ncol profiles of which only istartcol to
168  ! iendcol are to be processed, and there are nlev model levels.  The
169  ! output fluxes are written to the flux object, and all other
170  ! objects contain the input variables.  The variables may be defined
171  ! either in order of increasing or decreasing pressure, but if in
172  ! order of decreasing pressure then radiation_reverse will be called
173  ! to reverse the order for the computation and then reverse the
174  ! order of the output fluxes to match the inputs.
175  subroutine radiation(ncol, nlev, istartcol, iendcol, config, &
176       &  single_level, thermodynamics, gas, cloud, aerosol, flux)
177
178    use parkind1,                 only : jprb
179    use yomhook,                  only : lhook, dr_hook
180
181    use radiation_io,             only : nulout
182    use radiation_config,         only : config_type, &
183         &   IGasModelMonochromatic, IGasModelIFSRRTMG, &
184         &   ISolverMcICA, ISolverSpartacus, ISolverHomogeneous, &
185         &   ISolverTripleclouds
186    use radiation_single_level,   only : single_level_type
187    use radiation_thermodynamics, only : thermodynamics_type
188    use radiation_gas,            only : gas_type
189    use radiation_cloud,          only : cloud_type
190    use radiation_aerosol,        only : aerosol_type
191    use radiation_flux,           only : flux_type
192    use radiation_spartacus_sw,   only : solver_spartacus_sw
193    use radiation_spartacus_lw,   only : solver_spartacus_lw
194    use radiation_tripleclouds_sw,only : solver_tripleclouds_sw
195    use radiation_tripleclouds_lw,only : solver_tripleclouds_lw
196    use radiation_mcica_sw,       only : solver_mcica_sw
197    use radiation_mcica_lw,       only : solver_mcica_lw
198    use radiation_cloudless_sw,   only : solver_cloudless_sw
199    use radiation_cloudless_lw,   only : solver_cloudless_lw
200    use radiation_homogeneous_sw, only : solver_homogeneous_sw
201    use radiation_homogeneous_lw, only : solver_homogeneous_lw
202    use radiation_save,           only : save_radiative_properties
203
204    ! Treatment of gas and hydrometeor optics
205    use radiation_monochromatic,  only : &
206         &   gas_optics_mono         => gas_optics, &
207         &   cloud_optics_mono       => cloud_optics, &
208         &   add_aerosol_optics_mono => add_aerosol_optics
209    use radiation_ifs_rrtm,       only : gas_optics
210    use radiation_cloud_optics,   only : cloud_optics
211    use radiation_aerosol_optics, only : add_aerosol_optics
212
213    ! Inputs
214    integer, intent(in) :: ncol               ! number of columns
215    integer, intent(in) :: nlev               ! number of model levels
216    integer, intent(in) :: istartcol, iendcol ! range of columns to process
217    type(config_type),        intent(in)   :: config
218    type(single_level_type),  intent(in)   :: single_level
219    type(thermodynamics_type),intent(in)   :: thermodynamics
220    type(gas_type),           intent(in)   :: gas
221    type(cloud_type),         intent(inout):: cloud
222    type(aerosol_type),       intent(in)   :: aerosol
223    ! Output
224    type(flux_type),          intent(inout):: flux
225
226
227    ! Local variables
228
229    ! Layer optical depth, single scattering albedo and asymmetry factor of
230    ! gases and aerosols at each longwave g-point, where the latter
231    ! two variables are only defined if aerosol longwave scattering is
232    ! enabled (otherwise both are treated as zero).
233    real(jprb), dimension(config%n_g_lw,nlev,istartcol:iendcol) :: od_lw
234    real(jprb), dimension(config%n_g_lw_if_scattering,nlev,istartcol:iendcol) :: &
235         &  ssa_lw, g_lw
236
237    ! Layer in-cloud optical depth, single scattering albedo and
238    ! asymmetry factor of hydrometeors in each longwave band, where
239    ! the latter two variables are only defined if hydrometeor
240    ! longwave scattering is enabled (otherwise both are treated as
241    ! zero).
242    real(jprb), dimension(config%n_bands_lw,nlev,istartcol:iendcol) :: od_lw_cloud
243    real(jprb), dimension(config%n_bands_lw_if_scattering,nlev,istartcol:iendcol) :: &
244         &  ssa_lw_cloud, g_lw_cloud
245
246    ! Layer optical depth, single scattering albedo and asymmetry factor of
247    ! gases and aerosols at each shortwave g-point
248    real(jprb), dimension(config%n_g_sw,nlev,istartcol:iendcol) :: od_sw, ssa_sw, g_sw
249
250    ! Layer in-cloud optical depth, single scattering albedo and
251    ! asymmetry factor of hydrometeors in each shortwave band
252    real(jprb), dimension(config%n_bands_sw,nlev,istartcol:iendcol)   :: &
253         &  od_sw_cloud, ssa_sw_cloud, g_sw_cloud
254
255    ! The Planck function (emitted flux from a black body) at half
256    ! levels
257    real(jprb), dimension(config%n_g_lw,nlev+1,istartcol:iendcol) :: planck_hl
258
259    ! The longwave emission from and albedo of the surface in each
260    ! longwave g-point; note that these are weighted averages of the
261    ! values from individual tiles
262    real(jprb), dimension(config%n_g_lw, istartcol:iendcol) :: lw_emission
263    real(jprb), dimension(config%n_g_lw, istartcol:iendcol) :: lw_albedo
264
265    ! Direct and diffuse shortwave surface albedo in each shortwave
266    ! g-point; note that these are weighted averages of the values
267    ! from individual tiles
268    real(jprb), dimension(config%n_g_sw, istartcol:iendcol) :: sw_albedo_direct
269    real(jprb), dimension(config%n_g_sw, istartcol:iendcol) :: sw_albedo_diffuse
270
271    ! The incoming shortwave flux into a plane perpendicular to the
272    ! incoming radiation at top-of-atmosphere in each of the shortwave
273    ! g-points
274    real(jprb), dimension(config%n_g_sw,istartcol:iendcol) :: incoming_sw
275
276    character(len=100) :: rad_prop_file_name
277    character(*), parameter :: rad_prop_base_file_name = "radiative_properties"
278
279    real(jprb) :: hook_handle
280
281    if (lhook) call dr_hook('radiation_interface:radiation',0,hook_handle)
282
283    if (thermodynamics%pressure_hl(istartcol,2) &
284         &  < thermodynamics%pressure_hl(istartcol,1)) then
285      ! Input arrays are arranged in order of decreasing pressure /
286      ! increasing height: the following subroutine reverses them,
287      ! call the radiation scheme and then reverses the returned
288      ! fluxes
289      call radiation_reverse(ncol, nlev, istartcol, iendcol, config, &
290           &  single_level, thermodynamics, gas, cloud, aerosol, flux)
291    else
292
293      ! Input arrays arranged in order of increasing pressure /
294      ! decreasing height: progress normally
295
296      ! Extract surface albedos at each gridpoint
297      call single_level%get_albedos(istartcol, iendcol, config, &
298           &                        sw_albedo_direct, sw_albedo_diffuse, &
299           &                        lw_albedo)
300
301      ! Compute gas absorption optical depth in shortwave and
302      ! longwave, shortwave single scattering albedo (i.e. fraction of
303      ! extinction due to Rayleigh scattering), Planck functions and
304      ! incoming shortwave flux at each g-point, for the specified
305      ! range of atmospheric columns
306      if (config%i_gas_model == IGasModelMonochromatic) then
307        call gas_optics_mono(ncol,nlev,istartcol,iendcol, config, &
308             &  single_level, thermodynamics, gas, lw_albedo, &
309             &  od_lw, od_sw, ssa_sw, &
310             &  planck_hl, lw_emission, incoming_sw)
311      else
312        call gas_optics(ncol,nlev,istartcol,iendcol, config, &
313             &  single_level, thermodynamics, gas, &
314             &  od_lw, od_sw, ssa_sw, lw_albedo=lw_albedo, &
315             &  planck_hl=planck_hl, lw_emission=lw_emission, &
316             &  incoming_sw=incoming_sw)
317      end if
318
319      if (config%do_clouds) then
320        ! Crop the cloud fraction to remove clouds that have too small
321        ! a fraction or water content; after this, we can safely
322        ! assume that a cloud is present if cloud%fraction > 0.0.
323        call cloud%crop_cloud_fraction(istartcol, iendcol, &
324             &            config%cloud_fraction_threshold, &
325             &            config%cloud_mixing_ratio_threshold)
326
327        ! Compute hydrometeor absorption/scattering properties in each
328        ! shortwave and longwave band
329        if (config%i_gas_model == IGasModelMonochromatic) then
330          call cloud_optics_mono(nlev, istartcol, iendcol, &
331               &  config, thermodynamics, cloud, &
332               &  od_lw_cloud, ssa_lw_cloud, g_lw_cloud, &
333               &  od_sw_cloud, ssa_sw_cloud, g_sw_cloud)
334        else
335          call cloud_optics(nlev, istartcol, iendcol, &
336               &  config, thermodynamics, cloud, &
337               &  od_lw_cloud, ssa_lw_cloud, g_lw_cloud, &
338               &  od_sw_cloud, ssa_sw_cloud, g_sw_cloud)
339        end if
340      end if ! do_clouds
341
342      if (config%use_aerosols) then
343        if (config%i_gas_model == IGasModelMonochromatic) then
344!          call add_aerosol_optics_mono(nlev,istartcol,iendcol, &
345!               &  config, thermodynamics, gas, aerosol, &
346!               &  od_lw, ssa_lw, g_lw, od_sw, ssa_sw, g_sw)
347        else
348          call add_aerosol_optics(nlev,istartcol,iendcol, &
349               &  config, thermodynamics, gas, aerosol, &
350               &  od_lw, ssa_lw, g_lw, od_sw, ssa_sw, g_sw)
351        end if
352      else
353        g_sw = 0.0_jprb
354        if (config%do_lw_aerosol_scattering) then
355          ssa_lw = 0.0_jprb
356          g_lw   = 0.0_jprb
357        end if
358      end if
359
360      ! For diagnostic purposes, save these intermediate variables to
361      ! a NetCDF file
362      if (config%do_save_radiative_properties) then
363        if (istartcol == 1 .and. iendcol == ncol) then
364          rad_prop_file_name = rad_prop_base_file_name // ".nc"
365        else
366          write(rad_prop_file_name,'(a,a,i4.4,a,i4.4,a)') &
367               &  rad_prop_base_file_name, '_', istartcol, '-',iendcol,'.nc'
368        end if
369        call save_radiative_properties(trim(rad_prop_file_name), &
370             &  nlev, istartcol, iendcol, &
371             &  config, single_level, thermodynamics, cloud, &
372             &  planck_hl, lw_emission, lw_albedo, &
373             &  sw_albedo_direct, sw_albedo_diffuse, incoming_sw, &
374             &  od_lw, ssa_lw, g_lw, od_sw, ssa_sw, g_sw, &
375             &  od_lw_cloud, ssa_lw_cloud, g_lw_cloud, &
376             &  od_sw_cloud, ssa_sw_cloud, g_sw_cloud)
377      end if
378
379      if (config%do_lw) then
380        if (config%iverbose >= 2) then
381          write(nulout,'(a)') 'Computing longwave fluxes'
382        end if
383
384        if (config%i_solver_lw == ISolverMcICA) then
385          ! Compute fluxes using the McICA longwave solver
386          call solver_mcica_lw(nlev,istartcol,iendcol, &
387               &  config, single_level, cloud, &
388               &  od_lw, ssa_lw, g_lw, od_lw_cloud, ssa_lw_cloud, &
389               &  g_lw_cloud, planck_hl, lw_emission, lw_albedo, flux)
390        else if (config%i_solver_lw == ISolverSPARTACUS) then
391          ! Compute fluxes using the SPARTACUS longwave solver
392          call solver_spartacus_lw(nlev,istartcol,iendcol, &
393               &  config, thermodynamics, cloud, &
394               &  od_lw, ssa_lw, g_lw, od_lw_cloud, ssa_lw_cloud, g_lw_cloud, &
395               &  planck_hl, lw_emission, lw_albedo, flux)
396        else if (config%i_solver_lw == ISolverTripleclouds) then
397          ! Compute fluxes using the Tripleclouds longwave solver
398          call solver_tripleclouds_lw(nlev,istartcol,iendcol, &
399               &  config, cloud, &
400               &  od_lw, ssa_lw, g_lw, od_lw_cloud, ssa_lw_cloud, g_lw_cloud, &
401               &  planck_hl, lw_emission, lw_albedo, flux)
402        elseif (config%i_solver_lw == ISolverHomogeneous) then
403          ! Compute fluxes using the homogeneous solver
404          call solver_homogeneous_lw(nlev,istartcol,iendcol, &
405               &  config, cloud, &
406               &  od_lw, ssa_lw, g_lw, od_lw_cloud, ssa_lw_cloud, &
407               &  g_lw_cloud, planck_hl, lw_emission, lw_albedo, flux)
408        else
409          ! Compute fluxes using the cloudless solver
410          call solver_cloudless_lw(nlev,istartcol,iendcol, &
411               &  config, od_lw, ssa_lw, g_lw, &
412               &  planck_hl, lw_emission, lw_albedo, flux)
413        end if
414      end if
415
416      if (config%do_sw) then
417        if (config%iverbose >= 2) then
418          write(nulout,'(a)') 'Computing shortwave fluxes'
419        end if
420
421        if (config%i_solver_sw == ISolverMcICA) then
422          ! Compute fluxes using the McICA shortwave solver
423          call solver_mcica_sw(nlev,istartcol,iendcol, &
424               &  config, single_level, cloud, &
425               &  od_sw, ssa_sw, g_sw, od_sw_cloud, ssa_sw_cloud, &
426               &  g_sw_cloud, sw_albedo_direct, sw_albedo_diffuse, &
427               &  incoming_sw, flux)
428        else if (config%i_solver_sw == ISolverSPARTACUS) then
429          ! Compute fluxes using the SPARTACUS shortwave solver
430          call solver_spartacus_sw(nlev,istartcol,iendcol, &
431               &  config, single_level, thermodynamics, cloud, &
432               &  od_sw, ssa_sw, g_sw, od_sw_cloud, ssa_sw_cloud, &
433               &  g_sw_cloud, sw_albedo_direct, sw_albedo_diffuse, &
434               &  incoming_sw, flux)
435        else if (config%i_solver_sw == ISolverTripleclouds) then
436          ! Compute fluxes using the Tripleclouds shortwave solver
437          call solver_tripleclouds_sw(nlev,istartcol,iendcol, &
438               &  config, single_level, cloud, &
439               &  od_sw, ssa_sw, g_sw, od_sw_cloud, ssa_sw_cloud, &
440               &  g_sw_cloud, sw_albedo_direct, sw_albedo_diffuse, &
441               &  incoming_sw, flux)
442        elseif (config%i_solver_sw == ISolverHomogeneous) then
443          ! Compute fluxes using the homogeneous solver
444          call solver_homogeneous_sw(nlev,istartcol,iendcol, &
445               &  config, single_level, cloud, &
446               &  od_sw, ssa_sw, g_sw, od_sw_cloud, ssa_sw_cloud, &
447               &  g_sw_cloud, sw_albedo_direct, sw_albedo_diffuse, &
448               &  incoming_sw, flux)
449        else
450          ! Compute fluxes using the cloudless solver
451          call solver_cloudless_sw(nlev,istartcol,iendcol, &
452               &  config, single_level, od_sw, ssa_sw, g_sw, &
453               &  sw_albedo_direct, sw_albedo_diffuse, &
454               &  incoming_sw, flux)
455        end if
456      end if
457
458      ! Store surface downwelling fluxes in bands from fluxes in g
459      ! points
460      call flux%calc_surface_spectral(config, istartcol, iendcol)
461
462    end if
463   
464    if (lhook) call dr_hook('radiation_interface:radiation',1,hook_handle)
465
466  end subroutine radiation
467
468
469  !---------------------------------------------------------------------
470  ! If the input arrays are arranged in order of decreasing pressure /
471  ! increasing height then this subroutine reverses them, calls the
472  ! radiation scheme and then reverses the returned fluxes. Since this
473  ! subroutine calls, and is called by "radiation", it must be in this
474  ! module to avoid circular dependencies.
475  subroutine radiation_reverse(ncol, nlev, istartcol, iendcol, config, &
476       &  single_level, thermodynamics, gas, cloud, aerosol, flux)
477 
478    use parkind1, only : jprb
479
480    use radiation_io,             only : nulout
481    use radiation_config,         only : config_type
482    use radiation_single_level,   only : single_level_type
483    use radiation_thermodynamics, only : thermodynamics_type
484    use radiation_gas,            only : gas_type
485    use radiation_cloud,          only : cloud_type
486    use radiation_aerosol,        only : aerosol_type
487    use radiation_flux,           only : flux_type
488
489    ! Inputs
490    integer, intent(in) :: ncol               ! number of columns
491    integer, intent(in) :: nlev               ! number of model levels
492    integer, intent(in) :: istartcol, iendcol ! range of columns to process
493    type(config_type),        intent(in) :: config
494    type(single_level_type),  intent(in) :: single_level
495    type(thermodynamics_type),intent(in) :: thermodynamics
496    type(gas_type),           intent(in) :: gas
497    type(cloud_type),         intent(in) :: cloud
498    type(aerosol_type),       intent(in) :: aerosol
499    ! Output
500    type(flux_type),          intent(inout):: flux
501
502    ! Reversed data structures
503    type(thermodynamics_type) :: thermodynamics_rev
504    type(gas_type)            :: gas_rev
505    type(cloud_type)          :: cloud_rev
506    type(aerosol_type)        :: aerosol_rev
507    type(flux_type)           :: flux_rev
508
509    ! Start and end levels for aerosol data
510    integer :: istartlev, iendlev
511
512    if (config%iverbose >= 2) then
513      write(nulout,'(a)') 'Reversing arrays to be in order of increasing pressure'
514    end if
515
516    ! Allocate reversed arrays
517    call thermodynamics_rev%allocate(ncol, nlev)
518    call cloud_rev%allocate(ncol, nlev)
519    call flux_rev%allocate(config, istartcol, iendcol, nlev)
520    if (allocated(aerosol%mixing_ratio)) then
521      istartlev = nlev + 1 - aerosol%iendlev
522      iendlev   = nlev + 1 - aerosol%istartlev
523      call aerosol_rev%allocate(ncol, istartlev, iendlev, &
524           &                    config%n_aerosol_types)
525    end if
526
527    ! Fill reversed thermodynamic arrays
528    thermodynamics_rev%pressure_hl(istartcol:iendcol,:) &
529         &  = thermodynamics%pressure_hl(istartcol:iendcol, nlev+1:1:-1)
530    thermodynamics_rev%temperature_hl(istartcol:iendcol,:) &
531         &  = thermodynamics%temperature_hl(istartcol:iendcol, nlev+1:1:-1)
532
533    ! Fill reversed gas arrays
534    call gas%reverse(istartcol, iendcol, gas_rev)
535
536    if (config%do_clouds) then
537      ! Fill reversed cloud arrays
538      cloud_rev%q_liq(istartcol:iendcol,:) &
539           &  = cloud%q_liq(istartcol:iendcol,nlev:1:-1)
540      cloud_rev%re_liq(istartcol:iendcol,:) &
541           &  = cloud%re_liq(istartcol:iendcol,nlev:1:-1)
542      cloud_rev%q_ice(istartcol:iendcol,:) &
543           &  = cloud%q_ice(istartcol:iendcol,nlev:1:-1)
544      cloud_rev%re_ice(istartcol:iendcol,:) &
545           &  = cloud%re_ice(istartcol:iendcol,nlev:1:-1)
546      cloud_rev%fraction(istartcol:iendcol,:) &
547           &  = cloud%fraction(istartcol:iendcol,nlev:1:-1)
548      cloud_rev%overlap_param(istartcol:iendcol,:) &
549           &  = cloud%overlap_param(istartcol:iendcol,nlev-1:1:-1)
550      if (allocated(cloud%fractional_std)) then
551        cloud_rev%fractional_std(istartcol:iendcol,:) &
552             &  = cloud%fractional_std(istartcol:iendcol,nlev:1:-1)
553      else
554        cloud_rev%fractional_std(istartcol:iendcol,:) = 0.0_jprb       
555      end if
556      if (allocated(cloud%inv_cloud_effective_size)) then
557        cloud_rev%inv_cloud_effective_size(istartcol:iendcol,:) &
558             &  = cloud%inv_cloud_effective_size(istartcol:iendcol,nlev:1:-1)
559      else
560        cloud_rev%inv_cloud_effective_size(istartcol:iendcol,:) = 0.0_jprb
561      end if
562    end if
563
564    if (allocated(aerosol%mixing_ratio)) then
565      aerosol_rev%mixing_ratio(:,istartlev:iendlev,:) &
566           &  = aerosol%mixing_ratio(:,aerosol%iendlev:aerosol%istartlev:-1,:)
567    end if
568
569    ! Run radiation scheme on reversed profiles
570    call radiation(ncol, nlev,istartcol,iendcol, &
571         &  config, single_level, thermodynamics_rev, gas_rev, &
572         &  cloud_rev, aerosol_rev, flux_rev)
573
574    ! Reorder fluxes
575    if (allocated(flux%lw_up)) then
576      flux%lw_up(istartcol:iendcol,:) &
577           &  = flux_rev%lw_up(istartcol:iendcol,nlev+1:1:-1)
578      flux%lw_dn(istartcol:iendcol,:) &
579           &  = flux_rev%lw_dn(istartcol:iendcol,nlev+1:1:-1)
580      if (allocated(flux%lw_up_clear)) then
581        flux%lw_up_clear(istartcol:iendcol,:) &
582             &  = flux_rev%lw_up_clear(istartcol:iendcol,nlev+1:1:-1)
583        flux%lw_dn_clear(istartcol:iendcol,:) &
584             &  = flux_rev%lw_dn_clear(istartcol:iendcol,nlev+1:1:-1)
585      end if
586    end if
587    if (allocated(flux%sw_up)) then
588      flux%sw_up(istartcol:iendcol,:) &
589           &  = flux_rev%sw_up(istartcol:iendcol,nlev+1:1:-1)
590      flux%sw_dn(istartcol:iendcol,:) &
591           &  = flux_rev%sw_dn(istartcol:iendcol,nlev+1:1:-1)
592      if (allocated(flux%sw_dn_direct)) then
593        flux%sw_dn_direct(istartcol:iendcol,:) &
594             &  = flux_rev%sw_dn_direct(istartcol:iendcol,nlev+1:1:-1)
595      end if
596      if (allocated(flux%sw_up_clear)) then
597        flux%sw_up_clear(istartcol:iendcol,:) &
598             &  = flux_rev%sw_up_clear(istartcol:iendcol,nlev+1:1:-1)
599        flux%sw_dn_clear(istartcol:iendcol,:) &
600             &  = flux_rev%sw_dn_clear(istartcol:iendcol,nlev+1:1:-1)
601        if (allocated(flux%sw_dn_direct_clear)) then
602          flux%sw_dn_direct_clear(istartcol:iendcol,:) &
603               &  = flux_rev%sw_dn_direct_clear(istartcol:iendcol,nlev+1:1:-1)
604        end if
605      end if
606    end if
607
608    ! Deallocate reversed arrays
609    call thermodynamics_rev%deallocate
610    call gas_rev%deallocate
611    call cloud_rev%deallocate
612    call flux_rev%deallocate
613    if (allocated(aerosol%mixing_ratio)) then
614      call aerosol_rev%deallocate
615    end if
616
617  end subroutine radiation_reverse
618
619end module radiation_interface
Note: See TracBrowser for help on using the repository browser.