[3880] | 1 | ! radiation_single_level.F90 - Derived type for single-level fields |
---|
| 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 | ! 2019-01-14 R. Hogan Added out_of_physical_bounds routine |
---|
| 17 | ! 2019-01-18 R. Hogan Added weighted albedo mapping |
---|
| 18 | |
---|
| 19 | module radiation_single_level |
---|
| 20 | |
---|
| 21 | use parkind1, only : jprb |
---|
| 22 | |
---|
| 23 | implicit none |
---|
| 24 | public |
---|
| 25 | |
---|
| 26 | !--------------------------------------------------------------------- |
---|
| 27 | ! Derived type to contain variables that don't vary with height; |
---|
| 28 | ! mostly they are surface properties |
---|
| 29 | type single_level_type |
---|
| 30 | ! Note that skin temperature is only defined if |
---|
| 31 | ! is_simple_surface==.true. |
---|
| 32 | real(jprb), allocatable, dimension(:) :: & |
---|
| 33 | & cos_sza, & ! (ncol) Cosine of solar zenith angle |
---|
| 34 | & skin_temperature ! (ncol) Skin temperature (K) |
---|
| 35 | |
---|
| 36 | ! Shortwave albedo: if sw_albedo_direct is not allocated then |
---|
| 37 | ! sw_albedo will be used for both direct and diffuse solar |
---|
| 38 | ! radiation; otherwise, sw_albedo will be used for diffuse |
---|
| 39 | ! radiation and sw_albedo_direct for direct radiation. |
---|
| 40 | real(jprb), allocatable, dimension(:,:) :: & |
---|
| 41 | & sw_albedo, & ! (ncol,nalbedobands) |
---|
| 42 | & sw_albedo_direct ! (ncol,nalbedobands) |
---|
| 43 | |
---|
| 44 | ! If is_simple_surface==.true., we provide longwave emissivity |
---|
| 45 | ! coarser spectral intervals, while if |
---|
| 46 | ! is_simple_surface==.false. it will be generated by the surface |
---|
| 47 | ! radiation scheme and may optionally be on the full spectral grid |
---|
| 48 | ! of the atmospheric model. |
---|
| 49 | real(jprb), allocatable, dimension(:,:) :: & |
---|
| 50 | & lw_emissivity ! (ncol,nemissbands) If |
---|
| 51 | ! is_simple_surface==.false. then we specify the emission instead |
---|
| 52 | ! of the skin temperature, and again this may be on the emissivity |
---|
| 53 | ! bands or the full spectral grid |
---|
| 54 | real(jprb), allocatable, dimension(:,:) :: & |
---|
| 55 | & lw_emission ! (ncol,nemissbands) |
---|
| 56 | |
---|
| 57 | ! Incoming solar irradiance at the Earth is specified as the total |
---|
| 58 | ! across the shortwave spectrum. Note that this needs to be |
---|
| 59 | ! adjusted for Earth-Sun distance, solar cycle etc. To get |
---|
| 60 | ! normalized fluxes out, simply set this to 1.0. |
---|
| 61 | real(jprb) :: solar_irradiance = 1366.0_jprb ! W m-2 |
---|
| 62 | |
---|
| 63 | ! If config%use_spectral_solar_irradiance==true then this will be |
---|
| 64 | ! scaled by spectral_solar_scaling |
---|
| 65 | real(jprb), allocatable, dimension(:) :: & |
---|
| 66 | & spectral_solar_scaling ! (n_bands_sw) |
---|
| 67 | |
---|
| 68 | ! Seed for random number generator in McICA; it is expected that |
---|
| 69 | ! the host model generates this from the model time, longitude |
---|
| 70 | ! and latitude, in order that the model is deterministic |
---|
| 71 | integer, allocatable, dimension(:) :: iseed ! (ncol) |
---|
| 72 | |
---|
| 73 | ! Is the underlying surface a simple single "flat" tile? If so we |
---|
| 74 | ! describe it with skin_temperature and lw_emissivity. Otherwise |
---|
| 75 | ! we describe it with lw_emission and lw_emissivity coming out of |
---|
| 76 | ! the surface radiation scheme. |
---|
| 77 | logical :: is_simple_surface = .true. |
---|
| 78 | |
---|
| 79 | ! If is_simple_surface==.false., do we use the full number of g |
---|
| 80 | ! points for dimensioning sw_albedo, sw_albedo_direct and |
---|
| 81 | ! lw_emission? |
---|
| 82 | !logical :: use_full_spectrum_sw = .false. |
---|
| 83 | !logical :: use_full_spectrum_lw = .true. |
---|
| 84 | |
---|
| 85 | contains |
---|
| 86 | procedure :: allocate => allocate_single_level |
---|
| 87 | procedure :: deallocate => deallocate_single_level |
---|
| 88 | procedure :: init_seed_simple |
---|
| 89 | procedure :: get_albedos |
---|
| 90 | procedure :: out_of_physical_bounds |
---|
| 91 | |
---|
| 92 | end type single_level_type |
---|
| 93 | |
---|
| 94 | contains |
---|
| 95 | |
---|
| 96 | !--------------------------------------------------------------------- |
---|
| 97 | subroutine allocate_single_level(this, ncol, nalbedobands, nemisbands, & |
---|
| 98 | & use_sw_albedo_direct, is_simple_surface) |
---|
| 99 | |
---|
| 100 | use yomhook, only : lhook, dr_hook |
---|
| 101 | |
---|
| 102 | class(single_level_type), intent(inout) :: this |
---|
| 103 | integer, intent(in) :: ncol, nalbedobands, nemisbands |
---|
| 104 | logical, optional, intent(in) :: use_sw_albedo_direct |
---|
| 105 | logical, optional, intent(in) :: is_simple_surface |
---|
| 106 | |
---|
| 107 | real(jprb) :: hook_handle |
---|
| 108 | |
---|
| 109 | if (lhook) call dr_hook('radiation_single_level:allocate',0,hook_handle) |
---|
| 110 | |
---|
| 111 | call this%deallocate() |
---|
| 112 | |
---|
| 113 | if (present(is_simple_surface)) then |
---|
| 114 | this%is_simple_surface = is_simple_surface |
---|
| 115 | else |
---|
| 116 | this%is_simple_surface = .true. |
---|
| 117 | end if |
---|
| 118 | |
---|
| 119 | allocate(this%cos_sza(ncol)) |
---|
| 120 | |
---|
| 121 | if (this%is_simple_surface) then |
---|
| 122 | allocate(this%skin_temperature(ncol)) |
---|
| 123 | else |
---|
| 124 | allocate(this%lw_emission(ncol, nemisbands)) |
---|
| 125 | end if |
---|
| 126 | allocate(this%lw_emissivity(ncol, nemisbands)) |
---|
| 127 | |
---|
| 128 | allocate(this%sw_albedo(ncol, nalbedobands)) |
---|
| 129 | |
---|
| 130 | if (present(use_sw_albedo_direct)) then |
---|
| 131 | if (use_sw_albedo_direct) then |
---|
| 132 | allocate(this%sw_albedo_direct(ncol, nalbedobands)) |
---|
| 133 | end if |
---|
| 134 | end if |
---|
| 135 | |
---|
| 136 | allocate(this%iseed(ncol)) |
---|
| 137 | |
---|
| 138 | if (lhook) call dr_hook('radiation_single_level:allocate',1,hook_handle) |
---|
| 139 | |
---|
| 140 | end subroutine allocate_single_level |
---|
| 141 | |
---|
| 142 | |
---|
| 143 | !--------------------------------------------------------------------- |
---|
| 144 | subroutine deallocate_single_level(this) |
---|
| 145 | |
---|
| 146 | use yomhook, only : lhook, dr_hook |
---|
| 147 | |
---|
| 148 | class(single_level_type), intent(inout) :: this |
---|
| 149 | |
---|
| 150 | real(jprb) :: hook_handle |
---|
| 151 | |
---|
| 152 | if (lhook) call dr_hook('radiation_single_level:deallocate',0,hook_handle) |
---|
| 153 | |
---|
| 154 | if (allocated(this%cos_sza)) then |
---|
| 155 | deallocate(this%cos_sza) |
---|
| 156 | end if |
---|
| 157 | if (allocated(this%skin_temperature)) then |
---|
| 158 | deallocate(this%skin_temperature) |
---|
| 159 | end if |
---|
| 160 | if (allocated(this%sw_albedo)) then |
---|
| 161 | deallocate(this%sw_albedo) |
---|
| 162 | end if |
---|
| 163 | if (allocated(this%sw_albedo_direct)) then |
---|
| 164 | deallocate(this%sw_albedo_direct) |
---|
| 165 | end if |
---|
| 166 | if (allocated(this%lw_emissivity)) then |
---|
| 167 | deallocate(this%lw_emissivity) |
---|
| 168 | end if |
---|
| 169 | if (allocated(this%lw_emission)) then |
---|
| 170 | deallocate(this%lw_emission) |
---|
| 171 | end if |
---|
| 172 | if (allocated(this%spectral_solar_scaling)) then |
---|
| 173 | deallocate(this%spectral_solar_scaling) |
---|
| 174 | end if |
---|
| 175 | if (allocated(this%iseed)) then |
---|
| 176 | deallocate(this%iseed) |
---|
| 177 | end if |
---|
| 178 | |
---|
| 179 | if (lhook) call dr_hook('radiation_single_level:deallocate',1,hook_handle) |
---|
| 180 | |
---|
| 181 | end subroutine deallocate_single_level |
---|
| 182 | |
---|
| 183 | |
---|
| 184 | !--------------------------------------------------------------------- |
---|
| 185 | ! Unimaginative initialization of random-number seeds |
---|
| 186 | subroutine init_seed_simple(this, istartcol, iendcol) |
---|
| 187 | class(single_level_type), intent(inout) :: this |
---|
| 188 | integer, intent(in) :: istartcol, iendcol |
---|
| 189 | |
---|
| 190 | integer :: jcol |
---|
| 191 | |
---|
| 192 | if (.not. allocated(this%iseed)) then |
---|
| 193 | allocate(this%iseed(istartcol:iendcol)) |
---|
| 194 | end if |
---|
| 195 | |
---|
| 196 | do jcol = istartcol,iendcol |
---|
| 197 | this%iseed(jcol) = jcol |
---|
| 198 | end do |
---|
| 199 | end subroutine init_seed_simple |
---|
| 200 | |
---|
| 201 | |
---|
| 202 | !--------------------------------------------------------------------- |
---|
| 203 | ! Extract the shortwave and longwave surface albedos in each g-point |
---|
| 204 | subroutine get_albedos(this, istartcol, iendcol, config, & |
---|
| 205 | & sw_albedo_direct, sw_albedo_diffuse, lw_albedo) |
---|
| 206 | |
---|
| 207 | use radiation_config, only : config_type |
---|
| 208 | use radiation_io, only : nulerr, radiation_abort |
---|
| 209 | use yomhook, only : lhook, dr_hook |
---|
| 210 | |
---|
| 211 | class(single_level_type), intent(in) :: this |
---|
| 212 | type(config_type), intent(in) :: config |
---|
| 213 | integer, intent(in) :: istartcol, iendcol |
---|
| 214 | |
---|
| 215 | ! The longwave albedo of the surface in each longwave g-point; |
---|
| 216 | ! note that these are weighted averages of the values from |
---|
| 217 | ! individual tiles |
---|
| 218 | real(jprb), intent(out), optional & |
---|
| 219 | & :: lw_albedo(config%n_g_lw, istartcol:iendcol) |
---|
| 220 | |
---|
| 221 | ! Direct and diffuse shortwave surface albedo in each shortwave |
---|
| 222 | ! g-point; note that these are weighted averages of the values |
---|
| 223 | ! from individual tiles |
---|
| 224 | real(jprb), intent(out), dimension(config%n_g_sw, istartcol:iendcol) & |
---|
| 225 | & :: sw_albedo_direct, sw_albedo_diffuse |
---|
| 226 | |
---|
| 227 | ! Temporary storage of albedo in ecRad bands |
---|
| 228 | real(jprb) :: sw_albedo_band(istartcol:iendcol, config%n_bands_sw) |
---|
| 229 | real(jprb) :: lw_albedo_band (istartcol:iendcol, config%n_bands_lw) |
---|
| 230 | |
---|
| 231 | ! Number of albedo bands |
---|
| 232 | integer :: nalbedoband |
---|
| 233 | |
---|
| 234 | ! Loop indices for ecRad bands and albedo bands |
---|
| 235 | integer :: jband, jalbedoband |
---|
| 236 | |
---|
| 237 | real(jprb) :: hook_handle |
---|
| 238 | |
---|
| 239 | if (lhook) call dr_hook('radiation_single_level:get_albedos',0,hook_handle) |
---|
| 240 | |
---|
| 241 | ! Albedos/emissivities are stored in single_level in their own |
---|
| 242 | ! spectral intervals and with column as the first dimension |
---|
| 243 | if (config%use_canopy_full_spectrum_sw) then |
---|
| 244 | ! Albedos provided in each g point |
---|
| 245 | sw_albedo_diffuse = transpose(this%sw_albedo(istartcol:iendcol,:)) |
---|
| 246 | if (allocated(this%sw_albedo_direct)) then |
---|
| 247 | sw_albedo_direct = transpose(this%sw_albedo_direct(istartcol:iendcol,:)) |
---|
| 248 | end if |
---|
| 249 | elseif (.not. config%do_nearest_spectral_sw_albedo) then |
---|
| 250 | ! Albedos averaged accurately to ecRad spectral bands |
---|
| 251 | nalbedoband = size(config%sw_albedo_weights,1) |
---|
| 252 | sw_albedo_band = 0.0_jprb |
---|
| 253 | do jband = 1,config%n_bands_sw |
---|
| 254 | do jalbedoband = 1,nalbedoband |
---|
| 255 | if (config%sw_albedo_weights(jalbedoband,jband) /= 0.0_jprb) then |
---|
| 256 | sw_albedo_band(istartcol:iendcol,jband) & |
---|
| 257 | & = sw_albedo_band(istartcol:iendcol,jband) & |
---|
| 258 | & + config%sw_albedo_weights(jalbedoband,jband) & |
---|
| 259 | & * this%sw_albedo(istartcol:iendcol, jalbedoband) |
---|
| 260 | end if |
---|
| 261 | end do |
---|
| 262 | end do |
---|
| 263 | |
---|
| 264 | sw_albedo_diffuse = transpose(sw_albedo_band(istartcol:iendcol, & |
---|
| 265 | & config%i_band_from_reordered_g_sw)) |
---|
| 266 | if (allocated(this%sw_albedo_direct)) then |
---|
| 267 | sw_albedo_band = 0.0_jprb |
---|
| 268 | do jband = 1,config%n_bands_sw |
---|
| 269 | do jalbedoband = 1,nalbedoband |
---|
| 270 | if (config%sw_albedo_weights(jalbedoband,jband) /= 0.0_jprb) then |
---|
| 271 | sw_albedo_band(istartcol:iendcol,jband) & |
---|
| 272 | & = sw_albedo_band(istartcol:iendcol,jband) & |
---|
| 273 | & + config%sw_albedo_weights(jalbedoband,jband) & |
---|
| 274 | & * this%sw_albedo_direct(istartcol:iendcol, jalbedoband) |
---|
| 275 | end if |
---|
| 276 | end do |
---|
| 277 | end do |
---|
| 278 | sw_albedo_direct = transpose(sw_albedo_band(istartcol:iendcol, & |
---|
| 279 | & config%i_band_from_reordered_g_sw)) |
---|
| 280 | else |
---|
| 281 | sw_albedo_direct = sw_albedo_diffuse |
---|
| 282 | end if |
---|
| 283 | else |
---|
| 284 | ! Albedos mapped less accurately to ecRad spectral bands |
---|
| 285 | sw_albedo_diffuse = transpose(this%sw_albedo(istartcol:iendcol, & |
---|
| 286 | & config%i_albedo_from_band_sw(config%i_band_from_reordered_g_sw))) |
---|
| 287 | if (allocated(this%sw_albedo_direct)) then |
---|
| 288 | sw_albedo_direct = transpose(this%sw_albedo_direct(istartcol:iendcol, & |
---|
| 289 | & config%i_albedo_from_band_sw(config%i_band_from_reordered_g_sw))) |
---|
| 290 | else |
---|
| 291 | sw_albedo_direct = sw_albedo_diffuse |
---|
| 292 | end if |
---|
| 293 | end if |
---|
| 294 | |
---|
| 295 | if (present(lw_albedo)) then |
---|
| 296 | if (config%use_canopy_full_spectrum_lw) then |
---|
| 297 | if (config%n_g_lw /= size(this%lw_emissivity,2)) then |
---|
| 298 | write(nulerr,'(a)') '*** Error: single_level%lw_emissivity has the wrong number of spectral intervals' |
---|
| 299 | call radiation_abort() |
---|
| 300 | end if |
---|
| 301 | lw_albedo = 1.0_jprb - transpose(this%lw_emissivity(istartcol:iendcol,:)) |
---|
| 302 | else if (.not. config%do_nearest_spectral_lw_emiss) then |
---|
| 303 | ! Albedos averaged accurately to ecRad spectral bands |
---|
| 304 | nalbedoband = size(config%lw_emiss_weights,1) |
---|
| 305 | lw_albedo_band = 0.0_jprb |
---|
| 306 | do jband = 1,config%n_bands_lw |
---|
| 307 | do jalbedoband = 1,nalbedoband |
---|
| 308 | if (config%lw_emiss_weights(jalbedoband,jband) /= 0.0_jprb) then |
---|
| 309 | lw_albedo_band(istartcol:iendcol,jband) & |
---|
| 310 | & = lw_albedo_band(istartcol:iendcol,jband) & |
---|
| 311 | & + config%lw_emiss_weights(jalbedoband,jband) & |
---|
| 312 | & * (1.0_jprb-this%lw_emissivity(istartcol:iendcol, jalbedoband)) |
---|
| 313 | end if |
---|
| 314 | end do |
---|
| 315 | end do |
---|
| 316 | |
---|
| 317 | lw_albedo = transpose(lw_albedo_band(istartcol:iendcol, & |
---|
| 318 | & config%i_band_from_reordered_g_lw)) |
---|
| 319 | else |
---|
| 320 | lw_albedo = 1.0_jprb - transpose(this%lw_emissivity(istartcol:iendcol, & |
---|
| 321 | & config%i_emiss_from_band_lw(config%i_band_from_reordered_g_lw))) |
---|
| 322 | end if |
---|
| 323 | end if |
---|
| 324 | |
---|
| 325 | if (lhook) call dr_hook('radiation_single_level:get_albedos',1,hook_handle) |
---|
| 326 | |
---|
| 327 | end subroutine get_albedos |
---|
| 328 | |
---|
| 329 | |
---|
| 330 | !--------------------------------------------------------------------- |
---|
| 331 | ! Return .true. if the contents of a single_level structure are out |
---|
| 332 | ! of a physically sensible range, optionally only considering only |
---|
| 333 | ! columns between istartcol and iendcol |
---|
| 334 | function out_of_physical_bounds(this, istartcol, iendcol, do_fix) result(is_bad) |
---|
| 335 | |
---|
| 336 | use yomhook, only : lhook, dr_hook |
---|
| 337 | use radiation_config, only : out_of_bounds_1d, out_of_bounds_2d |
---|
| 338 | |
---|
| 339 | class(single_level_type), intent(inout) :: this |
---|
| 340 | integer, optional,intent(in) :: istartcol, iendcol |
---|
| 341 | logical, optional,intent(in) :: do_fix |
---|
| 342 | logical :: is_bad |
---|
| 343 | |
---|
| 344 | logical :: do_fix_local |
---|
| 345 | |
---|
| 346 | real(jprb) :: hook_handle |
---|
| 347 | |
---|
| 348 | if (lhook) call dr_hook('radiation_single_level:out_of_physical_bounds',0,hook_handle) |
---|
| 349 | |
---|
| 350 | if (present(do_fix)) then |
---|
| 351 | do_fix_local = do_fix |
---|
| 352 | else |
---|
| 353 | do_fix_local = .false. |
---|
| 354 | end if |
---|
| 355 | |
---|
| 356 | is_bad = out_of_bounds_1d(this%cos_sza, "cos_sza", -1.0_jprb, 1.0_jprb, & |
---|
| 357 | & do_fix_local, i1=istartcol, i2=iendcol) & |
---|
| 358 | & .or. out_of_bounds_1d(this%skin_temperature, "skin_temperature", 173.0_jprb, 373.0_jprb, & |
---|
| 359 | & do_fix_local, i1=istartcol, i2=iendcol) & |
---|
| 360 | & .or. out_of_bounds_2d(this%sw_albedo, "sw_albedo", 0.0_jprb, 1.0_jprb, & |
---|
| 361 | & do_fix_local, i1=istartcol, i2=iendcol) & |
---|
| 362 | & .or. out_of_bounds_2d(this%sw_albedo_direct, "sw_albedo", 0.0_jprb, 1.0_jprb, & |
---|
| 363 | & do_fix_local, i1=istartcol, i2=iendcol) & |
---|
| 364 | & .or. out_of_bounds_2d(this%lw_emissivity, "lw_emissivity", 0.0_jprb, 1.0_jprb, & |
---|
| 365 | & do_fix_local, i1=istartcol, i2=iendcol) |
---|
| 366 | |
---|
| 367 | if (lhook) call dr_hook('radiation_single_level:out_of_physical_bounds',1,hook_handle) |
---|
| 368 | |
---|
| 369 | end function out_of_physical_bounds |
---|
| 370 | |
---|
| 371 | end module radiation_single_level |
---|