GCC Code Coverage Report
Directory: ./ Exec Total Coverage
File: phylmd/sw_aeroAR4.F90 Lines: 0 207 0.0 %
Date: 2023-06-30 12:51:15 Branches: 0 470 0.0 %

Line Branch Exec Source
1
!
2
! $Id$
3
!
4
SUBROUTINE SW_AEROAR4(PSCT, PRMU0, PFRAC, &
5
     PPMB, PDP, &
6
     PPSOL, PALBD, PALBP,&
7
     PTAVE, PWV, PQS, POZON, PAER,&
8
     PCLDSW, PTAU, POMEGA, PCG,&
9
     PHEAT, PHEAT0,&
10
     PALBPLA,PTOPSW,PSOLSW,PTOPSW0,PSOLSW0,&
11
     ZFSUP,ZFSDN,ZFSUP0,ZFSDN0,&
12
     tauaero, pizaero, cgaero,&
13
     PTAUA, POMEGAA,&
14
     PTOPSWADAERO,PSOLSWADAERO,&
15
     PTOPSWAD0AERO,PSOLSWAD0AERO,&
16
     PTOPSWAIAERO,PSOLSWAIAERO,&
17
     PTOPSWAERO,PTOPSW0AERO,&
18
     PSOLSWAERO,PSOLSW0AERO,&
19
     PTOPSWCFAERO,PSOLSWCFAERO,&
20
     ok_ade, ok_aie, flag_aerosol, flag_aerosol_strat )
21
22
  USE dimphy
23
  USE phys_output_mod, ONLY : swaero_diag
24
  USE print_control_mod, ONLY: lunout
25
  USE aero_mod, ONLY : naero_grp
26
  IMPLICIT NONE
27
28
  INCLUDE "YOMCST.h"
29
  INCLUDE "clesphys.h"
30
  !
31
  !     ------------------------------------------------------------------
32
  !
33
  !     PURPOSE.
34
  !     --------
35
  !
36
  !          THIS ROUTINE COMPUTES THE SHORTWAVE RADIATION FLUXES IN TWO
37
  !     SPECTRAL INTERVALS FOLLOWING FOUQUART AND BONNEL (1980).
38
  !
39
  !     METHOD.
40
  !     -------
41
  !
42
  !          1. COMPUTES ABSORBER AMOUNTS                 (SWU)
43
  !          2. COMPUTES FLUXES IN 1ST SPECTRAL INTERVAL  (SW1S)
44
  !          3. COMPUTES FLUXES IN 2ND SPECTRAL INTERVAL  (SW2S)
45
  !
46
  !     REFERENCE.
47
  !     ----------
48
  !
49
  !        SEE RADIATION'S PART OF THE ECMWF RESEARCH DEPARTMENT
50
  !        DOCUMENTATION, AND FOUQUART AND BONNEL (1980)
51
  !
52
  !     AUTHOR.
53
  !     -------
54
  !        JEAN-JACQUES MORCRETTE  *ECMWF*
55
  !
56
  !     MODIFICATIONS.
57
  !     --------------
58
  !        ORIGINAL : 89-07-14
59
  !        1995-01-01  J.-J. MORCRETTE  Direct/Diffuse Albedo
60
  !        2003-11-27  J. QUAAS Introduce aerosol forcings (based on BOUCHER)
61
  !        2009-04     A. COZIC - C.DEANDREIS Indroduce NAT/BC/POM/DUST/SS aerosol forcing
62
  !        2012-09     O. BOUCHER - reorganise aerosol cases with ok_ade, ok_aie, flag_aerosol
63
  !     ------------------------------------------------------------------
64
  !
65
  !* ARGUMENTS:
66
  !
67
  REAL(KIND=8) PSCT  ! constante solaire (valeur conseillee: 1370)
68
69
  REAL(KIND=8) PPSOL(KDLON)        ! SURFACE PRESSURE (PA)
70
  REAL(KIND=8) PDP(KDLON,KFLEV)    ! LAYER THICKNESS (PA)
71
  REAL(KIND=8) PPMB(KDLON,KFLEV+1) ! HALF-LEVEL PRESSURE (MB)
72
73
  REAL(KIND=8) PRMU0(KDLON)  ! COSINE OF ZENITHAL ANGLE
74
  REAL(KIND=8) PFRAC(KDLON)  ! fraction de la journee
75
76
  REAL(KIND=8) PTAVE(KDLON,KFLEV)  ! LAYER TEMPERATURE (K)
77
  REAL(KIND=8) PWV(KDLON,KFLEV)    ! SPECIFI! HUMIDITY (KG/KG)
78
  REAL(KIND=8) PQS(KDLON,KFLEV)    ! SATURATED WATER VAPOUR (KG/KG)
79
  REAL(KIND=8) POZON(KDLON,KFLEV)  ! OZONE CONCENTRATION (KG/KG)
