Context Navigation

← Previous Revision
Next Revision →
Blame
Revision Log

sw_aeroAR4.f90

Last change on this file was 5274, checked in by abarral, 54 minutes ago
Replace yomcst.h by existing module
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: 24.6 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format