source: lmdz_wrf/trunk/WRFV3/lmdz/sw_aeroAR4.F90 @ 354

Last change on this file since 354 was 1, checked in by lfita, 10 years ago
  • -- --- Opening of the WRF+LMDZ coupling repository --- -- -

WRF: version v3.3
LMDZ: version v1818

More details in:

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