80
  REAL(KIND=8) PAER(KDLON,KFLEV,5) ! AEROSOLS' OPTICAL THICKNESS
81
82
  REAL(KIND=8) PALBD(KDLON,2)  ! albedo du sol (lumiere diffuse)
83
  REAL(KIND=8) PALBP(KDLON,2)  ! albedo du sol (lumiere parallele)
84
85
  REAL(KIND=8) PCLDSW(KDLON,KFLEV)    ! CLOUD FRACTION
86
  REAL(KIND=8) PTAU(KDLON,2,KFLEV)    ! CLOUD OPTICAL THICKNESS (pre-industrial value)
87
  REAL(KIND=8) PCG(KDLON,2,KFLEV)     ! ASYMETRY FACTOR
88
  REAL(KIND=8) POMEGA(KDLON,2,KFLEV)  ! SINGLE SCATTERING ALBEDO
89
90
  REAL(KIND=8) PHEAT(KDLON,KFLEV) ! SHORTWAVE HEATING (K/DAY)
91
  REAL(KIND=8) PHEAT0(KDLON,KFLEV)! SHORTWAVE HEATING (K/DAY) clear-sky
92
  REAL(KIND=8) PALBPLA(KDLON)     ! PLANETARY ALBEDO
93
  REAL(KIND=8) PTOPSW(KDLON)      ! SHORTWAVE FLUX AT T.O.A.
94
  REAL(KIND=8) PSOLSW(KDLON)      ! SHORTWAVE FLUX AT SURFACE
95
  REAL(KIND=8) PTOPSW0(KDLON)     ! SHORTWAVE FLUX AT T.O.A. (CLEAR-SKY)
96
  REAL(KIND=8) PSOLSW0(KDLON)     ! SHORTWAVE FLUX AT SURFACE (CLEAR-SKY)
97
  !
98
  !* LOCAL VARIABLES:
99
  !
100
  real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2
101
102
  REAL(KIND=8) ZOZ(KDLON,KFLEV)
103
  ! column-density of ozone in layer, in kilo-Dobsons
104
105
  REAL(KIND=8) ZAKI(KDLON,2)
106
  REAL(KIND=8) ZCLD(KDLON,KFLEV)
107
  REAL(KIND=8) ZCLEAR(KDLON)
108
  REAL(KIND=8) ZDSIG(KDLON,KFLEV)
109
  REAL(KIND=8) ZFACT(KDLON)
110
  REAL(KIND=8) ZFD(KDLON,KFLEV+1)
111
  REAL(KIND=8) ZFDOWN(KDLON,KFLEV+1)
112
  REAL(KIND=8) ZFU(KDLON,KFLEV+1)
113
  REAL(KIND=8) ZFUP(KDLON,KFLEV+1)
114
  REAL(KIND=8) ZRMU(KDLON)
115
  REAL(KIND=8) ZSEC(KDLON)
116
  REAL(KIND=8) ZUD(KDLON,5,KFLEV+1)
117
  REAL(KIND=8) ZCLDSW0(KDLON,KFLEV)
118
119
  REAL(KIND=8) ZFSUP(KDLON,KFLEV+1)
120
  REAL(KIND=8) ZFSDN(KDLON,KFLEV+1)
121
  REAL(KIND=8) ZFSUP0(KDLON,KFLEV+1)
122
  REAL(KIND=8) ZFSDN0(KDLON,KFLEV+1)
123
124
  INTEGER inu, jl, jk, i, k, kpl1
125
126
  INTEGER swpas  ! Every swpas steps, sw is calculated
127
  PARAMETER(swpas=1)
128
129
  INTEGER, SAVE :: itapsw = 0
130
  !$OMP THREADPRIVATE(itapsw)
131
  LOGICAL, SAVE :: appel1er = .TRUE.
132
  !$OMP THREADPRIVATE(appel1er)
133
  LOGICAL, SAVE :: initialized = .FALSE.
134
  !$OMP THREADPRIVATE(initialized)
135
136
  !jq-local flag introduced for aerosol forcings
137
  REAL(KIND=8), SAVE :: flag_aer
138
  !$OMP THREADPRIVATE(flag_aer)
139
140
  LOGICAL ok_ade, ok_aie    ! use aerosol forcings or not?
141
  INTEGER flag_aerosol_strat ! use stratospehric aerosols
142
  INTEGER flag_aerosol      ! global flag for aerosol 0 (no aerosol) or 1-5 (aerosols)
143
  REAL(KIND=8) tauaero(kdlon,kflev,naero_grp,2)  ! aerosol optical properties
144
  REAL(KIND=8) pizaero(kdlon,kflev,naero_grp,2)  ! (see aeropt.F)
145
  REAL(KIND=8) cgaero(kdlon,kflev,naero_grp,2)   ! -"-
146
  REAL(KIND=8) PTAUA(KDLON,2,KFLEV)    ! CLOUD OPTICAL THICKNESS (present-day value)
