Line data Source code
1 : ! path: $Source: /storm/rc1/cvsroot/rc/rrtmg_lw/src/rrtmg_lw_taumol.f90,v $
2 : ! author: $Author: mike $
3 : ! revision: $Revision: 1.7 $
4 : ! created: $Date: 2009/10/20 15:08:37 $
5 : !
6 : module rrtmg_lw_taumol
7 :
8 : ! --------------------------------------------------------------------------
9 : ! | |
10 : ! | Copyright 2002-2009, Atmospheric & Environmental Research, Inc. (AER). |
11 : ! | This software may be used, copied, or redistributed as long as it is |
12 : ! | not sold and this copyright notice is reproduced on each copy made. |
13 : ! | This model is provided as is without any express or implied warranties. |
14 : ! | (http://www.rtweb.aer.com/) |
15 : ! | |
16 : ! --------------------------------------------------------------------------
17 :
18 : ! ------- Modules -------
19 :
20 : use shr_kind_mod, only: r8 => shr_kind_r8
21 :
22 : use parrrtm, only: nbndlw, maxxsec, ngptlw
23 : use rrlw_con, only: oneminus
24 : use rrlw_wvn, only: nspa, nspb
25 :
26 : implicit none
27 :
28 : contains
29 :
30 : !----------------------------------------------------------------------------
31 486000 : subroutine taumol(nlayers, pavel, wx, coldry, &
32 486000 : laytrop, jp, jt, jt1, planklay, planklev, plankbnd, &
33 486000 : colh2o, colco2, colo3, coln2o, colco, colch4, colo2, &
34 486000 : colbrd, fac00, fac01, fac10, fac11, &
35 486000 : rat_h2oco2, rat_h2oco2_1, rat_h2oo3, rat_h2oo3_1, &
36 486000 : rat_h2on2o, rat_h2on2o_1, rat_h2och4, rat_h2och4_1, &
37 486000 : rat_n2oco2, rat_n2oco2_1, rat_o3co2, rat_o3co2_1, &
38 486000 : selffac, selffrac, indself, forfac, forfrac, indfor, &
39 486000 : minorfrac, scaleminor, scaleminorn2, indminor, &
40 486000 : fracs, taug)
41 : !----------------------------------------------------------------------------
42 :
43 : ! *******************************************************************************
44 : ! * *
45 : ! * Optical depths developed for the *
46 : ! * *
47 : ! * RAPID RADIATIVE TRANSFER MODEL (RRTM) *
48 : ! * *
49 : ! * *
50 : ! * ATMOSPHERIC AND ENVIRONMENTAL RESEARCH, INC. *
51 : ! * 131 HARTWELL AVENUE *
52 : ! * LEXINGTON, MA 02421 *
53 : ! * *
54 : ! * *
55 : ! * ELI J. MLAWER *
56 : ! * JENNIFER DELAMERE *
57 : ! * STEVEN J. TAUBMAN *
58 : ! * SHEPARD A. CLOUGH *
59 : ! * *
60 : ! * *
61 : ! * *
62 : ! * *
63 : ! * email: mlawer@aer.com *
64 : ! * email: jdelamer@aer.com *
65 : ! * *
66 : ! * The authors wish to acknowledge the contributions of the *
67 : ! * following people: Karen Cady-Pereira, Patrick D. Brown, *
68 : ! * Michael J. Iacono, Ronald E. Farren, Luke Chen, Robert Bergstrom. *
69 : ! * *
70 : ! *******************************************************************************
71 : ! * *
72 : ! * Revision for g-point reduction: Michael J. Iacono, AER, Inc. *
73 : ! * *
74 : ! *******************************************************************************
75 : ! * TAUMOL *
76 : ! * *
77 : ! * This file contains the subroutines TAUGBn (where n goes from *
78 : ! * 1 to 16). TAUGBn calculates the optical depths and Planck fractions *
79 : ! * per g-value and layer for band n. *
80 : ! * *
81 : ! * Output: optical depths (unitless) *
82 : ! * fractions needed to compute Planck functions at every layer *
83 : ! * and g-value *
84 : ! * *
85 : ! * COMMON /TAUGCOM/ TAUG(MXLAY,MG) *
86 : ! * COMMON /PLANKG/ FRACS(MXLAY,MG) *
87 : ! * *
88 : ! * Input *
89 : ! * *
90 : ! * COMMON /FEATURES/ NG(NBANDS),NSPA(NBANDS),NSPB(NBANDS) *
91 : ! * COMMON /PRECISE/ ONEMINUS *
92 : ! * COMMON /PROFILE/ NLAYERS,PAVEL(MXLAY),TAVEL(MXLAY), *
93 : ! * & PZ(0:MXLAY),TZ(0:MXLAY) *
94 : ! * COMMON /PROFDATA/ LAYTROP, *
95 : ! * & COLH2O(MXLAY),COLCO2(MXLAY),COLO3(MXLAY), *
96 : ! * & COLN2O(MXLAY),COLCO(MXLAY),COLCH4(MXLAY), *
97 : ! * & COLO2(MXLAY)
98 : ! * COMMON /INTFAC/ FAC00(MXLAY),FAC01(MXLAY), *
99 : ! * & FAC10(MXLAY),FAC11(MXLAY) *
100 : ! * COMMON /INTIND/ JP(MXLAY),JT(MXLAY),JT1(MXLAY) *
101 : ! * COMMON /SELF/ SELFFAC(MXLAY), SELFFRAC(MXLAY), INDSELF(MXLAY) *
102 : ! * *
103 : ! * Description: *
104 : ! * NG(IBAND) - number of g-values in band IBAND *
105 : ! * NSPA(IBAND) - for the lower atmosphere, the number of reference *
106 : ! * atmospheres that are stored for band IBAND per *
107 : ! * pressure level and temperature. Each of these *
108 : ! * atmospheres has different relative amounts of the *
109 : ! * key species for the band (i.e. different binary *
110 : ! * species parameters). *
111 : ! * NSPB(IBAND) - same for upper atmosphere *
112 : ! * ONEMINUS - since problems are caused in some cases by interpolation *
113 : ! * parameters equal to or greater than 1, for these cases *
114 : ! * these parameters are set to this value, slightly < 1. *
115 : ! * PAVEL - layer pressures (mb) *
116 : ! * TAVEL - layer temperatures (degrees K) *
117 : ! * PZ - level pressures (mb) *
118 : ! * TZ - level temperatures (degrees K) *
119 : ! * LAYTROP - layer at which switch is made from one combination of *
120 : ! * key species to another *
121 : ! * COLH2O, COLCO2, COLO3, COLN2O, COLCH4 - column amounts of water *
122 : ! * vapor,carbon dioxide, ozone, nitrous ozide, methane, *
123 : ! * respectively (molecules/cm**2) *
124 : ! * FACij(LAY) - for layer LAY, these are factors that are needed to *
125 : ! * compute the interpolation factors that multiply the *
126 : ! * appropriate reference k-values. A value of 0 (1) for *
127 : ! * i,j indicates that the corresponding factor multiplies *
128 : ! * reference k-value for the lower (higher) of the two *
129 : ! * appropriate temperatures, and altitudes, respectively. *
130 : ! * JP - the index of the lower (in altitude) of the two appropriate *
131 : ! * reference pressure levels needed for interpolation *
132 : ! * JT, JT1 - the indices of the lower of the two appropriate reference *
133 : ! * temperatures needed for interpolation (for pressure *
134 : ! * levels JP and JP+1, respectively) *
135 : ! * SELFFAC - scale factor needed for water vapor self-continuum, equals *
136 : ! * (water vapor density)/(atmospheric density at 296K and *
137 : ! * 1013 mb) *
138 : ! * SELFFRAC - factor needed for temperature interpolation of reference *
139 : ! * water vapor self-continuum data *
140 : ! * INDSELF - index of the lower of the two appropriate reference *
141 : ! * temperatures needed for the self-continuum interpolation *
142 : ! * FORFAC - scale factor needed for water vapor foreign-continuum. *
143 : ! * FORFRAC - factor needed for temperature interpolation of reference *
144 : ! * water vapor foreign-continuum data *
145 : ! * INDFOR - index of the lower of the two appropriate reference *
146 : ! * temperatures needed for the foreign-continuum interpolation *
147 : ! * *
148 : ! * Data input *
149 : ! * COMMON /Kn/ KA(NSPA(n),5,13,MG), KB(NSPB(n),5,13:59,MG), SELFREF(10,MG),*
150 : ! * FORREF(4,MG), KA_M'MGAS', KB_M'MGAS' *
151 : ! * (note: n is the band number,'MGAS' is the species name of the minor *
152 : ! * gas) *
153 : ! * *
154 : ! * Description: *
155 : ! * KA - k-values for low reference atmospheres (key-species only) *
156 : ! * (units: cm**2/molecule) *
157 : ! * KB - k-values for high reference atmospheres (key-species only) *
158 : ! * (units: cm**2/molecule) *
159 : ! * KA_M'MGAS' - k-values for low reference atmosphere minor species *
160 : ! * (units: cm**2/molecule) *
161 : ! * KB_M'MGAS' - k-values for high reference atmosphere minor species *
162 : ! * (units: cm**2/molecule) *
163 : ! * SELFREF - k-values for water vapor self-continuum for reference *
164 : ! * atmospheres (used below LAYTROP) *
165 : ! * (units: cm**2/molecule) *
166 : ! * FORREF - k-values for water vapor foreign-continuum for reference *
167 : ! * atmospheres (used below/above LAYTROP) *
168 : ! * (units: cm**2/molecule) *
169 : ! * *
170 : ! * DIMENSION ABSA(65*NSPA(n),MG), ABSB(235*NSPB(n),MG) *
171 : ! * EQUIVALENCE (KA,ABSA),(KB,ABSB) *
172 : ! * *
173 : !*******************************************************************************
174 :
175 : ! ------- Declarations -------
176 :
177 : ! ----- Input -----
178 : integer, intent(in) :: nlayers ! total number of layers
179 : real(kind=r8), intent(in) :: pavel(nlayers) ! layer pressures (mb)
180 : ! Dimensions: (nlayers)
181 : real(kind=r8), intent(in) :: wx(maxxsec,nlayers) ! cross-section amounts (mol/cm2)
182 : ! Dimensions: (maxxsec,nlayers)
183 : real(kind=r8), intent(in) :: coldry(nlayers) ! column amount (dry air)
184 : ! Dimensions: (nlayers)
185 :
186 : integer, intent(in) :: laytrop ! tropopause layer index
187 : integer, intent(in) :: jp(nlayers) !
188 : ! Dimensions: (nlayers)
189 : integer, intent(in) :: jt(nlayers) !
190 : ! Dimensions: (nlayers)
191 : integer, intent(in) :: jt1(nlayers) !
192 : ! Dimensions: (nlayers)
193 : real(kind=r8), intent(in) :: planklay(nlayers,nbndlw) !
194 : ! Dimensions: (nlayers,nbndlw)
195 : real(kind=r8), intent(in) :: planklev(0:nlayers,nbndlw) !
196 : ! Dimensions: (nlayers,nbndlw)
197 : real(kind=r8), intent(in) :: plankbnd(nbndlw) !
198 : ! Dimensions: (nbndlw)
199 :
200 : real(kind=r8), intent(in) :: colh2o(nlayers) ! column amount (h2o)
201 : ! Dimensions: (nlayers)
202 : real(kind=r8), intent(in) :: colco2(nlayers) ! column amount (co2)
203 : ! Dimensions: (nlayers)
204 : real(kind=r8), intent(in) :: colo3(nlayers) ! column amount (o3)
205 : ! Dimensions: (nlayers)
206 : real(kind=r8), intent(in) :: coln2o(nlayers) ! column amount (n2o)
207 : ! Dimensions: (nlayers)
208 : real(kind=r8), intent(in) :: colco(nlayers) ! column amount (co)
209 : ! Dimensions: (nlayers)
210 : real(kind=r8), intent(in) :: colch4(nlayers) ! column amount (ch4)
211 : ! Dimensions: (nlayers)
212 : real(kind=r8), intent(in) :: colo2(nlayers) ! column amount (o2)
213 : ! Dimensions: (nlayers)
214 : real(kind=r8), intent(in) :: colbrd(nlayers) ! column amount (broadening gases)
215 : ! Dimensions: (nlayers)
216 :
217 : integer, intent(in) :: indself(nlayers)
218 : ! Dimensions: (nlayers)
219 : integer, intent(in) :: indfor(nlayers)
220 : ! Dimensions: (nlayers)
221 : real(kind=r8), intent(in) :: selffac(nlayers)
222 : ! Dimensions: (nlayers)
223 : real(kind=r8), intent(in) :: selffrac(nlayers)
224 : ! Dimensions: (nlayers)
225 : real(kind=r8), intent(in) :: forfac(nlayers)
226 : ! Dimensions: (nlayers)
227 : real(kind=r8), intent(in) :: forfrac(nlayers)
228 : ! Dimensions: (nlayers)
229 :
230 : integer, intent(in) :: indminor(nlayers)
231 : ! Dimensions: (nlayers)
232 : real(kind=r8), intent(in) :: minorfrac(nlayers)
233 : ! Dimensions: (nlayers)
234 : real(kind=r8), intent(in) :: scaleminor(nlayers)
235 : ! Dimensions: (nlayers)
236 : real(kind=r8), intent(in) :: scaleminorn2(nlayers)
237 : ! Dimensions: (nlayers)
238 :
239 : real(kind=r8), dimension(nlayers), intent(in) :: & !
240 : fac00, fac01, & ! Dimensions: (nlayers)
241 : fac10, fac11
242 : real(kind=r8), dimension(nlayers), intent(in) :: & !
243 : rat_h2oco2,rat_h2oco2_1, &
244 : rat_h2oo3,rat_h2oo3_1, & ! Dimensions: (nlayers)
245 : rat_h2on2o,rat_h2on2o_1, &
246 : rat_h2och4,rat_h2och4_1, &
247 : rat_n2oco2,rat_n2oco2_1, &
248 : rat_o3co2,rat_o3co2_1
249 :
250 : ! ----- Output -----
251 : real(kind=r8), intent(out) :: fracs(nlayers,ngptlw) ! planck fractions
252 : ! Dimensions: (nlayers,ngptlw)
253 : real(kind=r8), intent(out) :: taug(nlayers,ngptlw) ! gaseous optical depth
254 : ! Dimensions: (nlayers,ngptlw)
255 :
256 : ! Calculate gaseous optical depth and planck fractions for each spectral band.
257 :
258 486000 : call taugb1
259 486000 : call taugb2
260 486000 : call taugb3
261 486000 : call taugb4
262 486000 : call taugb5
263 486000 : call taugb6
264 486000 : call taugb7
265 486000 : call taugb8
266 486000 : call taugb9
267 486000 : call taugb10
268 486000 : call taugb11
269 486000 : call taugb12
270 486000 : call taugb13
271 486000 : call taugb14
272 486000 : call taugb15
273 486000 : call taugb16
274 :
275 : contains
276 :
277 : !----------------------------------------------------------------------------
278 486000 : subroutine taugb1
279 : !----------------------------------------------------------------------------
280 :
281 : ! ------- Modifications -------
282 : ! Written by Eli J. Mlawer, Atmospheric & Environmental Research.
283 : ! Revised by Michael J. Iacono, Atmospheric & Environmental Research.
284 : !
285 : ! band 1: 10-350 cm-1 (low key - h2o; low minor - n2)
286 : ! (high key - h2o; high minor - n2)
287 : !
288 : ! note: previous versions of rrtm band 1:
289 : ! 10-250 cm-1 (low - h2o; high - h2o)
290 : !----------------------------------------------------------------------------
291 :
292 : ! ------- Modules -------
293 :
294 : use parrrtm, only : ng1
295 : use rrlw_kg01, only : fracrefa, fracrefb, absa, absb, &
296 : ka_mn2, kb_mn2, selfref, forref
297 :
298 : ! ------- Declarations -------
299 :
300 : ! Local
301 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
302 972000 : real(kind=r8) :: corradj(nlayers), scalen2(nlayers), tauself, taufor, taun2
303 :
304 :
305 : ! Minor gas mapping levels:
306 : ! lower - n2, p = 142.5490 mbar, t = 215.70 k
307 : ! upper - n2, p = 142.5490 mbar, t = 215.70 k
308 :
309 : ! Compute the optical depth by interpolating in ln(pressure) and
310 : ! temperature. Below laytrop, the water vapor self-continuum and
311 : ! foreign continuum is interpolated (in temperature) separately.
312 :
313 : ! Lower atmosphere loop
314 9219976 : do lay = 1, laytrop
315 :
316 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(1) + 1
317 8733976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(1) + 1
318 : !inds = indself(lay)
319 : !indf = indfor(lay)
320 : !indm = indminor(lay)
321 : !pp = pavel(lay)
322 8733976 : corradj(lay) = 1.
323 9219976 : scalen2(lay) = colbrd(lay) * scaleminorn2(lay)
324 : enddo
325 9219976 : do lay = 1, laytrop
326 9219976 : if (pavel(lay) .lt. 250._r8) then
327 2960573 : corradj(lay) = 1._r8 - 0.15_r8 * (250._r8-pavel(lay)) / 154.4_r8
328 : endif
329 : enddo
330 9219976 : do lay = 1, laytrop
331 96559736 : do ig = 1, ng1
332 436698800 : tauself = selffac(lay) * (selfref(indself(lay),ig) + selffrac(lay) * &
333 524038560 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
334 436698800 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
335 524038560 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
336 349359040 : taun2 = scalen2(lay)*(ka_mn2(indminor(lay),ig) + &
337 436698800 : minorfrac(lay) * (ka_mn2(indminor(lay)+1,ig) - ka_mn2(indminor(lay),ig)))
338 349359040 : taug(lay,ig) = corradj(lay) * (colh2o(lay) * &
339 349359040 : (fac00(lay) * absa(ind0(lay),ig) + &
340 349359040 : fac10(lay) * absa(ind0(lay)+1,ig) + &
341 349359040 : fac01(lay) * absa(ind1(lay),ig) + &
342 349359040 : fac11(lay) * absa(ind1(lay)+1,ig)) &
343 1746795200 : + tauself + taufor + taun2)
344 96073736 : fracs(lay,ig) = fracrefa(ig)
345 : enddo
346 : enddo
347 :
348 : ! Upper atmosphere loop
349 4874024 : do lay = laytrop+1, nlayers
350 :
351 4388024 : ind0(lay) = ((jp(lay)-13)*5+(jt(lay)-1))*nspb(1) + 1
352 4388024 : ind1(lay) = ((jp(lay)-12)*5+(jt1(lay)-1))*nspb(1) + 1
353 : !indf = indfor(lay)
354 : !indm = indminor(lay)
355 : !pp = pavel(lay)
356 4388024 : corradj(lay) = 1._r8 - 0.15_r8 * (pavel(lay) / 95.6_r8)
357 :
358 4874024 : scalen2(lay) = colbrd(lay) * scaleminorn2(lay)
359 : enddo
360 4874024 : do lay = laytrop+1, nlayers
361 48754264 : do ig = 1, ng1
362 175520960 : taufor = forfac(lay) * (forref(indfor(lay),ig) + &
363 219401200 : forfrac(lay) * (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
364 175520960 : taun2 = scalen2(lay)*(kb_mn2(indminor(lay),ig) + &
365 219401200 : minorfrac(lay) * (kb_mn2(indminor(lay)+1,ig) - kb_mn2(indminor(lay),ig)))
366 175520960 : taug(lay,ig) = corradj(lay) * (colh2o(lay) * &
367 175520960 : (fac00(lay) * absb(ind0(lay),ig) + &
368 175520960 : fac10(lay) * absb(ind0(lay)+1,ig) + &
369 175520960 : fac01(lay) * absb(ind1(lay),ig) + &
370 175520960 : fac11(lay) * absb(ind1(lay)+1,ig)) &
371 877604800 : + taufor + taun2)
372 48268264 : fracs(lay,ig) = fracrefb(ig)
373 : enddo
374 : enddo
375 :
376 486000 : end subroutine taugb1
377 :
378 : !----------------------------------------------------------------------------
379 486000 : subroutine taugb2
380 : !----------------------------------------------------------------------------
381 : !
382 : ! band 2: 350-500 cm-1 (low key - h2o; high key - h2o)
383 : !
384 : ! note: previous version of rrtm band 2:
385 : ! 250 - 500 cm-1 (low - h2o; high - h2o)
386 : !----------------------------------------------------------------------------
387 :
388 : ! ------- Modules -------
389 :
390 : use parrrtm, only : ng2, ngs1
391 : use rrlw_kg02, only : fracrefa, fracrefb, absa, absb, &
392 : selfref, forref
393 :
394 : ! ------- Declarations -------
395 :
396 : ! Local
397 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
398 972000 : real(kind=r8) :: corradj(nlayers), tauself, taufor
399 :
400 :
401 : ! Compute the optical depth by interpolating in ln(pressure) and
402 : ! temperature. Below laytrop, the water vapor self-continuum and
403 : ! foreign continuum is interpolated (in temperature) separately.
404 :
405 : ! Lower atmosphere loop
406 9219976 : do lay = 1, laytrop
407 :
408 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(2) + 1
409 8733976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(2) + 1
410 : !inds = indself(lay)
411 : !indf = indfor(lay)
412 : !pp = pavel(lay)
413 9219976 : corradj(lay) = 1._r8 - .05_r8 * (pavel(lay) - 100._r8) / 900._r8
414 : enddo
415 9219976 : do lay = 1, laytrop
416 114027688 : do ig = 1, ng2
417 524038560 : tauself = selffac(lay) * (selfref(indself(lay),ig) + selffrac(lay) * &
418 628846272 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
419 524038560 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
420 628846272 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
421 419230848 : taug(lay,ngs1+ig) = corradj(lay) * (colh2o(lay) * &
422 419230848 : (fac00(lay) * absa(ind0(lay),ig) + &
423 419230848 : fac10(lay) * absa(ind0(lay)+1,ig) + &
424 419230848 : fac01(lay) * absa(ind1(lay),ig) + &
425 419230848 : fac11(lay) * absa(ind1(lay)+1,ig)) &
426 2096154240 : + tauself + taufor)
427 113541688 : fracs(lay,ngs1+ig) = fracrefa(ig)
428 : enddo
429 : enddo
430 :
431 : ! Upper atmosphere loop
432 4874024 : do lay = laytrop+1, nlayers
433 :
434 4388024 : ind0(lay) = ((jp(lay)-13)*5+(jt(lay)-1))*nspb(2) + 1
435 4874024 : ind1(lay) = ((jp(lay)-12)*5+(jt1(lay)-1))*nspb(2) + 1
436 : !indf = indfor(lay)
437 : enddo
438 4874024 : do lay = laytrop+1, nlayers
439 57530312 : do ig = 1, ng2
440 210625152 : taufor = forfac(lay) * (forref(indfor(lay),ig) + &
441 263281440 : forfrac(lay) * (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
442 157968864 : taug(lay,ngs1+ig) = colh2o(lay) * &
443 210625152 : (fac00(lay) * absb(ind0(lay),ig) + &
444 210625152 : fac10(lay) * absb(ind0(lay)+1,ig) + &
445 210625152 : fac01(lay) * absb(ind1(lay),ig) + &
446 210625152 : fac11(lay) * absb(ind1(lay)+1,ig)) &
447 1000469472 : + taufor
448 57044312 : fracs(lay,ngs1+ig) = fracrefb(ig)
449 : enddo
450 : enddo
451 :
452 486000 : end subroutine taugb2
453 :
454 : !----------------------------------------------------------------------------
455 486000 : subroutine taugb3
456 : !----------------------------------------------------------------------------
457 : !
458 : ! band 3: 500-630 cm-1 (low key - h2o,co2; low minor - n2o)
459 : ! (high key - h2o,co2; high minor - n2o)
460 : !----------------------------------------------------------------------------
461 :
462 : ! ------- Modules -------
463 :
464 : use parrrtm, only : ng3, ngs2
465 : use rrlw_ref, only : chi_mls
466 : use rrlw_kg03, only : fracrefa, fracrefb, absa, absb,&
467 : ka_mn2o, kb_mn2o, selfref, forref
468 :
469 : ! ------- Declarations -------
470 :
471 : ! Local
472 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
473 972000 : integer, dimension(nlayers) :: js, js1, jmn2o, jpl
474 972000 : real(kind=r8), dimension(nlayers) :: speccomb, specparm, specmult, fs
475 972000 : real(kind=r8), dimension(nlayers) :: speccomb1, specparm1, specmult1, fs1
476 972000 : real(kind=r8), dimension(nlayers) :: speccomb_mn2o, specparm_mn2o, specmult_mn2o, &
477 972000 : fmn2o, fmn2omf, chi_n2o, ratn2o, adjfac, adjcoln2o
478 972000 : real(kind=r8), dimension(nlayers) :: speccomb_planck, specparm_planck, specmult_planck, fpl
479 : real(kind=r8) :: p, p4, fk0, fk1, fk2
480 972000 : real(kind=r8), dimension(nlayers) :: fac000, fac100, fac200, fac010, fac110, fac210
481 972000 : real(kind=r8), dimension(nlayers) :: fac001, fac101, fac201, fac011, fac111, fac211
482 : real(kind=r8) :: tauself, taufor, n2om1, n2om2, absn2o
483 : real(kind=r8) :: refrat_planck_a, refrat_planck_b, refrat_m_a, refrat_m_b
484 : real(kind=r8) :: tau_major, tau_major1
485 :
486 :
487 : ! Minor gas mapping levels:
488 : ! lower - n2o, p = 706.272 mbar, t = 278.94 k
489 : ! upper - n2o, p = 95.58 mbar, t = 215.7 k
490 :
491 : ! P = 212.725 mb
492 486000 : refrat_planck_a = chi_mls(1,9)/chi_mls(2,9)
493 :
494 : ! P = 95.58 mb
495 486000 : refrat_planck_b = chi_mls(1,13)/chi_mls(2,13)
496 :
497 : ! P = 706.270mb
498 486000 : refrat_m_a = chi_mls(1,3)/chi_mls(2,3)
499 :
500 : ! P = 95.58 mb
501 486000 : refrat_m_b = chi_mls(1,13)/chi_mls(2,13)
502 :
503 : ! Compute the optical depth by interpolating in ln(pressure) and
504 : ! temperature, and appropriate species. Below laytrop, the water vapor
505 : ! self-continuum and foreign continuum is interpolated (in temperature)
506 : ! separately.
507 :
508 : ! Lower atmosphere loop
509 9219976 : do lay = 1, laytrop
510 :
511 8733976 : speccomb(lay) = colh2o(lay) + rat_h2oco2(lay)*colco2(lay)
512 8733976 : specparm(lay) = colh2o(lay)/speccomb(lay)
513 8733976 : if (specparm(lay) .ge. oneminus) specparm(lay) = oneminus
514 8733976 : specmult(lay) = 8._r8*(specparm(lay))
515 8733976 : js(lay) = 1 + int(specmult(lay))
516 8733976 : fs(lay) = mod(specmult(lay),1.0_r8)
517 :
518 8733976 : speccomb1(lay) = colh2o(lay) + rat_h2oco2_1(lay)*colco2(lay)
519 8733976 : specparm1(lay) = colh2o(lay)/speccomb1(lay)
520 8733976 : if (specparm1(lay) .ge. oneminus) specparm1(lay) = oneminus
521 8733976 : specmult1(lay) = 8._r8*(specparm1(lay))
522 8733976 : js1(lay) = 1 + int(specmult1(lay))
523 8733976 : fs1(lay) = mod(specmult1(lay),1.0_r8)
524 :
525 8733976 : speccomb_mn2o(lay) = colh2o(lay) + refrat_m_a*colco2(lay)
526 8733976 : specparm_mn2o(lay) = colh2o(lay)/speccomb_mn2o(lay)
527 8733976 : if (specparm_mn2o(lay) .ge. oneminus) specparm_mn2o(lay) = oneminus
528 8733976 : specmult_mn2o(lay) = 8._r8*specparm_mn2o(lay)
529 8733976 : jmn2o(lay) = 1 + int(specmult_mn2o(lay))
530 8733976 : fmn2o(lay) = mod(specmult_mn2o(lay),1.0_r8)
531 8733976 : fmn2omf(lay) = minorfrac(lay)*fmn2o(lay)
532 : ! In atmospheres where the amount of N2O is too great to be considered
533 : ! a minor species, adjust the column amount of N2O by an empirical factor
534 : ! to obtain the proper contribution.
535 8733976 : chi_n2o(lay) = coln2o(lay)/coldry(lay)
536 9219976 : ratn2o(lay) = 1.e20_r8*chi_n2o(lay)/chi_mls(4,jp(lay)+1)
537 : enddo
538 9219976 : do lay = 1, laytrop
539 9219976 : if (ratn2o(lay) .gt. 1.5_r8) then
540 0 : adjfac(lay) = 0.5_r8+(ratn2o(lay)-0.5_r8)**0.65_r8
541 0 : adjcoln2o(lay) = adjfac(lay)*chi_mls(4,jp(lay)+1)*coldry(lay)*1.e-20_r8
542 : else
543 8733976 : adjcoln2o(lay) = coln2o(lay)
544 : endif
545 : enddo
546 9219976 : do lay = 1, laytrop
547 :
548 8733976 : speccomb_planck(lay) = colh2o(lay)+refrat_planck_a*colco2(lay)
549 8733976 : specparm_planck(lay) = colh2o(lay)/speccomb_planck(lay)
550 8733976 : if (specparm_planck(lay) .ge. oneminus) specparm_planck(lay)=oneminus
551 8733976 : specmult_planck(lay) = 8._r8*specparm_planck(lay)
552 8733976 : jpl(lay)= 1 + int(specmult_planck(lay))
553 8733976 : fpl(lay) = mod(specmult_planck(lay),1.0_r8)
554 :
555 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(3) + js(lay)
556 9219976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(3) + js1(lay)
557 : enddo
558 9219976 : do lay = 1, laytrop
559 :
560 8733976 : if (specparm(lay) .lt. 0.125_r8) then
561 770299 : p = fs(lay) - 1
562 770299 : p4 = p**4
563 770299 : fk0 = p4
564 770299 : fk1 = 1 - p - 2.0_r8*p4
565 770299 : fk2 = p + p4
566 770299 : fac000(lay) = fk0*fac00(lay)
567 770299 : fac100(lay) = fk1*fac00(lay)
568 770299 : fac200(lay) = fk2*fac00(lay)
569 770299 : fac010(lay) = fk0*fac10(lay)
570 770299 : fac110(lay) = fk1*fac10(lay)
571 770299 : fac210(lay) = fk2*fac10(lay)
572 7963677 : else if (specparm(lay) .gt. 0.875_r8) then
573 17165 : p = -fs(lay)
574 17165 : p4 = p**4
575 17165 : fk0 = p4
576 17165 : fk1 = 1 - p - 2.0_r8*p4
577 17165 : fk2 = p + p4
578 17165 : fac000(lay) = fk0*fac00(lay)
579 17165 : fac100(lay) = fk1*fac00(lay)
580 17165 : fac200(lay) = fk2*fac00(lay)
581 17165 : fac010(lay) = fk0*fac10(lay)
582 17165 : fac110(lay) = fk1*fac10(lay)
583 17165 : fac210(lay) = fk2*fac10(lay)
584 : else
585 7946512 : fac000(lay) = (1._r8 - fs(lay)) * fac00(lay)
586 7946512 : fac010(lay) = (1._r8 - fs(lay)) * fac10(lay)
587 7946512 : fac100(lay) = fs(lay) * fac00(lay)
588 7946512 : fac110(lay) = fs(lay) * fac10(lay)
589 : endif
590 9219976 : if (specparm1(lay) .lt. 0.125_r8) then
591 152313 : p = fs1(lay) - 1
592 152313 : p4 = p**4
593 152313 : fk0 = p4
594 152313 : fk1 = 1 - p - 2.0_r8*p4
595 152313 : fk2 = p + p4
596 152313 : fac001(lay) = fk0*fac01(lay)
597 152313 : fac101(lay) = fk1*fac01(lay)
598 152313 : fac201(lay) = fk2*fac01(lay)
599 152313 : fac011(lay) = fk0*fac11(lay)
600 152313 : fac111(lay) = fk1*fac11(lay)
601 152313 : fac211(lay) = fk2*fac11(lay)
602 8581663 : else if (specparm1(lay) .gt. 0.875_r8) then
603 351185 : p = -fs1(lay)
604 351185 : p4 = p**4
605 351185 : fk0 = p4
606 351185 : fk1 = 1 - p - 2.0_r8*p4
607 351185 : fk2 = p + p4
608 351185 : fac001(lay) = fk0*fac01(lay)
609 351185 : fac101(lay) = fk1*fac01(lay)
610 351185 : fac201(lay) = fk2*fac01(lay)
611 351185 : fac011(lay) = fk0*fac11(lay)
612 351185 : fac111(lay) = fk1*fac11(lay)
613 351185 : fac211(lay) = fk2*fac11(lay)
614 : else
615 8230478 : fac001(lay) = (1._r8 - fs1(lay)) * fac01(lay)
616 8230478 : fac011(lay) = (1._r8 - fs1(lay)) * fac11(lay)
617 8230478 : fac101(lay) = fs1(lay) * fac01(lay)
618 8230478 : fac111(lay) = fs1(lay) * fac11(lay)
619 : endif
620 : enddo
621 9219976 : do lay = 1, laytrop
622 :
623 148963592 : do ig = 1, ng3
624 698718080 : tauself = selffac(lay)* (selfref(indself(lay),ig) + selffrac(lay) * &
625 838461696 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
626 698718080 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
627 838461696 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
628 838461696 : n2om1 = ka_mn2o(jmn2o(lay),indminor(lay),ig) + fmn2o(lay) * &
629 978205312 : (ka_mn2o(jmn2o(lay)+1,indminor(lay),ig) - ka_mn2o(jmn2o(lay),indminor(lay),ig))
630 838461696 : n2om2 = ka_mn2o(jmn2o(lay),indminor(lay)+1,ig) + fmn2o(lay) * &
631 978205312 : (ka_mn2o(jmn2o(lay)+1,indminor(lay)+1,ig) - ka_mn2o(jmn2o(lay),indminor(lay)+1,ig))
632 139743616 : absn2o = n2om1 + minorfrac(lay) * (n2om2 - n2om1)
633 :
634 139743616 : if (specparm(lay) .lt. 0.125_r8) then
635 12324784 : tau_major = speccomb(lay) * &
636 49299136 : (fac000(lay) * absa(ind0(lay),ig) + &
637 49299136 : fac100(lay) * absa(ind0(lay)+1,ig) + &
638 49299136 : fac200(lay) * absa(ind0(lay)+2,ig) + &
639 49299136 : fac010(lay) * absa(ind0(lay)+9,ig) + &
640 49299136 : fac110(lay) * absa(ind0(lay)+10,ig) + &
641 271145248 : fac210(lay) * absa(ind0(lay)+11,ig))
642 127418832 : else if (specparm(lay) .gt. 0.875_r8) then
643 274640 : tau_major = speccomb(lay) * &
644 1098560 : (fac200(lay) * absa(ind0(lay)-1,ig) + &
645 1098560 : fac100(lay) * absa(ind0(lay),ig) + &
646 1098560 : fac000(lay) * absa(ind0(lay)+1,ig) + &
647 1098560 : fac210(lay) * absa(ind0(lay)+8,ig) + &
648 1098560 : fac110(lay) * absa(ind0(lay)+9,ig) + &
649 6042080 : fac010(lay) * absa(ind0(lay)+10,ig))
650 : else
651 127144192 : tau_major = speccomb(lay) * &
652 508576768 : (fac000(lay) * absa(ind0(lay),ig) + &
653 508576768 : fac100(lay) * absa(ind0(lay)+1,ig) + &
654 508576768 : fac010(lay) * absa(ind0(lay)+9,ig) + &
655 1780018688 : fac110(lay) * absa(ind0(lay)+10,ig))
656 : endif
657 :
658 139743616 : if (specparm1(lay) .lt. 0.125_r8) then
659 2437008 : tau_major1 = speccomb1(lay) * &
660 9748032 : (fac001(lay) * absa(ind1(lay),ig) + &
661 9748032 : fac101(lay) * absa(ind1(lay)+1,ig) + &
662 9748032 : fac201(lay) * absa(ind1(lay)+2,ig) + &
663 9748032 : fac011(lay) * absa(ind1(lay)+9,ig) + &
664 9748032 : fac111(lay) * absa(ind1(lay)+10,ig) + &
665 53614176 : fac211(lay) * absa(ind1(lay)+11,ig))
666 137306608 : else if (specparm1(lay) .gt. 0.875_r8) then
667 5618960 : tau_major1 = speccomb1(lay) * &
668 22475840 : (fac201(lay) * absa(ind1(lay)-1,ig) + &
669 22475840 : fac101(lay) * absa(ind1(lay),ig) + &
670 22475840 : fac001(lay) * absa(ind1(lay)+1,ig) + &
671 22475840 : fac211(lay) * absa(ind1(lay)+8,ig) + &
672 22475840 : fac111(lay) * absa(ind1(lay)+9,ig) + &
673 123617120 : fac011(lay) * absa(ind1(lay)+10,ig))
674 : else
675 131687648 : tau_major1 = speccomb1(lay) * &
676 526750592 : (fac001(lay) * absa(ind1(lay),ig) + &
677 526750592 : fac101(lay) * absa(ind1(lay)+1,ig) + &
678 526750592 : fac011(lay) * absa(ind1(lay)+9,ig) + &
679 1843627072 : fac111(lay) * absa(ind1(lay)+10,ig))
680 : endif
681 :
682 279487232 : taug(lay,ngs2+ig) = tau_major + tau_major1 &
683 : + tauself + taufor &
684 419230848 : + adjcoln2o(lay)*absn2o
685 838461696 : fracs(lay,ngs2+ig) = fracrefa(ig,jpl(lay)) + fpl(lay) * &
686 986939288 : (fracrefa(ig,jpl(lay)+1)-fracrefa(ig,jpl(lay)))
687 : enddo
688 : enddo
689 :
690 : ! Upper atmosphere loop
691 4874024 : do lay = laytrop+1, nlayers
692 :
693 4388024 : speccomb(lay) = colh2o(lay) + rat_h2oco2(lay)*colco2(lay)
694 4388024 : specparm(lay) = colh2o(lay)/speccomb(lay)
695 4388024 : if (specparm(lay) .ge. oneminus) specparm(lay) = oneminus
696 4388024 : specmult(lay) = 4._r8*(specparm(lay))
697 4388024 : js(lay) = 1 + int(specmult(lay))
698 4388024 : fs(lay) = mod(specmult(lay),1.0_r8)
699 :
700 4388024 : speccomb1(lay) = colh2o(lay) + rat_h2oco2_1(lay)*colco2(lay)
701 4388024 : specparm1(lay) = colh2o(lay)/speccomb1(lay)
702 4388024 : if (specparm1(lay) .ge. oneminus) specparm1(lay) = oneminus
703 4388024 : specmult1(lay) = 4._r8*(specparm1(lay))
704 4388024 : js1(lay) = 1 + int(specmult1(lay))
705 4388024 : fs1(lay) = mod(specmult1(lay),1.0_r8)
706 :
707 4388024 : fac000(lay) = (1._r8 - fs(lay)) * fac00(lay)
708 4388024 : fac010(lay) = (1._r8 - fs(lay)) * fac10(lay)
709 4388024 : fac100(lay) = fs(lay) * fac00(lay)
710 4388024 : fac110(lay) = fs(lay) * fac10(lay)
711 4388024 : fac001(lay) = (1._r8 - fs1(lay)) * fac01(lay)
712 4388024 : fac011(lay) = (1._r8 - fs1(lay)) * fac11(lay)
713 4388024 : fac101(lay) = fs1(lay) * fac01(lay)
714 4388024 : fac111(lay) = fs1(lay) * fac11(lay)
715 :
716 4388024 : speccomb_mn2o(lay) = colh2o(lay) + refrat_m_b*colco2(lay)
717 4388024 : specparm_mn2o(lay) = colh2o(lay)/speccomb_mn2o(lay)
718 4388024 : if (specparm_mn2o(lay) .ge. oneminus) specparm_mn2o(lay) = oneminus
719 4388024 : specmult_mn2o(lay) = 4._r8*specparm_mn2o(lay)
720 4388024 : jmn2o(lay) = 1 + int(specmult_mn2o(lay))
721 4388024 : fmn2o(lay) = mod(specmult_mn2o(lay),1.0_r8)
722 4388024 : fmn2omf(lay) = minorfrac(lay)*fmn2o(lay)
723 : ! In atmospheres where the amount of N2O is too great to be considered
724 : ! a minor species, adjust the column amount of N2O by an empirical factor
725 : ! to obtain the proper contribution.
726 4388024 : chi_n2o(lay) = coln2o(lay)/coldry(lay)
727 4874024 : ratn2o(lay) = 1.e20*chi_n2o(lay)/chi_mls(4,jp(lay)+1)
728 : enddo
729 4874024 : do lay = laytrop+1, nlayers
730 4874024 : if (ratn2o(lay) .gt. 1.5_r8) then
731 2441520 : adjfac(lay) = 0.5_r8+(ratn2o(lay)-0.5_r8)**0.65_r8
732 2441520 : adjcoln2o(lay) = adjfac(lay)*chi_mls(4,jp(lay)+1)*coldry(lay)*1.e-20_r8
733 : else
734 1946504 : adjcoln2o(lay) = coln2o(lay)
735 : endif
736 : enddo
737 4874024 : do lay = laytrop+1, nlayers
738 :
739 4388024 : speccomb_planck(lay) = colh2o(lay)+refrat_planck_b*colco2(lay)
740 4388024 : specparm_planck(lay) = colh2o(lay)/speccomb_planck(lay)
741 4388024 : if (specparm_planck(lay) .ge. oneminus) specparm_planck(lay)=oneminus
742 4388024 : specmult_planck(lay) = 4._r8*specparm_planck(lay)
743 4388024 : jpl(lay)= 1 + int(specmult_planck(lay))
744 4388024 : fpl(lay) = mod(specmult_planck(lay),1.0_r8)
745 :
746 4388024 : ind0(lay) = ((jp(lay)-13)*5+(jt(lay)-1))*nspb(3) + js(lay)
747 4874024 : ind1(lay) = ((jp(lay)-12)*5+(jt1(lay)-1))*nspb(3) + js1(lay)
748 : enddo
749 4874024 : do lay = laytrop+1, nlayers
750 :
751 75082408 : do ig = 1, ng3
752 280833536 : taufor = forfac(lay) * (forref(indfor(lay),ig) + &
753 351041920 : forfrac(lay) * (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
754 421250304 : n2om1 = kb_mn2o(jmn2o(lay),indminor(lay),ig) + fmn2o(lay) * &
755 491458688 : (kb_mn2o(jmn2o(lay)+1,indminor(lay),ig)-kb_mn2o(jmn2o(lay),indminor(lay),ig))
756 421250304 : n2om2 = kb_mn2o(jmn2o(lay),indminor(lay)+1,ig) + fmn2o(lay) * &
757 491458688 : (kb_mn2o(jmn2o(lay)+1,indminor(lay)+1,ig)-kb_mn2o(jmn2o(lay),indminor(lay)+1,ig))
758 70208384 : absn2o = n2om1 + minorfrac(lay) * (n2om2 - n2om1)
759 210625152 : taug(lay,ngs2+ig) = speccomb(lay) * &
760 280833536 : (fac000(lay) * absb(ind0(lay),ig) + &
761 280833536 : fac100(lay) * absb(ind0(lay)+1,ig) + &
762 280833536 : fac010(lay) * absb(ind0(lay)+5,ig) + &
763 280833536 : fac110(lay) * absb(ind0(lay)+6,ig)) &
764 70208384 : + speccomb1(lay) * &
765 280833536 : (fac001(lay) * absb(ind1(lay),ig) + &
766 280833536 : fac101(lay) * absb(ind1(lay)+1,ig) + &
767 280833536 : fac011(lay) * absb(ind1(lay)+5,ig) + &
768 280833536 : fac111(lay) * absb(ind1(lay)+6,ig)) &
769 : + taufor &
770 2597710208 : + adjcoln2o(lay)*absn2o
771 421250304 : fracs(lay,ngs2+ig) = fracrefb(ig,jpl(lay)) + fpl(lay) * &
772 495846712 : (fracrefb(ig,jpl(lay)+1)-fracrefb(ig,jpl(lay)))
773 : enddo
774 : enddo
775 :
776 486000 : end subroutine taugb3
777 :
778 : !----------------------------------------------------------------------------
779 486000 : subroutine taugb4
780 : !----------------------------------------------------------------------------
781 : !
782 : ! band 4: 630-700 cm-1 (low key - h2o,co2; high key - o3,co2)
783 : !----------------------------------------------------------------------------
784 :
785 : ! ------- Modules -------
786 :
787 : use parrrtm, only : ng4, ngs3
788 : use rrlw_ref, only : chi_mls
789 : use rrlw_kg04, only : fracrefa, fracrefb, absa, absb, &
790 : selfref, forref
791 :
792 : ! ------- Declarations -------
793 :
794 : ! Local
795 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
796 972000 : integer, dimension(nlayers) :: js, js1, jpl
797 972000 : real(kind=r8), dimension(nlayers) :: speccomb, specparm, specmult, fs
798 972000 : real(kind=r8), dimension(nlayers) :: speccomb1, specparm1, specmult1, fs1
799 972000 : real(kind=r8), dimension(nlayers) :: speccomb_planck, specparm_planck, specmult_planck, fpl
800 : real(kind=r8) :: p, p4, fk0, fk1, fk2
801 972000 : real(kind=r8), dimension(nlayers) :: fac000, fac100, fac200, fac010, fac110, fac210
802 972000 : real(kind=r8), dimension(nlayers) :: fac001, fac101, fac201, fac011, fac111, fac211
803 : real(kind=r8) :: tauself, taufor
804 : real(kind=r8) :: refrat_planck_a, refrat_planck_b
805 : real(kind=r8) :: tau_major, tau_major1
806 :
807 :
808 : ! P = 142.5940 mb
809 486000 : refrat_planck_a = chi_mls(1,11)/chi_mls(2,11)
810 :
811 : ! P = 95.58350 mb
812 486000 : refrat_planck_b = chi_mls(3,13)/chi_mls(2,13)
813 :
814 : ! Compute the optical depth by interpolating in ln(pressure) and
815 : ! temperature, and appropriate species. Below laytrop, the water
816 : ! vapor self-continuum and foreign continuum is interpolated (in temperature)
817 : ! separately.
818 :
819 : ! Lower atmosphere loop
820 9219976 : do lay = 1, laytrop
821 :
822 8733976 : speccomb(lay) = colh2o(lay) + rat_h2oco2(lay)*colco2(lay)
823 8733976 : specparm(lay) = colh2o(lay)/speccomb(lay)
824 8733976 : if (specparm(lay) .ge. oneminus) specparm(lay) = oneminus
825 8733976 : specmult(lay) = 8._r8*(specparm(lay))
826 8733976 : js(lay) = 1 + int(specmult(lay))
827 8733976 : fs(lay) = mod(specmult(lay),1.0_r8)
828 :
829 8733976 : speccomb1(lay) = colh2o(lay) + rat_h2oco2_1(lay)*colco2(lay)
830 8733976 : specparm1(lay) = colh2o(lay)/speccomb1(lay)
831 8733976 : if (specparm1(lay) .ge. oneminus) specparm1(lay) = oneminus
832 8733976 : specmult1(lay) = 8._r8*(specparm1(lay))
833 8733976 : js1(lay) = 1 + int(specmult1(lay))
834 8733976 : fs1(lay) = mod(specmult1(lay),1.0_r8)
835 :
836 8733976 : speccomb_planck(lay) = colh2o(lay)+refrat_planck_a*colco2(lay)
837 8733976 : specparm_planck(lay) = colh2o(lay)/speccomb_planck(lay)
838 8733976 : if (specparm_planck(lay) .ge. oneminus) specparm_planck(lay)=oneminus
839 8733976 : specmult_planck(lay) = 8._r8*specparm_planck(lay)
840 8733976 : jpl(lay)= 1 + int(specmult_planck(lay))
841 8733976 : fpl(lay) = mod(specmult_planck(lay),1.0_r8)
842 :
843 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(4) + js(lay)
844 9219976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(4) + js1(lay)
845 : enddo
846 9219976 : do lay = 1, laytrop
847 8733976 : if (specparm(lay) .lt. 0.125_r8) then
848 770299 : p = fs(lay) - 1
849 770299 : p4 = p**4
850 770299 : fk0 = p4
851 770299 : fk1 = 1 - p - 2.0_r8*p4
852 770299 : fk2 = p + p4
853 770299 : fac000(lay) = fk0*fac00(lay)
854 770299 : fac100(lay) = fk1*fac00(lay)
855 770299 : fac200(lay) = fk2*fac00(lay)
856 770299 : fac010(lay) = fk0*fac10(lay)
857 770299 : fac110(lay) = fk1*fac10(lay)
858 770299 : fac210(lay) = fk2*fac10(lay)
859 7963677 : else if (specparm(lay) .gt. 0.875_r8) then
860 17165 : p = -fs(lay)
861 17165 : p4 = p**4
862 17165 : fk0 = p4
863 17165 : fk1 = 1 - p - 2.0_r8*p4
864 17165 : fk2 = p + p4
865 17165 : fac000(lay) = fk0*fac00(lay)
866 17165 : fac100(lay) = fk1*fac00(lay)
867 17165 : fac200(lay) = fk2*fac00(lay)
868 17165 : fac010(lay) = fk0*fac10(lay)
869 17165 : fac110(lay) = fk1*fac10(lay)
870 17165 : fac210(lay) = fk2*fac10(lay)
871 : else
872 7946512 : fac000(lay) = (1._r8 - fs(lay)) * fac00(lay)
873 7946512 : fac010(lay) = (1._r8 - fs(lay)) * fac10(lay)
874 7946512 : fac100(lay) = fs(lay) * fac00(lay)
875 7946512 : fac110(lay) = fs(lay) * fac10(lay)
876 : endif
877 :
878 9219976 : if (specparm1(lay) .lt. 0.125_r8) then
879 152313 : p = fs1(lay) - 1
880 152313 : p4 = p**4
881 152313 : fk0 = p4
882 152313 : fk1 = 1 - p - 2.0_r8*p4
883 152313 : fk2 = p + p4
884 152313 : fac001(lay) = fk0*fac01(lay)
885 152313 : fac101(lay) = fk1*fac01(lay)
886 152313 : fac201(lay) = fk2*fac01(lay)
887 152313 : fac011(lay) = fk0*fac11(lay)
888 152313 : fac111(lay) = fk1*fac11(lay)
889 152313 : fac211(lay) = fk2*fac11(lay)
890 8581663 : else if (specparm1(lay) .gt. 0.875_r8) then
891 351185 : p = -fs1(lay)
892 351185 : p4 = p**4
893 351185 : fk0 = p4
894 351185 : fk1 = 1 - p - 2.0_r8*p4
895 351185 : fk2 = p + p4
896 351185 : fac001(lay) = fk0*fac01(lay)
897 351185 : fac101(lay) = fk1*fac01(lay)
898 351185 : fac201(lay) = fk2*fac01(lay)
899 351185 : fac011(lay) = fk0*fac11(lay)
900 351185 : fac111(lay) = fk1*fac11(lay)
901 351185 : fac211(lay) = fk2*fac11(lay)
902 : else
903 8230478 : fac001(lay) = (1._r8 - fs1(lay)) * fac01(lay)
904 8230478 : fac011(lay) = (1._r8 - fs1(lay)) * fac11(lay)
905 8230478 : fac101(lay) = fs1(lay) * fac01(lay)
906 8230478 : fac111(lay) = fs1(lay) * fac11(lay)
907 : endif
908 : enddo
909 9219976 : do lay = 1, laytrop
910 :
911 131495640 : do ig = 1, ng4
912 611378320 : tauself = selffac(lay)* (selfref(indself(lay),ig) + selffrac(lay) * &
913 733653984 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
914 611378320 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
915 733653984 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
916 :
917 122275664 : if (specparm(lay) .lt. 0.125_r8) then
918 10784186 : tau_major = speccomb(lay) * &
919 43136744 : (fac000(lay) * absa(ind0(lay),ig) + &
920 43136744 : fac100(lay) * absa(ind0(lay)+1,ig) + &
921 43136744 : fac200(lay) * absa(ind0(lay)+2,ig) + &
922 43136744 : fac010(lay) * absa(ind0(lay)+9,ig) + &
923 43136744 : fac110(lay) * absa(ind0(lay)+10,ig) + &
924 237252092 : fac210(lay) * absa(ind0(lay)+11,ig))
925 111491478 : else if (specparm(lay) .gt. 0.875_r8) then
926 240310 : tau_major = speccomb(lay) * &
927 961240 : (fac200(lay) * absa(ind0(lay)-1,ig) + &
928 961240 : fac100(lay) * absa(ind0(lay),ig) + &
929 961240 : fac000(lay) * absa(ind0(lay)+1,ig) + &
930 961240 : fac210(lay) * absa(ind0(lay)+8,ig) + &
931 961240 : fac110(lay) * absa(ind0(lay)+9,ig) + &
932 5286820 : fac010(lay) * absa(ind0(lay)+10,ig))
933 : else
934 111251168 : tau_major = speccomb(lay) * &
935 445004672 : (fac000(lay) * absa(ind0(lay),ig) + &
936 445004672 : fac100(lay) * absa(ind0(lay)+1,ig) + &
937 445004672 : fac010(lay) * absa(ind0(lay)+9,ig) + &
938 1557516352 : fac110(lay) * absa(ind0(lay)+10,ig))
939 : endif
940 :
941 122275664 : if (specparm1(lay) .lt. 0.125_r8) then
942 2132382 : tau_major1 = speccomb1(lay) * &
943 8529528 : (fac001(lay) * absa(ind1(lay),ig) + &
944 8529528 : fac101(lay) * absa(ind1(lay)+1,ig) + &
945 8529528 : fac201(lay) * absa(ind1(lay)+2,ig) + &
946 8529528 : fac011(lay) * absa(ind1(lay)+9,ig) + &
947 8529528 : fac111(lay) * absa(ind1(lay)+10,ig) + &
948 46912404 : fac211(lay) * absa(ind1(lay)+11,ig))
949 120143282 : else if (specparm1(lay) .gt. 0.875_r8) then
950 4916590 : tau_major1 = speccomb1(lay) * &
951 19666360 : (fac201(lay) * absa(ind1(lay)-1,ig) + &
952 19666360 : fac101(lay) * absa(ind1(lay),ig) + &
953 19666360 : fac001(lay) * absa(ind1(lay)+1,ig) + &
954 19666360 : fac211(lay) * absa(ind1(lay)+8,ig) + &
955 19666360 : fac111(lay) * absa(ind1(lay)+9,ig) + &
956 108164980 : fac011(lay) * absa(ind1(lay)+10,ig))
957 : else
958 115226692 : tau_major1 = speccomb1(lay) * &
959 460906768 : (fac001(lay) * absa(ind1(lay),ig) + &
960 460906768 : fac101(lay) * absa(ind1(lay)+1,ig) + &
961 460906768 : fac011(lay) * absa(ind1(lay)+9,ig) + &
962 1613173688 : fac111(lay) * absa(ind1(lay)+10,ig))
963 : endif
964 :
965 244551328 : taug(lay,ngs3+ig) = tau_major + tau_major1 &
966 244551328 : + tauself + taufor
967 733653984 : fracs(lay,ngs3+ig) = fracrefa(ig,jpl(lay)) + fpl(lay) * &
968 864663624 : (fracrefa(ig,jpl(lay)+1)-fracrefa(ig,jpl(lay)))
969 : enddo
970 : enddo
971 :
972 : ! Upper atmosphere loop
973 4874024 : do lay = laytrop+1, nlayers
974 :
975 4388024 : speccomb(lay) = colo3(lay) + rat_o3co2(lay)*colco2(lay)
976 4388024 : specparm(lay) = colo3(lay)/speccomb(lay)
977 4388024 : if (specparm(lay) .ge. oneminus) specparm(lay) = oneminus
978 4388024 : specmult(lay) = 4._r8*(specparm(lay))
979 4388024 : js(lay) = 1 + int(specmult(lay))
980 4388024 : fs(lay) = mod(specmult(lay),1.0_r8)
981 :
982 4388024 : speccomb1(lay) = colo3(lay) + rat_o3co2_1(lay)*colco2(lay)
983 4388024 : specparm1(lay) = colo3(lay)/speccomb1(lay)
984 4388024 : if (specparm1(lay) .ge. oneminus) specparm1(lay) = oneminus
985 4388024 : specmult1(lay) = 4._r8*(specparm1(lay))
986 4388024 : js1(lay) = 1 + int(specmult1(lay))
987 4388024 : fs1(lay) = mod(specmult1(lay),1.0_r8)
988 :
989 4388024 : fac000(lay) = (1._r8 - fs(lay)) * fac00(lay)
990 4388024 : fac010(lay) = (1._r8 - fs(lay)) * fac10(lay)
991 4388024 : fac100(lay) = fs(lay) * fac00(lay)
992 4388024 : fac110(lay) = fs(lay) * fac10(lay)
993 4388024 : fac001(lay) = (1._r8 - fs1(lay)) * fac01(lay)
994 4388024 : fac011(lay) = (1._r8 - fs1(lay)) * fac11(lay)
995 4388024 : fac101(lay) = fs1(lay) * fac01(lay)
996 4388024 : fac111(lay) = fs1(lay) * fac11(lay)
997 :
998 4388024 : speccomb_planck(lay) = colo3(lay)+refrat_planck_b*colco2(lay)
999 4388024 : specparm_planck(lay) = colo3(lay)/speccomb_planck(lay)
1000 4388024 : if (specparm_planck(lay) .ge. oneminus) specparm_planck(lay)=oneminus
1001 4388024 : specmult_planck(lay) = 4._r8*specparm_planck(lay)
1002 4388024 : jpl(lay)= 1 + int(specmult_planck(lay))
1003 4388024 : fpl(lay) = mod(specmult_planck(lay),1.0_r8)
1004 :
1005 4388024 : ind0(lay) = ((jp(lay)-13)*5+(jt(lay)-1))*nspb(4) + js(lay)
1006 4388024 : ind1(lay) = ((jp(lay)-12)*5+(jt1(lay)-1))*nspb(4) + js1(lay)
1007 65820360 : do ig = 1, ng4
1008 184297008 : taug(lay,ngs3+ig) = speccomb(lay) * &
1009 245729344 : (fac000(lay) * absb(ind0(lay),ig) + &
1010 245729344 : fac100(lay) * absb(ind0(lay)+1,ig) + &
1011 245729344 : fac010(lay) * absb(ind0(lay)+5,ig) + &
1012 245729344 : fac110(lay) * absb(ind0(lay)+6,ig)) &
1013 61432336 : + speccomb1(lay) * &
1014 245729344 : (fac001(lay) * absb(ind1(lay),ig) + &
1015 245729344 : fac101(lay) * absb(ind1(lay)+1,ig) + &
1016 245729344 : fac011(lay) * absb(ind1(lay)+5,ig) + &
1017 2027267088 : fac111(lay) * absb(ind1(lay)+6,ig))
1018 368594016 : fracs(lay,ngs3+ig) = fracrefb(ig,jpl(lay)) + fpl(lay) * &
1019 434414376 : (fracrefb(ig,jpl(lay)+1)-fracrefb(ig,jpl(lay)))
1020 : enddo
1021 :
1022 : ! Empirical modification to code to improve stratospheric cooling rates
1023 : ! for co2. Revised to apply weighting for g-point reduction in this band.
1024 :
1025 4388024 : taug(lay,ngs3+8)=taug(lay,ngs3+8)*0.92
1026 4388024 : taug(lay,ngs3+9)=taug(lay,ngs3+9)*0.88
1027 4388024 : taug(lay,ngs3+10)=taug(lay,ngs3+10)*1.07
1028 4388024 : taug(lay,ngs3+11)=taug(lay,ngs3+11)*1.1
1029 4388024 : taug(lay,ngs3+12)=taug(lay,ngs3+12)*0.99
1030 4388024 : taug(lay,ngs3+13)=taug(lay,ngs3+13)*0.88
1031 4874024 : taug(lay,ngs3+14)=taug(lay,ngs3+14)*0.943
1032 :
1033 : enddo
1034 :
1035 486000 : end subroutine taugb4
1036 :
1037 : !----------------------------------------------------------------------------
1038 486000 : subroutine taugb5
1039 : !----------------------------------------------------------------------------
1040 : !
1041 : ! band 5: 700-820 cm-1 (low key - h2o,co2; low minor - o3, ccl4)
1042 : ! (high key - o3,co2)
1043 : !----------------------------------------------------------------------------
1044 :
1045 : ! ------- Modules -------
1046 :
1047 : use parrrtm, only : ng5, ngs4
1048 : use rrlw_ref, only : chi_mls
1049 : use rrlw_kg05, only : fracrefa, fracrefb, absa, absb, &
1050 : ka_mo3, selfref, forref, ccl4
1051 :
1052 : ! ------- Declarations -------
1053 :
1054 : ! Local
1055 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
1056 972000 : integer, dimension(nlayers) :: js, js1, jmo3, jpl
1057 972000 : real(kind=r8), dimension(nlayers) :: speccomb, specparm, specmult, fs
1058 972000 : real(kind=r8), dimension(nlayers) :: speccomb1, specparm1, specmult1, fs1
1059 972000 : real(kind=r8), dimension(nlayers) :: speccomb_mo3, specparm_mo3, specmult_mo3, fmo3
1060 972000 : real(kind=r8), dimension(nlayers) :: speccomb_planck, specparm_planck, specmult_planck, fpl
1061 : real(kind=r8) :: p, p4, fk0, fk1, fk2
1062 972000 : real(kind=r8), dimension(nlayers) :: fac000, fac100, fac200, fac010, fac110, fac210
1063 972000 : real(kind=r8), dimension(nlayers) :: fac001, fac101, fac201, fac011, fac111, fac211
1064 : real(kind=r8) :: tauself, taufor, o3m1, o3m2, abso3
1065 : real(kind=r8) :: refrat_planck_a, refrat_planck_b, refrat_m_a
1066 : real(kind=r8) :: tau_major, tau_major1
1067 :
1068 :
1069 : ! Minor gas mapping level :
1070 : ! lower - o3, p = 317.34 mbar, t = 240.77 k
1071 : ! lower - ccl4
1072 :
1073 : ! Calculate reference ratio to be used in calculation of Planck
1074 : ! fraction in lower/upper atmosphere.
1075 :
1076 : ! P = 473.420 mb
1077 486000 : refrat_planck_a = chi_mls(1,5)/chi_mls(2,5)
1078 :
1079 : ! P = 0.2369 mb
1080 486000 : refrat_planck_b = chi_mls(3,43)/chi_mls(2,43)
1081 :
1082 : ! P = 317.3480
1083 486000 : refrat_m_a = chi_mls(1,7)/chi_mls(2,7)
1084 :
1085 : ! Compute the optical depth by interpolating in ln(pressure) and
1086 : ! temperature, and appropriate species. Below laytrop, the
1087 : ! water vapor self-continuum and foreign continuum is
1088 : ! interpolated (in temperature) separately.
1089 :
1090 : ! Lower atmosphere loop
1091 9219976 : do lay = 1, laytrop
1092 :
1093 8733976 : speccomb(lay) = colh2o(lay) + rat_h2oco2(lay)*colco2(lay)
1094 8733976 : specparm(lay) = colh2o(lay)/speccomb(lay)
1095 8733976 : if (specparm(lay) .ge. oneminus) specparm(lay) = oneminus
1096 8733976 : specmult(lay) = 8._r8*(specparm(lay))
1097 8733976 : js(lay) = 1 + int(specmult(lay))
1098 8733976 : fs(lay) = mod(specmult(lay),1.0_r8)
1099 :
1100 8733976 : speccomb1(lay) = colh2o(lay) + rat_h2oco2_1(lay)*colco2(lay)
1101 8733976 : specparm1(lay) = colh2o(lay)/speccomb1(lay)
1102 8733976 : if (specparm1(lay) .ge. oneminus) specparm1(lay) = oneminus
1103 8733976 : specmult1(lay) = 8._r8*(specparm1(lay))
1104 8733976 : js1(lay) = 1 + int(specmult1(lay))
1105 8733976 : fs1(lay) = mod(specmult1(lay),1.0_r8)
1106 :
1107 8733976 : speccomb_mo3(lay) = colh2o(lay) + refrat_m_a*colco2(lay)
1108 8733976 : specparm_mo3(lay) = colh2o(lay)/speccomb_mo3(lay)
1109 8733976 : if (specparm_mo3(lay) .ge. oneminus) specparm_mo3(lay) = oneminus
1110 8733976 : specmult_mo3(lay) = 8._r8*specparm_mo3(lay)
1111 8733976 : jmo3(lay) = 1 + int(specmult_mo3(lay))
1112 8733976 : fmo3(lay) = mod(specmult_mo3(lay),1.0_r8)
1113 :
1114 8733976 : speccomb_planck(lay) = colh2o(lay)+refrat_planck_a*colco2(lay)
1115 8733976 : specparm_planck(lay) = colh2o(lay)/speccomb_planck(lay)
1116 8733976 : if (specparm_planck(lay) .ge. oneminus) specparm_planck(lay)=oneminus
1117 8733976 : specmult_planck(lay) = 8._r8*specparm_planck(lay)
1118 8733976 : jpl(lay)= 1 + int(specmult_planck(lay))
1119 8733976 : fpl(lay) = mod(specmult_planck(lay),1.0_r8)
1120 :
1121 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(5) + js(lay)
1122 9219976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(5) + js1(lay)
1123 : enddo
1124 9219976 : do lay = 1, laytrop
1125 8733976 : if (specparm(lay) .lt. 0.125_r8) then
1126 770299 : p = fs(lay) - 1
1127 770299 : p4 = p**4
1128 770299 : fk0 = p4
1129 770299 : fk1 = 1 - p - 2.0_r8*p4
1130 770299 : fk2 = p + p4
1131 770299 : fac000(lay) = fk0*fac00(lay)
1132 770299 : fac100(lay) = fk1*fac00(lay)
1133 770299 : fac200(lay) = fk2*fac00(lay)
1134 770299 : fac010(lay) = fk0*fac10(lay)
1135 770299 : fac110(lay) = fk1*fac10(lay)
1136 770299 : fac210(lay) = fk2*fac10(lay)
1137 7963677 : else if (specparm(lay) .gt. 0.875_r8) then
1138 17165 : p = -fs(lay)
1139 17165 : p4 = p**4
1140 17165 : fk0 = p4
1141 17165 : fk1 = 1 - p - 2.0_r8*p4
1142 17165 : fk2 = p + p4
1143 17165 : fac000(lay) = fk0*fac00(lay)
1144 17165 : fac100(lay) = fk1*fac00(lay)
1145 17165 : fac200(lay) = fk2*fac00(lay)
1146 17165 : fac010(lay) = fk0*fac10(lay)
1147 17165 : fac110(lay) = fk1*fac10(lay)
1148 17165 : fac210(lay) = fk2*fac10(lay)
1149 : else
1150 7946512 : fac000(lay) = (1._r8 - fs(lay)) * fac00(lay)
1151 7946512 : fac010(lay) = (1._r8 - fs(lay)) * fac10(lay)
1152 7946512 : fac100(lay) = fs(lay) * fac00(lay)
1153 7946512 : fac110(lay) = fs(lay) * fac10(lay)
1154 : endif
1155 :
1156 9219976 : if (specparm1(lay) .lt. 0.125_r8) then
1157 152313 : p = fs1(lay) - 1
1158 152313 : p4 = p**4
1159 152313 : fk0 = p4
1160 152313 : fk1 = 1 - p - 2.0_r8*p4
1161 152313 : fk2 = p + p4
1162 152313 : fac001(lay) = fk0*fac01(lay)
1163 152313 : fac101(lay) = fk1*fac01(lay)
1164 152313 : fac201(lay) = fk2*fac01(lay)
1165 152313 : fac011(lay) = fk0*fac11(lay)
1166 152313 : fac111(lay) = fk1*fac11(lay)
1167 152313 : fac211(lay) = fk2*fac11(lay)
1168 8581663 : else if (specparm1(lay) .gt. 0.875_r8) then
1169 351185 : p = -fs1(lay)
1170 351185 : p4 = p**4
1171 351185 : fk0 = p4
1172 351185 : fk1 = 1 - p - 2.0_r8*p4
1173 351185 : fk2 = p + p4
1174 351185 : fac001(lay) = fk0*fac01(lay)
1175 351185 : fac101(lay) = fk1*fac01(lay)
1176 351185 : fac201(lay) = fk2*fac01(lay)
1177 351185 : fac011(lay) = fk0*fac11(lay)
1178 351185 : fac111(lay) = fk1*fac11(lay)
1179 351185 : fac211(lay) = fk2*fac11(lay)
1180 : else
1181 8230478 : fac001(lay) = (1._r8 - fs1(lay)) * fac01(lay)
1182 8230478 : fac011(lay) = (1._r8 - fs1(lay)) * fac11(lay)
1183 8230478 : fac101(lay) = fs1(lay) * fac01(lay)
1184 8230478 : fac111(lay) = fs1(lay) * fac11(lay)
1185 : endif
1186 : enddo
1187 9219976 : do lay = 1, laytrop
1188 :
1189 148963592 : do ig = 1, ng5
1190 698718080 : tauself = selffac(lay) * (selfref(indself(lay),ig) + selffrac(lay) * &
1191 838461696 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
1192 698718080 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
1193 838461696 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
1194 838461696 : o3m1 = ka_mo3(jmo3(lay),indminor(lay),ig) + fmo3(lay) * &
1195 978205312 : (ka_mo3(jmo3(lay)+1,indminor(lay),ig)-ka_mo3(jmo3(lay),indminor(lay),ig))
1196 838461696 : o3m2 = ka_mo3(jmo3(lay),indminor(lay)+1,ig) + fmo3(lay) * &
1197 978205312 : (ka_mo3(jmo3(lay)+1,indminor(lay)+1,ig)-ka_mo3(jmo3(lay),indminor(lay)+1,ig))
1198 139743616 : abso3 = o3m1 + minorfrac(lay)*(o3m2-o3m1)
1199 :
1200 139743616 : if (specparm(lay) .lt. 0.125_r8) then
1201 12324784 : tau_major = speccomb(lay) * &
1202 49299136 : (fac000(lay) * absa(ind0(lay),ig) + &
1203 49299136 : fac100(lay) * absa(ind0(lay)+1,ig) + &
1204 49299136 : fac200(lay) * absa(ind0(lay)+2,ig) + &
1205 49299136 : fac010(lay) * absa(ind0(lay)+9,ig) + &
1206 49299136 : fac110(lay) * absa(ind0(lay)+10,ig) + &
1207 271145248 : fac210(lay) * absa(ind0(lay)+11,ig))
1208 127418832 : else if (specparm(lay) .gt. 0.875_r8) then
1209 274640 : tau_major = speccomb(lay) * &
1210 1098560 : (fac200(lay) * absa(ind0(lay)-1,ig) + &
1211 1098560 : fac100(lay) * absa(ind0(lay),ig) + &
1212 1098560 : fac000(lay) * absa(ind0(lay)+1,ig) + &
1213 1098560 : fac210(lay) * absa(ind0(lay)+8,ig) + &
1214 1098560 : fac110(lay) * absa(ind0(lay)+9,ig) + &
1215 6042080 : fac010(lay) * absa(ind0(lay)+10,ig))
1216 : else
1217 127144192 : tau_major = speccomb(lay) * &
1218 508576768 : (fac000(lay) * absa(ind0(lay),ig) + &
1219 508576768 : fac100(lay) * absa(ind0(lay)+1,ig) + &
1220 508576768 : fac010(lay) * absa(ind0(lay)+9,ig) + &
1221 1780018688 : fac110(lay) * absa(ind0(lay)+10,ig))
1222 : endif
1223 :
1224 139743616 : if (specparm1(lay) .lt. 0.125_r8) then
1225 2437008 : tau_major1 = speccomb1(lay) * &
1226 9748032 : (fac001(lay) * absa(ind1(lay),ig) + &
1227 9748032 : fac101(lay) * absa(ind1(lay)+1,ig) + &
1228 9748032 : fac201(lay) * absa(ind1(lay)+2,ig) + &
1229 9748032 : fac011(lay) * absa(ind1(lay)+9,ig) + &
1230 9748032 : fac111(lay) * absa(ind1(lay)+10,ig) + &
1231 53614176 : fac211(lay) * absa(ind1(lay)+11,ig))
1232 137306608 : else if (specparm1(lay) .gt. 0.875_r8) then
1233 5618960 : tau_major1 = speccomb1(lay) * &
1234 22475840 : (fac201(lay) * absa(ind1(lay)-1,ig) + &
1235 22475840 : fac101(lay) * absa(ind1(lay),ig) + &
1236 22475840 : fac001(lay) * absa(ind1(lay)+1,ig) + &
1237 22475840 : fac211(lay) * absa(ind1(lay)+8,ig) + &
1238 22475840 : fac111(lay) * absa(ind1(lay)+9,ig) + &
1239 123617120 : fac011(lay) * absa(ind1(lay)+10,ig))
1240 : else
1241 131687648 : tau_major1 = speccomb1(lay) * &
1242 526750592 : (fac001(lay) * absa(ind1(lay),ig) + &
1243 526750592 : fac101(lay) * absa(ind1(lay)+1,ig) + &
1244 526750592 : fac011(lay) * absa(ind1(lay)+9,ig) + &
1245 1843627072 : fac111(lay) * absa(ind1(lay)+10,ig))
1246 : endif
1247 :
1248 279487232 : taug(lay,ngs4+ig) = tau_major + tau_major1 &
1249 : + tauself + taufor &
1250 139743616 : + abso3*colo3(lay) &
1251 558974464 : + wx(1,lay) * ccl4(ig)
1252 838461696 : fracs(lay,ngs4+ig) = fracrefa(ig,jpl(lay)) + fpl(lay) * &
1253 986939288 : (fracrefa(ig,jpl(lay)+1)-fracrefa(ig,jpl(lay)))
1254 : enddo
1255 : enddo
1256 :
1257 : ! Upper atmosphere loop
1258 4874024 : do lay = laytrop+1, nlayers
1259 :
1260 4388024 : speccomb(lay) = colo3(lay) + rat_o3co2(lay)*colco2(lay)
1261 4388024 : specparm(lay) = colo3(lay)/speccomb(lay)
1262 4388024 : if (specparm(lay) .ge. oneminus) specparm(lay) = oneminus
1263 4388024 : specmult(lay) = 4._r8*(specparm(lay))
1264 4388024 : js(lay) = 1 + int(specmult(lay))
1265 4388024 : fs(lay) = mod(specmult(lay),1.0_r8)
1266 :
1267 4388024 : speccomb1(lay) = colo3(lay) + rat_o3co2_1(lay)*colco2(lay)
1268 4388024 : specparm1(lay) = colo3(lay)/speccomb1(lay)
1269 4388024 : if (specparm1(lay) .ge. oneminus) specparm1(lay) = oneminus
1270 4388024 : specmult1(lay) = 4._r8*(specparm1(lay))
1271 4388024 : js1(lay) = 1 + int(specmult1(lay))
1272 4388024 : fs1(lay) = mod(specmult1(lay),1.0_r8)
1273 :
1274 4388024 : fac000(lay) = (1._r8 - fs(lay)) * fac00(lay)
1275 4388024 : fac010(lay) = (1._r8 - fs(lay)) * fac10(lay)
1276 4388024 : fac100(lay) = fs(lay) * fac00(lay)
1277 4388024 : fac110(lay) = fs(lay) * fac10(lay)
1278 4388024 : fac001(lay) = (1._r8 - fs1(lay)) * fac01(lay)
1279 4388024 : fac011(lay) = (1._r8 - fs1(lay)) * fac11(lay)
1280 4388024 : fac101(lay) = fs1(lay) * fac01(lay)
1281 4388024 : fac111(lay) = fs1(lay) * fac11(lay)
1282 :
1283 4388024 : speccomb_planck(lay) = colo3(lay)+refrat_planck_b*colco2(lay)
1284 4388024 : specparm_planck(lay) = colo3(lay)/speccomb_planck(lay)
1285 4388024 : if (specparm_planck(lay) .ge. oneminus) specparm_planck(lay)=oneminus
1286 4388024 : specmult_planck(lay) = 4._r8*specparm_planck(lay)
1287 4388024 : jpl(lay)= 1 + int(specmult_planck(lay))
1288 4388024 : fpl(lay) = mod(specmult_planck(lay),1.0_r8)
1289 :
1290 4388024 : ind0(lay) = ((jp(lay)-13)*5+(jt(lay)-1))*nspb(5) + js(lay)
1291 4874024 : ind1(lay) = ((jp(lay)-12)*5+(jt1(lay)-1))*nspb(5) + js1(lay)
1292 : enddo
1293 4874024 : do lay = laytrop+1, nlayers
1294 75082408 : do ig = 1, ng5
1295 210625152 : taug(lay,ngs4+ig) = speccomb(lay) * &
1296 280833536 : (fac000(lay) * absb(ind0(lay),ig) + &
1297 280833536 : fac100(lay) * absb(ind0(lay)+1,ig) + &
1298 280833536 : fac010(lay) * absb(ind0(lay)+5,ig) + &
1299 280833536 : fac110(lay) * absb(ind0(lay)+6,ig)) &
1300 70208384 : + speccomb1(lay) * &
1301 280833536 : (fac001(lay) * absb(ind1(lay),ig) + &
1302 280833536 : fac101(lay) * absb(ind1(lay)+1,ig) + &
1303 280833536 : fac011(lay) * absb(ind1(lay)+5,ig) + &
1304 280833536 : fac111(lay) * absb(ind1(lay)+6,ig)) &
1305 2597710208 : + wx(1,lay) * ccl4(ig)
1306 421250304 : fracs(lay,ngs4+ig) = fracrefb(ig,jpl(lay)) + fpl(lay) * &
1307 495846712 : (fracrefb(ig,jpl(lay)+1)-fracrefb(ig,jpl(lay)))
1308 : enddo
1309 : enddo
1310 :
1311 486000 : end subroutine taugb5
1312 :
1313 : !----------------------------------------------------------------------------
1314 486000 : subroutine taugb6
1315 : !----------------------------------------------------------------------------
1316 : !
1317 : ! band 6: 820-980 cm-1 (low key - h2o; low minor - co2)
1318 : ! (high key - nothing; high minor - cfc11, cfc12)
1319 : !----------------------------------------------------------------------------
1320 :
1321 : ! ------- Modules -------
1322 :
1323 : use parrrtm, only : ng6, ngs5
1324 : use rrlw_ref, only : chi_mls
1325 : use rrlw_kg06, only : fracrefa, absa, ka_mco2, &
1326 : selfref, forref, cfc11adj, cfc12
1327 :
1328 : ! ------- Declarations -------
1329 :
1330 : ! Local
1331 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
1332 972000 : real(kind=r8), dimension(nlayers) :: chi_co2, ratco2, adjfac, adjcolco2
1333 : real(kind=r8) :: tauself, taufor, absco2
1334 :
1335 :
1336 : ! Minor gas mapping level:
1337 : ! lower - co2, p = 706.2720 mb, t = 294.2 k
1338 : ! upper - cfc11, cfc12
1339 :
1340 : ! Compute the optical depth by interpolating in ln(pressure) and
1341 : ! temperature. The water vapor self-continuum and foreign continuum
1342 : ! is interpolated (in temperature) separately.
1343 :
1344 : ! Lower atmosphere loop
1345 9219976 : do lay = 1, laytrop
1346 :
1347 : ! In atmospheres where the amount of CO2 is too great to be considered
1348 : ! a minor species, adjust the column amount of CO2 by an empirical factor
1349 : ! to obtain the proper contribution.
1350 8733976 : chi_co2(lay) = colco2(lay)/(coldry(lay))
1351 9219976 : ratco2(lay) = 1.e20_r8*chi_co2(lay)/chi_mls(2,jp(lay)+1)
1352 : enddo
1353 9219976 : do lay = 1, laytrop
1354 9219976 : if (ratco2(lay) .gt. 3.0_r8) then
1355 0 : adjfac(lay) = 2.0_r8+(ratco2(lay)-2.0_r8)**0.77_r8
1356 0 : adjcolco2(lay) = adjfac(lay)*chi_mls(2,jp(lay)+1)*coldry(lay)*1.e-20_r8
1357 : else
1358 8733976 : adjcolco2(lay) = colco2(lay)
1359 : endif
1360 : enddo
1361 9219976 : do lay = 1, laytrop
1362 :
1363 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(6) + 1
1364 9219976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(6) + 1
1365 : enddo
1366 9219976 : do lay = 1, laytrop
1367 :
1368 79091784 : do ig = 1, ng6
1369 349359040 : tauself = selffac(lay) * (selfref(indself(lay),ig) + selffrac(lay) * &
1370 419230848 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
1371 349359040 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
1372 419230848 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
1373 279487232 : absco2 = (ka_mco2(indminor(lay),ig) + minorfrac(lay) * &
1374 349359040 : (ka_mco2(indminor(lay)+1,ig) - ka_mco2(indminor(lay),ig)))
1375 209615424 : taug(lay,ngs5+ig) = colh2o(lay) * &
1376 279487232 : (fac00(lay) * absa(ind0(lay),ig) + &
1377 279487232 : fac10(lay) * absa(ind0(lay)+1,ig) + &
1378 279487232 : fac01(lay) * absa(ind1(lay),ig) + &
1379 279487232 : fac11(lay) * absa(ind1(lay)+1,ig)) &
1380 : + tauself + taufor &
1381 69871808 : + adjcolco2(lay) * absco2 &
1382 139743616 : + wx(2,lay) * cfc11adj(ig) &
1383 1607051584 : + wx(3,lay) * cfc12(ig)
1384 78605784 : fracs(lay,ngs5+ig) = fracrefa(ig)
1385 : enddo
1386 : enddo
1387 :
1388 : ! Upper atmosphere loop
1389 : ! Nothing important goes on above laytrop in this band.
1390 4874024 : do lay = laytrop+1, nlayers
1391 :
1392 39978216 : do ig = 1, ng6
1393 70208384 : taug(lay,ngs5+ig) = 0.0_r8 &
1394 70208384 : + wx(2,lay) * cfc11adj(ig) &
1395 175520960 : + wx(3,lay) * cfc12(ig)
1396 39492216 : fracs(lay,ngs5+ig) = fracrefa(ig)
1397 : enddo
1398 : enddo
1399 :
1400 486000 : end subroutine taugb6
1401 :
1402 : !----------------------------------------------------------------------------
1403 486000 : subroutine taugb7
1404 : !----------------------------------------------------------------------------
1405 : !
1406 : ! band 7: 980-1080 cm-1 (low key - h2o,o3; low minor - co2)
1407 : ! (high key - o3; high minor - co2)
1408 : !----------------------------------------------------------------------------
1409 :
1410 : ! ------- Modules -------
1411 :
1412 : use parrrtm, only : ng7, ngs6
1413 : use rrlw_ref, only : chi_mls
1414 : use rrlw_kg07, only : fracrefa, fracrefb, absa, absb, &
1415 : ka_mco2, kb_mco2, selfref, forref
1416 :
1417 : ! ------- Declarations -------
1418 :
1419 : ! Local
1420 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
1421 972000 : integer, dimension(nlayers) :: js, js1, jmco2, jpl
1422 972000 : real(kind=r8), dimension(nlayers) :: speccomb, specparm, specmult, fs
1423 972000 : real(kind=r8), dimension(nlayers) :: speccomb1, specparm1, specmult1, fs1
1424 972000 : real(kind=r8), dimension(nlayers) :: speccomb_mco2, specparm_mco2, specmult_mco2, fmco2
1425 972000 : real(kind=r8), dimension(nlayers) :: speccomb_planck, specparm_planck, specmult_planck, fpl
1426 : real(kind=r8) :: p, p4, fk0, fk1, fk2
1427 972000 : real(kind=r8), dimension(nlayers) :: fac000, fac100, fac200, fac010, fac110, fac210
1428 972000 : real(kind=r8), dimension(nlayers) :: fac001, fac101, fac201, fac011, fac111, fac211
1429 : real(kind=r8) :: tauself, taufor, co2m1, co2m2, absco2
1430 972000 : real(kind=r8), dimension(nlayers) :: chi_co2, ratco2, adjfac, adjcolco2
1431 : real(kind=r8) :: refrat_planck_a, refrat_m_a
1432 : real(kind=r8) :: tau_major, tau_major1
1433 :
1434 :
1435 : ! Minor gas mapping level :
1436 : ! lower - co2, p = 706.2620 mbar, t= 278.94 k
1437 : ! upper - co2, p = 12.9350 mbar, t = 234.01 k
1438 :
1439 : ! Calculate reference ratio to be used in calculation of Planck
1440 : ! fraction in lower atmosphere.
1441 :
1442 : ! P = 706.2620 mb
1443 486000 : refrat_planck_a = chi_mls(1,3)/chi_mls(3,3)
1444 :
1445 : ! P = 706.2720 mb
1446 486000 : refrat_m_a = chi_mls(1,3)/chi_mls(3,3)
1447 :
1448 : ! Compute the optical depth by interpolating in ln(pressure),
1449 : ! temperature, and appropriate species. Below laytrop, the water
1450 : ! vapor self-continuum and foreign continuum is interpolated
1451 : ! (in temperature) separately.
1452 :
1453 : ! Lower atmosphere loop
1454 9219976 : do lay = 1, laytrop
1455 :
1456 8733976 : speccomb(lay) = colh2o(lay) + rat_h2oo3(lay)*colo3(lay)
1457 8733976 : specparm(lay) = colh2o(lay)/speccomb(lay)
1458 8733976 : if (specparm(lay) .ge. oneminus) specparm(lay) = oneminus
1459 8733976 : specmult(lay) = 8._r8*(specparm(lay))
1460 8733976 : js(lay) = 1 + int(specmult(lay))
1461 8733976 : fs(lay) = mod(specmult(lay),1.0_r8)
1462 :
1463 8733976 : speccomb1(lay) = colh2o(lay) + rat_h2oo3_1(lay)*colo3(lay)
1464 8733976 : specparm1(lay) = colh2o(lay)/speccomb1(lay)
1465 8733976 : if (specparm1(lay) .ge. oneminus) specparm1(lay) = oneminus
1466 8733976 : specmult1(lay) = 8._r8*(specparm1(lay))
1467 8733976 : js1(lay) = 1 + int(specmult1(lay))
1468 8733976 : fs1(lay) = mod(specmult1(lay),1.0_r8)
1469 :
1470 8733976 : speccomb_mco2(lay) = colh2o(lay) + refrat_m_a*colo3(lay)
1471 8733976 : specparm_mco2(lay) = colh2o(lay)/speccomb_mco2(lay)
1472 8733976 : if (specparm_mco2(lay) .ge. oneminus) specparm_mco2(lay) = oneminus
1473 8733976 : specmult_mco2(lay) = 8._r8*specparm_mco2(lay)
1474 :
1475 8733976 : jmco2(lay) = 1 + int(specmult_mco2(lay))
1476 8733976 : fmco2(lay) = mod(specmult_mco2(lay),1.0_r8)
1477 :
1478 : ! In atmospheres where the amount of CO2 is too great to be considered
1479 : ! a minor species, adjust the column amount of CO2 by an empirical factor
1480 : ! to obtain the proper contribution.
1481 8733976 : chi_co2(lay) = colco2(lay)/(coldry(lay))
1482 8733976 : ratco2(lay) = 1.e20*chi_co2(lay)/chi_mls(2,jp(lay)+1)
1483 8733976 : if (ratco2(lay) .gt. 3.0_r8) then
1484 0 : adjfac(lay) = 3.0_r8+(ratco2(lay)-3.0_r8)**0.79_r8
1485 0 : adjcolco2(lay) = adjfac(lay)*chi_mls(2,jp(lay)+1)*coldry(lay)*1.e-20_r8
1486 : else
1487 8733976 : adjcolco2(lay) = colco2(lay)
1488 : endif
1489 :
1490 8733976 : speccomb_planck(lay) = colh2o(lay)+refrat_planck_a*colo3(lay)
1491 8733976 : specparm_planck(lay) = colh2o(lay)/speccomb_planck(lay)
1492 8733976 : if (specparm_planck(lay) .ge. oneminus) specparm_planck(lay)=oneminus
1493 8733976 : specmult_planck(lay) = 8._r8*specparm_planck(lay)
1494 8733976 : jpl(lay)= 1 + int(specmult_planck(lay))
1495 8733976 : fpl(lay) = mod(specmult_planck(lay),1.0_r8)
1496 :
1497 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(7) + js(lay)
1498 8733976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(7) + js1(lay)
1499 :
1500 8733976 : if (specparm(lay) .lt. 0.125_r8) then
1501 567864 : p = fs(lay) - 1
1502 567864 : p4 = p**4
1503 567864 : fk0 = p4
1504 567864 : fk1 = 1 - p - 2.0_r8*p4
1505 567864 : fk2 = p + p4
1506 567864 : fac000(lay) = fk0*fac00(lay)
1507 567864 : fac100(lay) = fk1*fac00(lay)
1508 567864 : fac200(lay) = fk2*fac00(lay)
1509 567864 : fac010(lay) = fk0*fac10(lay)
1510 567864 : fac110(lay) = fk1*fac10(lay)
1511 567864 : fac210(lay) = fk2*fac10(lay)
1512 8166112 : else if (specparm(lay) .gt. 0.875_r8) then
1513 1359267 : p = -fs(lay)
1514 1359267 : p4 = p**4
1515 1359267 : fk0 = p4
1516 1359267 : fk1 = 1 - p - 2.0_r8*p4
1517 1359267 : fk2 = p + p4
1518 1359267 : fac000(lay) = fk0*fac00(lay)
1519 1359267 : fac100(lay) = fk1*fac00(lay)
1520 1359267 : fac200(lay) = fk2*fac00(lay)
1521 1359267 : fac010(lay) = fk0*fac10(lay)
1522 1359267 : fac110(lay) = fk1*fac10(lay)
1523 1359267 : fac210(lay) = fk2*fac10(lay)
1524 : else
1525 6806845 : fac000(lay) = (1._r8 - fs(lay)) * fac00(lay)
1526 6806845 : fac010(lay) = (1._r8 - fs(lay)) * fac10(lay)
1527 6806845 : fac100(lay) = fs(lay) * fac00(lay)
1528 6806845 : fac110(lay) = fs(lay) * fac10(lay)
1529 : endif
1530 8733976 : if (specparm(lay) .lt. 0.125_r8) then
1531 567864 : p = fs1(lay) - 1
1532 567864 : p4 = p**4
1533 567864 : fk0 = p4
1534 567864 : fk1 = 1 - p - 2.0_r8*p4
1535 567864 : fk2 = p + p4
1536 567864 : fac001(lay) = fk0*fac01(lay)
1537 567864 : fac101(lay) = fk1*fac01(lay)
1538 567864 : fac201(lay) = fk2*fac01(lay)
1539 567864 : fac011(lay) = fk0*fac11(lay)
1540 567864 : fac111(lay) = fk1*fac11(lay)
1541 567864 : fac211(lay) = fk2*fac11(lay)
1542 8166112 : else if (specparm1(lay) .gt. 0.875_r8) then
1543 2645821 : p = -fs1(lay)
1544 2645821 : p4 = p**4
1545 2645821 : fk0 = p4
1546 2645821 : fk1 = 1 - p - 2.0_r8*p4
1547 2645821 : fk2 = p + p4
1548 2645821 : fac001(lay) = fk0*fac01(lay)
1549 2645821 : fac101(lay) = fk1*fac01(lay)
1550 2645821 : fac201(lay) = fk2*fac01(lay)
1551 2645821 : fac011(lay) = fk0*fac11(lay)
1552 2645821 : fac111(lay) = fk1*fac11(lay)
1553 2645821 : fac211(lay) = fk2*fac11(lay)
1554 : else
1555 5520291 : fac001(lay) = (1._r8 - fs1(lay)) * fac01(lay)
1556 5520291 : fac011(lay) = (1._r8 - fs1(lay)) * fac11(lay)
1557 5520291 : fac101(lay) = fs1(lay) * fac01(lay)
1558 5520291 : fac111(lay) = fs1(lay) * fac11(lay)
1559 : endif
1560 :
1561 114027688 : do ig = 1, ng7
1562 524038560 : tauself = selffac(lay)* (selfref(indself(lay),ig) + selffrac(lay) * &
1563 628846272 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
1564 524038560 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
1565 628846272 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
1566 628846272 : co2m1 = ka_mco2(jmco2(lay),indminor(lay),ig) + fmco2(lay) * &
1567 733653984 : (ka_mco2(jmco2(lay)+1,indminor(lay),ig) - ka_mco2(jmco2(lay),indminor(lay),ig))
1568 628846272 : co2m2 = ka_mco2(jmco2(lay),indminor(lay)+1,ig) + fmco2(lay) * &
1569 733653984 : (ka_mco2(jmco2(lay)+1,indminor(lay)+1,ig) - ka_mco2(jmco2(lay),indminor(lay)+1,ig))
1570 104807712 : absco2 = co2m1 + minorfrac(lay) * (co2m2 - co2m1)
1571 :
1572 104807712 : if (specparm(lay) .lt. 0.125_r8) then
1573 6814368 : tau_major = speccomb(lay) * &
1574 27257472 : (fac000(lay) * absa(ind0(lay),ig) + &
1575 27257472 : fac100(lay) * absa(ind0(lay)+1,ig) + &
1576 27257472 : fac200(lay) * absa(ind0(lay)+2,ig) + &
1577 27257472 : fac010(lay) * absa(ind0(lay)+9,ig) + &
1578 27257472 : fac110(lay) * absa(ind0(lay)+10,ig) + &
1579 149916096 : fac210(lay) * absa(ind0(lay)+11,ig))
1580 97993344 : else if (specparm(lay) .gt. 0.875_r8) then
1581 16311204 : tau_major = speccomb(lay) * &
1582 65244816 : (fac200(lay) * absa(ind0(lay)-1,ig) + &
1583 65244816 : fac100(lay) * absa(ind0(lay),ig) + &
1584 65244816 : fac000(lay) * absa(ind0(lay)+1,ig) + &
1585 65244816 : fac210(lay) * absa(ind0(lay)+8,ig) + &
1586 65244816 : fac110(lay) * absa(ind0(lay)+9,ig) + &
1587 358846488 : fac010(lay) * absa(ind0(lay)+10,ig))
1588 : else
1589 81682140 : tau_major = speccomb(lay) * &
1590 326728560 : (fac000(lay) * absa(ind0(lay),ig) + &
1591 326728560 : fac100(lay) * absa(ind0(lay)+1,ig) + &
1592 326728560 : fac010(lay) * absa(ind0(lay)+9,ig) + &
1593 1143549960 : fac110(lay) * absa(ind0(lay)+10,ig))
1594 : endif
1595 :
1596 104807712 : if (specparm1(lay) .lt. 0.125_r8) then
1597 643848 : tau_major1 = speccomb1(lay) * &
1598 2575392 : (fac001(lay) * absa(ind1(lay),ig) + &
1599 2575392 : fac101(lay) * absa(ind1(lay)+1,ig) + &
1600 2575392 : fac201(lay) * absa(ind1(lay)+2,ig) + &
1601 2575392 : fac011(lay) * absa(ind1(lay)+9,ig) + &
1602 2575392 : fac111(lay) * absa(ind1(lay)+10,ig) + &
1603 14164656 : fac211(lay) * absa(ind1(lay)+11,ig))
1604 104163864 : else if (specparm1(lay) .gt. 0.875_r8) then
1605 31749852 : tau_major1 = speccomb1(lay) * &
1606 126999408 : (fac201(lay) * absa(ind1(lay)-1,ig) + &
1607 126999408 : fac101(lay) * absa(ind1(lay),ig) + &
1608 126999408 : fac001(lay) * absa(ind1(lay)+1,ig) + &
1609 126999408 : fac211(lay) * absa(ind1(lay)+8,ig) + &
1610 126999408 : fac111(lay) * absa(ind1(lay)+9,ig) + &
1611 698496744 : fac011(lay) * absa(ind1(lay)+10,ig))
1612 : else
1613 72414012 : tau_major1 = speccomb1(lay) * &
1614 289656048 : (fac001(lay) * absa(ind1(lay),ig) + &
1615 289656048 : fac101(lay) * absa(ind1(lay)+1,ig) + &
1616 289656048 : fac011(lay) * absa(ind1(lay)+9,ig) + &
1617 1013796168 : fac111(lay) * absa(ind1(lay)+10,ig))
1618 : endif
1619 :
1620 209615424 : taug(lay,ngs6+ig) = tau_major + tau_major1 &
1621 : + tauself + taufor &
1622 314423136 : + adjcolco2(lay)*absco2
1623 628846272 : fracs(lay,ngs6+ig) = fracrefa(ig,jpl(lay)) + fpl(lay) * &
1624 742387960 : (fracrefa(ig,jpl(lay)+1)-fracrefa(ig,jpl(lay)))
1625 : enddo
1626 : enddo
1627 :
1628 : ! Upper atmosphere loop
1629 4874024 : do lay = laytrop+1, nlayers
1630 :
1631 : ! In atmospheres where the amount of CO2 is too great to be considered
1632 : ! a minor species, adjust the column amount of CO2 by an empirical factor
1633 : ! to obtain the proper contribution.
1634 4388024 : chi_co2(lay) = colco2(lay)/(coldry(lay))
1635 4388024 : ratco2(lay) = 1.e20*chi_co2(lay)/chi_mls(2,jp(lay)+1)
1636 4388024 : if (ratco2(lay) .gt. 3.0_r8) then
1637 0 : adjfac(lay) = 2.0_r8+(ratco2(lay)-2.0_r8)**0.79_r8
1638 0 : adjcolco2(lay) = adjfac(lay)*chi_mls(2,jp(lay)+1)*coldry(lay)*1.e-20_r8
1639 : else
1640 4388024 : adjcolco2(lay) = colco2(lay)
1641 : endif
1642 :
1643 4388024 : ind0(lay) = ((jp(lay)-13)*5+(jt(lay)-1))*nspb(7) + 1
1644 4388024 : ind1(lay) = ((jp(lay)-12)*5+(jt1(lay)-1))*nspb(7) + 1
1645 :
1646 57044312 : do ig = 1, ng7
1647 210625152 : absco2 = kb_mco2(indminor(lay),ig) + minorfrac(lay) * &
1648 263281440 : (kb_mco2(indminor(lay)+1,ig) - kb_mco2(indminor(lay),ig))
1649 157968864 : taug(lay,ngs6+ig) = colo3(lay) * &
1650 210625152 : (fac00(lay) * absb(ind0(lay),ig) + &
1651 210625152 : fac10(lay) * absb(ind0(lay)+1,ig) + &
1652 210625152 : fac01(lay) * absb(ind1(lay),ig) + &
1653 210625152 : fac11(lay) * absb(ind1(lay)+1,ig)) &
1654 1053125760 : + adjcolco2(lay) * absco2
1655 57044312 : fracs(lay,ngs6+ig) = fracrefb(ig)
1656 : enddo
1657 :
1658 : ! Empirical modification to code to improve stratospheric cooling rates
1659 : ! for o3. Revised to apply weighting for g-point reduction in this band.
1660 :
1661 4388024 : taug(lay,ngs6+6)=taug(lay,ngs6+6)*0.92_r8
1662 4388024 : taug(lay,ngs6+7)=taug(lay,ngs6+7)*0.88_r8
1663 4388024 : taug(lay,ngs6+8)=taug(lay,ngs6+8)*1.07_r8
1664 4388024 : taug(lay,ngs6+9)=taug(lay,ngs6+9)*1.1_r8
1665 4388024 : taug(lay,ngs6+10)=taug(lay,ngs6+10)*0.99_r8
1666 4874024 : taug(lay,ngs6+11)=taug(lay,ngs6+11)*0.855_r8
1667 :
1668 : enddo
1669 :
1670 486000 : end subroutine taugb7
1671 :
1672 : !----------------------------------------------------------------------------
1673 486000 : subroutine taugb8
1674 : !----------------------------------------------------------------------------
1675 : !
1676 : ! band 8: 1080-1180 cm-1 (low key - h2o; low minor - co2,o3,n2o)
1677 : ! (high key - o3; high minor - co2, n2o)
1678 : !----------------------------------------------------------------------------
1679 :
1680 : ! ------- Modules -------
1681 :
1682 : use parrrtm, only : ng8, ngs7
1683 : use rrlw_ref, only : chi_mls
1684 : use rrlw_kg08, only : fracrefa, fracrefb, absa, absb, &
1685 : ka_mco2, ka_mn2o, ka_mo3, kb_mco2, kb_mn2o, &
1686 : selfref, forref, cfc12, cfc22adj
1687 :
1688 : ! ------- Declarations -------
1689 :
1690 : ! Local
1691 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
1692 : real(kind=r8) :: tauself, taufor, absco2, abso3, absn2o
1693 972000 : real(kind=r8), dimension(nlayers) :: chi_co2, ratco2, adjfac, adjcolco2
1694 :
1695 :
1696 : ! Minor gas mapping level:
1697 : ! lower - co2, p = 1053.63 mb, t = 294.2 k
1698 : ! lower - o3, p = 317.348 mb, t = 240.77 k
1699 : ! lower - n2o, p = 706.2720 mb, t= 278.94 k
1700 : ! lower - cfc12,cfc11
1701 : ! upper - co2, p = 35.1632 mb, t = 223.28 k
1702 : ! upper - n2o, p = 8.716e-2 mb, t = 226.03 k
1703 :
1704 : ! Compute the optical depth by interpolating in ln(pressure) and
1705 : ! temperature, and appropriate species. Below laytrop, the water vapor
1706 : ! self-continuum and foreign continuum is interpolated (in temperature)
1707 : ! separately.
1708 :
1709 : ! Lower atmosphere loop
1710 9219976 : do lay = 1, laytrop
1711 :
1712 : ! In atmospheres where the amount of CO2 is too great to be considered
1713 : ! a minor species, adjust the column amount of CO2 by an empirical factor
1714 : ! to obtain the proper contribution.
1715 8733976 : chi_co2(lay) = colco2(lay)/(coldry(lay))
1716 8733976 : ratco2(lay) = 1.e20_r8*chi_co2(lay)/chi_mls(2,jp(lay)+1)
1717 8733976 : if (ratco2(lay) .gt. 3.0_r8) then
1718 0 : adjfac(lay) = 2.0_r8+(ratco2(lay)-2.0_r8)**0.65_r8
1719 0 : adjcolco2(lay) = adjfac(lay)*chi_mls(2,jp(lay)+1)*coldry(lay)*1.e-20_r8
1720 : else
1721 8733976 : adjcolco2(lay) = colco2(lay)
1722 : endif
1723 :
1724 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(8) + 1
1725 8733976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(8) + 1
1726 :
1727 79091784 : do ig = 1, ng8
1728 349359040 : tauself = selffac(lay) * (selfref(indself(lay),ig) + selffrac(lay) * &
1729 419230848 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
1730 349359040 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
1731 419230848 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
1732 279487232 : absco2 = (ka_mco2(indminor(lay),ig) + minorfrac(lay) * &
1733 349359040 : (ka_mco2(indminor(lay)+1,ig) - ka_mco2(indminor(lay),ig)))
1734 279487232 : abso3 = (ka_mo3(indminor(lay),ig) + minorfrac(lay) * &
1735 349359040 : (ka_mo3(indminor(lay)+1,ig) - ka_mo3(indminor(lay),ig)))
1736 279487232 : absn2o = (ka_mn2o(indminor(lay),ig) + minorfrac(lay) * &
1737 349359040 : (ka_mn2o(indminor(lay)+1,ig) - ka_mn2o(indminor(lay),ig)))
1738 209615424 : taug(lay,ngs7+ig) = colh2o(lay) * &
1739 279487232 : (fac00(lay) * absa(ind0(lay),ig) + &
1740 279487232 : fac10(lay) * absa(ind0(lay)+1,ig) + &
1741 279487232 : fac01(lay) * absa(ind1(lay),ig) + &
1742 279487232 : fac11(lay) * absa(ind1(lay)+1,ig)) &
1743 : + tauself + taufor &
1744 69871808 : + adjcolco2(lay)*absco2 &
1745 69871808 : + colo3(lay) * abso3 &
1746 69871808 : + coln2o(lay) * absn2o &
1747 139743616 : + wx(3,lay) * cfc12(ig) &
1748 1746795200 : + wx(4,lay) * cfc22adj(ig)
1749 78605784 : fracs(lay,ngs7+ig) = fracrefa(ig)
1750 : enddo
1751 : enddo
1752 :
1753 : ! Upper atmosphere loop
1754 4874024 : do lay = laytrop+1, nlayers
1755 :
1756 : ! In atmospheres where the amount of CO2 is too great to be considered
1757 : ! a minor species, adjust the column amount of CO2 by an empirical factor
1758 : ! to obtain the proper contribution.
1759 4388024 : chi_co2(lay) = colco2(lay)/coldry(lay)
1760 4388024 : ratco2(lay) = 1.e20_r8*chi_co2(lay)/chi_mls(2,jp(lay)+1)
1761 4388024 : if (ratco2(lay) .gt. 3.0_r8) then
1762 0 : adjfac(lay) = 2.0_r8+(ratco2(lay)-2.0_r8)**0.65_r8
1763 0 : adjcolco2(lay) = adjfac(lay)*chi_mls(2,jp(lay)+1) * coldry(lay)*1.e-20_r8
1764 : else
1765 4388024 : adjcolco2(lay) = colco2(lay)
1766 : endif
1767 :
1768 4388024 : ind0(lay) = ((jp(lay)-13)*5+(jt(lay)-1))*nspb(8) + 1
1769 4388024 : ind1(lay) = ((jp(lay)-12)*5+(jt1(lay)-1))*nspb(8) + 1
1770 :
1771 39978216 : do ig = 1, ng8
1772 140416768 : absco2 = (kb_mco2(indminor(lay),ig) + minorfrac(lay) * &
1773 175520960 : (kb_mco2(indminor(lay)+1,ig) - kb_mco2(indminor(lay),ig)))
1774 140416768 : absn2o = (kb_mn2o(indminor(lay),ig) + minorfrac(lay) * &
1775 175520960 : (kb_mn2o(indminor(lay)+1,ig) - kb_mn2o(indminor(lay),ig)))
1776 105312576 : taug(lay,ngs7+ig) = colo3(lay) * &
1777 140416768 : (fac00(lay) * absb(ind0(lay),ig) + &
1778 140416768 : fac10(lay) * absb(ind0(lay)+1,ig) + &
1779 140416768 : fac01(lay) * absb(ind1(lay),ig) + &
1780 140416768 : fac11(lay) * absb(ind1(lay)+1,ig)) &
1781 35104192 : + adjcolco2(lay)*absco2 &
1782 35104192 : + coln2o(lay)*absn2o &
1783 70208384 : + wx(3,lay) * cfc12(ig) &
1784 842500608 : + wx(4,lay) * cfc22adj(ig)
1785 39492216 : fracs(lay,ngs7+ig) = fracrefb(ig)
1786 : enddo
1787 : enddo
1788 :
1789 486000 : end subroutine taugb8
1790 :
1791 : !----------------------------------------------------------------------------
1792 486000 : subroutine taugb9
1793 : !----------------------------------------------------------------------------
1794 : !
1795 : ! band 9: 1180-1390 cm-1 (low key - h2o,ch4; low minor - n2o)
1796 : ! (high key - ch4; high minor - n2o)
1797 : !----------------------------------------------------------------------------
1798 :
1799 : ! ------- Modules -------
1800 :
1801 : use parrrtm, only : ng9, ngs8
1802 : use rrlw_ref, only : chi_mls
1803 : use rrlw_kg09, only : fracrefa, fracrefb, absa, absb, &
1804 : ka_mn2o, kb_mn2o, selfref, forref
1805 :
1806 : ! ------- Declarations -------
1807 :
1808 : ! Local
1809 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
1810 972000 : integer, dimension(nlayers) :: js, js1, jmn2o, jpl
1811 972000 : real(kind=r8), dimension(nlayers) :: speccomb, specparm, specmult, fs
1812 972000 : real(kind=r8), dimension(nlayers) :: speccomb1, specparm1, specmult1, fs1
1813 972000 : real(kind=r8), dimension(nlayers) :: speccomb_mn2o, specparm_mn2o, specmult_mn2o, fmn2o
1814 972000 : real(kind=r8), dimension(nlayers) :: speccomb_planck, specparm_planck, specmult_planck, fpl
1815 : real(kind=r8) :: p, p4, fk0, fk1, fk2
1816 972000 : real(kind=r8), dimension(nlayers) :: fac000, fac100, fac200, fac010, fac110, fac210
1817 972000 : real(kind=r8), dimension(nlayers) :: fac001, fac101, fac201, fac011, fac111, fac211
1818 : real(kind=r8) :: tauself, taufor, n2om1, n2om2, absn2o
1819 972000 : real(kind=r8), dimension(nlayers) :: chi_n2o, ratn2o, adjfac, adjcoln2o
1820 : real(kind=r8) :: refrat_planck_a, refrat_m_a
1821 : real(kind=r8) :: tau_major, tau_major1
1822 :
1823 :
1824 : ! Minor gas mapping level :
1825 : ! lower - n2o, p = 706.272 mbar, t = 278.94 k
1826 : ! upper - n2o, p = 95.58 mbar, t = 215.7 k
1827 :
1828 : ! Calculate reference ratio to be used in calculation of Planck
1829 : ! fraction in lower/upper atmosphere.
1830 :
1831 : ! P = 212 mb
1832 486000 : refrat_planck_a = chi_mls(1,9)/chi_mls(6,9)
1833 :
1834 : ! P = 706.272 mb
1835 486000 : refrat_m_a = chi_mls(1,3)/chi_mls(6,3)
1836 :
1837 : ! Compute the optical depth by interpolating in ln(pressure),
1838 : ! temperature, and appropriate species. Below laytrop, the water
1839 : ! vapor self-continuum and foreign continuum is interpolated
1840 : ! (in temperature) separately.
1841 :
1842 : ! Lower atmosphere loop
1843 9219976 : do lay = 1, laytrop
1844 :
1845 8733976 : speccomb(lay) = colh2o(lay) + rat_h2och4(lay)*colch4(lay)
1846 8733976 : specparm(lay) = colh2o(lay)/speccomb(lay)
1847 8733976 : if (specparm(lay) .ge. oneminus) specparm(lay) = oneminus
1848 8733976 : specmult(lay) = 8._r8*(specparm(lay))
1849 8733976 : js(lay) = 1 + int(specmult(lay))
1850 8733976 : fs(lay) = mod(specmult(lay),1.0_r8)
1851 :
1852 8733976 : speccomb1(lay) = colh2o(lay) + rat_h2och4_1(lay)*colch4(lay)
1853 8733976 : specparm1(lay) = colh2o(lay)/speccomb1(lay)
1854 8733976 : if (specparm1(lay) .ge. oneminus) specparm1(lay) = oneminus
1855 8733976 : specmult1(lay) = 8._r8*(specparm1(lay))
1856 8733976 : js1(lay) = 1 + int(specmult1(lay))
1857 8733976 : fs1(lay) = mod(specmult1(lay),1.0_r8)
1858 :
1859 8733976 : speccomb_mn2o(lay) = colh2o(lay) + refrat_m_a*colch4(lay)
1860 8733976 : specparm_mn2o(lay) = colh2o(lay)/speccomb_mn2o(lay)
1861 8733976 : if (specparm_mn2o(lay) .ge. oneminus) specparm_mn2o(lay) = oneminus
1862 8733976 : specmult_mn2o(lay) = 8._r8*specparm_mn2o(lay)
1863 8733976 : jmn2o(lay) = 1 + int(specmult_mn2o(lay))
1864 8733976 : fmn2o(lay) = mod(specmult_mn2o(lay),1.0_r8)
1865 :
1866 : ! In atmospheres where the amount of N2O is too great to be considered
1867 : ! a minor species, adjust the column amount of N2O by an empirical factor
1868 : ! to obtain the proper contribution.
1869 8733976 : chi_n2o(lay) = coln2o(lay)/(coldry(lay))
1870 8733976 : ratn2o(lay) = 1.e20_r8*chi_n2o(lay)/chi_mls(4,jp(lay)+1)
1871 8733976 : if (ratn2o(lay) .gt. 1.5_r8) then
1872 0 : adjfac(lay) = 0.5_r8+(ratn2o(lay)-0.5_r8)**0.65_r8
1873 0 : adjcoln2o(lay) = adjfac(lay)*chi_mls(4,jp(lay)+1)*coldry(lay)*1.e-20_r8
1874 : else
1875 8733976 : adjcoln2o(lay) = coln2o(lay)
1876 : endif
1877 :
1878 8733976 : speccomb_planck(lay) = colh2o(lay)+refrat_planck_a*colch4(lay)
1879 8733976 : specparm_planck(lay) = colh2o(lay)/speccomb_planck(lay)
1880 8733976 : if (specparm_planck(lay) .ge. oneminus) specparm_planck(lay)=oneminus
1881 8733976 : specmult_planck(lay) = 8._r8*specparm_planck(lay)
1882 8733976 : jpl(lay)= 1 + int(specmult_planck(lay))
1883 8733976 : fpl(lay) = mod(specmult_planck(lay),1.0_r8)
1884 :
1885 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(9) + js(lay)
1886 8733976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(9) + js1(lay)
1887 :
1888 8733976 : if (specparm(lay) .lt. 0.125_r8) then
1889 292858 : p = fs(lay) - 1
1890 292858 : p4 = p**4
1891 292858 : fk0 = p4
1892 292858 : fk1 = 1 - p - 2.0_r8*p4
1893 292858 : fk2 = p + p4
1894 292858 : fac000(lay) = fk0*fac00(lay)
1895 292858 : fac100(lay) = fk1*fac00(lay)
1896 292858 : fac200(lay) = fk2*fac00(lay)
1897 292858 : fac010(lay) = fk0*fac10(lay)
1898 292858 : fac110(lay) = fk1*fac10(lay)
1899 292858 : fac210(lay) = fk2*fac10(lay)
1900 8441118 : else if (specparm(lay) .gt. 0.875_r8) then
1901 115205 : p = -fs(lay)
1902 115205 : p4 = p**4
1903 115205 : fk0 = p4
1904 115205 : fk1 = 1 - p - 2.0_r8*p4
1905 115205 : fk2 = p + p4
1906 115205 : fac000(lay) = fk0*fac00(lay)
1907 115205 : fac100(lay) = fk1*fac00(lay)
1908 115205 : fac200(lay) = fk2*fac00(lay)
1909 115205 : fac010(lay) = fk0*fac10(lay)
1910 115205 : fac110(lay) = fk1*fac10(lay)
1911 115205 : fac210(lay) = fk2*fac10(lay)
1912 : else
1913 8325913 : fac000(lay) = (1._r8 - fs(lay)) * fac00(lay)
1914 8325913 : fac010(lay) = (1._r8 - fs(lay)) * fac10(lay)
1915 8325913 : fac100(lay) = fs(lay) * fac00(lay)
1916 8325913 : fac110(lay) = fs(lay) * fac10(lay)
1917 : endif
1918 :
1919 8733976 : if (specparm1(lay) .lt. 0.125_r8) then
1920 26723 : p = fs1(lay) - 1
1921 26723 : p4 = p**4
1922 26723 : fk0 = p4
1923 26723 : fk1 = 1 - p - 2.0_r8*p4
1924 26723 : fk2 = p + p4
1925 26723 : fac001(lay) = fk0*fac01(lay)
1926 26723 : fac101(lay) = fk1*fac01(lay)
1927 26723 : fac201(lay) = fk2*fac01(lay)
1928 26723 : fac011(lay) = fk0*fac11(lay)
1929 26723 : fac111(lay) = fk1*fac11(lay)
1930 26723 : fac211(lay) = fk2*fac11(lay)
1931 8707253 : else if (specparm1(lay) .gt. 0.875_r8) then
1932 1046819 : p = -fs1(lay)
1933 1046819 : p4 = p**4
1934 1046819 : fk0 = p4
1935 1046819 : fk1 = 1 - p - 2.0_r8*p4
1936 1046819 : fk2 = p + p4
1937 1046819 : fac001(lay) = fk0*fac01(lay)
1938 1046819 : fac101(lay) = fk1*fac01(lay)
1939 1046819 : fac201(lay) = fk2*fac01(lay)
1940 1046819 : fac011(lay) = fk0*fac11(lay)
1941 1046819 : fac111(lay) = fk1*fac11(lay)
1942 1046819 : fac211(lay) = fk2*fac11(lay)
1943 : else
1944 7660434 : fac001(lay) = (1._r8 - fs1(lay)) * fac01(lay)
1945 7660434 : fac011(lay) = (1._r8 - fs1(lay)) * fac11(lay)
1946 7660434 : fac101(lay) = fs1(lay) * fac01(lay)
1947 7660434 : fac111(lay) = fs1(lay) * fac11(lay)
1948 : endif
1949 :
1950 114027688 : do ig = 1, ng9
1951 524038560 : tauself = selffac(lay)* (selfref(indself(lay),ig) + selffrac(lay) * &
1952 628846272 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
1953 524038560 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
1954 628846272 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
1955 628846272 : n2om1 = ka_mn2o(jmn2o(lay),indminor(lay),ig) + fmn2o(lay) * &
1956 733653984 : (ka_mn2o(jmn2o(lay)+1,indminor(lay),ig) - ka_mn2o(jmn2o(lay),indminor(lay),ig))
1957 628846272 : n2om2 = ka_mn2o(jmn2o(lay),indminor(lay)+1,ig) + fmn2o(lay) * &
1958 733653984 : (ka_mn2o(jmn2o(lay)+1,indminor(lay)+1,ig) - ka_mn2o(jmn2o(lay),indminor(lay)+1,ig))
1959 104807712 : absn2o = n2om1 + minorfrac(lay) * (n2om2 - n2om1)
1960 :
1961 104807712 : if (specparm(lay) .lt. 0.125_r8) then
1962 3514296 : tau_major = speccomb(lay) * &
1963 14057184 : (fac000(lay) * absa(ind0(lay),ig) + &
1964 14057184 : fac100(lay) * absa(ind0(lay)+1,ig) + &
1965 14057184 : fac200(lay) * absa(ind0(lay)+2,ig) + &
1966 14057184 : fac010(lay) * absa(ind0(lay)+9,ig) + &
1967 14057184 : fac110(lay) * absa(ind0(lay)+10,ig) + &
1968 77314512 : fac210(lay) * absa(ind0(lay)+11,ig))
1969 101293416 : else if (specparm(lay) .gt. 0.875_r8) then
1970 1382460 : tau_major = speccomb(lay) * &
1971 5529840 : (fac200(lay) * absa(ind0(lay)-1,ig) + &
1972 5529840 : fac100(lay) * absa(ind0(lay),ig) + &
1973 5529840 : fac000(lay) * absa(ind0(lay)+1,ig) + &
1974 5529840 : fac210(lay) * absa(ind0(lay)+8,ig) + &
1975 5529840 : fac110(lay) * absa(ind0(lay)+9,ig) + &
1976 30414120 : fac010(lay) * absa(ind0(lay)+10,ig))
1977 : else
1978 99910956 : tau_major = speccomb(lay) * &
1979 399643824 : (fac000(lay) * absa(ind0(lay),ig) + &
1980 399643824 : fac100(lay) * absa(ind0(lay)+1,ig) + &
1981 399643824 : fac010(lay) * absa(ind0(lay)+9,ig) + &
1982 1398753384 : fac110(lay) * absa(ind0(lay)+10,ig))
1983 : endif
1984 :
1985 104807712 : if (specparm1(lay) .lt. 0.125_r8) then
1986 320676 : tau_major1 = speccomb1(lay) * &
1987 1282704 : (fac001(lay) * absa(ind1(lay),ig) + &
1988 1282704 : fac101(lay) * absa(ind1(lay)+1,ig) + &
1989 1282704 : fac201(lay) * absa(ind1(lay)+2,ig) + &
1990 1282704 : fac011(lay) * absa(ind1(lay)+9,ig) + &
1991 1282704 : fac111(lay) * absa(ind1(lay)+10,ig) + &
1992 7054872 : fac211(lay) * absa(ind1(lay)+11,ig))
1993 104487036 : else if (specparm1(lay) .gt. 0.875_r8) then
1994 12561828 : tau_major1 = speccomb1(lay) * &
1995 50247312 : (fac201(lay) * absa(ind1(lay)-1,ig) + &
1996 50247312 : fac101(lay) * absa(ind1(lay),ig) + &
1997 50247312 : fac001(lay) * absa(ind1(lay)+1,ig) + &
1998 50247312 : fac211(lay) * absa(ind1(lay)+8,ig) + &
1999 50247312 : fac111(lay) * absa(ind1(lay)+9,ig) + &
2000 276360216 : fac011(lay) * absa(ind1(lay)+10,ig))
2001 : else
2002 91925208 : tau_major1 = speccomb1(lay) * &
2003 367700832 : (fac001(lay) * absa(ind1(lay),ig) + &
2004 367700832 : fac101(lay) * absa(ind1(lay)+1,ig) + &
2005 367700832 : fac011(lay) * absa(ind1(lay)+9,ig) + &
2006 1286952912 : fac111(lay) * absa(ind1(lay)+10,ig))
2007 : endif
2008 :
2009 209615424 : taug(lay,ngs8+ig) = tau_major + tau_major1 &
2010 : + tauself + taufor &
2011 314423136 : + adjcoln2o(lay)*absn2o
2012 628846272 : fracs(lay,ngs8+ig) = fracrefa(ig,jpl(lay)) + fpl(lay) * &
2013 742387960 : (fracrefa(ig,jpl(lay)+1)-fracrefa(ig,jpl(lay)))
2014 : enddo
2015 : enddo
2016 :
2017 : ! Upper atmosphere loop
2018 4874024 : do lay = laytrop+1, nlayers
2019 :
2020 : ! In atmospheres where the amount of N2O is too great to be considered
2021 : ! a minor species, adjust the column amount of N2O by an empirical factor
2022 : ! to obtain the proper contribution.
2023 4388024 : chi_n2o(lay) = coln2o(lay)/(coldry(lay))
2024 4388024 : ratn2o(lay) = 1.e20_r8*chi_n2o(lay)/chi_mls(4,jp(lay)+1)
2025 4388024 : if (ratn2o(lay) .gt. 1.5_r8) then
2026 2441520 : adjfac(lay) = 0.5_r8+(ratn2o(lay)-0.5_r8)**0.65_r8
2027 2441520 : adjcoln2o(lay) = adjfac(lay)*chi_mls(4,jp(lay)+1)*coldry(lay)*1.e-20_r8
2028 : else
2029 1946504 : adjcoln2o(lay) = coln2o(lay)
2030 : endif
2031 :
2032 4388024 : ind0(lay) = ((jp(lay)-13)*5+(jt(lay)-1))*nspb(9) + 1
2033 4388024 : ind1(lay) = ((jp(lay)-12)*5+(jt1(lay)-1))*nspb(9) + 1
2034 :
2035 57530312 : do ig = 1, ng9
2036 210625152 : absn2o = kb_mn2o(indminor(lay),ig) + minorfrac(lay) * &
2037 263281440 : (kb_mn2o(indminor(lay)+1,ig) - kb_mn2o(indminor(lay),ig))
2038 157968864 : taug(lay,ngs8+ig) = colch4(lay) * &
2039 210625152 : (fac00(lay) * absb(ind0(lay),ig) + &
2040 210625152 : fac10(lay) * absb(ind0(lay)+1,ig) + &
2041 210625152 : fac01(lay) * absb(ind1(lay),ig) + &
2042 210625152 : fac11(lay) * absb(ind1(lay)+1,ig)) &
2043 1053125760 : + adjcoln2o(lay)*absn2o
2044 57044312 : fracs(lay,ngs8+ig) = fracrefb(ig)
2045 : enddo
2046 : enddo
2047 :
2048 486000 : end subroutine taugb9
2049 :
2050 : !----------------------------------------------------------------------------
2051 486000 : subroutine taugb10
2052 : !----------------------------------------------------------------------------
2053 : !
2054 : ! band 10: 1390-1480 cm-1 (low key - h2o; high key - h2o)
2055 : !----------------------------------------------------------------------------
2056 :
2057 : ! ------- Modules -------
2058 :
2059 : use parrrtm, only : ng10, ngs9
2060 : use rrlw_kg10, only : fracrefa, fracrefb, absa, absb, &
2061 : selfref, forref
2062 :
2063 : ! ------- Declarations -------
2064 :
2065 : ! Local
2066 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
2067 : real(kind=r8) :: tauself, taufor
2068 :
2069 :
2070 : ! Compute the optical depth by interpolating in ln(pressure) and
2071 : ! temperature. Below laytrop, the water vapor self-continuum and
2072 : ! foreign continuum is interpolated (in temperature) separately.
2073 :
2074 : ! Lower atmosphere loop
2075 9219976 : do lay = 1, laytrop
2076 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(10) + 1
2077 8733976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(10) + 1
2078 :
2079 61623832 : do ig = 1, ng10
2080 262019280 : tauself = selffac(lay) * (selfref(indself(lay),ig) + selffrac(lay) * &
2081 314423136 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
2082 262019280 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
2083 314423136 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
2084 157211568 : taug(lay,ngs9+ig) = colh2o(lay) * &
2085 209615424 : (fac00(lay) * absa(ind0(lay),ig) + &
2086 209615424 : fac10(lay) * absa(ind0(lay)+1,ig) + &
2087 209615424 : fac01(lay) * absa(ind1(lay),ig) + &
2088 209615424 : fac11(lay) * absa(ind1(lay)+1,ig)) &
2089 995673264 : + tauself + taufor
2090 61137832 : fracs(lay,ngs9+ig) = fracrefa(ig)
2091 : enddo
2092 : enddo
2093 :
2094 : ! Upper atmosphere loop
2095 4874024 : do lay = laytrop+1, nlayers
2096 4388024 : ind0(lay) = ((jp(lay)-13)*5+(jt(lay)-1))*nspb(10) + 1
2097 4388024 : ind1(lay) = ((jp(lay)-12)*5+(jt1(lay)-1))*nspb(10) + 1
2098 :
2099 31202168 : do ig = 1, ng10
2100 131640720 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
2101 157968864 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
2102 78984432 : taug(lay,ngs9+ig) = colh2o(lay) * &
2103 105312576 : (fac00(lay) * absb(ind0(lay),ig) + &
2104 105312576 : fac10(lay) * absb(ind0(lay)+1,ig) + &
2105 105312576 : fac01(lay) * absb(ind1(lay),ig) + &
2106 105312576 : fac11(lay) * absb(ind1(lay)+1,ig)) &
2107 500234736 : + taufor
2108 30716168 : fracs(lay,ngs9+ig) = fracrefb(ig)
2109 : enddo
2110 : enddo
2111 :
2112 486000 : end subroutine taugb10
2113 :
2114 : !----------------------------------------------------------------------------
2115 486000 : subroutine taugb11
2116 : !----------------------------------------------------------------------------
2117 : !
2118 : ! band 11: 1480-1800 cm-1 (low - h2o; low minor - o2)
2119 : ! (high key - h2o; high minor - o2)
2120 : !----------------------------------------------------------------------------
2121 :
2122 : ! ------- Modules -------
2123 :
2124 : use parrrtm, only : ng11, ngs10
2125 : use rrlw_kg11, only : fracrefa, fracrefb, absa, absb, &
2126 : ka_mo2, kb_mo2, selfref, forref
2127 :
2128 : ! ------- Declarations -------
2129 :
2130 : ! Local
2131 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
2132 972000 : real(kind=r8) :: scaleo2(nlayers), tauself, taufor, tauo2
2133 :
2134 :
2135 : ! Minor gas mapping level :
2136 : ! lower - o2, p = 706.2720 mbar, t = 278.94 k
2137 : ! upper - o2, p = 4.758820 mbarm t = 250.85 k
2138 :
2139 : ! Compute the optical depth by interpolating in ln(pressure) and
2140 : ! temperature. Below laytrop, the water vapor self-continuum and
2141 : ! foreign continuum is interpolated (in temperature) separately.
2142 :
2143 : ! Lower atmosphere loop
2144 9219976 : do lay = 1, laytrop
2145 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(11) + 1
2146 8733976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(11) + 1
2147 8733976 : scaleo2(lay) = colo2(lay)*scaleminor(lay)
2148 79091784 : do ig = 1, ng11
2149 349359040 : tauself = selffac(lay) * (selfref(indself(lay),ig) + selffrac(lay) * &
2150 419230848 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
2151 349359040 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
2152 419230848 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
2153 349359040 : tauo2 = scaleo2(lay) * (ka_mo2(indminor(lay),ig) + minorfrac(lay) * &
2154 419230848 : (ka_mo2(indminor(lay)+1,ig) - ka_mo2(indminor(lay),ig)))
2155 209615424 : taug(lay,ngs10+ig) = colh2o(lay) * &
2156 279487232 : (fac00(lay) * absa(ind0(lay),ig) + &
2157 279487232 : fac10(lay) * absa(ind0(lay)+1,ig) + &
2158 279487232 : fac01(lay) * absa(ind1(lay),ig) + &
2159 279487232 : fac11(lay) * absa(ind1(lay)+1,ig)) &
2160 : + tauself + taufor &
2161 1327564352 : + tauo2
2162 78605784 : fracs(lay,ngs10+ig) = fracrefa(ig)
2163 : enddo
2164 : enddo
2165 :
2166 : ! Upper atmosphere loop
2167 4874024 : do lay = laytrop+1, nlayers
2168 4388024 : ind0(lay) = ((jp(lay)-13)*5+(jt(lay)-1))*nspb(11) + 1
2169 4388024 : ind1(lay) = ((jp(lay)-12)*5+(jt1(lay)-1))*nspb(11) + 1
2170 4388024 : scaleo2(lay) = colo2(lay)*scaleminor(lay)
2171 39978216 : do ig = 1, ng11
2172 175520960 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
2173 210625152 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
2174 175520960 : tauo2 = scaleo2(lay) * (kb_mo2(indminor(lay),ig) + minorfrac(lay) * &
2175 210625152 : (kb_mo2(indminor(lay)+1,ig) - kb_mo2(indminor(lay),ig)))
2176 105312576 : taug(lay,ngs10+ig) = colh2o(lay) * &
2177 140416768 : (fac00(lay) * absb(ind0(lay),ig) + &
2178 140416768 : fac10(lay) * absb(ind0(lay)+1,ig) + &
2179 140416768 : fac01(lay) * absb(ind1(lay),ig) + &
2180 140416768 : fac11(lay) * absb(ind1(lay)+1,ig)) &
2181 : + taufor &
2182 666979648 : + tauo2
2183 39492216 : fracs(lay,ngs10+ig) = fracrefb(ig)
2184 : enddo
2185 : enddo
2186 :
2187 486000 : end subroutine taugb11
2188 :
2189 : !----------------------------------------------------------------------------
2190 486000 : subroutine taugb12
2191 : !----------------------------------------------------------------------------
2192 : !
2193 : ! band 12: 1800-2080 cm-1 (low - h2o,co2; high - nothing)
2194 : !----------------------------------------------------------------------------
2195 :
2196 : ! ------- Modules -------
2197 :
2198 : use parrrtm, only : ng12, ngs11
2199 : use rrlw_ref, only : chi_mls
2200 : use rrlw_kg12, only : fracrefa, absa, &
2201 : selfref, forref
2202 :
2203 : ! ------- Declarations -------
2204 :
2205 : ! Local
2206 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
2207 972000 : integer, dimension(nlayers) :: js, js1, jpl
2208 972000 : real(kind=r8), dimension(nlayers) :: speccomb, specparm, specmult, fs
2209 972000 : real(kind=r8), dimension(nlayers) :: speccomb1, specparm1, specmult1, fs1
2210 972000 : real(kind=r8), dimension(nlayers) :: speccomb_planck, specparm_planck, specmult_planck, fpl
2211 : real(kind=r8) :: p, p4, fk0, fk1, fk2
2212 972000 : real(kind=r8), dimension(nlayers) :: fac000, fac100, fac200, fac010, fac110, fac210
2213 972000 : real(kind=r8), dimension(nlayers) :: fac001, fac101, fac201, fac011, fac111, fac211
2214 : real(kind=r8) :: tauself, taufor
2215 : real(kind=r8) :: refrat_planck_a
2216 : real(kind=r8) :: tau_major, tau_major1
2217 :
2218 :
2219 : ! Calculate reference ratio to be used in calculation of Planck
2220 : ! fraction in lower/upper atmosphere.
2221 :
2222 : ! P = 174.164 mb
2223 486000 : refrat_planck_a = chi_mls(1,10)/chi_mls(2,10)
2224 :
2225 : ! Compute the optical depth by interpolating in ln(pressure),
2226 : ! temperature, and appropriate species. Below laytrop, the water
2227 : ! vapor self-continuum adn foreign continuum is interpolated
2228 : ! (in temperature) separately.
2229 :
2230 : ! Lower atmosphere loop
2231 9219976 : do lay = 1, laytrop
2232 :
2233 8733976 : speccomb(lay) = colh2o(lay) + rat_h2oco2(lay)*colco2(lay)
2234 8733976 : specparm(lay) = colh2o(lay)/speccomb(lay)
2235 8733976 : if (specparm(lay) .ge. oneminus) specparm(lay) = oneminus
2236 8733976 : specmult(lay) = 8._r8*(specparm(lay))
2237 8733976 : js(lay) = 1 + int(specmult(lay))
2238 8733976 : fs(lay) = mod(specmult(lay),1.0_r8)
2239 :
2240 8733976 : speccomb1(lay) = colh2o(lay) + rat_h2oco2_1(lay)*colco2(lay)
2241 8733976 : specparm1(lay) = colh2o(lay)/speccomb1(lay)
2242 8733976 : if (specparm1(lay) .ge. oneminus) specparm1(lay) = oneminus
2243 8733976 : specmult1(lay) = 8._r8*(specparm1(lay))
2244 8733976 : js1(lay) = 1 + int(specmult1(lay))
2245 8733976 : fs1(lay) = mod(specmult1(lay),1.0_r8)
2246 :
2247 8733976 : speccomb_planck(lay) = colh2o(lay)+refrat_planck_a*colco2(lay)
2248 8733976 : specparm_planck(lay) = colh2o(lay)/speccomb_planck(lay)
2249 8733976 : if (specparm_planck(lay) .ge. oneminus) specparm_planck(lay)=oneminus
2250 8733976 : specmult_planck(lay) = 8._r8*specparm_planck(lay)
2251 8733976 : jpl(lay)= 1 + int(specmult_planck(lay))
2252 8733976 : fpl(lay) = mod(specmult_planck(lay),1.0_r8)
2253 :
2254 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(12) + js(lay)
2255 9219976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(12) + js1(lay)
2256 : enddo
2257 9219976 : do lay = 1, laytrop
2258 9219976 : if (specparm(lay) .lt. 0.125_r8) then
2259 770299 : p = fs(lay) - 1
2260 770299 : p4 = p**4
2261 770299 : fk0 = p4
2262 770299 : fk1 = 1 - p - 2.0_r8*p4
2263 770299 : fk2 = p + p4
2264 770299 : fac000(lay) = fk0*fac00(lay)
2265 770299 : fac100(lay) = fk1*fac00(lay)
2266 770299 : fac200(lay) = fk2*fac00(lay)
2267 770299 : fac010(lay) = fk0*fac10(lay)
2268 770299 : fac110(lay) = fk1*fac10(lay)
2269 770299 : fac210(lay) = fk2*fac10(lay)
2270 7963677 : else if (specparm(lay) .gt. 0.875_r8) then
2271 17165 : p = -fs(lay)
2272 17165 : p4 = p**4
2273 17165 : fk0 = p4
2274 17165 : fk1 = 1 - p - 2.0_r8*p4
2275 17165 : fk2 = p + p4
2276 17165 : fac000(lay) = fk0*fac00(lay)
2277 17165 : fac100(lay) = fk1*fac00(lay)
2278 17165 : fac200(lay) = fk2*fac00(lay)
2279 17165 : fac010(lay) = fk0*fac10(lay)
2280 17165 : fac110(lay) = fk1*fac10(lay)
2281 17165 : fac210(lay) = fk2*fac10(lay)
2282 : else
2283 7946512 : fac000(lay) = (1._r8 - fs(lay)) * fac00(lay)
2284 7946512 : fac010(lay) = (1._r8 - fs(lay)) * fac10(lay)
2285 7946512 : fac100(lay) = fs(lay) * fac00(lay)
2286 7946512 : fac110(lay) = fs(lay) * fac10(lay)
2287 : endif
2288 :
2289 : enddo
2290 9219976 : do lay = 1, laytrop
2291 9219976 : if (specparm1(lay) .lt. 0.125_r8) then
2292 152313 : p = fs1(lay) - 1
2293 152313 : p4 = p**4
2294 152313 : fk0 = p4
2295 152313 : fk1 = 1 - p - 2.0_r8*p4
2296 152313 : fk2 = p + p4
2297 152313 : fac001(lay) = fk0*fac01(lay)
2298 152313 : fac101(lay) = fk1*fac01(lay)
2299 152313 : fac201(lay) = fk2*fac01(lay)
2300 152313 : fac011(lay) = fk0*fac11(lay)
2301 152313 : fac111(lay) = fk1*fac11(lay)
2302 152313 : fac211(lay) = fk2*fac11(lay)
2303 8581663 : else if (specparm1(lay) .gt. 0.875_r8) then
2304 351185 : p = -fs1(lay)
2305 351185 : p4 = p**4
2306 351185 : fk0 = p4
2307 351185 : fk1 = 1 - p - 2.0_r8*p4
2308 351185 : fk2 = p + p4
2309 351185 : fac001(lay) = fk0*fac01(lay)
2310 351185 : fac101(lay) = fk1*fac01(lay)
2311 351185 : fac201(lay) = fk2*fac01(lay)
2312 351185 : fac011(lay) = fk0*fac11(lay)
2313 351185 : fac111(lay) = fk1*fac11(lay)
2314 351185 : fac211(lay) = fk2*fac11(lay)
2315 : else
2316 8230478 : fac001(lay) = (1._r8 - fs1(lay)) * fac01(lay)
2317 8230478 : fac011(lay) = (1._r8 - fs1(lay)) * fac11(lay)
2318 8230478 : fac101(lay) = fs1(lay) * fac01(lay)
2319 8230478 : fac111(lay) = fs1(lay) * fac11(lay)
2320 : endif
2321 : enddo
2322 9219976 : do lay = 1, laytrop
2323 :
2324 79091784 : do ig = 1, ng12
2325 349359040 : tauself = selffac(lay)* (selfref(indself(lay),ig) + selffrac(lay) * &
2326 419230848 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
2327 349359040 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
2328 419230848 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
2329 :
2330 69871808 : if (specparm(lay) .lt. 0.125_r8) then
2331 6162392 : tau_major = speccomb(lay) * &
2332 24649568 : (fac000(lay) * absa(ind0(lay),ig) + &
2333 24649568 : fac100(lay) * absa(ind0(lay)+1,ig) + &
2334 24649568 : fac200(lay) * absa(ind0(lay)+2,ig) + &
2335 24649568 : fac010(lay) * absa(ind0(lay)+9,ig) + &
2336 24649568 : fac110(lay) * absa(ind0(lay)+10,ig) + &
2337 135572624 : fac210(lay) * absa(ind0(lay)+11,ig))
2338 63709416 : else if (specparm(lay) .gt. 0.875_r8) then
2339 137320 : tau_major = speccomb(lay) * &
2340 549280 : (fac200(lay) * absa(ind0(lay)-1,ig) + &
2341 549280 : fac100(lay) * absa(ind0(lay),ig) + &
2342 549280 : fac000(lay) * absa(ind0(lay)+1,ig) + &
2343 549280 : fac210(lay) * absa(ind0(lay)+8,ig) + &
2344 549280 : fac110(lay) * absa(ind0(lay)+9,ig) + &
2345 3021040 : fac010(lay) * absa(ind0(lay)+10,ig))
2346 : else
2347 63572096 : tau_major = speccomb(lay) * &
2348 254288384 : (fac000(lay) * absa(ind0(lay),ig) + &
2349 254288384 : fac100(lay) * absa(ind0(lay)+1,ig) + &
2350 254288384 : fac010(lay) * absa(ind0(lay)+9,ig) + &
2351 890009344 : fac110(lay) * absa(ind0(lay)+10,ig))
2352 : endif
2353 :
2354 69871808 : if (specparm1(lay) .lt. 0.125_r8) then
2355 1218504 : tau_major1 = speccomb1(lay) * &
2356 4874016 : (fac001(lay) * absa(ind1(lay),ig) + &
2357 4874016 : fac101(lay) * absa(ind1(lay)+1,ig) + &
2358 4874016 : fac201(lay) * absa(ind1(lay)+2,ig) + &
2359 4874016 : fac011(lay) * absa(ind1(lay)+9,ig) + &
2360 4874016 : fac111(lay) * absa(ind1(lay)+10,ig) + &
2361 26807088 : fac211(lay) * absa(ind1(lay)+11,ig))
2362 68653304 : else if (specparm1(lay) .gt. 0.875_r8) then
2363 2809480 : tau_major1 = speccomb1(lay) * &
2364 11237920 : (fac201(lay) * absa(ind1(lay)-1,ig) + &
2365 11237920 : fac101(lay) * absa(ind1(lay),ig) + &
2366 11237920 : fac001(lay) * absa(ind1(lay)+1,ig) + &
2367 11237920 : fac211(lay) * absa(ind1(lay)+8,ig) + &
2368 11237920 : fac111(lay) * absa(ind1(lay)+9,ig) + &
2369 61808560 : fac011(lay) * absa(ind1(lay)+10,ig))
2370 : else
2371 65843824 : tau_major1 = speccomb1(lay) * &
2372 263375296 : (fac001(lay) * absa(ind1(lay),ig) + &
2373 263375296 : fac101(lay) * absa(ind1(lay)+1,ig) + &
2374 263375296 : fac011(lay) * absa(ind1(lay)+9,ig) + &
2375 921813536 : fac111(lay) * absa(ind1(lay)+10,ig))
2376 : endif
2377 :
2378 139743616 : taug(lay,ngs11+ig) = tau_major + tau_major1 &
2379 139743616 : + tauself + taufor
2380 419230848 : fracs(lay,ngs11+ig) = fracrefa(ig,jpl(lay)) + fpl(lay) * &
2381 497836632 : (fracrefa(ig,jpl(lay)+1)-fracrefa(ig,jpl(lay)))
2382 : enddo
2383 : enddo
2384 :
2385 : ! Upper atmosphere loop
2386 4874024 : do lay = laytrop+1, nlayers
2387 39978216 : do ig = 1, ng12
2388 35104192 : taug(lay,ngs11+ig) = 0.0_r8
2389 39492216 : fracs(lay,ngs11+ig) = 0.0_r8
2390 : enddo
2391 : enddo
2392 :
2393 486000 : end subroutine taugb12
2394 :
2395 : !----------------------------------------------------------------------------
2396 486000 : subroutine taugb13
2397 : !----------------------------------------------------------------------------
2398 : !
2399 : ! band 13: 2080-2250 cm-1 (low key - h2o,n2o; high minor - o3 minor)
2400 : !----------------------------------------------------------------------------
2401 :
2402 : ! ------- Modules -------
2403 :
2404 : use parrrtm, only : ng13, ngs12
2405 : use rrlw_ref, only : chi_mls
2406 : use rrlw_kg13, only : fracrefa, fracrefb, absa, &
2407 : ka_mco2, ka_mco, kb_mo3, selfref, forref
2408 :
2409 : ! ------- Declarations -------
2410 :
2411 : ! Local
2412 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
2413 972000 : integer, dimension(nlayers) :: js, js1, jmco2, jmco, jpl
2414 972000 : real(kind=r8), dimension(nlayers) :: speccomb, specparm, specmult, fs
2415 972000 : real(kind=r8), dimension(nlayers) :: speccomb1, specparm1, specmult1, fs1
2416 972000 : real(kind=r8), dimension(nlayers) :: speccomb_mco2, specparm_mco2, specmult_mco2, fmco2
2417 972000 : real(kind=r8), dimension(nlayers) :: speccomb_mco, specparm_mco, specmult_mco, fmco
2418 972000 : real(kind=r8), dimension(nlayers) :: speccomb_planck, specparm_planck, specmult_planck, fpl
2419 : real(kind=r8) :: p, p4, fk0, fk1, fk2
2420 972000 : real(kind=r8), dimension(nlayers) :: fac000, fac100, fac200, fac010, fac110, fac210
2421 972000 : real(kind=r8), dimension(nlayers) :: fac001, fac101, fac201, fac011, fac111, fac211
2422 : real(kind=r8) :: tauself, taufor, co2m1, co2m2, absco2
2423 : real(kind=r8) :: com1, com2, absco, abso3
2424 972000 : real(kind=r8), dimension(nlayers) :: chi_co2, ratco2, adjfac, adjcolco2
2425 : real(kind=r8) :: refrat_planck_a, refrat_m_a, refrat_m_a3
2426 : real(kind=r8) :: tau_major, tau_major1
2427 :
2428 :
2429 : ! Minor gas mapping levels :
2430 : ! lower - co2, p = 1053.63 mb, t = 294.2 k
2431 : ! lower - co, p = 706 mb, t = 278.94 k
2432 : ! upper - o3, p = 95.5835 mb, t = 215.7 k
2433 :
2434 : ! Calculate reference ratio to be used in calculation of Planck
2435 : ! fraction in lower/upper atmosphere.
2436 :
2437 : ! P = 473.420 mb (Level 5)
2438 486000 : refrat_planck_a = chi_mls(1,5)/chi_mls(4,5)
2439 :
2440 : ! P = 1053. (Level 1)
2441 486000 : refrat_m_a = chi_mls(1,1)/chi_mls(4,1)
2442 :
2443 : ! P = 706. (Level 3)
2444 486000 : refrat_m_a3 = chi_mls(1,3)/chi_mls(4,3)
2445 :
2446 : ! Compute the optical depth by interpolating in ln(pressure),
2447 : ! temperature, and appropriate species. Below laytrop, the water
2448 : ! vapor self-continuum and foreign continuum is interpolated
2449 : ! (in temperature) separately.
2450 :
2451 : ! Lower atmosphere loop
2452 9219976 : do lay = 1, laytrop
2453 :
2454 8733976 : speccomb(lay) = colh2o(lay) + rat_h2on2o(lay)*coln2o(lay)
2455 8733976 : specparm(lay) = colh2o(lay)/speccomb(lay)
2456 8733976 : if (specparm(lay) .ge. oneminus) specparm(lay) = oneminus
2457 8733976 : specmult(lay) = 8._r8*(specparm(lay))
2458 8733976 : js(lay) = 1 + int(specmult(lay))
2459 8733976 : fs(lay) = mod(specmult(lay),1.0_r8)
2460 :
2461 8733976 : speccomb1(lay) = colh2o(lay) + rat_h2on2o_1(lay)*coln2o(lay)
2462 8733976 : specparm1(lay) = colh2o(lay)/speccomb1(lay)
2463 8733976 : if (specparm1(lay) .ge. oneminus) specparm1(lay) = oneminus
2464 8733976 : specmult1(lay) = 8._r8*(specparm1(lay))
2465 8733976 : js1(lay) = 1 + int(specmult1(lay))
2466 8733976 : fs1(lay) = mod(specmult1(lay),1.0_r8)
2467 :
2468 8733976 : speccomb_mco2(lay) = colh2o(lay) + refrat_m_a*coln2o(lay)
2469 8733976 : specparm_mco2(lay) = colh2o(lay)/speccomb_mco2(lay)
2470 8733976 : if (specparm_mco2(lay) .ge. oneminus) specparm_mco2(lay) = oneminus
2471 8733976 : specmult_mco2(lay) = 8._r8*specparm_mco2(lay)
2472 8733976 : jmco2(lay) = 1 + int(specmult_mco2(lay))
2473 8733976 : fmco2(lay) = mod(specmult_mco2(lay),1.0_r8)
2474 :
2475 : ! In atmospheres where the amount of CO2 is too great to be considered
2476 : ! a minor species, adjust the column amount of CO2 by an empirical factor
2477 : ! to obtain the proper contribution.
2478 8733976 : chi_co2(lay) = colco2(lay)/(coldry(lay))
2479 9219976 : ratco2(lay) = 1.e20_r8*chi_co2(lay)/3.55e-4_r8
2480 : enddo
2481 9219976 : do lay = 1, laytrop
2482 9219976 : if (ratco2(lay) .gt. 3.0_r8) then
2483 0 : adjfac(lay) = 2.0_r8+(ratco2(lay)-2.0_r8)**0.68_r8
2484 0 : adjcolco2(lay) = adjfac(lay)*3.55e-4*coldry(lay)*1.e-20_r8
2485 : else
2486 8733976 : adjcolco2(lay) = colco2(lay)
2487 : endif
2488 :
2489 : enddo
2490 9219976 : do lay = 1, laytrop
2491 8733976 : speccomb_mco(lay) = colh2o(lay) + refrat_m_a3*coln2o(lay)
2492 8733976 : specparm_mco(lay) = colh2o(lay)/speccomb_mco(lay)
2493 8733976 : if (specparm_mco(lay) .ge. oneminus) specparm_mco(lay) = oneminus
2494 8733976 : specmult_mco(lay) = 8._r8*specparm_mco(lay)
2495 8733976 : jmco(lay) = 1 + int(specmult_mco(lay))
2496 8733976 : fmco(lay) = mod(specmult_mco(lay),1.0_r8)
2497 :
2498 8733976 : speccomb_planck(lay) = colh2o(lay)+refrat_planck_a*coln2o(lay)
2499 8733976 : specparm_planck(lay) = colh2o(lay)/speccomb_planck(lay)
2500 8733976 : if (specparm_planck(lay) .ge. oneminus) specparm_planck(lay)=oneminus
2501 8733976 : specmult_planck(lay) = 8._r8*specparm_planck(lay)
2502 8733976 : jpl(lay)= 1 + int(specmult_planck(lay))
2503 8733976 : fpl(lay) = mod(specmult_planck(lay),1.0_r8)
2504 :
2505 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(13) + js(lay)
2506 9219976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(13) + js1(lay)
2507 : enddo
2508 9219976 : do lay = 1, laytrop
2509 :
2510 9219976 : if (specparm(lay) .lt. 0.125_r8) then
2511 858585 : p = fs(lay) - 1
2512 858585 : p4 = p**4
2513 858585 : fk0 = p4
2514 858585 : fk1 = 1 - p - 2.0_r8*p4
2515 858585 : fk2 = p + p4
2516 858585 : fac000(lay) = fk0*fac00(lay)
2517 858585 : fac100(lay) = fk1*fac00(lay)
2518 858585 : fac200(lay) = fk2*fac00(lay)
2519 858585 : fac010(lay) = fk0*fac10(lay)
2520 858585 : fac110(lay) = fk1*fac10(lay)
2521 858585 : fac210(lay) = fk2*fac10(lay)
2522 7875391 : else if (specparm(lay) .gt. 0.875_r8) then
2523 6119 : p = -fs(lay)
2524 6119 : p4 = p**4
2525 6119 : fk0 = p4
2526 6119 : fk1 = 1 - p - 2.0_r8*p4
2527 6119 : fk2 = p + p4
2528 6119 : fac000(lay) = fk0*fac00(lay)
2529 6119 : fac100(lay) = fk1*fac00(lay)
2530 6119 : fac200(lay) = fk2*fac00(lay)
2531 6119 : fac010(lay) = fk0*fac10(lay)
2532 6119 : fac110(lay) = fk1*fac10(lay)
2533 6119 : fac210(lay) = fk2*fac10(lay)
2534 : else
2535 7869272 : fac000(lay) = (1._r8 - fs(lay)) * fac00(lay)
2536 7869272 : fac010(lay) = (1._r8 - fs(lay)) * fac10(lay)
2537 7869272 : fac100(lay) = fs(lay) * fac00(lay)
2538 7869272 : fac110(lay) = fs(lay) * fac10(lay)
2539 : endif
2540 :
2541 : enddo
2542 9219976 : do lay = 1, laytrop
2543 9219976 : if (specparm1(lay) .lt. 0.125_r8) then
2544 197310 : p = fs1(lay) - 1
2545 197310 : p4 = p**4
2546 197310 : fk0 = p4
2547 197310 : fk1 = 1 - p - 2.0_r8*p4
2548 197310 : fk2 = p + p4
2549 197310 : fac001(lay) = fk0*fac01(lay)
2550 197310 : fac101(lay) = fk1*fac01(lay)
2551 197310 : fac201(lay) = fk2*fac01(lay)
2552 197310 : fac011(lay) = fk0*fac11(lay)
2553 197310 : fac111(lay) = fk1*fac11(lay)
2554 197310 : fac211(lay) = fk2*fac11(lay)
2555 8536666 : else if (specparm1(lay) .gt. 0.875_r8) then
2556 226059 : p = -fs1(lay)
2557 226059 : p4 = p**4
2558 226059 : fk0 = p4
2559 226059 : fk1 = 1 - p - 2.0_r8*p4
2560 226059 : fk2 = p + p4
2561 226059 : fac001(lay) = fk0*fac01(lay)
2562 226059 : fac101(lay) = fk1*fac01(lay)
2563 226059 : fac201(lay) = fk2*fac01(lay)
2564 226059 : fac011(lay) = fk0*fac11(lay)
2565 226059 : fac111(lay) = fk1*fac11(lay)
2566 226059 : fac211(lay) = fk2*fac11(lay)
2567 : else
2568 8310607 : fac001(lay) = (1._r8 - fs1(lay)) * fac01(lay)
2569 8310607 : fac011(lay) = (1._r8 - fs1(lay)) * fac11(lay)
2570 8310607 : fac101(lay) = fs1(lay) * fac01(lay)
2571 8310607 : fac111(lay) = fs1(lay) * fac11(lay)
2572 : endif
2573 : enddo
2574 9219976 : do lay = 1, laytrop
2575 :
2576 44155880 : do ig = 1, ng13
2577 174679520 : tauself = selffac(lay)* (selfref(indself(lay),ig) + selffrac(lay) * &
2578 209615424 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
2579 174679520 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
2580 209615424 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
2581 209615424 : co2m1 = ka_mco2(jmco2(lay),indminor(lay),ig) + fmco2(lay) * &
2582 244551328 : (ka_mco2(jmco2(lay)+1,indminor(lay),ig) - ka_mco2(jmco2(lay),indminor(lay),ig))
2583 209615424 : co2m2 = ka_mco2(jmco2(lay),indminor(lay)+1,ig) + fmco2(lay) * &
2584 244551328 : (ka_mco2(jmco2(lay)+1,indminor(lay)+1,ig) - ka_mco2(jmco2(lay),indminor(lay)+1,ig))
2585 34935904 : absco2 = co2m1 + minorfrac(lay) * (co2m2 - co2m1)
2586 209615424 : com1 = ka_mco(jmco(lay),indminor(lay),ig) + fmco(lay) * &
2587 244551328 : (ka_mco(jmco(lay)+1,indminor(lay),ig) - ka_mco(jmco(lay),indminor(lay),ig))
2588 209615424 : com2 = ka_mco(jmco(lay),indminor(lay)+1,ig) + fmco(lay) * &
2589 244551328 : (ka_mco(jmco(lay)+1,indminor(lay)+1,ig) - ka_mco(jmco(lay),indminor(lay)+1,ig))
2590 34935904 : absco = com1 + minorfrac(lay) * (com2 - com1)
2591 :
2592 34935904 : if (specparm(lay) .lt. 0.125_r8) then
2593 3434340 : tau_major = speccomb(lay) * &
2594 13737360 : (fac000(lay) * absa(ind0(lay),ig) + &
2595 13737360 : fac100(lay) * absa(ind0(lay)+1,ig) + &
2596 13737360 : fac200(lay) * absa(ind0(lay)+2,ig) + &
2597 13737360 : fac010(lay) * absa(ind0(lay)+9,ig) + &
2598 13737360 : fac110(lay) * absa(ind0(lay)+10,ig) + &
2599 75555480 : fac210(lay) * absa(ind0(lay)+11,ig))
2600 31501564 : else if (specparm(lay) .gt. 0.875_r8) then
2601 24476 : tau_major = speccomb(lay) * &
2602 97904 : (fac200(lay) * absa(ind0(lay)-1,ig) + &
2603 97904 : fac100(lay) * absa(ind0(lay),ig) + &
2604 97904 : fac000(lay) * absa(ind0(lay)+1,ig) + &
2605 97904 : fac210(lay) * absa(ind0(lay)+8,ig) + &
2606 97904 : fac110(lay) * absa(ind0(lay)+9,ig) + &
2607 538472 : fac010(lay) * absa(ind0(lay)+10,ig))
2608 : else
2609 31477088 : tau_major = speccomb(lay) * &
2610 125908352 : (fac000(lay) * absa(ind0(lay),ig) + &
2611 125908352 : fac100(lay) * absa(ind0(lay)+1,ig) + &
2612 125908352 : fac010(lay) * absa(ind0(lay)+9,ig) + &
2613 440679232 : fac110(lay) * absa(ind0(lay)+10,ig))
2614 : endif
2615 :
2616 34935904 : if (specparm1(lay) .lt. 0.125_r8) then
2617 789240 : tau_major1 = speccomb1(lay) * &
2618 3156960 : (fac001(lay) * absa(ind1(lay),ig) + &
2619 3156960 : fac101(lay) * absa(ind1(lay)+1,ig) + &
2620 3156960 : fac201(lay) * absa(ind1(lay)+2,ig) + &
2621 3156960 : fac011(lay) * absa(ind1(lay)+9,ig) + &
2622 3156960 : fac111(lay) * absa(ind1(lay)+10,ig) + &
2623 17363280 : fac211(lay) * absa(ind1(lay)+11,ig))
2624 34146664 : else if (specparm1(lay) .gt. 0.875_r8) then
2625 904236 : tau_major1 = speccomb1(lay) * &
2626 3616944 : (fac201(lay) * absa(ind1(lay)-1,ig) + &
2627 3616944 : fac101(lay) * absa(ind1(lay),ig) + &
2628 3616944 : fac001(lay) * absa(ind1(lay)+1,ig) + &
2629 3616944 : fac211(lay) * absa(ind1(lay)+8,ig) + &
2630 3616944 : fac111(lay) * absa(ind1(lay)+9,ig) + &
2631 19893192 : fac011(lay) * absa(ind1(lay)+10,ig))
2632 : else
2633 33242428 : tau_major1 = speccomb1(lay) * &
2634 132969712 : (fac001(lay) * absa(ind1(lay),ig) + &
2635 132969712 : fac101(lay) * absa(ind1(lay)+1,ig) + &
2636 132969712 : fac011(lay) * absa(ind1(lay)+9,ig) + &
2637 465393992 : fac111(lay) * absa(ind1(lay)+10,ig))
2638 : endif
2639 :
2640 69871808 : taug(lay,ngs12+ig) = tau_major + tau_major1 &
2641 : + tauself + taufor &
2642 34935904 : + adjcolco2(lay)*absco2 &
2643 139743616 : + colco(lay)*absco
2644 209615424 : fracs(lay,ngs12+ig) = fracrefa(ig,jpl(lay)) + fpl(lay) * &
2645 253285304 : (fracrefa(ig,jpl(lay)+1)-fracrefa(ig,jpl(lay)))
2646 : enddo
2647 : enddo
2648 :
2649 : ! Upper atmosphere loop
2650 4874024 : do lay = laytrop+1, nlayers
2651 22426120 : do ig = 1, ng13
2652 70208384 : abso3 = kb_mo3(indminor(lay),ig) + minorfrac(lay) * &
2653 87760480 : (kb_mo3(indminor(lay)+1,ig) - kb_mo3(indminor(lay),ig))
2654 17552096 : taug(lay,ngs12+ig) = colo3(lay)*abso3
2655 21940120 : fracs(lay,ngs12+ig) = fracrefb(ig)
2656 : enddo
2657 : enddo
2658 :
2659 486000 : end subroutine taugb13
2660 :
2661 : !----------------------------------------------------------------------------
2662 486000 : subroutine taugb14
2663 : !----------------------------------------------------------------------------
2664 : !
2665 : ! band 14: 2250-2380 cm-1 (low - co2; high - co2)
2666 : !----------------------------------------------------------------------------
2667 :
2668 : ! ------- Modules -------
2669 :
2670 : use parrrtm, only : ng14, ngs13
2671 : use rrlw_kg14, only : fracrefa, fracrefb, absa, absb, &
2672 : selfref, forref
2673 :
2674 : ! ------- Declarations -------
2675 :
2676 : ! Local
2677 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
2678 : real(kind=r8) :: tauself, taufor
2679 :
2680 :
2681 : ! Compute the optical depth by interpolating in ln(pressure) and
2682 : ! temperature. Below laytrop, the water vapor self-continuum
2683 : ! and foreign continuum is interpolated (in temperature) separately.
2684 :
2685 : ! Lower atmosphere loop
2686 9219976 : do lay = 1, laytrop
2687 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(14) + 1
2688 8733976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(14) + 1
2689 26687928 : do ig = 1, ng14
2690 87339760 : tauself = selffac(lay) * (selfref(indself(lay),ig) + selffrac(lay) * &
2691 104807712 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
2692 87339760 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
2693 104807712 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
2694 52403856 : taug(lay,ngs13+ig) = colco2(lay) * &
2695 69871808 : (fac00(lay) * absa(ind0(lay),ig) + &
2696 69871808 : fac10(lay) * absa(ind0(lay)+1,ig) + &
2697 69871808 : fac01(lay) * absa(ind1(lay),ig) + &
2698 69871808 : fac11(lay) * absa(ind1(lay)+1,ig)) &
2699 331891088 : + tauself + taufor
2700 26201928 : fracs(lay,ngs13+ig) = fracrefa(ig)
2701 : enddo
2702 : enddo
2703 :
2704 : ! Upper atmosphere loop
2705 4874024 : do lay = laytrop+1, nlayers
2706 4388024 : ind0(lay) = ((jp(lay)-13)*5+(jt(lay)-1))*nspb(14) + 1
2707 4388024 : ind1(lay) = ((jp(lay)-12)*5+(jt1(lay)-1))*nspb(14) + 1
2708 13650072 : do ig = 1, ng14
2709 26328144 : taug(lay,ngs13+ig) = colco2(lay) * &
2710 35104192 : (fac00(lay) * absb(ind0(lay),ig) + &
2711 35104192 : fac10(lay) * absb(ind0(lay)+1,ig) + &
2712 35104192 : fac01(lay) * absb(ind1(lay),ig) + &
2713 140416768 : fac11(lay) * absb(ind1(lay)+1,ig))
2714 13164072 : fracs(lay,ngs13+ig) = fracrefb(ig)
2715 : enddo
2716 : enddo
2717 :
2718 486000 : end subroutine taugb14
2719 :
2720 : !----------------------------------------------------------------------------
2721 486000 : subroutine taugb15
2722 : !----------------------------------------------------------------------------
2723 : !
2724 : ! band 15: 2380-2600 cm-1 (low - n2o,co2; low minor - n2)
2725 : ! (high - nothing)
2726 : !----------------------------------------------------------------------------
2727 :
2728 : ! ------- Modules -------
2729 :
2730 : use parrrtm, only : ng15, ngs14
2731 : use rrlw_ref, only : chi_mls
2732 : use rrlw_kg15, only : fracrefa, absa, &
2733 : ka_mn2, selfref, forref
2734 :
2735 : ! ------- Declarations -------
2736 :
2737 : ! Local
2738 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
2739 972000 : integer, dimension(nlayers) :: js, js1, jmn2, jpl
2740 972000 : real(kind=r8), dimension(nlayers) :: speccomb, specparm, specmult, fs
2741 972000 : real(kind=r8), dimension(nlayers) :: speccomb1, specparm1, specmult1, fs1
2742 972000 : real(kind=r8), dimension(nlayers) :: speccomb_mn2, specparm_mn2, specmult_mn2, fmn2
2743 972000 : real(kind=r8), dimension(nlayers) :: speccomb_planck, specparm_planck, specmult_planck, fpl
2744 : real(kind=r8) :: p, p4, fk0, fk1, fk2
2745 972000 : real(kind=r8), dimension(nlayers) :: fac000, fac100, fac200, fac010, fac110, fac210
2746 972000 : real(kind=r8), dimension(nlayers) :: fac001, fac101, fac201, fac011, fac111, fac211
2747 972000 : real(kind=r8) :: scalen2(nlayers), tauself, taufor, n2m1, n2m2, taun2
2748 : real(kind=r8) :: refrat_planck_a, refrat_m_a
2749 : real(kind=r8) :: tau_major, tau_major1
2750 :
2751 :
2752 : ! Minor gas mapping level :
2753 : ! Lower - Nitrogen Continuum, P = 1053., T = 294.
2754 :
2755 : ! Calculate reference ratio to be used in calculation of Planck
2756 : ! fraction in lower atmosphere.
2757 : ! P = 1053. mb (Level 1)
2758 486000 : refrat_planck_a = chi_mls(4,1)/chi_mls(2,1)
2759 :
2760 : ! P = 1053.
2761 486000 : refrat_m_a = chi_mls(4,1)/chi_mls(2,1)
2762 :
2763 : ! Compute the optical depth by interpolating in ln(pressure),
2764 : ! temperature, and appropriate species. Below laytrop, the water
2765 : ! vapor self-continuum and foreign continuum is interpolated
2766 : ! (in temperature) separately.
2767 :
2768 : ! Lower atmosphere loop
2769 9219976 : do lay = 1, laytrop
2770 :
2771 8733976 : speccomb(lay) = coln2o(lay) + rat_n2oco2(lay)*colco2(lay)
2772 8733976 : specparm(lay) = coln2o(lay)/speccomb(lay)
2773 8733976 : if (specparm(lay) .ge. oneminus) specparm(lay) = oneminus
2774 8733976 : specmult(lay) = 8._r8*(specparm(lay))
2775 8733976 : js(lay) = 1 + int(specmult(lay))
2776 8733976 : fs(lay) = mod(specmult(lay),1.0_r8)
2777 :
2778 8733976 : speccomb1(lay) = coln2o(lay) + rat_n2oco2_1(lay)*colco2(lay)
2779 8733976 : specparm1(lay) = coln2o(lay)/speccomb1(lay)
2780 8733976 : if (specparm1(lay) .ge. oneminus) specparm1(lay) = oneminus
2781 8733976 : specmult1(lay) = 8._r8*(specparm1(lay))
2782 8733976 : js1(lay) = 1 + int(specmult1(lay))
2783 8733976 : fs1(lay) = mod(specmult1(lay),1.0_r8)
2784 :
2785 8733976 : speccomb_mn2(lay) = coln2o(lay) + refrat_m_a*colco2(lay)
2786 8733976 : specparm_mn2(lay) = coln2o(lay)/speccomb_mn2(lay)
2787 8733976 : if (specparm_mn2(lay) .ge. oneminus) specparm_mn2(lay) = oneminus
2788 8733976 : specmult_mn2(lay) = 8._r8*specparm_mn2(lay)
2789 8733976 : jmn2(lay) = 1 + int(specmult_mn2(lay))
2790 8733976 : fmn2(lay) = mod(specmult_mn2(lay),1.0_r8)
2791 :
2792 8733976 : speccomb_planck(lay) = coln2o(lay)+refrat_planck_a*colco2(lay)
2793 8733976 : specparm_planck(lay) = coln2o(lay)/speccomb_planck(lay)
2794 8733976 : if (specparm_planck(lay) .ge. oneminus) specparm_planck(lay)=oneminus
2795 8733976 : specmult_planck(lay) = 8._r8*specparm_planck(lay)
2796 8733976 : jpl(lay)= 1 + int(specmult_planck(lay))
2797 8733976 : fpl(lay) = mod(specmult_planck(lay),1.0_r8)
2798 :
2799 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(15) + js(lay)
2800 8733976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(15) + js1(lay)
2801 :
2802 9219976 : scalen2(lay) = colbrd(lay)*scaleminor(lay)
2803 : enddo
2804 9219976 : do lay = 1, laytrop
2805 9219976 : if (specparm(lay) .lt. 0.125_r8) then
2806 0 : p = fs(lay) - 1
2807 0 : p4 = p**4
2808 0 : fk0 = p4
2809 0 : fk1 = 1 - p - 2.0_r8*p4
2810 0 : fk2 = p + p4
2811 0 : fac000(lay) = fk0*fac00(lay)
2812 0 : fac100(lay) = fk1*fac00(lay)
2813 0 : fac200(lay) = fk2*fac00(lay)
2814 0 : fac010(lay) = fk0*fac10(lay)
2815 0 : fac110(lay) = fk1*fac10(lay)
2816 0 : fac210(lay) = fk2*fac10(lay)
2817 8733976 : else if (specparm(lay) .gt. 0.875_r8) then
2818 0 : p = -fs(lay)
2819 0 : p4 = p**4
2820 0 : fk0 = p4
2821 0 : fk1 = 1 - p - 2.0_r8*p4
2822 0 : fk2 = p + p4
2823 0 : fac000(lay) = fk0*fac00(lay)
2824 0 : fac100(lay) = fk1*fac00(lay)
2825 0 : fac200(lay) = fk2*fac00(lay)
2826 0 : fac010(lay) = fk0*fac10(lay)
2827 0 : fac110(lay) = fk1*fac10(lay)
2828 0 : fac210(lay) = fk2*fac10(lay)
2829 : else
2830 8733976 : fac000(lay) = (1._r8 - fs(lay)) * fac00(lay)
2831 8733976 : fac010(lay) = (1._r8 - fs(lay)) * fac10(lay)
2832 8733976 : fac100(lay) = fs(lay) * fac00(lay)
2833 8733976 : fac110(lay) = fs(lay) * fac10(lay)
2834 : endif
2835 : enddo
2836 9219976 : do lay = 1, laytrop
2837 9219976 : if (specparm1(lay) .lt. 0.125_r8) then
2838 0 : p = fs1(lay) - 1
2839 0 : p4 = p**4
2840 0 : fk0 = p4
2841 0 : fk1 = 1 - p - 2.0_r8*p4
2842 0 : fk2 = p + p4
2843 0 : fac001(lay) = fk0*fac01(lay)
2844 0 : fac101(lay) = fk1*fac01(lay)
2845 0 : fac201(lay) = fk2*fac01(lay)
2846 0 : fac011(lay) = fk0*fac11(lay)
2847 0 : fac111(lay) = fk1*fac11(lay)
2848 0 : fac211(lay) = fk2*fac11(lay)
2849 8733976 : else if (specparm1(lay) .gt. 0.875_r8) then
2850 0 : p = -fs1(lay)
2851 0 : p4 = p**4
2852 0 : fk0 = p4
2853 0 : fk1 = 1 - p - 2.0_r8*p4
2854 0 : fk2 = p + p4
2855 0 : fac001(lay) = fk0*fac01(lay)
2856 0 : fac101(lay) = fk1*fac01(lay)
2857 0 : fac201(lay) = fk2*fac01(lay)
2858 0 : fac011(lay) = fk0*fac11(lay)
2859 0 : fac111(lay) = fk1*fac11(lay)
2860 0 : fac211(lay) = fk2*fac11(lay)
2861 : else
2862 8733976 : fac001(lay) = (1._r8 - fs1(lay)) * fac01(lay)
2863 8733976 : fac011(lay) = (1._r8 - fs1(lay)) * fac11(lay)
2864 8733976 : fac101(lay) = fs1(lay) * fac01(lay)
2865 8733976 : fac111(lay) = fs1(lay) * fac11(lay)
2866 : endif
2867 :
2868 : enddo
2869 9219976 : do lay = 1, laytrop
2870 26687928 : do ig = 1, ng15
2871 87339760 : tauself = selffac(lay)* (selfref(indself(lay),ig) + selffrac(lay) * &
2872 104807712 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
2873 87339760 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
2874 104807712 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
2875 104807712 : n2m1 = ka_mn2(jmn2(lay),indminor(lay),ig) + fmn2(lay) * &
2876 122275664 : (ka_mn2(jmn2(lay)+1,indminor(lay),ig) - ka_mn2(jmn2(lay),indminor(lay),ig))
2877 104807712 : n2m2 = ka_mn2(jmn2(lay),indminor(lay)+1,ig) + fmn2(lay) * &
2878 122275664 : (ka_mn2(jmn2(lay)+1,indminor(lay)+1,ig) - ka_mn2(jmn2(lay),indminor(lay)+1,ig))
2879 17467952 : taun2 = scalen2(lay) * (n2m1 + minorfrac(lay) * (n2m2 - n2m1))
2880 :
2881 17467952 : if (specparm(lay) .lt. 0.125_r8) then
2882 0 : tau_major = speccomb(lay) * &
2883 0 : (fac000(lay) * absa(ind0(lay),ig) + &
2884 0 : fac100(lay) * absa(ind0(lay)+1,ig) + &
2885 0 : fac200(lay) * absa(ind0(lay)+2,ig) + &
2886 0 : fac010(lay) * absa(ind0(lay)+9,ig) + &
2887 0 : fac110(lay) * absa(ind0(lay)+10,ig) + &
2888 0 : fac210(lay) * absa(ind0(lay)+11,ig))
2889 17467952 : else if (specparm(lay) .gt. 0.875_r8) then
2890 0 : tau_major = speccomb(lay) * &
2891 0 : (fac200(lay) * absa(ind0(lay)-1,ig) + &
2892 0 : fac100(lay) * absa(ind0(lay),ig) + &
2893 0 : fac000(lay) * absa(ind0(lay)+1,ig) + &
2894 0 : fac210(lay) * absa(ind0(lay)+8,ig) + &
2895 0 : fac110(lay) * absa(ind0(lay)+9,ig) + &
2896 0 : fac010(lay) * absa(ind0(lay)+10,ig))
2897 : else
2898 17467952 : tau_major = speccomb(lay) * &
2899 69871808 : (fac000(lay) * absa(ind0(lay),ig) + &
2900 69871808 : fac100(lay) * absa(ind0(lay)+1,ig) + &
2901 69871808 : fac010(lay) * absa(ind0(lay)+9,ig) + &
2902 244551328 : fac110(lay) * absa(ind0(lay)+10,ig))
2903 : endif
2904 :
2905 17467952 : if (specparm1(lay) .lt. 0.125_r8) then
2906 0 : tau_major1 = speccomb1(lay) * &
2907 0 : (fac001(lay) * absa(ind1(lay),ig) + &
2908 0 : fac101(lay) * absa(ind1(lay)+1,ig) + &
2909 0 : fac201(lay) * absa(ind1(lay)+2,ig) + &
2910 0 : fac011(lay) * absa(ind1(lay)+9,ig) + &
2911 0 : fac111(lay) * absa(ind1(lay)+10,ig) + &
2912 0 : fac211(lay) * absa(ind1(lay)+11,ig))
2913 17467952 : else if (specparm1(lay) .gt. 0.875_r8) then
2914 0 : tau_major1 = speccomb1(lay) * &
2915 0 : (fac201(lay) * absa(ind1(lay)-1,ig) + &
2916 0 : fac101(lay) * absa(ind1(lay),ig) + &
2917 0 : fac001(lay) * absa(ind1(lay)+1,ig) + &
2918 0 : fac211(lay) * absa(ind1(lay)+8,ig) + &
2919 0 : fac111(lay) * absa(ind1(lay)+9,ig) + &
2920 0 : fac011(lay) * absa(ind1(lay)+10,ig))
2921 : else
2922 17467952 : tau_major1 = speccomb1(lay) * &
2923 69871808 : (fac001(lay) * absa(ind1(lay),ig) + &
2924 69871808 : fac101(lay) * absa(ind1(lay)+1,ig) + &
2925 69871808 : fac011(lay) * absa(ind1(lay)+9,ig) + &
2926 244551328 : fac111(lay) * absa(ind1(lay)+10,ig))
2927 : endif
2928 :
2929 34935904 : taug(lay,ngs14+ig) = tau_major + tau_major1 &
2930 : + tauself + taufor &
2931 34935904 : + taun2
2932 104807712 : fracs(lay,ngs14+ig) = fracrefa(ig,jpl(lay)) + fpl(lay) * &
2933 131009640 : (fracrefa(ig,jpl(lay)+1)-fracrefa(ig,jpl(lay)))
2934 : enddo
2935 : enddo
2936 :
2937 : ! Upper atmosphere loop
2938 4874024 : do lay = laytrop+1, nlayers
2939 13650072 : do ig = 1, ng15
2940 8776048 : taug(lay,ngs14+ig) = 0.0_r8
2941 13164072 : fracs(lay,ngs14+ig) = 0.0_r8
2942 : enddo
2943 : enddo
2944 :
2945 486000 : end subroutine taugb15
2946 :
2947 : !----------------------------------------------------------------------------
2948 486000 : subroutine taugb16
2949 : !----------------------------------------------------------------------------
2950 : !
2951 : ! band 16: 2600-3250 cm-1 (low key- h2o,ch4; high key - ch4)
2952 : !----------------------------------------------------------------------------
2953 :
2954 : ! ------- Modules -------
2955 :
2956 : use parrrtm, only : ng16, ngs15
2957 : use rrlw_ref, only : chi_mls
2958 : use rrlw_kg16, only : fracrefa, fracrefb, absa, absb, &
2959 : selfref, forref
2960 :
2961 : ! ------- Declarations -------
2962 :
2963 : ! Local
2964 972000 : integer :: lay, ind0(nlayers), ind1(nlayers), ig
2965 972000 : integer, dimension(nlayers) :: js, js1, jpl
2966 972000 : real(kind=r8), dimension(nlayers) :: speccomb, specparm, specmult, fs
2967 972000 : real(kind=r8), dimension(nlayers) :: speccomb1, specparm1, specmult1, fs1
2968 972000 : real(kind=r8), dimension(nlayers) :: speccomb_planck, specparm_planck, specmult_planck, fpl
2969 : real(kind=r8) :: p, p4, fk0, fk1, fk2
2970 972000 : real(kind=r8), dimension(nlayers) :: fac000, fac100, fac200, fac010, fac110, fac210
2971 972000 : real(kind=r8), dimension(nlayers) :: fac001, fac101, fac201, fac011, fac111, fac211
2972 : real(kind=r8) :: tauself, taufor
2973 : real(kind=r8) :: refrat_planck_a
2974 : real(kind=r8) :: tau_major, tau_major1
2975 :
2976 :
2977 : ! Calculate reference ratio to be used in calculation of Planck
2978 : ! fraction in lower atmosphere.
2979 :
2980 : ! P = 387. mb (Level 6)
2981 486000 : refrat_planck_a = chi_mls(1,6)/chi_mls(6,6)
2982 :
2983 : ! Compute the optical depth by interpolating in ln(pressure),
2984 : ! temperature,and appropriate species. Below laytrop, the water
2985 : ! vapor self-continuum and foreign continuum is interpolated
2986 : ! (in temperature) separately.
2987 :
2988 : ! Lower atmosphere loop
2989 9219976 : do lay = 1, laytrop
2990 :
2991 8733976 : speccomb(lay) = colh2o(lay) + rat_h2och4(lay)*colch4(lay)
2992 8733976 : specparm(lay) = colh2o(lay)/speccomb(lay)
2993 8733976 : if (specparm(lay) .ge. oneminus) specparm(lay) = oneminus
2994 8733976 : specmult(lay) = 8._r8*(specparm(lay))
2995 8733976 : js(lay) = 1 + int(specmult(lay))
2996 8733976 : fs(lay) = mod(specmult(lay),1.0_r8)
2997 :
2998 8733976 : speccomb1(lay) = colh2o(lay) + rat_h2och4_1(lay)*colch4(lay)
2999 8733976 : specparm1(lay) = colh2o(lay)/speccomb1(lay)
3000 8733976 : if (specparm1(lay) .ge. oneminus) specparm1(lay) = oneminus
3001 8733976 : specmult1(lay) = 8._r8*(specparm1(lay))
3002 8733976 : js1(lay) = 1 + int(specmult1(lay))
3003 8733976 : fs1(lay) = mod(specmult1(lay),1.0_r8)
3004 :
3005 8733976 : speccomb_planck(lay) = colh2o(lay)+refrat_planck_a*colch4(lay)
3006 8733976 : specparm_planck(lay) = colh2o(lay)/speccomb_planck(lay)
3007 8733976 : if (specparm_planck(lay) .ge. oneminus) specparm_planck(lay)=oneminus
3008 8733976 : specmult_planck(lay) = 8._r8*specparm_planck(lay)
3009 8733976 : jpl(lay)= 1 + int(specmult_planck(lay))
3010 8733976 : fpl(lay) = mod(specmult_planck(lay),1.0_r8)
3011 :
3012 8733976 : ind0(lay) = ((jp(lay)-1)*5+(jt(lay)-1))*nspa(16) + js(lay)
3013 9219976 : ind1(lay) = (jp(lay)*5+(jt1(lay)-1))*nspa(16) + js1(lay)
3014 : enddo
3015 9219976 : do lay = 1, laytrop
3016 9219976 : if (specparm(lay) .lt. 0.125_r8) then
3017 292858 : p = fs(lay) - 1
3018 292858 : p4 = p**4
3019 292858 : fk0 = p4
3020 292858 : fk1 = 1 - p - 2.0_r8*p4
3021 292858 : fk2 = p + p4
3022 292858 : fac000(lay) = fk0*fac00(lay)
3023 292858 : fac100(lay) = fk1*fac00(lay)
3024 292858 : fac200(lay) = fk2*fac00(lay)
3025 292858 : fac010(lay) = fk0*fac10(lay)
3026 292858 : fac110(lay) = fk1*fac10(lay)
3027 292858 : fac210(lay) = fk2*fac10(lay)
3028 8441118 : else if (specparm(lay) .gt. 0.875_r8) then
3029 115205 : p = -fs(lay)
3030 115205 : p4 = p**4
3031 115205 : fk0 = p4
3032 115205 : fk1 = 1 - p - 2.0_r8*p4
3033 115205 : fk2 = p + p4
3034 115205 : fac000(lay) = fk0*fac00(lay)
3035 115205 : fac100(lay) = fk1*fac00(lay)
3036 115205 : fac200(lay) = fk2*fac00(lay)
3037 115205 : fac010(lay) = fk0*fac10(lay)
3038 115205 : fac110(lay) = fk1*fac10(lay)
3039 115205 : fac210(lay) = fk2*fac10(lay)
3040 : else
3041 8325913 : fac000(lay) = (1._r8 - fs(lay)) * fac00(lay)
3042 8325913 : fac010(lay) = (1._r8 - fs(lay)) * fac10(lay)
3043 8325913 : fac100(lay) = fs(lay) * fac00(lay)
3044 8325913 : fac110(lay) = fs(lay) * fac10(lay)
3045 : endif
3046 : enddo
3047 9219976 : do lay = 1, laytrop
3048 :
3049 9219976 : if (specparm1(lay) .lt. 0.125_r8) then
3050 26723 : p = fs1(lay) - 1
3051 26723 : p4 = p**4
3052 26723 : fk0 = p4
3053 26723 : fk1 = 1 - p - 2.0_r8*p4
3054 26723 : fk2 = p + p4
3055 26723 : fac001(lay) = fk0*fac01(lay)
3056 26723 : fac101(lay) = fk1*fac01(lay)
3057 26723 : fac201(lay) = fk2*fac01(lay)
3058 26723 : fac011(lay) = fk0*fac11(lay)
3059 26723 : fac111(lay) = fk1*fac11(lay)
3060 26723 : fac211(lay) = fk2*fac11(lay)
3061 8707253 : else if (specparm1(lay) .gt. 0.875_r8) then
3062 1046819 : p = -fs1(lay)
3063 1046819 : p4 = p**4
3064 1046819 : fk0 = p4
3065 1046819 : fk1 = 1 - p - 2.0_r8*p4
3066 1046819 : fk2 = p + p4
3067 1046819 : fac001(lay) = fk0*fac01(lay)
3068 1046819 : fac101(lay) = fk1*fac01(lay)
3069 1046819 : fac201(lay) = fk2*fac01(lay)
3070 1046819 : fac011(lay) = fk0*fac11(lay)
3071 1046819 : fac111(lay) = fk1*fac11(lay)
3072 1046819 : fac211(lay) = fk2*fac11(lay)
3073 : else
3074 7660434 : fac001(lay) = (1._r8 - fs1(lay)) * fac01(lay)
3075 7660434 : fac011(lay) = (1._r8 - fs1(lay)) * fac11(lay)
3076 7660434 : fac101(lay) = fs1(lay) * fac01(lay)
3077 7660434 : fac111(lay) = fs1(lay) * fac11(lay)
3078 : endif
3079 : enddo
3080 9219976 : do lay = 1, laytrop
3081 :
3082 26687928 : do ig = 1, ng16
3083 87339760 : tauself = selffac(lay)* (selfref(indself(lay),ig) + selffrac(lay) * &
3084 104807712 : (selfref(indself(lay)+1,ig) - selfref(indself(lay),ig)))
3085 87339760 : taufor = forfac(lay) * (forref(indfor(lay),ig) + forfrac(lay) * &
3086 104807712 : (forref(indfor(lay)+1,ig) - forref(indfor(lay),ig)))
3087 :
3088 17467952 : if (specparm(lay) .lt. 0.125_r8) then
3089 585716 : tau_major = speccomb(lay) * &
3090 2342864 : (fac000(lay) * absa(ind0(lay),ig) + &
3091 2342864 : fac100(lay) * absa(ind0(lay)+1,ig) + &
3092 2342864 : fac200(lay) * absa(ind0(lay)+2,ig) + &
3093 2342864 : fac010(lay) * absa(ind0(lay)+9,ig) + &
3094 2342864 : fac110(lay) * absa(ind0(lay)+10,ig) + &
3095 12885752 : fac210(lay) * absa(ind0(lay)+11,ig))
3096 16882236 : else if (specparm(lay) .gt. 0.875_r8) then
3097 230410 : tau_major = speccomb(lay) * &
3098 921640 : (fac200(lay) * absa(ind0(lay)-1,ig) + &
3099 921640 : fac100(lay) * absa(ind0(lay),ig) + &
3100 921640 : fac000(lay) * absa(ind0(lay)+1,ig) + &
3101 921640 : fac210(lay) * absa(ind0(lay)+8,ig) + &
3102 921640 : fac110(lay) * absa(ind0(lay)+9,ig) + &
3103 5069020 : fac010(lay) * absa(ind0(lay)+10,ig))
3104 : else
3105 16651826 : tau_major = speccomb(lay) * &
3106 66607304 : (fac000(lay) * absa(ind0(lay),ig) + &
3107 66607304 : fac100(lay) * absa(ind0(lay)+1,ig) + &
3108 66607304 : fac010(lay) * absa(ind0(lay)+9,ig) + &
3109 233125564 : fac110(lay) * absa(ind0(lay)+10,ig))
3110 : endif
3111 :
3112 17467952 : if (specparm1(lay) .lt. 0.125_r8) then
3113 53446 : tau_major1 = speccomb1(lay) * &
3114 213784 : (fac001(lay) * absa(ind1(lay),ig) + &
3115 213784 : fac101(lay) * absa(ind1(lay)+1,ig) + &
3116 213784 : fac201(lay) * absa(ind1(lay)+2,ig) + &
3117 213784 : fac011(lay) * absa(ind1(lay)+9,ig) + &
3118 213784 : fac111(lay) * absa(ind1(lay)+10,ig) + &
3119 1175812 : fac211(lay) * absa(ind1(lay)+11,ig))
3120 17414506 : else if (specparm1(lay) .gt. 0.875_r8) then
3121 2093638 : tau_major1 = speccomb1(lay) * &
3122 8374552 : (fac201(lay) * absa(ind1(lay)-1,ig) + &
3123 8374552 : fac101(lay) * absa(ind1(lay),ig) + &
3124 8374552 : fac001(lay) * absa(ind1(lay)+1,ig) + &
3125 8374552 : fac211(lay) * absa(ind1(lay)+8,ig) + &
3126 8374552 : fac111(lay) * absa(ind1(lay)+9,ig) + &
3127 46060036 : fac011(lay) * absa(ind1(lay)+10,ig))
3128 : else
3129 15320868 : tau_major1 = speccomb1(lay) * &
3130 61283472 : (fac001(lay) * absa(ind1(lay),ig) + &
3131 61283472 : fac101(lay) * absa(ind1(lay)+1,ig) + &
3132 61283472 : fac011(lay) * absa(ind1(lay)+9,ig) + &
3133 214492152 : fac111(lay) * absa(ind1(lay)+10,ig))
3134 : endif
3135 :
3136 34935904 : taug(lay,ngs15+ig) = tau_major + tau_major1 &
3137 34935904 : + tauself + taufor
3138 104807712 : fracs(lay,ngs15+ig) = fracrefa(ig,jpl(lay)) + fpl(lay) * &
3139 131009640 : (fracrefa(ig,jpl(lay)+1)-fracrefa(ig,jpl(lay)))
3140 : enddo
3141 : enddo
3142 :
3143 : ! Upper atmosphere loop
3144 4874024 : do lay = laytrop+1, nlayers
3145 4388024 : ind0(lay) = ((jp(lay)-13)*5+(jt(lay)-1))*nspb(16) + 1
3146 4874024 : ind1(lay) = ((jp(lay)-12)*5+(jt1(lay)-1))*nspb(16) + 1
3147 : enddo
3148 4874024 : do lay = laytrop+1, nlayers
3149 13650072 : do ig = 1, ng16
3150 26328144 : taug(lay,ngs15+ig) = colch4(lay) * &
3151 35104192 : (fac00(lay) * absb(ind0(lay),ig) + &
3152 35104192 : fac10(lay) * absb(ind0(lay)+1,ig) + &
3153 35104192 : fac01(lay) * absb(ind1(lay),ig) + &
3154 140416768 : fac11(lay) * absb(ind1(lay)+1,ig))
3155 13164072 : fracs(lay,ngs15+ig) = fracrefb(ig)
3156 : enddo
3157 : enddo
3158 :
3159 486000 : end subroutine taugb16
3160 :
3161 : end subroutine taumol
3162 :
3163 : end module rrtmg_lw_taumol
3164 :
|
