1 2! KGEN-generated Fortran source file 3! 4! Filename : rrtmg_lw_rtrnmc.f90 5! Generated at: 2015-07-06 23:28:45 6! KGEN version: 0.4.13 7 8 MODULE rrtmg_lw_rtrnmc 9 USE kgen_utils_mod, ONLY : kgen_dp, check_t, kgen_init_check, kgen_print_check 10 ! -------------------------------------------------------------------------- 11 ! | | 12 ! | Copyright 2002-2007, Atmospheric & Environmental Research, Inc. (AER). | 13 ! | This software may be used, copied, or redistributed as long as it is | 14 ! | not sold and this copyright notice is reproduced on each copy made. | 15 ! | This model is provided as is without any express or implied warranties. | 16 ! | (http://www.rtweb.aer.com/) | 17 ! | | 18 ! -------------------------------------------------------------------------- 19 ! --------- Modules ---------- 20 USE shr_kind_mod, ONLY: r8 => shr_kind_r8 21 ! use parkind, only : jpim, jprb 22 USE parrrtm, ONLY: ngptlw 23 USE parrrtm, ONLY: nbndlw 24 USE rrlw_con, ONLY: fluxfac 25 USE rrlw_con, ONLY: heatfac 26 USE rrlw_wvn, ONLY: ngb 27 USE rrlw_wvn, ONLY: ngs 28 USE rrlw_wvn, ONLY: delwave 29 USE rrlw_tbl, ONLY: bpade 30 USE rrlw_tbl, ONLY: tblint 31 USE rrlw_tbl, ONLY: tfn_tbl 32 USE rrlw_tbl, ONLY: exp_tbl 33 USE rrlw_tbl, ONLY: tau_tbl 34 USE rrlw_vsn, ONLY: hvrrtc 35 IMPLICIT NONE 36 37#ifdef OLD_RTRNMC 38 public :: rtrnmc_old 39#else 40 public :: rtrnmc 41#endif 42 CONTAINS 43 44 ! write subroutines 45 ! No subroutines 46 ! No module extern variables 47 !----------------------------------------------------------------------------- 48 49#ifdef OLD_RTRNMC 50 SUBROUTINE rtrnmc_old(nlayers, istart, iend, iout, pz, semiss, ncbands, cldfmc, taucmc, planklay, planklev, plankbnd, pwvcm, & 51 fracs, taut, totuflux, totdflux, fnet, htr, totuclfl, totdclfl, fnetc, htrc, totufluxs, totdfluxs) 52 !----------------------------------------------------------------------------- 53 ! 54 ! Original version: E. J. Mlawer, et al. RRTM_V3.0 55 ! Revision for GCMs: Michael J. Iacono; October, 2002 56 ! Revision for F90: Michael J. Iacono; June, 2006 57 ! 58 ! This program calculates the upward fluxes, downward fluxes, and 59 ! heating rates for an arbitrary clear or cloudy atmosphere. The input 60 ! to this program is the atmospheric profile, all Planck function 61 ! information, and the cloud fraction by layer. A variable diffusivity 62 ! angle (SECDIFF) is used for the angle integration. Bands 2-3 and 5-9 63 ! use a value for SECDIFF that varies from 1.50 to 1.80 as a function of 64 ! the column water vapor, and other bands use a value of 1.66. The Gaussian 65 ! weight appropriate to this angle (WTDIFF=0.5) is applied here. Note that 66 ! use of the emissivity angle for the flux integration can cause errors of 67 ! 1 to 4 W/m2 within cloudy layers. 68 ! Clouds are treated with the McICA stochastic approach and maximum-random 69 ! cloud overlap. 70 !*************************************************************************** 71 ! ------- Declarations ------- 72 ! ----- Input ----- 73 INTEGER, intent(in) :: nlayers ! total number of layers 74 INTEGER, intent(in) :: istart ! beginning band of calculation 75 INTEGER, intent(in) :: iend ! ending band of calculation 76 INTEGER, intent(in) :: iout ! output option flag 77 ! Atmosphere 78 REAL(KIND=r8), intent(in) :: pz(0:) ! level (interface) pressures (hPa, mb) 79 ! Dimensions: (0:nlayers) 80 REAL(KIND=r8), intent(in) :: pwvcm ! precipitable water vapor (cm) 81 REAL(KIND=r8), intent(in) :: semiss(:) ! lw surface emissivity 82 ! Dimensions: (nbndlw) 83 REAL(KIND=r8), intent(in) :: planklay(:,:) ! 84 ! Dimensions: (nlayers,nbndlw) 85 REAL(KIND=r8), intent(in) :: planklev(0:,:) ! 86 ! Dimensions: (0:nlayers,nbndlw) 87 REAL(KIND=r8), intent(in) :: plankbnd(:) ! 88 ! Dimensions: (nbndlw) 89 REAL(KIND=r8), intent(in) :: fracs(:,:) ! 90 ! Dimensions: (nlayers,ngptw) 91 REAL(KIND=r8), intent(in) :: taut(:,:) ! gaseous + aerosol optical depths 92 ! Dimensions: (nlayers,ngptlw) 93 ! Clouds 94 INTEGER, intent(in) :: ncbands ! number of cloud spectral bands 95 REAL(KIND=r8), intent(in) :: cldfmc(:,:) ! layer cloud fraction [mcica] 96 ! Dimensions: (ngptlw,nlayers) 97 REAL(KIND=r8), intent(in) :: taucmc(:,:) ! layer cloud optical depth [mcica] 98 ! Dimensions: (ngptlw,nlayers) 99 ! ----- Output ----- 100 REAL(KIND=r8), intent(out) :: totuflux(0:) ! upward longwave flux (w/m2) 101 ! Dimensions: (0:nlayers) 102 REAL(KIND=r8), intent(out) :: totdflux(0:) ! downward longwave flux (w/m2) 103 ! Dimensions: (0:nlayers) 104 REAL(KIND=r8), intent(out) :: fnet(0:) ! net longwave flux (w/m2) 105 ! Dimensions: (0:nlayers) 106 REAL(KIND=r8), intent(out) :: htr(0:) ! longwave heating rate (k/day) 107 ! Dimensions: (0:nlayers) 108 REAL(KIND=r8), intent(out) :: totuclfl(0:) ! clear sky upward longwave flux (w/m2) 109 ! Dimensions: (0:nlayers) 110 REAL(KIND=r8), intent(out) :: totdclfl(0:) ! clear sky downward longwave flux (w/m2) 111 ! Dimensions: (0:nlayers) 112 REAL(KIND=r8), intent(out) :: fnetc(0:) ! clear sky net longwave flux (w/m2) 113 ! Dimensions: (0:nlayers) 114 REAL(KIND=r8), intent(out) :: htrc(0:) ! clear sky longwave heating rate (k/day) 115 ! Dimensions: (0:nlayers) 116 REAL(KIND=r8), intent(out) :: totufluxs(:,0:) ! upward longwave flux spectral (w/m2) 117 ! Dimensions: (nbndlw, 0:nlayers) 118 REAL(KIND=r8), intent(out) :: totdfluxs(:,0:) ! downward longwave flux spectral (w/m2) 119 ! Dimensions: (nbndlw, 0:nlayers) 120 ! ----- Local ----- 121 ! Declarations for radiative transfer 122 REAL(KIND=r8) :: abscld(nlayers,ngptlw) 123 REAL(KIND=r8) :: atot(nlayers) 124 REAL(KIND=r8) :: atrans(nlayers) 125 REAL(KIND=r8) :: bbugas(nlayers) 126 REAL(KIND=r8) :: bbutot(nlayers) 127 REAL(KIND=r8) :: clrurad(0:nlayers) 128 REAL(KIND=r8) :: clrdrad(0:nlayers) 129 REAL(KIND=r8) :: efclfrac(nlayers,ngptlw) 130 REAL(KIND=r8) :: uflux(0:nlayers) 131 REAL(KIND=r8) :: dflux(0:nlayers) 132 REAL(KIND=r8) :: urad(0:nlayers) 133 REAL(KIND=r8) :: drad(0:nlayers) 134 REAL(KIND=r8) :: uclfl(0:nlayers) 135 REAL(KIND=r8) :: dclfl(0:nlayers) 136 REAL(KIND=r8) :: odcld(nlayers,ngptlw) 137 REAL(KIND=r8) :: secdiff(nbndlw) ! secant of diffusivity angle 138 REAL(KIND=r8) :: a0(nbndlw) 139 REAL(KIND=r8) :: a1(nbndlw) 140 REAL(KIND=r8) :: a2(nbndlw) ! diffusivity angle adjustment coefficients 141 REAL(KIND=r8) :: wtdiff 142 REAL(KIND=r8) :: rec_6 143 REAL(KIND=r8) :: transcld 144 REAL(KIND=r8) :: radld 145 REAL(KIND=r8) :: radclrd 146 REAL(KIND=r8) :: plfrac 147 REAL(KIND=r8) :: blay 148 REAL(KIND=r8) :: dplankup 149 REAL(KIND=r8) :: dplankdn 150 REAL(KIND=r8) :: odepth 151 REAL(KIND=r8) :: odtot 152 REAL(KIND=r8) :: odepth_rec 153 REAL(KIND=r8) :: gassrc 154 REAL(KIND=r8) :: odtot_rec 155 REAL(KIND=r8) :: bbdtot 156 REAL(KIND=r8) :: bbd 157 REAL(KIND=r8) :: tblind 158 REAL(KIND=r8) :: tfactot 159 REAL(KIND=r8) :: tfacgas 160 REAL(KIND=r8) :: transc 161 REAL(KIND=r8) :: tausfac 162 REAL(KIND=r8) :: rad0 163 REAL(KIND=r8) :: reflect 164 REAL(KIND=r8) :: radlu 165 REAL(KIND=r8) :: radclru 166 INTEGER :: icldlyr(nlayers) ! flag for cloud in layer 167 INTEGER :: ibnd 168 INTEGER :: lay 169 INTEGER :: ig 170 INTEGER :: ib 171 INTEGER :: iband 172 INTEGER :: lev 173 INTEGER :: l ! loop indices 174 INTEGER :: igc ! g-point interval counter 175 INTEGER :: iclddn ! flag for cloud in down path 176 INTEGER :: ittot 177 INTEGER :: itgas 178 INTEGER :: itr ! lookup table indices 179 ! ------- Definitions ------- 180 ! input 181 ! nlayers ! number of model layers 182 ! ngptlw ! total number of g-point subintervals 183 ! nbndlw ! number of longwave spectral bands 184 ! ncbands ! number of spectral bands for clouds 185 ! secdiff ! diffusivity angle 186 ! wtdiff ! weight for radiance to flux conversion 187 ! pavel ! layer pressures (mb) 188 ! pz ! level (interface) pressures (mb) 189 ! tavel ! layer temperatures (k) 190 ! tz ! level (interface) temperatures(mb) 191 ! tbound ! surface temperature (k) 192 ! cldfrac ! layer cloud fraction 193 ! taucloud ! layer cloud optical depth 194 ! itr ! integer look-up table index 195 ! icldlyr ! flag for cloudy layers 196 ! iclddn ! flag for cloud in column at any layer 197 ! semiss ! surface emissivities for each band 198 ! reflect ! surface reflectance 199 ! bpade ! 1/(pade constant) 200 ! tau_tbl ! clear sky optical depth look-up table 201 ! exp_tbl ! exponential look-up table for transmittance 202 ! tfn_tbl ! tau transition function look-up table 203 ! local 204 ! atrans ! gaseous absorptivity 205 ! abscld ! cloud absorptivity 206 ! atot ! combined gaseous and cloud absorptivity 207 ! odclr ! clear sky (gaseous) optical depth 208 ! odcld ! cloud optical depth 209 ! odtot ! optical depth of gas and cloud 210 ! tfacgas ! gas-only pade factor, used for planck fn 211 ! tfactot ! gas and cloud pade factor, used for planck fn 212 ! bbdgas ! gas-only planck function for downward rt 213 ! bbugas ! gas-only planck function for upward rt 214 ! bbdtot ! gas and cloud planck function for downward rt 215 ! bbutot ! gas and cloud planck function for upward calc. 216 ! gassrc ! source radiance due to gas only 217 ! efclfrac ! effective cloud fraction 218 ! radlu ! spectrally summed upward radiance 219 ! radclru ! spectrally summed clear sky upward radiance 220 ! urad ! upward radiance by layer 221 ! clrurad ! clear sky upward radiance by layer 222 ! radld ! spectrally summed downward radiance 223 ! radclrd ! spectrally summed clear sky downward radiance 224 ! drad ! downward radiance by layer 225 ! clrdrad ! clear sky downward radiance by layer 226 ! output 227 ! totuflux ! upward longwave flux (w/m2) 228 ! totdflux ! downward longwave flux (w/m2) 229 ! fnet ! net longwave flux (w/m2) 230 ! htr ! longwave heating rate (k/day) 231 ! totuclfl ! clear sky upward longwave flux (w/m2) 232 ! totdclfl ! clear sky downward longwave flux (w/m2) 233 ! fnetc ! clear sky net longwave flux (w/m2) 234 ! htrc ! clear sky longwave heating rate (k/day) 235 ! This secant and weight corresponds to the standard diffusivity 236 ! angle. This initial value is redefined below for some bands. 237 data wtdiff /0.5_r8/ 238 data rec_6 /0.166667_r8/ 239 ! Reset diffusivity angle for Bands 2-3 and 5-9 to vary (between 1.50 240 ! and 1.80) as a function of total column water vapor. The function 241 ! has been defined to minimize flux and cooling rate errors in these bands 242 ! over a wide range of precipitable water values. 243 data a0 / 1.66_r8, 1.55_r8, 1.58_r8, 1.66_r8, & 244 1.54_r8, 1.454_r8, 1.89_r8, 1.33_r8, & 245 1.668_r8, 1.66_r8, 1.66_r8, 1.66_r8, & 246 1.66_r8, 1.66_r8, 1.66_r8, 1.66_r8 / 247 data a1 / 0.00_r8, 0.25_r8, 0.22_r8, 0.00_r8, & 248 0.13_r8, 0.446_r8, -0.10_r8, 0.40_r8, & 249 -0.006_r8, 0.00_r8, 0.00_r8, 0.00_r8, & 250 0.00_r8, 0.00_r8, 0.00_r8, 0.00_r8 / 251 data a2 / 0.00_r8, -12.0_r8, -11.7_r8, 0.00_r8, & 252 -0.72_r8,-0.243_r8, 0.19_r8,-0.062_r8, & 253 0.414_r8, 0.00_r8, 0.00_r8, 0.00_r8, & 254 0.00_r8, 0.00_r8, 0.00_r8, 0.00_r8 / 255 hvrrtc = '$Revision: 1.3 $' 256 do ibnd = 1,nbndlw 257 if (ibnd.eq.1 .or. ibnd.eq.4 .or. ibnd.ge.10) then 258 secdiff(ibnd) = 1.66_r8 259 else 260 secdiff(ibnd) = a0(ibnd) + a1(ibnd)*exp(a2(ibnd)*pwvcm) 261 endif 262 enddo 263 if (pwvcm.lt.1.0) secdiff(6) = 1.80_r8 264 if (pwvcm.gt.7.1) secdiff(7) = 1.50_r8 265 urad(0) = 0.0_r8 266 drad(0) = 0.0_r8 267 totuflux(0) = 0.0_r8 268 totdflux(0) = 0.0_r8 269 clrurad(0) = 0.0_r8 270 clrdrad(0) = 0.0_r8 271 totuclfl(0) = 0.0_r8 272 totdclfl(0) = 0.0_r8 273 do lay = 1, nlayers 274 urad(lay) = 0.0_r8 275 drad(lay) = 0.0_r8 276 totuflux(lay) = 0.0_r8 277 totdflux(lay) = 0.0_r8 278 clrurad(lay) = 0.0_r8 279 clrdrad(lay) = 0.0_r8 280 totuclfl(lay) = 0.0_r8 281 totdclfl(lay) = 0.0_r8 282 icldlyr(lay) = 0 283 ! Change to band loop? 284 do ig = 1, ngptlw 285 if (cldfmc(ig,lay) .eq. 1._r8) then 286 ib = ngb(ig) 287 odcld(lay,ig) = secdiff(ib) * taucmc(ig,lay) 288 transcld = exp(-odcld(lay,ig)) 289 abscld(lay,ig) = 1._r8 - transcld 290 efclfrac(lay,ig) = abscld(lay,ig) * cldfmc(ig,lay) 291 icldlyr(lay) = 1 292 else 293 odcld(lay,ig) = 0.0_r8 294 abscld(lay,ig) = 0.0_r8 295 efclfrac(lay,ig) = 0.0_r8 296 endif 297 enddo 298 enddo 299 igc = 1 300 ! Loop over frequency bands. 301 do iband = istart, iend 302 ! Reinitialize g-point counter for each band if output for each band is requested. 303 if (iout.gt.0.and.iband.ge.2) igc = ngs(iband-1)+1 304 ! Loop over g-channels. 305 1000 continue 306 ! Radiative transfer starts here. 307 radld = 0._r8 308 radclrd = 0._r8 309 iclddn = 0 310 ! Downward radiative transfer loop. 311 do lev = nlayers, 1, -1 312 plfrac = fracs(lev,igc) 313 blay = planklay(lev,iband) 314 dplankup = planklev(lev,iband) - blay 315 dplankdn = planklev(lev-1,iband) - blay 316 odepth = secdiff(iband) * taut(lev,igc) 317 if (odepth .lt. 0.0_r8) odepth = 0.0_r8 318 ! Cloudy layer 319 if (icldlyr(lev).eq.1) then 320 iclddn = 1 321 odtot = odepth + odcld(lev,igc) 322 if (odtot .lt. 0.06_r8) then 323 atrans(lev) = odepth - 0.5_r8*odepth*odepth 324 odepth_rec = rec_6*odepth 325 gassrc = plfrac*(blay+dplankdn*odepth_rec)*atrans(lev) 326 atot(lev) = odtot - 0.5_r8*odtot*odtot 327 odtot_rec = rec_6*odtot 328 bbdtot = plfrac * (blay+dplankdn*odtot_rec) 329 bbd = plfrac*(blay+dplankdn*odepth_rec) 330 radld = radld - radld * (atrans(lev) + & 331 efclfrac(lev,igc) * (1. - atrans(lev))) + & 332 gassrc + cldfmc(igc,lev) * & 333 (bbdtot * atot(lev) - gassrc) 334 drad(lev-1) = drad(lev-1) + radld 335 bbugas(lev) = plfrac * (blay+dplankup*odepth_rec) 336 bbutot(lev) = plfrac * (blay+dplankup*odtot_rec) 337 elseif (odepth .le. 0.06_r8) then 338 atrans(lev) = odepth - 0.5_r8*odepth*odepth 339 odepth_rec = rec_6*odepth 340 gassrc = plfrac*(blay+dplankdn*odepth_rec)*atrans(lev) 341 odtot = odepth + odcld(lev,igc) 342 tblind = odtot/(bpade+odtot) 343 ittot = tblint*tblind + 0.5_r8 344 tfactot = tfn_tbl(ittot) 345 bbdtot = plfrac * (blay + tfactot*dplankdn) 346 bbd = plfrac*(blay+dplankdn*odepth_rec) 347 atot(lev) = 1. - exp_tbl(ittot) 348 radld = radld - radld * (atrans(lev) + & 349 efclfrac(lev,igc) * (1._r8 - atrans(lev))) + & 350 gassrc + cldfmc(igc,lev) * & 351 (bbdtot * atot(lev) - gassrc) 352 drad(lev-1) = drad(lev-1) + radld 353 bbugas(lev) = plfrac * (blay + dplankup*odepth_rec) 354 bbutot(lev) = plfrac * (blay + tfactot * dplankup) 355 else 356 tblind = odepth/(bpade+odepth) 357 itgas = tblint*tblind+0.5_r8 358 odepth = tau_tbl(itgas) 359 atrans(lev) = 1._r8 - exp_tbl(itgas) 360 tfacgas = tfn_tbl(itgas) 361 gassrc = atrans(lev) * plfrac * (blay + tfacgas*dplankdn) 362 odtot = odepth + odcld(lev,igc) 363 tblind = odtot/(bpade+odtot) 364 ittot = tblint*tblind + 0.5_r8 365 tfactot = tfn_tbl(ittot) 366 bbdtot = plfrac * (blay + tfactot*dplankdn) 367 bbd = plfrac*(blay+tfacgas*dplankdn) 368 atot(lev) = 1._r8 - exp_tbl(ittot) 369 radld = radld - radld * (atrans(lev) + & 370 efclfrac(lev,igc) * (1._r8 - atrans(lev))) + & 371 gassrc + cldfmc(igc,lev) * & 372 (bbdtot * atot(lev) - gassrc) 373 drad(lev-1) = drad(lev-1) + radld 374 bbugas(lev) = plfrac * (blay + tfacgas * dplankup) 375 bbutot(lev) = plfrac * (blay + tfactot * dplankup) 376 endif 377 ! Clear layer 378 else 379 if (odepth .le. 0.06_r8) then 380 atrans(lev) = odepth-0.5_r8*odepth*odepth 381 odepth = rec_6*odepth 382 bbd = plfrac*(blay+dplankdn*odepth) 383 bbugas(lev) = plfrac*(blay+dplankup*odepth) 384 else 385 tblind = odepth/(bpade+odepth) 386 itr = tblint*tblind+0.5_r8 387 transc = exp_tbl(itr) 388 atrans(lev) = 1._r8-transc 389 tausfac = tfn_tbl(itr) 390 bbd = plfrac*(blay+tausfac*dplankdn) 391 bbugas(lev) = plfrac * (blay + tausfac * dplankup) 392 endif 393 radld = radld + (bbd-radld)*atrans(lev) 394 drad(lev-1) = drad(lev-1) + radld 395 endif 396 ! Set clear sky stream to total sky stream as long as layers 397 ! remain clear. Streams diverge when a cloud is reached (iclddn=1), 398 ! and clear sky stream must be computed separately from that point. 399 if (iclddn.eq.1) then 400 radclrd = radclrd + (bbd-radclrd) * atrans(lev) 401 clrdrad(lev-1) = clrdrad(lev-1) + radclrd 402 else 403 radclrd = radld 404 clrdrad(lev-1) = drad(lev-1) 405 endif 406 enddo 407 ! Spectral emissivity & reflectance 408 ! Include the contribution of spectrally varying longwave emissivity 409 ! and reflection from the surface to the upward radiative transfer. 410 ! Note: Spectral and Lambertian reflection are identical for the 411 ! diffusivity angle flux integration used here. 412 rad0 = fracs(1,igc) * plankbnd(iband) 413 ! Add in specular reflection of surface downward radiance. 414 reflect = 1._r8 - semiss(iband) 415 radlu = rad0 + reflect * radld 416 radclru = rad0 + reflect * radclrd 417 ! Upward radiative transfer loop. 418 urad(0) = urad(0) + radlu 419 clrurad(0) = clrurad(0) + radclru 420 do lev = 1, nlayers 421 ! Cloudy layer 422 if (icldlyr(lev) .eq. 1) then 423 gassrc = bbugas(lev) * atrans(lev) 424 radlu = radlu - radlu * (atrans(lev) + & 425 efclfrac(lev,igc) * (1._r8 - atrans(lev))) + & 426 gassrc + cldfmc(igc,lev) * & 427 (bbutot(lev) * atot(lev) - gassrc) 428 urad(lev) = urad(lev) + radlu 429 ! Clear layer 430 else 431 radlu = radlu + (bbugas(lev)-radlu)*atrans(lev) 432 urad(lev) = urad(lev) + radlu 433 endif 434 ! Set clear sky stream to total sky stream as long as all layers 435 ! are clear (iclddn=0). Streams must be calculated separately at 436 ! all layers when a cloud is present (ICLDDN=1), because surface 437 ! reflectance is different for each stream. 438 if (iclddn.eq.1) then 439 radclru = radclru + (bbugas(lev)-radclru)*atrans(lev) 440 clrurad(lev) = clrurad(lev) + radclru 441 else 442 radclru = radlu 443 clrurad(lev) = urad(lev) 444 endif 445 enddo 446 ! Increment g-point counter 447 igc = igc + 1 448 ! Return to continue radiative transfer for all g-channels in present band 449 if (igc .le. ngs(iband)) go to 1000 450 ! Process longwave output from band for total and clear streams. 451 ! Calculate upward, downward, and net flux. 452 do lev = nlayers, 0, -1 453 uflux(lev) = urad(lev)*wtdiff 454 dflux(lev) = drad(lev)*wtdiff 455 urad(lev) = 0.0_r8 456 drad(lev) = 0.0_r8 457 totuflux(lev) = totuflux(lev) + uflux(lev) * delwave(iband) 458 totdflux(lev) = totdflux(lev) + dflux(lev) * delwave(iband) 459 uclfl(lev) = clrurad(lev)*wtdiff 460 dclfl(lev) = clrdrad(lev)*wtdiff 461 clrurad(lev) = 0.0_r8 462 clrdrad(lev) = 0.0_r8 463 totuclfl(lev) = totuclfl(lev) + uclfl(lev) * delwave(iband) 464 totdclfl(lev) = totdclfl(lev) + dclfl(lev) * delwave(iband) 465 totufluxs(iband,lev) = uflux(lev) * delwave(iband) 466 totdfluxs(iband,lev) = dflux(lev) * delwave(iband) 467 enddo 468 ! End spectral band loop 469 enddo 470 ! Calculate fluxes at surface 471 totuflux(0) = totuflux(0) * fluxfac 472 totdflux(0) = totdflux(0) * fluxfac 473 totufluxs(:,0) = totufluxs(:,0) * fluxfac 474 totdfluxs(:,0) = totdfluxs(:,0) * fluxfac 475 fnet(0) = totuflux(0) - totdflux(0) 476 totuclfl(0) = totuclfl(0) * fluxfac 477 totdclfl(0) = totdclfl(0) * fluxfac 478 fnetc(0) = totuclfl(0) - totdclfl(0) 479 ! Calculate fluxes at model levels 480 do lev = 1, nlayers 481 totuflux(lev) = totuflux(lev) * fluxfac 482 totdflux(lev) = totdflux(lev) * fluxfac 483 totufluxs(:,lev) = totufluxs(:,lev) * fluxfac 484 totdfluxs(:,lev) = totdfluxs(:,lev) * fluxfac 485 fnet(lev) = totuflux(lev) - totdflux(lev) 486 totuclfl(lev) = totuclfl(lev) * fluxfac 487 totdclfl(lev) = totdclfl(lev) * fluxfac 488 fnetc(lev) = totuclfl(lev) - totdclfl(lev) 489 l = lev - 1 490 ! Calculate heating rates at model layers 491 htr(l)=heatfac*(fnet(l)-fnet(lev))/(pz(l)-pz(lev)) 492 htrc(l)=heatfac*(fnetc(l)-fnetc(lev))/(pz(l)-pz(lev)) 493 enddo 494 ! Set heating rate to zero in top layer 495 htr(nlayers) = 0.0_r8 496 htrc(nlayers) = 0.0_r8 497 END SUBROUTINE rtrnmc_old 498#else 499 SUBROUTINE rtrnmc(ncol,nlayers, istart, iend, iout, pz, semiss, ncbands, cldfmc, taucmc, planklay, planklev, plankbnd, pwvcm, & 500 fracs, taut, totuflux, totdflux, fnet, htr, totuclfl, totdclfl, fnetc, htrc, totufluxs, totdfluxs) 501 !----------------------------------------------------------------------------- 502 ! 503 ! Original version: E. J. Mlawer, et al. RRTM_V3.0 504 ! Revision for GCMs: Michael J. Iacono; October, 2002 505 ! Revision for F90: Michael J. Iacono; June, 2006 506 ! 507 ! This program calculates the upward fluxes, downward fluxes, and 508 ! heating rates for an arbitrary clear or cloudy atmosphere. The input 509 ! to this program is the atmospheric profile, all Planck function 510 ! information, and the cloud fraction by layer. A variable diffusivity 511 ! angle (SECDIFF) is used for the angle integration. Bands 2-3 and 5-9 512 ! use a value for SECDIFF that varies from 1.50 to 1.80 as a function of 513 ! the column water vapor, and other bands use a value of 1.66. The Gaussian 514 ! weight appropriate to this angle (WTDIFF=0.5) is applied here. Note that 515 ! use of the emissivity angle for the flux integration can cause errors of 516 ! 1 to 4 W/m2 within cloudy layers. 517 ! Clouds are treated with the McICA stochastic approach and maximum-random 518 ! cloud overlap. 519 !*************************************************************************** 520 ! ------- Declarations ------- 521 ! ----- Input ----- 522 INTEGER, intent(in) :: ncol ! total number of columns 523 INTEGER, intent(in) :: nlayers ! total number of layers 524 INTEGER, intent(in) :: istart ! beginning band of calculation 525 INTEGER, intent(in) :: iend ! ending band of calculation 526 INTEGER, intent(in) :: iout ! output option flag 527 ! Atmosphere 528 REAL(KIND=r8), intent(in) :: pz(:,0:) ! level (interface) pressures (hPa, mb) 529 ! Dimensions: (ncol,0:nlayers) 530 REAL(KIND=r8), intent(in) :: pwvcm(:) ! precipitable water vapor (cm) 531 ! Dimensions: (ncol) 532 REAL(KIND=r8), intent(in) :: semiss(:,:) ! lw surface emissivity 533 ! Dimensions: (ncol,nbndlw) 534 REAL(KIND=r8), intent(in) :: planklay(:,:,:) ! 535 ! Dimensions: (ncol,nlayers,nbndlw) 536 REAL(KIND=r8), intent(in) :: planklev(:,0:,:) ! 537 ! Dimensions: (ncol,0:nlayers,nbndlw) 538 REAL(KIND=r8), intent(in) :: plankbnd(:,:) ! 539 ! Dimensions: (ncol,nbndlw) 540 REAL(KIND=r8), intent(in) :: fracs(:,:,:) ! 541 ! Dimensions: (ncol,nlayers,ngptw) 542 REAL(KIND=r8), intent(in) :: taut(:,:,:) ! gaseous + aerosol optical depths 543 ! Dimensions: (ncol,nlayers,ngptlw) 544 ! Clouds 545 INTEGER, intent(in) :: ncbands(:) ! number of cloud spectral bands 546 ! Dimensions: (ncol) 547 REAL(KIND=r8), intent(in) :: cldfmc(:,:,:) ! layer cloud fraction [mcica] 548 ! Dimensions: (ncol,ngptlw,nlayers) 549 REAL(KIND=r8), intent(in) :: taucmc(:,:,:) ! layer cloud optical depth [mcica] 550 ! Dimensions: (ncol,ngptlw,nlayers) 551 ! ----- Output ----- 552 REAL(KIND=r8), intent(out) :: totuflux(:,0:) ! upward longwave flux (w/m2) 553 ! Dimensions: (ncol,0:nlayers) 554 REAL(KIND=r8), intent(out) :: totdflux(:,0:) ! downward longwave flux (w/m2) 555 ! Dimensions: (ncol,0:nlayers) 556 REAL(KIND=r8), intent(out) :: fnet(:,0:) ! net longwave flux (w/m2) 557 ! Dimensions: (ncol,0:nlayers) 558 REAL(KIND=r8), intent(out) :: htr(:,0:) ! longwave heating rate (k/day) 559 ! Dimensions: (ncol,0:nlayers) 560 REAL(KIND=r8), intent(out) :: totuclfl(:,0:) ! clear sky upward longwave flux (w/m2) 561 ! Dimensions: (ncol,0:nlayers) 562 REAL(KIND=r8), intent(out) :: totdclfl(:,0:) ! clear sky downward longwave flux (w/m2) 563 ! Dimensions: (ncol,0:nlayers) 564 REAL(KIND=r8), intent(out) :: fnetc(:,0:) ! clear sky net longwave flux (w/m2) 565 ! Dimensions: (ncol,0:nlayers) 566 REAL(KIND=r8), intent(out) :: htrc(:,0:) ! clear sky longwave heating rate (k/day) 567 ! Dimensions: (ncol,0:nlayers) 568 REAL(KIND=r8), intent(out) :: totufluxs(:,:,0:) ! upward longwave flux spectral (w/m2) 569 ! Dimensions: (ncol,nbndlw, 0:nlayers) 570 REAL(KIND=r8), intent(out) :: totdfluxs(:,:,0:) ! downward longwave flux spectral (w/m2) 571 ! Dimensions: (ncol,nbndlw, 0:nlayers) 572 ! ----- Local ----- 573 ! Declarations for radiative transfer 574 REAL(KIND=r8) :: abscld(nlayers,ngptlw) 575 REAL(KIND=r8) :: atot(nlayers) 576 REAL(KIND=r8) :: atrans(nlayers) 577 REAL(KIND=r8) :: bbugas(nlayers) 578 REAL(KIND=r8) :: bbutot(nlayers) 579 REAL(KIND=r8) :: clrurad(0:nlayers) 580 REAL(KIND=r8) :: clrdrad(0:nlayers) 581 REAL(KIND=r8) :: efclfrac(nlayers,ngptlw) 582 REAL(KIND=r8) :: uflux(0:nlayers) 583 REAL(KIND=r8) :: dflux(0:nlayers) 584 REAL(KIND=r8) :: urad(0:nlayers) 585 REAL(KIND=r8) :: drad(0:nlayers) 586 REAL(KIND=r8) :: uclfl(0:nlayers) 587 REAL(KIND=r8) :: dclfl(0:nlayers) 588 REAL(KIND=r8) :: odcld(nlayers,ngptlw) 589 REAL(KIND=r8) :: secdiff(nbndlw) ! secant of diffusivity angle 590 REAL(KIND=r8) :: a0(nbndlw) 591 REAL(KIND=r8) :: a1(nbndlw) 592 REAL(KIND=r8) :: a2(nbndlw) ! diffusivity angle adjustment coefficients 593 REAL(KIND=r8) :: wtdiff 594 REAL(KIND=r8) :: rec_6 595 REAL(KIND=r8) :: transcld 596 REAL(KIND=r8) :: radld 597 REAL(KIND=r8) :: radclrd 598 REAL(KIND=r8) :: plfrac 599 REAL(KIND=r8) :: blay 600 REAL(KIND=r8) :: dplankup 601 REAL(KIND=r8) :: dplankdn 602 REAL(KIND=r8) :: odepth 603 REAL(KIND=r8) :: odtot 604 REAL(KIND=r8) :: odepth_rec 605 REAL(KIND=r8) :: gassrc 606 REAL(KIND=r8) :: odtot_rec 607 REAL(KIND=r8) :: bbdtot 608 REAL(KIND=r8) :: bbd 609 REAL(KIND=r8) :: tblind 610 REAL(KIND=r8) :: tfactot 611 REAL(KIND=r8) :: tfacgas 612 REAL(KIND=r8) :: transc 613 REAL(KIND=r8) :: tausfac 614 REAL(KIND=r8) :: rad0 615 REAL(KIND=r8) :: reflect 616 REAL(KIND=r8) :: radlu 617 REAL(KIND=r8) :: radclru 618 INTEGER :: icldlyr(nlayers) ! flag for cloud in layer 619 INTEGER :: ibnd 620 INTEGER :: lay 621 INTEGER :: ig 622 INTEGER :: ib 623 INTEGER :: iband 624 INTEGER :: lev 625 INTEGER :: l ! loop indices 626 INTEGER :: igc ! g-point interval counter 627 INTEGER :: iclddn ! flag for cloud in down path 628 INTEGER :: ittot 629 INTEGER :: itgas 630 INTEGER :: itr ! lookup table indices 631 ! ------- Definitions ------- 632 ! input 633 ! nlayers ! number of model layers 634 ! ngptlw ! total number of g-point subintervals 635 ! nbndlw ! number of longwave spectral bands 636 ! ncbands ! number of spectral bands for clouds 637 ! secdiff ! diffusivity angle 638 ! wtdiff ! weight for radiance to flux conversion 639 ! pavel ! layer pressures (mb) 640 ! pz ! level (interface) pressures (mb) 641 ! tavel ! layer temperatures (k) 642 ! tz ! level (interface) temperatures(mb) 643 ! tbound ! surface temperature (k) 644 ! cldfrac ! layer cloud fraction 645 ! taucloud ! layer cloud optical depth 646 ! itr ! integer look-up table index 647 ! icldlyr ! flag for cloudy layers 648 ! iclddn ! flag for cloud in column at any layer 649 ! semiss ! surface emissivities for each band 650 ! reflect ! surface reflectance 651 ! bpade ! 1/(pade constant) 652 ! tau_tbl ! clear sky optical depth look-up table 653 ! exp_tbl ! exponential look-up table for transmittance 654 ! tfn_tbl ! tau transition function look-up table 655 ! local 656 ! atrans ! gaseous absorptivity 657 ! abscld ! cloud absorptivity 658 ! atot ! combined gaseous and cloud absorptivity 659 ! odclr ! clear sky (gaseous) optical depth 660 ! odcld ! cloud optical depth 661 ! odtot ! optical depth of gas and cloud 662 ! tfacgas ! gas-only pade factor, used for planck fn 663 ! tfactot ! gas and cloud pade factor, used for planck fn 664 ! bbdgas ! gas-only planck function for downward rt 665 ! bbugas ! gas-only planck function for upward rt 666 ! bbdtot ! gas and cloud planck function for downward rt 667 ! bbutot ! gas and cloud planck function for upward calc. 668 ! gassrc ! source radiance due to gas only 669 ! efclfrac ! effective cloud fraction 670 ! radlu ! spectrally summed upward radiance 671 ! radclru ! spectrally summed clear sky upward radiance 672 ! urad ! upward radiance by layer 673 ! clrurad ! clear sky upward radiance by layer 674 ! radld ! spectrally summed downward radiance 675 ! radclrd ! spectrally summed clear sky downward radiance 676 ! drad ! downward radiance by layer 677 ! clrdrad ! clear sky downward radiance by layer 678 ! output 679 ! totuflux ! upward longwave flux (w/m2) 680 ! totdflux ! downward longwave flux (w/m2) 681 ! fnet ! net longwave flux (w/m2) 682 ! htr ! longwave heating rate (k/day) 683 ! totuclfl ! clear sky upward longwave flux (w/m2) 684 ! totdclfl ! clear sky downward longwave flux (w/m2) 685 ! fnetc ! clear sky net longwave flux (w/m2) 686 ! htrc ! clear sky longwave heating rate (k/day) 687 ! This secant and weight corresponds to the standard diffusivity 688 ! angle. This initial value is redefined below for some bands. 689 data wtdiff /0.5_r8/ 690 data rec_6 /0.166667_r8/ 691 ! Reset diffusivity angle for Bands 2-3 and 5-9 to vary (between 1.50 692 ! and 1.80) as a function of total column water vapor. The function 693 ! has been defined to minimize flux and cooling rate errors in these bands 694 ! over a wide range of precipitable water values. 695 data a0 / 1.66_r8, 1.55_r8, 1.58_r8, 1.66_r8, & 696 1.54_r8, 1.454_r8, 1.89_r8, 1.33_r8, & 697 1.668_r8, 1.66_r8, 1.66_r8, 1.66_r8, & 698 1.66_r8, 1.66_r8, 1.66_r8, 1.66_r8 / 699 data a1 / 0.00_r8, 0.25_r8, 0.22_r8, 0.00_r8, & 700 0.13_r8, 0.446_r8, -0.10_r8, 0.40_r8, & 701 -0.006_r8, 0.00_r8, 0.00_r8, 0.00_r8, & 702 0.00_r8, 0.00_r8, 0.00_r8, 0.00_r8 / 703 data a2 / 0.00_r8, -12.0_r8, -11.7_r8, 0.00_r8, & 704 -0.72_r8,-0.243_r8, 0.19_r8,-0.062_r8, & 705 0.414_r8, 0.00_r8, 0.00_r8, 0.00_r8, & 706 0.00_r8, 0.00_r8, 0.00_r8, 0.00_r8 / 707 integer iplon 708 hvrrtc = '$Revision: 1.3 $' 709 do iplon=1,ncol 710 do ibnd = 1,nbndlw 711 if (ibnd.eq.1 .or. ibnd.eq.4 .or. ibnd.ge.10) then 712 secdiff(ibnd) = 1.66_r8 713 else 714 secdiff(ibnd) = a0(ibnd) + a1(ibnd)*exp(a2(ibnd)*pwvcm(iplon)) 715 endif 716 enddo 717 if (pwvcm(iplon).lt.1.0) secdiff(6) = 1.80_r8 718 if (pwvcm(iplon).gt.7.1) secdiff(7) = 1.50_r8 719 urad(0) = 0.0_r8 720 drad(0) = 0.0_r8 721 totuflux(iplon,0) = 0.0_r8 722 totdflux(iplon,0) = 0.0_r8 723 clrurad(0) = 0.0_r8 724 clrdrad(0) = 0.0_r8 725 totuclfl(iplon,0) = 0.0_r8 726 totdclfl(iplon,0) = 0.0_r8 727 do lay = 1, nlayers 728 urad(lay) = 0.0_r8 729 drad(lay) = 0.0_r8 730 totuflux(iplon,lay) = 0.0_r8 731 totdflux(iplon,lay) = 0.0_r8 732 clrurad(lay) = 0.0_r8 733 clrdrad(lay) = 0.0_r8 734 totuclfl(iplon,lay) = 0.0_r8 735 totdclfl(iplon,lay) = 0.0_r8 736 icldlyr(lay) = 0 737 ! Change to band loop? 738 do ig = 1, ngptlw 739 if (cldfmc(iplon,ig,lay) .eq. 1._r8) then 740 ib = ngb(ig) 741 odcld(lay,ig) = secdiff(ib) * taucmc(iplon,ig,lay) 742 transcld = exp(-odcld(lay,ig)) 743 abscld(lay,ig) = 1._r8 - transcld 744 efclfrac(lay,ig) = abscld(lay,ig) * cldfmc(iplon,ig,lay) 745 icldlyr(lay) = 1 746 else 747 odcld(lay,ig) = 0.0_r8 748 abscld(lay,ig) = 0.0_r8 749 efclfrac(lay,ig) = 0.0_r8 750 endif 751 enddo 752 enddo 753 igc = 1 754 ! Loop over frequency bands. 755 do iband = istart, iend 756 ! Reinitialize g-point counter for each band if output for each band is requested. 757 if (iout.gt.0.and.iband.ge.2) igc = ngs(iband-1)+1 758 ! Loop over g-channels. 759 1000 continue 760 ! Radiative transfer starts here. 761 radld = 0._r8 762 radclrd = 0._r8 763 iclddn = 0 764 ! Downward radiative transfer loop. 765 do lev = nlayers, 1, -1 766 plfrac = fracs(iplon,lev,igc) 767 blay = planklay(iplon,lev,iband) 768 dplankup = planklev(iplon,lev,iband) - blay 769 dplankdn = planklev(iplon,lev-1,iband) - blay 770 odepth = secdiff(iband) * taut(iplon,lev,igc) 771 if (odepth .lt. 0.0_r8) odepth = 0.0_r8 772 ! Cloudy layer 773 if (icldlyr(lev).eq.1) then 774 iclddn = 1 775 odtot = odepth + odcld(lev,igc) 776 if (odtot .lt. 0.06_r8) then 777 atrans(lev) = odepth - 0.5_r8*odepth*odepth 778 odepth_rec = rec_6*odepth 779 gassrc = plfrac*(blay+dplankdn*odepth_rec)*atrans(lev) 780 atot(lev) = odtot - 0.5_r8*odtot*odtot 781 odtot_rec = rec_6*odtot 782 bbdtot = plfrac * (blay+dplankdn*odtot_rec) 783 bbd = plfrac*(blay+dplankdn*odepth_rec) 784 radld = radld - radld * (atrans(lev) + & 785 efclfrac(lev,igc) * (1. - atrans(lev))) + & 786 gassrc + cldfmc(iplon,igc,lev) * & 787 (bbdtot * atot(lev) - gassrc) 788 drad( lev-1) = drad(lev-1) + radld 789 bbugas(lev) = plfrac * (blay+dplankup*odepth_rec) 790 bbutot(lev) = plfrac * (blay+dplankup*odtot_rec) 791 elseif (odepth .le. 0.06_r8) then 792 atrans(lev) = odepth - 0.5_r8*odepth*odepth 793 odepth_rec = rec_6*odepth 794 gassrc = plfrac*(blay+dplankdn*odepth_rec)*atrans(lev) 795 odtot = odepth + odcld(lev,igc) 796 tblind = odtot/(bpade+odtot) 797 ittot = tblint*tblind + 0.5_r8 798 tfactot = tfn_tbl(ittot) 799 bbdtot = plfrac * (blay + tfactot*dplankdn) 800 bbd = plfrac*(blay+dplankdn*odepth_rec) 801 atot(lev) = 1. - exp_tbl(ittot) 802 radld = radld - radld * (atrans(lev) + & 803 efclfrac(lev,igc) * (1._r8 - atrans(lev))) + & 804 gassrc + cldfmc(iplon,igc,lev) * & 805 (bbdtot * atot(lev) - gassrc) 806 drad(lev-1) = drad(lev-1) + radld 807 bbugas(lev) = plfrac * (blay + dplankup*odepth_rec) 808 bbutot(lev) = plfrac * (blay + tfactot * dplankup) 809 else 810 tblind = odepth/(bpade+odepth) 811 itgas = tblint*tblind+0.5_r8 812 odepth = tau_tbl(itgas) 813 atrans(lev) = 1._r8 - exp_tbl(itgas) 814 tfacgas = tfn_tbl(itgas) 815 gassrc = atrans(lev) * plfrac * (blay + tfacgas*dplankdn) 816 odtot = odepth + odcld(lev,igc) 817 tblind = odtot/(bpade+odtot) 818 ittot = tblint*tblind + 0.5_r8 819 tfactot = tfn_tbl(ittot) 820 bbdtot = plfrac * (blay + tfactot*dplankdn) 821 bbd = plfrac*(blay+tfacgas*dplankdn) 822 atot(lev) = 1._r8 - exp_tbl(ittot) 823 radld = radld - radld * (atrans(lev) + & 824 efclfrac(lev,igc) * (1._r8 - atrans(lev))) + & 825 gassrc + cldfmc(iplon,igc,lev) * & 826 (bbdtot * atot(lev) - gassrc) 827 drad(lev-1) = drad(lev-1) + radld 828 bbugas(lev) = plfrac * (blay + tfacgas * dplankup) 829 bbutot(lev) = plfrac * (blay + tfactot * dplankup) 830 endif 831 ! Clear layer 832 else 833 if (odepth .le. 0.06_r8) then 834 atrans(lev) = odepth-0.5_r8*odepth*odepth 835 odepth = rec_6*odepth 836 bbd = plfrac*(blay+dplankdn*odepth) 837 bbugas(lev) = plfrac*(blay+dplankup*odepth) 838 else 839 tblind = odepth/(bpade+odepth) 840 itr = tblint*tblind+0.5_r8 841 transc = exp_tbl(itr) 842 atrans(lev) = 1._r8-transc 843 tausfac = tfn_tbl(itr) 844 bbd = plfrac*(blay+tausfac*dplankdn) 845 bbugas(lev) = plfrac * (blay + tausfac * dplankup) 846 endif 847 radld = radld + (bbd-radld)*atrans(lev) 848 drad(lev-1) = drad(lev-1) + radld 849 endif 850 ! Set clear sky stream to total sky stream as long as layers 851 ! remain clear. Streams diverge when a cloud is reached (iclddn=1), 852 ! and clear sky stream must be computed separately from that point. 853 if (iclddn.eq.1) then 854 radclrd = radclrd + (bbd-radclrd) * atrans(lev) 855 clrdrad(lev-1) = clrdrad(lev-1) + radclrd 856 else 857 radclrd = radld 858 clrdrad(lev-1) = drad(lev-1) 859 endif 860 enddo 861 ! Spectral emissivity & reflectance 862 ! Include the contribution of spectrally varying longwave emissivity 863 ! and reflection from the surface to the upward radiative transfer. 864 ! Note: Spectral and Lambertian reflection are identical for the 865 ! diffusivity angle flux integration used here. 866 rad0 = fracs(iplon,1,igc) * plankbnd(iplon,iband) 867 ! Add in specular reflection of surface downward radiance. 868 reflect = 1._r8 - semiss(iplon,iband) 869 radlu = rad0 + reflect * radld 870 radclru = rad0 + reflect * radclrd 871 ! Upward radiative transfer loop. 872 urad(0) = urad(0) + radlu 873 clrurad(0) = clrurad(0) + radclru 874 do lev = 1, nlayers 875 ! Cloudy layer 876 if (icldlyr(lev) .eq. 1) then 877 gassrc = bbugas(lev) * atrans(lev) 878 radlu = radlu - radlu * (atrans(lev) + & 879 efclfrac(lev,igc) * (1._r8 - atrans(lev))) + & 880 gassrc + cldfmc(iplon,igc,lev) * & 881 (bbutot(lev) * atot(lev) - gassrc) 882 urad(lev) = urad(lev) + radlu 883 ! Clear layer 884 else 885 radlu = radlu + (bbugas(lev)-radlu)*atrans(lev) 886 urad(lev) = urad(lev) + radlu 887 endif 888 ! Set clear sky stream to total sky stream as long as all layers 889 ! are clear (iclddn=0). Streams must be calculated separately at 890 ! all layers when a cloud is present (ICLDDN=1), because surface 891 ! reflectance is different for each stream. 892 if (iclddn.eq.1) then 893 radclru = radclru + (bbugas(lev)-radclru)*atrans(lev) 894 clrurad(lev) = clrurad(lev) + radclru 895 else 896 radclru = radlu 897 clrurad(lev) = urad(lev) 898 endif 899 enddo 900 ! Increment g-point counter 901 igc = igc + 1 902 ! Return to continue radiative transfer for all g-channels in present band 903 if (igc .le. ngs(iband)) go to 1000 904 ! Process longwave output from band for total and clear streams. 905 ! Calculate upward, downward, and net flux. 906 do lev = nlayers, 0, -1 907 uflux(lev) = urad(lev)*wtdiff 908 dflux(lev) = drad(lev)*wtdiff 909 urad(lev) = 0.0_r8 910 drad(lev) = 0.0_r8 911 totuflux(iplon,lev) = totuflux(iplon,lev) + uflux(lev) * delwave(iband) 912 totdflux(iplon,lev) = totdflux(iplon,lev) + dflux(lev) * delwave(iband) 913 uclfl(lev) = clrurad(lev)*wtdiff 914 dclfl(lev) = clrdrad(lev)*wtdiff 915 clrurad(lev) = 0.0_r8 916 clrdrad(lev) = 0.0_r8 917 totuclfl(iplon,lev) = totuclfl(iplon,lev) + uclfl(lev) * delwave(iband) 918 totdclfl(iplon,lev) = totdclfl(iplon,lev) + dclfl(lev) * delwave(iband) 919 totufluxs(iplon,iband,lev) = uflux(lev) * delwave(iband) 920 totdfluxs(iplon,iband,lev) = dflux(lev) * delwave(iband) 921 enddo 922 ! End spectral band loop 923 enddo 924 enddo 925 do iplon=1,ncol 926 ! Calculate fluxes at surface 927 totuflux(iplon,0) = totuflux(iplon,0) * fluxfac 928 totdflux(iplon,0) = totdflux(iplon,0) * fluxfac 929 totufluxs(iplon,:,0) = totufluxs(iplon,:,0) * fluxfac 930 totdfluxs(iplon,:,0) = totdfluxs(iplon,:,0) * fluxfac 931 fnet(iplon,0) = totuflux(iplon,0) - totdflux(iplon,0) 932 totuclfl(iplon,0) = totuclfl(iplon,0) * fluxfac 933 totdclfl(iplon,0) = totdclfl(iplon,0) * fluxfac 934 fnetc(iplon,0) = totuclfl(iplon,0) - totdclfl(iplon,0) 935 enddo 936 ! Calculate fluxes at model levels 937 do lev = 1, nlayers 938 do iplon=1,ncol 939 totuflux(iplon,lev) = totuflux(iplon,lev) * fluxfac 940 totdflux(iplon,lev) = totdflux(iplon,lev) * fluxfac 941 totufluxs(iplon,:,lev) = totufluxs(iplon,:,lev) * fluxfac 942 totdfluxs(iplon,:,lev) = totdfluxs(iplon,:,lev) * fluxfac 943 fnet(iplon,lev) = totuflux(iplon,lev) - totdflux(iplon,lev) 944 totuclfl(iplon,lev) = totuclfl(iplon,lev) * fluxfac 945 totdclfl(iplon,lev) = totdclfl(iplon,lev) * fluxfac 946 fnetc(iplon,lev) = totuclfl(iplon,lev) - totdclfl(iplon,lev) 947 l = lev - 1 948 ! Calculate heating rates at model layers 949 htr(iplon,l)=heatfac*(fnet(iplon,l)-fnet(iplon,lev))/(pz(iplon,l)-pz(iplon,lev)) 950 htrc(iplon,l)=heatfac*(fnetc(iplon,l)-fnetc(iplon,lev))/(pz(iplon,l)-pz(iplon,lev)) 951 enddo 952 enddo 953 ! Set heating rate to zero in top layer 954 do iplon=1,ncol 955 htr(iplon,nlayers) = 0.0_r8 956 htrc(iplon,nlayers) = 0.0_r8 957 enddo 958 END SUBROUTINE rtrnmc 959#endif 960 961 END MODULE rrtmg_lw_rtrnmc 962