147
  REAL(KIND=8) POMEGAA(KDLON,2,KFLEV)  ! SINGLE SCATTERING ALBEDO
148
  REAL(KIND=8) PTOPSWADAERO(KDLON)     ! SHORTWAVE FLUX AT T.O.A.(+AEROSOL DIR)
149
  REAL(KIND=8) PSOLSWADAERO(KDLON)     ! SHORTWAVE FLUX AT SURFACE(+AEROSOL DIR)
150
  REAL(KIND=8) PTOPSWAD0AERO(KDLON)    ! SHORTWAVE FLUX AT T.O.A.(+AEROSOL DIR)
151
  REAL(KIND=8) PSOLSWAD0AERO(KDLON)    ! SHORTWAVE FLUX AT SURFACE(+AEROSOL DIR)
152
  REAL(KIND=8) PTOPSWAIAERO(KDLON)     ! SHORTWAVE FLUX AT T.O.A.(+AEROSOL IND)
153
  REAL(KIND=8) PSOLSWAIAERO(KDLON)     ! SHORTWAVE FLUX AT SURFACE(+AEROSOL IND)
154
  REAL(KIND=8) PTOPSWAERO(KDLON,9)     ! SW TOA AS DRF nat & ant
155
  REAL(KIND=8) PTOPSW0AERO(KDLON,9)    ! SW SRF AS DRF nat & ant
156
  REAL(KIND=8) PSOLSWAERO(KDLON,9)     ! SW TOA CS DRF nat & ant
157
  REAL(KIND=8) PSOLSW0AERO(KDLON,9)    ! SW SRF CS DRF nat & ant
158
  REAL(KIND=8) PTOPSWCFAERO(KDLON,3)   !  SW TOA AS cloudRF nat & ant
159
  REAL(KIND=8) PSOLSWCFAERO(KDLON,3)   !  SW SRF AS cloudRF nat & ant
160
161
  !jq - Fluxes including aerosol effects
162
  REAL(KIND=8),ALLOCATABLE,SAVE :: ZFSUPAD_AERO(:,:)
163
  !$OMP THREADPRIVATE(ZFSUPAD_AERO)
164
  REAL(KIND=8),ALLOCATABLE,SAVE :: ZFSDNAD_AERO(:,:)
165
  !$OMP THREADPRIVATE(ZFSDNAD_AERO)
166
  !jq - Fluxes including aerosol effects
167
  REAL(KIND=8),ALLOCATABLE,SAVE :: ZFSUPAD0_AERO(:,:)
168
  !$OMP THREADPRIVATE(ZFSUPAD0_AERO)
169
  REAL(KIND=8),ALLOCATABLE,SAVE :: ZFSDNAD0_AERO(:,:)
170
  !$OMP THREADPRIVATE(ZFSDNAD0_AERO)
171
  REAL(KIND=8),ALLOCATABLE,SAVE :: ZFSUPAI_AERO(:,:)
172
  !$OMP THREADPRIVATE(ZFSUPAI_AERO)
173
  REAL(KIND=8),ALLOCATABLE,SAVE :: ZFSDNAI_AERO(:,:)
174
  !$OMP THREADPRIVATE(ZFSDNAI_AERO)
175
  REAL(KIND=8),ALLOCATABLE,SAVE ::  ZFSUP_AERO(:,:,:)
176
  !$OMP THREADPRIVATE(ZFSUP_AERO)
177
  REAL(KIND=8),ALLOCATABLE,SAVE ::  ZFSDN_AERO(:,:,:)
178
  !$OMP THREADPRIVATE(ZFSDN_AERO)
179
  REAL(KIND=8),ALLOCATABLE,SAVE ::  ZFSUP0_AERO(:,:,:)
180
  !$OMP THREADPRIVATE(ZFSUP0_AERO)
181
  REAL(KIND=8),ALLOCATABLE,SAVE ::  ZFSDN0_AERO(:,:,:)
182
  !$OMP THREADPRIVATE(ZFSDN0_AERO)
183
184
! Key to define the aerosol effect acting on climate
185
! OB: AEROSOLFEEDBACK_ACTIVE is now a LOGICAL
186
! TRUE: fluxes use natural and/or anthropogenic aerosols according to ok_ade and ok_aie, DEFAULT
187
! FALSE: fluxes use no aerosols (case 1)
188
189
  LOGICAL,SAVE :: AEROSOLFEEDBACK_ACTIVE = .TRUE.
190
!$OMP THREADPRIVATE(AEROSOLFEEDBACK_ACTIVE)
191
192
      CHARACTER (LEN=20) :: modname='sw_aeroAR4'
193
      CHARACTER (LEN=80) :: abort_message
194
195
  IF(.NOT.initialized) THEN
196
     flag_aer=0.
197
     initialized=.TRUE.
