rrtm_ecrt_140gp.F90

Go to the documentation of this file.

!****************** SUBROUTINE RRTM_ECRT_140GP **************************

SUBROUTINE rrtm_ecrt_140gp &
 &( iplon, klon , klev, kcld &
 &, paer , paph , pap &
 &, pts  , pth  , pt &
 &, zemis, zemiw &
 &, pq   , pcco2, pozn, pcldf, ptaucld, ptclear &
 &, cldfrac,taucld,coldry,wkl,wx &
 &, tauaerl,pavel,tavel,pz,tz,tbound,nlayers,semiss,ireflect)

!     Reformatted for F90 by JJMorcrette, ECMWF, 980714

!     Read in atmospheric profile from ECMWF radiation code, and prepare it
!     for use in RRTM.  Set other RRTM input parameters.  Values are passed
!     back through existing RRTM arrays and commons.

!- Modifications

!     2000-05-15 Deborah Salmond  Speed-up


#include "tsmbkind.h"

USE parrrtm  , ONLY : jpband   ,jpg      ,jpxsec   ,jpgpt    ,jplay   ,&
            &jpinpx
USE yoerad   , ONLY : novlp
USE yoerdi   , ONLY : rcardi   ,rch4     ,rn2o    ,rcfc11  ,rcfc12
USE yoesw    , ONLY : raer

!------------------------------Arguments--------------------------------


IMPLICIT NONE


!     DUMMY INTEGER SCALARS
integer_m :: iplon
integer_m :: kcld

!     DUMMY REAL SCALARS
real_b :: ptclear

integer_m :: kidia                   ! First atmosphere index
integer_m :: kfdia                   ! Last atmosphere index
integer_m :: klon                    ! Number of atmospheres (longitudes)
integer_m :: klev                    ! Number of atmospheric layers
real_b :: paer(klon,6,klev)          ! Aerosol optical thickness
real_b :: pap(klon,klev)             ! Layer pressures (Pa)
real_b :: paph(klon,klev+1)          ! Interface pressures (Pa)
real_b :: pts(klon)                  ! Surface temperature (K)
real_b :: pth(klon,klev+1)           ! Interface temperatures (K)
real_b :: pt(klon,klev)              ! Layer temperature (K)
real_b :: zemis(klon)                ! Non-window surface emissivity
real_b :: zemiw(klon)                ! Window surface emissivity
real_b :: pq(klon,klev)              ! H2O specific humidity (mmr)
real_b :: pozn(klon,klev)            ! O3 mass mixing ratio
real_b :: pcco2                      ! CO2 mass mixing ratio
!      real rch4                       ! CH4 mass mixing ratio
!      real rn2o                       ! N2O mass mixing ratio
!      real rcfc11                     ! CFC11 mass mixing ratio
!      real rcfc12                     ! CFC12 mass mixing ratio
real_b :: pcldf(klon,klev)           ! Cloud fraction
real_b :: ptaucld(klon,klev,jpband)  ! Cloud optical depth
real_b :: cldfrac(jplay)             ! Cloud fraction
real_b :: taucld(jplay,jpband)       ! Spectral optical thickness
real_b :: coldry(jplay)
real_b :: wkl(jpinpx,jplay)
real_b :: wx(jpxsec,jplay)           ! Amount of trace gases

!- from AER
real_b :: tauaerl(jplay,jpband)

!- from PROFILE             
real_b :: pavel(jplay)
real_b :: tavel(jplay)
real_b :: pz(0:jplay)
real_b :: tz(0:jplay)
real_b :: tbound
integer_m :: nlayers

!- from SURFACE             
real_b :: semiss(jpband)
integer_m :: ireflect

real_b :: ztauaer(5)
real_b :: zc1j(0:klev)               ! total cloud from top and level k
integer_m ::  ixindx(jpinpx)         ! Indices of trace gases accounted for

real_b :: amd                  ! Effective molecular weight of dry air (g/mol)
real_b :: amw                  ! Molecular weight of water vapor (g/mol)
real_b :: amco2                ! Molecular weight of carbon dioxide (g/mol)
real_b :: amo                  ! Molecular weight of ozone (g/mol)
real_b :: amch4                ! Molecular weight of methane (g/mol)
real_b :: amn2o                ! Molecular weight of nitrous oxide (g/mol)
real_b :: amc11                ! Molecular weight of CFC11 (g/mol) - CFCL3
real_b :: amc12                ! Molecular weight of CFC12 (g/mol) - CF2CL2
real_b :: avgdro               ! Avogadro's number (molecules/mole)
real_b :: gravit               ! Gravitational acceleration (cm/sec2)

