source: LMDZ5/trunk/libf/phylmd/sw_aeroAR4.F90 @ 2183

Last change on this file since 2183 was 1907, checked in by lguez, 11 years ago

Added a copyright property to every file of the distribution, except
for the fcm files (which have their own copyright). Use svn propget on
a file to see the copyright. For instance:

$ svn propget copyright libf/phylmd/physiq.F90
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

Also added the files defining the CeCILL version 2 license, in French
and English, at the top of the LMDZ tree.

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