198
     ALLOCATE(ZFSUPAD_AERO(KDLON,KFLEV+1))
199
     ALLOCATE(ZFSDNAD_AERO(KDLON,KFLEV+1))
200
     ALLOCATE(ZFSUPAD0_AERO(KDLON,KFLEV+1))
201
     ALLOCATE(ZFSDNAD0_AERO(KDLON,KFLEV+1))
202
     ALLOCATE(ZFSUPAI_AERO(KDLON,KFLEV+1))
203
     ALLOCATE(ZFSDNAI_AERO(KDLON,KFLEV+1))
204
!-OB decrease size of these arrays to what is needed
205
!                | direct effect
206
!ind effect      | no aerosol   natural  total
207
!natural (PTAU)  |   1            3       2     --ZFSUP/ZFSDN
208
!total (PTAUA)   |                5       4     --ZFSUP/ZFSDN
209
!no cloud        |   1            3       2     --ZFSUP0/ZFSDN0
210
! so we need which case when ?
211
! ok_ade and ok_aie = 4-5, 4-2 and 2
212
! ok_ade and not ok_aie = 2-3 and 2
213
! not ok_ade and ok_aie = 5-3 and 5
214
! not ok_ade and not ok_aie = 3
215
! therefore the cases have the folliwng switches
216
! 3 = not ok_ade or not ok_aie
217
! 4 = ok_ade and ok_aie
218
! 2 = ok_ade
219
! 5 = ok_aie
220
     ALLOCATE(ZFSUP_AERO (KDLON,KFLEV+1,5))
221
     ALLOCATE(ZFSDN_AERO (KDLON,KFLEV+1,5))
222
     ALLOCATE(ZFSUP0_AERO(KDLON,KFLEV+1,3))
223
     ALLOCATE(ZFSDN0_AERO(KDLON,KFLEV+1,3))
224
! end OB modif
225
     ZFSUPAD_AERO(:,:)=0.
226
     ZFSDNAD_AERO(:,:)=0.
227
     ZFSUPAD0_AERO(:,:)=0.
228
     ZFSDNAD0_AERO(:,:)=0.
229
     ZFSUPAI_AERO(:,:)=0.
230
     ZFSDNAI_AERO(:,:)=0.
231
     ZFSUP_AERO (:,:,:)=0.
232
     ZFSDN_AERO (:,:,:)=0.
233
     ZFSUP0_AERO(:,:,:)=0.
234
     ZFSDN0_AERO(:,:,:)=0.
235
  ENDIF
236
237
  IF (appel1er) THEN
238
     WRITE(lunout,*)'SW calling frequency : ', swpas
239
     WRITE(lunout,*) "   In general, it should be 1"
240
     appel1er = .FALSE.
241
  ENDIF
242
  !     ------------------------------------------------------------------
243
  IF (MOD(itapsw,swpas).EQ.0) THEN
244
245
     DO JK = 1 , KFLEV
246
        DO JL = 1, KDLON
247
           ZCLDSW0(JL,JK) = 0.0
248
           ZOZ(JL,JK) = POZON(JL,JK)*46.6968/RG &
249
                *PDP(JL,JK)*(101325.0/PPSOL(JL))
250
        ENDDO
251
     ENDDO
252
253
! clear sky with no aerosols at all is computed IF ACTIVEFEEDBACK_ACTIVE is false or for extended diag
254
     IF ( swaero_diag .or. .not. AEROSOLFEEDBACK_ACTIVE .OR. flag_aerosol .EQ. 0 ) THEN
255
256
     ! clear-sky: zero aerosol effect
257
     flag_aer=0.0
