phylmd/Radlwsw/sw2s.f

module sw2s_m

  IMPLICIT NONE

contains

  SUBROUTINE sw2s(knu, flag_aer, tauae, pizae, cgae, paki, palbd, palbp, &
       pcg, pcld, pclear, pdsig, pomega, poz, prmu, psec, ptau, pud, &
       pwv, pqs, pfdown, pfup)
    USE dimens_m
    USE dimphy
    USE raddim
    USE radepsi
    use swclr_m, only: swclr
    use swde_m, only: swde
    use swr_m, only: swr

    ! ------------------------------------------------------------------
    ! PURPOSE.
    ! --------

    ! THIS ROUTINE COMPUTES THE SHORTWAVE RADIATION FLUXES IN THE
    ! SECOND SPECTRAL INTERVAL FOLLOWING FOUQUART AND BONNEL (1980).

    ! METHOD.
    ! -------

    ! 1. COMPUTES REFLECTIVITY/TRANSMISSIVITY CORRESPONDING TO
    ! CONTINUUM SCATTERING
    ! 2. COMPUTES REFLECTIVITY/TRANSMISSIVITY CORRESPONDING FOR
    ! A GREY MOLECULAR ABSORPTION
    ! 3. LAPLACE TRANSFORM ON THE PREVIOUS TO GET EFFECTIVE AMOUNTS
    ! OF ABSORBERS
    ! 4. APPLY H2O AND U.M.G. TRANSMISSION FUNCTIONS
    ! 5. MULTIPLY BY OZONE TRANSMISSION FUNCTION

    ! REFERENCE.
    ! ----------

    ! SEE RADIATION'S PART OF THE ECMWF RESEARCH DEPARTMENT
    ! DOCUMENTATION, AND FOUQUART AND BONNEL (1980)

    ! AUTHOR.
    ! -------
    ! JEAN-JACQUES MORCRETTE  *ECMWF*

    ! MODIFICATIONS.
    ! --------------
    ! ORIGINAL : 89-07-14
    ! 94-11-15   J.-J. MORCRETTE    DIRECT/DIFFUSE ALBEDO
    ! ------------------------------------------------------------------
    ! * ARGUMENTS:

    INTEGER knu
    ! -OB
    DOUBLE PRECISION flag_aer
    DOUBLE PRECISION tauae(kdlon, kflev, 2)
    DOUBLE PRECISION pizae(kdlon, kflev, 2)
    DOUBLE PRECISION cgae(kdlon, kflev, 2)
    DOUBLE PRECISION paki(kdlon, 2)
    DOUBLE PRECISION palbd(kdlon, 2)
    DOUBLE PRECISION palbp(kdlon, 2)
    DOUBLE PRECISION pcg(kdlon, 2, kflev)
    DOUBLE PRECISION pcld(kdlon, kflev)
    DOUBLE PRECISION pclear(kdlon)
    DOUBLE PRECISION pdsig(kdlon, kflev)
    DOUBLE PRECISION pomega(kdlon, 2, kflev)
    DOUBLE PRECISION poz(kdlon, kflev)
    DOUBLE PRECISION pqs(kdlon, kflev)
    DOUBLE PRECISION prmu(kdlon)
    DOUBLE PRECISION psec(kdlon)
    DOUBLE PRECISION ptau(kdlon, 2, kflev)
    DOUBLE PRECISION pud(kdlon, 5, kflev+1)
    DOUBLE PRECISION pwv(kdlon, kflev)

    DOUBLE PRECISION pfdown(kdlon, kflev+1)
    DOUBLE PRECISION pfup(kdlon, kflev+1)

    ! * LOCAL VARIABLES:

    INTEGER iind2(2), iind3(3)
    DOUBLE PRECISION zcgaz(kdlon, kflev)
    DOUBLE PRECISION zfd(kdlon, kflev+1)
    DOUBLE PRECISION zfu(kdlon, kflev+1)
    DOUBLE PRECISION zg(kdlon)
    DOUBLE PRECISION zgg(kdlon)
    DOUBLE PRECISION zpizaz(kdlon, kflev)
    DOUBLE PRECISION zrayl(kdlon)
    DOUBLE PRECISION zray1(kdlon, kflev+1)
    DOUBLE PRECISION zray2(kdlon, kflev+1)
    DOUBLE PRECISION zref(kdlon)
    DOUBLE PRECISION zrefz(kdlon, 2, kflev+1)
    DOUBLE PRECISION zre1(kdlon)
    DOUBLE PRECISION zre2(kdlon)
    DOUBLE PRECISION zrj(kdlon, 6, kflev+1)
    DOUBLE PRECISION zrj0(kdlon, 6, kflev+1)
    DOUBLE PRECISION zrk(kdlon, 6, kflev+1)
    DOUBLE PRECISION zrk0(kdlon, 6, kflev+1)
    DOUBLE PRECISION zrl(kdlon, 8)
    DOUBLE PRECISION zrmue(kdlon, kflev+1)
    DOUBLE PRECISION zrmu0(kdlon, kflev+1)
    DOUBLE PRECISION zrmuz(kdlon)
    DOUBLE PRECISION zrneb(kdlon)
    DOUBLE PRECISION zr1(kdlon)
    DOUBLE PRECISION zr2(kdlon, 2)
    DOUBLE PRECISION zr3(kdlon, 3)
    DOUBLE PRECISION zr4(kdlon)
    DOUBLE PRECISION zr21(kdlon)
    DOUBLE PRECISION zr22(kdlon)
    DOUBLE PRECISION zs(kdlon)
    DOUBLE PRECISION ztauaz(kdlon, kflev)
    DOUBLE PRECISION zto1(kdlon)
    DOUBLE PRECISION ztr(kdlon, 2, kflev+1)
    DOUBLE PRECISION ztra1(kdlon, kflev+1)
    DOUBLE PRECISION ztra2(kdlon, kflev+1)
    DOUBLE PRECISION ztr1(kdlon)
    DOUBLE PRECISION ztr2(kdlon)
    DOUBLE PRECISION zw(kdlon)
    DOUBLE PRECISION zw1(kdlon)
    DOUBLE PRECISION zw2(kdlon, 2)
    DOUBLE PRECISION zw3(kdlon, 3)
    DOUBLE PRECISION zw4(kdlon)
    DOUBLE PRECISION zw5(kdlon)

    INTEGER jl, jk, k, jaj, ikm1, ikl, jn, jabs, jkm1
    INTEGER jref, jkl, jklp1, jajp, jkki, jkkp4, jn2j, iabs
    DOUBLE PRECISION zrmum1, zwh2o, zcneb, zaa, zbb, zrki, zre11

    ! * Prescribed Data:

    DOUBLE PRECISION rsun(2)
    SAVE rsun
    DOUBLE PRECISION rray(2, 6)
    SAVE rray
    DATA rsun(1)/0.441676d0/
    DATA rsun(2)/0.558324d0/
    DATA (rray(1,k), k=1, 6)/.428937d-01, .890743d+00, -.288555d+01, &
         .522744d+01, -.469173d+01, .161645d+01/
    DATA (rray(2,k), k=1, 6)/.697200d-02, .173297d-01, -.850903d-01, &
         .248261d+00, -.302031d+00, .129662d+00/

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

    ! *         1.     SECOND SPECTRAL INTERVAL (0.68-4.00 MICRON)
    ! -------------------------------------------


    ! *         1.1    OPTICAL THICKNESS FOR RAYLEIGH SCATTERING
    ! -----------------------------------------


    DO jl = 1, kdlon
       zrmum1 = 1. - prmu(jl)
       zrayl(jl) = rray(knu, 1) + zrmum1*(rray(knu,2)+zrmum1*(rray(knu, &
            3)+zrmum1*(rray(knu,4)+zrmum1*(rray(knu,5)+zrmum1*rray(knu,6)))))
    END DO


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

    ! *         2.    CONTINUUM SCATTERING CALCULATIONS
    ! ---------------------------------


    ! *         2.1   CLEAR-SKY FRACTION OF THE COLUMN
    ! --------------------------------


    CALL swclr(knu, flag_aer, tauae, pizae, cgae, palbp, pdsig, zrayl, &
         psec, zcgaz, zpizaz, zray1, zray2, zrefz, zrj0, zrk0, zrmu0, ztauaz, &
         ztra1, ztra2)


    ! *         2.2   CLOUDY FRACTION OF THE COLUMN
    ! -----------------------------


    CALL swr(knu, palbd, pcg, pcld, pomega, psec, ptau, zcgaz, &
         zpizaz, zray1, zray2, zrefz, zrj, zrk, zrmue, ztauaz, ztra1, ztra2)


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

    ! *         3.    SCATTERING CALCULATIONS WITH GREY MOLECULAR ABSORPTION
    ! ------------------------------------------------------


    jn = 2

    DO jabs = 1, 2


       ! *         3.1  SURFACE CONDITIONS
       ! ------------------


       DO jl = 1, kdlon
          zrefz(jl, 2, 1) = palbd(jl, knu)
          zrefz(jl, 1, 1) = palbd(jl, knu)
       END DO


       ! *         3.2  INTRODUCING CLOUD EFFECTS
       ! -------------------------


       DO jk = 2, kflev + 1
          jkm1 = jk - 1
          ikl = kflev + 1 - jkm1
          DO jl = 1, kdlon
             zrneb(jl) = pcld(jl, jkm1)
             IF (jabs==1 .AND. zrneb(jl)>2.*zeelog) THEN
                zwh2o = max(pwv(jl,jkm1), zeelog)
                zcneb = max(zeelog, min(zrneb(jl),1.-zeelog))
                zbb = pud(jl, jabs, jkm1)*pqs(jl, jkm1)/zwh2o
                zaa = max((pud(jl,jabs,jkm1)-zcneb*zbb)/(1.-zcneb), zeelog)
             ELSE
                zaa = pud(jl, jabs, jkm1)
                zbb = zaa
             END IF
             zrki = paki(jl, jabs)
             zs(jl) = exp(-zrki*zaa*1.66)
             zg(jl) = exp(-zrki*zaa/zrmue(jl,jk))
             ztr1(jl) = 0.
             zre1(jl) = 0.
             ztr2(jl) = 0.
             zre2(jl) = 0.

             zw(jl) = pomega(jl, knu, jkm1)
             zto1(jl) = ptau(jl, knu, jkm1)/zw(jl) + ztauaz(jl, jkm1)/zpizaz(jl, &
                  jkm1) + zbb*zrki

             zr21(jl) = ptau(jl, knu, jkm1) + ztauaz(jl, jkm1)
             zr22(jl) = ptau(jl, knu, jkm1)/zr21(jl)
             zgg(jl) = zr22(jl)*pcg(jl, knu, jkm1) + (1.-zr22(jl))*zcgaz(jl, jkm1)
             zw(jl) = zr21(jl)/zto1(jl)
             zref(jl) = zrefz(jl, 1, jkm1)
             zrmuz(jl) = zrmue(jl, jk)
          END DO

          CALL swde(zgg, zref, zrmuz, zto1, zw, zre1, zre2, ztr1, ztr2)

          DO jl = 1, kdlon

             zrefz(jl, 2, jk) = (1.-zrneb(jl))*(zray1(jl,jkm1)+zrefz(jl,2,jkm1)* &
                  ztra1(jl,jkm1)*ztra2(jl,jkm1))*zg(jl)*zs(jl) + zrneb(jl)*zre1(jl)

             ztr(jl, 2, jkm1) = zrneb(jl)*ztr1(jl) + (ztra1(jl,jkm1))*zg(jl)*(1.- &
                  zrneb(jl))

             zrefz(jl, 1, jk) = (1.-zrneb(jl))*(zray1(jl,jkm1)+zrefz(jl,1,jkm1)* &
                  ztra1(jl,jkm1)*ztra2(jl,jkm1)/(1.-zray2(jl,jkm1)*zrefz(jl,1, &
                  jkm1)))*zg(jl)*zs(jl) + zrneb(jl)*zre2(jl)

             ztr(jl, 1, jkm1) = zrneb(jl)*ztr2(jl) + (ztra1(jl,jkm1)/(1.-zray2(jl, &
                  jkm1)*zrefz(jl,1,jkm1)))*zg(jl)*(1.-zrneb(jl))

          END DO
       END DO

       ! *         3.3  REFLECT./TRANSMISSIVITY BETWEEN SURFACE AND LEVEL
       ! -------------------------------------------------


       DO jref = 1, 2

          jn = jn + 1

          DO jl = 1, kdlon
             zrj(jl, jn, kflev+1) = 1.
             zrk(jl, jn, kflev+1) = zrefz(jl, jref, kflev+1)
          END DO

          DO jk = 1, kflev
             jkl = kflev + 1 - jk
             jklp1 = jkl + 1
             DO jl = 1, kdlon
                zre11 = zrj(jl, jn, jklp1)*ztr(jl, jref, jkl)
                zrj(jl, jn, jkl) = zre11
                zrk(jl, jn, jkl) = zre11*zrefz(jl, jref, jkl)
             END DO
          END DO
       END DO
    END DO


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

    ! *         4.    INVERT GREY AND CONTINUUM FLUXES
    ! --------------------------------


    ! *         4.1   UPWARD (ZRK) AND DOWNWARD (ZRJ) PSEUDO-FLUXES
    ! ---------------------------------------------


    DO jk = 1, kflev + 1
       DO jaj = 1, 5, 2
          jajp = jaj + 1
          DO jl = 1, kdlon
             zrj(jl, jaj, jk) = zrj(jl, jaj, jk) - zrj(jl, jajp, jk)
             zrk(jl, jaj, jk) = zrk(jl, jaj, jk) - zrk(jl, jajp, jk)
             zrj(jl, jaj, jk) = max(zrj(jl,jaj,jk), zeelog)
             zrk(jl, jaj, jk) = max(zrk(jl,jaj,jk), zeelog)
          END DO
       END DO
    END DO

    DO jk = 1, kflev + 1
       DO jaj = 2, 6, 2
          DO jl = 1, kdlon
             zrj(jl, jaj, jk) = max(zrj(jl,jaj,jk), zeelog)
             zrk(jl, jaj, jk) = max(zrk(jl,jaj,jk), zeelog)
          END DO
       END DO
    END DO

    ! *         4.2    EFFECTIVE ABSORBER AMOUNTS BY INVERSE LAPLACE
    ! ---------------------------------------------


    DO jk = 1, kflev + 1
       jkki = 1
       DO jaj = 1, 2
          iind2(1) = jaj
          iind2(2) = jaj
          DO jn = 1, 2
             jn2j = jn + 2*jaj
             jkkp4 = jkki + 4

             ! *         4.2.1  EFFECTIVE ABSORBER AMOUNTS
             ! --------------------------


             DO jl = 1, kdlon
                zw2(jl, 1) = log(zrj(jl,jn,jk)/zrj(jl,jn2j,jk))/paki(jl, jaj)
                zw2(jl, 2) = log(zrk(jl,jn,jk)/zrk(jl,jn2j,jk))/paki(jl, jaj)
             END DO

             ! *         4.2.2  TRANSMISSION FUNCTION
             ! ---------------------


             CALL swtt1(knu, 2, iind2, zw2, zr2)

             DO jl = 1, kdlon
                zrl(jl, jkki) = zr2(jl, 1)
                zrl(jl, jkkp4) = zr2(jl, 2)
             END DO

             jkki = jkki + 1
          END DO
       END DO

       ! *         4.3    UPWARD AND DOWNWARD FLUXES WITH H2O AND UMG ABSORPTION
       ! ------------------------------------------------------


       DO jl = 1, kdlon
          pfdown(jl, jk) = zrj(jl, 1, jk)*zrl(jl, 1)*zrl(jl, 3) + &
               zrj(jl, 2, jk)*zrl(jl, 2)*zrl(jl, 4)
          pfup(jl, jk) = zrk(jl, 1, jk)*zrl(jl, 5)*zrl(jl, 7) + &
               zrk(jl, 2, jk)*zrl(jl, 6)*zrl(jl, 8)
       END DO
    END DO


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

    ! *         5.    MOLECULAR ABSORPTION ON CLEAR-SKY FLUXES
    ! ----------------------------------------


    ! *         5.1   DOWNWARD FLUXES
    ! ---------------


    jaj = 2
    iind3(1) = 1
    iind3(2) = 2
    iind3(3) = 3

    DO jl = 1, kdlon
       zw3(jl, 1) = 0.
       zw3(jl, 2) = 0.
       zw3(jl, 3) = 0.
       zw4(jl) = 0.
       zw5(jl) = 0.
       zr4(jl) = 1.
       zfd(jl, kflev+1) = zrj0(jl, jaj, kflev+1)
    END DO
    DO jk = 1, kflev
       ikl = kflev + 1 - jk
       DO jl = 1, kdlon
          zw3(jl, 1) = zw3(jl, 1) + pud(jl, 1, ikl)/zrmu0(jl, ikl)
          zw3(jl, 2) = zw3(jl, 2) + pud(jl, 2, ikl)/zrmu0(jl, ikl)
          zw3(jl, 3) = zw3(jl, 3) + poz(jl, ikl)/zrmu0(jl, ikl)
          zw4(jl) = zw4(jl) + pud(jl, 4, ikl)/zrmu0(jl, ikl)
          zw5(jl) = zw5(jl) + pud(jl, 5, ikl)/zrmu0(jl, ikl)
       END DO

       CALL swtt1(knu, 3, iind3, zw3, zr3)

       DO jl = 1, kdlon
          ! ZR4(JL) = EXP(-RSWCE*ZW4(JL)-RSWCP*ZW5(JL))
          zfd(jl, ikl) = zr3(jl, 1)*zr3(jl, 2)*zr3(jl, 3)*zr4(jl)* &
               zrj0(jl, jaj, ikl)
       END DO
    END DO


    ! *         5.2   UPWARD FLUXES
    ! -------------


    DO jl = 1, kdlon
       zfu(jl, 1) = zfd(jl, 1)*palbp(jl, knu)
    END DO

    DO jk = 2, kflev + 1
       ikm1 = jk - 1
       DO jl = 1, kdlon
          zw3(jl, 1) = zw3(jl, 1) + pud(jl, 1, ikm1)*1.66
          zw3(jl, 2) = zw3(jl, 2) + pud(jl, 2, ikm1)*1.66
          zw3(jl, 3) = zw3(jl, 3) + poz(jl, ikm1)*1.66
          zw4(jl) = zw4(jl) + pud(jl, 4, ikm1)*1.66
          zw5(jl) = zw5(jl) + pud(jl, 5, ikm1)*1.66
       END DO

       CALL swtt1(knu, 3, iind3, zw3, zr3)

       DO jl = 1, kdlon
          ! ZR4(JL) = EXP(-RSWCE*ZW4(JL)-RSWCP*ZW5(JL))
          zfu(jl, jk) = zr3(jl, 1)*zr3(jl, 2)*zr3(jl, 3)*zr4(jl)* &
               zrk0(jl, jaj, jk)
       END DO
    END DO


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

    ! *         6.     INTRODUCTION OF OZONE AND H2O CONTINUUM ABSORPTION
    ! --------------------------------------------------

    iabs = 3

    ! *         6.1    DOWNWARD FLUXES
    ! ---------------

    DO jl = 1, kdlon
       zw1(jl) = 0.
       zw4(jl) = 0.
       zw5(jl) = 0.
       zr1(jl) = 0.
       pfdown(jl, kflev+1) = ((1.-pclear(jl))*pfdown(jl,kflev+1)+pclear(jl)*zfd( &
            jl,kflev+1))*rsun(knu)
    END DO

    DO jk = 1, kflev
       ikl = kflev + 1 - jk
       DO jl = 1, kdlon
          zw1(jl) = zw1(jl) + poz(jl, ikl)/zrmue(jl, ikl)
          zw4(jl) = zw4(jl) + pud(jl, 4, ikl)/zrmue(jl, ikl)
          zw5(jl) = zw5(jl) + pud(jl, 5, ikl)/zrmue(jl, ikl)
          ! ZR4(JL) = EXP(-RSWCE*ZW4(JL)-RSWCP*ZW5(JL))
       END DO

       CALL swtt(knu, iabs, zw1, zr1)

       DO jl = 1, kdlon
          pfdown(jl, ikl) = ((1.-pclear(jl))*zr1(jl)*zr4(jl)*pfdown(jl,ikl)+ &
               pclear(jl)*zfd(jl,ikl))*rsun(knu)
       END DO
    END DO


    ! *         6.2    UPWARD FLUXES
    ! -------------

    DO jl = 1, kdlon
       pfup(jl, 1) = ((1.-pclear(jl))*zr1(jl)*zr4(jl)*pfup(jl,1)+pclear(jl)*zfu( &
            jl,1))*rsun(knu)
    END DO

    DO jk = 2, kflev + 1
       ikm1 = jk - 1
       DO jl = 1, kdlon
          zw1(jl) = zw1(jl) + poz(jl, ikm1)*1.66
          zw4(jl) = zw4(jl) + pud(jl, 4, ikm1)*1.66
          zw5(jl) = zw5(jl) + pud(jl, 5, ikm1)*1.66
          ! ZR4(JL) = EXP(-RSWCE*ZW4(JL)-RSWCP*ZW5(JL))
       END DO

       CALL swtt(knu, iabs, zw1, zr1)

       DO jl = 1, kdlon
          pfup(jl, jk) = ((1.-pclear(jl))*zr1(jl)*zr4(jl)*pfup(jl,jk)+pclear(jl)* &
               zfu(jl,jk))*rsun(knu)
       END DO
    END DO

  END SUBROUTINE sw2s

