298 module module_radiation_driver
302 use physcons
, only : eps => con_eps, &
303 & epsm1 => con_epsm1, &
304 & fvirt => con_fvirt &
306 use funcphys
, only : fpvs
311 use module_radiation_aerosols
,only : nf_aesw,
nf_aelw,
setaer, &
333 character(40),
parameter :: &
334 & VTAGRAD=
'NCEP-Radiation_driver v5.2 Jan 2013 '
359 integer :: month0=0, iyear0=0, monthd=0
365 integer,
parameter ::
ltp = 0
367 logical,
parameter :: lextop = (
ltp > 0)
385 subroutine radinit( si, NLAY, me )
496 integer,
intent(in) :: NLAY, me
498 real (kind=kind_phys),
intent(in) :: si(:)
517 print *,
' NEW RADIATION PROGRAM STRUCTURES BECAME OPER. ',
520 print *,
' - Selected Control Flag settings: ICTMflg=',
ictmflg,
529 print *,
' LTP =',
ltp,
', add extra top layer =',lextop
532 print *,
' Data usage is limited by initial condition!'
533 print *,
' No volcanic aerosols'
537 print *,
' - ISUBCLW=',
isubclw,
' No McICA, use grid ',
538 'averaged cloud in LW radiation'
540 print *,
' - ISUBCLW=',
isubclw,
' Use McICA with fixed ',
541 'permutation seeds for LW random number generator'
543 print *,
' - ISUBCLW=',
isubclw,
' Use McICA with random ',
544 'permutation seeds for LW random number generator'
546 print *,
' - ERROR!!! ISUBCLW=',
isubclw,
' is not a ',
552 print *,
' - ISUBCSW=',
isubcsw,
' No McICA, use grid ',
553 'averaged cloud in SW radiation'
555 print *,
' - ISUBCSW=',
isubcsw,
' Use McICA with fixed ',
556 'permutation seeds for SW random number generator'
558 print *,
' - ISUBCSW=',
isubcsw,
' Use McICA with random ',
559 'permutation seeds for SW random number generator'
561 print *,
' - ERROR!!! ISUBCSW=',
isubcsw,
' is not a ',
567 print *,
' - *** Notice *** ISUBCSW /= ISUBCLW !!!',
625 subroutine radupdate( idate,jdate,deltsw,deltim,lsswr, me, &
626 & slag,sdec,cdec,solcon)
690 integer,
intent(in) :: idate(:), jdate(:), me
691 logical,
intent(in) :: lsswr
693 real (kind=kind_phys),
intent(in) :: deltsw, deltim
696 real (kind=kind_phys),
intent(out) :: slag, sdec, cdec, solcon
699 integer :: iyear, imon, iday, ihour
700 integer :: kyear, kmon, kday, khour
731 if ( month0 /= imon )
then
744 elseif ( iyear0 /= iyear )
then
747 lsol_chg = (
isolar==4 .and. lmon_chg )
753 & ( jdate,kyear,deltsw,deltim,lsol_chg, me,
755 & slag,sdec,cdec,solcon
768 if ( monthd /= kmon )
then
945 & xlon,xlat,tsfc,snowd,sncovr,snoalb,zorl,hprim,
953 & htrlw,topflw,sfcflw,tsflw,semis,cldcov,
955 &, htrlw0,htrsw0,htrswb,htrlwb
1245 integer,
intent(in) :: IX,IM, LM, NTRAC, NFXR, me,
1246 & ntoz, ntcw, ncld, ipt, kdt
1247 integer,
intent(in) :: icsdsw(im), icsdlw(im), jdate(8)
1249 logical,
intent(in) :: lsswr, lslwr, lssav, lprnt, uni_cld, &
1252 real (kind=kind_phys),
dimension(IX,LM+1),
intent(in) :: prsi
1254 real (kind=kind_phys),
dimension(IX,LM),
intent(in) :: prsl,
1255 & prslk, tgrs, qgrs, vvl, fcice, frain, rrime, deltaq, cnvw,
1257 real (kind=kind_phys),
dimension(IM),
intent(in) :: flgmin
1258 real(kind=kind_phys),
intent(in) :: sup
1260 real (kind=kind_phys),
dimension(IM),
intent(in) :: slmsk,
1261 & xlon, xlat, tsfc, snowd, zorl, hprim, alvsf, alnsf, alvwf,
1262 & alnwf, facsf, facwf, cv, cvt, cvb, fice, tisfc,
1263 & sncovr, snoalb, sinlat, coslat
1265 real (kind=kind_phys),
intent(in) :: solcon, dtlw, dtsw, solhr,
1266 & tracer(ix,lm,ntrac)
1268 real (kind=kind_phys),
dimension(IX,LM),
intent(inout):: cldcov
1271 real (kind=kind_phys),
dimension(IX,LM),
intent(out):: htrsw,htrlw
1273 real (kind=kind_phys),
dimension(IX,4),
intent(out) :: dswcmp,
1276 real (kind=kind_phys),
dimension(IM),
intent(out):: tsflw,
1277 & sfalb, semis, coszen, coszdg
1279 type(
topfsw_type),
dimension(IM),
intent(out) :: topfsw
1280 type(
sfcfsw_type),
dimension(IM),
intent(out) :: sfcfsw
1282 type(
topflw_type),
dimension(IM),
intent(out) :: topflw
1283 type(
sfcflw_type),
dimension(IM),
intent(out) :: sfcflw
1286 real (kind=kind_phys),
intent(inout) :: fluxr(ix,nfxr)
1289 real (kind=kind_phys),
dimension(IX,LM,NBDSW),
optional,
1290 &
intent(out) :: htrswb
1291 real (kind=kind_phys),
dimension(IX,LM,NBDLW),
optional,
1292 &
intent(out) :: htrlwb
1293 real (kind=kind_phys),
dimension(ix,lm),
optional, &
1294 & intent(out) :: htrlw0
1295 real (kind=kind_phys),
dimension(ix,lm),
optional, &
1296 & intent(out) :: htrsw0
1299 real (kind=kind_phys),
dimension(IM,LM+1+LTP):: plvl, tlvl
1301 real (kind=kind_phys),
dimension(IM,LM+LTP) :: plyr, tlyr, qlyr, &
1302 & olyr, rhly, qstl, vvel, clw, prslk1, tem2da, tem2db, tvly
1304 real (kind=kind_phys),
dimension(IM) :: tsfa, cvt1, cvb1, tem1d, &
1305 & sfcemis, tsfg, tskn
1307 real (kind=kind_phys),
dimension(IM,LM+LTP,NF_CLDS) :: clouds
1308 real (kind=kind_phys),
dimension(IM,LM+LTP,NF_VGAS) :: gasvmr
1309 real (kind=kind_phys),
dimension(IM, NF_ALBD) :: sfcalb
1310 real (kind=kind_phys),
dimension(IM, NSPC1) :: aerodp
1311 real (kind=kind_phys),
dimension(IM,LM+LTP,NTRAC) :: tracer1
1313 real (kind=kind_phys),
dimension(IM,LM+LTP,NBDSW,NF_AESW)::faersw
1314 real (kind=kind_phys),
dimension(IM,LM+LTP,NBDLW,NF_AELW)::faerlw
1316 real (kind=kind_phys),
dimension(IM,LM+LTP) :: htswc
1317 real (kind=kind_phys),
dimension(IM,LM+LTP) :: htlwc
1319 real (kind=kind_phys),
dimension(IM,LM+LTP) :: gcice, grain, grime
1323 real (kind=kind_phys),
dimension(IM,LM+LTP) :: htsw0
1326 real (kind=kind_phys),
dimension(IM,LM+LTP,NBDSW) :: htswb
1328 real (kind=kind_phys),
dimension(IM,LM+LTP) :: htlw0
1330 real (kind=kind_phys),
dimension(IM,LM+LTP,NBDLW) :: htlwb
1332 real (kind=kind_phys) :: raddt, es, qs, delt, tem0d, cldsa(im,5)
1334 integer :: i, j, k, k1, lv, icec, itop, ibtc, nday, idxday(im), &
1335 & mbota(IM,3), mtopa(IM,3), LP1, nb, LMK, LMP, kd, lla, llb,
1382 raddt = min(dtsw, dtlw)
1404 if (
itsfc == 0 )
then
1426 plvl(i,k1) = 0.01 * prsi(i,k)
1427 plyr(i,k1) = 0.01 * prsl(i,k)
1428 tlyr(i,k1) = tgrs(i,k)
1429 prslk1(i,k1) = prslk(i,k)
1433 es = min( prsl(i,k), fpvs( tgrs(i,k) ) )
1434 qs = max(
qmin, eps * es / (prsl(i,k) + epsm1*es) )
1435 rhly(i,k1) = max( 0.0, min( 1.0, max(
qmin, qgrs(i,k))/qs ) )
1443 tracer1(i,k1,j) = tracer(i,k,j)
1450 plvl(i,lp1+kd) = 0.01 * prsi(i,lp1)
1457 plyr(i,lyb) = 0.5 * plvl(i,lla)
1458 tlyr(i,lyb) = tlyr(i,lya)
1460 prslk1(i,lyb) = (plyr(i,lyb)*0.00001) ** rocp
1462 rhly(i,lyb) = rhly(i,lya)
1463 qstl(i,lyb) = qstl(i,lya)
1469 tracer1(i,lyb,j) = tracer1(i,lya,j)
1481 gcice(i,k1) = fcice(i,k)
1482 grain(i,k1) = frain(i,k)
1483 grime(i,k1) = rrime(i,k)
1489 gcice(i,lyb) = fcice(i,lya)
1490 grain(i,lyb) = frain(i,lya)
1491 grime(i,lyb) = rrime(i,lya)
1502 olyr(i,k) = max(
qmin, tracer1(i,k,ntoz) )
1523 & ( xlon,sinlat,coslat,solhr, im, me,
1544 & ( plvl, xlon, xlat,
1554 tem2da(i,k) = log( plyr(i,k) )
1555 tem2db(i,k) = log( plvl(i,k) )
1563 tem2da(i,1) = log( plyr(i,1) )
1565 tsfa(i) = tlyr(i,lmk)
1566 tlvl(i,1) = tlyr(i,1)
1567 tlvl(i,lmp) = tskn(i)
1574 qlyr(i,k1) = max( tem1d(i), qgrs(i,k) )
1575 tem1d(i) = min(
qme5, qlyr(i,k1) )
1576 tvly(i,k1) = tgrs(i,k) * (1.0 + fvirt*qlyr(i,k1))
1582 qlyr(i,lyb) = qlyr(i,lya)
1583 tvly(i,lyb) = tvly(i,lya)
1589 tlvl(i,k) = tlyr(i,k) + (tlyr(i,k-1) - tlyr(i,k))
1599 tem2da(i,1) = log( plyr(i,1) )
1600 tem2db(i,1) = log( plvl(i,1) )
1603 tlvl(i,lmp) = tlyr(i,lmk)
1608 qlyr(i,k) = max( tem1d(i), qgrs(i,k) )
1609 tem1d(i) = min(
qme5, qlyr(i,k) )
1610 tvly(i,k) = tgrs(i,k) * (1.0 + fvirt*qlyr(i,k))
1616 qlyr(i,lyb) = qlyr(i,lya)
1617 tvly(i,lyb) = tvly(i,lya)
1623 tlvl(i,k+1) = tlyr(i,k) + (tlyr(i,k+1) - tlyr(i,k))
1635 if (coszen(i) >= 0.0001)
then
1654 & ( plvl,plyr,prslk1,tvly,rhly,slmsk,tracer1,xlon,xlat,
1657 & faersw,faerlw,aerodp
1684 clw(i,k) = clw(i,k) + tracer1(i,k,lv)
1691 if ( clw(i,k) <
epsq ) clw(i,k) = 0.0
1704 clw(i,k) = clw(i,k) + cnvw(i,k)
1714 & ( plyr,plvl,tlyr,tvly,qlyr,qstl,rhly,clw,
1718 & clouds,cldsa,mtopa,mbota
1726 & ( plyr,plvl,tlyr,tvly,qlyr,qstl,rhly,clw,
1730 & clouds,cldsa,mtopa,mbota
1737 & ( plyr,plvl,tlyr,tvly,qlyr,qstl,rhly,clw,cnvw,cnvc,
1742 & clouds,cldsa,mtopa,mbota
1752 cvt1(i) = 0.01 * cvt(i)
1753 cvb1(i) = 0.01 * cvb(i)
1762 vvel(i,k1) = 0.01 * vvl(i,k)
1768 vvel(i,lyb) = vvel(i,lya)
1776 & ( plyr,plvl,tlyr,rhly,vvel,cv,cvt1,cvb1,
1780 & clouds,cldsa,mtopa,mbota
1796 & ( slmsk,snowd,sncovr,snoalb,zorl,coszen,tsfg,tsfa,hprim,
1806 sfalb(i) = max(0.01, 0.5 * (sfcalb(i,2) + sfcalb(i,4)))
1814 if (
present(htrswb) .and.
present(htrsw0))
then
1818 & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr,
1823 & htswc,topfsw,sfcfsw
1826 &, hsw0=htsw0,hswb=htswb,fdncmp=scmpsw
1834 htrswb(i,k,j) = htswb(i,k1,j)
1839 else if (
present(htrswb) .and. .not.
present(htrsw0) )
then
1843 & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr,
1848 & htswc,topfsw,sfcfsw
1851 &, hswb=htswb,fdncmp=scmpsw
1858 htrswb(i,k,j) = htswb(i,k1,j)
1863 else if (
present(htrsw0) .and. .not.
present(htrswb) )
then
1867 & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr,
1872 & htswc,topfsw,sfcfsw
1875 &, hsw0=htsw0,fdncmp=scmpsw
1882 & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr,
1887 & htswc,topfsw,sfcfsw
1899 htrsw(i,k) = htswc(i,k1)
1902 if (
present(htrsw0))
then
1906 htrsw0(i,k) = htsw0(i,k1)
1914 dswcmp(i,1) = scmpsw(i)%nirbm
1915 dswcmp(i,2) = scmpsw(i)%nirdf
1916 dswcmp(i,3) = scmpsw(i)%visbm
1917 dswcmp(i,4) = scmpsw(i)%visdf
1919 uswcmp(i,1) = scmpsw(i)%nirbm * sfcalb(i,1)
1920 uswcmp(i,2) = scmpsw(i)%nirdf * sfcalb(i,2)
1921 uswcmp(i,3) = scmpsw(i)%visbm * sfcalb(i,3)
1922 uswcmp(i,4) = scmpsw(i)%visdf * sfcalb(i,4)
1935 scmpsw =
cmpfsw_type( 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 )
1947 if (
present(htrswb) )
then
1956 if (
present(htrsw0) )
then
1976 & ( xlon,xlat,slmsk,snowd,sncovr,zorl,tsfg,tsfa,hprim,
1985 if (
present(htrlwb) .and.
present(htrlw0) )
then
1989 & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr,
1993 & htlwc,topflw,sfcflw
2005 htrlwb(i,k,j) = htlwb(i,k1,j)
2010 else if (
present(htrlwb) .and. .not.
present(htrlw0) )
then
2014 & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr,
2018 & htlwc,topflw,sfcflw
2029 htrlwb(i,k,j) = htlwb(i,k1,j)
2033 else if (
present(htrlw0) .and. .not.
present(htrlwb) )
then
2038 & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr,
2042 & htlwc,topflw,sfcflw
2052 & ( plyr,plvl,tlyr,tlvl,qlyr,olyr,gasvmr,
2056 & htlwc,topflw,sfcflw
2064 semis(i) = sfcemis(i)
2073 htrlw(i,k) = htlwc(i,k1)
2077 if (
present(htrlw0))
then
2081 htrlw0(i,k) = htlw0(i,k1)
2095 fluxr(i,34) = fluxr(i,34) + dtsw*aerodp(i,1)
2096 fluxr(i,35) = fluxr(i,35) + dtsw*aerodp(i,2)
2097 fluxr(i,36) = fluxr(i,36) + dtsw*aerodp(i,3)
2098 fluxr(i,37) = fluxr(i,37) + dtsw*aerodp(i,4)
2099 fluxr(i,38) = fluxr(i,38) + dtsw*aerodp(i,5)
2100 fluxr(i,39) = fluxr(i,39) + dtsw*aerodp(i,6)
2108 fluxr(i,1 ) = fluxr(i,1 ) + dtlw * topflw(i)%upfxc
2109 fluxr(i,19) = fluxr(i,19) + dtlw * sfcflw(i)%dnfxc
2110 fluxr(i,20) = fluxr(i,20) + dtlw * sfcflw(i)%upfxc
2112 fluxr(i,28) = fluxr(i,28) + dtlw * topflw(i)%upfx0
2113 fluxr(i,30) = fluxr(i,30) + dtlw * sfcflw(i)%dnfx0
2114 fluxr(i,33) = fluxr(i,33) + dtlw * sfcflw(i)%upfx0
2123 if (coszen(i) > 0.)
then
2126 tem0d = dtsw * coszdg(i) / coszen(i)
2127 fluxr(i,2 ) = fluxr(i,2) + topfsw(i)%upfxc * tem0d
2128 fluxr(i,3 ) = fluxr(i,3) + sfcfsw(i)%upfxc * tem0d
2129 fluxr(i,4 ) = fluxr(i,4) + sfcfsw(i)%dnfxc * tem0d
2132 fluxr(i,21) = fluxr(i,21) + scmpsw(i)%uvbfc * tem0d
2133 fluxr(i,22) = fluxr(i,22) + scmpsw(i)%uvbf0 * tem0d
2136 fluxr(i,23) = fluxr(i,23) + topfsw(i)%dnfxc * tem0d
2139 fluxr(i,24) = fluxr(i,24) + scmpsw(i)%visbm * tem0d
2140 fluxr(i,25) = fluxr(i,25) + scmpsw(i)%visdf * tem0d
2141 fluxr(i,26) = fluxr(i,26) + scmpsw(i)%nirbm * tem0d
2142 fluxr(i,27) = fluxr(i,27) + scmpsw(i)%nirdf * tem0d
2145 fluxr(i,29) = fluxr(i,29) + topfsw(i)%upfx0 * tem0d
2146 fluxr(i,31) = fluxr(i,31) + sfcfsw(i)%upfx0 * tem0d
2147 fluxr(i,32) = fluxr(i,32) + sfcfsw(i)%dnfx0 * tem0d
2154 if (lsswr .or. lslwr)
then
2156 fluxr(i,17) = fluxr(i,17) + raddt * cldsa(i,4)
2157 fluxr(i,18) = fluxr(i,18) + raddt * cldsa(i,5)
2166 tem0d = raddt * cldsa(i,j)
2167 itop = mtopa(i,j) - kd
2168 ibtc = mbota(i,j) - kd
2169 fluxr(i, 8-j) = fluxr(i, 8-j) + tem0d
2170 fluxr(i,11-j) = fluxr(i,11-j) + tem0d * prsi(i,itop+kt)
2171 fluxr(i,14-j) = fluxr(i,14-j) + tem0d * prsi(i,ibtc+kb)
2172 fluxr(i,17-j) = fluxr(i,17-j) + tem0d * tgrs(i,itop)
2177 if (.not. uni_cld)
then
2181 cldcov(i,k) = clouds(i,k1,1)
2216 end subroutine grrad
2223 end module module_radiation_driver
subroutine, public sol_update
This subroutine computes solar parameters at forecast time.
integer, save iovrsw
cloud overlapping control flag for SW
Define type construct for radiation fluxes at surface.
integer, save isubcsw
sub-column cloud approx flag in SW radiation
real(kind=kind_phys) epsq
EPSQ=1.0e-12.
Define type construct for optional component downward fluxes at surface.
subroutine, public diagcld1
This subroutine computes cloud fractions for radiation calculations.
integer, save iaerflg
aerosol effect control flag
real(kind=kind_phys) qmin
QMIN=1.0e-10.
logical, save lcnorm
in-cld condensate control flag
logical, save lcrick
eliminating CRICK control flag
subroutine, public gas_update
This subroutine reads in 2-d monthly co2 data set for a specified year. Data are in a 15 degree lat/l...
subroutine, public radupdate(idate, jdate, deltsw, deltim, lsswr, me,
This subroutine calls many update subroutines to check and update radiation required but time varying...
real(kind=kind_phys) qme5
QME5=1.0e-7.
real, parameter prsmin
toa pressure minimum value in mb (hPa)
subroutine, public gas_init
This subroutine sets up ozone, co2, etc. parameters. If climatology ozone then read in monthly ozone ...
subroutine, public progcld2
This subroutine computes cloud related quantities using ferrier's prognostic cloud microphysics schem...
subroutine, public coszmn
This subroutine computes mean cos solar zenith angle over SW calling interval.
integer, save ialbflg
surface albedo scheme control flag
subroutine, public sfc_init
This subroutine is the initialization program for surface radiation related quantities (albedo...
subroutine, public swrad
This subroutine is the main SW radiation routine.
logical loz1st
first-time climatological ozone data read flag
integer, parameter, public nf_aelw
num of output fields for LW rad
Define type construct for optional radiation flux profiles.
integer, save ico2flg
co2 data source control flag
integer, save icmphys
cloud micorphysics scheme control flag
integer, save icldflg
cloud optical property scheme control flag
Define type construct for optional radiation flux profiles.
integer, save iemsflg
surface emissivity scheme control flag
subroutine, public grrad
This subroutine is the driver of radiation calculation subroutines. It sets up profile variables for ...
This module contains SW band parameters set up.
subroutine, public setalb
This subroutine computes four components of surface albedos (i.e., vis-nir, direct-diffused) accordin...
logical, save lnoprec
precip effect on radiation flag (Ferrier microphysics)
subroutine, public progcld1
This subroutine computes cloud related quantities using zhao/moorthi's prognostic cloud microphysics ...
subroutine, public rswinit
This subroutine initializes non-varying module variables, conversion factors, and look-up tables...
integer, parameter ltp
optional extra top layer on top of low ceiling models LTP=0: no extra top layer ...
subroutine, public radinit(si, NLAY, me)
This subroutine is the initialization of radiation calculations.
This module contains LW band parameters set up.
Define type construct for radiation fluxes at toa.
subroutine, public aer_init
The initialization program is to set up necessary parameters and working arrays.
subroutine, public sol_init
This subroutine initializes astronomy process, and set up module constants.
subroutine, public progcld3
This subroutine computes cloud related quantities using zhao/moorthi's prognostic cloud microphysics ...
integer, save ioznflg
ozone data source control flag
subroutine, public cld_init
This subroutine is an initialization program for cloud-radiation calculations and sets up boundary la...
subroutine, public setaer
This subroutine computes aerosols optical properties.
subroutine, public getozn
This subroutine sets up climatological ozone profile for radiation calculation. This code is original...
integer, parameter, public nf_clds
number of fields in cloud array
subroutine, public rlwinit
This subroutine performs calculations necessary for the initialization of the longwave model...
integer, save ictmflg
external data time/date control flag
real(kind=kind_phys) qme6
QME6=1.0e-7.
integer, save isolar
solar constant scheme control flag
integer, parameter, public nspc1
total+species
subroutine, public getgases
This subroutine sets up global distribution of radiation absorbing gases in volume mixing ratio...
integer itsfc
control flag for lw sfc air/ground interface temp setting
Define type construct for radiation fluxes at toa.
integer, save isubclw
sub-column cloud approx flag in LW radiation
integer, parameter, public nf_vgas
number of gas species
integer, save ivflip
vertical profile indexing flag
subroutine, public aer_update
This subroutine checks and updates time varying climatology aerosol data sets.
subroutine, public setemis
This subroutine computes surface emissivity for LW radiation.
subroutine, public lwrad
This subroutine is the main LW radiation routine.
Define type construct for radiation fluxes at surface.
integer, save iovrlw
cloud overlapping control flag for LW