[4773] | 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 | |
---|
| 17 | module radiation_ecckd_interface |
---|
| 18 | |
---|
| 19 | implicit none |
---|
| 20 | |
---|
| 21 | public :: setup_gas_optics, set_gas_units, gas_optics !, planck_function |
---|
| 22 | |
---|
| 23 | contains |
---|
| 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 | use yomhook, only : lhook, dr_hook, jphook |
---|
| 32 | |
---|
| 33 | type(config_type), intent(inout), target :: config |
---|
| 34 | |
---|
| 35 | integer :: jj |
---|
| 36 | |
---|
| 37 | real(jphook) :: hook_handle |
---|
| 38 | |
---|
| 39 | if (lhook) call dr_hook('radiation_ecckd_interface:setup_gas_optics',0,hook_handle) |
---|
| 40 | |
---|
[4853] | 41 | if (config%do_sw .and. config%i_gas_model_sw == IGasModelECCKD) then |
---|
[4773] | 42 | |
---|
| 43 | ! Read shortwave ecCKD gas optics NetCDF file |
---|
| 44 | call config%gas_optics_sw%read(trim(config%gas_optics_sw_file_name), & |
---|
| 45 | & config%iverbosesetup) |
---|
| 46 | |
---|
| 47 | ! Copy over relevant properties |
---|
| 48 | config%n_g_sw = config%gas_optics_sw%ng |
---|
| 49 | |
---|
| 50 | if (config%do_cloud_aerosol_per_sw_g_point) then |
---|
| 51 | ! Bands and g points are the same |
---|
| 52 | config%n_bands_sw = config%n_g_sw |
---|
| 53 | else |
---|
| 54 | ! Bands are groups of g points and span a continuous region of |
---|
| 55 | ! wavenumber space |
---|
| 56 | config%n_bands_sw = config%gas_optics_sw%spectral_def%nband |
---|
| 57 | end if |
---|
| 58 | |
---|
[4853] | 59 | allocate(config%i_band_from_g_sw (config%n_g_sw)) |
---|
| 60 | allocate(config%i_band_from_reordered_g_sw(config%n_g_sw)) |
---|
| 61 | allocate(config%i_g_from_reordered_g_sw (config%n_g_sw)) |
---|
[4773] | 62 | |
---|
| 63 | if (config%do_cloud_aerosol_per_sw_g_point) then |
---|
| 64 | config%i_band_from_g_sw = [ (jj, jj = 1,config%n_g_sw) ] |
---|
| 65 | config%i_band_from_reordered_g_sw = [ (jj, jj = 1,config%n_g_sw) ] |
---|
| 66 | else |
---|
| 67 | config%i_band_from_g_sw & |
---|
| 68 | & = config%gas_optics_sw%spectral_def%i_band_number |
---|
| 69 | config%i_band_from_reordered_g_sw & |
---|
| 70 | & = config%gas_optics_sw%spectral_def%i_band_number |
---|
| 71 | end if |
---|
| 72 | config%i_g_from_reordered_g_sw = [ (jj, jj = 1,config%n_g_sw) ] |
---|
| 73 | |
---|
| 74 | if (config%iverbosesetup >= 2) then |
---|
| 75 | call config%gas_optics_sw%print() |
---|
| 76 | end if |
---|
| 77 | |
---|
| 78 | ! Override solar spectral irradiance, if filename provided |
---|
| 79 | if (config%use_spectral_solar_cycle) then |
---|
| 80 | call config%gas_optics_sw%read_spectral_solar_cycle(config%ssi_file_name, & |
---|
| 81 | & config%iverbosesetup, config%use_updated_solar_spectrum) |
---|
| 82 | end if |
---|
| 83 | |
---|
| 84 | end if |
---|
| 85 | |
---|
[4853] | 86 | if (config%do_lw .and. config%i_gas_model_lw == IGasModelECCKD) then |
---|
[4773] | 87 | |
---|
| 88 | ! Read longwave ecCKD gas optics NetCDF file |
---|
| 89 | call config%gas_optics_lw%read(trim(config%gas_optics_lw_file_name), & |
---|
| 90 | & config%iverbosesetup) |
---|
| 91 | |
---|
| 92 | ! Copy over relevant properties |
---|
| 93 | config%n_g_lw = config%gas_optics_lw%ng |
---|
| 94 | |
---|
| 95 | if (config%do_cloud_aerosol_per_lw_g_point) then |
---|
| 96 | ! Bands and g points are the same |
---|
| 97 | config%n_bands_lw = config%n_g_lw |
---|
| 98 | else |
---|
| 99 | ! Bands are groups of g points and span a continuous region of |
---|
| 100 | ! wavenumber space |
---|
| 101 | config%n_bands_lw = config%gas_optics_lw%spectral_def%nband |
---|
| 102 | end if |
---|
| 103 | |
---|
[4853] | 104 | allocate(config%i_band_from_g_lw (config%n_g_lw)) |
---|
| 105 | allocate(config%i_band_from_reordered_g_lw(config%n_g_lw)) |
---|
| 106 | allocate(config%i_g_from_reordered_g_lw (config%n_g_lw)) |
---|
[4773] | 107 | |
---|
| 108 | if (config%do_cloud_aerosol_per_lw_g_point) then |
---|
| 109 | config%i_band_from_g_lw = [ (jj, jj = 1,config%n_g_lw) ] |
---|
| 110 | config%i_band_from_reordered_g_lw = [ (jj, jj = 1,config%n_g_lw) ] |
---|
| 111 | else |
---|
| 112 | config%i_band_from_g_lw & |
---|
| 113 | & = config%gas_optics_lw%spectral_def%i_band_number |
---|
| 114 | config%i_band_from_reordered_g_lw & |
---|
| 115 | & = config%gas_optics_lw%spectral_def%i_band_number |
---|
| 116 | end if |
---|
| 117 | config%i_g_from_reordered_g_lw = [ (jj, jj = 1,config%n_g_lw) ] |
---|
| 118 | |
---|
| 119 | if (config%iverbosesetup >= 2) then |
---|
| 120 | call config%gas_optics_lw%print() |
---|
| 121 | end if |
---|
| 122 | |
---|
| 123 | end if |
---|
| 124 | |
---|
| 125 | ! The i_spec_* variables are used solely for storing spectral |
---|
| 126 | ! data, and this can either be by band or by g-point |
---|
| 127 | if (config%do_save_spectral_flux) then |
---|
| 128 | if (config%do_save_gpoint_flux) then |
---|
| 129 | config%n_spec_sw = config%n_g_sw |
---|
| 130 | config%n_spec_lw = config%n_g_lw |
---|
| 131 | config%i_spec_from_reordered_g_sw => config%i_g_from_reordered_g_sw |
---|
| 132 | config%i_spec_from_reordered_g_lw => config%i_g_from_reordered_g_lw |
---|
| 133 | else |
---|
| 134 | config%n_spec_sw = config%n_bands_sw |
---|
| 135 | config%n_spec_lw = config%n_bands_lw |
---|
| 136 | config%i_spec_from_reordered_g_sw => config%i_band_from_reordered_g_sw |
---|
| 137 | config%i_spec_from_reordered_g_lw => config%i_band_from_reordered_g_lw |
---|
| 138 | end if |
---|
| 139 | else |
---|
| 140 | config%n_spec_sw = 0 |
---|
| 141 | config%n_spec_lw = 0 |
---|
| 142 | nullify(config%i_spec_from_reordered_g_sw) |
---|
| 143 | nullify(config%i_spec_from_reordered_g_lw) |
---|
| 144 | end if |
---|
| 145 | |
---|
| 146 | if (lhook) call dr_hook('radiation_ecckd_interface:setup_gas_optics',1,hook_handle) |
---|
| 147 | |
---|
| 148 | end subroutine setup_gas_optics |
---|
| 149 | |
---|
| 150 | |
---|
| 151 | !--------------------------------------------------------------------- |
---|
| 152 | ! Scale gas mixing ratios according to required units |
---|
| 153 | subroutine set_gas_units(gas) |
---|
| 154 | |
---|
| 155 | use radiation_gas, only : gas_type, IVolumeMixingRatio |
---|
| 156 | use yomhook, only : lhook, dr_hook, jphook |
---|
| 157 | |
---|
| 158 | type(gas_type), intent(inout) :: gas |
---|
| 159 | |
---|
| 160 | real(jphook) :: hook_handle |
---|
| 161 | |
---|
| 162 | if (lhook) call dr_hook('radiation_ecckd_interface:set_gas_units',0,hook_handle) |
---|
| 163 | |
---|
| 164 | call gas%set_units(IVolumeMixingRatio) |
---|
| 165 | |
---|
| 166 | if (lhook) call dr_hook('radiation_ecckd_interface:set_gas_units',1,hook_handle) |
---|
| 167 | |
---|
| 168 | end subroutine set_gas_units |
---|
| 169 | |
---|
| 170 | |
---|
| 171 | !--------------------------------------------------------------------- |
---|
| 172 | ! Compute gas optical depths, shortwave scattering, Planck function |
---|
| 173 | ! and incoming shortwave radiation at top-of-atmosphere |
---|
| 174 | subroutine gas_optics(ncol,nlev,istartcol,iendcol, & |
---|
| 175 | & config, single_level, thermodynamics, gas, & |
---|
| 176 | & od_lw, od_sw, ssa_sw, lw_albedo, planck_hl, lw_emission, & |
---|
| 177 | & incoming_sw) |
---|
| 178 | |
---|
| 179 | use parkind1, only : jprb |
---|
| 180 | use yomhook, only : lhook, dr_hook, jphook |
---|
| 181 | |
---|
[4853] | 182 | use radiation_config, only : config_type, IGasModelECCKD |
---|
[4773] | 183 | use radiation_thermodynamics, only : thermodynamics_type |
---|
| 184 | use radiation_single_level, only : single_level_type |
---|
| 185 | use radiation_gas_constants, only : NMaxGases |
---|
| 186 | use radiation_gas |
---|
| 187 | |
---|
| 188 | integer, intent(in) :: ncol ! number of columns |
---|
| 189 | integer, intent(in) :: nlev ! number of levels |
---|
| 190 | integer, intent(in) :: istartcol, iendcol ! range of cols to process |
---|
| 191 | type(config_type), intent(in) :: config |
---|
| 192 | type(single_level_type), intent(in) :: single_level |
---|
| 193 | type(thermodynamics_type),intent(in) :: thermodynamics |
---|
| 194 | type(gas_type), intent(in) :: gas |
---|
| 195 | |
---|
| 196 | ! Longwave albedo of the surface |
---|
| 197 | real(jprb), dimension(config%n_g_lw,istartcol:iendcol), & |
---|
| 198 | & intent(in), optional :: lw_albedo |
---|
| 199 | |
---|
| 200 | ! Gaseous layer optical depth in longwave and shortwave, and |
---|
| 201 | ! shortwave single scattering albedo (i.e. fraction of extinction |
---|
| 202 | ! due to Rayleigh scattering) at each g-point |
---|
| 203 | real(jprb), dimension(config%n_g_lw,nlev,istartcol:iendcol), intent(out) :: & |
---|
| 204 | & od_lw |
---|
| 205 | real(jprb), dimension(config%n_g_sw,nlev,istartcol:iendcol), intent(out) :: & |
---|
| 206 | & od_sw, ssa_sw |
---|
| 207 | |
---|
| 208 | ! The Planck function (emitted flux from a black body) at half |
---|
| 209 | ! levels at each longwave g-point |
---|
| 210 | real(jprb), dimension(config%n_g_lw,nlev+1,istartcol:iendcol), & |
---|
| 211 | & intent(out), optional :: planck_hl |
---|
| 212 | ! Planck function for the surface (W m-2) |
---|
| 213 | real(jprb), dimension(config%n_g_lw,istartcol:iendcol), & |
---|
| 214 | & intent(out), optional :: lw_emission |
---|
| 215 | |
---|
| 216 | ! The incoming shortwave flux into a plane perpendicular to the |
---|
| 217 | ! incoming radiation at top-of-atmosphere in each of the shortwave |
---|
| 218 | ! g-points |
---|
| 219 | real(jprb), dimension(config%n_g_sw,istartcol:iendcol), & |
---|
| 220 | & intent(out), optional :: incoming_sw |
---|
| 221 | |
---|
| 222 | ! Temperature at full levels (K) |
---|
| 223 | real(jprb) :: temperature_fl(istartcol:iendcol,nlev) |
---|
| 224 | |
---|
[4853] | 225 | real(jprb) :: concentration_scaling(NMaxGases) |
---|
| 226 | |
---|
| 227 | logical :: is_volume_mixing_ratio |
---|
| 228 | |
---|
| 229 | integer :: jcol, jlev, jg |
---|
[4773] | 230 | |
---|
| 231 | real(jphook) :: hook_handle |
---|
| 232 | |
---|
| 233 | if (lhook) call dr_hook('radiation_ecckd_interface:gas_optics',0,hook_handle) |
---|
| 234 | |
---|
| 235 | !temperature_fl(istartcol:iendcol,:) & |
---|
| 236 | ! & = 0.5_jprb * (thermodynamics%temperature_hl(istartcol:iendcol,1:nlev) & |
---|
| 237 | ! & +thermodynamics%temperature_hl(istartcol:iendcol,2:nlev+1)) |
---|
| 238 | |
---|
| 239 | temperature_fl(istartcol:iendcol,:) & |
---|
| 240 | & = (thermodynamics%temperature_hl(istartcol:iendcol,1:nlev) & |
---|
| 241 | & *thermodynamics%pressure_hl(istartcol:iendcol,1:nlev) & |
---|
| 242 | & +thermodynamics%temperature_hl(istartcol:iendcol,2:nlev+1) & |
---|
| 243 | & *thermodynamics%pressure_hl(istartcol:iendcol,2:nlev+1)) & |
---|
| 244 | & / (thermodynamics%pressure_hl(istartcol:iendcol,1:nlev) & |
---|
| 245 | & +thermodynamics%pressure_hl(istartcol:iendcol,2:nlev+1)) |
---|
| 246 | |
---|
[4853] | 247 | ! Check that the gas concentrations are stored in volume mixing |
---|
| 248 | ! ratio with no scaling; if not, return a vector of scalings |
---|
| 249 | call gas%assert_units(IVolumeMixingRatio, scale_factor=1.0_jprb, & |
---|
| 250 | & istatus=is_volume_mixing_ratio) |
---|
| 251 | if (.not. is_volume_mixing_ratio) then |
---|
| 252 | call gas%get_scaling(IVolumeMixingRatio, concentration_scaling) |
---|
| 253 | else |
---|
| 254 | concentration_scaling = 1.0_jprb |
---|
| 255 | end if |
---|
| 256 | |
---|
| 257 | if (config%do_sw .and. config%i_gas_model_sw == IGasModelECCKD) then |
---|
| 258 | |
---|
| 259 | if (is_volume_mixing_ratio) then |
---|
| 260 | call config%gas_optics_sw%calc_optical_depth(ncol,nlev,istartcol,iendcol, & |
---|
| 261 | & NMaxGases, thermodynamics%pressure_hl, & |
---|
| 262 | & temperature_fl, gas%mixing_ratio, & |
---|
| 263 | & od_sw, rayleigh_od_fl=ssa_sw) |
---|
| 264 | else |
---|
| 265 | call config%gas_optics_sw%calc_optical_depth(ncol,nlev,istartcol,iendcol, & |
---|
| 266 | & NMaxGases, thermodynamics%pressure_hl, & |
---|
| 267 | & temperature_fl, gas%mixing_ratio, & |
---|
| 268 | & od_sw, rayleigh_od_fl=ssa_sw, concentration_scaling=concentration_scaling) |
---|
| 269 | end if |
---|
| 270 | |
---|
[4773] | 271 | ! At this point od_sw = absorption optical depth and ssa_sw = |
---|
| 272 | ! rayleigh optical depth: convert to total optical depth and |
---|
| 273 | ! single-scattering albedo |
---|
[4853] | 274 | do jcol = istartcol,iendcol |
---|
| 275 | do jlev = 1, nlev |
---|
| 276 | do jg = 1, config%n_g_sw |
---|
| 277 | od_sw(jg,jlev,jcol) = od_sw(jg,jlev,jcol) + ssa_sw(jg,jlev,jcol) |
---|
| 278 | ssa_sw(jg,jlev,jcol) = ssa_sw(jg,jlev,jcol) / od_sw(jg,jlev,jcol) |
---|
| 279 | end do |
---|
| 280 | end do |
---|
| 281 | end do |
---|
[4773] | 282 | |
---|
| 283 | if (present(incoming_sw)) then |
---|
| 284 | if (single_level%spectral_solar_cycle_multiplier == 0.0_jprb) then |
---|
| 285 | call config%gas_optics_sw%calc_incoming_sw(single_level%solar_irradiance, & |
---|
| 286 | & incoming_sw) |
---|
| 287 | else |
---|
| 288 | call config%gas_optics_sw%calc_incoming_sw(single_level%solar_irradiance, & |
---|
| 289 | & incoming_sw, single_level%spectral_solar_cycle_multiplier) |
---|
| 290 | end if |
---|
| 291 | end if |
---|
| 292 | |
---|
| 293 | end if |
---|
| 294 | |
---|
[4853] | 295 | if (config%do_lw .and. config%i_gas_model_lw == IGasModelECCKD) then |
---|
[4773] | 296 | |
---|
[4853] | 297 | if (is_volume_mixing_ratio) then |
---|
| 298 | call config%gas_optics_lw%calc_optical_depth(ncol,nlev,istartcol,iendcol, & |
---|
| 299 | & NMaxGases, thermodynamics%pressure_hl, & |
---|
| 300 | & temperature_fl, gas%mixing_ratio, & |
---|
| 301 | & od_lw) |
---|
| 302 | else |
---|
| 303 | call config%gas_optics_lw%calc_optical_depth(ncol,nlev,istartcol,iendcol, & |
---|
| 304 | & NMaxGases, thermodynamics%pressure_hl, & |
---|
| 305 | & temperature_fl, gas%mixing_ratio, & |
---|
| 306 | & od_lw, concentration_scaling=concentration_scaling) |
---|
| 307 | end if |
---|
[4773] | 308 | |
---|
| 309 | ! Calculate the Planck function for each g point |
---|
| 310 | do jcol = istartcol,iendcol |
---|
| 311 | call config%gas_optics_lw%calc_planck_function(nlev+1, & |
---|
| 312 | & thermodynamics%temperature_hl(jcol,:), planck_hl(:,:,jcol)) |
---|
| 313 | end do |
---|
| 314 | call config%gas_optics_lw%calc_planck_function(iendcol+1-istartcol, & |
---|
| 315 | & single_level%skin_temperature(istartcol:iendcol), & |
---|
| 316 | & lw_emission(:,:)) |
---|
[4853] | 317 | !NEC$ forced_collapse |
---|
[4773] | 318 | lw_emission = lw_emission * (1.0_jprb - lw_albedo) |
---|
| 319 | |
---|
| 320 | end if |
---|
| 321 | |
---|
| 322 | if (lhook) call dr_hook('radiation_ecckd_interface:gas_optics',1,hook_handle) |
---|
| 323 | |
---|
| 324 | end subroutine gas_optics |
---|
| 325 | |
---|
| 326 | ! !--------------------------------------------------------------------- |
---|
| 327 | ! ! Externally facing function for computing the Planck function |
---|
| 328 | ! ! without reference to any gas profile; typically this would be used |
---|
| 329 | ! ! for computing the emission by a surface. |
---|
| 330 | ! subroutine planck_function(config, temperature, planck_surf) |
---|
| 331 | |
---|
| 332 | ! use parkind1, only : jprb |
---|
| 333 | ! use radiation_config, only : config_type |
---|
| 334 | |
---|
| 335 | ! type(config_type), intent(in) :: config |
---|
| 336 | ! real(jprb), intent(in) :: temperature |
---|
| 337 | |
---|
| 338 | ! ! Planck function of the surface (W m-2) |
---|
| 339 | ! real(jprb), dimension(config%n_g_lw), intent(out) :: planck_surf |
---|
| 340 | |
---|
| 341 | ! end subroutine planck_function |
---|
| 342 | |
---|
| 343 | end module radiation_ecckd_interface |
---|