#include subroutine radctl(hbuf ,clat ,coslat ,lat ,lwup , $ pmid ,pint ,pmln ,piln ,t , $ h2ommr ,cld ,effcld ,clwp ,coszrs , $ asdir ,asdif ,aldir ,aldif ,fsns , $ qrs ,qrl ,flwds ,rel ,rei , $ fice ,sols ,soll ,solsd ,solld , $ absnxt ,abstot ,emstot ,diag) !----------------------------------------------------------------------- ! ! Driver for radiation computation. ! ! Radiation uses cgs units, so conversions must be done from ! model fields to radiation fields. ! ! Calling sequence: ! ! radinp Converts units of model fields and computes ozone ! mixing ratio for solar scheme ! ! radcsw Performs solar computation ! radalb Computes surface albedos ! radded Computes delta-Eddington solution ! radclr Computes diagnostic clear sky fluxes ! ! radclw Performs longwave computation ! ! radtpl Computes path quantities ! radems Computes emissivity ! radabs Computes absorptivity ! !---------------------------Code history-------------------------------- ! ! Original version: CCM1 ! Standardized: J. Rosinski, June 1992 ! Reviewed: J. Kiehl, B. Briegleb, August 1992 ! ! Modified: B. Briegleb, March 1995 to add aerosol ! to shortwave code ! ! Reviewed: J. Kiehl, April 1996 ! Reviewed: B. Briegleb, May 1996 ! !----------------------------------------------------------------------- ! ! $Id: radctl.F,v 1.1 2001/08/21 16:38:31 bshen Exp $ ! $Author: bshen $ ! use mod_diag #include !------------------------------Parameters------------------------------- #include #include #include !------------------------------Commons---------------------------------- ! SJL #include !----------------------------------------------------------------------- #include !------------------------------Arguments-------------------------------- ! ! Input arguments ! ! SJL #if (defined R4BYTE) ! 4-byte real*4 absnxt(plond,plev,4) ! Nearest layer absorptivities real*4 abstot(plond,plevp,plevp) ! Non-adjacent layer absorptivites real*4 emstot(plond,plevp) ! Total emissivity real*4 n2o(plond,plev) ! nitrous oxide mass mixing ratio real*4 ch4(plond,plev) ! methane mass mixing ratio real*4 cfc11(plond,plev) ! cfc11 mass mixing ratio real*4 cfc12(plond,plev) ! cfc12 mass mixing ratio #else real absnxt(plond,plev,4) ! Nearest layer absorptivities real abstot(plond,plevp,plevp) ! Non-adjacent layer absorptivites real emstot(plond,plevp) ! Total emissivity real n2o(plond,plev) ! nitrous oxide mass mixing ratio real ch4(plond,plev) ! methane mass mixing ratio real cfc11(plond,plev) ! cfc11 mass mixing ratio real cfc12(plond,plev) ! cfc12 mass mixing ratio #endif ! SJL type (diag_type) diag(pdiag) ! History field attributes real*4 hbuf(*) ! History buffer integer lat ! Latitude row index real lwup(plond) ! Longwave up flux at surface real pmid(plond,plev) ! Model level pressures real pint(plond,plevp) ! Model interface pressures real pmln(plond,plev) ! Natural log of pmid real rel(plond,plev) ! liquid effective drop size (microns) real rei(plond,plev) ! ice effective drop size (microns) real fice(plond,plev) ! fractional ice content within cloud real piln(plond,plevp) ! Natural log of pint real t(plond,plev) ! Model level temperatures real h2ommr(plond,plev) ! Model level specific humidity real cld(plond,plevp) ! Fractional cloud cover real effcld(plond,plevp) ! Effective fractional cloud cover real clwp(plond,plev) ! Cloud liquid water path real coszrs(plond) ! Cosine solar zenith angle real asdir(plond) ! albedo shortwave direct real asdif(plond) ! albedo shortwave diffuse real aldir(plond) ! albedo longwave direct real aldif(plond) ! albedo longwave diffuse real clat ! current latitude(radians) real coslat ! cosine latitude ! ! Output solar arguments ! real fsns(plond) ! Surface absorbed solar flux real sols(plond) ! Downward solar rad onto surface (sw direct) real soll(plond) ! Downward solar rad onto surface (lw direct) real solsd(plond) ! Downward solar rad onto surface (sw diffuse) real solld(plond) ! Downward solar rad onto surface (lw diffuse) real qrs(plond,plev) ! Solar heating rate ! ! Output longwave arguments ! real qrl(plond,plev) ! Longwave cooling rate real flwds(plond) ! Surface down longwave flux ! !---------------------------Local variables----------------------------- ! integer i ! index real solin(plond) ! Solar incident flux real fsds(plond) ! Flux Shortwave Downwelling Surface real fsnt(plond) ! Net column abs solar flux at model top real fsntc(plond) ! Clear sky total column abs solar flux real fsnsc(plond) ! Clear sky surface abs solar flux real flnt(plond) ! Net outgoing lw flux at model top real flns(plond) ! Srf longwave cooling (up-down) flux real flntc(plond) ! Clear sky lw flux at model top real flnsc(plond) ! Clear sky lw flux at srf (up-down) real pbr(plond,plevr) ! Model mid-level pressures (dynes/cm2) real pnm(plond,plevrp) ! Model interface pressures (dynes/cm2) real o3vmr(plond,plevr) ! Ozone volume mixing ratio real o3mmr(plond,plevr) ! Ozone mass mixing ratio real plco2(plond,plevrp) ! Prs weighted CO2 path real plh2o(plond,plevrp) ! Prs weighted H2O path real tclrsf(plond,plevrp) ! Total clear sky fraction level to space real aermmr(plond,plevr) ! level aerosol mass mixing ratio real rh(plond,plevr) ! level relative humidity (fraction) real lwupcgs(plond) ! Upward longwave flux in cgs units ! JDC ADD real albedo(plond) ! surface albedo real pardif(plond) ! diffuse photosythetically active radiation (0.35-0.70 um) real pardir(plond) ! direct photosythetically active radiation (0.35-0.70 um) real wk(plond) ! temporary space real taucl(plond,plev) ! water cloud extinction optical depth real tauci(plond,plev) ! ice cloud extinction optical depth real eccf ! Earth/sun distance factor ! ! Declare local arrays to which model input arrays are interpolated here. ! Current default is none since radiation grid = model grid. ! !-------------------------------------------------------------------------- ! ! Interpolate model input arrays to radiation vertical grid. Currently this ! is a do-nothing routine because radiation grid = model grid. ! ! SJL ! call torgrid(pmid ,pint ,pmln ,piln ,t , & ! h2ommr ,cld ,effcld ,clwp , & ! pmid ,pint ,pmln ,piln ,t , & ! h2ommr ,cld ,effcld ,clwp ) ! SJL ! Interpolate ozone volume mixing ratio to model levels call radozn(lat ,pmid ,o3vmr ) call outfld(diag,iO3VMR,o3vmr ,plond, lat, hbuf) ! ! Set latitude dependent radiation input ! call radinp(pmid ,pint ,h2ommr ,cld ,o3vmr , $ pbr ,pnm ,plco2 ,plh2o ,tclrsf , $ eccf ,o3mmr ) ! ! Solar radiation computation ! if (dosw) then ! ! Specify aerosol mass mixing ratio ! call aermix(pnm ,aermmr ,rh ) call t_startf('radcsw') taucl = undef tauci = undef call radcsw(pnm ,h2ommr ,o3mmr ,aermmr ,rh , $ cld ,clwp ,rel ,rei ,fice , $ eccf ,coszrs ,asdir ,asdif ,aldir , $ aldif ,solin ,qrs ,fsns ,fsnt , $ fsds ,fsnsc ,fsntc ,sols ,soll , $ solsd ,solld ,pardif ,pardir ,taucl, tauci) call t_stopf('radcsw') ! ! Convert units of shortwave fields needed by rest of model from CGS to MKS ! do i=1,plon solin(i) = solin(i)*1.e-3 fsds(i) = fsds(i)*1.e-3 fsnt(i) = fsnt(i) *1.e-3 fsns(i) = fsns(i) *1.e-3 fsntc(i) = fsntc(i)*1.e-3 fsnsc(i) = fsnsc(i)*1.e-3 if( fsds(i) > 0.0 ) then albedo(i) = ( fsds(i)-fsns(i) ) / fsds(i) albedo(i) = min( 1.0, albedo(i) ) albedo(i) = max( 0.0, albedo(i) ) else albedo(i) = undef end if end do ! ! Dump shortwave radiation information to history tape buffer (diagnostics) ! call outfld(diag,iSOLIN, solin ,plond,lat,hbuf) call outfld(diag,iFSDS, fsds ,plond,lat,hbuf) call outfld(diag,iFSNT, fsnt ,plond,lat,hbuf) call outfld(diag,iFSNS, fsns ,plond,lat,hbuf) call outfld(diag,iFSNTC,fsntc ,plond,lat,hbuf) call outfld(diag,iFSNSC,fsnsc ,plond,lat,hbuf) ! JDC ADD call outfld(diag,iALBEDO, albedo ,plond,lat,hbuf) call outfld(diag,iPARDF, pardif ,plond,lat,hbuf) call outfld(diag,iPARDR, pardir ,plond,lat,hbuf) call outfld(diag,iTAUCLI, tauci ,plond,lat,hbuf) call outfld(diag,iTAUCLW, taucl ,plond,lat,hbuf) if (diag(iOSR)%pick .gt. 0) then do i = 1, plond wk(i) = solin(i) - fsnt(i) enddo call outfld(diag,iOSR, wk ,plond,lat,hbuf) endif if (diag(iOSRCLR)%pick .gt. 0) then do i = 1, plond wk(i) = fsntc(i) - fsnt(i) enddo call outfld(diag,iOSRCLR, wk ,plond,lat,hbuf) endif end if ! ! Longwave radiation computation ! if (dolw) then ! ! Convert upward longwave flux units to CGS ! do i=1,plon lwupcgs(i) = lwup(i)*1000. end do ! ! Specify trace gas mixing ratios ! call trcmix(pmid, clat, coslat, n2o, ch4, cfc11, cfc12) ! call radclw(lat ,lwupcgs ,t ,h2ommr ,o3vmr , $ pbr ,pnm ,pmln ,piln ,plco2 , $ plh2o ,n2o ,ch4 ,cfc11 ,cfc12 , $ effcld ,tclrsf ,qrl ,flns ,flnt , $ flnsc ,flntc ,flwds , $ absnxt ,abstot ,emstot ) ! ! Convert units of longwave fields needed by rest of model from CGS to MKS ! do i=1,plon flnt(i) = flnt(i)*1.e-3 flns(i) = flns(i)*1.e-3 flntc(i) = flntc(i)*1.e-3 flnsc(i) = flnsc(i)*1.e-3 flwds(i) = flwds(i)*1.e-3 end do ! ! Dump longwave radiation information to history tape buffer (diagnostics) ! call outfld(diag,iFLNT, flnt ,plond,lat,hbuf) call outfld(diag,iFLNTC,flntc ,plond,lat,hbuf) call outfld(diag,iFLNS, flns ,plond,lat,hbuf) call outfld(diag,iFLNSC,flnsc ,plond,lat,hbuf) end if ! ! Interpolate radiation output arrays to model vertical grid. Currently this ! is a do-nothing routine because radiation grid = model grid. ! SJL ! call fmrgrid(qrs ,qrl , & ! qrs ,qrl ) return end