258
     CALL SWU_LMDAR4(PSCT,ZCLDSW0,PPMB,PPSOL,&
259
          PRMU0,PFRAC,PTAVE,PWV,&
260
          ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
261
     INU = 1
262
     CALL SW1S_LMDAR4(INU,PAER, flag_aer, &
263
          tauaero(:,:,1,:), pizaero(:,:,1,:), cgaero(:,:,1,:),&
264
          PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZCLDSW0,&
265
          ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
266
          ZFD, ZFU)
267
     INU = 2
268
     CALL SW2S_LMDAR4(INU, PAER, flag_aer, &
269
          tauaero(:,:,1,:), pizaero(:,:,1,:), cgaero(:,:,1,:),&
270
          ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, ZCLDSW0,&
271
          ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
272
          PWV, PQS,&
273
          ZFDOWN, ZFUP)
274
     DO JK = 1 , KFLEV+1
275
        DO JL = 1, KDLON
276
           ZFSUP0_AERO(JL,JK,1) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
277
           ZFSDN0_AERO(JL,JK,1) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
278
        ENDDO
279
     ENDDO
280
     ENDIF ! swaero_diag .or. .not. AEROSOLFEEDBACK_ACTIVE
281
282
! cloudy sky with no aerosols at all is either computed IF no indirect effect is asked for, or for extended diag
283
     IF ( swaero_diag .or. .not. AEROSOLFEEDBACK_ACTIVE .OR. flag_aerosol .EQ. 0 ) THEN
284
     ! cloudy-sky: zero aerosol effect
285
     flag_aer=0.0
286
     CALL SWU_LMDAR4(PSCT,PCLDSW,PPMB,PPSOL,&
287
          PRMU0,PFRAC,PTAVE,PWV,&
288
          ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
289
     INU = 1
290
     CALL SW1S_LMDAR4(INU, PAER, flag_aer, &
291
          tauaero(:,:,1,:), pizaero(:,:,1,:), cgaero(:,:,1,:),&
292
          PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
293
          ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
294
          ZFD, ZFU)
295
     INU = 2
296
     CALL SW2S_LMDAR4(INU, PAER, flag_aer, &
297
          tauaero(:,:,1,:), pizaero(:,:,1,:), cgaero(:,:,1,:),&
298
          ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
299
          ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
300
          PWV, PQS,&
301
          ZFDOWN, ZFUP)
302
303
     DO JK = 1 , KFLEV+1
304
        DO JL = 1, KDLON
305
           ZFSUP_AERO(JL,JK,1) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
306
           ZFSDN_AERO(JL,JK,1) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
307
        ENDDO
308
     ENDDO
309
     ENDIF ! swaero_diag .or. .not. AEROSOLFEEDBACK_ACTIVE
310
311
     IF (flag_aerosol.GT.0 .OR. flag_aerosol_strat.GT.0) THEN
312
313
     IF (ok_ade.and.swaero_diag .or. .not. ok_ade) THEN
314
315
        ! clear sky direct effect natural aerosol
316
        ! CAS AER (3)
317
        flag_aer=1.0
318
        CALL SWU_LMDAR4(PSCT,ZCLDSW0,PPMB,PPSOL,&
319
             PRMU0,PFRAC,PTAVE,PWV,&
320
             ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
321
        INU = 1
322
        CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
323
             tauaero(:,:,3,:), pizaero(:,:,3,:), cgaero(:,:,3,:),&
324
             PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
325
             ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
326
             ZFD, ZFU)
327
        INU = 2
328
        CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
329
             tauaero(:,:,3,:), pizaero(:,:,3,:), cgaero(:,:,3,:),&
330
             ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
331
             ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
332
             PWV, PQS,&
333
             ZFDOWN, ZFUP)
334
335
        DO JK = 1 , KFLEV+1
336
           DO JL = 1, KDLON
337
              ZFSUP0_AERO(JL,JK,3) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
338
              ZFSDN0_AERO(JL,JK,3) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
339
           ENDDO
340
        ENDDO
341
     ENDIF !--end not swaero_diag or not ok_ade
342
343
     IF (ok_ade) THEN
344
345
        ! clear sky direct effect of total aerosol
346
        ! CAS AER (2)
347
        flag_aer=1.0
348
        CALL SWU_LMDAR4(PSCT,ZCLDSW0,PPMB,PPSOL,&
349
             PRMU0,PFRAC,PTAVE,PWV,&
350
             ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
351
        INU = 1
352
        CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
353
             tauaero(:,:,2,:), pizaero(:,:,2,:), cgaero(:,:,2,:),&
354
             PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
355
             ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
356
             ZFD, ZFU)
357
        INU = 2
358
        CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
359
             tauaero(:,:,2,:), pizaero(:,:,2,:), cgaero(:,:,2,:),&
360
             ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
361
             ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
362
             PWV, PQS,&
363
             ZFDOWN, ZFUP)
364
365
        DO JK = 1 , KFLEV+1
366
           DO JL = 1, KDLON
367
              ZFSUP0_AERO(JL,JK,2) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
368
              ZFSDN0_AERO(JL,JK,2) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
369
           ENDDO
370
        ENDDO
371
372
        ! cloudy-sky with natural aerosols for indirect effect
373
        ! but total aerosols for direct effect
374
        ! PTAU
375
        ! CAS AER (2)
376
        flag_aer=1.0
377
        CALL SWU_LMDAR4(PSCT,PCLDSW,PPMB,PPSOL,&
378
             PRMU0,PFRAC,PTAVE,PWV,&
379
             ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
380
        INU = 1
381
        CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
382
             tauaero(:,:,2,:), pizaero(:,:,2,:), cgaero(:,:,2,:),&
383
             PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
384
             ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
385
             ZFD, ZFU)
386
        INU = 2
387
        CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
388
             tauaero(:,:,2,:), pizaero(:,:,2,:), cgaero(:,:,2,:),&
389
             ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
390
             ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
391
             PWV, PQS,&
392
             ZFDOWN, ZFUP)
393
394
        DO JK = 1 , KFLEV+1
395
           DO JL = 1, KDLON
396
              ZFSUP_AERO(JL,JK,2) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