end module sw2s_m
1	guez	178	module sw2s_m
2
3	guez	81	IMPLICIT NONE
4	guez	24
5	guez	178	contains
6	guez	81
7	guez	178	SUBROUTINE sw2s(knu, flag_aer, tauae, pizae, cgae, paki, palbd, palbp, &
8			pcg, pcld, pclear, pdsig, pomega, poz, prmu, psec, ptau, pud, &
9			pwv, pqs, pfdown, pfup)
10			USE dimens_m
11			USE dimphy
12			USE raddim
13			USE radepsi
14			use swclr_m, only: swclr
15	guez	207	use swde_m, only: swde
16	guez	178	use swr_m, only: swr
17	guez	81
18	guez	178	! ------------------------------------------------------------------
19			! PURPOSE.
20			! --------
21	guez	81
22	guez	178	! THIS ROUTINE COMPUTES THE SHORTWAVE RADIATION FLUXES IN THE
23			! SECOND SPECTRAL INTERVAL FOLLOWING FOUQUART AND BONNEL (1980).
24	guez	81
25	guez	178	! METHOD.
26			! -------
27	guez	81
28	guez	178	! 1. COMPUTES REFLECTIVITY/TRANSMISSIVITY CORRESPONDING TO
29			! CONTINUUM SCATTERING
30			! 2. COMPUTES REFLECTIVITY/TRANSMISSIVITY CORRESPONDING FOR
31			! A GREY MOLECULAR ABSORPTION
32			! 3. LAPLACE TRANSFORM ON THE PREVIOUS TO GET EFFECTIVE AMOUNTS
33			! OF ABSORBERS
34			! 4. APPLY H2O AND U.M.G. TRANSMISSION FUNCTIONS
35			! 5. MULTIPLY BY OZONE TRANSMISSION FUNCTION
36	guez	81
37	guez	178	! REFERENCE.
38			! ----------
39	guez	81
40	guez	178	! SEE RADIATION'S PART OF THE ECMWF RESEARCH DEPARTMENT
41			! DOCUMENTATION, AND FOUQUART AND BONNEL (1980)
42	guez	81
43	guez	178	! AUTHOR.
44			! -------
45			! JEAN-JACQUES MORCRETTE ECMWF
46	guez	81
47	guez	178	! MODIFICATIONS.
48			! --------------
49			! ORIGINAL : 89-07-14
50			! 94-11-15 J.-J. MORCRETTE DIRECT/DIFFUSE ALBEDO
51			! ------------------------------------------------------------------
52			! * ARGUMENTS:
53	guez	81
54	guez	178	INTEGER knu
55			! -OB
56			DOUBLE PRECISION flag_aer
57			DOUBLE PRECISION tauae(kdlon, kflev, 2)
58			DOUBLE PRECISION pizae(kdlon, kflev, 2)
59			DOUBLE PRECISION cgae(kdlon, kflev, 2)
60			DOUBLE PRECISION paki(kdlon, 2)
61			DOUBLE PRECISION palbd(kdlon, 2)
62			DOUBLE PRECISION palbp(kdlon, 2)
63			DOUBLE PRECISION pcg(kdlon, 2, kflev)
64			DOUBLE PRECISION pcld(kdlon, kflev)
65			DOUBLE PRECISION pclear(kdlon)
66			DOUBLE PRECISION pdsig(kdlon, kflev)
67			DOUBLE PRECISION pomega(kdlon, 2, kflev)
68			DOUBLE PRECISION poz(kdlon, kflev)
69			DOUBLE PRECISION pqs(kdlon, kflev)
70			DOUBLE PRECISION prmu(kdlon)
71			DOUBLE PRECISION psec(kdlon)
72			DOUBLE PRECISION ptau(kdlon, 2, kflev)
73			DOUBLE PRECISION pud(kdlon, 5, kflev+1)
74			DOUBLE PRECISION pwv(kdlon, kflev)
75	guez	81
76	guez	178	DOUBLE PRECISION pfdown(kdlon, kflev+1)
77			DOUBLE PRECISION pfup(kdlon, kflev+1)
78	guez	81
79	guez	178	! * LOCAL VARIABLES:
80	guez	81
81	guez	178	INTEGER iind2(2), iind3(3)
82			DOUBLE PRECISION zcgaz(kdlon, kflev)
83			DOUBLE PRECISION zfd(kdlon, kflev+1)
84			DOUBLE PRECISION zfu(kdlon, kflev+1)
85			DOUBLE PRECISION zg(kdlon)
86			DOUBLE PRECISION zgg(kdlon)
87			DOUBLE PRECISION zpizaz(kdlon, kflev)
88			DOUBLE PRECISION zrayl(kdlon)
89			DOUBLE PRECISION zray1(kdlon, kflev+1)
90			DOUBLE PRECISION zray2(kdlon, kflev+1)
91			DOUBLE PRECISION zref(kdlon)
92			DOUBLE PRECISION zrefz(kdlon, 2, kflev+1)
93			DOUBLE PRECISION zre1(kdlon)
94			DOUBLE PRECISION zre2(kdlon)
95			DOUBLE PRECISION zrj(kdlon, 6, kflev+1)
96			DOUBLE PRECISION zrj0(kdlon, 6, kflev+1)
97			DOUBLE PRECISION zrk(kdlon, 6, kflev+1)
98			DOUBLE PRECISION zrk0(kdlon, 6, kflev+1)
99			DOUBLE PRECISION zrl(kdlon, 8)
100			DOUBLE PRECISION zrmue(kdlon, kflev+1)
101			DOUBLE PRECISION zrmu0(kdlon, kflev+1)
102			DOUBLE PRECISION zrmuz(kdlon)
103			DOUBLE PRECISION zrneb(kdlon)
104			DOUBLE PRECISION zr1(kdlon)
105			DOUBLE PRECISION zr2(kdlon, 2)
106			DOUBLE PRECISION zr3(kdlon, 3)
107			DOUBLE PRECISION zr4(kdlon)
108			DOUBLE PRECISION zr21(kdlon)
109			DOUBLE PRECISION zr22(kdlon)
110			DOUBLE PRECISION zs(kdlon)
111			DOUBLE PRECISION ztauaz(kdlon, kflev)
112			DOUBLE PRECISION zto1(kdlon)
113			DOUBLE PRECISION ztr(kdlon, 2, kflev+1)
114			DOUBLE PRECISION ztra1(kdlon, kflev+1)
115			DOUBLE PRECISION ztra2(kdlon, kflev+1)
116			DOUBLE PRECISION ztr1(kdlon)
117			DOUBLE PRECISION ztr2(kdlon)
118			DOUBLE PRECISION zw(kdlon)
119			DOUBLE PRECISION zw1(kdlon)
120			DOUBLE PRECISION zw2(kdlon, 2)
121			DOUBLE PRECISION zw3(kdlon, 3)
122			DOUBLE PRECISION zw4(kdlon)
123			DOUBLE PRECISION zw5(kdlon)
124	guez	81
125	guez	178	INTEGER jl, jk, k, jaj, ikm1, ikl, jn, jabs, jkm1
126			INTEGER jref, jkl, jklp1, jajp, jkki, jkkp4, jn2j, iabs
127			DOUBLE PRECISION zrmum1, zwh2o, zcneb, zaa, zbb, zrki, zre11
128	guez	81
129	guez	178	! * Prescribed Data:
130	guez	81
131	guez	178	DOUBLE PRECISION rsun(2)
132			SAVE rsun
133			DOUBLE PRECISION rray(2, 6)
134			SAVE rray
135			DATA rsun(1)/0.441676d0/
136			DATA rsun(2)/0.558324d0/
137			DATA (rray(1,k), k=1, 6)/.428937d-01, .890743d+00, -.288555d+01, &
138			.522744d+01, -.469173d+01, .161645d+01/
139			DATA (rray(2,k), k=1, 6)/.697200d-02, .173297d-01, -.850903d-01, &
140			.248261d+00, -.302031d+00, .129662d+00/
141	guez	81
142	guez	178	! ------------------------------------------------------------------
143	guez	81
144	guez	178	! * 1. SECOND SPECTRAL INTERVAL (0.68-4.00 MICRON)
145			! -------------------------------------------
146	guez	81
147
148
149	guez	178	! * 1.1 OPTICAL THICKNESS FOR RAYLEIGH SCATTERING
150			! -----------------------------------------
151	guez	81
152
153	guez	178	DO jl = 1, kdlon
154			zrmum1 = 1. - prmu(jl)
155			zrayl(jl) = rray(knu, 1) + zrmum1(rray(knu,2)+zrmum1(rray(knu, &
156			3)+zrmum1(rray(knu,4)+zrmum1(rray(knu,5)+zrmum1*rray(knu,6)))))
157			END DO
158	guez	81
159
160	guez	178	! ------------------------------------------------------------------
161	guez	81
162	guez	178	! * 2. CONTINUUM SCATTERING CALCULATIONS
163			! ---------------------------------
164	guez	81
165
166	guez	178	! * 2.1 CLEAR-SKY FRACTION OF THE COLUMN
167			! --------------------------------
168	guez	81
169
170	guez	178	CALL swclr(knu, flag_aer, tauae, pizae, cgae, palbp, pdsig, zrayl, &
171			psec, zcgaz, zpizaz, zray1, zray2, zrefz, zrj0, zrk0, zrmu0, ztauaz, &
172			ztra1, ztra2)
173	guez	81
174
175	guez	178	! * 2.2 CLOUDY FRACTION OF THE COLUMN
176			! -----------------------------
177	guez	81
178
179	guez	178	CALL swr(knu, palbd, pcg, pcld, pomega, psec, ptau, zcgaz, &
180			zpizaz, zray1, zray2, zrefz, zrj, zrk, zrmue, ztauaz, ztra1, ztra2)
181	guez	81
182
183	guez	178	! ------------------------------------------------------------------
184	guez	81
185	guez	178	! * 3. SCATTERING CALCULATIONS WITH GREY MOLECULAR ABSORPTION
186			! ------------------------------------------------------
187	guez	81
188
189	guez	178	jn = 2
190	guez	81
191	guez	178	DO jabs = 1, 2
192	guez	81
193
194	guez	178	! * 3.1 SURFACE CONDITIONS
195			! ------------------
196	guez	81
197
198	guez	178	DO jl = 1, kdlon
199			zrefz(jl, 2, 1) = palbd(jl, knu)
200			zrefz(jl, 1, 1) = palbd(jl, knu)
201			END DO
202	guez	81
203
204	guez	178	! * 3.2 INTRODUCING CLOUD EFFECTS
205			! -------------------------
206	guez	81
207
208	guez	178	DO jk = 2, kflev + 1
209			jkm1 = jk - 1
210			ikl = kflev + 1 - jkm1
211			DO jl = 1, kdlon
212			zrneb(jl) = pcld(jl, jkm1)
213			IF (jabs==1 .AND. zrneb(jl)>2.*zeelog) THEN
214			zwh2o = max(pwv(jl,jkm1), zeelog)
215			zcneb = max(zeelog, min(zrneb(jl),1.-zeelog))
216			zbb = pud(jl, jabs, jkm1)*pqs(jl, jkm1)/zwh2o
217			zaa = max((pud(jl,jabs,jkm1)-zcneb*zbb)/(1.-zcneb), zeelog)
218			ELSE
219			zaa = pud(jl, jabs, jkm1)
220			zbb = zaa
221			END IF
222			zrki = paki(jl, jabs)
223			zs(jl) = exp(-zrkizaa1.66)
224			zg(jl) = exp(-zrki*zaa/zrmue(jl,jk))
225			ztr1(jl) = 0.
226			zre1(jl) = 0.
227			ztr2(jl) = 0.
228			zre2(jl) = 0.
229	guez	81
230	guez	178	zw(jl) = pomega(jl, knu, jkm1)
231			zto1(jl) = ptau(jl, knu, jkm1)/zw(jl) + ztauaz(jl, jkm1)/zpizaz(jl, &
232			jkm1) + zbb*zrki
233	guez	81
234	guez	178	zr21(jl) = ptau(jl, knu, jkm1) + ztauaz(jl, jkm1)
235			zr22(jl) = ptau(jl, knu, jkm1)/zr21(jl)
236			zgg(jl) = zr22(jl)pcg(jl, knu, jkm1) + (1.-zr22(jl))zcgaz(jl, jkm1)
237			zw(jl) = zr21(jl)/zto1(jl)
238			zref(jl) = zrefz(jl, 1, jkm1)
239			zrmuz(jl) = zrmue(jl, jk)
240			END DO
241	guez	81
242	guez	178	CALL swde(zgg, zref, zrmuz, zto1, zw, zre1, zre2, ztr1, ztr2)
243	guez	81
244	guez	178	DO jl = 1, kdlon
245	guez	81
246	guez	178	zrefz(jl, 2, jk) = (1.-zrneb(jl))(zray1(jl,jkm1)+zrefz(jl,2,jkm1) &
247			ztra1(jl,jkm1)ztra2(jl,jkm1))zg(jl)zs(jl) + zrneb(jl)zre1(jl)
248	guez	81
249	guez	178	ztr(jl, 2, jkm1) = zrneb(jl)ztr1(jl) + (ztra1(jl,jkm1))zg(jl)*(1.- &
250			zrneb(jl))
251	guez	81
252	guez	178	zrefz(jl, 1, jk) = (1.-zrneb(jl))(zray1(jl,jkm1)+zrefz(jl,1,jkm1) &
253			ztra1(jl,jkm1)ztra2(jl,jkm1)/(1.-zray2(jl,jkm1)zrefz(jl,1, &
254			jkm1)))zg(jl)zs(jl) + zrneb(jl)*zre2(jl)
255	guez	81
256	guez	178	ztr(jl, 1, jkm1) = zrneb(jl)*ztr2(jl) + (ztra1(jl,jkm1)/(1.-zray2(jl, &
257			jkm1)zrefz(jl,1,jkm1)))zg(jl)*(1.-zrneb(jl))
258	guez	81
259	guez	178	END DO
260			END DO
261	guez	81
262	guez	178	! * 3.3 REFLECT./TRANSMISSIVITY BETWEEN SURFACE AND LEVEL
263			! -------------------------------------------------
264	guez	81
265
266	guez	178	DO jref = 1, 2
267	guez	81
268	guez	178	jn = jn + 1
269
270			DO jl = 1, kdlon
271			zrj(jl, jn, kflev+1) = 1.
272			zrk(jl, jn, kflev+1) = zrefz(jl, jref, kflev+1)
273			END DO
274
275			DO jk = 1, kflev
276			jkl = kflev + 1 - jk
277			jklp1 = jkl + 1
278			DO jl = 1, kdlon
279			zre11 = zrj(jl, jn, jklp1)*ztr(jl, jref, jkl)
280			zrj(jl, jn, jkl) = zre11
281			zrk(jl, jn, jkl) = zre11*zrefz(jl, jref, jkl)
282			END DO
283			END DO
284			END DO
285	guez	81	END DO
286
287
288	guez	178	! ------------------------------------------------------------------
289	guez	81
290	guez	178	! * 4. INVERT GREY AND CONTINUUM FLUXES
291			! --------------------------------
292	guez	81
293
294
295	guez	178	! * 4.1 UPWARD (ZRK) AND DOWNWARD (ZRJ) PSEUDO-FLUXES
296			! ---------------------------------------------
297	guez	81
298
299	guez	178	DO jk = 1, kflev + 1
300			DO jaj = 1, 5, 2
301			jajp = jaj + 1
302			DO jl = 1, kdlon
303			zrj(jl, jaj, jk) = zrj(jl, jaj, jk) - zrj(jl, jajp, jk)
304			zrk(jl, jaj, jk) = zrk(jl, jaj, jk) - zrk(jl, jajp, jk)
305			zrj(jl, jaj, jk) = max(zrj(jl,jaj,jk), zeelog)
306			zrk(jl, jaj, jk) = max(zrk(jl,jaj,jk), zeelog)
307			END DO
308			END DO
309	guez	81	END DO
310
311	guez	178	DO jk = 1, kflev + 1
312			DO jaj = 2, 6, 2
313			DO jl = 1, kdlon
314			zrj(jl, jaj, jk) = max(zrj(jl,jaj,jk), zeelog)
315			zrk(jl, jaj, jk) = max(zrk(jl,jaj,jk), zeelog)
316			END DO
317			END DO
318	guez	81	END DO
319
320	guez	178	! * 4.2 EFFECTIVE ABSORBER AMOUNTS BY INVERSE LAPLACE
321			! ---------------------------------------------
322	guez	81
323
324	guez	178	DO jk = 1, kflev + 1
325			jkki = 1
326			DO jaj = 1, 2
327			iind2(1) = jaj
328			iind2(2) = jaj
329			DO jn = 1, 2
330			jn2j = jn + 2*jaj
331			jkkp4 = jkki + 4
332	guez	81
333	guez	178	! * 4.2.1 EFFECTIVE ABSORBER AMOUNTS
334			! --------------------------
335	guez	81
336
337	guez	178	DO jl = 1, kdlon
338			zw2(jl, 1) = log(zrj(jl,jn,jk)/zrj(jl,jn2j,jk))/paki(jl, jaj)
339			zw2(jl, 2) = log(zrk(jl,jn,jk)/zrk(jl,jn2j,jk))/paki(jl, jaj)
340			END DO
341	guez	81
342	guez	178	! * 4.2.2 TRANSMISSION FUNCTION
343			! ---------------------
344	guez	81
345
346	guez	178	CALL swtt1(knu, 2, iind2, zw2, zr2)
347	guez	81
348	guez	178	DO jl = 1, kdlon
349			zrl(jl, jkki) = zr2(jl, 1)
350			zrl(jl, jkkp4) = zr2(jl, 2)
351			END DO
352	guez	81
353	guez	178	jkki = jkki + 1
354			END DO
355			END DO
356	guez	81
357	guez	178	! * 4.3 UPWARD AND DOWNWARD FLUXES WITH H2O AND UMG ABSORPTION
358			! ------------------------------------------------------
359	guez	81
360
361	guez	178	DO jl = 1, kdlon
362			pfdown(jl, jk) = zrj(jl, 1, jk)zrl(jl, 1)zrl(jl, 3) + &
363			zrj(jl, 2, jk)zrl(jl, 2)zrl(jl, 4)
364			pfup(jl, jk) = zrk(jl, 1, jk)zrl(jl, 5)zrl(jl, 7) + &
365			zrk(jl, 2, jk)zrl(jl, 6)zrl(jl, 8)
366			END DO
367	guez	81	END DO
368
369
370	guez	178	! ------------------------------------------------------------------
371	guez	81
372	guez	178	! * 5. MOLECULAR ABSORPTION ON CLEAR-SKY FLUXES
373			! ----------------------------------------
374	guez	81
375
376
377	guez	178	! * 5.1 DOWNWARD FLUXES
378			! ---------------
379	guez	81
380
381	guez	178	jaj = 2
382			iind3(1) = 1
383			iind3(2) = 2
384			iind3(3) = 3
385	guez	81
386			DO jl = 1, kdlon
387	guez	178	zw3(jl, 1) = 0.
388			zw3(jl, 2) = 0.
389			zw3(jl, 3) = 0.
390			zw4(jl) = 0.
391			zw5(jl) = 0.
392			zr4(jl) = 1.
393			zfd(jl, kflev+1) = zrj0(jl, jaj, kflev+1)
394	guez	81	END DO
395	guez	178	DO jk = 1, kflev
396			ikl = kflev + 1 - jk
397			DO jl = 1, kdlon
398			zw3(jl, 1) = zw3(jl, 1) + pud(jl, 1, ikl)/zrmu0(jl, ikl)
399			zw3(jl, 2) = zw3(jl, 2) + pud(jl, 2, ikl)/zrmu0(jl, ikl)
400			zw3(jl, 3) = zw3(jl, 3) + poz(jl, ikl)/zrmu0(jl, ikl)
401			zw4(jl) = zw4(jl) + pud(jl, 4, ikl)/zrmu0(jl, ikl)
402			zw5(jl) = zw5(jl) + pud(jl, 5, ikl)/zrmu0(jl, ikl)
403			END DO
404	guez	81
405	guez	178	CALL swtt1(knu, 3, iind3, zw3, zr3)
406	guez	81
407	guez	178	DO jl = 1, kdlon
408			! ZR4(JL) = EXP(-RSWCEZW4(JL)-RSWCPZW5(JL))
409			zfd(jl, ikl) = zr3(jl, 1)zr3(jl, 2)zr3(jl, 3)zr4(jl) &
410			zrj0(jl, jaj, ikl)
411			END DO
412	guez	81	END DO
413
414
415	guez	178	! * 5.2 UPWARD FLUXES
416			! -------------
417	guez	81
418
419			DO jl = 1, kdlon
420	guez	178	zfu(jl, 1) = zfd(jl, 1)*palbp(jl, knu)
421	guez	81	END DO
422
423	guez	178	DO jk = 2, kflev + 1
424			ikm1 = jk - 1
425			DO jl = 1, kdlon
426			zw3(jl, 1) = zw3(jl, 1) + pud(jl, 1, ikm1)*1.66
427			zw3(jl, 2) = zw3(jl, 2) + pud(jl, 2, ikm1)*1.66
428			zw3(jl, 3) = zw3(jl, 3) + poz(jl, ikm1)*1.66
429			zw4(jl) = zw4(jl) + pud(jl, 4, ikm1)*1.66
430			zw5(jl) = zw5(jl) + pud(jl, 5, ikm1)*1.66
431			END DO
432	guez	81
433	guez	178	CALL swtt1(knu, 3, iind3, zw3, zr3)
434
435			DO jl = 1, kdlon
436			! ZR4(JL) = EXP(-RSWCEZW4(JL)-RSWCPZW5(JL))
437			zfu(jl, jk) = zr3(jl, 1)zr3(jl, 2)zr3(jl, 3)zr4(jl) &
438			zrk0(jl, jaj, jk)
439			END DO
440	guez	81	END DO
441
442
443	guez	178	! ------------------------------------------------------------------
444	guez	81
445	guez	178	! * 6. INTRODUCTION OF OZONE AND H2O CONTINUUM ABSORPTION
446			! --------------------------------------------------
447	guez	81
448	guez	178	iabs = 3
449	guez	81
450	guez	178	! * 6.1 DOWNWARD FLUXES
451			! ---------------
452	guez	81
453			DO jl = 1, kdlon
454	guez	178	zw1(jl) = 0.
455			zw4(jl) = 0.
456			zw5(jl) = 0.
457			zr1(jl) = 0.
458			pfdown(jl, kflev+1) = ((1.-pclear(jl))pfdown(jl,kflev+1)+pclear(jl)zfd( &
459			jl,kflev+1))*rsun(knu)
460	guez	81	END DO
461
462	guez	178	DO jk = 1, kflev
463			ikl = kflev + 1 - jk
464			DO jl = 1, kdlon
465			zw1(jl) = zw1(jl) + poz(jl, ikl)/zrmue(jl, ikl)
466			zw4(jl) = zw4(jl) + pud(jl, 4, ikl)/zrmue(jl, ikl)
467			zw5(jl) = zw5(jl) + pud(jl, 5, ikl)/zrmue(jl, ikl)
468			! ZR4(JL) = EXP(-RSWCEZW4(JL)-RSWCPZW5(JL))
469			END DO
470	guez	81
471	guez	178	CALL swtt(knu, iabs, zw1, zr1)
472
473			DO jl = 1, kdlon
474			pfdown(jl, ikl) = ((1.-pclear(jl))zr1(jl)zr4(jl)*pfdown(jl,ikl)+ &
475			pclear(jl)zfd(jl,ikl))rsun(knu)
476			END DO
477	guez	81	END DO
478
479
480	guez	178	! * 6.2 UPWARD FLUXES
481			! -------------
482	guez	81
483			DO jl = 1, kdlon
484	guez	178	pfup(jl, 1) = ((1.-pclear(jl))zr1(jl)zr4(jl)pfup(jl,1)+pclear(jl)zfu( &
485			jl,1))*rsun(knu)
486	guez	81	END DO
487
488	guez	178	DO jk = 2, kflev + 1
489			ikm1 = jk - 1
490			DO jl = 1, kdlon
491			zw1(jl) = zw1(jl) + poz(jl, ikm1)*1.66
492			zw4(jl) = zw4(jl) + pud(jl, 4, ikm1)*1.66
493			zw5(jl) = zw5(jl) + pud(jl, 5, ikm1)*1.66
494			! ZR4(JL) = EXP(-RSWCEZW4(JL)-RSWCPZW5(JL))
495			END DO
496	guez	81
497	guez	178	CALL swtt(knu, iabs, zw1, zr1)
498
499			DO jl = 1, kdlon
500			pfup(jl, jk) = ((1.-pclear(jl))zr1(jl)zr4(jl)pfup(jl,jk)+pclear(jl) &
501			zfu(jl,jk))*rsun(knu)
502			END DO
503	guez	81	END DO
504
505	guez	178	END SUBROUTINE sw2s
506	guez	81
507	guez	178	end module sw2s_m