source: LMDZ4/trunk/libf/phylmd/sw_aeroAR4.F90 @ 1664

Last change on this file since 1664 was 1307, checked in by Laurent Fairhead, 15 years ago

Changes in output options for aerosols


Changements dans les options de sorties pour les aerosols

ACo

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