source: trunk/libf/phylmd/sw_aeroAR4.F90 @ 1

Last change on this file since 1 was 1, checked in by emillour, 14 years ago

Import initial LMDZ5
File size: 21.4 KB
RevLine 
[1]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.