phylmd/Radlwsw/sw.f

module sw_m

  IMPLICIT none

contains

  SUBROUTINE SW(PSCT, PRMU0, PFRAC, PPMB, PDP, PPSOL, PALBD, PALBP, PTAVE, &
       PWV, PQS, POZON, PCLDSW, PTAU, POMEGA, PCG, PHEAT, PHEAT0, PALBPLA, &
       PTOPSW, PSOLSW, PTOPSW0, PSOLSW0, ZFSUP, ZFSDN, ZFSUP0, ZFSDN0, &
       PTOPSWAD, PSOLSWAD, PTOPSWAI, PSOLSWAI, ok_ade)

    ! Purpose.
    ! This routine computes the shortwave radiation fluxes in two
    ! spectral intervals following Fouquart and Bonnel (1980).

    ! Method.
    ! 1. Computes absorber amounts (swu)
    ! 2. Computes fluxes in 1st spectral interval (SW1S)
    ! 3. Computes fluxes in 2nd spectral interval (SW2S)

    ! Reference.
    ! See radiation part of the ECMWF research department
    ! documentation, and Fouquart and Bonnel (1980)

    ! Author.
    ! Jean-Jacques Morcrette *ecmwf*

    ! Modifications.
    ! Original: 89-07-14
    ! 95-01-01 J.-J. Morcrette direct/diffuse albedo
    ! 03-11-27 J. Quaas Introduce aerosol forcings (based on Boucher)

    USE raddim, ONLY: kdlon, kflev
    USE suphec_m, ONLY: rcpd, rday, rg
    use sw1s_m, only: sw1s
    use sw2s_m, only: sw2s

    ! ARGUMENTS:

    DOUBLE PRECISION PSCT ! constante solaire (valeur conseillee: 1370)
    DOUBLE PRECISION PRMU0(KDLON) ! COSINE OF ZENITHAL ANGLE
    DOUBLE PRECISION PFRAC(KDLON) ! fraction de la journee
    DOUBLE PRECISION PPMB(KDLON, KFLEV+1) ! HALF-LEVEL PRESSURE (MB)
    DOUBLE PRECISION PDP(KDLON, KFLEV) ! LAYER THICKNESS (PA)
    DOUBLE PRECISION PPSOL(KDLON) ! SURFACE PRESSURE (PA)
    DOUBLE PRECISION PALBD(KDLON, 2) ! albedo du sol (lumiere diffuse)
    DOUBLE PRECISION PALBP(KDLON, 2) ! albedo du sol (lumiere parallele)
    DOUBLE PRECISION PTAVE(KDLON, KFLEV) ! LAYER TEMPERATURE (K)
    DOUBLE PRECISION PWV(KDLON, KFLEV) ! SPECIFIC HUMIDITY (KG/KG)
    DOUBLE PRECISION PQS(KDLON, KFLEV) ! SATURATED WATER VAPOUR (KG/KG)
    DOUBLE PRECISION POZON(KDLON, KFLEV) ! OZONE CONCENTRATION (KG/KG)
    DOUBLE PRECISION PCLDSW(KDLON, KFLEV) ! CLOUD FRACTION
    DOUBLE PRECISION PTAU(KDLON, 2, KFLEV) ! CLOUD OPTICAL THICKNESS
    DOUBLE PRECISION POMEGA(KDLON, 2, KFLEV) ! SINGLE SCATTERING ALBEDO
    DOUBLE PRECISION PCG(KDLON, 2, KFLEV) ! ASYMETRY FACTOR
    DOUBLE PRECISION PHEAT(KDLON, KFLEV) ! SHORTWAVE HEATING (K/DAY)
    DOUBLE PRECISION PHEAT0(KDLON, KFLEV)! SHORTWAVE HEATING (K/DAY) clear-sky
    DOUBLE PRECISION PALBPLA(KDLON) ! PLANETARY ALBEDO
    DOUBLE PRECISION PTOPSW(KDLON) ! SHORTWAVE FLUX AT T.O.A.
    DOUBLE PRECISION PSOLSW(KDLON) ! SHORTWAVE FLUX AT SURFACE
    DOUBLE PRECISION PTOPSW0(KDLON) ! SHORTWAVE FLUX AT T.O.A. (CLEAR-SKY)
    DOUBLE PRECISION PSOLSW0(KDLON) ! SHORTWAVE FLUX AT SURFACE (CLEAR-SKY)
    DOUBLE PRECISION ZFSUP(KDLON, KFLEV+1)
    DOUBLE PRECISION ZFSDN(KDLON, KFLEV+1)
    DOUBLE PRECISION ZFSUP0(KDLON, KFLEV+1)
    DOUBLE PRECISION ZFSDN0(KDLON, KFLEV+1)

    DOUBLE PRECISION, intent(out):: PTOPSWAD(KDLON) 
    ! (diagnosed aerosol forcing)SHORTWAVE FLUX AT T.O.A.(+AEROSOL DIR)

    DOUBLE PRECISION, intent(out):: PSOLSWAD(KDLON) 
    ! (diagnosed aerosol forcing)SHORTWAVE FLUX AT SURFACE(+AEROSOL DIR)

    DOUBLE PRECISION, intent(out):: PTOPSWAI(KDLON) 
    ! (diagnosed aerosol forcing)SHORTWAVE FLUX AT T.O.A.(+AEROSOL IND)

    DOUBLE PRECISION, intent(out):: PSOLSWAI(KDLON) 
    ! (diagnosed aerosol forcing)SHORTWAVE FLUX AT SURFACE(+AEROSOL IND)

    logical, intent(in):: ok_ade ! use aerosol forcings or not?

    ! Local:

    DOUBLE PRECISION ZOZ(KDLON, KFLEV)
    DOUBLE PRECISION ZAKI(KDLON, 2) 
    DOUBLE PRECISION ZCLD(KDLON, KFLEV)
    DOUBLE PRECISION ZCLEAR(KDLON) 
    DOUBLE PRECISION ZDSIG(KDLON, KFLEV)
    DOUBLE PRECISION ZFACT(KDLON)
    DOUBLE PRECISION ZFD(KDLON, KFLEV+1)
    DOUBLE PRECISION ZFDOWN(KDLON, KFLEV+1)
    DOUBLE PRECISION ZFU(KDLON, KFLEV+1)
    DOUBLE PRECISION ZFUP(KDLON, KFLEV+1)
    DOUBLE PRECISION ZRMU(KDLON)
    DOUBLE PRECISION ZSEC(KDLON)
    DOUBLE PRECISION ZUD(KDLON, 5, KFLEV+1)
    DOUBLE PRECISION ZCLDSW0(KDLON, KFLEV)

    INTEGER inu, jl, jk, i, k, kpl1

    INTEGER, PARAMETER:: swpas = 1 ! Every swpas steps, sw is calculated

    INTEGER:: itapsw = 0
    LOGICAL:: appel1er = .TRUE.
    !jq-Introduced for aerosol forcings
    logical, save:: flag_aer

    !jq - Fluxes including aerosol effects
    DOUBLE PRECISION, save:: ZFSUPAD(KDLON, KFLEV+1)
    DOUBLE PRECISION, save:: ZFSDNAD(KDLON, KFLEV+1)
    DOUBLE PRECISION, save:: ZFSUPAI(KDLON, KFLEV+1)
    DOUBLE PRECISION, save:: ZFSDNAI(KDLON, KFLEV+1)

    logical:: initialized = .false.
    REAL, PARAMETER :: dobson_u = 2.1415E-05 ! Dobson unit, in kg m-2

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

    if(.not.initialized) then
       flag_aer=.false.
       initialized=.TRUE.
       ZFSUPAD = 0.
       ZFSDNAD = 0.
       ZFSUPAI = 0.
       ZFSDNAI = 0.
    endif
    !rv

    IF (appel1er) THEN
       PRINT*, 'SW calling frequency: ', swpas
       PRINT*, " In general, it should be 1"
       appel1er = .FALSE.
    ENDIF

    IF (MOD(itapsw, swpas).EQ.0) THEN
       DO JK = 1 , KFLEV
          DO JL = 1, KDLON
             ZCLDSW0(JL, JK) = 0.0
             ZOZ(JL, JK) = POZON(JL, JK) / (dobson_u * 1E3 * rg) * PDP(JL, JK)
          ENDDO
       ENDDO

       ! clear-sky:
       CALL SWU(PSCT, ZCLDSW0, PPMB, PPSOL, &
            PRMU0, PFRAC, PTAVE, PWV, &
            ZAKI, ZCLD, ZCLEAR, ZDSIG, ZFACT, ZRMU, ZSEC, ZUD)
       INU = 1
       CALL SW1S(INU, flag_aer, PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZDSIG, &
            POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD, ZFD, ZFU)
       INU = 2
       CALL SW2S(INU, flag_aer, ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZDSIG, &
            POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD, PWV, PQS, ZFDOWN, ZFUP)
       DO JK = 1 , KFLEV+1
          DO JL = 1, KDLON
             ZFSUP0(JL, JK) = (ZFUP(JL, JK) + ZFU(JL, JK)) * ZFACT(JL)
             ZFSDN0(JL, JK) = (ZFDOWN(JL, JK) + ZFD(JL, JK)) * ZFACT(JL)
          ENDDO
       ENDDO

       flag_aer= .false.
       CALL SWU(PSCT, PCLDSW, PPMB, PPSOL, &
            PRMU0, PFRAC, PTAVE, PWV, &
            ZAKI, ZCLD, ZCLEAR, ZDSIG, ZFACT, ZRMU, ZSEC, ZUD)
       INU = 1
       CALL SW1S(INU, .false., PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZDSIG, &
            POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD, ZFD, ZFU)
       INU = 2
       CALL SW2S(INU, .false., ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZDSIG, &
            POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD, PWV, PQS, ZFDOWN, ZFUP)

       ! cloudy-sky:

       DO JK = 1 , KFLEV+1
          DO JL = 1, KDLON
             ZFSUP(JL, JK) = (ZFUP(JL, JK) + ZFU(JL, JK)) * ZFACT(JL)
             ZFSDN(JL, JK) = (ZFDOWN(JL, JK) + ZFD(JL, JK)) * ZFACT(JL)
          ENDDO
       ENDDO

       IF (ok_ade) THEN
          ! cloudy-sky + aerosol dir OB
          flag_aer= .true.
          CALL SWU(PSCT, PCLDSW, PPMB, PPSOL, PRMU0, PFRAC, PTAVE, PWV, ZAKI, &
               ZCLD, ZCLEAR, ZDSIG, ZFACT, ZRMU, ZSEC, ZUD)
          INU = 1
          CALL SW1S(INU, .true., PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZDSIG, &
               POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD, ZFD, ZFU)
          INU = 2
          CALL SW2S(INU, .true., ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, &
               ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD, PWV, PQS, ZFDOWN, &
               ZFUP)
          DO JK = 1 , KFLEV+1
             DO JL = 1, KDLON
                ZFSUPAD(JL, JK) = ZFSUP(JL, JK) 
                ZFSDNAD(JL, JK) = ZFSDN(JL, JK) 
                ZFSUP(JL, JK) = (ZFUP(JL, JK) + ZFU(JL, JK)) * ZFACT(JL)
                ZFSDN(JL, JK) = (ZFDOWN(JL, JK) + ZFD(JL, JK)) * ZFACT(JL)
             ENDDO
          ENDDO
       ENDIF

       itapsw = 0
    ENDIF
    itapsw = itapsw + 1

    DO k = 1, KFLEV
       kpl1 = k+1
       DO i = 1, KDLON
          PHEAT(i, k) = -(ZFSUP(i, kpl1)-ZFSUP(i, k)) &
               -(ZFSDN(i, k)-ZFSDN(i, kpl1))
          PHEAT(i, k) = PHEAT(i, k) * RDAY*RG/RCPD / PDP(i, k)
          PHEAT0(i, k) = -(ZFSUP0(i, kpl1)-ZFSUP0(i, k)) &
               -(ZFSDN0(i, k)-ZFSDN0(i, kpl1))
          PHEAT0(i, k) = PHEAT0(i, k) * RDAY*RG/RCPD / PDP(i, k)
       ENDDO
    ENDDO
    DO i = 1, KDLON
       PALBPLA(i) = ZFSUP(i, KFLEV+1)/(ZFSDN(i, KFLEV+1)+1.0e-20)

       PSOLSW(i) = ZFSDN(i, 1) - ZFSUP(i, 1)
       PTOPSW(i) = ZFSDN(i, KFLEV+1) - ZFSUP(i, KFLEV+1)

       PSOLSW0(i) = ZFSDN0(i, 1) - ZFSUP0(i, 1)
       PTOPSW0(i) = ZFSDN0(i, KFLEV+1) - ZFSUP0(i, KFLEV+1)

       PSOLSWAD(i) = ZFSDNAD(i, 1) - ZFSUPAD(i, 1)
       PTOPSWAD(i) = ZFSDNAD(i, KFLEV+1) - ZFSUPAD(i, KFLEV+1)

       PSOLSWAI(i) = ZFSDNAI(i, 1) - ZFSUPAI(i, 1)
       PTOPSWAI(i) = ZFSDNAI(i, KFLEV+1) - ZFSUPAI(i, KFLEV+1)
    ENDDO

  END SUBROUTINE SW

