rrtmg_sw_spcvrt.f90 coverage: 0.00 %func 0.00 %block
1) ! path: $Source: /storm/rc1/cvsroot/rc/rrtmg_sw/src/rrtmg_sw_spcvrt.f90,v $
2) ! author: $Author: mike $
3) ! revision: $Revision: 1.2 $
4) ! created: $Date: 2007/08/23 20:40:14 $
5)
6) module rrtmg_sw_spcvrt
7)
8) ! --------------------------------------------------------------------------
9) ! | |
10) ! | Copyright 2002-2007, 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) use ppgrid, only: pcols, begchunk, endchunk
22)
23) ! use parkind, only : jpim, jprb
24) use parrrsw, only : nbndsw, ngptsw, mxmol, jpband
25) use rrsw_tbl, only : tblint, bpade, od_lo, exp_tbl
26) use rrsw_vsn, only : hvrspv, hnamspv
27) use rrsw_wvn, only : ngc, ngs
28) use rrtmg_sw_reftra, only: reftra_sw
29) use rrtmg_sw_taumol, only: taumol_sw
30) use rrtmg_sw_vrtqdr, only: vrtqdr_sw
31)
32) implicit none
33)
34) contains
35)
36) ! ---------------------------------------------------------------------------
37) subroutine spcvrt_sw &
38) (lchnk, iplon, nlayers, istart, iend, icpr, idelm, iout, dotau, &
39) pavel, tavel, pz, tz, tbound, palbd, palbp, &
40) pclfr, ptauc, pasyc, pomgc, ptaucorig, &
41) ptaua, pasya, pomga, prmu0, coldry, wkl, adjflux, &
42) laytrop, layswtch, laylow, jp, jt, jt1, &
43) co2mult, colch4, colco2, colh2o, colmol, coln2o, colo2, colo3, &
44) fac00, fac01, fac10, fac11, &
45) selffac, selffrac, indself, forfac, forfrac, indfor, &
46) pbbfd, pbbfu, pbbcd, pbbcu, puvfd, puvcd, pnifd, pnicd, pnifu, pnicu, &
47) pbbfddir, pbbcddir, puvfddir, puvcddir, pnifddir, pnicddir)
48) ! ---------------------------------------------------------------------------
49) !
50) ! Purpose: Contains spectral loop to compute the shortwave radiative fluxes,
51) ! using the two-stream method of H. Barker.
52) !
53) ! Interface: *spcvrt_sw* is called from *rrtmg_sw.F90* or rrtmg_sw.1col.F90*
54) !
55) ! Method:
56) ! Adapted from two-stream model of H. Barker;
57) ! Two-stream model options (selected with kmodts in rrtmg_sw_reftra.F90):
58) ! 1: Eddington, 2: PIFM, Zdunkowski et al., 3: discret ordinates
59) !
60) ! Modifications:
61) !
62) ! Original: H. Barker
63) ! Revision: Merge with RRTMG_SW: J.-J.Morcrette, ECMWF, Feb 2003
64) ! Revision: Add adjustment for Earth/Sun distance : MJIacono, AER, Oct 2003
65) ! Revision: Bug fix for use of PALBP and PALBD: MJIacono, AER, Nov 2003
66) ! Revision: Bug fix to apply delta scaling to clear sky: AER, Dec 2004
67) ! Revision: Code modified so that delta scaling is not done in cloudy profiles
68) ! if routine cldprop is used; delta scaling can be applied by swithcing
69) ! code below if cldprop is not used to get cloud properties.
70) ! AER, Jan 2005
71) ! Revision: Uniform formatting for RRTMG: MJIacono, AER, Jul 2006
72) ! Revision: Use exponential lookup table for transmittance: MJIacono, AER,
73) ! Aug 2007
74) !
75) ! ------------------------------------------------------------------
76)
77) ! ------- Declarations ------
78)
79) ! -------- Input -------
80)
81) integer, intent(in) :: lchnk
82) integer, intent(in) :: iplon ! column loop index
83) integer, intent(in) :: nlayers
84) integer, intent(in) :: istart
85) integer, intent(in) :: iend
86) integer, intent(in) :: icpr
87) integer, intent(in) :: idelm ! delta-m scaling flag
88) ! [0 = direct and diffuse fluxes are unscaled]
89) ! [1 = direct and diffuse fluxes are scaled]
90) integer, intent(in) :: iout
91) logical, intent(in) :: dotau ! True -> do tau calculation
92)
93) integer, intent(in) :: laytrop
94) integer, intent(in) :: layswtch
95) integer, intent(in) :: laylow
96)
97) integer, intent(in) :: indfor(:)
98) ! Dimensions: (nlayers)
99) integer, intent(in) :: indself(:)
100) ! Dimensions: (nlayers)
101) integer, intent(in) :: jp(:)
102) ! Dimensions: (nlayers)
103) integer, intent(in) :: jt(:)
104) ! Dimensions: (nlayers)
105) integer, intent(in) :: jt1(:)
106) ! Dimensions: (nlayers)
107)
108) real(kind=r8), intent(in) :: pavel(:) ! layer pressure (hPa, mb)
109) ! Dimensions: (nlayers)
110) real(kind=r8), intent(in) :: tavel(:) ! layer temperature (K)
111) ! Dimensions: (nlayers)
112) real(kind=r8), intent(in) :: pz(0:) ! level (interface) pressure (hPa, mb)
113) ! Dimensions: (0:nlayers)
114) real(kind=r8), intent(in) :: tz(0:) ! level temperatures (hPa, mb)
115) ! Dimensions: (0:nlayers)
116) real(kind=r8), intent(in) :: tbound ! surface temperature (K)
117) real(kind=r8), intent(in) :: wkl(:,:) ! molecular amounts (mol/cm2)
118) ! Dimensions: (mxmol,nlayers)
119) real(kind=r8), intent(in) :: coldry(:) ! dry air column density (mol/cm2)
120) ! Dimensions: (nlayers)
121) real(kind=r8), intent(in) :: colmol(:)
122) ! Dimensions: (nlayers)
123) real(kind=r8), intent(in) :: adjflux(:) ! Earth/Sun distance adjustment
124) ! Dimensions: (jpband)
125)
126) real(kind=r8), intent(in) :: palbd(:) ! surface albedo (diffuse)
127) ! Dimensions: (nbndsw)
128) real(kind=r8), intent(in) :: palbp(:) ! surface albedo (direct)
129) ! Dimensions: (nbndsw)
130) real(kind=r8), intent(in) :: prmu0 ! cosine of solar zenith angle
131) real(kind=r8), intent(in) :: pclfr(:) ! cloud fraction
132) ! Dimensions: (nlayers)
133) real(kind=r8), intent(in) :: ptauc(:,:) ! cloud optical depth
134) ! Dimensions: (nlayers,nbndsw)
135) real(kind=r8), intent(in) :: pasyc(:,:) ! cloud asymmetry parameter
136) ! Dimensions: (nlayers,nbndsw)
137) real(kind=r8), intent(in) :: pomgc(:,:) ! cloud single scattering albedo
138) ! Dimensions: (nlayers,nbndsw)
139) real(kind=r8), intent(in) :: ptaucorig(:,:) ! cloud optical depth, non-delta scaled
140) ! Dimensions: (nlayers,nbndsw)
141) real(kind=r8), intent(in) :: ptaua(:,:) ! aerosol optical depth
142) ! Dimensions: (nlayers,nbndsw)
143) real(kind=r8), intent(in) :: pasya(:,:) ! aerosol asymmetry parameter
144) ! Dimensions: (nlayers,nbndsw)
145) real(kind=r8), intent(in) :: pomga(:,:) ! aerosol single scattering albedo
146) ! Dimensions: (nlayers,nbndsw)
147)
148) real(kind=r8), intent(in) :: colh2o(:)
149) ! Dimensions: (nlayers)
150) real(kind=r8), intent(in) :: colco2(:)
151) ! Dimensions: (nlayers)
152) real(kind=r8), intent(in) :: colch4(:)
153) ! Dimensions: (nlayers)
154) real(kind=r8), intent(in) :: co2mult(:)
155) ! Dimensions: (nlayers)
156) real(kind=r8), intent(in) :: colo3(:)
157) ! Dimensions: (nlayers)
158) real(kind=r8), intent(in) :: colo2(:)
159) ! Dimensions: (nlayers)
160) real(kind=r8), intent(in) :: coln2o(:)
161) ! Dimensions: (nlayers)
162)
163) real(kind=r8), intent(in) :: forfac(:)
164) ! Dimensions: (nlayers)
165) real(kind=r8), intent(in) :: forfrac(:)
166) ! Dimensions: (nlayers)
167) real(kind=r8), intent(in) :: selffac(:)
168) ! Dimensions: (nlayers)
169) real(kind=r8), intent(in) :: selffrac(:)
170) ! Dimensions: (nlayers)
171) real(kind=r8), intent(in) :: fac00(:)
172) ! Dimensions: (nlayers)
173) real(kind=r8), intent(in) :: fac01(:)
174) ! Dimensions: (nlayers)
175) real(kind=r8), intent(in) :: fac10(:)
176) ! Dimensions: (nlayers)
177) real(kind=r8), intent(in) :: fac11(:)
178) ! Dimensions: (nlayers)
179)
180) ! ------- Output -------
181) ! All Dimensions: (nlayers+1)
182) real(kind=r8), intent(out) :: pbbcd(:)
183) real(kind=r8), intent(out) :: pbbcu(:)
184) real(kind=r8), intent(out) :: pbbfd(:)
185) real(kind=r8), intent(out) :: pbbfu(:)
186) real(kind=r8), intent(out) :: pbbfddir(:)
187) real(kind=r8), intent(out) :: pbbcddir(:)
188)
189) real(kind=r8), intent(out) :: puvcd(:)
190) real(kind=r8), intent(out) :: puvfd(:)
191) real(kind=r8), intent(out) :: puvcddir(:)
192) real(kind=r8), intent(out) :: puvfddir(:)
193)
194) real(kind=r8), intent(out) :: pnicd(:)
195) real(kind=r8), intent(out) :: pnifd(:)
196) real(kind=r8), intent(out) :: pnicddir(:)
197) real(kind=r8), intent(out) :: pnifddir(:)
198)
199) ! Added for net near-IR flux diagnostic
200) real(kind=r8), intent(out) :: pnicu(:)
201) real(kind=r8), intent(out) :: pnifu(:)
202)
203) ! Output - inactive ! All Dimensions: (nlayers+1)
204) ! real(kind=r8), intent(out) :: puvcu(:)
205) ! real(kind=r8), intent(out) :: puvfu(:)
206) ! real(kind=r8), intent(out) :: pnicu(:)
207) ! real(kind=r8), intent(out) :: pnifu(:)
208) ! real(kind=r8), intent(out) :: pvscd(:)
209) ! real(kind=r8), intent(out) :: pvscu(:)
210) ! real(kind=r8), intent(out) :: pvsfd(:)
211) ! real(kind=r8), intent(out) :: pvsfu(:)
212)
213)
214) ! ------- Local -------
215)
216) logical :: lrtchkclr(nlayers),lrtchkcld(nlayers)
217)
218) integer :: klev
219) integer :: ib1, ib2, ibm, igt, ikl, ikp, ikx
220) integer :: iw, jb, jg, jl, jk
221) ! integer, parameter :: nuv = ??
222) ! integer, parameter :: nvs = ??
223) integer :: itind
224)
225) real(kind=r8) :: tblind, ze1
226) real(kind=r8) :: zclear, zcloud
227) real(kind=r8) :: zdbt(nlayers+1), zdbt_nodel(nlayers+1)
228) real(kind=r8) :: zgc(nlayers), zgcc(nlayers), zgco(nlayers)
229) real(kind=r8) :: zomc(nlayers), zomcc(nlayers), zomco(nlayers)
230) real(kind=r8) :: zrdnd(nlayers+1), zrdndc(nlayers+1)
231) real(kind=r8) :: zref(nlayers+1), zrefc(nlayers+1), zrefo(nlayers+1)
232) real(kind=r8) :: zrefd(nlayers+1), zrefdc(nlayers+1), zrefdo(nlayers+1)
233) real(kind=r8) :: zrup(nlayers+1), zrupd(nlayers+1)
234) real(kind=r8) :: zrupc(nlayers+1), zrupdc(nlayers+1)
235) real(kind=r8) :: zs1(nlayers+1)
236) real(kind=r8) :: ztauc(nlayers), ztauo(nlayers)
237) real(kind=r8) :: ztdn(nlayers+1), ztdnd(nlayers+1), ztdbt(nlayers+1)
238) real(kind=r8) :: ztoc(nlayers), ztor(nlayers)
239) real(kind=r8) :: ztra(nlayers+1), ztrac(nlayers+1), ztrao(nlayers+1)
240) real(kind=r8) :: ztrad(nlayers+1), ztradc(nlayers+1), ztrado(nlayers+1)
241) real(kind=r8) :: zdbtc(nlayers+1), ztdbtc(nlayers+1)
242) real(kind=r8) :: zincflx(ngptsw), zdbtc_nodel(nlayers+1)
243) real(kind=r8) :: ztdbt_nodel(nlayers+1), ztdbtc_nodel(nlayers+1)
244)
245) real(kind=r8) :: zdbtmc, zdbtmo, zf, zgw, zreflect
246) real(kind=r8) :: zwf, tauorig, repclc
247) ! real(kind=r8) :: zincflux ! inactive
248)
249) ! Arrays from rrtmg_sw_taumoln routines
250)
251) ! real(kind=r8) :: ztaug(nlayers,16), ztaur(nlayers,16)
252) ! real(kind=r8) :: zsflxzen(16)
253) real(kind=r8) :: ztaug(nlayers,ngptsw), ztaur(nlayers,ngptsw)
254) real(kind=r8) :: zsflxzen(ngptsw)
255)
256) ! Arrays from rrtmg_sw_vrtqdr routine
257)
258) real(kind=r8) :: zcd(nlayers+1,ngptsw), zcu(nlayers+1,ngptsw)
259) real(kind=r8) :: zfd(nlayers+1,ngptsw), zfu(nlayers+1,ngptsw)
260)
261) ! Inactive arrays
262) ! real(kind=r8) :: zbbcd(nlayers+1), zbbcu(nlayers+1)
263) ! real(kind=r8) :: zbbfd(nlayers+1), zbbfu(nlayers+1)
264) ! real(kind=r8) :: zbbfddir(nlayers+1), zbbcddir(nlayers+1)
265)
266) ! Local storage arrays for zsflxzen, ztaug, and ztaur
267)
268) real(kind=r8), allocatable, save :: ztaugst(:,:,:,:) ! Gaseous optical depth storage array
269) real(kind=r8), allocatable, save :: ztaurst(:,:,:,:) ! Rayleigh optical depth storage array
270) real(kind=r8), allocatable, save :: zsflxzenst(:,:,:) ! Solar source function storage array
271)
272) ! Allocate storage arrays
273) if (.not.allocated(ztaugst)) allocate (ztaugst(pcols,nlayers,ngptsw,begchunk:endchunk))
274) if (.not.allocated(ztaurst)) allocate (ztaurst(pcols,nlayers,ngptsw,begchunk:endchunk))
275) if (.not.allocated(zsflxzenst)) allocate (zsflxzenst(pcols,ngptsw,begchunk:endchunk))
276)
277) ! ------------------------------------------------------------------
278)
279) ! Initializations
280)
281) ib1 = istart
282) ib2 = iend
283) klev = nlayers
284) iw = 0
285) repclc = 1.e-12_r8
286) ! zincflux = 0.0_r8
287)
288) do jk=1,klev+1
289) pbbcd(jk)=0._r8
290) pbbcu(jk)=0._r8
291) pbbfd(jk)=0._r8
292) pbbfu(jk)=0._r8
293) pbbcddir(jk)=0._r8
294) pbbfddir(jk)=0._r8
295) puvcd(jk)=0._r8
296) puvfd(jk)=0._r8
297) puvcddir(jk)=0._r8
298) puvfddir(jk)=0._r8
299) pnicd(jk)=0._r8
300) pnifd(jk)=0._r8
301) pnicddir(jk)=0._r8
302) pnifddir(jk)=0._r8
303) pnicu(jk)=0._r8
304) pnifu(jk)=0._r8
305) enddo
306)
307)
308) ! Perform calculations of optical depth and solar source function at interval
309) ! specified by dotau and store output arrays zsflxzen, ztaug, and ztaur in memory.
310) ! In intervening time steps, obtain zsflxzen, ztaug, and ztaur from stored arrays.
311) ! Memory storage should be replaced with I/O to netCDF file.
312)
313) if (dotau) then
314)
315) ! Calculate the optical depths for gaseous absorption and Rayleigh scattering
316)
317) call taumol_sw(klev, &
318) colh2o, colco2, colch4, colo2, colo3, colmol, &
319) laytrop, jp, jt, jt1, &
320) fac00, fac01, fac10, fac11, &
321) selffac, selffrac, indself, forfac, forfrac, indfor, &
322) zsflxzen, ztaug, ztaur)
323)
324) ! Store ztaug, ztaur, and zsflxzen for use when dotau is false
325) ztaugst(iplon,:nlayers,:ngptsw,lchnk) = ztaug(:nlayers,:ngptsw)
326) ztaurst(iplon,:nlayers,:ngptsw,lchnk) = ztaur(:nlayers,:ngptsw)
327) zsflxzenst(iplon,:ngptsw,lchnk) = zsflxzen(:ngptsw)
328)
329) else
330)
331) ! Restore ztaug, ztaur, and zsflxzen from storage when dotau is false
332) ztaug(:nlayers,:ngptsw) = ztaugst(iplon,:nlayers,:ngptsw,lchnk)
333) ztaur(:nlayers,:ngptsw) = ztaurst(iplon,:nlayers,:ngptsw,lchnk)
334) zsflxzen(:ngptsw) = zsflxzenst(iplon,:ngptsw,lchnk)
335)
336) endif
337)
338) ! Top of shortwave spectral band loop, jb = 16 -> 29; ibm = 1 -> 14
339)
340) jb = ib1-1 ! ???
341) do jb = ib1, ib2
342) ibm = jb-15
343) igt = ngc(ibm)
344)
345) ! Reinitialize g-point counter for each band if output for each band is requested.
346) if (iout.gt.0.and.ibm.ge.2) iw = ngs(ibm-1)
347)
348) ! do jk=1,klev+1
349) ! zbbcd(jk)=0.0_r8
350) ! zbbcu(jk)=0.0_r8
351) ! zbbfd(jk)=0.0_r8
352) ! zbbfu(jk)=0.0_r8
353) ! enddo
354)
355) ! Top of g-point interval loop within each band (iw is cumulative counter)
356) do jg = 1,igt
357) iw = iw+1
358)
359) ! Apply adjustments for correct Earth/Sun distance and zenith angle to incoming solar flux
360) zincflx(iw) = adjflux(jb) * zsflxzen(iw) * prmu0
361) ! zincflux = zincflux + adjflux(jb) * zsflxzen(iw) * prmu0 ! inactive
362)
363) ! Compute layer reflectances and transmittances for direct and diffuse sources,
364) ! first clear then cloudy
365)
366) ! zrefc(jk) direct albedo for clear
367) ! zrefo(jk) direct albedo for cloud
368) ! zrefdc(jk) diffuse albedo for clear
369) ! zrefdo(jk) diffuse albedo for cloud
370) ! ztrac(jk) direct transmittance for clear
371) ! ztrao(jk) direct transmittance for cloudy
372) ! ztradc(jk) diffuse transmittance for clear
373) ! ztrado(jk) diffuse transmittance for cloudy
374) !
375) ! zref(jk) direct reflectance
376) ! zrefd(jk) diffuse reflectance
377) ! ztra(jk) direct transmittance
378) ! ztrad(jk) diffuse transmittance
379) !
380) ! zdbtc(jk) clear direct beam transmittance
381) ! zdbto(jk) cloudy direct beam transmittance
382) ! zdbt(jk) layer mean direct beam transmittance
383) ! ztdbt(jk) total direct beam transmittance at levels
384)
385) ! Clear-sky
386) ! TOA direct beam
387) ztdbtc(1)=1.0_r8
388) ztdbtc_nodel(1)=1.0_r8
389) ! Surface values
390) zdbtc(klev+1) =0.0_r8
391) ztrac(klev+1) =0.0_r8
392) ztradc(klev+1)=0.0_r8
393) zrefc(klev+1) =palbp(ibm)
394) zrefdc(klev+1)=palbd(ibm)
395) zrupc(klev+1) =palbp(ibm)
396) zrupdc(klev+1)=palbd(ibm)
397)
398) ! Total sky
399) ! TOA direct beam
400) ztdbt(1)=1.0_r8
401) ztdbt_nodel(1)=1.0_r8
402) ! Surface values
403) zdbt(klev+1) =0.0_r8
404) ztra(klev+1) =0.0_r8
405) ztrad(klev+1)=0.0_r8
406) zref(klev+1) =palbp(ibm)
407) zrefd(klev+1)=palbd(ibm)
408) zrup(klev+1) =palbp(ibm)
409) zrupd(klev+1)=palbd(ibm)
410)
411)
412) ! Top of layer loop
413) do jk=1,klev
414)
415) ! Note: two-stream calculations proceed from top to bottom;
416) ! RRTMG_SW quantities are given bottom to top and are reversed here
417)
418) ikl=klev+1-jk
419)
420) ! Set logical flag to do REFTRA calculation
421) ! Do REFTRA for all clear layers
422) lrtchkclr(jk)=.true.
423)
424) ! Do REFTRA only for cloudy layers in profile, since already done for clear layers
425) lrtchkcld(jk)=.false.
426) lrtchkcld(jk)=(pclfr(ikl) > repclc)
427)
428) ! Clear-sky optical parameters - this section inactive
429) ! Original
430) ! ztauc(jk) = ztaur(ikl,iw) + ztaug(ikl,iw)
431) ! zomcc(jk) = ztaur(ikl,iw) / ztauc(jk)
432) ! zgcc(jk) = 0.0001_r8
433) ! Total sky optical parameters
434) ! ztauo(jk) = ztaur(ikl,iw) + ztaug(ikl,iw) + ptauc(ikl,ibm)
435) ! zomco(jk) = ptauc(ikl,ibm) * pomgc(ikl,ibm) + ztaur(ikl,iw)
436) ! zgco (jk) = (ptauc(ikl,ibm) * pomgc(ikl,ibm) * pasyc(ikl,ibm) + &
437) ! ztaur(ikl,iw) * 0.0001_r8) / zomco(jk)
438) ! zomco(jk) = zomco(jk) / ztauo(jk)
439)
440) ! Clear-sky optical parameters including aerosols
441) ztauc(jk) = ztaur(ikl,iw) + ztaug(ikl,iw) + ptaua(ikl,ibm)
442) zomcc(jk) = ztaur(ikl,iw) * 1.0_r8 + ptaua(ikl,ibm) * pomga(ikl,ibm)
443) zgcc(jk) = pasya(ikl,ibm) * pomga(ikl,ibm) * ptaua(ikl,ibm) / zomcc(jk)
444) zomcc(jk) = zomcc(jk) / ztauc(jk)
445)
446) ! Pre-delta-scaling clear and cloudy direct beam transmittance (must use 'orig', unscaled cloud OD)
447) ! \/\/\/ This block of code is only needed for unscaled direct beam calculation
448) if (idelm .eq. 0) then
449) !
450) zclear = 1.0_r8 - pclfr(ikl)
451) zcloud = pclfr(ikl)
452)
453) ! Clear
454) ! zdbtmc = exp(-ztauc(jk) / prmu0)
455)
456) ! Use exponential lookup table for transmittance, or expansion of exponential for low tau
457) ze1 = ztauc(jk) / prmu0
458) if (ze1 .le. od_lo) then
459) zdbtmc = 1._r8 - ze1 + 0.5_r8 * ze1 * ze1
460) else
461) tblind = ze1 / (bpade + ze1)
462) itind = tblint * tblind + 0.5_r8
463) zdbtmc = exp_tbl(itind)
464) endif
465)
466) zdbtc_nodel(jk) = zdbtmc
467) ztdbtc_nodel(jk+1) = zdbtc_nodel(jk) * ztdbtc_nodel(jk)
468)
469) ! Clear + Cloud
470) tauorig = ztauc(jk) + ptaucorig(ikl,ibm)
471) ! zdbtmo = exp(-tauorig / prmu0)
472)
473) ! Use exponential lookup table for transmittance, or expansion of exponential for low tau
474) ze1 = tauorig / prmu0
475) if (ze1 .le. od_lo) then
476) zdbtmo = 1._r8 - ze1 + 0.5_r8 * ze1 * ze1
477) else
478) tblind = ze1 / (bpade + ze1)
479) itind = tblint * tblind + 0.5_r8
480) zdbtmo = exp_tbl(itind)
481) endif
482)
483) zdbt_nodel(jk) = zclear * zdbtmc + zcloud * zdbtmo
484) ztdbt_nodel(jk+1) = zdbt_nodel(jk) * ztdbt_nodel(jk)
485)
486) endif
487) ! /\/\/\ Above code only needed for unscaled direct beam calculation
488)
489)
490) ! Delta scaling - clear
491) zf = zgcc(jk) * zgcc(jk)
492) zwf = zomcc(jk) * zf
493) ztauc(jk) = (1.0_r8 - zwf) * ztauc(jk)
494) zomcc(jk) = (zomcc(jk) - zwf) / (1.0_r8 - zwf)
495) zgcc (jk) = (zgcc(jk) - zf) / (1.0_r8 - zf)
496)
497) ! Total sky optical parameters (cloud properties already delta-scaled)
498) ! Use this code if cloud properties are derived in rrtmg_sw_cldprop
499) if (icpr .ge. 1) then
500) ztauo(jk) = ztauc(jk) + ptauc(ikl,ibm)
501) zomco(jk) = ztauc(jk) * zomcc(jk) + ptauc(ikl,ibm) * pomgc(ikl,ibm)
502) zgco (jk) = (ptauc(ikl,ibm) * pomgc(ikl,ibm) * pasyc(ikl,ibm) + &
503) ztauc(jk) * zomcc(jk) * zgcc(jk)) / zomco(jk)
504) zomco(jk) = zomco(jk) / ztauo(jk)
505)
506) ! Total sky optical parameters (if cloud properties not delta scaled)
507) ! Use this code if cloud properties are not derived in rrtmg_sw_cldprop
508) elseif (icpr .eq. 0) then
509) ztauo(jk) = ztaur(ikl,iw) + ztaug(ikl,iw) + ptaua(ikl,ibm) + ptauc(ikl,ibm)
510) zomco(jk) = ptaua(ikl,ibm) * pomga(ikl,ibm) + ptauc(ikl,ibm) * pomgc(ikl,ibm) + &
511) ztaur(ikl,iw) * 1.0_r8
512) zgco (jk) = (ptauc(ikl,ibm) * pomgc(ikl,ibm) * pasyc(ikl,ibm) + &
513) ptaua(ikl,ibm)*pomga(ikl,ibm)*pasya(ikl,ibm)) / zomco(jk)
514) zomco(jk) = zomco(jk) / ztauo(jk)
515)
516) ! Delta scaling - clouds
517) ! Use only if subroutine rrtmg_sw_cldprop is not used to get cloud properties and to apply delta scaling
518) zf = zgco(jk) * zgco(jk)
519) zwf = zomco(jk) * zf
520) ztauo(jk) = (1._r8 - zwf) * ztauo(jk)
521) zomco(jk) = (zomco(jk) - zwf) / (1.0_r8 - zwf)
522) zgco (jk) = (zgco(jk) - zf) / (1.0_r8 - zf)
523) endif
524)
525) ! End of layer loop
526) enddo
527)
528)
529) ! Clear sky reflectivities
530) call reftra_sw (klev, &
531) lrtchkclr, zgcc, prmu0, ztauc, zomcc, &
532) zrefc, zrefdc, ztrac, ztradc)
533)
534) ! Total sky reflectivities
535) call reftra_sw (klev, &
536) lrtchkcld, zgco, prmu0, ztauo, zomco, &
537) zrefo, zrefdo, ztrao, ztrado)
538)
539)
540) do jk=1,klev
541)
542) ! Combine clear and cloudy contributions for total sky
543) ikl = klev+1-jk
544) zclear = 1.0_r8 - pclfr(ikl)
545) zcloud = pclfr(ikl)
546)
547) zref(jk) = zclear*zrefc(jk) + zcloud*zrefo(jk)
548) zrefd(jk)= zclear*zrefdc(jk) + zcloud*zrefdo(jk)
549) ztra(jk) = zclear*ztrac(jk) + zcloud*ztrao(jk)
550) ztrad(jk)= zclear*ztradc(jk) + zcloud*ztrado(jk)
551)
552) ! Direct beam transmittance
553)
554) ! Clear
555) ! zdbtmc = exp(-ztauc(jk) / prmu0)
556)
557) ! Use exponential lookup table for transmittance, or expansion of
558) ! exponential for low tau
559) ze1 = ztauc(jk) / prmu0
560) if (ze1 .le. od_lo) then
561) zdbtmc = 1._r8 - ze1 + 0.5_r8 * ze1 * ze1
562) else
563) tblind = ze1 / (bpade + ze1)
564) itind = tblint * tblind + 0.5_r8
565) zdbtmc = exp_tbl(itind)
566) endif
567)
568) zdbtc(jk) = zdbtmc
569) ztdbtc(jk+1) = zdbtc(jk)*ztdbtc(jk)
570)
571) ! Clear + Cloud
572) ! zdbtmo = exp(-ztauo(jk) / prmu0)
573)
574) ! Use exponential lookup table for transmittance, or expansion of
575) ! exponential for low tau
576) ze1 = ztauo(jk) / prmu0
577) if (ze1 .le. od_lo) then
578) zdbtmo = 1._r8 - ze1 + 0.5_r8 * ze1 * ze1
579) else
580) tblind = ze1 / (bpade + ze1)
581) itind = tblint * tblind + 0.5_r8
582) zdbtmo = exp_tbl(itind)
583) endif
584)
585) zdbt(jk) = zclear*zdbtmc + zcloud*zdbtmo
586) ztdbt(jk+1) = zdbt(jk)*ztdbt(jk)
587)
588) enddo
589)
590) ! Vertical quadrature for clear-sky fluxes
591)
592) call vrtqdr_sw (klev, iw, &
593) zrefc, zrefdc, ztrac, ztradc, &
594) zdbtc, zrdndc, zrupc, zrupdc, ztdbtc, &
595) zcd, zcu)
596)
597) ! Vertical quadrature for cloudy fluxes
598)
599) call vrtqdr_sw (klev, iw, &
600) zref, zrefd, ztra, ztrad, &
601) zdbt, zrdnd, zrup, zrupd, ztdbt, &
602) zfd, zfu)
603)
604) ! Upwelling and downwelling fluxes at levels
605) ! Two-stream calculations go from top to bottom;
606) ! layer indexing is reversed to go bottom to top for output arrays
607)
608) do jk=1,klev+1
609) ikl=klev+2-jk
610)
611) ! Accumulate spectral fluxes over bands - inactive
612) ! zbbfu(ikl) = zbbfu(ikl) + zincflx(iw)*zfu(jk,iw)
613) ! zbbfd(ikl) = zbbfd(ikl) + zincflx(iw)*zfd(jk,iw)
614) ! zbbcu(ikl) = zbbcu(ikl) + zincflx(iw)*zcu(jk,iw)
615) ! zbbcd(ikl) = zbbcd(ikl) + zincflx(iw)*zcd(jk,iw)
616) ! zbbfddir(ikl) = zbbfddir(ikl) + zincflx(iw)*ztdbt_nodel(jk)
617) ! zbbcddir(ikl) = zbbcddir(ikl) + zincflx(iw)*ztdbtc_nodel(jk)
618)
619) ! Accumulate spectral fluxes over whole spectrum
620) pbbfu(ikl) = pbbfu(ikl) + zincflx(iw)*zfu(jk,iw)
621) pbbfd(ikl) = pbbfd(ikl) + zincflx(iw)*zfd(jk,iw)
622) pbbcu(ikl) = pbbcu(ikl) + zincflx(iw)*zcu(jk,iw)
623) pbbcd(ikl) = pbbcd(ikl) + zincflx(iw)*zcd(jk,iw)
624) if (idelm .eq. 0) then
625) pbbfddir(ikl) = pbbfddir(ikl) + zincflx(iw)*ztdbt_nodel(jk)
626) pbbcddir(ikl) = pbbcddir(ikl) + zincflx(iw)*ztdbtc_nodel(jk)
627) elseif (idelm .eq. 1) then
628) pbbfddir(ikl) = pbbfddir(ikl) + zincflx(iw)*ztdbt(jk)
629) pbbcddir(ikl) = pbbcddir(ikl) + zincflx(iw)*ztdbtc(jk)
630) endif
631)
632) ! Accumulate direct fluxes for UV/visible bands
633) if (ibm >= 10 .and. ibm <= 13) then
634) puvcd(ikl) = puvcd(ikl) + zincflx(iw)*zcd(jk,iw)
635) puvfd(ikl) = puvfd(ikl) + zincflx(iw)*zfd(jk,iw)
636) if (idelm .eq. 0) then
637) puvfddir(ikl) = puvfddir(ikl) + zincflx(iw)*ztdbt_nodel(jk)
638) puvcddir(ikl) = puvcddir(ikl) + zincflx(iw)*ztdbtc_nodel(jk)
639) elseif (idelm .eq. 1) then
640) puvfddir(ikl) = puvfddir(ikl) + zincflx(iw)*ztdbt(jk)
641) puvcddir(ikl) = puvcddir(ikl) + zincflx(iw)*ztdbtc(jk)
642) endif
643) ! Accumulate direct fluxes for near-IR bands
644) else if (ibm == 14 .or. ibm <= 9) then
645) pnicd(ikl) = pnicd(ikl) + zincflx(iw)*zcd(jk,iw)
646) pnifd(ikl) = pnifd(ikl) + zincflx(iw)*zfd(jk,iw)
647) if (idelm .eq. 0) then
648) pnifddir(ikl) = pnifddir(ikl) + zincflx(iw)*ztdbt_nodel(jk)
649) pnicddir(ikl) = pnicddir(ikl) + zincflx(iw)*ztdbtc_nodel(jk)
650) elseif (idelm .eq. 1) then
651) pnifddir(ikl) = pnifddir(ikl) + zincflx(iw)*ztdbt(jk)
652) pnicddir(ikl) = pnicddir(ikl) + zincflx(iw)*ztdbtc(jk)
653) endif
654) ! Added for net near-IR flux diagnostic
655) pnicu(ikl) = pnicu(ikl) + zincflx(iw)*zcu(jk,iw)
656) pnifu(ikl) = pnifu(ikl) + zincflx(iw)*zfu(jk,iw)
657) endif
658)
659) enddo
660)
661) ! End loop on jg, g-point interval
662) enddo
663)
664) ! End loop on jb, spectral band
665) enddo
666)
667) end subroutine spcvrt_sw
668)
669) end module rrtmg_sw_spcvrt
670)
671)