! Atomic weights for conversion from mass to volume mixing ratios; these
!  are the same values used in ECRT to assure accurate conversion to vmr
data amd   /  28.970_jprb    /
data amw   /  18.0154_jprb   /
data amco2 /  44.011_jprb    /
data amo   /  47.9982_jprb   /
data amch4 /  16.043_jprb    /
data amn2o /  44.013_jprb    /
data amc11 / 137.3686_jprb   /
data amc12 / 120.9140_jprb   /
data avgdro/ 6.02214e23_jprb /
data gravit/ 9.80665e02_jprb /

!     LOCAL INTEGER SCALARS
integer_m :: iatm, imol, ix, ixmax, j1, j2, jae, jb, jk, jl, l, jis
integer_m :: nmol, nxmol

!     LOCAL REAL SCALARS
real_b :: amm, zcldly, zclear, zcloud, zepsec

! ***

! *** mji
! Initialize all molecular amounts and aerosol optical depths to zero here, 
! then pass ECRT amounts into RRTM arrays below.

!      DATA ZWKL /MAXPRDW*0.0/
!      DATA ZWX  /MAXPROD*0.0/
!      DATA KREFLECT /0/

! Activate cross section molecules:
!     NXMOL     - number of cross-sections input by user
!     IXINDX(I) - index of cross-section molecule corresponding to Ith
!                 cross-section specified by user
!                 = 0 -- not allowed in RRTM
!                 = 1 -- CCL4
!                 = 2 -- CFC11
!                 = 3 -- CFC12
!                 = 4 -- CFC22
!      DATA KXMOL  /2/
!      DATA KXINDX /0,2,3,0,31*0/

!      IREFLECT=KREFLECT
!      NXMOL=KXMOL

!print *,'Just entering RRTM_ECRT_140GP KLEV=',KLEV,' IPLON=',IPLON

ireflect=0
nxmol=2

DO j1=1,35
  ixindx(j1)=0
  DO j2=1,klev
    wkl(j1,j2)=_zero_
  ENDDO
ENDDO
ixindx(2)=2
ixindx(3)=3
DO j1=1,jpxsec
  DO j2=1,klev
    wx(j1,j2)=_zero_
  ENDDO
ENDDO

!     Set parameters needed for RRTM execution:
iatm    = 0
!      IXSECT  = 1
!      NUMANGS = 0
!      IOUT    = -1
ixmax   = 4

!     Bands 6,7,8 are considered the 'window' and allowed to have a
!     different surface emissivity (as in ECMWF).  Eli wrote this part....
semiss(1)  = zemis(iplon)
semiss(2)  = zemis(iplon)
semiss(3)  = zemis(iplon)
semiss(4)  = zemis(iplon)
semiss(5)  = zemis(iplon)
semiss(6)  = zemiw(iplon)
semiss(7)  = zemiw(iplon)
semiss(8)  = zemiw(iplon)
semiss(9)  = zemis(iplon)
semiss(10) = zemis(iplon)
semiss(11) = zemis(iplon)
semiss(12) = zemis(iplon)
semiss(13) = zemis(iplon)
semiss(14) = zemis(iplon)
semiss(15) = zemis(iplon)
semiss(16) = zemis(iplon)

!print *,'after SEMISS'

!     Set surface temperature.  

tbound = pts(iplon)
!print *,'after TBOUND=',TBOUND

!     Install ECRT arrays into RRTM arrays for pressure, temperature,
!     and molecular amounts.  Pressures are converted from Pascals
!     (ECRT) to mb (RRTM).  H2O, CO2, O3 and trace gas amounts are 
!     converted from mass mixing ratio to volume mixing ratio.  CO2
!     converted with same dry air and CO2 molecular weights used in 
!     ECRT to assure correct conversion back to the proper CO2 vmr.
!     The dry air column COLDRY (in molec/cm2) is calculated from 
!     the level pressures PZ (in mb) based on the hydrostatic equation
!     and includes a correction to account for H2O in the layer.  The
!     molecular weight of moist air (amm) is calculated for each layer.
!     Note: RRTM levels count from bottom to top, while the ECRT input
!     variables count from the top down and must be reversed here.

