source: LMDZ6/trunk/libf/phylmd/ecrad.v1.5.1/radiation_ecckd_interface.F90 @ 5319

Last change on this file since 5319 was 4677, checked in by idelkadi, 14 months ago

Implementation in the LMDZ code of the double call of the ECRAD radiative transfer code to estimate the 3D radiative effect of clouds.

  • This double call of Ecrad is controlled by the ok_3Deffect logic key.
  • If this key is enabled, 2 files of parameter configuration "namelists" for ECRAD are required at runtime: namelist_ecrad and namelist_ecrad_s2.
  • If this key is deactivated, the configuration and initialization part (reading namelist and netcdf files) is performed only once during simulation (1st call to ECRAD). Otherwise, configuration and initialization are performed each time Ecrad is called.
File size: 11.1 KB
Line 
1! radiation_ecckd_interface.F90 - Interface to ecCKD gas optics model
2!
3! (C) Copyright 2020- 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! License: see the COPYING file for details
15!
16
17module radiation_ecckd_interface
18
19  implicit none
20
21  public  :: setup_gas_optics, set_gas_units, gas_optics !, planck_function
22
23contains
24
25  !---------------------------------------------------------------------
26  ! Setup the ecCKD generalized gas optics model
27  subroutine setup_gas_optics(config)
28
29    use parkind1, only : jprb
30    use radiation_config
31
32    type(config_type), intent(inout), target :: config
33
34    integer :: jj
35   
36    if (config%do_sw) then
37
38      ! Read shortwave ecCKD gas optics NetCDF file
39      call config%gas_optics_sw%read(trim(config%gas_optics_sw_file_name), &
40           &                         config%iverbosesetup)
41
42      ! Copy over relevant properties
43      config%n_g_sw     = config%gas_optics_sw%ng
44
45      if (config%do_cloud_aerosol_per_sw_g_point) then
46        ! Bands and g points are the same
47        config%n_bands_sw = config%n_g_sw
48      else
49        ! Bands are groups of g points and span a continuous region of
50        ! wavenumber space
51        config%n_bands_sw = config%gas_optics_sw%spectral_def%nband
52      end if
53
54      if (allocated(config%i_band_from_g_sw)) deallocate(config%i_band_from_g_sw)
55      allocate(config%i_band_from_g_sw(config%n_g_sw))
56      if (allocated(config%i_band_from_reordered_g_sw)) deallocate(config%i_band_from_reordered_g_sw)
57      allocate(config%i_band_from_reordered_g_sw(config%n_g_sw))
58      if (allocated(config%i_g_from_reordered_g_sw)) deallocate(config%i_g_from_reordered_g_sw)
59      allocate(config%i_g_from_reordered_g_sw(config%n_g_sw))
60       
61      if (config%do_cloud_aerosol_per_sw_g_point) then
62        config%i_band_from_g_sw           = [ (jj, jj = 1,config%n_g_sw) ]
63        config%i_band_from_reordered_g_sw = [ (jj, jj = 1,config%n_g_sw) ]
64      else
65        config%i_band_from_g_sw &
66             &  = config%gas_optics_sw%spectral_def%i_band_number
67        config%i_band_from_reordered_g_sw &
68             &  = config%gas_optics_sw%spectral_def%i_band_number
69      end if
70      config%i_g_from_reordered_g_sw      = [ (jj, jj = 1,config%n_g_sw) ]
71
72      if (config%iverbosesetup >= 2) then
73        call config%gas_optics_sw%print()
74      end if
75
76    end if
77
78    if (config%do_lw) then
79
80      ! Read longwave ecCKD gas optics NetCDF file
81      call config%gas_optics_lw%read(trim(config%gas_optics_lw_file_name), &
82           &                         config%iverbosesetup)
83
84      ! Copy over relevant properties
85      config%n_g_lw     = config%gas_optics_lw%ng
86
87      if (config%do_cloud_aerosol_per_lw_g_point) then
88        ! Bands and g points are the same
89        config%n_bands_lw = config%n_g_lw
90      else
91        ! Bands are groups of g points and span a continuous region of
92        ! wavenumber space
93        config%n_bands_lw = config%gas_optics_lw%spectral_def%nband
94      end if
95
96      if (allocated(config%i_band_from_g_lw)) deallocate(config%i_band_from_g_lw)
97      allocate(config%i_band_from_g_lw          (config%n_g_lw))
98      if (allocated(config%i_band_from_reordered_g_lw)) deallocate(config%i_band_from_reordered_g_lw)
99      allocate(config%i_band_from_reordered_g_lw(config%n_g_lw))
100      if (allocated(config%i_g_from_reordered_g_lw)) deallocate(config%i_g_from_reordered_g_lw)
101      allocate(config%i_g_from_reordered_g_lw   (config%n_g_lw))
102
103      if (config%do_cloud_aerosol_per_lw_g_point) then
104        config%i_band_from_g_lw           = [ (jj, jj = 1,config%n_g_lw) ]
105        config%i_band_from_reordered_g_lw = [ (jj, jj = 1,config%n_g_lw) ]
106      else
107        config%i_band_from_g_lw &
108             &  = config%gas_optics_lw%spectral_def%i_band_number
109        config%i_band_from_reordered_g_lw &
110             &  = config%gas_optics_lw%spectral_def%i_band_number
111      end if
112      config%i_g_from_reordered_g_lw      = [ (jj, jj = 1,config%n_g_lw) ]
113
114      if (config%iverbosesetup >= 2) then
115        call config%gas_optics_lw%print()
116      end if
117
118    end if
119
120    ! The i_spec_* variables are used solely for storing spectral
121    ! data, and this can either be by band or by g-point
122    if (config%do_save_spectral_flux) then
123      if (config%do_save_gpoint_flux) then
124        config%n_spec_sw = config%n_g_sw
125        config%n_spec_lw = config%n_g_lw
126        config%i_spec_from_reordered_g_sw => config%i_g_from_reordered_g_sw
127        config%i_spec_from_reordered_g_lw => config%i_g_from_reordered_g_lw
128      else
129        config%n_spec_sw = config%n_bands_sw
130        config%n_spec_lw = config%n_bands_lw
131        config%i_spec_from_reordered_g_sw => config%i_band_from_reordered_g_sw
132        config%i_spec_from_reordered_g_lw => config%i_band_from_reordered_g_lw
133      end if
134    else
135      config%n_spec_sw = 0
136      config%n_spec_lw = 0
137      nullify(config%i_spec_from_reordered_g_sw)
138      nullify(config%i_spec_from_reordered_g_lw)
139    end if
140
141  end subroutine setup_gas_optics
142
143
144  !---------------------------------------------------------------------
145  ! Scale gas mixing ratios according to required units
146  subroutine set_gas_units(gas)
147
148    use radiation_gas,           only : gas_type, IVolumeMixingRatio
149    type(gas_type),    intent(inout) :: gas
150
151    call gas%set_units(IVolumeMixingRatio)
152
153  end subroutine set_gas_units
154
155
156  !---------------------------------------------------------------------
157  ! Compute gas optical depths, shortwave scattering, Planck function
158  ! and incoming shortwave radiation at top-of-atmosphere
159  subroutine gas_optics(ncol,nlev,istartcol,iendcol, &
160       &  config, single_level, thermodynamics, gas, &
161       &  od_lw, od_sw, ssa_sw, lw_albedo, planck_hl, lw_emission, &
162       &  incoming_sw)
163
164    use parkind1, only : jprb
165    use radiation_config,         only : config_type
166    use radiation_thermodynamics, only : thermodynamics_type
167    use radiation_single_level,   only : single_level_type
168    use radiation_gas_constants,  only : NMaxGases
169    use radiation_gas
170
171    integer,                  intent(in) :: ncol               ! number of columns
172    integer,                  intent(in) :: nlev               ! number of levels
173    integer,                  intent(in) :: istartcol, iendcol ! range of cols to process
174    type(config_type),        intent(in) :: config
175    type(single_level_type),  intent(in) :: single_level
176    type(thermodynamics_type),intent(in) :: thermodynamics
177    type(gas_type),           intent(in) :: gas
178
179    ! Longwave albedo of the surface
180    real(jprb), dimension(config%n_g_lw,istartcol:iendcol), &
181         &  intent(in), optional :: lw_albedo
182
183    ! Gaseous layer optical depth in longwave and shortwave, and
184    ! shortwave single scattering albedo (i.e. fraction of extinction
185    ! due to Rayleigh scattering) at each g-point
186    real(jprb), dimension(config%n_g_lw,nlev,istartcol:iendcol), intent(out) :: &
187         &   od_lw
188    real(jprb), dimension(config%n_g_sw,nlev,istartcol:iendcol), intent(out) :: &
189         &   od_sw, ssa_sw
190
191    ! The Planck function (emitted flux from a black body) at half
192    ! levels at each longwave g-point
193    real(jprb), dimension(config%n_g_lw,nlev+1,istartcol:iendcol), &
194         &   intent(out), optional :: planck_hl
195    ! Planck function for the surface (W m-2)
196    real(jprb), dimension(config%n_g_lw,istartcol:iendcol), &
197         &   intent(out), optional :: lw_emission
198
199    ! The incoming shortwave flux into a plane perpendicular to the
200    ! incoming radiation at top-of-atmosphere in each of the shortwave
201    ! g-points
202    real(jprb), dimension(config%n_g_sw,istartcol:iendcol), &
203         &   intent(out), optional :: incoming_sw
204
205    ! Temperature at full levels (K)
206    real(jprb) :: temperature_fl(istartcol:iendcol,nlev)
207
208    integer :: jcol
209
210    !temperature_fl(istartcol:iendcol,:) &
211    !     &  = 0.5_jprb * (thermodynamics%temperature_hl(istartcol:iendcol,1:nlev) &
212    !     &               +thermodynamics%temperature_hl(istartcol:iendcol,2:nlev+1))
213 
214    temperature_fl(istartcol:iendcol,:) &
215         &  = (thermodynamics%temperature_hl(istartcol:iendcol,1:nlev) &
216         &     *thermodynamics%pressure_hl(istartcol:iendcol,1:nlev) &
217         &    +thermodynamics%temperature_hl(istartcol:iendcol,2:nlev+1) &
218         &     *thermodynamics%pressure_hl(istartcol:iendcol,2:nlev+1)) &
219         &  / (thermodynamics%pressure_hl(istartcol:iendcol,1:nlev) &
220         &    +thermodynamics%pressure_hl(istartcol:iendcol,2:nlev+1))
221 
222    if (config%do_sw) then
223
224      call config%gas_optics_sw%calc_optical_depth(ncol,nlev,istartcol,iendcol, &
225           &  NMaxGases, thermodynamics%pressure_hl, &
226           &  temperature_fl, &
227           &  gas%mixing_ratio, &
228!           &  reshape(gas%mixing_ratio(istartcol:iendcol,:,:), &
229!           &          [nlev,iendcol-istartcol+1,NMaxGases],order=[2,1,3]), &
230           &  od_sw, rayleigh_od_fl=ssa_sw)
231      ! At this point od_sw = absorption optical depth and ssa_sw =
232      ! rayleigh optical depth: convert to total optical depth and
233      ! single-scattering albedo
234      od_sw = od_sw + ssa_sw
235      ssa_sw = ssa_sw / od_sw
236
237      if (present(incoming_sw)) then
238        call config%gas_optics_sw%calc_incoming_sw(single_level%solar_irradiance, incoming_sw)
239      end if
240
241    end if
242
243    if (config%do_lw) then
244
245      call config%gas_optics_lw%calc_optical_depth(ncol,nlev,istartcol,iendcol, &
246           &  NMaxGases, thermodynamics%pressure_hl, &
247           &  temperature_fl, &
248           &  gas%mixing_ratio, &
249!           &  reshape(gas%mixing_ratio(istartcol:iendcol,:,:), &
250!           &          [nlev,iendcol-istartcol+1,NMaxGases],order=[2,1,3]), &
251           &  od_lw)
252
253      ! Calculate the Planck function for each g point
254      do jcol = istartcol,iendcol
255        call config%gas_optics_lw%calc_planck_function(nlev+1, &
256             &  thermodynamics%temperature_hl(jcol,:), planck_hl(:,:,jcol))
257      end do
258      call config%gas_optics_lw%calc_planck_function(iendcol+1-istartcol, &
259           &  single_level%skin_temperature(istartcol:iendcol), &
260           &  lw_emission(:,:))
261      lw_emission = lw_emission * (1.0_jprb - lw_albedo)
262
263    end if
264
265  end subroutine gas_optics
266
267  ! !---------------------------------------------------------------------
268  ! ! Externally facing function for computing the Planck function
269  ! ! without reference to any gas profile; typically this would be used
270  ! ! for computing the emission by a surface.
271  ! subroutine planck_function(config, temperature, planck_surf)
272
273  !   use parkind1,                 only : jprb
274  !   use radiation_config,         only : config_type
275
276  !   type(config_type), intent(in) :: config
277  !   real(jprb),        intent(in) :: temperature
278
279  !   ! Planck function of the surface (W m-2)
280  !   real(jprb), dimension(config%n_g_lw), intent(out) :: planck_surf
281
282  ! end subroutine planck_function
283
284end module radiation_ecckd_interface
Note: See TracBrowser for help on using the repository browser.