407 module module_radsw_main
410 use physparam
, only : iswrate, iswrgas, iswcliq, iswcice, &
411 & isubcsw, icldflg, iovrsw, ivflip, &
413 use physcons
, only : con_g, con_cp, con_avgd, con_amd, &
417 use mersenne_twister
, only : random_setseed, random_number, &
427 character(40),
parameter :: &
428 & VTAGSW=
'NCEP SW v5.1 Nov 2012 -RRTMG-SW v3.8 ' 440 real (kind=kind_phys),
parameter :: eps = 1.0e-6
441 real (kind=kind_phys),
parameter :: oneminus= 1.0 - eps
443 real (kind=kind_phys),
parameter ::
bpade = 1.0/0.278
444 real (kind=kind_phys),
parameter :: stpfac = 296.0/1013.0
445 real (kind=kind_phys),
parameter :: ftiny = 1.0e-12
446 real (kind=kind_phys),
parameter :: flimit = 1.0e-20
448 real (kind=kind_phys),
parameter ::
s0 = 1368.22
450 real (kind=kind_phys),
parameter :: f_zero = 0.0
451 real (kind=kind_phys),
parameter :: f_one = 1.0
454 real (kind=kind_phys),
parameter :: amdw = con_amd/con_amw
455 real (kind=kind_phys),
parameter :: amdo3 = con_amd/con_amo3
458 integer,
dimension(nblow:nbhgh) :: nspa, nspb
460 integer,
dimension(nblow:nbhgh) ::
idxsfc 462 integer,
dimension(nblow:nbhgh) ::
idxebc 464 data nspa(:) / 9, 9, 9, 9, 1, 9, 9, 1, 9, 1, 0, 1, 9, 1 /
465 data nspb(:) / 1, 5, 1, 1, 1, 5, 1, 0, 1, 0, 0, 1, 5, 1 /
468 data idxsfc(:) / 1, 1, 1, 1, 1, 1, 1, 1, 0, 2, 2, 2, 2, 1 /
469 data idxebc(:) / 5, 5, 4, 4, 3, 3, 2, 2, 1, 1, 1, 1, 1, 5 /
485 integer,
parameter ::
nuvb = 27
488 logical :: lhswb = .false.
489 logical :: lhsw0 = .false.
490 logical :: lflxprf= .false.
491 logical :: lfdncmp= .false.
587 & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr, &
588 & clouds,icseed,aerosols,sfcalb, &
589 & cosz,solcon,nday,idxday, &
590 & npts, nlay, nlp1, lprnt, &
591 & hswc,topflx,sfcflx, &
592 & hsw0,hswb,flxprf,fdncmp &
771 integer,
intent(in) :: npts, nlay, nlp1, NDAY
773 integer,
dimension(:),
intent(in) :: idxday, icseed
775 logical,
intent(in) :: lprnt
777 real (kind=kind_phys),
dimension(npts,nlp1),
intent(in) :: &
779 real (kind=kind_phys),
dimension(npts,nlay),
intent(in) :: &
780 & plyr, tlyr, qlyr, olyr
781 real (kind=kind_phys),
dimension(npts,4),
intent(in) :: sfcalb
783 real (kind=kind_phys),
dimension(npts,nlay,9),
intent(in):: gasvmr
784 real (kind=kind_phys),
dimension(npts,nlay,11):: clouds
785 real (kind=kind_phys),
dimension(npts,nlay,nbdsw,3),
intent(in):: &
788 real (kind=kind_phys),
intent(in) :: cosz(npts), solcon
791 real (kind=kind_phys),
dimension(npts,nlay),
intent(out) :: hswc
793 type(
topfsw_type),
dimension(npts),
intent(out) :: topflx
794 type(
sfcfsw_type),
dimension(npts),
intent(out) :: sfcflx
797 real (kind=kind_phys),
dimension(npts,nlay,nbdsw),
optional, &
798 & intent(out) :: hswb
800 real (kind=kind_phys),
dimension(npts,nlay),
optional, &
801 & intent(out) :: hsw0
802 type (
profsw_type),
dimension(npts,nlp1),
optional, &
803 & intent(out) :: flxprf
805 & intent(out) :: fdncmp
808 real (kind=kind_phys),
dimension(nlay,ngptsw) :: cldfmc, &
810 real (kind=kind_phys),
dimension(nlp1,nbdsw):: fxupc, fxdnc, &
813 real (kind=kind_phys),
dimension(nlay,nbdsw) :: &
814 & tauae, ssaae, asyae, taucw, ssacw, asycw
816 real (kind=kind_phys),
dimension(ngptsw) :: sfluxzen
818 real (kind=kind_phys),
dimension(nlay) :: cldfrc, delp, &
819 & pavel, tavel, coldry, colmol, h2ovmr, o3vmr, temcol, &
820 & cliqp, reliq, cicep, reice, cdat1, cdat2, cdat3, cdat4, &
821 & cfrac, fac00, fac01, fac10, fac11, forfac, forfrac, &
822 & selffac, selffrac, rfdelp
824 real (kind=kind_phys),
dimension(nlp1) :: fnet, flxdc, flxuc, &
827 real (kind=kind_phys),
dimension(2) :: albbm, albdf, sfbmc, &
828 & sfbm0, sfdfc, sfdf0
830 real (kind=kind_phys) :: cosz1, sntz1, tem0, tem1, tem2, s0fac, &
831 & ssolar, zcf0, zcf1, ftoau0, ftoauc, ftoadc, &
832 & fsfcu0, fsfcuc, fsfcd0, fsfcdc, suvbfc, suvbf0
836 real (kind=kind_phys) :: colamt(nlay,
maxgas)
838 integer,
dimension(npts) :: ipseed
839 integer,
dimension(nlay) :: indfor, indself, jp, jt, jt1
841 integer :: i, ib, ipt, j1, k, kk, laytrop, mb
846 lhswb =
present ( hswb )
847 lhsw0 =
present ( hsw0 )
848 lflxprf=
present ( flxprf )
849 lfdncmp=
present ( fdncmp )
861 sfcflx =
sfcfsw_type( f_zero, f_zero, f_zero, f_zero )
865 flxprf =
profsw_type( f_zero, f_zero, f_zero, f_zero )
869 fdncmp =
cmpfsw_type(f_zero,f_zero,f_zero,f_zero,f_zero,f_zero)
883 if ( isubcsw == 1 )
then 887 elseif ( isubcsw == 2 )
then 889 ipseed(i) = icseed(i)
894 write(0,*)
' In radsw, isubcsw, ipsdsw0,ipseed =', &
900 lab_do_ipt :
do ipt = 1, nday
905 sntz1 = f_one / cosz(j1)
906 ssolar = s0fac * cosz(j1)
909 albbm(1) = sfcalb(j1,1)
910 albdf(1) = sfcalb(j1,2)
911 albbm(2) = sfcalb(j1,3)
912 albdf(2) = sfcalb(j1,4)
917 if (ivflip == 0)
then 920 tem2 = 1.0e-20 * 1.0e3 * con_avgd
924 pavel(k) = plyr(j1,kk)
925 tavel(k) = tlyr(j1,kk)
926 delp(k) = plvl(j1,kk+1) - plvl(j1,kk)
938 h2ovmr(k)= max(f_zero,qlyr(j1,kk)*amdw/(f_one-qlyr(j1,kk)))
939 o3vmr(k)= max(f_zero,olyr(j1,kk)*amdo3)
941 tem0 = (f_one - h2ovmr(k))*con_amd + h2ovmr(k)*con_amw
942 coldry(k) = tem2 * delp(k) / (tem1*tem0*(f_one + h2ovmr(k)))
943 temcol(k) = 1.0e-12 * coldry(k)
945 colamt(k,1) = max(f_zero, coldry(k)*h2ovmr(k))
946 colamt(k,2) = max(temcol(k), coldry(k)*gasvmr(j1,kk,1))
947 colamt(k,3) = max(f_zero, coldry(k)*o3vmr(k))
948 colmol(k) = coldry(k) + colamt(k,1)
954 if (iswrgas > 0)
then 957 colamt(k,4) = max(temcol(k), coldry(k)*gasvmr(j1,kk,2))
958 colamt(k,5) = max(temcol(k), coldry(k)*gasvmr(j1,kk,3))
959 colamt(k,6) = max(temcol(k), coldry(k)*gasvmr(j1,kk,4))
964 colamt(k,4) = temcol(k)
965 colamt(k,5) = temcol(k)
966 colamt(k,6) = temcol(k)
976 tauae(k,ib) = aerosols(j1,kk,ib,1)
977 ssaae(k,ib) = aerosols(j1,kk,ib,2)
978 asyae(k,ib) = aerosols(j1,kk,ib,3)
983 if (iswcliq > 0)
then 986 cfrac(k) = clouds(j1,kk,1)
987 cliqp(k) = clouds(j1,kk,2)
988 reliq(k) = clouds(j1,kk,3)
989 cicep(k) = clouds(j1,kk,4)
990 reice(k) = clouds(j1,kk,5)
991 cdat1(k) = clouds(j1,kk,6)
992 cdat2(k) = clouds(j1,kk,7)
993 cdat3(k) = clouds(j1,kk,8)
994 cdat4(k) = clouds(j1,kk,9)
999 cfrac(k) = clouds(j1,kk,1)
1000 cdat1(k) = clouds(j1,kk,2)
1001 cdat2(k) = clouds(j1,kk,3)
1002 cdat3(k) = clouds(j1,kk,4)
1008 tem1 = 100.0 * con_g
1009 tem2 = 1.0e-20 * 1.0e3 * con_avgd
1012 pavel(k) = plyr(j1,k)
1013 tavel(k) = tlyr(j1,k)
1014 delp(k) = plvl(j1,k) - plvl(j1,k+1)
1022 h2ovmr(k)= max(f_zero,qlyr(j1,k)*amdw/(f_one-qlyr(j1,k)))
1023 o3vmr(k)= max(f_zero,olyr(j1,k)*amdo3)
1025 tem0 = (f_one - h2ovmr(k))*con_amd + h2ovmr(k)*con_amw
1026 coldry(k) = tem2 * delp(k) / (tem1*tem0*(f_one + h2ovmr(k)))
1027 temcol(k) = 1.0e-12 * coldry(k)
1029 colamt(k,1) = max(f_zero, coldry(k)*h2ovmr(k))
1030 colamt(k,2) = max(temcol(k), coldry(k)*gasvmr(j1,k,1))
1031 colamt(k,3) = max(f_zero, coldry(k)*o3vmr(k))
1032 colmol(k) = coldry(k) + colamt(k,1)
1038 write(0,*)
' pavel=',pavel
1039 write(0,*)
' tavel=',tavel
1040 write(0,*)
' delp=',delp
1041 write(0,*)
' h2ovmr=',h2ovmr*1000
1042 write(0,*)
' o3vmr=',o3vmr*1000000
1049 if (iswrgas > 0)
then 1051 colamt(k,4) = max(temcol(k), coldry(k)*gasvmr(j1,k,2))
1052 colamt(k,5) = max(temcol(k), coldry(k)*gasvmr(j1,k,3))
1053 colamt(k,6) = max(temcol(k), coldry(k)*gasvmr(j1,k,4))
1058 colamt(k,4) = temcol(k)
1059 colamt(k,5) = temcol(k)
1060 colamt(k,6) = temcol(k)
1069 tauae(k,ib) = aerosols(j1,k,ib,1)
1070 ssaae(k,ib) = aerosols(j1,k,ib,2)
1071 asyae(k,ib) = aerosols(j1,k,ib,3)
1075 if (iswcliq > 0)
then 1077 cfrac(k) = clouds(j1,k,1)
1078 cliqp(k) = clouds(j1,k,2)
1079 reliq(k) = clouds(j1,k,3)
1080 cicep(k) = clouds(j1,k,4)
1081 reice(k) = clouds(j1,k,5)
1082 cdat1(k) = clouds(j1,k,6)
1083 cdat2(k) = clouds(j1,k,7)
1084 cdat3(k) = clouds(j1,k,8)
1085 cdat4(k) = clouds(j1,k,9)
1089 cfrac(k) = clouds(j1,k,1)
1090 cdat1(k) = clouds(j1,k,2)
1091 cdat2(k) = clouds(j1,k,3)
1092 cdat3(k) = clouds(j1,k,4)
1105 if (iovrsw == 0)
then 1107 zcf0 = zcf0 * (f_one - cfrac(k))
1109 else if (iovrsw == 1)
then 1111 if (cfrac(k) > ftiny)
then 1112 zcf1 = min( zcf1, f_one-cfrac(k) )
1113 elseif (zcf1 < f_one)
then 1119 else if (iovrsw == 2)
then 1121 zcf0 = min( zcf0, f_one-cfrac(k) )
1125 if (zcf0 <= ftiny) zcf0 = f_zero
1126 if (zcf0 > oneminus) zcf0 = f_one
1132 if (zcf1 > f_zero)
then 1136 & ( cfrac,cliqp,reliq,cicep,reice,cdat1,cdat2,cdat3,cdat4, &
1137 & zcf1, nlay, ipseed(j1), &
1139 & taucw, ssacw, asycw, cldfrc, cldfmc &
1154 clouds(j1,k,10) = taucw(k,10)
1160 & ( pavel,tavel,h2ovmr, nlay,nlp1, &
1162 & laytrop,jp,jt,jt1,fac00,fac01,fac10,fac11, &
1163 & selffac,selffrac,indself,forfac,forfrac,indfor &
1170 & ( colamt,colmol,fac00,fac01,fac10,fac11,jp,jt,jt1,laytrop, &
1171 & forfac,forfrac,indfor,selffac,selffrac,indself, nlay, &
1173 & sfluxzen, taug, taur &
1182 if ( isubcsw <= 0 )
then 1186 & ( ssolar,cosz1,sntz1,albbm,albdf,sfluxzen,cldfrc, &
1187 & zcf1,zcf0,taug,taur,tauae,ssaae,asyae,taucw,ssacw,asycw, &
1190 & fxupc,fxdnc,fxup0,fxdn0, &
1191 & ftoauc,ftoau0,ftoadc,fsfcuc,fsfcu0,fsfcdc,fsfcd0, &
1192 & sfbmc,sfdfc,sfbm0,sfdf0,suvbfc,suvbf0 &
1199 & ( ssolar,cosz1,sntz1,albbm,albdf,sfluxzen,cldfmc, &
1200 & zcf1,zcf0,taug,taur,tauae,ssaae,asyae,taucw,ssacw,asycw, &
1203 & fxupc,fxdnc,fxup0,fxdn0, &
1204 & ftoauc,ftoau0,ftoadc,fsfcuc,fsfcu0,fsfcdc,fsfcd0, &
1205 & sfbmc,sfdfc,sfbm0,sfdf0,suvbfc,suvbf0 &
1218 flxuc(k) = flxuc(k) + fxupc(k,ib)
1219 flxdc(k) = flxdc(k) + fxdnc(k,ib)
1225 if ( lhsw0 .or. lflxprf )
then 1231 flxu0(k) = flxu0(k) + fxup0(k,ib)
1232 flxd0(k) = flxd0(k) + fxdn0(k,ib)
1245 fdncmp(j1)%uvbf0 = suvbf0
1246 fdncmp(j1)%uvbfc = suvbfc
1249 fdncmp(j1)%nirbm = sfbmc(1)
1250 fdncmp(j1)%nirdf = sfdfc(1)
1251 fdncmp(j1)%visbm = sfbmc(2)
1252 fdncmp(j1)%visdf = sfdfc(2)
1257 topflx(j1)%upfxc = ftoauc
1258 topflx(j1)%dnfxc = ftoadc
1259 topflx(j1)%upfx0 = ftoau0
1261 sfcflx(j1)%upfxc = fsfcuc
1262 sfcflx(j1)%dnfxc = fsfcdc
1263 sfcflx(j1)%upfx0 = fsfcu0
1264 sfcflx(j1)%dnfx0 = fsfcd0
1266 if (ivflip == 0)
then 1270 fnet(1) = flxdc(1) - flxuc(1)
1274 fnet(k) = flxdc(k) - flxuc(k)
1275 hswc(j1,kk) = (fnet(k)-fnet(k-1)) * rfdelp(k-1)
1283 flxprf(j1,kk)%upfxc = flxuc(k)
1284 flxprf(j1,kk)%dnfxc = flxdc(k)
1285 flxprf(j1,kk)%upfx0 = flxu0(k)
1286 flxprf(j1,kk)%dnfx0 = flxd0(k)
1293 fnet(1) = flxd0(1) - flxu0(1)
1297 fnet(k) = flxd0(k) - flxu0(k)
1298 hsw0(j1,kk) = (fnet(k)-fnet(k-1)) * rfdelp(k-1)
1306 fnet(1) = fxdnc(1,mb) - fxupc(1,mb)
1310 fnet(k) = fxdnc(k,mb) - fxupc(k,mb)
1311 hswb(j1,kk,mb) = (fnet(k) - fnet(k-1)) * rfdelp(k-1)
1320 fnet(1) = flxdc(1) - flxuc(1)
1323 fnet(k) = flxdc(k) - flxuc(k)
1324 hswc(j1,k-1) = (fnet(k)-fnet(k-1)) * rfdelp(k-1)
1331 flxprf(j1,k)%upfxc = flxuc(k)
1332 flxprf(j1,k)%dnfxc = flxdc(k)
1333 flxprf(j1,k)%upfx0 = flxu0(k)
1334 flxprf(j1,k)%dnfx0 = flxd0(k)
1341 fnet(1) = flxd0(1) - flxu0(1)
1344 fnet(k) = flxd0(k) - flxu0(k)
1345 hsw0(j1,k-1) = (fnet(k)-fnet(k-1)) * rfdelp(k-1)
1353 fnet(1) = fxdnc(1,mb) - fxupc(1,mb)
1356 fnet(k+1) = fxdnc(k+1,mb) - fxupc(k+1,mb)
1357 hswb(j1,k,mb) = (fnet(k+1) - fnet(k)) * rfdelp(k)
1368 end subroutine swrad 1433 integer,
intent(in) :: me
1438 real (kind=kind_phys),
parameter :: expeps = 1.e-20
1442 real (kind=kind_phys) :: tfn, tau
1447 if ( iovrsw<0 .or. iovrsw>2 )
then 1448 print *,
' *** Error in specification of cloud overlap flag', &
1449 &
' IOVRSW=',iovrsw,
' in RSWINIT !!' 1454 print *,
' - Using AER Shortwave Radiation, Version: ',vtagsw
1456 if (iswmode == 1)
then 1457 print *,
' --- Delta-eddington 2-stream transfer scheme' 1458 else if (iswmode == 2)
then 1459 print *,
' --- PIFM 2-stream transfer scheme' 1460 else if (iswmode == 3)
then 1461 print *,
' --- Discrete ordinates 2-stream transfer scheme' 1464 if (iswrgas <= 0)
then 1465 print *,
' --- Rare gases absorption is NOT included in SW' 1467 print *,
' --- Include rare gases N2O, CH4, O2, absorptions',&
1471 if ( isubcsw == 0 )
then 1472 print *,
' --- Using standard grid average clouds, no ', &
1473 &
'sub-column clouds approximation applied' 1474 elseif ( isubcsw == 1 )
then 1475 print *,
' --- Using MCICA sub-colum clouds approximation ', &
1476 &
'with a prescribed sequence of permutation seeds' 1477 elseif ( isubcsw == 2 )
then 1478 print *,
' --- Using MCICA sub-colum clouds approximation ', &
1479 &
'with provided input array of permutation seeds' 1481 print *,
' *** Error in specification of sub-column cloud ', &
1482 &
' control flag isubcsw =',isubcsw,
' !!' 1489 if ((icldflg == 0 .and. iswcliq /= 0) .or. &
1490 & (icldflg == 1 .and. iswcliq == 0))
then 1491 print *,
' *** Model cloud scheme inconsistent with SW', &
1492 &
' radiation cloud radiative property setup !!' 1499 if (iswrate == 1)
then 1502 heatfac = con_g * 864.0 / con_cp
1504 heatfac = con_g * 1.0e-2 / con_cp
1518 tfn = float(i) / float(
ntbmx-i)
1564 & ( cfrac,cliqp,reliq,cicep,reice,cdat1,cdat2,cdat3,cdat4, &
1565 & cf1, nlay, ipseed, &
1566 & taucw, ssacw, asycw, cldfrc, cldfmc &
1647 integer,
intent(in) :: nlay, ipseed
1648 real (kind=kind_phys),
intent(in) :: cf1
1650 real (kind=kind_phys),
dimension(nlay),
intent(in) :: cliqp, &
1651 & reliq, cicep, reice, cdat1, cdat2, cdat3, cdat4, cfrac
1654 real (kind=kind_phys),
dimension(nlay,ngptsw),
intent(out) :: &
1656 real (kind=kind_phys),
dimension(nlay,nbdsw),
intent(out) :: &
1657 & taucw, ssacw, asycw
1658 real (kind=kind_phys),
dimension(nlay),
intent(out) :: cldfrc
1661 real (kind=kind_phys),
dimension(nblow:nbhgh) :: tauliq, tauice, &
1662 & ssaliq, ssaice, ssaran, ssasnw, asyliq, asyice, &
1664 real (kind=kind_phys),
dimension(nlay) :: cldf
1666 real (kind=kind_phys) :: dgeice, factor, fint, tauran, tausnw, &
1667 & cldliq, refliq, cldice, refice, cldran, cldsnw, refsnw, &
1668 & extcoliq, ssacoliq, asycoliq, extcoice, ssacoice, asycoice,&
1671 logical :: lcloudy(nlay,
ngptsw)
1672 integer :: ia, ib, ig, jb, k, index
1679 taucw(k,ib) = f_zero
1681 asycw(k,ib) = f_zero
1687 lab_if_iswcliq :
if (iswcliq > 0)
then 1689 lab_do_k :
do k = 1, nlay
1690 lab_if_cld :
if (cfrac(k) > ftiny)
then 1708 dgesnw = 1.0315 * refsnw
1710 tauran = cldran *
a0r 1718 if (cldsnw>f_zero .and. refsnw>10.0_kind_phys)
then 1720 tausnw = cldsnw*1.09087*(
a0s + a1s/dgesnw)
1726 ssaran(ib) = tauran * (f_one -
b0r(ib))
1727 ssasnw(ib) = tausnw * (f_one - (
b0s(ib)+b1s(ib)*dgesnw))
1728 asyran(ib) = ssaran(ib) *
c0r(ib)
1729 asysnw(ib) = ssasnw(ib) *
c0s(ib)
1739 if ( cldliq <= f_zero )
then 1746 if ( iswcliq == 1 )
then 1747 factor = refliq - 1.5
1748 index = max( 1, min( 57, int( factor ) ))
1749 fint = factor - float(index)
1752 extcoliq = max(f_zero,
extliq1(index,ib) &
1754 ssacoliq = max(f_zero, min(f_one,
ssaliq1(index,ib) &
1757 asycoliq = max(f_zero, min(f_one,
asyliq1(index,ib) &
1761 tauliq(ib) = cldliq * extcoliq
1762 ssaliq(ib) = tauliq(ib) * ssacoliq
1763 asyliq(ib) = ssaliq(ib) * asycoliq
1770 if ( cldice <= f_zero )
then 1781 if ( iswcice == 1 )
then 1782 refice = min(130.0_kind_phys,max(13.0_kind_phys,refice))
1787 extcoice = max(f_zero,
abari(ia)+
bbari(ia)/refice )
1788 ssacoice = max(f_zero, min(f_one, &
1790 asycoice = max(f_zero, min(f_one, &
1794 tauice(ib) = cldice * extcoice
1795 ssaice(ib) = tauice(ib) * ssacoice
1796 asyice(ib) = ssaice(ib) * asycoice
1801 elseif ( iswcice == 2 )
then 1802 refice = min(131.0_kind_phys,max(5.0_kind_phys,refice))
1804 factor = (refice - 2.0) / 3.0
1805 index = max( 1, min( 42, int( factor ) ))
1806 fint = factor - float(index)
1809 extcoice = max(f_zero,
extice2(index,ib) &
1811 ssacoice = max(f_zero, min(f_one,
ssaice2(index,ib) &
1813 asycoice = max(f_zero, min(f_one,
asyice2(index,ib) &
1817 tauice(ib) = cldice * extcoice
1818 ssaice(ib) = tauice(ib) * ssacoice
1819 asyice(ib) = ssaice(ib) * asycoice
1825 elseif ( iswcice == 3 )
then 1826 dgeice = max( 5.0, min( 140.0, 1.0315*refice ))
1828 factor = (dgeice - 2.0) / 3.0
1829 index = max( 1, min( 45, int( factor ) ))
1830 fint = factor - float(index)
1833 extcoice = max(f_zero,
extice3(index,ib) &
1835 ssacoice = max(f_zero, min(f_one,
ssaice3(index,ib) &
1837 asycoice = max(f_zero, min(f_one,
asyice3(index,ib) &
1843 tauice(ib) = cldice * extcoice
1844 ssaice(ib) = tauice(ib) * ssacoice
1845 asyice(ib) = ssaice(ib) * asycoice
1853 taucw(k,ib) = tauliq(jb)+tauice(jb)+tauran+tausnw
1854 ssacw(k,ib) = ssaliq(jb)+ssaice(jb)+ssaran(jb)+ssasnw(jb)
1855 asycw(k,ib) = asyliq(jb)+asyice(jb)+asyran(jb)+asysnw(jb)
1864 if (cfrac(k) > ftiny)
then 1866 taucw(k,ib) = cdat1(k)
1867 ssacw(k,ib) = cdat1(k) * cdat2(k)
1868 asycw(k,ib) = ssacw(k,ib) * cdat3(k)
1873 endif lab_if_iswcliq
1878 if ( isubcsw > 0 )
then 1881 where (cldf(:) < ftiny)
1889 & ( cldf, nlay, ipseed, &
1896 if ( lcloudy(k,ig) )
then 1897 cldfmc(k,ig) = f_one
1899 cldfmc(k,ig) = f_zero
1907 cldfrc(k) = cfrac(k) / cf1
1925 & ( cldf, nlay, ipseed, &
1952 integer,
intent(in) :: nlay, ipseed
1954 real (kind=kind_phys),
dimension(nlay),
intent(in) :: cldf
1957 logical,
dimension(nlay,ngptsw),
intent(out):: lcloudy
1960 real (kind=kind_phys) :: cdfunc(nlay,
ngptsw), tem1, &
1961 & rand2d(nlay*ngptsw), rand1d(ngptsw)
1963 type(random_stat) :: stat
1971 call random_setseed &
1980 select case ( iovrsw )
1984 call random_number &
1993 cdfunc(k,n) = rand2d(k1)
1999 call random_number &
2008 cdfunc(k,n) = rand2d(k1)
2020 tem1 = f_one - cldf(k1)
2023 if ( cdfunc(k1,n) > tem1 )
then 2024 cdfunc(k,n) = cdfunc(k1,n)
2026 cdfunc(k,n) = cdfunc(k,n) * tem1
2051 call random_number &
2069 tem1 = f_one - cldf(k)
2072 lcloudy(k,n) = cdfunc(k,n) >= tem1
2107 & ( pavel,tavel,h2ovmr, nlay,nlp1, &
2108 & laytrop,jp,jt,jt1,fac00,fac01,fac10,fac11, &
2109 & selffac,selffrac,indself,forfac,forfrac,indfor &
2149 integer,
intent(in) :: nlay, nlp1
2151 real (kind=kind_phys),
dimension(:),
intent(in) :: pavel, tavel, &
2155 integer,
dimension(nlay),
intent(out) :: indself, indfor, &
2157 integer,
intent(out) :: laytrop
2159 real (kind=kind_phys),
dimension(nlay),
intent(out) :: fac00, &
2160 & fac01, fac10, fac11, selffac, selffrac, forfac, forfrac
2163 real (kind=kind_phys) :: plog, fp, fp1, ft, ft1, tem1, tem2
2165 integer :: i, k, jp1
2173 forfac(k) = pavel(k)*stpfac / (tavel(k)*(f_one + h2ovmr(k)))
2180 plog = log(pavel(k))
2181 jp(k) = max(1, min(58, int(36.0 - 5.0*(plog+0.04)) ))
2183 fp = 5.0 * (
preflog(jp(k)) - plog)
2192 tem1 = (tavel(k) -
tref(jp(k))) / 15.0
2193 tem2 = (tavel(k) -
tref(jp1 )) / 15.0
2194 jt(k) = max(1, min(4, int(3.0 + tem1) ))
2195 jt1(k) = max(1, min(4, int(3.0 + tem2) ))
2196 ft = tem1 - float(jt(k) - 3)
2197 ft1 = tem2 - float(jt1(k) - 3)
2208 fac00(k) = fp1 * (f_one - ft)
2210 fac01(k) = fp * (f_one - ft1)
2215 if ( plog > 4.56 )
then 2222 tem1 = (332.0 - tavel(k)) / 36.0
2223 indfor(k) = min(2, max(1, int(tem1)))
2224 forfrac(k) = tem1 - float(indfor(k))
2229 tem2 = (tavel(k) - 188.0) / 7.2
2230 indself(k) = min(9, max(1, int(tem2)-7))
2231 selffrac(k) = tem2 - float(indself(k) + 7)
2232 selffac(k) = h2ovmr(k) * forfac(k)
2239 tem1 = (tavel(k) - 188.0) / 36.0
2241 forfrac(k) = tem1 - f_one
2244 selffrac(k) = f_zero
2297 & ( ssolar,cosz,sntz,albbm,albdf,sfluxzen,cldfrc, &
2298 & cf1,cf0,taug,taur,tauae,ssaae,asyae,taucw,ssacw,asycw, &
2300 & fxupc,fxdnc,fxup0,fxdn0, &
2301 & ftoauc,ftoau0,ftoadc,fsfcuc,fsfcu0,fsfcdc,fsfcd0, &
2302 & sfbmc,sfdfc,sfbm0,sfdf0,suvbfc,suvbf0 &
2396 real (kind=kind_phys),
parameter :: zcrit = 0.9999995
2397 real (kind=kind_phys),
parameter :: zsr3 = sqrt(3.0)
2398 real (kind=kind_phys),
parameter :: od_lo = 0.06
2399 real (kind=kind_phys),
parameter :: eps1 = 1.0e-8
2402 integer,
intent(in) :: nlay, nlp1
2404 real (kind=kind_phys),
dimension(nlay,ngptsw),
intent(in) :: &
2406 real (kind=kind_phys),
dimension(nlay,nbdsw),
intent(in) :: &
2407 & taucw, ssacw, asycw, tauae, ssaae, asyae
2409 real (kind=kind_phys),
dimension(ngptsw),
intent(in) :: sfluxzen
2410 real (kind=kind_phys),
dimension(nlay),
intent(in) :: cldfrc
2412 real (kind=kind_phys),
dimension(2),
intent(in) :: albbm, albdf
2414 real (kind=kind_phys),
intent(in) :: cosz, sntz, cf1, cf0, ssolar
2417 real (kind=kind_phys),
dimension(nlp1,nbdsw),
intent(out) :: &
2418 & fxupc, fxdnc, fxup0, fxdn0
2420 real (kind=kind_phys),
dimension(2),
intent(out) :: sfbmc, sfdfc, &
2423 real (kind=kind_phys),
intent(out) :: suvbfc, suvbf0, ftoadc, &
2424 & ftoauc, ftoau0, fsfcuc, fsfcu0, fsfcdc, fsfcd0
2427 real (kind=kind_phys),
dimension(nlay) :: ztaus, zssas, zasys, &
2430 real (kind=kind_phys),
dimension(nlp1) :: zrefb, zrefd, ztrab, &
2431 & ztrad, ztdbt, zldbt, zfu, zfd
2433 real (kind=kind_phys) :: ztau1, zssa1, zasy1, ztau0, zssa0, &
2434 & zasy0, zasy3, zssaw, zasyw, zgam1, zgam2, zgam3, zgam4, &
2435 & zc0, zc1, za1, za2, zb1, zb2, zrk, zrk2, zrp, zrp1, zrm1, &
2436 & zrpp, zrkg1, zrkg3, zrkg4, zexp1, zexm1, zexp2, zexm2, &
2437 & zexp3, zexp4, zden1, ze1r45, ftind, zsolar, zrefb1, &
2438 & zrefd1, ztrab1, ztrad1, ztdbt0, zr1, zr2, zr3, zr4, zr5, &
2439 & zt1, zt2, zt3, zf1, zf2, zrpp1
2441 integer :: ib, ibd, jb, jg, k, kp, itind
2448 fxdnc(k,ib) = f_zero
2449 fxupc(k,ib) = f_zero
2450 fxdn0(k,ib) = f_zero
2451 fxup0(k,ib) = f_zero
2479 lab_do_jg :
do jg = 1,
ngptsw 2485 zsolar = ssolar * sfluxzen(jg)
2494 zrefb(1) = albbm(ibd)
2495 zrefd(1) = albdf(ibd)
2497 zrefb(1) = 0.5 * (albbm(1) + albbm(2))
2498 zrefd(1) = 0.5 * (albdf(1) + albdf(2))
2517 ztau0 = max( ftiny, taur(k,jg)+taug(k,jg)+tauae(k,ib) )
2518 zssa0 = taur(k,jg) + tauae(k,ib)*ssaae(k,ib)
2519 zasy0 = asyae(k,ib)*ssaae(k,ib)*tauae(k,ib)
2520 zssaw = min( oneminus, zssa0 / ztau0 )
2521 zasyw = zasy0 / max( ftiny, zssa0 )
2532 ztau1 = (f_one - za2) * ztau0
2533 zssa1 = (zssaw - za2) / (f_one - za2)
2535 zasy1 = zasyw / (f_one + zasyw)
2536 zasy3 = 0.75 * zasy1
2539 if ( iswmode == 1 )
then 2540 zgam1 = 1.75 - zssa1 * (f_one + zasy3)
2541 zgam2 =-0.25 + zssa1 * (f_one - zasy3)
2542 zgam3 = 0.5 - zasy3 * cosz
2543 elseif ( iswmode == 2 )
then 2544 zgam1 = 2.0 - zssa1 * (1.25 + zasy3)
2545 zgam2 = 0.75* zssa1 * (f_one- zasy1)
2546 zgam3 = 0.5 - zasy3 * cosz
2547 elseif ( iswmode == 3 )
then 2548 zgam1 = zsr3 * (2.0 - zssa1 * (1.0 + zasy1)) * 0.5
2549 zgam2 = zsr3 * zssa1 * (1.0 - zasy1) * 0.5
2550 zgam3 = (1.0 - zsr3 * zasy1 * cosz) * 0.5
2552 zgam4 = f_one - zgam3
2556 if ( zssaw >= zcrit )
then 2557 za1 = zgam1 * cosz - zgam3
2563 zb1 = min( ztau1*sntz , 500.0 )
2564 if ( zb1 <= od_lo )
then 2565 zb2 = f_one - zb1 + 0.5*zb1*zb1
2567 ftind = zb1 / (
bpade + zb1)
2568 itind = ftind*
ntbmx + 0.5
2573 zrefb(kp) = max(f_zero, min(f_one, &
2574 & (za2 - za1*(f_one - zb2))/(f_one + za2) ))
2575 ztrab(kp) = max(f_zero, min(f_one, f_one-zrefb(kp) ))
2578 zrefd(kp) = max(f_zero, min(f_one, za2/(f_one + za2) ))
2579 ztrad(kp) = max(f_zero, min(f_one, f_one-zrefd(kp) ))
2582 za1 = zgam1*zgam4 + zgam2*zgam3
2583 za2 = zgam1*zgam3 + zgam2*zgam4
2584 zrk = sqrt( (zgam1 - zgam2) * (zgam1 + zgam2) )
2590 zrpp1= f_one - zrp*zrp
2591 zrpp = sign( max(flimit, abs(zrpp1)), zrpp1 )
2596 zr1 = zrm1 * (za2 + zrkg3)
2597 zr2 = zrp1 * (za2 - zrkg3)
2598 zr3 = zrk2 * (zgam3 - za2*cosz)
2600 zr5 = zrpp * (zrk - zgam1)
2602 zt1 = zrp1 * (za1 + zrkg4)
2603 zt2 = zrm1 * (za1 - zrkg4)
2604 zt3 = zrk2 * (zgam4 + za1*cosz)
2609 zb1 = min( zrk*ztau1, 500.0 )
2610 if ( zb1 <= od_lo )
then 2611 zexm1 = f_one - zb1 + 0.5*zb1*zb1
2613 ftind = zb1 / (
bpade + zb1)
2614 itind = ftind*
ntbmx + 0.5
2617 zexp1 = f_one / zexm1
2619 zb2 = min( sntz*ztau1, 500.0 )
2620 if ( zb2 <= od_lo )
then 2621 zexm2 = f_one - zb2 + 0.5*zb2*zb2
2623 ftind = zb2 / (
bpade + zb2)
2624 itind = ftind*
ntbmx + 0.5
2627 zexp2 = f_one / zexm2
2628 ze1r45 = zr4*zexp1 + zr5*zexm1
2631 if (ze1r45>=-eps1 .and. ze1r45<=eps1)
then 2635 zden1 = zssa1 / ze1r45
2636 zrefb(kp) = max(f_zero, min(f_one, &
2637 & (zr1*zexp1 - zr2*zexm1 - zr3*zexm2)*zden1 ))
2638 ztrab(kp) = max(f_zero, min(f_one, zexm2*(f_one &
2639 & - (zt1*zexp1 - zt2*zexm1 - zt3*zexp2)*zden1) ))
2643 zden1 = zr4 / (ze1r45 * zrkg1)
2644 zrefd(kp) = max(f_zero, min(f_one, &
2645 & zgam2*(zexp1 - zexm1)*zden1 ))
2646 ztrad(kp) = max(f_zero, min(f_one, zrk2*zden1 ))
2654 if ( zr1 <= od_lo )
then 2655 zexp3 = f_one - zr1 + 0.5*zr1*zr1
2657 ftind = zr1 / (
bpade + zr1)
2658 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
2662 ztdbt(k) = zexp3 * ztdbt(kp)
2669 if ( zr1 <= od_lo )
then 2670 zexp4 = f_one - zr1 + 0.5*zr1*zr1
2672 ftind = zr1 / (
bpade + zr1)
2673 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
2678 ztdbt0 = zexp4 * ztdbt0
2683 & ( zrefb,zrefd,ztrab,ztrad,zldbt,ztdbt, &
2691 fxup0(k,ib) = fxup0(k,ib) + zsolar*zfu(k)
2692 fxdn0(k,ib) = fxdn0(k,ib) + zsolar*zfd(k)
2697 zb2 = zsolar*(zfd(1) - ztdbt0)
2700 sfbm0(ibd) = sfbm0(ibd) + zb1
2701 sfdf0(ibd) = sfdf0(ibd) + zb2
2705 sfbm0(1) = sfbm0(1) + zf1
2706 sfdf0(1) = sfdf0(1) + zf2
2707 sfbm0(2) = sfbm0(2) + zf1
2708 sfdf0(2) = sfdf0(2) + zf2
2723 if ( cf1 > eps )
then 2731 zc0 = f_one - cldfrc(k)
2733 if ( zc1 > ftiny )
then 2735 ztau0 = ztaus(k) + taucw(k,ib)
2736 zssa0 = zssas(k) + ssacw(k,ib)
2737 zasy0 = zasys(k) + asycw(k,ib)
2738 zssaw = min(oneminus, zssa0 / ztau0)
2739 zasyw = zasy0 / max(ftiny, zssa0)
2745 ztau1 = (f_one - za2) * ztau0
2746 zssa1 = (zssaw - za2) / (f_one - za2)
2748 zasy1 = zasyw / (f_one + zasyw)
2749 zasy3 = 0.75 * zasy1
2752 if ( iswmode == 1 )
then 2753 zgam1 = 1.75 - zssa1 * (f_one + zasy3)
2754 zgam2 =-0.25 + zssa1 * (f_one - zasy3)
2755 zgam3 = 0.5 - zasy3 * cosz
2756 elseif ( iswmode == 2 )
then 2757 zgam1 = 2.0 - zssa1 * (1.25 + zasy3)
2758 zgam2 = 0.75* zssa1 * (f_one- zasy1)
2759 zgam3 = 0.5 - zasy3 * cosz
2760 elseif ( iswmode == 3 )
then 2761 zgam1 = zsr3 * (2.0 - zssa1 * (1.0 + zasy1)) * 0.5
2762 zgam2 = zsr3 * zssa1 * (1.0 - zasy1) * 0.5
2763 zgam3 = (1.0 - zsr3 * zasy1 * cosz) * 0.5
2765 zgam4 = f_one - zgam3
2774 if ( zssaw >= zcrit )
then 2775 za1 = zgam1 * cosz - zgam3
2781 zb1 = min( ztau1*sntz , 500.0 )
2782 if ( zb1 <= od_lo )
then 2783 zb2 = f_one - zb1 + 0.5*zb1*zb1
2785 ftind = zb1 / (
bpade + zb1)
2786 itind = ftind*
ntbmx + 0.5
2791 zrefb(kp) = max(f_zero, min(f_one, &
2792 & (za2 - za1*(f_one - zb2))/(f_one + za2) ))
2793 ztrab(kp) = max(f_zero, min(f_one, f_one-zrefb(kp)))
2796 zrefd(kp) = max(f_zero, min(f_one, za2 / (f_one+za2) ))
2797 ztrad(kp) = max(f_zero, min(f_one, f_one - zrefd(kp) ))
2800 za1 = zgam1*zgam4 + zgam2*zgam3
2801 za2 = zgam1*zgam3 + zgam2*zgam4
2802 zrk = sqrt( (zgam1 - zgam2) * (zgam1 + zgam2) )
2808 zrpp1= f_one - zrp*zrp
2809 zrpp = sign( max(flimit, abs(zrpp1)), zrpp1 )
2814 zr1 = zrm1 * (za2 + zrkg3)
2815 zr2 = zrp1 * (za2 - zrkg3)
2816 zr3 = zrk2 * (zgam3 - za2*cosz)
2818 zr5 = zrpp * (zrk - zgam1)
2820 zt1 = zrp1 * (za1 + zrkg4)
2821 zt2 = zrm1 * (za1 - zrkg4)
2822 zt3 = zrk2 * (zgam4 + za1*cosz)
2827 zb1 = min( zrk*ztau1, 500.0 )
2828 if ( zb1 <= od_lo )
then 2829 zexm1 = f_one - zb1 + 0.5*zb1*zb1
2831 ftind = zb1 / (
bpade + zb1)
2832 itind = ftind*
ntbmx + 0.5
2835 zexp1 = f_one / zexm1
2837 zb2 = min( ztau1*sntz, 500.0 )
2838 if ( zb2 <= od_lo )
then 2839 zexm2 = f_one - zb2 + 0.5*zb2*zb2
2841 ftind = zb2 / (
bpade + zb2)
2842 itind = ftind*
ntbmx + 0.5
2845 zexp2 = f_one / zexm2
2846 ze1r45 = zr4*zexp1 + zr5*zexm1
2849 if ( ze1r45>=-eps1 .and. ze1r45<=eps1 )
then 2853 zden1 = zssa1 / ze1r45
2854 zrefb(kp) = max(f_zero, min(f_one, &
2855 & (zr1*zexp1-zr2*zexm1-zr3*zexm2)*zden1 ))
2856 ztrab(kp) = max(f_zero, min(f_one, zexm2*(f_one - &
2857 & (zt1*zexp1-zt2*zexm1-zt3*zexp2)*zden1) ))
2861 zden1 = zr4 / (ze1r45 * zrkg1)
2862 zrefd(kp) = max(f_zero, min(f_one, &
2863 & zgam2*(zexp1 - zexm1)*zden1 ))
2864 ztrad(kp) = max(f_zero, min(f_one, zrk2*zden1 ))
2870 zrefb(kp) = zc0*zrefb1 + zc1*zrefb(kp)
2871 zrefd(kp) = zc0*zrefd1 + zc1*zrefd(kp)
2872 ztrab(kp) = zc0*ztrab1 + zc1*ztrab(kp)
2873 ztrad(kp) = zc0*ztrad1 + zc1*ztrad(kp)
2880 if ( zr1 <= od_lo )
then 2881 zexp3 = f_one - zr1 + 0.5*zr1*zr1
2883 ftind = zr1 / (
bpade + zr1)
2884 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
2888 zldbt(kp) = zc0*zldbt(kp) + zc1*zexp3
2889 ztdbt(k) = zldbt(kp) * ztdbt(kp)
2895 if ( zr1 <= od_lo )
then 2896 zexp4 = f_one - zr1 + 0.5*zr1*zr1
2898 ftind = zr1 / (
bpade + zr1)
2899 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
2903 ztdbt0 = (zc0*zldbt0(k) + zc1*zexp4) * ztdbt0
2908 ztdbt(k) = zldbt(kp) * ztdbt(kp)
2911 ztdbt0 = zldbt0(k) * ztdbt0
2920 & ( zrefb,zrefd,ztrab,ztrad,zldbt,ztdbt, &
2928 fxupc(k,ib) = fxupc(k,ib) + zsolar*zfu(k)
2929 fxdnc(k,ib) = fxdnc(k,ib) + zsolar*zfd(k)
2935 zb2 = zsolar*(zfd(1) - ztdbt0)
2938 sfbmc(ibd) = sfbmc(ibd) + zb1
2939 sfdfc(ibd) = sfdfc(ibd) + zb2
2943 sfbmc(1) = sfbmc(1) + zf1
2944 sfdfc(1) = sfdfc(1) + zf2
2945 sfbmc(2) = sfbmc(2) + zf1
2946 sfdfc(2) = sfdfc(2) + zf2
2958 ftoadc = ftoadc + fxdn0(nlp1,ib)
2959 ftoau0 = ftoau0 + fxup0(nlp1,ib)
2960 fsfcu0 = fsfcu0 + fxup0(1,ib)
2961 fsfcd0 = fsfcd0 + fxdn0(1,ib)
2966 suvbf0 = fxdn0(1,ibd)
2968 if ( cf1 <= eps )
then 2971 fxupc(k,ib) = fxup0(k,ib)
2972 fxdnc(k,ib) = fxdn0(k,ib)
2991 fxupc(k,ib) = cf1*fxupc(k,ib) + cf0*fxup0(k,ib)
2992 fxdnc(k,ib) = cf1*fxdnc(k,ib) + cf0*fxdn0(k,ib)
2997 ftoauc = ftoauc + fxupc(nlp1,ib)
2998 fsfcuc = fsfcuc + fxupc(1,ib)
2999 fsfcdc = fsfcdc + fxdnc(1,ib)
3003 suvbfc = fxdnc(1,ibd)
3006 sfbmc(1) = cf1*sfbmc(1) + cf0*sfbm0(1)
3007 sfbmc(2) = cf1*sfbmc(2) + cf0*sfbm0(2)
3008 sfdfc(1) = cf1*sfdfc(1) + cf0*sfdf0(1)
3009 sfdfc(2) = cf1*sfdfc(2) + cf0*sfdf0(2)
3059 & ( ssolar,cosz,sntz,albbm,albdf,sfluxzen,cldfmc, &
3060 & cf1,cf0,taug,taur,tauae,ssaae,asyae,taucw,ssacw,asycw, &
3062 & fxupc,fxdnc,fxup0,fxdn0, &
3063 & ftoauc,ftoau0,ftoadc,fsfcuc,fsfcu0,fsfcdc,fsfcd0, &
3064 & sfbmc,sfdfc,sfbm0,sfdf0,suvbfc,suvbf0 &
3160 real (kind=kind_phys),
parameter :: zcrit = 0.9999995
3161 real (kind=kind_phys),
parameter :: zsr3 = sqrt(3.0)
3162 real (kind=kind_phys),
parameter :: od_lo = 0.06
3163 real (kind=kind_phys),
parameter :: eps1 = 1.0e-8
3166 integer,
intent(in) :: nlay, nlp1
3168 real (kind=kind_phys),
dimension(nlay,ngptsw),
intent(in) :: &
3169 & taug, taur, cldfmc
3170 real (kind=kind_phys),
dimension(nlay,nbdsw),
intent(in) :: &
3171 & taucw, ssacw, asycw, tauae, ssaae, asyae
3173 real (kind=kind_phys),
dimension(ngptsw),
intent(in) :: sfluxzen
3175 real (kind=kind_phys),
dimension(2),
intent(in) :: albbm, albdf
3177 real (kind=kind_phys),
intent(in) :: cosz, sntz, cf1, cf0, ssolar
3180 real (kind=kind_phys),
dimension(nlp1,nbdsw),
intent(out) :: &
3181 & fxupc, fxdnc, fxup0, fxdn0
3183 real (kind=kind_phys),
dimension(2),
intent(out) :: sfbmc, sfdfc, &
3186 real (kind=kind_phys),
intent(out) :: suvbfc, suvbf0, ftoadc, &
3187 & ftoauc, ftoau0, fsfcuc, fsfcu0, fsfcdc, fsfcd0
3190 real (kind=kind_phys),
dimension(nlay) :: ztaus, zssas, zasys, &
3193 real (kind=kind_phys),
dimension(nlp1) :: zrefb, zrefd, ztrab, &
3194 & ztrad, ztdbt, zldbt, zfu, zfd
3196 real (kind=kind_phys) :: ztau1, zssa1, zasy1, ztau0, zssa0, &
3197 & zasy0, zasy3, zssaw, zasyw, zgam1, zgam2, zgam3, zgam4, &
3198 & za1, za2, zb1, zb2, zrk, zrk2, zrp, zrp1, zrm1, zrpp, &
3199 & zrkg1, zrkg3, zrkg4, zexp1, zexm1, zexp2, zexm2, zden1, &
3200 & zexp3, zexp4, ze1r45, ftind, zsolar, ztdbt0, zr1, zr2, &
3201 & zr3, zr4, zr5, zt1, zt2, zt3, zf1, zf2, zrpp1
3203 integer :: ib, ibd, jb, jg, k, kp, itind
3211 fxdnc(k,ib) = f_zero
3212 fxupc(k,ib) = f_zero
3213 fxdn0(k,ib) = f_zero
3214 fxup0(k,ib) = f_zero
3242 lab_do_jg :
do jg = 1,
ngptsw 3248 zsolar = ssolar * sfluxzen(jg)
3257 zrefb(1) = albbm(ibd)
3258 zrefd(1) = albdf(ibd)
3260 zrefb(1) = 0.5 * (albbm(1) + albbm(2))
3261 zrefd(1) = 0.5 * (albdf(1) + albdf(2))
3279 ztau0 = max( ftiny, taur(k,jg)+taug(k,jg)+tauae(k,ib) )
3280 zssa0 = taur(k,jg) + tauae(k,ib)*ssaae(k,ib)
3281 zasy0 = asyae(k,ib)*ssaae(k,ib)*tauae(k,ib)
3282 zssaw = min( oneminus, zssa0 / ztau0 )
3283 zasyw = zasy0 / max( ftiny, zssa0 )
3294 ztau1 = (f_one - za2) * ztau0
3295 zssa1 = (zssaw - za2) / (f_one - za2)
3297 zasy1 = zasyw / (f_one + zasyw)
3298 zasy3 = 0.75 * zasy1
3301 if ( iswmode == 1 )
then 3302 zgam1 = 1.75 - zssa1 * (f_one + zasy3)
3303 zgam2 =-0.25 + zssa1 * (f_one - zasy3)
3304 zgam3 = 0.5 - zasy3 * cosz
3305 elseif ( iswmode == 2 )
then 3306 zgam1 = 2.0 - zssa1 * (1.25 + zasy3)
3307 zgam2 = 0.75* zssa1 * (f_one- zasy1)
3308 zgam3 = 0.5 - zasy3 * cosz
3309 elseif ( iswmode == 3 )
then 3310 zgam1 = zsr3 * (2.0 - zssa1 * (1.0 + zasy1)) * 0.5
3311 zgam2 = zsr3 * zssa1 * (1.0 - zasy1) * 0.5
3312 zgam3 = (1.0 - zsr3 * zasy1 * cosz) * 0.5
3314 zgam4 = f_one - zgam3
3318 if ( zssaw >= zcrit )
then 3319 za1 = zgam1 * cosz - zgam3
3325 zb1 = min( ztau1*sntz , 500.0 )
3326 if ( zb1 <= od_lo )
then 3327 zb2 = f_one - zb1 + 0.5*zb1*zb1
3329 ftind = zb1 / (
bpade + zb1)
3330 itind = ftind*
ntbmx + 0.5
3335 zrefb(kp) = max(f_zero, min(f_one, &
3336 & (za2 - za1*(f_one - zb2))/(f_one + za2) ))
3337 ztrab(kp) = max(f_zero, min(f_one, f_one-zrefb(kp) ))
3340 zrefd(kp) = max(f_zero, min(f_one, za2/(f_one + za2) ))
3341 ztrad(kp) = max(f_zero, min(f_one, f_one-zrefd(kp) ))
3344 za1 = zgam1*zgam4 + zgam2*zgam3
3345 za2 = zgam1*zgam3 + zgam2*zgam4
3346 zrk = sqrt( (zgam1 - zgam2) * (zgam1 + zgam2) )
3352 zrpp1= f_one - zrp*zrp
3353 zrpp = sign( max(flimit, abs(zrpp1)), zrpp1 )
3358 zr1 = zrm1 * (za2 + zrkg3)
3359 zr2 = zrp1 * (za2 - zrkg3)
3360 zr3 = zrk2 * (zgam3 - za2*cosz)
3362 zr5 = zrpp * (zrk - zgam1)
3364 zt1 = zrp1 * (za1 + zrkg4)
3365 zt2 = zrm1 * (za1 - zrkg4)
3366 zt3 = zrk2 * (zgam4 + za1*cosz)
3371 zb1 = min( zrk*ztau1, 500.0 )
3372 if ( zb1 <= od_lo )
then 3373 zexm1 = f_one - zb1 + 0.5*zb1*zb1
3375 ftind = zb1 / (
bpade + zb1)
3376 itind = ftind*
ntbmx + 0.5
3379 zexp1 = f_one / zexm1
3381 zb2 = min( sntz*ztau1, 500.0 )
3382 if ( zb2 <= od_lo )
then 3383 zexm2 = f_one - zb2 + 0.5*zb2*zb2
3385 ftind = zb2 / (
bpade + zb2)
3386 itind = ftind*
ntbmx + 0.5
3389 zexp2 = f_one / zexm2
3390 ze1r45 = zr4*zexp1 + zr5*zexm1
3393 if (ze1r45>=-eps1 .and. ze1r45<=eps1)
then 3397 zden1 = zssa1 / ze1r45
3398 zrefb(kp) = max(f_zero, min(f_one, &
3399 & (zr1*zexp1 - zr2*zexm1 - zr3*zexm2)*zden1 ))
3400 ztrab(kp) = max(f_zero, min(f_one, zexm2*(f_one &
3401 & - (zt1*zexp1 - zt2*zexm1 - zt3*zexp2)*zden1) ))
3405 zden1 = zr4 / (ze1r45 * zrkg1)
3406 zrefd(kp) = max(f_zero, min(f_one, &
3407 & zgam2*(zexp1 - zexm1)*zden1 ))
3408 ztrad(kp) = max(f_zero, min(f_one, zrk2*zden1 ))
3416 if ( zr1 <= od_lo )
then 3417 zexp3 = f_one - zr1 + 0.5*zr1*zr1
3419 ftind = zr1 / (
bpade + zr1)
3420 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
3424 ztdbt(k) = zexp3 * ztdbt(kp)
3431 if ( zr1 <= od_lo )
then 3432 zexp4 = f_one - zr1 + 0.5*zr1*zr1
3434 ftind = zr1 / (
bpade + zr1)
3435 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
3440 ztdbt0 = zexp4 * ztdbt0
3445 & ( zrefb,zrefd,ztrab,ztrad,zldbt,ztdbt, &
3453 fxup0(k,ib) = fxup0(k,ib) + zsolar*zfu(k)
3454 fxdn0(k,ib) = fxdn0(k,ib) + zsolar*zfd(k)
3459 zb2 = zsolar*(zfd(1) - ztdbt0)
3462 sfbm0(ibd) = sfbm0(ibd) + zb1
3463 sfdf0(ibd) = sfdf0(ibd) + zb2
3467 sfbm0(1) = sfbm0(1) + zf1
3468 sfdf0(1) = sfdf0(1) + zf2
3469 sfbm0(2) = sfbm0(2) + zf1
3470 sfdf0(2) = sfdf0(2) + zf2
3485 if ( cf1 > eps )
then 3493 if ( cldfmc(k,jg) > ftiny )
then 3495 ztau0 = ztaus(k) + taucw(k,ib)
3496 zssa0 = zssas(k) + ssacw(k,ib)
3497 zasy0 = zasys(k) + asycw(k,ib)
3498 zssaw = min(oneminus, zssa0 / ztau0)
3499 zasyw = zasy0 / max(ftiny, zssa0)
3505 ztau1 = (f_one - za2) * ztau0
3506 zssa1 = (zssaw - za2) / (f_one - za2)
3508 zasy1 = zasyw / (f_one + zasyw)
3509 zasy3 = 0.75 * zasy1
3512 if ( iswmode == 1 )
then 3513 zgam1 = 1.75 - zssa1 * (f_one + zasy3)
3514 zgam2 =-0.25 + zssa1 * (f_one - zasy3)
3515 zgam3 = 0.5 - zasy3 * cosz
3516 elseif ( iswmode == 2 )
then 3517 zgam1 = 2.0 - zssa1 * (1.25 + zasy3)
3518 zgam2 = 0.75* zssa1 * (f_one- zasy1)
3519 zgam3 = 0.5 - zasy3 * cosz
3520 elseif ( iswmode == 3 )
then 3521 zgam1 = zsr3 * (2.0 - zssa1 * (1.0 + zasy1)) * 0.5
3522 zgam2 = zsr3 * zssa1 * (1.0 - zasy1) * 0.5
3523 zgam3 = (1.0 - zsr3 * zasy1 * cosz) * 0.5
3525 zgam4 = f_one - zgam3
3529 if ( zssaw >= zcrit )
then 3530 za1 = zgam1 * cosz - zgam3
3536 zb1 = min( ztau1*sntz , 500.0 )
3537 if ( zb1 <= od_lo )
then 3538 zb2 = f_one - zb1 + 0.5*zb1*zb1
3540 ftind = zb1 / (
bpade + zb1)
3541 itind = ftind*
ntbmx + 0.5
3546 zrefb(kp) = max(f_zero, min(f_one, &
3547 & (za2 - za1*(f_one - zb2))/(f_one + za2) ))
3548 ztrab(kp) = max(f_zero, min(f_one, f_one-zrefb(kp)))
3551 zrefd(kp) = max(f_zero, min(f_one, za2 / (f_one+za2) ))
3552 ztrad(kp) = max(f_zero, min(f_one, f_one - zrefd(kp) ))
3555 za1 = zgam1*zgam4 + zgam2*zgam3
3556 za2 = zgam1*zgam3 + zgam2*zgam4
3557 zrk = sqrt( (zgam1 - zgam2) * (zgam1 + zgam2) )
3563 zrpp1= f_one - zrp*zrp
3564 zrpp = sign( max(flimit, abs(zrpp1)), zrpp1 )
3569 zr1 = zrm1 * (za2 + zrkg3)
3570 zr2 = zrp1 * (za2 - zrkg3)
3571 zr3 = zrk2 * (zgam3 - za2*cosz)
3573 zr5 = zrpp * (zrk - zgam1)
3575 zt1 = zrp1 * (za1 + zrkg4)
3576 zt2 = zrm1 * (za1 - zrkg4)
3577 zt3 = zrk2 * (zgam4 + za1*cosz)
3582 zb1 = min( zrk*ztau1, 500.0 )
3583 if ( zb1 <= od_lo )
then 3584 zexm1 = f_one - zb1 + 0.5*zb1*zb1
3586 ftind = zb1 / (
bpade + zb1)
3587 itind = ftind*
ntbmx + 0.5
3590 zexp1 = f_one / zexm1
3592 zb2 = min( ztau1*sntz, 500.0 )
3593 if ( zb2 <= od_lo )
then 3594 zexm2 = f_one - zb2 + 0.5*zb2*zb2
3596 ftind = zb2 / (
bpade + zb2)
3597 itind = ftind*
ntbmx + 0.5
3600 zexp2 = f_one / zexm2
3601 ze1r45 = zr4*zexp1 + zr5*zexm1
3604 if ( ze1r45>=-eps1 .and. ze1r45<=eps1 )
then 3608 zden1 = zssa1 / ze1r45
3609 zrefb(kp) = max(f_zero, min(f_one, &
3610 & (zr1*zexp1-zr2*zexm1-zr3*zexm2)*zden1 ))
3611 ztrab(kp) = max(f_zero, min(f_one, zexm2*(f_one - &
3612 & (zt1*zexp1-zt2*zexm1-zt3*zexp2)*zden1) ))
3616 zden1 = zr4 / (ze1r45 * zrkg1)
3617 zrefd(kp) = max(f_zero, min(f_one, &
3618 & zgam2*(zexp1 - zexm1)*zden1 ))
3619 ztrad(kp) = max(f_zero, min(f_one, zrk2*zden1 ))
3627 if ( zr1 <= od_lo )
then 3628 zexp3 = f_one - zr1 + 0.5*zr1*zr1
3630 ftind = zr1 / (
bpade + zr1)
3631 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
3636 ztdbt(k) = zexp3 * ztdbt(kp)
3642 if ( zr1 <= od_lo )
then 3643 zexp4 = f_one - zr1 + 0.5*zr1*zr1
3645 ftind = zr1 / (
bpade + zr1)
3646 itind = max(0, min(
ntbmx, int(0.5+
ntbmx*ftind) ))
3650 ztdbt0 = zexp4 * ztdbt0
3655 ztdbt(k) = zldbt(kp) * ztdbt(kp)
3658 ztdbt0 = zldbt0(k) * ztdbt0
3667 & ( zrefb,zrefd,ztrab,ztrad,zldbt,ztdbt, &
3675 fxupc(k,ib) = fxupc(k,ib) + zsolar*zfu(k)
3676 fxdnc(k,ib) = fxdnc(k,ib) + zsolar*zfd(k)
3682 zb2 = zsolar*(zfd(1) - ztdbt0)
3685 sfbmc(ibd) = sfbmc(ibd) + zb1
3686 sfdfc(ibd) = sfdfc(ibd) + zb2
3690 sfbmc(1) = sfbmc(1) + zf1
3691 sfdfc(1) = sfdfc(1) + zf2
3692 sfbmc(2) = sfbmc(2) + zf1
3693 sfdfc(2) = sfdfc(2) + zf2
3705 ftoadc = ftoadc + fxdn0(nlp1,ib)
3706 ftoau0 = ftoau0 + fxup0(nlp1,ib)
3707 fsfcu0 = fsfcu0 + fxup0(1,ib)
3708 fsfcd0 = fsfcd0 + fxdn0(1,ib)
3713 suvbf0 = fxdn0(1,ibd)
3715 if ( cf1 <= eps )
then 3718 fxupc(k,ib) = fxup0(k,ib)
3719 fxdnc(k,ib) = fxdn0(k,ib)
3737 ftoauc = ftoauc + fxupc(nlp1,ib)
3738 fsfcuc = fsfcuc + fxupc(1,ib)
3739 fsfcdc = fsfcdc + fxdnc(1,ib)
3743 suvbfc = fxdnc(1,ibd)
3766 & ( zrefb,zrefd,ztrab,ztrad,zldbt,ztdbt, &
3798 integer,
intent(in) :: nlay, nlp1
3800 real (kind=kind_phys),
dimension(nlp1),
intent(in) :: zrefb, &
3801 & zrefd, ztrab, ztrad, ztdbt, zldbt
3804 real (kind=kind_phys),
dimension(nlp1),
intent(out) :: zfu, zfd
3807 real (kind=kind_phys),
dimension(nlp1) :: zrupb,zrupd,zrdnd,ztdn
3809 real (kind=kind_phys) :: zden1
3824 zden1 = f_one / ( f_one - zrupd(k)*zrefd(kp) )
3825 zrupb(kp) = zrefb(kp) + ( ztrad(kp) * &
3826 & ( (ztrab(kp) - zldbt(kp))*zrupd(k) + &
3827 & zldbt(kp)*zrupb(k)) ) * zden1
3828 zrupd(kp) = zrefd(kp) + ztrad(kp)*ztrad(kp)*zrupd(k)*zden1
3833 zrdnd(nlp1) = f_zero
3834 ztdn(nlay) = ztrab(nlp1)
3835 zrdnd(nlay) = zrefd(nlp1)
3839 zden1 = f_one / (f_one - zrefd(k)*zrdnd(k))
3840 ztdn(k-1) = ztdbt(k)*ztrab(k) + ( ztrad(k) * &
3841 & ( (ztdn(k) - ztdbt(k)) + ztdbt(k) * &
3842 & zrefb(k)*zrdnd(k) )) * zden1
3843 zrdnd(k-1) = zrefd(k) + ztrad(k)*ztrad(k)*zrdnd(k)*zden1
3848 zden1 = f_one / (f_one - zrdnd(k)*zrupd(k))
3849 zfu(k) = ( ztdbt(k)*zrupb(k) + &
3850 & (ztdn(k) - ztdbt(k))*zrupd(k) ) * zden1
3851 zfd(k) = ztdbt(k) + ( ztdn(k) - ztdbt(k) + &
3852 & ztdbt(k)*zrupb(k)*zrdnd(k) ) * zden1
3902 & ( colamt,colmol,fac00,fac01,fac10,fac11,jp,jt,jt1,laytrop, &
3903 & forfac,forfrac,indfor,selffac,selffrac,indself, nlay, &
3904 & sfluxzen, taug, taur &
4004 integer,
intent(in) :: nlay, laytrop
4006 integer,
dimension(nlay),
intent(in) :: indfor, indself, &
4009 real (kind=kind_phys),
dimension(nlay),
intent(in) :: colmol, &
4010 & fac00, fac01, fac10, fac11, forfac, forfrac, selffac, &
4013 real (kind=kind_phys),
dimension(nlay,maxgas),
intent(in) :: colamt
4016 real (kind=kind_phys),
dimension(ngptsw),
intent(out) :: sfluxzen
4018 real (kind=kind_phys),
dimension(nlay,ngptsw),
intent(out) :: &
4022 real (kind=kind_phys) :: fs, speccomb, specmult, colm1, colm2
4024 integer,
dimension(nlay,nblow:nbhgh) :: id0, id1
4026 integer :: ibd, j, jb, js, k, klow, khgh, klim, ks, njb, ns
4037 id0(k,jb) = ((jp(k)-1)*5 + (jt(k)-1)) * nspa(jb)
4038 id1(k,jb) = ( jp(k) *5 + (jt1(k)-1)) * nspa(jb)
4041 do k = laytrop+1, nlay
4042 id0(k,jb) = ((jp(k)-13)*5 + (jt(k)-1)) * nspb(jb)
4043 id1(k,jb) = ((jp(k)-12)*5 + (jt1(k)-1)) * nspb(jb)
4054 case (16, 20, 23, 25, 26, 29)
4063 sfluxzen(ns+j) = scalekur *
sfluxref01(j,1,ibd)
4068 if (jb==17 .or. jb==28)
then 4071 lab_do_k1 :
do k = laytrop, nlay-1
4078 colm1 = colamt(ks,
ix1(jb))
4079 colm2 = colamt(ks,
ix2(jb))
4080 speccomb = colm1 +
strrat(jb)*colm2
4081 specmult =
specwt(jb) * min( oneminus, colm1/speccomb )
4082 js = 1 + int( specmult )
4083 fs = mod(specmult, f_one)
4093 lab_do_k2 :
do k = 1, laytrop-1
4100 colm1 = colamt(ks,
ix1(jb))
4101 colm2 = colamt(ks,
ix2(jb))
4102 speccomb = colm1 +
strrat(jb)*colm2
4103 specmult =
specwt(jb) * min( oneminus, colm1/speccomb )
4104 js = 1 + int( specmult )
4105 fs = mod(specmult, f_one)
4168 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
4169 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
4171 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4172 integer :: inds, indf, indsp, indfp, j, js, k
4183 tauray = colmol(k) *
rayl 4186 taur(k,ns16+j) = tauray
4191 speccomb = colamt(k,1) +
strrat(16)*colamt(k,5)
4192 specmult = 8.0 * min( oneminus, colamt(k,1)/speccomb )
4194 js = 1 + int( specmult )
4195 fs = mod( specmult, f_one )
4197 fac000 = fs1 * fac00(k)
4198 fac010 = fs1 * fac10(k)
4199 fac100 = fs * fac00(k)
4200 fac110 = fs * fac10(k)
4201 fac001 = fs1 * fac01(k)
4202 fac011 = fs1 * fac11(k)
4203 fac101 = fs * fac01(k)
4204 fac111 = fs * fac11(k)
4206 ind01 = id0(k,16) + js
4210 ind11 = id1(k,16) + js
4220 taug(k,ns16+j) = speccomb &
4221 & *( fac000 *
absa(ind01,j) + fac100 *
absa(ind02,j) &
4222 & + fac010 *
absa(ind03,j) + fac110 *
absa(ind04,j) &
4223 & + fac001 *
absa(ind11,j) + fac101 *
absa(ind12,j) &
4224 & + fac011 *
absa(ind13,j) + fac111 *
absa(ind14,j) ) &
4225 & + colamt(k,1) * (selffac(k) * (
selfref(inds,j) &
4227 & + forfac(k) * (forref(indf,j) + forfrac(k) &
4228 & * (forref(indfp,j) - forref(indf,j))))
4232 do k = laytrop+1, nlay
4233 ind01 = id0(k,16) + 1
4235 ind11 = id1(k,16) + 1
4239 taug(k,ns16+j) = colamt(k,5) &
4240 & * ( fac00(k)*
absb(ind01,j) + fac10(k)*
absb(ind02,j) &
4241 & + fac01(k)*
absb(ind11,j) + fac11(k)*
absb(ind12,j) )
4264 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
4265 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
4267 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4268 integer :: inds, indf, indsp, indfp, j, js, k
4279 tauray = colmol(k) *
rayl 4282 taur(k,ns17+j) = tauray
4287 speccomb = colamt(k,1) +
strrat(17)*colamt(k,2)
4288 specmult = 8.0 * min(oneminus, colamt(k,1) / speccomb)
4290 js = 1 + int(specmult)
4291 fs = mod(specmult, f_one)
4293 fac000 = fs1 * fac00(k)
4294 fac010 = fs1 * fac10(k)
4295 fac100 = fs * fac00(k)
4296 fac110 = fs * fac10(k)
4297 fac001 = fs1 * fac01(k)
4298 fac011 = fs1 * fac11(k)
4299 fac101 = fs * fac01(k)
4300 fac111 = fs * fac11(k)
4302 ind01 = id0(k,17) + js
4306 ind11 = id1(k,17) + js
4317 taug(k,ns17+j) = speccomb &
4318 & * ( fac000 *
absa(ind01,j) + fac100 *
absa(ind02,j) &
4319 & + fac010 *
absa(ind03,j) + fac110 *
absa(ind04,j) &
4320 & + fac001 *
absa(ind11,j) + fac101 *
absa(ind12,j) &
4321 & + fac011 *
absa(ind13,j) + fac111 *
absa(ind14,j) ) &
4322 & + colamt(k,1) * (selffac(k) * (
selfref(inds,j) &
4324 & + forfac(k) * (forref(indf,j) + forfrac(k) &
4325 & * (forref(indfp,j) - forref(indf,j))))
4329 do k = laytrop+1, nlay
4330 speccomb = colamt(k,1) +
strrat(17)*colamt(k,2)
4331 specmult = 4.0 * min(oneminus, colamt(k,1) / speccomb)
4333 js = 1 + int(specmult)
4334 fs = mod(specmult, f_one)
4336 fac000 = fs1 * fac00(k)
4337 fac010 = fs1 * fac10(k)
4338 fac100 = fs * fac00(k)
4339 fac110 = fs * fac10(k)
4340 fac001 = fs1 * fac01(k)
4341 fac011 = fs1 * fac11(k)
4342 fac101 = fs * fac01(k)
4343 fac111 = fs * fac11(k)
4345 ind01 = id0(k,17) + js
4349 ind11 = id1(k,17) + js
4358 taug(k,ns17+j) = speccomb &
4359 & * ( fac000 *
absb(ind01,j) + fac100 *
absb(ind02,j) &
4360 & + fac010 *
absb(ind03,j) + fac110 *
absb(ind04,j) &
4361 & + fac001 *
absb(ind11,j) + fac101 *
absb(ind12,j) &
4362 & + fac011 *
absb(ind13,j) + fac111 *
absb(ind14,j) ) &
4363 & + colamt(k,1) * forfac(k) * (forref(indf,j) &
4364 & + forfrac(k) * (forref(indfp,j) - forref(indf,j)))
4387 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
4388 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
4390 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4391 integer :: inds, indf, indsp, indfp, j, js, k
4402 tauray = colmol(k) *
rayl 4405 taur(k,ns18+j) = tauray
4410 speccomb = colamt(k,1) +
strrat(18)*colamt(k,5)
4411 specmult = 8.0 * min(oneminus, colamt(k,1) / speccomb)
4413 js = 1 + int(specmult)
4414 fs = mod(specmult, f_one)
4416 fac000 = fs1 * fac00(k)
4417 fac010 = fs1 * fac10(k)
4418 fac100 = fs * fac00(k)
4419 fac110 = fs * fac10(k)
4420 fac001 = fs1 * fac01(k)
4421 fac011 = fs1 * fac11(k)
4422 fac101 = fs * fac01(k)
4423 fac111 = fs * fac11(k)
4425 ind01 = id0(k,18) + js
4429 ind11 = id1(k,18) + js
4440 taug(k,ns18+j) = speccomb &
4441 & * ( fac000 *
absa(ind01,j) + fac100 *
absa(ind02,j) &
4442 & + fac010 *
absa(ind03,j) + fac110 *
absa(ind04,j) &
4443 & + fac001 *
absa(ind11,j) + fac101 *
absa(ind12,j) &
4444 & + fac011 *
absa(ind13,j) + fac111 *
absa(ind14,j) ) &
4445 & + colamt(k,1) * (selffac(k) * (
selfref(inds,j) &
4447 & + forfac(k) * (forref(indf,j) + forfrac(k) &
4448 & * (forref(indfp,j) - forref(indf,j))))
4452 do k = laytrop+1, nlay
4453 ind01 = id0(k,18) + 1
4455 ind11 = id1(k,18) + 1
4459 taug(k,ns18+j) = colamt(k,5) &
4460 & * ( fac00(k)*
absb(ind01,j) + fac10(k)*
absb(ind02,j) &
4461 & + fac01(k)*
absb(ind11,j) + fac11(k)*
absb(ind12,j) )
4483 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
4484 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
4486 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4487 integer :: inds, indf, indsp, indfp, j, js, k
4498 tauray = colmol(k) *
rayl 4501 taur(k,ns19+j) = tauray
4506 speccomb = colamt(k,1) +
strrat(19)*colamt(k,2)
4507 specmult = 8.0 * min(oneminus, colamt(k,1) / speccomb)
4509 js = 1 + int(specmult)
4510 fs = mod(specmult, f_one)
4512 fac000 = fs1 * fac00(k)
4513 fac010 = fs1 * fac10(k)
4514 fac100 = fs * fac00(k)
4515 fac110 = fs * fac10(k)
4516 fac001 = fs1 * fac01(k)
4517 fac011 = fs1 * fac11(k)
4518 fac101 = fs * fac01(k)
4519 fac111 = fs * fac11(k)
4521 ind01 = id0(k,19) + js
4525 ind11 = id1(k,19) + js
4536 taug(k,ns19+j) = speccomb &
4537 & * ( fac000 *
absa(ind01,j) + fac100 *
absa(ind02,j) &
4538 & + fac010 *
absa(ind03,j) + fac110 *
absa(ind04,j) &
4539 & + fac001 *
absa(ind11,j) + fac101 *
absa(ind12,j) &
4540 & + fac011 *
absa(ind13,j) + fac111 *
absa(ind14,j) ) &
4541 & + colamt(k,1) * (selffac(k) * (
selfref(inds,j) &
4543 & + forfac(k) * (forref(indf,j) + forfrac(k) &
4544 & * (forref(indfp,j) - forref(indf,j))))
4548 do k = laytrop+1, nlay
4549 ind01 = id0(k,19) + 1
4551 ind11 = id1(k,19) + 1
4555 taug(k,ns19+j) = colamt(k,2) &
4556 & * ( fac00(k)*
absb(ind01,j) + fac10(k)*
absb(ind02,j) &
4557 & + fac01(k)*
absb(ind11,j) + fac11(k)*
absb(ind12,j) )
4579 real (kind=kind_phys) :: tauray
4581 integer :: ind01, ind02, ind11, ind12
4582 integer :: inds, indf, indsp, indfp, j, k
4593 tauray = colmol(k) *
rayl 4596 taur(k,ns20+j) = tauray
4601 ind01 = id0(k,20) + 1
4603 ind11 = id1(k,20) + 1
4612 taug(k,ns20+j) = colamt(k,1) &
4613 & * ( (fac00(k)*
absa(ind01,j) + fac10(k)*
absa(ind02,j) &
4614 & + fac01(k)*
absa(ind11,j) + fac11(k)*
absa(ind12,j)) &
4615 & + selffac(k) * (
selfref(inds,j) + selffrac(k) &
4617 & + forfac(k) * (forref(indf,j) + forfrac(k) &
4618 & * (forref(indfp,j) - forref(indf,j))) ) &
4619 & + colamt(k,5) *
absch4(j)
4623 do k = laytrop+1, nlay
4624 ind01 = id0(k,20) + 1
4626 ind11 = id1(k,20) + 1
4633 taug(k,ns20+j) = colamt(k,1) &
4634 & * ( fac00(k)*
absb(ind01,j) + fac10(k)*
absb(ind02,j) &
4635 & + fac01(k)*
absb(ind11,j) + fac11(k)*
absb(ind12,j) &
4636 & + forfac(k) * (forref(indf,j) + forfrac(k) &
4637 & * (forref(indfp,j) - forref(indf,j))) ) &
4638 & + colamt(k,5) *
absch4(j)
4661 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
4662 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
4664 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4665 integer :: inds, indf, indsp, indfp, j, js, k
4676 tauray = colmol(k) *
rayl 4679 taur(k,ns21+j) = tauray
4684 speccomb = colamt(k,1) +
strrat(21)*colamt(k,2)
4685 specmult = 8.0 * min(oneminus, colamt(k,1) / speccomb)
4687 js = 1 + int(specmult)
4688 fs = mod(specmult, f_one)
4690 fac000 = fs1 * fac00(k)
4691 fac010 = fs1 * fac10(k)
4692 fac100 = fs * fac00(k)
4693 fac110 = fs * fac10(k)
4694 fac001 = fs1 * fac01(k)
4695 fac011 = fs1 * fac11(k)
4696 fac101 = fs * fac01(k)
4697 fac111 = fs * fac11(k)
4699 ind01 = id0(k,21) + js
4703 ind11 = id1(k,21) + js
4714 taug(k,ns21+j) = speccomb &
4715 & * ( fac000 *
absa(ind01,j) + fac100 *
absa(ind02,j) &
4716 & + fac010 *
absa(ind03,j) + fac110 *
absa(ind04,j) &
4717 & + fac001 *
absa(ind11,j) + fac101 *
absa(ind12,j) &
4718 & + fac011 *
absa(ind13,j) + fac111 *
absa(ind14,j) ) &
4719 & + colamt(k,1) * (selffac(k) * (
selfref(inds,j) &
4721 & + forfac(k) * (forref(indf,j) + forfrac(k) &
4722 & * (forref(indfp,j) - forref(indf,j))))
4726 do k = laytrop+1, nlay
4727 speccomb = colamt(k,1) +
strrat(21)*colamt(k,2)
4728 specmult = 4.0 * min(oneminus, colamt(k,1) / speccomb)
4730 js = 1 + int(specmult)
4731 fs = mod(specmult, f_one)
4733 fac000 = fs1 * fac00(k)
4734 fac010 = fs1 * fac10(k)
4735 fac100 = fs * fac00(k)
4736 fac110 = fs * fac10(k)
4737 fac001 = fs1 * fac01(k)
4738 fac011 = fs1 * fac11(k)
4739 fac101 = fs * fac01(k)
4740 fac111 = fs * fac11(k)
4742 ind01 = id0(k,21) + js
4746 ind11 = id1(k,21) + js
4755 taug(k,ns21+j) = speccomb &
4756 & * ( fac000 *
absb(ind01,j) + fac100 *
absb(ind02,j) &
4757 & + fac010 *
absb(ind03,j) + fac110 *
absb(ind04,j) &
4758 & + fac001 *
absb(ind11,j) + fac101 *
absb(ind12,j) &
4759 & + fac011 *
absb(ind13,j) + fac111 *
absb(ind14,j) ) &
4760 & + colamt(k,1) * forfac(k) * (forref(indf,j) &
4761 & + forfrac(k) * (forref(indfp,j) - forref(indf,j)))
4783 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
4784 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111, &
4785 & o2adj, o2cont, o2tem
4787 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4788 integer :: inds, indf, indsp, indfp, j, js, k
4798 o2tem = 4.35e-4 / (350.0*2.0)
4806 tauray = colmol(k) *
rayl 4809 taur(k,ns22+j) = tauray
4814 o2cont = o2tem * colamt(k,6)
4815 speccomb = colamt(k,1) +
strrat(22)*colamt(k,6)
4816 specmult = 8.0 * min(oneminus, colamt(k,1) / speccomb)
4818 js = 1 + int(specmult)
4819 fs = mod(specmult, f_one)
4821 fac000 = fs1 * fac00(k)
4822 fac010 = fs1 * fac10(k)
4823 fac100 = fs * fac00(k)
4824 fac110 = fs * fac10(k)
4825 fac001 = fs1 * fac01(k)
4826 fac011 = fs1 * fac11(k)
4827 fac101 = fs * fac01(k)
4828 fac111 = fs * fac11(k)
4830 ind01 = id0(k,22) + js
4834 ind11 = id1(k,22) + js
4845 taug(k,ns22+j) = speccomb &
4846 & * ( fac000 *
absa(ind01,j) + fac100 *
absa(ind02,j) &
4847 & + fac010 *
absa(ind03,j) + fac110 *
absa(ind04,j) &
4848 & + fac001 *
absa(ind11,j) + fac101 *
absa(ind12,j) &
4849 & + fac011 *
absa(ind13,j) + fac111 *
absa(ind14,j) ) &
4850 & + colamt(k,1) * (selffac(k) * (
selfref(inds,j) &
4852 & + forfac(k) * (forref(indf,j) + forfrac(k) &
4853 & * (forref(indfp,j) - forref(indf,j)))) + o2cont
4857 do k = laytrop+1, nlay
4858 o2cont = o2tem * colamt(k,6)
4860 ind01 = id0(k,22) + 1
4862 ind11 = id1(k,22) + 1
4866 taug(k,ns22+j) = colamt(k,6) * o2adj &
4867 & * ( fac00(k)*
absb(ind01,j) + fac10(k)*
absb(ind02,j) &
4868 & + fac01(k)*
absb(ind11,j) + fac11(k)*
absb(ind12,j) ) &
4892 integer :: ind01, ind02, ind11, ind12
4893 integer :: inds, indf, indsp, indfp, j, k
4905 taur(k,ns23+j) = colmol(k) *
rayl(j)
4910 ind01 = id0(k,23) + 1
4912 ind11 = id1(k,23) + 1
4921 taug(k,ns23+j) = colamt(k,1) * (
givfac &
4922 & * ( fac00(k)*
absa(ind01,j) + fac10(k)*
absa(ind02,j) &
4923 & + fac01(k)*
absa(ind11,j) + fac11(k)*
absa(ind12,j) ) &
4924 & + selffac(k) * (
selfref(inds,j) + selffrac(k) &
4926 & + forfac(k) * (forref(indf,j) + forfrac(k) &
4927 & * (forref(indfp,j) - forref(indf,j))))
4931 do k = laytrop+1, nlay
4933 taug(k,ns23+j) = f_zero
4955 real (kind=kind_phys) :: speccomb, specmult, fs, fs1, &
4956 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
4958 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
4959 integer :: inds, indf, indsp, indfp, j, js, k
4970 speccomb = colamt(k,1) +
strrat(24)*colamt(k,6)
4971 specmult = 8.0 * min(oneminus, colamt(k,1) / speccomb)
4973 js = 1 + int(specmult)
4974 fs = mod(specmult, f_one)
4976 fac000 = fs1 * fac00(k)
4977 fac010 = fs1 * fac10(k)
4978 fac100 = fs * fac00(k)
4979 fac110 = fs * fac10(k)
4980 fac001 = fs1 * fac01(k)
4981 fac011 = fs1 * fac11(k)
4982 fac101 = fs * fac01(k)
4983 fac111 = fs * fac11(k)
4985 ind01 = id0(k,24) + js
4989 ind11 = id1(k,24) + js
5000 taug(k,ns24+j) = speccomb &
5001 & * ( fac000 *
absa(ind01,j) + fac100 *
absa(ind02,j) &
5002 & + fac010 *
absa(ind03,j) + fac110 *
absa(ind04,j) &
5003 & + fac001 *
absa(ind11,j) + fac101 *
absa(ind12,j) &
5004 & + fac011 *
absa(ind13,j) + fac111 *
absa(ind14,j) ) &
5005 & + colamt(k,3) *
abso3a(j) + colamt(k,1) &
5006 & * (selffac(k) * (
selfref(inds,j) + selffrac(k) &
5008 & + forfac(k) * (forref(indf,j) + forfrac(k) &
5009 & * (forref(indfp,j) - forref(indf,j))))
5011 taur(k,ns24+j) = colmol(k) &
5016 do k = laytrop+1, nlay
5017 ind01 = id0(k,24) + 1
5019 ind11 = id1(k,24) + 1
5023 taug(k,ns24+j) = colamt(k,6) &
5024 & * ( fac00(k)*
absb(ind01,j) + fac10(k)*
absb(ind02,j) &
5025 & + fac01(k)*
absb(ind11,j) + fac11(k)*
absb(ind12,j) ) &
5026 & + colamt(k,3) *
abso3b(j)
5028 taur(k,ns24+j) = colmol(k) * raylb(j)
5051 integer :: ind01, ind02, ind11, ind12
5064 taur(k,ns25+j) = colmol(k) *
rayl(j)
5069 ind01 = id0(k,25) + 1
5071 ind11 = id1(k,25) + 1
5075 taug(k,ns25+j) = colamt(k,1) &
5076 & * ( fac00(k)*
absa(ind01,j) + fac10(k)*
absa(ind02,j) &
5077 & + fac01(k)*
absa(ind11,j) + fac11(k)*
absa(ind12,j) ) &
5078 & + colamt(k,3) *
abso3a(j)
5082 do k = laytrop+1, nlay
5084 taug(k,ns25+j) = colamt(k,3) *
abso3b(j)
5119 taug(k,ns26+j) = f_zero
5120 taur(k,ns26+j) = colmol(k) *
rayl(j)
5143 integer :: ind01, ind02, ind11, ind12
5156 taur(k,ns27+j) = colmol(k) *
rayl(j)
5161 ind01 = id0(k,27) + 1
5163 ind11 = id1(k,27) + 1
5167 taug(k,ns27+j) = colamt(k,3) &
5168 & * ( fac00(k)*
absa(ind01,j) + fac10(k)*
absa(ind02,j) &
5169 & + fac01(k)*
absa(ind11,j) + fac11(k)*
absa(ind12,j) )
5173 do k = laytrop+1, nlay
5174 ind01 = id0(k,27) + 1
5176 ind11 = id1(k,27) + 1
5180 taug(k,ns27+j) = colamt(k,3) &
5181 & * ( fac00(k)*
absb(ind01,j) + fac10(k)*
absb(ind02,j) &
5182 & + fac01(k)*
absb(ind11,j) + fac11(k)*
absb(ind12,j) )
5205 real (kind=kind_phys) :: speccomb, specmult, tauray, fs, fs1, &
5206 & fac000,fac001,fac010,fac011, fac100,fac101,fac110,fac111
5208 integer :: ind01, ind02, ind03, ind04, ind11, ind12, ind13, ind14
5220 tauray = colmol(k) *
rayl 5223 taur(k,ns28+j) = tauray
5228 speccomb = colamt(k,3) +
strrat(28)*colamt(k,6)
5229 specmult = 8.0 * min(oneminus, colamt(k,3) / speccomb)
5231 js = 1 + int(specmult)
5232 fs = mod(specmult, f_one)
5234 fac000 = fs1 * fac00(k)
5235 fac010 = fs1 * fac10(k)
5236 fac100 = fs * fac00(k)
5237 fac110 = fs * fac10(k)
5238 fac001 = fs1 * fac01(k)
5239 fac011 = fs1 * fac11(k)
5240 fac101 = fs * fac01(k)
5241 fac111 = fs * fac11(k)
5243 ind01 = id0(k,28) + js
5247 ind11 = id1(k,28) + js
5253 taug(k,ns28+j) = speccomb &
5254 & * ( fac000 *
absa(ind01,j) + fac100 *
absa(ind02,j) &
5255 & + fac010 *
absa(ind03,j) + fac110 *
absa(ind04,j) &
5256 & + fac001 *
absa(ind11,j) + fac101 *
absa(ind12,j) &
5257 & + fac011 *
absa(ind13,j) + fac111 *
absa(ind14,j) )
5261 do k = laytrop+1, nlay
5262 speccomb = colamt(k,3) +
strrat(28)*colamt(k,6)
5263 specmult = 4.0 * min(oneminus, colamt(k,3) / speccomb)
5265 js = 1 + int(specmult)
5266 fs = mod(specmult, f_one)
5268 fac000 = fs1 * fac00(k)
5269 fac010 = fs1 * fac10(k)
5270 fac100 = fs * fac00(k)
5271 fac110 = fs * fac10(k)
5272 fac001 = fs1 * fac01(k)
5273 fac011 = fs1 * fac11(k)
5274 fac101 = fs * fac01(k)
5275 fac111 = fs * fac11(k)
5277 ind01 = id0(k,28) + js
5281 ind11 = id1(k,28) + js
5287 taug(k,ns28+j) = speccomb &
5288 & * ( fac000 *
absb(ind01,j) + fac100 *
absb(ind02,j) &
5289 & + fac010 *
absb(ind03,j) + fac110 *
absb(ind04,j) &
5290 & + fac001 *
absb(ind11,j) + fac101 *
absb(ind12,j) &
5291 & + fac011 *
absb(ind13,j) + fac111 *
absb(ind14,j) )
5314 real (kind=kind_phys) :: tauray
5316 integer :: ind01, ind02, ind11, ind12
5317 integer :: inds, indf, indsp, indfp, j, k
5328 tauray = colmol(k) *
rayl 5331 taur(k,ns29+j) = tauray
5336 ind01 = id0(k,29) + 1
5338 ind11 = id1(k,29) + 1
5347 taug(k,ns29+j) = colamt(k,1) &
5348 & * ( (fac00(k)*
absa(ind01,j) + fac10(k)*
absa(ind02,j) &
5349 & + fac01(k)*
absa(ind11,j) + fac11(k)*
absa(ind12,j) ) &
5350 & + selffac(k) * (
selfref(inds,j) + selffrac(k) &
5352 & + forfac(k) * (forref(indf,j) + forfrac(k) &
5353 & * (forref(indfp,j) - forref(indf,j)))) &
5354 & + colamt(k,2) *
absco2(j)
5358 do k = laytrop+1, nlay
5359 ind01 = id0(k,29) + 1
5361 ind11 = id1(k,29) + 1
5365 taug(k,ns29+j) = colamt(k,2) &
5366 & * ( fac00(k)*
absb(ind01,j) + fac10(k)*
absb(ind02,j) &
5367 & + fac01(k)*
absb(ind11,j) + fac11(k)*
absb(ind12,j) ) &
5368 & + colamt(k,1) *
absh2o(j)
5384 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...
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) ...
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(5), public ebari
asymmetry coefficients
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...
integer, dimension(nblow:nbhgh), public ix2
indexes for 2nd entries of the two key species for each of the spectral bands
real(kind=kind_phys), dimension(0:ntbmx) exp_tbl
those data will be set up only once by "rswinit"
integer, dimension(nblow:nbhgh), public ix1
indexes for 1st entries of the two key species for each of the spectral bands
derived type for special components of surface SW fluxes
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...
real(kind=kind_phys), dimension(nblow:nbhgh), public c0r
asymmetry coefficients
subroutine vrtqdr(zrefb, zrefd, ztrab, ztrad, zldbt, ztdbt, NLAY, NLP1, zfu, zfd)
This subroutine is called by spcvrtc() and spcvrtm(), and computes the upward and downward radiation ...
This module sets up absorption coeffients for band 27: 29000-38000 cm-1 (low - o3; high - o3) ...
real(kind=kind_phys), dimension(43, nblow:nbhgh), public asyice2
asymmetry coefficients
real(kind=kind_phys), dimension(43, nblow:nbhgh), public extice2
extinction coefficients
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 index
This module sets up absorption coeffients for band 20: 5150-6150 cm-1 (low - h2o; high - h2o) ...
integer, parameter ntbmx
index upper limit of optical depth and transmittance tables
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...
real(kind=kind_phys), public a0s
optical depth coefficients
real(kind=kind_phys), dimension(ngmax, mfs01, mfb01), target, public sfluxref01
spectral solar fluxes, j=1,2,...,7 for SW band number of 16,20,23,25,26,27,29
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
real(kind=kind_phys), dimension(58, nblow:nbhgh), public extliq1
extinction coefficients
This module contains cloud radiative property coefficients.
real(kind=kind_phys), dimension(nblow:nbhgh), public c0s
asymmetry coefficients
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...
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 ...
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(ngmax, mfs02, mfb02), target, public sfluxref02
spectral solar fluxes, j=1,2 for SW band number of 17 and 28
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 ...
subroutine, public swrad(plyr, plvl, tlyr, tlvl, qlyr, olyr, gasvmr, clouds, icseed, aerosols, sfcalb, cosz, solcon, NDAY, idxday, npts, nlay, nlp1, lprnt, hswc, topflx, sfcflx, HSW0, HSWB, FLXPRF, FDNCMP)
This subroutine is the main SW radiation routine.
subroutine spcvrtc(ssolar, cosz, sntz, albbm, albdf, sfluxzen, cldfrc, cf1, cf0, taug, taur, tauae, ssaae, asyae, taucw, ssacw, asycw, nlay, nlp1, fxupc, fxdnc, fxup0, fxdn0, ftoauc, ftoau0, ftoadc, fsfcuc, fsfcu0, fsfcdc, fsfcd0, sfbmc, sfdfc, sfbm0, sfdf0, suvbfc, suvbf0)
This subroutine computes the shortwave radiative fluxes using two-stream method.
This module sets up absorption coeffients for band 18: 4000-4650 cm-1 (low - h2o,ch4; high - ch4) ...
subroutine taumol22
The subroutine computes the optical depth in band 22: 7700-8050 cm-1 (low - h2o,o2; high - o2) ...
real(kind=kind_phys), dimension(nblow:nbhgh), public b0s
single scattering albedo coefficients
derived type for SW fluxes' column profiles (at layer interfaces)
integer, dimension(nblow:nbhgh) idxebc
band index for cld prop
subroutine cldprop(cfrac, cliqp, reliq, cicep, reice, cdat1, cdat2, cdat3, cdat4, cf1, nlay, ipseed, taucw, ssacw, asycw, cldfrc, cldfmc)
This subroutine computes the cloud optical properties for each cloudy layer and g-point interval...
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), dimension(43, nblow:nbhgh), public ssaice2
single scattering albedo coefficients
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
logarithm,ln(p), of a 59-level standard pressure profile
subroutine taumol26
The subroutine computes the optical depth in band 26: 22650-29000 cm-1 (low - nothing; high - nothing...
real(kind=kind_phys), dimension(ng25), public abso3b
o3
This module contains the reference pressures (in logarithm form) at 59 vertical levels (TOA is omitte...
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) ...
subroutine, public rswinit(me)
This subroutine initializes non-varying module variables, conversion factors, and look-up tables...
This module is for specifying the band structures and program parameters used by the RRTMG-SW scheme...
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 index
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 setcoef(pavel, tavel, h2ovmr, nlay, nlp1, laytrop, jp, jt, jt1, fac00, fac01, fac10, fac11, selffac, selffrac, indself, forfac, forfrac, indfor)
This subroutine computes various coefficients needed in radiative transfer calculation.
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) ...
real(kind=kind_phys), dimension(nblow:nbhgh), public strrat
coefficients for computing binary absorbing species
real(kind=kind_phys), dimension(46, nblow:nbhgh), public extice3
extinction coefficients
real(kind=kind_phys), dimension(5), public fbari
asymmetry coefficients
real(kind=kind_phys), dimension(5), public dbari
single scattering albedo coefficients
subroutine spcvrtm(ssolar, cosz, sntz, albbm, albdf, sfluxzen, cldfmc, cf1, cf0, taug, taur, tauae, ssaae, asyae, taucw, ssacw, asycw, nlay, nlp1, fxupc, fxdnc, fxup0, fxdn0, ftoauc, ftoau0, ftoadc, fsfcuc, fsfcu0, fsfcdc, fsfcd0, sfbmc, sfdfc, sfbm0, sfdf0, suvbfc, suvbf0)
This subroutine computes the shortwave radiative fluxes using two-stream method of h...
This module sets up absorption coeffients for band 29: 820-2600 cm-1 (low - h2o; high - co2) ...
This module sets up absorption coeffients for band 17: 3250-4000 cm-1 (low - h2o,co2; high - h2o...
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...
real(kind=kind_phys), dimension(58, nblow:nbhgh), public asyliq1
asymmetry coefficients
real(kind=kind_phys), dimension(46, nblow:nbhgh), public asyice3
asymmetry coefficients
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
real(kind=kind_phys), dimension(ngmax, mfs03, mfb03), target, public sfluxref03
spectral solar fluxes, j=1,2,...,5 for SW band number of 18,19,21,22,24
integer, parameter ipsdsw0
initial permutation seed used for sub-column cloud scheme
This module contains various indexes and coefficients for SW spectral bands, as well as the spectral ...
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), 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(58, nblow:nbhgh), public ssaliq1
single scattering albedo coefficients
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
integer, dimension(nblow:nbhgh) ngs
array contains values of NS16-NS29
integer, dimension(nblow:nbhgh), public ibx
band index (3rd index in array sfluxref described below)
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
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...
real(kind=kind_phys), dimension(46, nblow:nbhgh), public ssaice3
single scattering albedo coefficients
subroutine mcica_subcol(cldf, nlay, ipseed, lcloudy)
This subroutine computes the sub-colum cloud profile flag array.
real(kind=kind_phys), public a0r
optical depth coefficients
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), dimension(5), public cbari
single scattering albedo coefficients
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.
derived type for SW 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
reverse checking of band index for each g-point
real(kind=kind_phys), dimension(ng29), public absh2o
h2o
integer, parameter maxgas
maximum number 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...
integer, dimension(nblow:nbhgh), public layreffr
reference pressure level for each of the spectral bands
subroutine taumol(colamt, colmol, fac00, fac01, fac10, fac11, jp, jt, jt1, laytrop, forfac, forfrac, indfor, selffac, selffrac, indself, nlay, sfluxzen, taug, taur)
This subroutine calculates optical depths for gaseous absorption and rayleigh scattering subroutine...
real(kind=kind_phys), dimension(nblow:nbhgh), public b0r
single scattering albedo coefficients
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...
real(kind=kind_phys), dimension(5), public abari
extinction coefficients
This module sets up absorption coefficients for band 16: 2600-3250 cm-1 (low - h2o, ch4; high - ch4)
derived type for SW fluxes at surface
real(kind=kind_phys), dimension(ng24), public abso3b
o3
real(kind=kind_phys), dimension(59) tref
MLS standard atmosphere temperature at the standard pressure levels.
real(kind=kind_phys), dimension(5), public bbari
extinction coefficients
real(kind=kind_phys), dimension(nblow:nbhgh), public specwt
weighting parameters for major absorbers in each band
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), dimension(ng25), public rayl
rayleigh extinction coefficient
real(kind=kind_phys), parameter, public rayl
rayleigh extinction coefficient at