397
              ZFSDN_AERO(JL,JK,2) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
398
           ENDDO
399
        ENDDO
400
401
     ENDIF !-end ok_ade
402
403
     IF ( .not. ok_ade .or. .not. ok_aie ) THEN
404
405
        ! cloudy-sky with natural aerosols for indirect effect
406
        ! and natural aerosols for direct effect
407
        ! PTAU
408
        ! CAS AER (3)
409
        ! cloudy-sky direct effect natural aerosol
410
        flag_aer=1.0
411
        CALL SWU_LMDAR4(PSCT,PCLDSW,PPMB,PPSOL,&
412
             PRMU0,PFRAC,PTAVE,PWV,&
413
             ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
414
        INU = 1
415
        CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
416
             tauaero(:,:,3,:), pizaero(:,:,3,:), cgaero(:,:,3,:),&
417
             PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
418
             ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
419
             ZFD, ZFU)
420
        INU = 2
421
        CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
422
             tauaero(:,:,3,:), pizaero(:,:,3,:), cgaero(:,:,3,:),&
423
             ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
424
             ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
425
             PWV, PQS,&
426
             ZFDOWN, ZFUP)
427
428
        DO JK = 1 , KFLEV+1
429
           DO JL = 1, KDLON
430
              ZFSUP_AERO(JL,JK,3) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
431
              ZFSDN_AERO(JL,JK,3) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
432
           ENDDO
433
        ENDDO
434
435
     ENDIF  !--true/false or false/true
436
437
     IF (ok_ade .and. ok_aie) THEN
438
439
        ! cloudy-sky with total aerosols for indirect effect
440
        ! and total aerosols for direct effect
441
        ! PTAUA
442
        ! CAS AER (2)
443
        flag_aer=1.0
444
        CALL SWU_LMDAR4(PSCT,PCLDSW,PPMB,PPSOL,&
445
             PRMU0,PFRAC,PTAVE,PWV,&
446
             ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
447
        INU = 1
448
        CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
449
             tauaero(:,:,2,:), pizaero(:,:,2,:), cgaero(:,:,2,:),&
450
             PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
451
             ZDSIG, POMEGAA, ZOZ, ZRMU, ZSEC, PTAUA, ZUD,&
452
             ZFD, ZFU)
453
        INU = 2
454
        CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
455
             tauaero(:,:,2,:), pizaero(:,:,2,:), cgaero(:,:,2,:),&
456
             ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
457
             ZDSIG, POMEGAA, ZOZ, ZRMU, ZSEC, PTAUA, ZUD,&
458
             PWV, PQS,&
459
             ZFDOWN, ZFUP)
460
461
        DO JK = 1 , KFLEV+1
462
           DO JL = 1, KDLON
463
              ZFSUP_AERO(JL,JK,4) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
464
              ZFSDN_AERO(JL,JK,4) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
465
           ENDDO
466
        ENDDO
467
468
      ENDIF ! ok_ade .and. ok_aie
469
470
     IF (ok_aie) THEN
471
        ! cloudy-sky with total aerosols for indirect effect
472
        ! and natural aerosols for direct effect
473
        ! PTAUA
474
        ! CAS AER (3)
475
        flag_aer=1.0
476
        CALL SWU_LMDAR4(PSCT,PCLDSW,PPMB,PPSOL,&
477
             PRMU0,PFRAC,PTAVE,PWV,&
478
             ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
479
        INU = 1
480
        CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
481
             tauaero(:,:,3,:), pizaero(:,:,3,:), cgaero(:,:,3,:),&
482
             PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
483
             ZDSIG, POMEGAA, ZOZ, ZRMU, ZSEC, PTAUA, ZUD,&
484
             ZFD, ZFU)
485
        INU = 2
486
        CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
487
             tauaero(:,:,3,:), pizaero(:,:,3,:), cgaero(:,:,3,:),&
488
             ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
489
             ZDSIG, POMEGAA, ZOZ, ZRMU, ZSEC, PTAUA, ZUD,&
490
             PWV, PQS,&
491
             ZFDOWN, ZFUP)
492
493
        DO JK = 1 , KFLEV+1
494
           DO JL = 1, KDLON
495
              ZFSUP_AERO(JL,JK,5) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
496
              ZFSDN_AERO(JL,JK,5) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
497
           ENDDO
498
        ENDDO
499
500
     ENDIF ! ok_aie
501
502
     ENDIF !--if flag_aerosol GT 0 OR flag_aerosol_strat GT 0
503
504
     itapsw = 0
505
  ENDIF
506
  itapsw = itapsw + 1
507
508
  IF  ( AEROSOLFEEDBACK_ACTIVE .AND. (flag_aerosol.GT.0 .OR. flag_aerosol_strat.GT.0) ) THEN
509
  IF ( ok_ade .and. ok_aie  ) THEN
510
    ZFSUP(:,:) =    ZFSUP_AERO(:,:,4)