nlayers = klev
nmol = 6
pz(0) = paph(iplon,klev+1)/100._jprb
tz(0) = pth(iplon,klev+1)
DO l = 1, klev
  pavel(l) = pap(iplon,klev-l+1)/100._jprb
  tavel(l) = pt(iplon,klev-l+1)
  pz(l) = paph(iplon,klev-l+1)/100._jprb
  tz(l) = pth(iplon,klev-l+1)
  wkl(1,l) = pq(iplon,klev-l+1)*amd/amw
  wkl(2,l) = pcco2*amd/amco2
  wkl(3,l) = pozn(iplon,klev-l+1)*amd/amo
  wkl(4,l) = rn2o*amd/amn2o
  wkl(6,l) = rch4*amd/amch4
  amm = (1-wkl(1,l))*amd + wkl(1,l)*amw
  coldry(l) = (pz(l-1)-pz(l))*1.e3_jprb*avgdro/(gravit*amm*(1+wkl(1,l)))
ENDDO

!print *,'after WKL'
!print 9001,((RAER(JIS,JAE),JAE=1,6),JIS=1,5)
9001 format(1x,6e12.5)


!- Fill RRTM aerosol arrays with operational ECMWF aerosols,
!  do the mixing and distribute over the 16 spectral intervals

DO l=1,klev
  jk=klev-l+1
!  print 9002,JK,(PAER(IPLON,JK,JAE),JAE=1,6)
9002 format(1x,i3,6e12.5)
  
  
!       DO JAE=1,5
  jae=1
  ztauaer(jae) =&
   &(raer(jae,1)*paer(iplon,1,jk)+raer(jae,2)*paer(iplon,2,jk)&
   &+raer(jae,3)*paer(iplon,3,jk)+raer(jae,4)*paer(iplon,4,jk)&
   &+raer(jae,5)*paer(iplon,5,jk)+raer(jae,6)*paer(iplon,6,jk))
!   &/(PAPH(IPLON,JK+1)-PAPH(IPLON,JK))
  tauaerl(l, 1)=ztauaer(1)
  tauaerl(l, 2)=ztauaer(1)
  jae=2
  ztauaer(jae) =&
   &(raer(jae,1)*paer(iplon,1,jk)+raer(jae,2)*paer(iplon,2,jk)&
   &+raer(jae,3)*paer(iplon,3,jk)+raer(jae,4)*paer(iplon,4,jk)&
   &+raer(jae,5)*paer(iplon,5,jk)+raer(jae,6)*paer(iplon,6,jk))
!   &/(PAPH(IPLON,JK+1)-PAPH(IPLON,JK))
  tauaerl(l, 3)=ztauaer(2)
  tauaerl(l, 4)=ztauaer(2)
  tauaerl(l, 5)=ztauaer(2)
  jae=3
  ztauaer(jae) =&
   &(raer(jae,1)*paer(iplon,1,jk)+raer(jae,2)*paer(iplon,2,jk)&
   &+raer(jae,3)*paer(iplon,3,jk)+raer(jae,4)*paer(iplon,4,jk)&
   &+raer(jae,5)*paer(iplon,5,jk)+raer(jae,6)*paer(iplon,6,jk))
!   &/(PAPH(IPLON,JK+1)-PAPH(IPLON,JK))
  tauaerl(l, 6)=ztauaer(3)
  tauaerl(l, 8)=ztauaer(3)
  tauaerl(l, 9)=ztauaer(3)
  jae=4
  ztauaer(jae) =&
   &(raer(jae,1)*paer(iplon,1,jk)+raer(jae,2)*paer(iplon,2,jk)&
   &+raer(jae,3)*paer(iplon,3,jk)+raer(jae,4)*paer(iplon,4,jk)&
   &+raer(jae,5)*paer(iplon,5,jk)+raer(jae,6)*paer(iplon,6,jk))
!   &/(PAPH(IPLON,JK+1)-PAPH(IPLON,JK))
  tauaerl(l, 7)=ztauaer(4)
  jae=5
  ztauaer(jae) =&
   &(raer(jae,1)*paer(iplon,1,jk)+raer(jae,2)*paer(iplon,2,jk)&
   &+raer(jae,3)*paer(iplon,3,jk)+raer(jae,4)*paer(iplon,4,jk)&
   &+raer(jae,5)*paer(iplon,5,jk)+raer(jae,6)*paer(iplon,6,jk))