end module sw_m
1	module sw_m
2
3	IMPLICIT none
4
5	contains
6
7	SUBROUTINE SW(PSCT, PRMU0, PFRAC, PPMB, PDP, PPSOL, PALBD, PALBP, PTAVE, &
8	PWV, PQS, POZON, PCLDSW, PTAU, POMEGA, PCG, PHEAT, PHEAT0, PALBPLA, &
9	PTOPSW, PSOLSW, PTOPSW0, PSOLSW0, ZFSUP, ZFSDN, ZFSUP0, ZFSDN0, &
10	PTOPSWAD, PSOLSWAD, PTOPSWAI, PSOLSWAI, ok_ade)
11
12	! Purpose.
13	! This routine computes the shortwave radiation fluxes in two
14	! spectral intervals following Fouquart and Bonnel (1980).
15
16	! Method.
17	! 1. Computes absorber amounts (swu)
18	! 2. Computes fluxes in 1st spectral interval (SW1S)
19	! 3. Computes fluxes in 2nd spectral interval (SW2S)
20
21	! Reference.
22	! See radiation part of the ECMWF research department
23	! documentation, and Fouquart and Bonnel (1980)
24
25	! Author.
26	! Jean-Jacques Morcrette ecmwf
27
28	! Modifications.
29	! Original: 89-07-14
30	! 95-01-01 J.-J. Morcrette direct/diffuse albedo
31	! 03-11-27 J. Quaas Introduce aerosol forcings (based on Boucher)
32
33	USE raddim, ONLY: kdlon, kflev
34	USE suphec_m, ONLY: rcpd, rday, rg
35	use sw1s_m, only: sw1s
36	use sw2s_m, only: sw2s
37
38	! ARGUMENTS:
39
40	DOUBLE PRECISION PSCT ! constante solaire (valeur conseillee: 1370)
41	DOUBLE PRECISION PRMU0(KDLON) ! COSINE OF ZENITHAL ANGLE
42	DOUBLE PRECISION PFRAC(KDLON) ! fraction de la journee
43	DOUBLE PRECISION PPMB(KDLON, KFLEV+1) ! HALF-LEVEL PRESSURE (MB)
44	DOUBLE PRECISION PDP(KDLON, KFLEV) ! LAYER THICKNESS (PA)
45	DOUBLE PRECISION PPSOL(KDLON) ! SURFACE PRESSURE (PA)
46	DOUBLE PRECISION PALBD(KDLON, 2) ! albedo du sol (lumiere diffuse)
47	DOUBLE PRECISION PALBP(KDLON, 2) ! albedo du sol (lumiere parallele)
48	DOUBLE PRECISION PTAVE(KDLON, KFLEV) ! LAYER TEMPERATURE (K)
49	DOUBLE PRECISION PWV(KDLON, KFLEV) ! SPECIFIC HUMIDITY (KG/KG)
50	DOUBLE PRECISION PQS(KDLON, KFLEV) ! SATURATED WATER VAPOUR (KG/KG)
51	DOUBLE PRECISION POZON(KDLON, KFLEV) ! OZONE CONCENTRATION (KG/KG)
52	DOUBLE PRECISION PCLDSW(KDLON, KFLEV) ! CLOUD FRACTION
53	DOUBLE PRECISION PTAU(KDLON, 2, KFLEV) ! CLOUD OPTICAL THICKNESS
54	DOUBLE PRECISION POMEGA(KDLON, 2, KFLEV) ! SINGLE SCATTERING ALBEDO
55	DOUBLE PRECISION PCG(KDLON, 2, KFLEV) ! ASYMETRY FACTOR
56	DOUBLE PRECISION PHEAT(KDLON, KFLEV) ! SHORTWAVE HEATING (K/DAY)
57	DOUBLE PRECISION PHEAT0(KDLON, KFLEV)! SHORTWAVE HEATING (K/DAY) clear-sky
58	DOUBLE PRECISION PALBPLA(KDLON) ! PLANETARY ALBEDO
59	DOUBLE PRECISION PTOPSW(KDLON) ! SHORTWAVE FLUX AT T.O.A.
60	DOUBLE PRECISION PSOLSW(KDLON) ! SHORTWAVE FLUX AT SURFACE
61	DOUBLE PRECISION PTOPSW0(KDLON) ! SHORTWAVE FLUX AT T.O.A. (CLEAR-SKY)
62	DOUBLE PRECISION PSOLSW0(KDLON) ! SHORTWAVE FLUX AT SURFACE (CLEAR-SKY)
63	DOUBLE PRECISION ZFSUP(KDLON, KFLEV+1)
64	DOUBLE PRECISION ZFSDN(KDLON, KFLEV+1)
65	DOUBLE PRECISION ZFSUP0(KDLON, KFLEV+1)
66	DOUBLE PRECISION ZFSDN0(KDLON, KFLEV+1)
67
68	DOUBLE PRECISION, intent(out):: PTOPSWAD(KDLON)
69	! (diagnosed aerosol forcing)SHORTWAVE FLUX AT T.O.A.(+AEROSOL DIR)
70
71	DOUBLE PRECISION, intent(out):: PSOLSWAD(KDLON)
72	! (diagnosed aerosol forcing)SHORTWAVE FLUX AT SURFACE(+AEROSOL DIR)
73
74	DOUBLE PRECISION, intent(out):: PTOPSWAI(KDLON)
75	! (diagnosed aerosol forcing)SHORTWAVE FLUX AT T.O.A.(+AEROSOL IND)
76
77	DOUBLE PRECISION, intent(out):: PSOLSWAI(KDLON)
78	! (diagnosed aerosol forcing)SHORTWAVE FLUX AT SURFACE(+AEROSOL IND)
79
80	logical, intent(in):: ok_ade ! use aerosol forcings or not?
81
82	! Local:
83
84	DOUBLE PRECISION ZOZ(KDLON, KFLEV)
85	DOUBLE PRECISION ZAKI(KDLON, 2)
86	DOUBLE PRECISION ZCLD(KDLON, KFLEV)
87	DOUBLE PRECISION ZCLEAR(KDLON)
88	DOUBLE PRECISION ZDSIG(KDLON, KFLEV)
89	DOUBLE PRECISION ZFACT(KDLON)
90	DOUBLE PRECISION ZFD(KDLON, KFLEV+1)
91	DOUBLE PRECISION ZFDOWN(KDLON, KFLEV+1)
92	DOUBLE PRECISION ZFU(KDLON, KFLEV+1)
93	DOUBLE PRECISION ZFUP(KDLON, KFLEV+1)
94	DOUBLE PRECISION ZRMU(KDLON)
95	DOUBLE PRECISION ZSEC(KDLON)
96	DOUBLE PRECISION ZUD(KDLON, 5, KFLEV+1)
97	DOUBLE PRECISION ZCLDSW0(KDLON, KFLEV)
98
99	INTEGER inu, jl, jk, i, k, kpl1
100
101	INTEGER, PARAMETER:: swpas = 1 ! Every swpas steps, sw is calculated
102
103	INTEGER:: itapsw = 0
104	LOGICAL:: appel1er = .TRUE.
105	!jq-Introduced for aerosol forcings
106	logical, save:: flag_aer
107
108	!jq - Fluxes including aerosol effects
109	DOUBLE PRECISION, save:: ZFSUPAD(KDLON, KFLEV+1)
110	DOUBLE PRECISION, save:: ZFSDNAD(KDLON, KFLEV+1)
111	DOUBLE PRECISION, save:: ZFSUPAI(KDLON, KFLEV+1)
112	DOUBLE PRECISION, save:: ZFSDNAI(KDLON, KFLEV+1)
113
114	logical:: initialized = .false.
115	REAL, PARAMETER :: dobson_u = 2.1415E-05 ! Dobson unit, in kg m-2
116
117	!-------------------------------------------------------------------
118
119	if(.not.initialized) then
120	flag_aer=.false.
121	initialized=.TRUE.
122	ZFSUPAD = 0.
123	ZFSDNAD = 0.
124	ZFSUPAI = 0.
125	ZFSDNAI = 0.
126	endif
127	!rv
128
129	IF (appel1er) THEN
130	PRINT*, 'SW calling frequency: ', swpas
131	PRINT*, " In general, it should be 1"
132	appel1er = .FALSE.
133	ENDIF
134
135	IF (MOD(itapsw, swpas).EQ.0) THEN
136	DO JK = 1 , KFLEV
137	DO JL = 1, KDLON
138	ZCLDSW0(JL, JK) = 0.0
139	ZOZ(JL, JK) = POZON(JL, JK) / (dobson_u * 1E3 * rg) * PDP(JL, JK)
140	ENDDO
141	ENDDO
142
143	! clear-sky:
144	CALL SWU(PSCT, ZCLDSW0, PPMB, PPSOL, &
145	PRMU0, PFRAC, PTAVE, PWV, &
146	ZAKI, ZCLD, ZCLEAR, ZDSIG, ZFACT, ZRMU, ZSEC, ZUD)
147	INU = 1
148	CALL SW1S(INU, flag_aer, PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZDSIG, &
149	POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD, ZFD, ZFU)
150	INU = 2
151	CALL SW2S(INU, flag_aer, ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZDSIG, &
152	POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD, PWV, PQS, ZFDOWN, ZFUP)
153	DO JK = 1 , KFLEV+1
154	DO JL = 1, KDLON
155	ZFSUP0(JL, JK) = (ZFUP(JL, JK) + ZFU(JL, JK)) * ZFACT(JL)
156	ZFSDN0(JL, JK) = (ZFDOWN(JL, JK) + ZFD(JL, JK)) * ZFACT(JL)
157	ENDDO
158	ENDDO
159
160	flag_aer= .false.
161	CALL SWU(PSCT, PCLDSW, PPMB, PPSOL, &
162	PRMU0, PFRAC, PTAVE, PWV, &
163	ZAKI, ZCLD, ZCLEAR, ZDSIG, ZFACT, ZRMU, ZSEC, ZUD)
164	INU = 1
165	CALL SW1S(INU, .false., PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZDSIG, &
166	POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD, ZFD, ZFU)
167	INU = 2
168	CALL SW2S(INU, .false., ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZDSIG, &
169	POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD, PWV, PQS, ZFDOWN, ZFUP)
170
171	! cloudy-sky:
172
173	DO JK = 1 , KFLEV+1
174	DO JL = 1, KDLON
175	ZFSUP(JL, JK) = (ZFUP(JL, JK) + ZFU(JL, JK)) * ZFACT(JL)
176	ZFSDN(JL, JK) = (ZFDOWN(JL, JK) + ZFD(JL, JK)) * ZFACT(JL)
177	ENDDO
178	ENDDO
179
180	IF (ok_ade) THEN
181	! cloudy-sky + aerosol dir OB
182	flag_aer= .true.
183	CALL SWU(PSCT, PCLDSW, PPMB, PPSOL, PRMU0, PFRAC, PTAVE, PWV, ZAKI, &
184	ZCLD, ZCLEAR, ZDSIG, ZFACT, ZRMU, ZSEC, ZUD)
185	INU = 1
186	CALL SW1S(INU, .true., PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZDSIG, &
187	POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD, ZFD, ZFU)
188	INU = 2
189	CALL SW2S(INU, .true., ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, &
190	ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD, PWV, PQS, ZFDOWN, &
191	ZFUP)
192	DO JK = 1 , KFLEV+1
193	DO JL = 1, KDLON
194	ZFSUPAD(JL, JK) = ZFSUP(JL, JK)
195	ZFSDNAD(JL, JK) = ZFSDN(JL, JK)
196	ZFSUP(JL, JK) = (ZFUP(JL, JK) + ZFU(JL, JK)) * ZFACT(JL)
197	ZFSDN(JL, JK) = (ZFDOWN(JL, JK) + ZFD(JL, JK)) * ZFACT(JL)
198	ENDDO
199	ENDDO
200	ENDIF
201
202	itapsw = 0
203	ENDIF
204	itapsw = itapsw + 1
205
206	DO k = 1, KFLEV
207	kpl1 = k+1
208	DO i = 1, KDLON
209	PHEAT(i, k) = -(ZFSUP(i, kpl1)-ZFSUP(i, k)) &
210	-(ZFSDN(i, k)-ZFSDN(i, kpl1))
211	PHEAT(i, k) = PHEAT(i, k) * RDAY*RG/RCPD / PDP(i, k)
212	PHEAT0(i, k) = -(ZFSUP0(i, kpl1)-ZFSUP0(i, k)) &
213	-(ZFSDN0(i, k)-ZFSDN0(i, kpl1))
214	PHEAT0(i, k) = PHEAT0(i, k) * RDAY*RG/RCPD / PDP(i, k)
215	ENDDO
216	ENDDO
217	DO i = 1, KDLON
218	PALBPLA(i) = ZFSUP(i, KFLEV+1)/(ZFSDN(i, KFLEV+1)+1.0e-20)
219
220	PSOLSW(i) = ZFSDN(i, 1) - ZFSUP(i, 1)
221	PTOPSW(i) = ZFSDN(i, KFLEV+1) - ZFSUP(i, KFLEV+1)
222
223	PSOLSW0(i) = ZFSDN0(i, 1) - ZFSUP0(i, 1)
224	PTOPSW0(i) = ZFSDN0(i, KFLEV+1) - ZFSUP0(i, KFLEV+1)
225
226	PSOLSWAD(i) = ZFSDNAD(i, 1) - ZFSUPAD(i, 1)
227	PTOPSWAD(i) = ZFSDNAD(i, KFLEV+1) - ZFSUPAD(i, KFLEV+1)
228
229	PSOLSWAI(i) = ZFSDNAI(i, 1) - ZFSUPAI(i, 1)
230	PTOPSWAI(i) = ZFSDNAI(i, KFLEV+1) - ZFSUPAI(i, KFLEV+1)
231	ENDDO
232
233	END SUBROUTINE SW
234
235	end module sw_m