511
    ZFSDN(:,:) =    ZFSDN_AERO(:,:,4)
512
    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,2)
513
    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,2)
514
  ENDIF
515
516
  IF ( ok_ade .and. (.not. ok_aie) )  THEN
517
    ZFSUP(:,:) =    ZFSUP_AERO(:,:,2)
518
    ZFSDN(:,:) =    ZFSDN_AERO(:,:,2)
519
    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,2)
520
    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,2)
521
  ENDIF
522
523
  IF ( (.not. ok_ade) .and. ok_aie  )  THEN
524
    ZFSUP(:,:) =    ZFSUP_AERO(:,:,5)
525
    ZFSDN(:,:) =    ZFSDN_AERO(:,:,5)
526
    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,3)
527
    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,3)
528
  ENDIF
529
530
  IF ((.not. ok_ade) .and. (.not. ok_aie)) THEN
531
    ZFSUP(:,:) =    ZFSUP_AERO(:,:,3)
532
    ZFSDN(:,:) =    ZFSDN_AERO(:,:,3)
533
    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,3)
534
    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,3)
535
  ENDIF
536
537
! MS the following allows to compute the forcing diagostics without
538
! letting the aerosol forcing act on the meteorology
539
! SEE logic above
540
  ELSE
541
    ZFSUP(:,:) =    ZFSUP_AERO(:,:,1)
542
    ZFSDN(:,:) =    ZFSDN_AERO(:,:,1)
543
    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,1)
544
    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,1)
545
  ENDIF
546
547
! Now computes heating rates
548
  DO k = 1, KFLEV
549
     kpl1 = k+1
550
     DO i = 1, KDLON
551
        PHEAT(i,k) = -(ZFSUP(i,kpl1)-ZFSUP(i,k))-(ZFSDN(i,k)-ZFSDN(i,kpl1))
552
        PHEAT(i,k) = PHEAT(i,k) * RDAY*RG/RCPD / PDP(i,k)
553
        PHEAT0(i,k) = -(ZFSUP0(i,kpl1)-ZFSUP0(i,k))-(ZFSDN0(i,k)-ZFSDN0(i,kpl1))
554
        PHEAT0(i,k) = PHEAT0(i,k) * RDAY*RG/RCPD / PDP(i,k)
555
     ENDDO
556
  ENDDO
557
558
  DO i = 1, KDLON
559
! effective SW surface albedo calculation
560
     PALBPLA(i) = ZFSUP(i,KFLEV+1)/(ZFSDN(i,KFLEV+1)+1.0e-20)
561
562
! clear sky net fluxes at TOA and SRF
563
     PSOLSW0(i) = ZFSDN0(i,1) - ZFSUP0(i,1)
564
     PTOPSW0(i) = ZFSDN0(i,KFLEV+1) - ZFSUP0(i,KFLEV+1)
565
566
! cloudy sky net fluxes at TOA and SRF
567
     PSOLSW(i) = ZFSDN(i,1) - ZFSUP(i,1)
568
     PTOPSW(i) = ZFSDN(i,KFLEV+1) - ZFSUP(i,KFLEV+1)
569
570
! net anthropogenic forcing direct and 1st indirect effect diagnostics
571
! requires a natural aerosol field read and used
572
! Difference of net fluxes from double call to radiation
573
574
IF (ok_ade) THEN
575
576
! indices 1: natural; 2 anthropogenic
577
578
! TOA/SRF all sky natural forcing
579
     PSOLSWAERO(i,1) = (ZFSDN_AERO(i,1,3) - ZFSUP_AERO(i,1,3))-(ZFSDN_AERO(i,1,1) - ZFSUP_AERO(i,1,1))
580
     PTOPSWAERO(i,1) = (ZFSDN_AERO(i,KFLEV+1,3) - ZFSUP_AERO(i,KFLEV+1,3))- (ZFSDN_AERO(i,KFLEV+1,1) - ZFSUP_AERO(i,KFLEV+1,1))
581
582
! TOA/SRF clear sky natural forcing
583
     PSOLSW0AERO(i,1) = (ZFSDN0_AERO(i,1,3) - ZFSUP0_AERO(i,1,3))-(ZFSDN0_AERO(i,1,1) - ZFSUP0_AERO(i,1,1))
584
     PTOPSW0AERO(i,1) = (ZFSDN0_AERO(i,KFLEV+1,3) - ZFSUP0_AERO(i,KFLEV+1,3))-(ZFSDN0_AERO(i,KFLEV+1,1) - ZFSUP0_AERO(i,KFLEV+1,1))
585
586
   IF (ok_aie) THEN
587
588
! TOA/SRF all sky anthropogenic forcing
589
     PSOLSWAERO(i,2) = (ZFSDN_AERO(i,1,4) - ZFSUP_AERO(i,1,4))-(ZFSDN_AERO(i,1,5) - ZFSUP_AERO(i,1,5))
