source: LMDZ5/branches/testing/libf/phylmd/sw_aeroAR4.F90 @ 3946

Last change on this file since 3946 was 2542, checked in by Laurent Fairhead, 8 years ago

Merged trunk changes r2487:2541 into testing branch

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
File size: 23.9 KB
RevLine 
[1237]1!
2! $Id$
3!
[1153]4SUBROUTINE SW_AEROAR4(PSCT, PRMU0, PFRAC, &
[1150]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,&
[1246]19     PTOPSWCFAERO,PSOLSWCFAERO,&
[1795]20     ok_ade, ok_aie, flag_aerosol, flag_aerosol_strat )
[1150]21
22  USE dimphy
[1669]23  USE phys_output_mod, ONLY : swaero_diag
[2408]24  USE print_control_mod, ONLY: lunout
25  USE aero_mod, ONLY : naero_grp
[1150]26  IMPLICIT NONE
27
28#include "YOMCST.h"
[1246]29#include "clesphys.h"
[1150]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
[1669]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
[1150]63  !     ------------------------------------------------------------------
64  !
65  !* ARGUMENTS:
66  !
[1220]67  REAL(KIND=8) PSCT  ! constante solaire (valeur conseillee: 1370)
[1150]68
[1220]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)
[1150]72
[1220]73  REAL(KIND=8) PRMU0(KDLON)  ! COSINE OF ZENITHAL ANGLE
74  REAL(KIND=8) PFRAC(KDLON)  ! fraction de la journee
[1150]75
[1220]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
[1150]81
[1220]82  REAL(KIND=8) PALBD(KDLON,2)  ! albedo du sol (lumiere diffuse)
83  REAL(KIND=8) PALBP(KDLON,2)  ! albedo du sol (lumiere parallele)
[1150]84
[1220]85  REAL(KIND=8) PCLDSW(KDLON,KFLEV)    ! CLOUD FRACTION
[1669]86  REAL(KIND=8) PTAU(KDLON,2,KFLEV)    ! CLOUD OPTICAL THICKNESS (pre-industrial value)
[1220]87  REAL(KIND=8) PCG(KDLON,2,KFLEV)     ! ASYMETRY FACTOR
88  REAL(KIND=8) POMEGA(KDLON,2,KFLEV)  ! SINGLE SCATTERING ALBEDO
[1150]89
[1220]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)
[1150]97  !
98  !* LOCAL VARIABLES:
99  !
[1215]100  real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2
101
[1220]102  REAL(KIND=8) ZOZ(KDLON,KFLEV)
[1215]103  ! column-density of ozone in layer, in kilo-Dobsons
104
[1220]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)
[1150]118
[1220]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)
[1150]123
124  INTEGER inu, jl, jk, i, k, kpl1
125
126  INTEGER swpas  ! Every swpas steps, sw is calculated
127  PARAMETER(swpas=1)
128
[1159]129  INTEGER, SAVE :: itapsw = 0
130  !$OMP THREADPRIVATE(itapsw)
131  LOGICAL, SAVE :: appel1er = .TRUE.
[1150]132  !$OMP THREADPRIVATE(appel1er)
[1159]133  LOGICAL, SAVE :: initialized = .FALSE.
134  !$OMP THREADPRIVATE(initialized)
135
[1669]136  !jq-local flag introduced for aerosol forcings
[1220]137  REAL(KIND=8), SAVE :: flag_aer
[1159]138  !$OMP THREADPRIVATE(flag_aer)
139
[1150]140  LOGICAL ok_ade, ok_aie    ! use aerosol forcings or not?
[2542]141  INTEGER flag_aerosol_strat ! use stratospehric aerosols
[1669]142  INTEGER flag_aerosol      ! global flag for aerosol 0 (no aerosol) or 1-5 (aerosols)
[2408]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)   ! -"-
[1669]146  REAL(KIND=8) PTAUA(KDLON,2,KFLEV)    ! CLOUD OPTICAL THICKNESS (present-day value)
[1220]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)
[1669]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)
[1220]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)
[1669]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
[1246]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
[1150]160
161  !jq - Fluxes including aerosol effects
[1220]162  REAL(KIND=8),ALLOCATABLE,SAVE :: ZFSUPAD_AERO(:,:)
[1150]163  !$OMP THREADPRIVATE(ZFSUPAD_AERO)
[1220]164  REAL(KIND=8),ALLOCATABLE,SAVE :: ZFSDNAD_AERO(:,:)
[1150]165  !$OMP THREADPRIVATE(ZFSDNAD_AERO)
166  !jq - Fluxes including aerosol effects
[1220]167  REAL(KIND=8),ALLOCATABLE,SAVE :: ZFSUPAD0_AERO(:,:)
[1150]168  !$OMP THREADPRIVATE(ZFSUPAD0_AERO)
[1220]169  REAL(KIND=8),ALLOCATABLE,SAVE :: ZFSDNAD0_AERO(:,:)
[1150]170  !$OMP THREADPRIVATE(ZFSDNAD0_AERO)
[1220]171  REAL(KIND=8),ALLOCATABLE,SAVE :: ZFSUPAI_AERO(:,:)
[1150]172  !$OMP THREADPRIVATE(ZFSUPAI_AERO)
[1220]173  REAL(KIND=8),ALLOCATABLE,SAVE :: ZFSDNAI_AERO(:,:)
[1150]174  !$OMP THREADPRIVATE(ZFSDNAI_AERO)
[1220]175  REAL(KIND=8),ALLOCATABLE,SAVE ::  ZFSUP_AERO(:,:,:)
[1150]176  !$OMP THREADPRIVATE(ZFSUP_AERO)
[1220]177  REAL(KIND=8),ALLOCATABLE,SAVE ::  ZFSDN_AERO(:,:,:)
[1150]178  !$OMP THREADPRIVATE(ZFSDN_AERO)
[1220]179  REAL(KIND=8),ALLOCATABLE,SAVE ::  ZFSUP0_AERO(:,:,:)
[1150]180  !$OMP THREADPRIVATE(ZFSUP0_AERO)
[1220]181  REAL(KIND=8),ALLOCATABLE,SAVE ::  ZFSDN0_AERO(:,:,:)
[1150]182  !$OMP THREADPRIVATE(ZFSDN0_AERO)
183
[1246]184! Key to define the aerosol effect acting on climate
[1669]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)
[1150]188
[1669]189  LOGICAL,SAVE :: AEROSOLFEEDBACK_ACTIVE = .TRUE.
[1249]190!$OMP THREADPRIVATE(AEROSOLFEEDBACK_ACTIVE) 
[1246]191
[1307]192      CHARACTER (LEN=20) :: modname='sw_aeroAR4'
193      CHARACTER (LEN=80) :: abort_message
194
[1150]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))
[1669]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
[1150]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
[1669]238     WRITE(lunout,*)'SW calling frequency : ', swpas
[1664]239     WRITE(lunout,*) "   In general, it should be 1"
[1150]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
[1237]248           ZOZ(JL,JK) = POZON(JL,JK)*46.6968/RG &
249                *PDP(JL,JK)*(101325.0/PPSOL(JL))
[1150]250        ENDDO
251     ENDDO
252
[1669]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   
[1150]255
[1246]256     ! clear-sky: zero aerosol effect
[1150]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
[1237]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)
[1150]278        ENDDO
279     ENDDO
[1669]280     ENDIF ! swaero_diag .or. .not. AEROSOLFEEDBACK_ACTIVE
[1150]281
[1669]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   
[1246]284     ! cloudy-sky: zero aerosol effect
[1150]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
[1237]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)
[1150]307        ENDDO
308     ENDDO
[1669]309     ENDIF ! swaero_diag .or. .not. AEROSOLFEEDBACK_ACTIVE
[1150]310
[2542]311     IF (flag_aerosol.GT.0 .OR. flag_aerosol_strat.GT.0) THEN
[1150]312
[1669]313     IF (ok_ade.and.swaero_diag .or. .not. ok_ade) THEN
[1150]314
[1669]315        ! clear sky direct effect natural aerosol
316        ! CAS AER (3)
[1150]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,&
[1669]323             tauaero(:,:,3,:), pizaero(:,:,3,:), cgaero(:,:,3,:),&
[1150]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,&
[1669]329             tauaero(:,:,3,:), pizaero(:,:,3,:), cgaero(:,:,3,:),&
[1150]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
[1669]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)
[1150]339           ENDDO
340        ENDDO
[1669]341     ENDIF !--end not swaero_diag or not ok_ade
[1150]342
[1669]343     IF (ok_ade) THEN
344
345        ! clear sky direct effect of total aerosol
346        ! CAS AER (2)
[1150]347        flag_aer=1.0
[1669]348        CALL SWU_LMDAR4(PSCT,ZCLDSW0,PPMB,PPSOL,&
[1150]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
[1669]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)
[1150]369           ENDDO
370        ENDDO
371
[1669]372        ! cloudy-sky with natural aerosols for indirect effect
373        ! but total aerosols for direct effect
374        ! PTAU
375        ! CAS AER (2)
[1150]376        flag_aer=1.0
[1669]377        CALL SWU_LMDAR4(PSCT,PCLDSW,PPMB,PPSOL,&
[1150]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,&
[1669]382             tauaero(:,:,2,:), pizaero(:,:,2,:), cgaero(:,:,2,:),&
[1150]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,&
[1669]388             tauaero(:,:,2,:), pizaero(:,:,2,:), cgaero(:,:,2,:),&
[1150]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
[1669]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)
[1150]398           ENDDO
399        ENDDO
400
[1669]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)
[1246]409        ! cloudy-sky direct effect natural aerosol
[1150]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
[1669]435     ENDIF  !--true/false or false/true
[1150]436
[1669]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)
[1150]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)
[1669]460
[1150]461        DO JK = 1 , KFLEV+1
462           DO JL = 1, KDLON
[1237]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)
[1150]465           ENDDO
466        ENDDO
[1669]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
[1150]500     ENDIF ! ok_aie     
501
[2542]502     ENDIF !--if flag_aerosol GT 0 OR flag_aerosol_strat GT 0
[1669]503
[1150]504     itapsw = 0
505  ENDIF
506  itapsw = itapsw + 1
507
[2542]508  IF  ( AEROSOLFEEDBACK_ACTIVE .AND. (flag_aerosol.GT.0 .OR. flag_aerosol_strat.GT.0) ) THEN
[1237]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
[1669]515
[1237]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
[1246]522
[1237]523  IF ( (.not. ok_ade) .and. ok_aie  )  THEN
[1669]524    ZFSUP(:,:) =    ZFSUP_AERO(:,:,5)
525    ZFSDN(:,:) =    ZFSDN_AERO(:,:,5)
526    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,3)
527    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,3)
[1237]528  ENDIF
[1669]529
[1237]530  IF ((.not. ok_ade) .and. (.not. ok_aie)) THEN
[1669]531    ZFSUP(:,:) =    ZFSUP_AERO(:,:,3)
532    ZFSDN(:,:) =    ZFSDN_AERO(:,:,3)
533    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,3)
534    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,3)
[1237]535  ENDIF
536
537! MS the following allows to compute the forcing diagostics without
538! letting the aerosol forcing act on the meteorology
[1246]539! SEE logic above
[1669]540  ELSE
541    ZFSUP(:,:) =    ZFSUP_AERO(:,:,1)
542    ZFSDN(:,:) =    ZFSDN_AERO(:,:,1)
543    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,1)
544    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,1)
[1237]545  ENDIF
546
[1669]547! Now computes heating rates
[1150]548  DO k = 1, KFLEV
549     kpl1 = k+1
550     DO i = 1, KDLON
[1237]551        PHEAT(i,k) = -(ZFSUP(i,kpl1)-ZFSUP(i,k))-(ZFSDN(i,k)-ZFSDN(i,kpl1))
[1150]552        PHEAT(i,k) = PHEAT(i,k) * RDAY*RG/RCPD / PDP(i,k)
[1237]553        PHEAT0(i,k) = -(ZFSUP0(i,kpl1)-ZFSUP0(i,k))-(ZFSDN0(i,k)-ZFSDN0(i,kpl1))
[1150]554        PHEAT0(i,k) = PHEAT0(i,k) * RDAY*RG/RCPD / PDP(i,k)
555     ENDDO
556  ENDDO
[1237]557
[1150]558  DO i = 1, KDLON
[1246]559! effective SW surface albedo calculation
[1150]560     PALBPLA(i) = ZFSUP(i,KFLEV+1)/(ZFSDN(i,KFLEV+1)+1.0e-20)
[1246]561     
562! clear sky net fluxes at TOA and SRF
[1150]563     PSOLSW0(i) = ZFSDN0(i,1) - ZFSUP0(i,1)
564     PTOPSW0(i) = ZFSDN0(i,KFLEV+1) - ZFSUP0(i,KFLEV+1)
565
[1246]566! cloudy sky net fluxes at TOA and SRF
[1150]567     PSOLSW(i) = ZFSDN(i,1) - ZFSUP(i,1)
568     PTOPSW(i) = ZFSDN(i,KFLEV+1) - ZFSUP(i,KFLEV+1)
569
[1246]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
[1150]573
[1246]574IF (ok_ade) THEN
[1150]575
[1246]576! indices 1: natural; 2 anthropogenic
[1669]577
[1246]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))
[1150]581
[1669]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
[1246]588! TOA/SRF all sky anthropogenic forcing
[1669]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
[1246]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
[1669]598   ENDIF
[1246]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
[1669]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
[1246]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
627ENDIF
628
629IF (ok_aie) THEN
[1669]630   IF (ok_ade) THEN
[1237]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))
[1669]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
[1246]637ENDIF
[1237]638
[1669]639ENDDO
640
[1153]641END SUBROUTINE SW_AEROAR4
Note: See TracBrowser for help on using the repository browser.