!   &/(PAPH(IPLON,JK+1)-PAPH(IPLON,JK))
!       END DO
  tauaerl(l,10)=ztauaer(5)
  tauaerl(l,11)=ztauaer(5)
  tauaerl(l,12)=ztauaer(5)
  tauaerl(l,13)=ztauaer(5)
  tauaerl(l,14)=ztauaer(5)
  tauaerl(l,15)=ztauaer(5)
  tauaerl(l,16)=ztauaer(5)
!  print 9003,L,(ZTAUAER(JAE),JAE=1,5)
9003 format(1x,'rrtm_ecrt ZTAUAER:',i3,5e13.6)  
ENDDO

DO l = 1, klev
!- Set cross section molecule amounts from ECRT; convert to vmr
  wx(2,l) = rcfc11*amd/amc11
  wx(3,l) = rcfc12*amd/amc12
!-- DS_000515  
END DO

!- Here, all molecules in WKL and WX are in volume mixing ratio; convert to
!  molec/cm2 based on COLDRY for use in RRTM
DO imol = 1, nmol
  DO l = 1, klev
!-- DS_000515  
    wkl(imol,l) = coldry(l) * wkl(imol,l)
  END DO  
ENDDO
  
DO ix = 1,jpxsec
  IF (ixindx(ix)  /=  0) THEN
!-- DS_000515  
    DO l=1 , klev
      wx(ixindx(ix),l) = coldry(l) * wx(ix,l) * 1.e-20_jprb
    END DO  
  ENDIF
ENDDO


!- Approximate treatment for various cloud overlaps
zclear=_one_
zcloud=_zero_
zc1j(0)=_zero_
zepsec=1.e-03_jprb
jl=iplon

IF (novlp == 1) THEN

  DO jk=1,klev
    IF (pcldf(jl,jk) > zepsec) THEN
      zcldly=pcldf(jl,jk)
      zclear=zclear &
       &*(_one_-max( zcldly , zcloud ))&
       &/(_one_-min( zcloud , _one_-zepsec ))
      zcloud = zcldly
      zc1j(jk)= _one_ - zclear
    ELSE
      zcldly=_zero_
      zclear=zclear &
       &*(_one_-max( zcldly , zcloud ))&
       &/(_one_-min( zcloud , _one_-zepsec ))
      zcloud = zcldly
      zc1j(jk)= _one_ - zclear
    ENDIF
  ENDDO

ELSEIF (novlp == 2) THEN

  DO jk=1,klev
    IF (pcldf(jl,jk) > zepsec) THEN
      zcldly=pcldf(jl,jk)
      zcloud = max( zcldly , zcloud )
      zc1j(jk) = zcloud
    ELSE
      zcldly=_zero_
      zcloud = max( zcldly , zcloud )
      zc1j(jk) = zcloud
    ENDIF
  ENDDO

ELSEIF (novlp == 3) THEN

  DO jk=1,klev
    IF (pcldf(jl,jk) > zepsec) THEN
      zcldly=pcldf(jl,jk)
      zclear = zclear * (_one_-zcldly)
      zcloud = _one_ - zclear
      zc1j(jk) = zcloud
    ELSE
      zcldly=_zero_
      zclear = zclear * (_one_-zcldly)
      zcloud = _one_ - zclear
      zc1j(jk) = zcloud
    ENDIF
  ENDDO

ENDIF
ptclear=_one_-zc1j(klev)

! Transfer cloud fraction and cloud optical depth to RRTM arrays; 
! invert array index for pcldf to go from bottom to top for RRTM

!- clear-sky column
IF (ptclear  >  _one_-zepsec) THEN
  kcld=0
  DO l = 1, klev
    cldfrac(l) = _zero_
  ENDDO
  DO jb=1,jpband
    DO l=1,klev
      taucld(l,jb) = _zero_
    ENDDO
  ENDDO

ELSE

!- cloudy column
!   The diffusivity factor (Savijarvi, 1997) on the cloud optical 
!   thickness TAUCLD has already been applied in RADLSW

  kcld=1
  DO l=1,klev
    cldfrac(l) = pcldf(iplon,l)
  ENDDO
  DO jb=1,jpband
    DO l=1,klev
      taucld(l,jb) = ptaucld(iplon,l,jb)
    ENDDO
  ENDDO

ENDIF

!     ------------------------------------------------------------------

RETURN
END SUBROUTINE rrtm_ecrt_140gp