590
     PTOPSWAERO(i,2) = (ZFSDN_AERO(i,KFLEV+1,4) - ZFSUP_AERO(i,KFLEV+1,4))- (ZFSDN_AERO(i,KFLEV+1,5) - ZFSUP_AERO(i,KFLEV+1,5))
591
592
   ELSE
593
594
! TOA/SRF all sky anthropogenic forcing
595
     PSOLSWAERO(i,2) = (ZFSDN_AERO(i,1,2) - ZFSUP_AERO(i,1,2))-(ZFSDN_AERO(i,1,3) - ZFSUP_AERO(i,1,3))
596
     PTOPSWAERO(i,2) = (ZFSDN_AERO(i,KFLEV+1,2) - ZFSUP_AERO(i,KFLEV+1,2))- (ZFSDN_AERO(i,KFLEV+1,3) - ZFSUP_AERO(i,KFLEV+1,3))
597
598
   ENDIF
599
600
! TOA/SRF clear sky anthropogenic forcing
601
     PSOLSW0AERO(i,2) = (ZFSDN0_AERO(i,1,2) - ZFSUP0_AERO(i,1,2))-(ZFSDN0_AERO(i,1,3) - ZFSUP0_AERO(i,1,3))
602
     PTOPSW0AERO(i,2) = (ZFSDN0_AERO(i,KFLEV+1,2) - ZFSUP0_AERO(i,KFLEV+1,2))-(ZFSDN0_AERO(i,KFLEV+1,3) - ZFSUP0_AERO(i,KFLEV+1,3))
603
604
! direct anthropogenic forcing , as in old LMDzT, however differences of net fluxes
605
     PSOLSWADAERO(i) = PSOLSWAERO(i,2)
606
     PTOPSWADAERO(i) = PTOPSWAERO(i,2)
607
     PSOLSWAD0AERO(i) = PSOLSW0AERO(i,2)
608
     PTOPSWAD0AERO(i) = PTOPSW0AERO(i,2)
609
610
! OB: these diagnostics may not always work but who need them
611
! Cloud forcing indices 1: natural; 2 anthropogenic; 3: zero aerosol direct effect
612
! Instantaneously computed cloudy sky direct aerosol effect, cloud forcing due to aerosols above clouds
613
! natural
614
     PSOLSWCFAERO(i,1) = PSOLSWAERO(i,1) - PSOLSW0AERO(i,1)
615
     PTOPSWCFAERO(i,1) = PTOPSWAERO(i,1) - PTOPSW0AERO(i,1)
616
617
! Instantaneously computed cloudy SKY DIRECT aerosol effect, cloud forcing due to aerosols above clouds
618
! anthropogenic
619
     PSOLSWCFAERO(i,2) = PSOLSWAERO(i,2) - PSOLSW0AERO(i,2)
620
     PTOPSWCFAERO(i,2) = PTOPSWAERO(i,2) - PTOPSW0AERO(i,2)
621
622
! Cloudforcing without aerosol
623
! zero
624
     PSOLSWCFAERO(i,3) = (ZFSDN_AERO(i,1,1) - ZFSUP_AERO(i,1,1))-(ZFSDN0_AERO(i,1,1) - ZFSUP0_AERO(i,1,1))
625
     PTOPSWCFAERO(i,3) = (ZFSDN_AERO(i,KFLEV+1,1) - ZFSUP_AERO(i,KFLEV+1,1))- (ZFSDN0_AERO(i,KFLEV+1,1) - ZFSUP0_AERO(i,KFLEV+1,1))
626
627
ENDIF
628
629
IF (ok_aie) THEN
630
   IF (ok_ade) THEN
631
     PSOLSWAIAERO(i) = (ZFSDN_AERO(i,1,4) - ZFSUP_AERO(i,1,4))-(ZFSDN_AERO(i,1,2) - ZFSUP_AERO(i,1,2))
632
     PTOPSWAIAERO(i) = (ZFSDN_AERO(i,KFLEV+1,4) - ZFSUP_AERO(i,KFLEV+1,4))-(ZFSDN_AERO(i,KFLEV+1,2) - ZFSUP_AERO(i,KFLEV+1,2))
633
   ELSE
634
     PSOLSWAIAERO(i) = (ZFSDN_AERO(i,1,5) - ZFSUP_AERO(i,1,5))-(ZFSDN_AERO(i,1,3) - ZFSUP_AERO(i,1,3))
635
     PTOPSWAIAERO(i) = (ZFSDN_AERO(i,KFLEV+1,5) - ZFSUP_AERO(i,KFLEV+1,5))-(ZFSDN_AERO(i,KFLEV+1,3) - ZFSUP_AERO(i,KFLEV+1,3))
636
   ENDIF
637
ENDIF
638
639
ENDDO
640
641
END SUBROUTINE SW_AEROAR4