304 module module_radsw_main
314 use mersenne_twister
, only : random_setseed, random_number, &
324 character(40),
parameter :: &
325 & VTAGSW=
'NCEP SW v5.1 Nov 2012 -RRTMG-SW v3.8 '
337 real (kind=kind_phys),
parameter :: eps = 1.0e-6
338 real (kind=kind_phys),
parameter :: oneminus= 1.0 - eps
340 real (kind=kind_phys),
parameter ::
bpade = 1.0/0.278
341 real (kind=kind_phys),
parameter :: stpfac = 296.0/1013.0
342 real (kind=kind_phys),
parameter :: ftiny = 1.0e-12
344 real (kind=kind_phys),
parameter ::
s0 = 1368.22
346 real (kind=kind_phys),
parameter :: f_zero = 0.0
347 real (kind=kind_phys),
parameter :: f_one = 1.0
354 integer,
dimension(nblow:nbhgh) :: nspa, nspb
356 integer,
dimension(nblow:nbhgh) ::
idxsfc
358 integer,
dimension(nblow:nbhgh) ::
idxebc
360 data nspa(:) / 9, 9, 9, 9, 1, 9, 9, 1, 9, 1, 0, 1, 9, 1 /
361 data nspb(:) / 1, 5, 1, 1, 1, 5, 1, 0, 1, 0, 0, 1, 5, 1 /
364 data idxsfc(:) / 1, 1, 1, 1, 1, 1, 1, 1, 0, 2, 2, 2, 2, 1 /
365 data idxebc(:) / 5, 5, 4, 4, 3, 3, 2, 2, 1, 1, 1, 1, 1, 5 /
381 integer,
parameter ::
nuvb = 27
384 logical :: lhswb = .false.
385 logical :: lhsw0 = .false.
386 logical :: lflxprf= .false.
387 logical :: lfdncmp= .false.
484 & clouds,icseed,aerosols,sfcalb,
488 & hsw0,hswb,flxprf,fdncmp
667 integer,
intent(in) :: npts, nlay, nlp1, NDAY
669 integer,
dimension(:),
intent(in) :: idxday, icseed
671 logical,
intent(in) :: lprnt
673 real (kind=kind_phys),
dimension(npts,nlp1),
intent(in) :: &
675 real (kind=kind_phys),
dimension(npts,nlay),
intent(in) :: &
676 & plyr, tlyr, qlyr, olyr
677 real (kind=kind_phys),
dimension(npts,4),
intent(in) :: sfcalb
679 real (kind=kind_phys),
dimension(npts,nlay,9),
intent(in):: gasvmr
680 real (kind=kind_phys),
dimension(npts,nlay,9),
intent(in):: clouds
681 real (kind=kind_phys),
dimension(npts,nlay,nbdsw,3),
intent(in):: &
684 real (kind=kind_phys),
intent(in) :: cosz(npts), solcon
687 real (kind=kind_phys),
dimension(npts,nlay),
intent(out) :: hswc
689 type(
topfsw_type),
dimension(npts),
intent(out) :: topflx
690 type(
sfcfsw_type),
dimension(npts),
intent(out) :: sfcflx
693 real (kind=kind_phys),
dimension(npts,nlay,nbdsw),
optional, &
694 & intent(out) :: hswb
696 real (kind=kind_phys),
dimension(npts,nlay),
optional, &
697 & intent(out) :: hsw0
698 type (
profsw_type),
dimension(npts,nlp1),
optional, &
699 & intent(out) :: flxprf
701 & intent(out) :: fdncmp
704 real (kind=kind_phys),
dimension(nlay,ngptsw) :: cldfmc, &
706 real (kind=kind_phys),
dimension(nlp1,nbdsw):: fxupc, fxdnc, &
709 real (kind=kind_phys),
dimension(nlay,nbdsw) :: &
710 & tauae, ssaae, asyae, taucw, ssacw, asycw
712 real (kind=kind_phys),
dimension(ngptsw) :: sfluxzen
714 real (kind=kind_phys),
dimension(nlay) :: cldfrc, delp, &
715 & pavel, tavel, coldry, colmol, h2ovmr, o3vmr, temcol,
720 real (kind=kind_phys),
dimension(nlp1) :: fnet, flxdc, flxuc, &
723 real (kind=kind_phys),
dimension(2) :: albbm, albdf, sfbmc, &
724 & sfbm0, sfdfc, sfdf0
726 real (kind=kind_phys) :: cosz1, sntz1, tem0, tem1, tem2, s0fac, &
727 & ssolar, zcf0, zcf1, ftoau0, ftoauc, ftoadc,
732 real (kind=kind_phys) :: colamt(nlay,
maxgas)
734 integer,
dimension(npts) :: ipseed
735 integer,
dimension(nlay) :: indfor, indself, jp, jt, jt1
737 integer :: i, ib, ipt, j1, k, kk, laytrop, mb
742 lhswb =
present ( hswb )
743 lhsw0 =
present ( hsw0 )
744 lflxprf=
present ( flxprf )
745 lfdncmp=
present ( fdncmp )
757 sfcflx =
sfcfsw_type( f_zero, f_zero, f_zero, f_zero )
761 flxprf =
profsw_type( f_zero, f_zero, f_zero, f_zero )
765 fdncmp =
cmpfsw_type(f_zero,f_zero,f_zero,f_zero,f_zero,f_zero)
785 ipseed(i) = icseed(i)
790 write(0,*)
' In radsw, isubcsw, ipsdsw0,ipseed =',
796 lab_do_ipt :
do ipt = 1, nday
801 sntz1 = f_one / cosz(j1)
802 ssolar = s0fac * cosz(j1)
805 albbm(1) = sfcalb(j1,1)
806 albdf(1) = sfcalb(j1,2)
807 albbm(2) = sfcalb(j1,3)
808 albdf(2) = sfcalb(j1,4)
820 pavel(k) = plyr(j1,kk)
821 tavel(k) = tlyr(j1,kk)
822 delp(k) = plvl(j1,kk+1) - plvl(j1,kk)
834 h2ovmr(k)= max(f_zero,qlyr(j1,kk)*amdw/(f_one-qlyr(j1,kk)))
835 o3vmr(k)= max(f_zero,olyr(j1,kk)*amdo3)
838 coldry(k) = tem2 * delp(k) / (tem1*tem0*(f_one + h2ovmr(k)))
839 temcol(k) = 1.0e-12 * coldry(k)
841 colamt(k,1) = max(f_zero, coldry(k)*h2ovmr(k))
842 colamt(k,2) = max(temcol(k), coldry(k)*gasvmr(j1,kk,1))
843 colamt(k,3) = max(f_zero, coldry(k)*o3vmr(k))
844 colmol(k) = coldry(k) + colamt(k,1)
853 colamt(k,4) = max(temcol(k), coldry(k)*gasvmr(j1,kk,2))
854 colamt(k,5) = max(temcol(k), coldry(k)*gasvmr(j1,kk,3))
855 colamt(k,6) = max(temcol(k), coldry(k)*gasvmr(j1,kk,4))
860 colamt(k,4) = temcol(k)
861 colamt(k,5) = temcol(k)
862 colamt(k,6) = temcol(k)
872 tauae(k,ib) = aerosols(j1,kk,ib,1)
873 ssaae(k,ib) = aerosols(j1,kk,ib,2)
874 asyae(k,ib) = aerosols(j1,kk,ib,3)
882 cfrac(k) = clouds(j1,kk,1)
883 cliqp(k) = clouds(j1,kk,2)
884 reliq(k) = clouds(j1,kk,3)
885 cicep(k) = clouds(j1,kk,4)
886 reice(k) = clouds(j1,kk,5)
887 cdat1(k) = clouds(j1,kk,6)
888 cdat2(k) = clouds(j1,kk,7)
889 cdat3(k) = clouds(j1,kk,8)
890 cdat4(k) = clouds(j1,kk,9)
895 cfrac(k) = clouds(j1,kk,1)
896 cdat1(k) = clouds(j1,kk,2)
897 cdat2(k) = clouds(j1,kk,3)
898 cdat3(k) = clouds(j1,kk,4)
908 pavel(k) = plyr(j1,k)
909 tavel(k) = tlyr(j1,k)
910 delp(k) = plvl(j1,k) - plvl(j1,k+1)
918 h2ovmr(k)= max(f_zero,qlyr(j1,k)*amdw/(f_one-qlyr(j1,k)))
919 o3vmr(k)= max(f_zero,olyr(j1,k)*amdo3)
922 coldry(k) = tem2 * delp(k) / (tem1*tem0*(f_one + h2ovmr(k)))
923 temcol(k) = 1.0e-12 * coldry(k)
925 colamt(k,1) = max(f_zero, coldry(k)*h2ovmr(k))
926 colamt(k,2) = max(temcol(k), coldry(k)*gasvmr(j1,k,1))
927 colamt(k,3) = max(f_zero, coldry(k)*o3vmr(k))
928 colmol(k) = coldry(k) + colamt(k,1)
934 write(0,*)
' pavel=',pavel
935 write(0,*)
' tavel=',tavel
936 write(0,*)
' delp=',delp
937 write(0,*)
' h2ovmr=',h2ovmr*1000
938 write(0,*)
' o3vmr=',o3vmr*1000000
947 colamt(k,4) = max(temcol(k), coldry(k)*gasvmr(j1,k,2))
948 colamt(k,5) = max(temcol(k), coldry(k)*gasvmr(j1,k,3))
949 colamt(k,6) = max(temcol(k), coldry(k)*gasvmr(j1,k,4))
954 colamt(k,4) = temcol(k)
955 colamt(k,5) = temcol(k)
956 colamt(k,6) = temcol(k)
965 tauae(k,ib) = aerosols(j1,k,ib,1)
966 ssaae(k,ib) = aerosols(j1,k,ib,2)
967 asyae(k,ib) = aerosols(j1,k,ib,3)
973 cfrac(k) = clouds(j1,k,1)
974 cliqp(k) = clouds(j1,k,2)
975 reliq(k) = clouds(j1,k,3)
976 cicep(k) = clouds(j1,k,4)
977 reice(k) = clouds(j1,k,5)
978 cdat1(k) = clouds(j1,k,6)
979 cdat2(k) = clouds(j1,k,7)
980 cdat3(k) = clouds(j1,k,8)
981 cdat4(k) = clouds(j1,k,9)
985 cfrac(k) = clouds(j1,k,1)
986 cdat1(k) = clouds(j1,k,2)
987 cdat2(k) = clouds(j1,k,3)
988 cdat3(k) = clouds(j1,k,4)
1003 zcf0 = zcf0 * (f_one - cfrac(k))
1005 else if (
iovrsw == 1)
then
1007 if (cfrac(k) > ftiny)
then
1008 zcf1 = min( zcf1, f_one-cfrac(k) )
1009 elseif (zcf1 < f_one)
then
1015 else if (
iovrsw == 2)
then
1017 zcf0 = min( zcf0, f_one-cfrac(k) )
1021 if (zcf0 <= ftiny) zcf0 = f_zero
1022 if (zcf0 > oneminus) zcf0 = f_one
1028 if (zcf1 > f_zero)
then
1032 & ( cfrac,cliqp,reliq,cicep,reice,cdat1,cdat2,cdat3,cdat4,
1035 & taucw, ssacw, asycw, cldfrc, cldfmc
1054 & ( pavel,tavel,h2ovmr, nlay,nlp1,
1056 & laytrop,jp,jt,jt1,fac00,fac01,fac10,fac11,
1064 & ( colamt,colmol,fac00,fac01,fac10,fac11,jp,jt,jt1,laytrop,
1067 & sfluxzen, taug, taur
1080 & ( ssolar,cosz1,sntz1,albbm,albdf,sfluxzen,cldfrc,
1084 & fxupc,fxdnc,fxup0,fxdn0,
1093 & ( ssolar,cosz1,sntz1,albbm,albdf,sfluxzen,cldfmc,
1097 & fxupc,fxdnc,fxup0,fxdn0,
1112 flxuc(k) = flxuc(k) + fxupc(k,ib)
1113 flxdc(k) = flxdc(k) + fxdnc(k,ib)
1119 if ( lhsw0 .or. lflxprf )
then
1125 flxu0(k) = flxu0(k) + fxup0(k,ib)
1126 flxd0(k) = flxd0(k) + fxdn0(k,ib)
1139 fdncmp(j1)%uvbf0 = suvbf0
1140 fdncmp(j1)%uvbfc = suvbfc
1143 fdncmp(j1)%nirbm = sfbmc(1)
1144 fdncmp(j1)%nirdf = sfdfc(1)
1145 fdncmp(j1)%visbm = sfbmc(2)
1146 fdncmp(j1)%visdf = sfdfc(2)
1151 topflx(j1)%upfxc = ftoauc
1152 topflx(j1)%dnfxc = ftoadc
1153 topflx(j1)%upfx0 = ftoau0
1155 sfcflx(j1)%upfxc = fsfcuc
1156 sfcflx(j1)%dnfxc = fsfcdc
1157 sfcflx(j1)%upfx0 = fsfcu0
1158 sfcflx(j1)%dnfx0 = fsfcd0
1164 fnet(1) = flxdc(1) - flxuc(1)
1168 fnet(k) = flxdc(k) - flxuc(k)
1169 hswc(j1,kk) = (fnet(k)-fnet(k-1)) * rfdelp(k-1)
1177 flxprf(j1,kk)%upfxc = flxuc(k)
1178 flxprf(j1,kk)%dnfxc = flxdc(k)
1179 flxprf(j1,kk)%upfx0 = flxu0(k)
1180 flxprf(j1,kk)%dnfx0 = flxd0(k)
1187 fnet(1) = flxd0(1) - flxu0(1)
1191 fnet(k) = flxd0(k) - flxu0(k)
1192 hsw0(j1,kk) = (fnet(k)-fnet(k-1)) * rfdelp(k-1)
1200 fnet(1) = fxdnc(1,mb) - fxupc(1,mb)
1204 fnet(k) = fxdnc(k,mb) - fxupc(k,mb)
1205 hswb(j1,kk,mb) = (fnet(k) - fnet(k-1)) * rfdelp(k-1)
1214 fnet(1) = flxdc(1) - flxuc(1)
1217 fnet(k) = flxdc(k) - flxuc(k)
1218 hswc(j1,k-1) = (fnet(k)-fnet(k-1)) * rfdelp(k-1)
1225 flxprf(j1,k)%upfxc = flxuc(k)
1226 flxprf(j1,k)%dnfxc = flxdc(k)
1227 flxprf(j1,k)%upfx0 = flxu0(k)
1228 flxprf(j1,k)%dnfx0 = flxd0(k)
1235 fnet(1) = flxd0(1) - flxu0(1)
1238 fnet(k) = flxd0(k) - flxu0(k)
1239 hsw0(j1,k-1) = (fnet(k)-fnet(k-1)) * rfdelp(k-1)
1247 fnet(1) = fxdnc(1,mb) - fxupc(1,mb)
1250 fnet(k+1) = fxdnc(k+1,mb) - fxupc(k+1,mb)
1251 hswb(j1,k,mb) = (fnet(k+1) - fnet(k)) * rfdelp(k)
1262 end subroutine swrad
1327 integer,
intent(in) :: me
1332 real (kind=kind_phys),
parameter :: expeps = 1.e-20
1336 real (kind=kind_phys) :: tfn, tau
1342 print *,
' *** Error in specification of cloud overlap flag',
1343 ' IOVRSW=',
iovrsw,
' in RSWINIT !!'
1348 print *,
' - Using AER Shortwave Radiation, Version: ',vtagsw
1351 print *,
' --- Delta-eddington 2-stream transfer scheme'
1353 print *,
' --- PIFM 2-stream transfer scheme'
1355 print *,
' --- Discrete ordinates 2-stream transfer scheme'
1359 print *,
' --- Rare gases absorption is NOT included in SW'
1361 print *,
' --- Include rare gases N2O, CH4, O2, absorptions',
1366 print *,
' --- Using standard grid average clouds, no ',
1367 'sub-column clouds approximation applied'
1369 print *,
' --- Using MCICA sub-colum clouds approximation ',
1370 'with a prescribed sequence of permutation seeds'
1372 print *,
' --- Using MCICA sub-colum clouds approximation ',
1373 'with provided input array of permutation seeds'
1375 print *,
' *** Error in specification of sub-column cloud ',
1376 ' control flag isubcsw =',
isubcsw,
' !!'
1385 print *,
' *** Model cloud scheme inconsistent with SW',
1386 ' radiation cloud radiative property setup !!'
1412 tfn = float(i) / float(
ntbmx-i)
1459 & cf1, nlay, ipseed,
1541 integer,
intent(in) :: nlay, ipseed
1542 real (kind=kind_phys),
intent(in) :: cf1
1544 real (kind=kind_phys),
dimension(nlay),
intent(in) :: cliqp, &
1545 & reliq, cicep, reice, cdat1, cdat2, cdat3, cdat4, cfrac
1548 real (kind=kind_phys),
dimension(nlay,ngptsw),
intent(out) :: &
1550 real (kind=kind_phys),
dimension(nlay,nbdsw),
intent(out) :: &
1551 & taucw, ssacw, asycw
1552 real (kind=kind_phys),
dimension(nlay),
intent(out) :: cldfrc
1555 real (kind=kind_phys),
dimension(nblow:nbhgh) :: tauliq, tauice, &
1556 & ssaliq, ssaice, ssaran, ssasnw, asyliq, asyice,
1558 real (kind=kind_phys),
dimension(nlay) :: cldf
1560 real (kind=kind_phys) :: dgeice, factor, fint, tauran, tausnw, &
1561 & cldliq, refliq, cldice, refice, cldran, cldsnw, refsnw,
1565 logical :: lcloudy(nlay,
ngptsw)
1566 integer :: ia, ib, ig, jb, k, index
1573 taucw(k,ib) = f_zero
1575 asycw(k,ib) = f_zero
1581 lab_if_iswcliq :
if (
iswcliq > 0)
then
1583 lab_do_k :
do k = 1, nlay
1584 lab_if_cld :
if (cfrac(k) > ftiny)
then
1602 dgesnw = 1.0315 * refsnw
1604 tauran = cldran * a0r
1612 if (cldsnw>f_zero .and. refsnw>10.0_kind_phys)
then
1614 tausnw = cldsnw*1.09087*(a0s + a1s/dgesnw)
1620 ssaran(ib) = tauran * (f_one - b0r(ib))
1621 ssasnw(ib) = tausnw * (f_one - (b0s(ib)+b1s(ib)*dgesnw))
1622 asyran(ib) = ssaran(ib) * c0r(ib)
1623 asysnw(ib) = ssasnw(ib) * c0s(ib)
1633 if ( cldliq <= f_zero )
then
1641 factor = refliq - 1.5
1642 index = max( 1, min( 57, int( factor ) ))
1643 fint = factor - float(index)
1646 extcoliq = max(f_zero, extliq1(index,ib)
1655 tauliq(ib) = cldliq * extcoliq
1656 ssaliq(ib) = tauliq(ib) * ssacoliq
1657 asyliq(ib) = ssaliq(ib) * asycoliq
1664 if ( cldice <= f_zero )
then
1676 refice = min(130.0_kind_phys,max(13.0_kind_phys,refice))
1681 extcoice = max(f_zero, abari(ia)+bbari(ia)/refice )
1682 ssacoice = max(f_zero, min(f_one,
1688 tauice(ib) = cldice * extcoice
1689 ssaice(ib) = tauice(ib) * ssacoice
1690 asyice(ib) = ssaice(ib) * asycoice
1696 refice = min(131.0_kind_phys,max(5.0_kind_phys,refice))
1698 factor = (refice - 2.0) / 3.0
1699 index = max( 1, min( 42, int( factor ) ))
1700 fint = factor - float(index)
1703 extcoice = max(f_zero, extice2(index,ib)
1711 tauice(ib) = cldice * extcoice
1712 ssaice(ib) = tauice(ib) * ssacoice
1713 asyice(ib) = ssaice(ib) * asycoice
1720 dgeice = max( 5.0, min( 140.0, 1.0315*refice ))
1722 factor = (dgeice - 2.0) / 3.0
1723 index = max( 1, min( 45, int( factor ) ))
1724 fint = factor - float(index)
1727 extcoice = max(f_zero, extice3(index,ib)
1737 tauice(ib) = cldice * extcoice
1738 ssaice(ib) = tauice(ib) * ssacoice
1739 asyice(ib) = ssaice(ib) * asycoice
1747 taucw(k,ib) = tauliq(jb)+tauice(jb)+tauran+tausnw
1748 ssacw(k,ib) = ssaliq(jb)+ssaice(jb)+ssaran(jb)+ssasnw(jb)
1749 asycw(k,ib) = asyliq(jb)+asyice(jb)+asyran(jb)+asysnw(jb)
1758 if (cfrac(k) > ftiny)
then
1760 taucw(k,ib) = cdat1(k)
1761 ssacw(k,ib) = cdat1(k) * cdat2(k)
1762 asycw(k,ib) = ssacw(k,ib) * cdat3(k)
1767 endif lab_if_iswcliq
1775 where (cldf(:) < ftiny)
1783 & ( cldf, nlay, ipseed,
1790 if ( lcloudy(k,ig) )
then
1791 cldfmc(k,ig) = f_one
1793 cldfmc(k,ig) = f_zero
1801 cldfrc(k) = cfrac(k) / cf1
1846 integer,
intent(in) :: nlay, ipseed
1848 real (kind=kind_phys),
dimension(nlay),
intent(in) :: cldf
1851 logical,
dimension(nlay,ngptsw),
intent(out):: lcloudy
1854 real (kind=kind_phys) :: cdfunc(nlay,
ngptsw), tem1, &
1855 & rand2d(nlay*ngptsw), rand1d(ngptsw)
1857 type(random_stat) :: stat
1887 cdfunc(k,n) = rand2d(k1)
1902 cdfunc(k,n) = rand2d(k1)
1914 tem1 = f_one - cldf(k1)
1917 if ( cdfunc(k1,n) > tem1 )
then
1918 cdfunc(k,n) = cdfunc(k1,n)
1920 cdfunc(k,n) = cdfunc(k,n) * tem1
1963 tem1 = f_one - cldf(k)
1966 lcloudy(k,n) = cdfunc(k,n) >= tem1
2002 & laytrop,jp,jt,jt1,fac00,fac01,fac10,fac11,
2003 & selffac,selffrac,indself,forfac,forfrac,indfor
2043 integer,
intent(in) :: nlay, nlp1
2045 real (kind=kind_phys),
dimension(:),
intent(in) :: pavel, tavel, &
2049 integer,
dimension(nlay),
intent(out) :: indself, indfor, &
2051 integer,
intent(out) :: laytrop
2053 real (kind=kind_phys),
dimension(nlay),
intent(out) :: fac00, &
2054 & fac01, fac10, fac11, selffac, selffrac, forfac, forfrac
2057 real (kind=kind_phys) :: plog, fp, fp1, ft, ft1, tem1, tem2
2059 integer :: i, k, jp1
2067 forfac(k) = pavel(k)*stpfac / (tavel(k)*(f_one + h2ovmr(k)))
2074 plog = log(pavel(k))
2075 jp(k) = max(1, min(58, int(36.0 - 5.0*(plog+0.04)) ))
2077 fp = 5.0 * (
preflog(jp(k)) - plog)
2086 tem1 = (tavel(k) - tref(jp(k))) / 15.0
2087 tem2 = (tavel(k) - tref(jp1 )) / 15.0
2088 jt(k) = max(1, min(4, int(3.0 + tem1) ))
2089 jt1(k) = max(1, min(4, int(3.0 + tem2) ))
2090 ft = tem1 - float(jt(k) - 3)
2091 ft1 = tem2 - float(jt1(k) - 3)
2102 fac00(k) = fp1 * (f_one - ft)
2104 fac01(k) = fp * (f_one - ft1)
2109 if ( plog > 4.56 )
then
2116 tem1 = (332.0 - tavel(k)) / 36.0
2117 indfor(k) = min(2, max(1, int(tem1)))
2118 forfrac(k) = tem1 - float(indfor(k))
2123 tem2 = (tavel(k) - 188.0) / 7.2
2124 indself(k) = min(9, max(1, int(tem2)-7))
2125 selffrac(k) = tem2 - float(indself(k) + 7)
2126 selffac(k) = h2ovmr(k) * forfac(k)
2133 tem1 = (tavel(k) - 188.0) / 36.0
2135 forfrac(k) = tem1 - f_one
2138 selffrac(k) = f_zero
2192 & cf1,cf0,taug,taur,tauae,ssaae,asyae,taucw,ssacw,asycw,
2195 & ftoauc,ftoau0,ftoadc,fsfcuc,fsfcu0,fsfcdc,fsfcd0,
2290 real (kind=kind_phys),
parameter :: zcrit = 0.9999995
2291 real (kind=kind_phys),
parameter :: zsr3 = sqrt(3.0)
2292 real (kind=kind_phys),
parameter :: od_lo = 0.06
2293 real (kind=kind_phys),
parameter :: eps1 = 1.0e-8
2296 integer,
intent(in) :: nlay, nlp1
2298 real (kind=kind_phys),
dimension(nlay,ngptsw),
intent(in) :: &
2300 real (kind=kind_phys),
dimension(nlay,nbdsw),
intent(in) :: &
2301 & taucw, ssacw, asycw, tauae, ssaae, asyae
2303 real (kind=kind_phys),
dimension(ngptsw),
intent(in) :: sfluxzen
2304 real (kind=kind_phys),
dimension(nlay),
intent(in) :: cldfrc
2306 real (kind=kind_phys),
dimension(2),
intent(in) :: albbm, albdf
2308 real (kind=kind_phys),
intent(in) :: cosz, sntz, cf1, cf0, ssolar
2311 real (kind=kind_phys),
dimension(nlp1,nbdsw),
intent(out) :: &
2312 & fxupc, fxdnc, fxup0, fxdn0
2314 real (kind=kind_phys),
dimension(2),
intent(out) :: sfbmc, sfdfc, &
2317 real (kind=kind_phys),
intent(out) :: suvbfc, suvbf0, ftoadc, &
2318 & ftoauc, ftoau0, fsfcuc, fsfcu0, fsfcdc, fsfcd0
2321 real (kind=kind_phys),
dimension(nlay) :: ztaus, zssas, zasys, &
2324 real (kind=kind_phys),
dimension(nlp1) :: zrefb, zrefd, ztrab, &
2325 & ztrad, ztdbt, zldbt, zfu, zfd
2327 real (kind=kind_phys) :: ztau1, zssa1, zasy1, ztau0, zssa0, &
2328 & zasy0, zasy3, zssaw, zasyw, zgam1, zgam2, zgam3, zgam4,
2335 integer :: ib, ibd, jb, jg, k, kp, itind
2342 fxdnc(k,ib) = f_zero
2343 fxupc(k,ib) = f_zero
2344 fxdn0(k,ib) = f_zero
2345 fxup0(k,ib) = f_zero
2373 lab_do_jg :
do jg = 1,
ngptsw
2379 zsolar = ssolar * sfluxzen(jg)
2388 zrefb(1) = albbm(ibd)
2389 zrefd(1) = albdf(ibd)
2391 zrefb(1) = 0.5 * (albbm(1) + albbm(2))
2392 zrefd(1) = 0.5 * (albdf(1) + albdf(2))
2411 ztau0 = max( ftiny, taur(k,jg)+taug(k,jg)+tauae(k,ib) )
2412 zssa0 = taur(k,jg) + tauae(k,ib)*ssaae(k,ib)
2413 zasy0 = asyae(k,ib)*ssaae(k,ib)*tauae(k,ib)
2414 zssaw = min( oneminus, zssa0 / ztau0 )
2415 zasyw = zasy0 / max( ftiny, zssa0 )
2426 ztau1 = (f_one - za2) * ztau0
2427 zssa1 = (zssaw - za2) / (f_one - za2)
2429 zasy1 = zasyw / (f_one + zasyw)
2430 zasy3 = 0.75 * zasy1
2434 zgam1 = 1.75 - zssa1 * (f_one + zasy3)
2435 zgam2 =-0.25 + zssa1 * (f_one - zasy3)
2436 zgam3 = 0.5 - zasy3 * cosz
2438 zgam1 = 2.0 - zssa1 * (1.25 + zasy3)
2439 zgam2 = 0.75* zssa1 * (f_one- zasy1)
2440 zgam3 = 0.5 - zasy3 * cosz
2442 zgam1 = zsr3 * (2.0 - zssa1 * (1.0 + zasy1)) * 0.5
2443 zgam2 = zsr3 * zssa1 * (1.0 - zasy1) * 0.5
2444 zgam3 = (1.0 - zsr3 * zasy1 * cosz) * 0.5
2446 zgam4 = f_one - zgam3
2450 if ( zssaw >= zcrit )
then
2451 za1 = zgam1 * cosz - zgam3
2457 zb1 = min( ztau1*sntz , 500.0 )
2458 if ( zb1 <= od_lo )
then
2459 zb2 = f_one - zb1 + 0.5*zb1*zb1
2461 ftind = zb1 / (
bpade + zb1)
2462 itind = ftind*
ntbmx + 0.5
2467 zrefb(kp) = max(f_zero, min(f_one,
2472 zrefd(kp) = max(f_zero, min(f_one, za2/(f_one + za2) ))
2473 ztrad(kp) = max(f_zero, min(f_one, f_one-zrefd(kp) ))
2476 za1 = zgam1*zgam4 + zgam2*zgam3
2477 za2 = zgam1*zgam3 + zgam2*zgam4
2478 zrk = sqrt( (zgam1 - zgam2) * (zgam1 + zgam2) )
2484 zrpp = f_one - zrp*zrp
2489 zr1 = zrm1 * (za2 + zrkg3)
2490 zr2 = zrp1 * (za2 - zrkg3)
2491 zr3 = zrk2 * (zgam3 - za2*cosz)
2493 zr5 = zrpp * (zrk - zgam1)
2495 zt1 = zrp1 * (za1 + zrkg4)
2496 zt2 = zrm1 * (za1 - zrkg4)
2497 zt3 = zrk2 * (zgam4 + za1*cosz)
2502 zb1 = min( zrk*ztau1, 500.0 )
2503 if ( zb1 <= od_lo )
then
2504 zexm1 = f_one - zb1 + 0.5*zb1*zb1
2506 ftind = zb1 / (
bpade + zb1)
2507 itind = ftind*
ntbmx + 0.5
2510 zexp1 = f_one / zexm1
2512 zb2 = min( sntz*ztau1, 500.0 )
2513 if ( zb2 <= od_lo )
then
2514 zexm2 = f_one - zb2 + 0.5*zb2*zb2
2516 ftind = zb2 / (
bpade + zb2)
2517 itind = ftind*
ntbmx + 0.5
2520 zexp2 = f_one / zexm2
2521 ze1r45 = zr4*zexp1 + zr5*zexm1
2524 if (ze1r45>=-eps1 .and. ze1r45<=eps1)
then
2528 zden1 = zssa1 / ze1r45
2529 zrefb(kp) = max(f_zero, min(f_one,
2536 zden1 = zr4 / (ze1r45 * zrkg1)
2537 zrefd(kp) = max(f_zero, min(f_one,
2547 if ( zr1 <= od_lo )
then
2548 zexp3 = f_one - zr1 + 0.5*zr1*zr1
2550 ftind = zr1 / (
bpade + zr1)
2551 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
2555 ztdbt(k) = zexp3 * ztdbt(kp)
2562 if ( zr1 <= od_lo )
then
2563 zexp4 = f_one - zr1 + 0.5*zr1*zr1
2565 ftind = zr1 / (
bpade + zr1)
2566 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
2571 ztdbt0 = zexp4 * ztdbt0
2576 & ( zrefb,zrefd,ztrab,ztrad,zldbt,ztdbt,
2584 fxup0(k,ib) = fxup0(k,ib) + zsolar*zfu(k)
2585 fxdn0(k,ib) = fxdn0(k,ib) + zsolar*zfd(k)
2590 zb2 = zsolar*(zfd(1) - ztdbt0)
2593 sfbm0(ibd) = sfbm0(ibd) + zb1
2594 sfdf0(ibd) = sfdf0(ibd) + zb2
2598 sfbm0(1) = sfbm0(1) + zf1
2599 sfdf0(1) = sfdf0(1) + zf2
2600 sfbm0(2) = sfbm0(2) + zf1
2601 sfdf0(2) = sfdf0(2) + zf2
2616 if ( cf1 > eps )
then
2624 zc0 = f_one - cldfrc(k)
2626 if ( zc1 > ftiny )
then
2628 ztau0 = ztaus(k) + taucw(k,ib)
2629 zssa0 = zssas(k) + ssacw(k,ib)
2630 zasy0 = zasys(k) + asycw(k,ib)
2631 zssaw = min(oneminus, zssa0 / ztau0)
2632 zasyw = zasy0 / max(ftiny, zssa0)
2638 ztau1 = (f_one - za2) * ztau0
2639 zssa1 = (zssaw - za2) / (f_one - za2)
2641 zasy1 = zasyw / (f_one + zasyw)
2642 zasy3 = 0.75 * zasy1
2646 zgam1 = 1.75 - zssa1 * (f_one + zasy3)
2647 zgam2 =-0.25 + zssa1 * (f_one - zasy3)
2648 zgam3 = 0.5 - zasy3 * cosz
2650 zgam1 = 2.0 - zssa1 * (1.25 + zasy3)
2651 zgam2 = 0.75* zssa1 * (f_one- zasy1)
2652 zgam3 = 0.5 - zasy3 * cosz
2654 zgam1 = zsr3 * (2.0 - zssa1 * (1.0 + zasy1)) * 0.5
2655 zgam2 = zsr3 * zssa1 * (1.0 - zasy1) * 0.5
2656 zgam3 = (1.0 - zsr3 * zasy1 * cosz) * 0.5
2658 zgam4 = f_one - zgam3
2667 if ( zssaw >= zcrit )
then
2668 za1 = zgam1 * cosz - zgam3
2674 zb1 = min( ztau1*sntz , 500.0 )
2675 if ( zb1 <= od_lo )
then
2676 zb2 = f_one - zb1 + 0.5*zb1*zb1
2678 ftind = zb1 / (
bpade + zb1)
2679 itind = ftind*
ntbmx + 0.5
2684 zrefb(kp) = max(f_zero, min(f_one,
2689 zrefd(kp) = max(f_zero, min(f_one, za2 / (f_one+za2) ))
2690 ztrad(kp) = max(f_zero, min(f_one, f_one - zrefd(kp) ))
2693 za1 = zgam1*zgam4 + zgam2*zgam3
2694 za2 = zgam1*zgam3 + zgam2*zgam4
2695 zrk = sqrt( (zgam1 - zgam2) * (zgam1 + zgam2) )
2701 zrpp = f_one - zrp*zrp
2706 zr1 = zrm1 * (za2 + zrkg3)
2707 zr2 = zrp1 * (za2 - zrkg3)
2708 zr3 = zrk2 * (zgam3 - za2*cosz)
2710 zr5 = zrpp * (zrk - zgam1)
2712 zt1 = zrp1 * (za1 + zrkg4)
2713 zt2 = zrm1 * (za1 - zrkg4)
2714 zt3 = zrk2 * (zgam4 + za1*cosz)
2719 zb1 = min( zrk*ztau1, 500.0 )
2720 if ( zb1 <= od_lo )
then
2721 zexm1 = f_one - zb1 + 0.5*zb1*zb1
2723 ftind = zb1 / (
bpade + zb1)
2724 itind = ftind*
ntbmx + 0.5
2727 zexp1 = f_one / zexm1
2729 zb2 = min( ztau1*sntz, 500.0 )
2730 if ( zb2 <= od_lo )
then
2731 zexm2 = f_one - zb2 + 0.5*zb2*zb2
2733 ftind = zb2 / (
bpade + zb2)
2734 itind = ftind*
ntbmx + 0.5
2737 zexp2 = f_one / zexm2
2738 ze1r45 = zr4*zexp1 + zr5*zexm1
2741 if ( ze1r45>=-eps1 .and. ze1r45<=eps1 )
then
2745 zden1 = zssa1 / ze1r45
2746 zrefb(kp) = max(f_zero, min(f_one,
2753 zden1 = zr4 / (ze1r45 * zrkg1)
2754 zrefd(kp) = max(f_zero, min(f_one,
2762 zrefb(kp) = zc0*zrefb1 + zc1*zrefb(kp)
2763 zrefd(kp) = zc0*zrefd1 + zc1*zrefd(kp)
2764 ztrab(kp) = zc0*ztrab1 + zc1*ztrab(kp)
2765 ztrad(kp) = zc0*ztrad1 + zc1*ztrad(kp)
2772 if ( zr1 <= od_lo )
then
2773 zexp3 = f_one - zr1 + 0.5*zr1*zr1
2775 ftind = zr1 / (
bpade + zr1)
2776 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
2780 zldbt(kp) = zc0*zldbt(kp) + zc1*zexp3
2781 ztdbt(k) = zldbt(kp) * ztdbt(kp)
2787 if ( zr1 <= od_lo )
then
2788 zexp4 = f_one - zr1 + 0.5*zr1*zr1
2790 ftind = zr1 / (
bpade + zr1)
2791 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
2795 ztdbt0 = (zc0*zldbt0(k) + zc1*zexp4) * ztdbt0
2800 ztdbt(k) = zldbt(kp) * ztdbt(kp)
2803 ztdbt0 = zldbt0(k) * ztdbt0
2812 & ( zrefb,zrefd,ztrab,ztrad,zldbt,ztdbt,
2820 fxupc(k,ib) = fxupc(k,ib) + zsolar*zfu(k)
2821 fxdnc(k,ib) = fxdnc(k,ib) + zsolar*zfd(k)
2827 zb2 = zsolar*(zfd(1) - ztdbt0)
2830 sfbmc(ibd) = sfbmc(ibd) + zb1
2831 sfdfc(ibd) = sfdfc(ibd) + zb2
2835 sfbmc(1) = sfbmc(1) + zf1
2836 sfdfc(1) = sfdfc(1) + zf2
2837 sfbmc(2) = sfbmc(2) + zf1
2838 sfdfc(2) = sfdfc(2) + zf2
2850 ftoadc = ftoadc + fxdn0(nlp1,ib)
2851 ftoau0 = ftoau0 + fxup0(nlp1,ib)
2852 fsfcu0 = fsfcu0 + fxup0(1,ib)
2853 fsfcd0 = fsfcd0 + fxdn0(1,ib)
2858 suvbf0 = fxdn0(1,ibd)
2860 if ( cf1 <= eps )
then
2863 fxupc(k,ib) = fxup0(k,ib)
2864 fxdnc(k,ib) = fxdn0(k,ib)
2883 fxupc(k,ib) = cf1*fxupc(k,ib) + cf0*fxup0(k,ib)
2884 fxdnc(k,ib) = cf1*fxdnc(k,ib) + cf0*fxdn0(k,ib)
2889 ftoauc = ftoauc + fxupc(nlp1,ib)
2890 fsfcuc = fsfcuc + fxupc(1,ib)
2891 fsfcdc = fsfcdc + fxdnc(1,ib)
2895 suvbfc = fxdnc(1,ibd)
2898 sfbmc(1) = cf1*sfbmc(1) + cf0*sfbm0(1)
2899 sfbmc(2) = cf1*sfbmc(2) + cf0*sfbm0(2)
2900 sfdfc(1) = cf1*sfdfc(1) + cf0*sfdf0(1)
2901 sfdfc(2) = cf1*sfdfc(2) + cf0*sfdf0(2)
2952 & cf1,cf0,taug,taur,tauae,ssaae,asyae,taucw,ssacw,asycw,
2955 & ftoauc,ftoau0,ftoadc,fsfcuc,fsfcu0,fsfcdc,fsfcd0,
3052 real (kind=kind_phys),
parameter :: zcrit = 0.9999995
3053 real (kind=kind_phys),
parameter :: zsr3 = sqrt(3.0)
3054 real (kind=kind_phys),
parameter :: od_lo = 0.06
3055 real (kind=kind_phys),
parameter :: eps1 = 1.0e-8
3058 integer,
intent(in) :: nlay, nlp1
3060 real (kind=kind_phys),
dimension(nlay,ngptsw),
intent(in) :: &
3061 & taug, taur, cldfmc
3062 real (kind=kind_phys),
dimension(nlay,nbdsw),
intent(in) :: &
3063 & taucw, ssacw, asycw, tauae, ssaae, asyae
3065 real (kind=kind_phys),
dimension(ngptsw),
intent(in) :: sfluxzen
3067 real (kind=kind_phys),
dimension(2),
intent(in) :: albbm, albdf
3069 real (kind=kind_phys),
intent(in) :: cosz, sntz, cf1, cf0, ssolar
3072 real (kind=kind_phys),
dimension(nlp1,nbdsw),
intent(out) :: &
3073 & fxupc, fxdnc, fxup0, fxdn0
3075 real (kind=kind_phys),
dimension(2),
intent(out) :: sfbmc, sfdfc, &
3078 real (kind=kind_phys),
intent(out) :: suvbfc, suvbf0, ftoadc, &
3079 & ftoauc, ftoau0, fsfcuc, fsfcu0, fsfcdc, fsfcd0
3082 real (kind=kind_phys),
dimension(nlay) :: ztaus, zssas, zasys, &
3085 real (kind=kind_phys),
dimension(nlp1) :: zrefb, zrefd, ztrab, &
3086 & ztrad, ztdbt, zldbt, zfu, zfd
3088 real (kind=kind_phys) :: ztau1, zssa1, zasy1, ztau0, zssa0, &
3089 & zasy0, zasy3, zssaw, zasyw, zgam1, zgam2, zgam3, zgam4,
3095 integer :: ib, ibd, jb, jg, k, kp, itind
3103 fxdnc(k,ib) = f_zero
3104 fxupc(k,ib) = f_zero
3105 fxdn0(k,ib) = f_zero
3106 fxup0(k,ib) = f_zero
3134 lab_do_jg :
do jg = 1,
ngptsw
3140 zsolar = ssolar * sfluxzen(jg)
3149 zrefb(1) = albbm(ibd)
3150 zrefd(1) = albdf(ibd)
3152 zrefb(1) = 0.5 * (albbm(1) + albbm(2))
3153 zrefd(1) = 0.5 * (albdf(1) + albdf(2))
3171 ztau0 = max( ftiny, taur(k,jg)+taug(k,jg)+tauae(k,ib) )
3172 zssa0 = taur(k,jg) + tauae(k,ib)*ssaae(k,ib)
3173 zasy0 = asyae(k,ib)*ssaae(k,ib)*tauae(k,ib)
3174 zssaw = min( oneminus, zssa0 / ztau0 )
3175 zasyw = zasy0 / max( ftiny, zssa0 )
3186 ztau1 = (f_one - za2) * ztau0
3187 zssa1 = (zssaw - za2) / (f_one - za2)
3189 zasy1 = zasyw / (f_one + zasyw)
3190 zasy3 = 0.75 * zasy1
3194 zgam1 = 1.75 - zssa1 * (f_one + zasy3)
3195 zgam2 =-0.25 + zssa1 * (f_one - zasy3)
3196 zgam3 = 0.5 - zasy3 * cosz
3198 zgam1 = 2.0 - zssa1 * (1.25 + zasy3)
3199 zgam2 = 0.75* zssa1 * (f_one- zasy1)
3200 zgam3 = 0.5 - zasy3 * cosz
3202 zgam1 = zsr3 * (2.0 - zssa1 * (1.0 + zasy1)) * 0.5
3203 zgam2 = zsr3 * zssa1 * (1.0 - zasy1) * 0.5
3204 zgam3 = (1.0 - zsr3 * zasy1 * cosz) * 0.5
3206 zgam4 = f_one - zgam3
3210 if ( zssaw >= zcrit )
then
3211 za1 = zgam1 * cosz - zgam3
3217 zb1 = min( ztau1*sntz , 500.0 )
3218 if ( zb1 <= od_lo )
then
3219 zb2 = f_one - zb1 + 0.5*zb1*zb1
3221 ftind = zb1 / (
bpade + zb1)
3222 itind = ftind*
ntbmx + 0.5
3227 zrefb(kp) = max(f_zero, min(f_one,
3232 zrefd(kp) = max(f_zero, min(f_one, za2/(f_one + za2) ))
3233 ztrad(kp) = max(f_zero, min(f_one, f_one-zrefd(kp) ))
3236 za1 = zgam1*zgam4 + zgam2*zgam3
3237 za2 = zgam1*zgam3 + zgam2*zgam4
3238 zrk = sqrt( (zgam1 - zgam2) * (zgam1 + zgam2) )
3244 zrpp = f_one - zrp*zrp
3249 zr1 = zrm1 * (za2 + zrkg3)
3250 zr2 = zrp1 * (za2 - zrkg3)
3251 zr3 = zrk2 * (zgam3 - za2*cosz)
3253 zr5 = zrpp * (zrk - zgam1)
3255 zt1 = zrp1 * (za1 + zrkg4)
3256 zt2 = zrm1 * (za1 - zrkg4)
3257 zt3 = zrk2 * (zgam4 + za1*cosz)
3262 zb1 = min( zrk*ztau1, 500.0 )
3263 if ( zb1 <= od_lo )
then
3264 zexm1 = f_one - zb1 + 0.5*zb1*zb1
3266 ftind = zb1 / (
bpade + zb1)
3267 itind = ftind*
ntbmx + 0.5
3270 zexp1 = f_one / zexm1
3272 zb2 = min( sntz*ztau1, 500.0 )
3273 if ( zb2 <= od_lo )
then
3274 zexm2 = f_one - zb2 + 0.5*zb2*zb2
3276 ftind = zb2 / (
bpade + zb2)
3277 itind = ftind*
ntbmx + 0.5
3280 zexp2 = f_one / zexm2
3281 ze1r45 = zr4*zexp1 + zr5*zexm1
3284 if (ze1r45>=-eps1 .and. ze1r45<=eps1)
then
3288 zden1 = zssa1 / ze1r45
3289 zrefb(kp) = max(f_zero, min(f_one,
3296 zden1 = zr4 / (ze1r45 * zrkg1)
3297 zrefd(kp) = max(f_zero, min(f_one,
3307 if ( zr1 <= od_lo )
then
3308 zexp3 = f_one - zr1 + 0.5*zr1*zr1
3310 ftind = zr1 / (
bpade + zr1)
3311 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
3315 ztdbt(k) = zexp3 * ztdbt(kp)
3322 if ( zr1 <= od_lo )
then
3323 zexp4 = f_one - zr1 + 0.5*zr1*zr1
3325 ftind = zr1 / (
bpade + zr1)
3326 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
3331 ztdbt0 = zexp4 * ztdbt0
3336 & ( zrefb,zrefd,ztrab,ztrad,zldbt,ztdbt,
3344 fxup0(k,ib) = fxup0(k,ib) + zsolar*zfu(k)
3345 fxdn0(k,ib) = fxdn0(k,ib) + zsolar*zfd(k)
3350 zb2 = zsolar*(zfd(1) - ztdbt0)
3353 sfbm0(ibd) = sfbm0(ibd) + zb1
3354 sfdf0(ibd) = sfdf0(ibd) + zb2
3358 sfbm0(1) = sfbm0(1) + zf1
3359 sfdf0(1) = sfdf0(1) + zf2
3360 sfbm0(2) = sfbm0(2) + zf1
3361 sfdf0(2) = sfdf0(2) + zf2
3376 if ( cf1 > eps )
then
3384 if ( cldfmc(k,jg) > ftiny )
then
3386 ztau0 = ztaus(k) + taucw(k,ib)
3387 zssa0 = zssas(k) + ssacw(k,ib)
3388 zasy0 = zasys(k) + asycw(k,ib)
3389 zssaw = min(oneminus, zssa0 / ztau0)
3390 zasyw = zasy0 / max(ftiny, zssa0)
3396 ztau1 = (f_one - za2) * ztau0
3397 zssa1 = (zssaw - za2) / (f_one - za2)
3399 zasy1 = zasyw / (f_one + zasyw)
3400 zasy3 = 0.75 * zasy1
3404 zgam1 = 1.75 - zssa1 * (f_one + zasy3)
3405 zgam2 =-0.25 + zssa1 * (f_one - zasy3)
3406 zgam3 = 0.5 - zasy3 * cosz
3408 zgam1 = 2.0 - zssa1 * (1.25 + zasy3)
3409 zgam2 = 0.75* zssa1 * (f_one- zasy1)
3410 zgam3 = 0.5 - zasy3 * cosz
3412 zgam1 = zsr3 * (2.0 - zssa1 * (1.0 + zasy1)) * 0.5
3413 zgam2 = zsr3 * zssa1 * (1.0 - zasy1) * 0.5
3414 zgam3 = (1.0 - zsr3 * zasy1 * cosz) * 0.5
3416 zgam4 = f_one - zgam3
3420 if ( zssaw >= zcrit )
then
3421 za1 = zgam1 * cosz - zgam3
3427 zb1 = min( ztau1*sntz , 500.0 )
3428 if ( zb1 <= od_lo )
then
3429 zb2 = f_one - zb1 + 0.5*zb1*zb1
3431 ftind = zb1 / (
bpade + zb1)
3432 itind = ftind*
ntbmx + 0.5
3437 zrefb(kp) = max(f_zero, min(f_one,
3442 zrefd(kp) = max(f_zero, min(f_one, za2 / (f_one+za2) ))
3443 ztrad(kp) = max(f_zero, min(f_one, f_one - zrefd(kp) ))
3446 za1 = zgam1*zgam4 + zgam2*zgam3
3447 za2 = zgam1*zgam3 + zgam2*zgam4
3448 zrk = sqrt( (zgam1 - zgam2) * (zgam1 + zgam2) )
3454 zrpp = f_one - zrp*zrp
3459 zr1 = zrm1 * (za2 + zrkg3)
3460 zr2 = zrp1 * (za2 - zrkg3)
3461 zr3 = zrk2 * (zgam3 - za2*cosz)
3463 zr5 = zrpp * (zrk - zgam1)
3465 zt1 = zrp1 * (za1 + zrkg4)
3466 zt2 = zrm1 * (za1 - zrkg4)
3467 zt3 = zrk2 * (zgam4 + za1*cosz)
3472 zb1 = min( zrk*ztau1, 500.0 )
3473 if ( zb1 <= od_lo )
then
3474 zexm1 = f_one - zb1 + 0.5*zb1*zb1
3476 ftind = zb1 / (
bpade + zb1)
3477 itind = ftind*
ntbmx + 0.5
3480 zexp1 = f_one / zexm1
3482 zb2 = min( ztau1*sntz, 500.0 )
3483 if ( zb2 <= od_lo )
then
3484 zexm2 = f_one - zb2 + 0.5*zb2*zb2
3486 ftind = zb2 / (
bpade + zb2)
3487 itind = ftind*
ntbmx + 0.5
3490 zexp2 = f_one / zexm2
3491 ze1r45 = zr4*zexp1 + zr5*zexm1
3494 if ( ze1r45>=-eps1 .and. ze1r45<=eps1 )
then
3498 zden1 = zssa1 / ze1r45
3499 zrefb(kp) = max(f_zero, min(f_one,
3506 zden1 = zr4 / (ze1r45 * zrkg1)
3507 zrefd(kp) = max(f_zero, min(f_one,
3517 if ( zr1 <= od_lo )
then
3518 zexp3 = f_one - zr1 + 0.5*zr1*zr1
3520 ftind = zr1 / (
bpade + zr1)
3521 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
3526 ztdbt(k) = zexp3 * ztdbt(kp)
3532 if ( zr1 <= od_lo )
then
3533 zexp4 = f_one - zr1 + 0.5*zr1*zr1
3535 ftind = zr1 / (
bpade + zr1)
3536 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
3540 ztdbt0 = zexp4 * ztdbt0
3545 ztdbt(k) = zldbt(kp) * ztdbt(kp)
3548 ztdbt0 = zldbt0(k) * ztdbt0
3557 & ( zrefb,zrefd,ztrab,ztrad,zldbt,ztdbt,
3565 fxupc(k,ib) = fxupc(k,ib) + zsolar*zfu(k)
3566 fxdnc(k,ib) = fxdnc(k,ib) + zsolar*zfd(k)
3572 zb2 = zsolar*(zfd(1) - ztdbt0)
3575 sfbmc(ibd) = sfbmc(ibd) + zb1
3576 sfdfc(ibd) = sfdfc(ibd) + zb2
3580 sfbmc(1) = sfbmc(1) + zf1
3581 sfdfc(1) = sfdfc(1) + zf2
3582 sfbmc(2) = sfbmc(2) + zf1
3583 sfdfc(2) = sfdfc(2) + zf2
3595 ftoadc = ftoadc + fxdn0(nlp1,ib)
3596 ftoau0 = ftoau0 + fxup0(nlp1,ib)
3597 fsfcu0 = fsfcu0 + fxup0(1,ib)
3598 fsfcd0 = fsfcd0 + fxdn0(1,ib)
3603 suvbf0 = fxdn0(1,ibd)
3605 if ( cf1 <= eps )
then
3608 fxupc(k,ib) = fxup0(k,ib)
3609 fxdnc(k,ib) = fxdn0(k,ib)
3627 ftoauc = ftoauc + fxupc(nlp1,ib)
3628 fsfcuc = fsfcuc + fxupc(1,ib)
3629 fsfcdc = fsfcdc + fxdnc(1,ib)
3633 suvbfc = fxdnc(1,ibd)
3688 integer,
intent(in) :: nlay, nlp1
3690 real (kind=kind_phys),
dimension(nlp1),
intent(in) :: zrefb,
3691 & zrefd, ztrab, ztrad, ztdbt, zldbt
3694 real (kind=kind_phys),
dimension(nlp1),
intent(out) :: zfu, zfd
3697 real (kind=kind_phys),
dimension(nlp1) :: zrupb,zrupd,zrdnd,ztdn
3699 real (kind=kind_phys) :: zden1
3714 zden1 = f_one / ( f_one - zrupd(k)*zrefd(kp) )
3715 zrupb(kp) = zrefb(kp) + ( ztrad(kp) *
3723 zrdnd(nlp1) = f_zero
3724 ztdn(nlay) = ztrab(nlp1)
3725 zrdnd(nlay) = zrefd(nlp1)
3729 zden1 = f_one / (f_one - zrefd(k)*zrdnd(k))
3730 ztdn(k-1) = ztdbt(k)*ztrab(k) + ( ztrad(k) *
3738 zden1 = f_one / (f_one - zrdnd(k)*zrupd(k))
3739 zfu(k) = ( ztdbt(k)*zrupb(k) +
3793 & forfac,forfrac,indfor,selffac,selffrac,indself, nlay,
3894 integer,
intent(in) :: nlay, laytrop
3896 integer,
dimension(nlay),
intent(in) :: indfor, indself, &
3899 real (kind=kind_phys),
dimension(nlay),
intent(in) :: colmol, &
3900 & fac00, fac01, fac10, fac11, forfac, forfrac, selffac,
3903 real (kind=kind_phys),
dimension(nlay,maxgas),
intent(in) :: colamt
3906 real (kind=kind_phys),
dimension(ngptsw),
intent(out) :: sfluxzen
3908 real (kind=kind_phys),
dimension(nlay,ngptsw),
intent(out) :: &
3912 real (kind=kind_phys) :: fs, speccomb, specmult, colm1, colm2
3914 integer,
dimension(nlay,nblow:nbhgh) :: id0, id1
3916 integer :: ibd, j, jb, js, k, klow, khgh, klim, ks, njb, ns
3927 id0(k,jb) = ((jp(k)-1)*5 + (jt(k)-1)) * nspa(jb)
3928 id1(k,jb) = ( jp(k) *5 + (jt1(k)-1)) * nspa(jb)
3931 do k = laytrop+1, nlay
3932 id0(k,jb) = ((jp(k)-13)*5 + (jt(k)-1)) * nspb(jb)
3933 id1(k,jb) = ((jp(k)-12)*5 + (jt1(k)-1)) * nspb(jb)
3944 case (16, 20, 23, 25, 26, 29)
3947 sfluxzen(ns+j) = sfluxref01(j,1,ibd)
3953 sfluxzen(ns+j) = scalekur * sfluxref01(j,1,ibd)
3958 if (jb==17 .or. jb==28)
then
3961 lab_do_k1 :
do k = laytrop, nlay-1
3962 if (jp(k)<layreffr(jb) .and. jp(k+1)>=layreffr(jb))
then
3968 colm1 = colamt(ks,ix1(jb))
3969 colm2 = colamt(ks,ix2(jb))
3970 speccomb = colm1 + strrat(jb)*colm2
3971 specmult = specwt(jb) * min( oneminus, colm1/speccomb )
3972 js = 1 + int( specmult )
3973 fs = mod(specmult, f_one)
3976 sfluxzen(ns+j) = sfluxref02(j,js,ibd)
3983 lab_do_k2 :
do k = 1, laytrop-1
3984 if (jp(k)<layreffr(jb) .and. jp(k+1)>=layreffr(jb))
then
3990 colm1 = colamt(ks,ix1(jb))
3991 colm2 = colamt(ks,ix2(jb))
3992 speccomb = colm1 + strrat(jb)*colm2
3993 specmult = specwt(jb) * min( oneminus, colm1/speccomb )
3994 js = 1 + int( specmult )
3995 fs = mod(specmult, f_one)
3998 sfluxzen(ns+j) = sfluxref03(j,js,ibd)
4058 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
4059 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
4061 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4062 integer :: inds, indf, indsp, indfp, j, js, k
4073 tauray = colmol(k) *
rayl
4076 taur(k,ns16+j) = tauray
4081 speccomb = colamt(k,1) + strrat(16)*colamt(k,5)
4082 specmult = 8.0 * min( oneminus, colamt(k,1)/speccomb )
4084 js = 1 + int( specmult )
4085 fs = mod( specmult, f_one )
4087 fac000 = fs1 * fac00(k)
4088 fac010 = fs1 * fac10(k)
4089 fac100 = fs * fac00(k)
4090 fac110 = fs * fac10(k)
4091 fac001 = fs1 * fac01(k)
4092 fac011 = fs1 * fac11(k)
4093 fac101 = fs * fac01(k)
4094 fac111 = fs * fac11(k)
4096 ind01 = id0(k,16) + js
4100 ind11 = id1(k,16) + js
4110 taug(k,ns16+j) = speccomb
4122 do k = laytrop+1, nlay
4123 ind01 = id0(k,16) + 1
4125 ind11 = id1(k,16) + 1
4129 taug(k,ns16+j) = colamt(k,5)
4154 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
4155 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
4157 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4158 integer :: inds, indf, indsp, indfp, j, js, k
4169 tauray = colmol(k) *
rayl
4172 taur(k,ns17+j) = tauray
4177 speccomb = colamt(k,1) + strrat(17)*colamt(k,2)
4178 specmult = 8.0 * min(oneminus, colamt(k,1) / speccomb)
4180 js = 1 + int(specmult)
4181 fs = mod(specmult, f_one)
4183 fac000 = fs1 * fac00(k)
4184 fac010 = fs1 * fac10(k)
4185 fac100 = fs * fac00(k)
4186 fac110 = fs * fac10(k)
4187 fac001 = fs1 * fac01(k)
4188 fac011 = fs1 * fac11(k)
4189 fac101 = fs * fac01(k)
4190 fac111 = fs * fac11(k)
4192 ind01 = id0(k,17) + js
4196 ind11 = id1(k,17) + js
4207 taug(k,ns17+j) = speccomb
4219 do k = laytrop+1, nlay
4220 speccomb = colamt(k,1) + strrat(17)*colamt(k,2)
4221 specmult = 4.0 * min(oneminus, colamt(k,1) / speccomb)
4223 js = 1 + int(specmult)
4224 fs = mod(specmult, f_one)
4226 fac000 = fs1 * fac00(k)
4227 fac010 = fs1 * fac10(k)
4228 fac100 = fs * fac00(k)
4229 fac110 = fs * fac10(k)
4230 fac001 = fs1 * fac01(k)
4231 fac011 = fs1 * fac11(k)
4232 fac101 = fs * fac01(k)
4233 fac111 = fs * fac11(k)
4235 ind01 = id0(k,17) + js
4239 ind11 = id1(k,17) + js
4248 taug(k,ns17+j) = speccomb
4277 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
4278 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
4280 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4281 integer :: inds, indf, indsp, indfp, j, js, k
4292 tauray = colmol(k) *
rayl
4295 taur(k,ns18+j) = tauray
4300 speccomb = colamt(k,1) + strrat(18)*colamt(k,5)
4301 specmult = 8.0 * min(oneminus, colamt(k,1) / speccomb)
4303 js = 1 + int(specmult)
4304 fs = mod(specmult, f_one)
4306 fac000 = fs1 * fac00(k)
4307 fac010 = fs1 * fac10(k)
4308 fac100 = fs * fac00(k)
4309 fac110 = fs * fac10(k)
4310 fac001 = fs1 * fac01(k)
4311 fac011 = fs1 * fac11(k)
4312 fac101 = fs * fac01(k)
4313 fac111 = fs * fac11(k)
4315 ind01 = id0(k,18) + js
4319 ind11 = id1(k,18) + js
4330 taug(k,ns18+j) = speccomb
4342 do k = laytrop+1, nlay
4343 ind01 = id0(k,18) + 1
4345 ind11 = id1(k,18) + 1
4349 taug(k,ns18+j) = colamt(k,5)
4373 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
4374 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
4376 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4377 integer :: inds, indf, indsp, indfp, j, js, k
4388 tauray = colmol(k) *
rayl
4391 taur(k,ns19+j) = tauray
4396 speccomb = colamt(k,1) + strrat(19)*colamt(k,2)
4397 specmult = 8.0 * min(oneminus, colamt(k,1) / speccomb)
4399 js = 1 + int(specmult)
4400 fs = mod(specmult, f_one)
4402 fac000 = fs1 * fac00(k)
4403 fac010 = fs1 * fac10(k)
4404 fac100 = fs * fac00(k)
4405 fac110 = fs * fac10(k)
4406 fac001 = fs1 * fac01(k)
4407 fac011 = fs1 * fac11(k)
4408 fac101 = fs * fac01(k)
4409 fac111 = fs * fac11(k)
4411 ind01 = id0(k,19) + js
4415 ind11 = id1(k,19) + js
4426 taug(k,ns19+j) = speccomb
4438 do k = laytrop+1, nlay
4439 ind01 = id0(k,19) + 1
4441 ind11 = id1(k,19) + 1
4445 taug(k,ns19+j) = colamt(k,2)
4469 real (kind=kind_phys) :: tauray
4471 integer :: ind01, ind02, ind11, ind12
4472 integer :: inds, indf, indsp, indfp, j, k
4483 tauray = colmol(k) *
rayl
4486 taur(k,ns20+j) = tauray
4491 ind01 = id0(k,20) + 1
4493 ind11 = id1(k,20) + 1
4502 taug(k,ns20+j) = colamt(k,1)
4513 do k = laytrop+1, nlay
4514 ind01 = id0(k,20) + 1
4516 ind11 = id1(k,20) + 1
4523 taug(k,ns20+j) = colamt(k,1)
4551 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
4552 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
4554 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4555 integer :: inds, indf, indsp, indfp, j, js, k
4566 tauray = colmol(k) *
rayl
4569 taur(k,ns21+j) = tauray
4574 speccomb = colamt(k,1) + strrat(21)*colamt(k,2)
4575 specmult = 8.0 * min(oneminus, colamt(k,1) / speccomb)
4577 js = 1 + int(specmult)
4578 fs = mod(specmult, f_one)
4580 fac000 = fs1 * fac00(k)
4581 fac010 = fs1 * fac10(k)
4582 fac100 = fs * fac00(k)
4583 fac110 = fs * fac10(k)
4584 fac001 = fs1 * fac01(k)
4585 fac011 = fs1 * fac11(k)
4586 fac101 = fs * fac01(k)
4587 fac111 = fs * fac11(k)
4589 ind01 = id0(k,21) + js
4593 ind11 = id1(k,21) + js
4604 taug(k,ns21+j) = speccomb
4616 do k = laytrop+1, nlay
4617 speccomb = colamt(k,1) + strrat(21)*colamt(k,2)
4618 specmult = 4.0 * min(oneminus, colamt(k,1) / speccomb)
4620 js = 1 + int(specmult)
4621 fs = mod(specmult, f_one)
4623 fac000 = fs1 * fac00(k)
4624 fac010 = fs1 * fac10(k)
4625 fac100 = fs * fac00(k)
4626 fac110 = fs * fac10(k)
4627 fac001 = fs1 * fac01(k)
4628 fac011 = fs1 * fac11(k)
4629 fac101 = fs * fac01(k)
4630 fac111 = fs * fac11(k)
4632 ind01 = id0(k,21) + js
4636 ind11 = id1(k,21) + js
4645 taug(k,ns21+j) = speccomb
4673 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
4674 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111,
4677 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4678 integer :: inds, indf, indsp, indfp, j, js, k
4688 o2tem = 4.35e-4 / (350.0*2.0)
4696 tauray = colmol(k) *
rayl
4699 taur(k,ns22+j) = tauray
4704 o2cont = o2tem * colamt(k,6)
4705 speccomb = colamt(k,1) + strrat(22)*colamt(k,6)
4706 specmult = 8.0 * min(oneminus, colamt(k,1) / speccomb)
4708 js = 1 + int(specmult)
4709 fs = mod(specmult, f_one)
4711 fac000 = fs1 * fac00(k)
4712 fac010 = fs1 * fac10(k)
4713 fac100 = fs * fac00(k)
4714 fac110 = fs * fac10(k)
4715 fac001 = fs1 * fac01(k)
4716 fac011 = fs1 * fac11(k)
4717 fac101 = fs * fac01(k)
4718 fac111 = fs * fac11(k)
4720 ind01 = id0(k,22) + js
4724 ind11 = id1(k,22) + js
4735 taug(k,ns22+j) = speccomb
4747 do k = laytrop+1, nlay
4748 o2cont = o2tem * colamt(k,6)
4750 ind01 = id0(k,22) + 1
4752 ind11 = id1(k,22) + 1
4756 taug(k,ns22+j) = colamt(k,6) * o2adj
4782 integer :: ind01, ind02, ind11, ind12
4783 integer :: inds, indf, indsp, indfp, j, k
4795 taur(k,ns23+j) = colmol(k) *
rayl(j)
4800 ind01 = id0(k,23) + 1
4802 ind11 = id1(k,23) + 1
4811 taug(k,ns23+j) = colamt(k,1) * (
givfac
4821 do k = laytrop+1, nlay
4823 taug(k,ns23+j) = f_zero
4845 real (kind=kind_phys) :: speccomb, specmult, fs, fs1, &
4846 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
4848 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4849 integer :: inds, indf, indsp, indfp, j, js, k
4860 speccomb = colamt(k,1) + strrat(24)*colamt(k,6)
4861 specmult = 8.0 * min(oneminus, colamt(k,1) / speccomb)
4863 js = 1 + int(specmult)
4864 fs = mod(specmult, f_one)
4866 fac000 = fs1 * fac00(k)
4867 fac010 = fs1 * fac10(k)
4868 fac100 = fs * fac00(k)
4869 fac110 = fs * fac10(k)
4870 fac001 = fs1 * fac01(k)
4871 fac011 = fs1 * fac11(k)
4872 fac101 = fs * fac01(k)
4873 fac111 = fs * fac11(k)
4875 ind01 = id0(k,24) + js
4879 ind11 = id1(k,24) + js
4890 taug(k,ns24+j) = speccomb
4906 do k = laytrop+1, nlay
4907 ind01 = id0(k,24) + 1
4909 ind11 = id1(k,24) + 1
4913 taug(k,ns24+j) = colamt(k,6)
4941 integer :: ind01, ind02, ind11, ind12
4954 taur(k,ns25+j) = colmol(k) *
rayl(j)
4959 ind01 = id0(k,25) + 1
4961 ind11 = id1(k,25) + 1
4965 taug(k,ns25+j) = colamt(k,1)
4972 do k = laytrop+1, nlay
4974 taug(k,ns25+j) = colamt(k,3) *
abso3b(j)
5009 taug(k,ns26+j) = f_zero
5010 taur(k,ns26+j) = colmol(k) *
rayl(j)
5033 integer :: ind01, ind02, ind11, ind12
5046 taur(k,ns27+j) = colmol(k) *
rayl(j)
5051 ind01 = id0(k,27) + 1
5053 ind11 = id1(k,27) + 1
5057 taug(k,ns27+j) = colamt(k,3)
5063 do k = laytrop+1, nlay
5064 ind01 = id0(k,27) + 1
5066 ind11 = id1(k,27) + 1
5070 taug(k,ns27+j) = colamt(k,3)
5095 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
5096 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
5098 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
5110 tauray = colmol(k) *
rayl
5113 taur(k,ns28+j) = tauray
5118 speccomb = colamt(k,3) + strrat(28)*colamt(k,6)
5119 specmult = 8.0 * min(oneminus, colamt(k,3) / speccomb)
5121 js = 1 + int(specmult)
5122 fs = mod(specmult, f_one)
5124 fac000 = fs1 * fac00(k)
5125 fac010 = fs1 * fac10(k)
5126 fac100 = fs * fac00(k)
5127 fac110 = fs * fac10(k)
5128 fac001 = fs1 * fac01(k)
5129 fac011 = fs1 * fac11(k)
5130 fac101 = fs * fac01(k)
5131 fac111 = fs * fac11(k)
5133 ind01 = id0(k,28) + js
5137 ind11 = id1(k,28) + js
5143 taug(k,ns28+j) = speccomb
5151 do k = laytrop+1, nlay
5152 speccomb = colamt(k,3) + strrat(28)*colamt(k,6)
5153 specmult = 4.0 * min(oneminus, colamt(k,3) / speccomb)
5155 js = 1 + int(specmult)
5156 fs = mod(specmult, f_one)
5158 fac000 = fs1 * fac00(k)
5159 fac010 = fs1 * fac10(k)
5160 fac100 = fs * fac00(k)
5161 fac110 = fs * fac10(k)
5162 fac001 = fs1 * fac01(k)
5163 fac011 = fs1 * fac11(k)
5164 fac101 = fs * fac01(k)
5165 fac111 = fs * fac11(k)
5167 ind01 = id0(k,28) + js
5171 ind11 = id1(k,28) + js
5177 taug(k,ns28+j) = speccomb
5204 real (kind=kind_phys) :: tauray
5206 integer :: ind01, ind02, ind11, ind12
5207 integer :: inds, indf, indsp, indfp, j, k
5218 tauray = colmol(k) *
rayl
5221 taur(k,ns29+j) = tauray
5226 ind01 = id0(k,29) + 1
5228 ind11 = id1(k,29) + 1
5237 taug(k,ns29+j) = colamt(k,1)
5248 do k = laytrop+1, nlay
5249 ind01 = id0(k,29) + 1
5251 ind11 = id1(k,29) + 1
5255 taug(k,ns29+j) = colamt(k,2)
5274 end module module_radsw_main
real(kind=kind_phys), dimension(msf22, ng22), public selfref
the array selfref contains the coefficient of the water vapor self-continuum (including the energy te...
integer, save iovrsw
cloud overlapping control flag for SW
real(kind=kind_phys), parameter, public rayl
rayleigh extinction coefficient at
This module sets up absorption coeffients for band 19: 4650-5150 cm-1 (low - h2o,co2; high - co2) ...
integer, save isubcsw
sub-column cloud approx flag in SW radiation
subroutine spcvrtm
This subroutine computes the shortwave radiative fluxes using two-stream method of h...
real(kind=kind_phys), parameter con_amw
molecular wght of water vapor ( )
subroutine spcvrtc
This subroutine computes the shortwave radiative fluxes using two-stream method.
real(kind=kind_phys), dimension(msf24, ng24), public selfref
the array selfref contains the coefficient of the water vapor self-continuum (including the energy te...
real(kind=kind_phys), parameter, public rayl
rayleigh extinction coefficient at
real(kind=kind_phys), dimension(msf18, ng18), public selfref
the array selfref contains the coefficient of the water vapor self-continuum (including the energy te...
This module sets up absorption coeffients for band 24: 12850-16000 cm-1 (low - h2o, o2; high - o2)
real(kind=kind_phys), dimension(msa19, ng19), public absa
the array absa(585,NG19) (ka(9,5,13,NG19)) contains absorption coefs at the 16 chosen g-values for a ...
subroutine taumol17
The subroutine computes the optical depth in band 17: 3250-4000 cm-1 (low - h2o,co2; high - h2o...
real(kind=kind_phys) heatfac
the factor for heating rates (in k/day, or k/sec set by subroutine 'rswinit')
real(kind=kind_phys), dimension(msa23, ng23), public absa
the array absa(65,NG23) (ka(5,13,NG23)) contains absorption coefs at the 16 chosen g-values for a ran...
real(kind=kind_phys), parameter con_g
gravity ( )
real(kind=kind_phys), dimension(0:ntbmx) exp_tbl
those data will be set up only once by "rswinit"
Define type construct for optional component downward fluxes at surface.
real(kind=kind_phys), dimension(ng24), public abso3a
o3
real(kind=kind_phys), dimension(msb28, ng28), public absb
the array absb(1175,6) (kb(5,5,13:59,6)) contains absorption coefs at the 16 chosen g-values for a ra...
This module sets up absorption coeffients for band 27: 29000-38000 cm-1 (low - o3; high - o3) ...
integer, parameter iswrgas
SW rare gases effect control flag (ch4,n2o,o2,...): =0:no; =1:yes.
real(kind=kind_phys), dimension(msf16, ng16), public selfref
the array selfref contains the coefficient of the water vapor self-continuum (including the energy te...
integer, parameter nbhgh
band range upper limit
This module sets up absorption coeffients for band 20: 5150-6150 cm-1 (low - h2o; high - h2o) ...
integer, parameter ntbmx
index upper limit of trans table
real(kind=kind_phys), dimension(msb24, ng24), public absb
the array absb(235,8) (kb(5,13:59,8)) contains absorption coefs at the 16 chosen g-values for a range...
integer, parameter nuvb
uv-b band index
real(kind=kind_phys), dimension(ng26), public rayl
rayleigh extinction coefficient at all v
integer, parameter ngptsw
total number of g-point in all bands
real(kind=kind_phys), dimension(msf23, ng23), public selfref
the array selfref contains the coefficient of the water vapor self-continuum (including the energy te...
real(kind=kind_phys), parameter, public rayl
rayleigh extinction coefficient at
subroutine, public swrad
This subroutine is the main SW radiation routine.
integer, parameter iswmode
SW control flag for 2-stream transfer scheme =1:delta-eddington (Joseph et al. 1976 ) =2:pifm (Zd...
This module contains coefficients of cloud optical properties for each of the spectral bands...
subroutine mcica_subcol
This subroutine computes the sub-colum cloud profile flag array.
real(kind=kind_phys), dimension(msb22, ng22), public absb
the array absb(235,2) (kb(5,13:59,2)) contains absorption coefs at the 16 chosen g-values for a range...
real(kind=kind_phys), dimension(ng24, mfx24), public rayla
rayleigh extinction coefficient at all v
This module sets up absorption coeffients for band 23: 8050-12850 cm-1 (low - h2o; high - nothing) ...
real(kind=kind_phys), dimension(msf19, ng19), public selfref
the array selfref contains the coefficient of the water vapor self-continuum (including the energy te...
subroutine cldprop
This subroutine computes the cloud optical properties for each cloudy layer and g-point interval...
real(kind=kind_phys), dimension(msa29, ng29), public absa
the array absa(65,NG29) (ka(5,13,NG29)) contains absorption coefs at the 16 chosen g-values for a ran...
real(kind=kind_phys), dimension(msa22, ng22), public absa
the array absa(585,NG22) (ka(9,5,13,NG22)) contains absorption coefs at the 16 chosen g-values for a ...
subroutine vrtqdr
This subroutine is called by spcvrtc() and spcvrtm(), and computes the upward and downward radiation ...
real(kind=kind_phys), dimension(msa27, ng27), public absa
the array absa(65,NG27) (ka(5,13,NG27)) contains absorption coefs at the 16 chosen g-values for a ran...
real(kind=kind_phys), dimension(msb27, ng27), public absb
the array absb(235,8) (kb(5,13:59,8)) contains absorption coefs at the 16 chosen g-values for a range...
real(kind=kind_phys), dimension(msf20, ng20), public selfref
the array selfref contains the coefficient of the water vapor self-continuum (including the energy te...
real(kind=kind_phys), dimension(msa16, ng16), public absa
the array absa(585,NG16) (ka(9,5,13,NG16)) contains absorption coefs at the 16 chosen g-values for a ...
real(kind=kind_phys), dimension(msa24, ng24), public absa
the array absa(585,NG24) (ka(9,5,13,NG24)) contains absorption coefs at the 16 chosen g-values for a ...
This module sets up absorption coeffients for band 18: 4000-4650 cm-1 (low - h2o,ch4; high - ch4) ...
integer, save icldflg
cloud optical property scheme control flag
subroutine taumol22
The subroutine computes the optical depth in band 22: 7700-8050 cm-1 (low - h2o,o2; high - o2) ...
Define type construct for optional radiation flux profiles.
integer, dimension(nblow:nbhgh) idxebc
band index for cld prop
real(kind=kind_phys), dimension(msb29, ng29), public absb
the array absb(235,12) (kb(5,13:59,12)) contains absorption coefs at the 16 chosen g-values for a ran...
real(kind=kind_phys), dimension(ng20), public absch4
ch4
real(kind=kind_phys), parameter con_cp
spec heat air at p ( )
subroutine taumol23
The subroutine computes the optical depth in band 23: 8050-12850 cm-1 (low - h2o; high - nothing) ...
real(kind=kind_phys), dimension(59) preflog
Reference pressure and temperature.
subroutine taumol26
The subroutine computes the optical depth in band 26: 22650-29000 cm-1 (low - nothing; high - nothing...
integer, save iswcice
SW optical property for ice clouds (only iswcliq>0) =0:not defined yet =1:input cip...
real(kind=kind_phys), dimension(ng25), public abso3b
o3
This module contains reference temperature and pressure.
integer, dimension(nblow:nbhgh) idxsfc
band index for sfc flux
subroutine taumol16
The subroutine computes the optical depth in band 16: 2600-3250 cm-1 (low - h2o,ch4; high - ch4) ...
This module contains SW band parameters set up.
real(kind=kind_phys), dimension(msb21, ng21), public absb
the array absb(1175,10) (kb(5,5,13:59,10)) contains absorption coefs at the 16 chosen g-values for a ...
integer, parameter nblow
band range lower limit
subroutine taumol19
The subroutine computes the optical depth in band 19: 4650-5150 cm-1 (low - h2o,co2; high - co2) ...
This module sets up absorption coeffients for band 26: 22650-29000 cm-1 (low - nothing; high - nothin...
subroutine taumol28
The subroutine computes the optical depth in band 28: 38000-50000 cm-1 (low - o3,o2; high - o3...
subroutine taumol18
The subroutine computes the optical depth in band 18: 4000-4650 cm-1 (low - h2o,ch4; high - ch4) ...
integer, parameter iswrate
SW heating rate unit control flag: =1:k/day; =2:k/second.
subroutine setcoef
This subroutine computes various coefficients needed in radiative transfer calculation.
This module sets up absorption coeffients for band 29: 820-2600 cm-1 (low - h2o; high - co2) ...
integer, save iswcliq
SW optical property for liquid clouds =0:input cld opt depth, ignoring iswcice setting =1:input c...
This module sets up absorption coeffients for band 17: 3250-4000 cm-1 (low - h2o,co2; high - h2o...
subroutine, public rswinit
This subroutine initializes non-varying module variables, conversion factors, and look-up tables...
subroutine taumol29
The subroutine computes the optical depth in band 29: 820-2600 cm-1 (low - h2o; high - co2) ...
real(kind=kind_phys), dimension(msa25, ng25), public absa
the array absa(65,NG25) (ka(5,13,NG25)) contains absorption coefs at the 16 chosen g-values for a ran...
subroutine taumol27
The subroutine computes the optical depth in band 27: 29000-38000 cm-1 (low - o3; high - o3) ...
This module sets up absorption coeffients for band 22: 7700-8050 cm-1 (low - h2o, o2; high - o2) ...
real(kind=kind_phys), dimension(msf29, ng29), public selfref
the array selfref contains the coefficient of the water vapor self-continuum (including the energy te...
real(kind=kind_phys), parameter s0
internal solar constant
integer, parameter ipsdsw0
initial permutation seed used for sub-column cloud scheme
This module contains spectral distribution of solar radiation flux used to obtain the incoming solar ...
real(kind=kind_phys), dimension(msa17, ng17), public absa
the array absa(585,NG17) (ka((9,5,13,NG17)) contains absorption coefs at the 16 chosen g-values for a...
real(kind=kind_phys), dimension(ng23), public rayl
rayleigh extinction coefficient at all v
real(kind=kind_phys), parameter con_avgd
avogadro constant ( )
real(kind=kind_phys), dimension(ng25), public abso3a
o3
subroutine taumol20
The subroutine computes the optical depth in band 20: 5150-6150 cm-1 (low - h2o; high - h2o) ...
subroutine taumol21
The subroutine computes the optical depth in band 21: 6150-7700 cm-1 (low - h2o,co2; high - h2o...
real(kind=kind_phys), dimension(msa21, ng21), public absa
the array absa(585,NG21) (ka((9,5,13,NG21)) contains absorption coefs at the 16 chosen g-values for a...
This module sets up absorption coeffients for band 25: 16000-22650 cm-1 (low - h2o; high - nothing) ...
real(kind=kind_phys), dimension(msb20, ng20), public absb
the array absb(235,10) (kb(5,13:59,10)) contains absorption coefs at the 16 chosen g-values for a ran...
real(kind=kind_phys), dimension(msa18, ng18), public absa
the array absa(585,NG18) (ka(9,5,13,NG18)) contains absorption coefs at the 16 chosen g-values for a ...
real(kind=kind_phys), dimension(msb16, ng16), public absb
the array absb(235,6) (kb(5,13:59,6)) contains absorption coefs at the 16 chosen g-values for a range...
real(kind=kind_phys), parameter bpade
pade approx constant
real(kind=kind_phys), parameter, public rayl
rayleigh extinction coefficient at
real(kind=kind_phys), dimension(msb17, ng17), public absb
the array absb(1175,12) (kb(5,5,13:59,12)) contains absorption coefs at the 16 chosen g-values for a ...
real(kind=kind_phys), parameter, public rayl
rayleigh extinction coefficient
real(kind=kind_phys), parameter con_amo3
molecular wght of o3 ( )
subroutine taumol25
The subroutine computes the optical depth in band 25: 16000-22650 cm-1 (low - h2o; high - nothing) ...
real(kind=kind_phys), dimension(msa20, ng20), public absa
the array absa(65,NG20) (ka(5,13,NG20)) contains absorption coefs at the 16 chosen g-values for a ran...
This module sets up absorption coeffients for band 21: 6150-7700 cm-1 (low - h2o,co2; high - h2o...
real(kind=kind_phys), dimension(msb18, ng18), public absb
the array absb(235,8) (kb(5,13:59,8)) contains absorption coefs at the 16 chosen g-values for a range...
real(kind=kind_phys), parameter, public rayl
rayleigh extinction coefficient at v =
real(kind=kind_phys), parameter, public givfac
average giver et al. correction factor for this band.
Define type construct for radiation fluxes at toa.
real(kind=kind_phys), dimension(msf21, ng21), public selfref
the array selfref contains the coefficient of the water vapor self-continuum (including the energy te...
real(kind=kind_phys), dimension(ng27), public rayl
rayleigh extinction coefficient
real(kind=kind_phys), dimension(msf17, ng17), public selfref
the array selfref contains the coefficient of the water vapor self-continuum (including the energy te...
real(kind=kind_phys), parameter, public rayl
rayleigh extinction coefficient at
integer, dimension(ngptsw) ngb
band index for each g-point
integer, save ivflip
vertical profile indexing flag
real(kind=kind_phys), dimension(ng29), public absh2o
h2o
integer, parameter maxgas
max num of absorbing gases
real(kind=kind_phys), dimension(msb19, ng19), public absb
the array absb(235,8) (kb(5,13:59,8)) contains absorption coefs at the 16 chosen g-values for a range...
This module sets up absorption coeffients for band 28: 38000-50000 cm-1 (low - o3,o2; high - o3,o2)
real(kind=kind_phys), dimension(msa28, ng28), public absa
the array absa(585,NG28) (ka((9,5,13,NG28)) contains absorption coefs at the 16 chosen g-values for a...
This module sets up absorption coefficients for band 16: 2600-3250 cm-1 (low - h2o, ch4; high - ch4)
Define type construct for radiation fluxes at surface.
subroutine taumol
This subroutine calculates optical depths for gaseous absorption and rayleigh scattering subroutine...
real(kind=kind_phys), dimension(ng24), public abso3b
o3
real(kind=kind_phys), parameter, public rayl
rayleigh extinction coefficient at v = 3625
subroutine taumol24
The subroutine computes the optical depth in band 24: 12850-16000 cm-1 (low - h2o,o2; high - o2)
real(kind=kind_phys), dimension(ng29), public absco2
co2
real(kind=kind_phys), parameter con_amd
molecular wght of dry air ( )
real(kind=kind_phys), dimension(ng25), public rayl
rayleigh extinction coefficient
real(kind=kind_phys), parameter, public rayl
rayleigh extinction coefficient at