Context Navigation

radlwsw_m.F90 @ 4773

Last change on this file since 4773 was 4773, checked in by idelkadi, 6 months ago

Update of Ecrad in LMDZ The same organization of the Ecrad offline version is retained in order to facilitate the updating of Ecrad in LMDZ and the comparison between online and offline results. version 1.6.1 of Ecrad (https://github.com/lguez/ecrad.git)

Implementation of the double call of Ecrad in LMDZ

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
Property svn:keywords set to Author Date Id Revi

File size: 63.8 KB

Line
1	!
2	! $Id: radlwsw_m.F90 4773 2023-12-14 22:00:52Z idelkadi $
3	!
4	module radlwsw_m
5
6	IMPLICIT NONE
7
8	contains
9
10	SUBROUTINE radlwsw( &
11	dist, rmu0, fract, &
12	!albedo SB >>>
13	! paprs, pplay,tsol,alb1, alb2, &
14	paprs, pplay,tsol,SFRWL,alb_dir, alb_dif, &
15	!albedo SB <<<
16	t,q,wo,&
17	cldfra, cldemi, cldtaupd,&
18	ok_ade, ok_aie, ok_volcan, flag_volc_surfstrat, flag_aerosol,&
19	flag_aerosol_strat, flag_aer_feedback, &
20	tau_aero, piz_aero, cg_aero,&
21	tau_aero_sw_rrtm, piz_aero_sw_rrtm, cg_aero_sw_rrtm,& ! rajoute par OB RRTM
22	tau_aero_lw_rrtm, & ! rajoute par C.Kleinschmitt pour RRTM
23	cldtaupi, &
24	qsat, flwc, fiwc, &
25	ref_liq, ref_ice, ref_liq_pi, ref_ice_pi, &
26	namelist_ecrad_file, &
27	heat,heat0,cool,cool0,albpla,&
28	heat_volc, cool_volc,&
29	topsw,toplw,solsw,solswfdiff,sollw,&
30	sollwdown,&
31	topsw0,toplw0,solsw0,sollw0,&
32	lwdnc0, lwdn0, lwdn, lwupc0, lwup0, lwup,&
33	swdnc0, swdn0, swdn, swupc0, swup0, swup,&
34	topswad_aero, solswad_aero,&
35	topswai_aero, solswai_aero, &
36	topswad0_aero, solswad0_aero,&
37	topsw_aero, topsw0_aero,&
38	solsw_aero, solsw0_aero, &
39	topswcf_aero, solswcf_aero,&
40	!-C. Kleinschmitt for LW diagnostics
41	toplwad_aero, sollwad_aero,&
42	toplwai_aero, sollwai_aero, &
43	toplwad0_aero, sollwad0_aero, &
44	!-end
45	ZLWFT0_i, ZFLDN0, ZFLUP0, &
46	ZSWFT0_i, ZFSDN0, ZFSUP0, &
47	cloud_cover_sw)
48
49	! Modules necessaires
50	USE DIMPHY
51	USE assert_m, ONLY : assert
52	USE infotrac_phy, ONLY : type_trac
53	USE write_field_phy
54
55	#ifdef REPROBUS
56	USE CHEM_REP, ONLY : solaireTIME, ok_SUNTIME, ndimozon
57	#endif
58
59	#ifdef CPP_RRTM
60	! modules necessaires au rayonnement
61	! -----------------------------------------
62	USE YOERAD , ONLY : NLW, LRRTM ,LCCNL ,LCCNO ,&
63	NRADIP , NRADLP , NICEOPT, NLIQOPT ,RCCNLND , RCCNSEA
64	USE YOELW , ONLY : NSIL ,NTRA ,NUA ,TSTAND ,XP
65	USE YOESW , ONLY : RYFWCA ,RYFWCB ,RYFWCC ,RYFWCD,&
66	RYFWCE ,RYFWCF ,REBCUA ,REBCUB ,REBCUC,&
67	REBCUD ,REBCUE ,REBCUF ,REBCUI ,REBCUJ,&
68	REBCUG ,REBCUH ,RHSAVI ,RFULIO ,RFLAA0,&
69	RFLAA1 ,RFLBB0 ,RFLBB1 ,RFLBB2 ,RFLBB3,&
70	RFLCC0 ,RFLCC1 ,RFLCC2 ,RFLCC3 ,RFLDD0,&
71	RFLDD1 ,RFLDD2 ,RFLDD3 ,RFUETA ,RASWCA,&
72	RASWCB ,RASWCC ,RASWCD ,RASWCE ,RASWCF
73	USE YOERDU , ONLY : NUAER ,NTRAER ,REPLOG ,REPSC ,REPSCW ,DIFF
74	USE YOERRTWN , ONLY : DELWAVE ,TOTPLNK
75	USE YOMPHY3 , ONLY : RII0
76	#endif
77	USE aero_mod
78
79	! AI 02.2021
80	! Besoin pour ECRAD de pctsrf, zmasq, longitude, altitude
81	#ifdef CPP_ECRAD
82	USE phys_local_var_mod, ONLY: rhcl, m_allaer
83	USE geometry_mod, ONLY: latitude, longitude
84	USE phys_state_var_mod, ONLY: pctsrf
85	USE indice_sol_mod
86	USE time_phylmdz_mod, only: current_time
87	USE phys_cal_mod, only: day_cur
88	USE interface_lmdz_ecrad
89	#endif
90
91	!======================================================================
92	! Auteur(s): Z.X. Li (LMD/CNRS) date: 19960719
93	! Objet: interface entre le modele et les rayonnements
94	! Arguments:
95	! INPUTS
96	! dist----- input-R- distance astronomique terre-soleil
97	! rmu0----- input-R- cosinus de l'angle zenithal
98	! fract---- input-R- duree d'ensoleillement normalisee
99	! co2_ppm-- input-R- concentration du gaz carbonique (en ppm)
100	! paprs---- input-R- pression a inter-couche (Pa)
101	! pplay---- input-R- pression au milieu de couche (Pa)
102	! tsol----- input-R- temperature du sol (en K)
103	! alb1----- input-R- albedo du sol(entre 0 et 1) dans l'interval visible
104	! alb2----- input-R- albedo du sol(entre 0 et 1) dans l'interval proche infra-rouge
105	! t-------- input-R- temperature (K)
106	! q-------- input-R- vapeur d'eau (en kg/kg)
107	! cldfra--- input-R- fraction nuageuse (entre 0 et 1)
108	! cldtaupd- input-R- epaisseur optique des nuages dans le visible (present-day value)
109	! cldemi--- input-R- emissivite des nuages dans l'IR (entre 0 et 1)
110	! ok_ade--- input-L- apply the Aerosol Direct Effect or not?
111	! ok_aie--- input-L- apply the Aerosol Indirect Effect or not?
112	! ok_volcan input-L- activate volcanic diags (SW heat & LW cool rate, SW & LW flux)
113	! flag_volc_surfstrat input-I- activate volcanic surf cooling or strato heating (or nothing)
114	! flag_aerosol input-I- aerosol flag from 0 to 6
115	! flag_aerosol_strat input-I- use stratospheric aerosols flag (0, 1, 2)
116	! flag_aer_feedback input-I- activate aerosol radiative feedback (T, F)
117	! tau_ae, piz_ae, cg_ae input-R- aerosol optical properties (calculated in aeropt.F)
118	! cldtaupi input-R- epaisseur optique des nuages dans le visible
119	! calculated for pre-industrial (pi) aerosol concentrations, i.e. with smaller
120	! droplet concentration, thus larger droplets, thus generally cdltaupi cldtaupd
121	! it is needed for the diagnostics of the aerosol indirect radiative forcing
122	!
123	! OUTPUTS
124	! heat-----output-R- echauffement atmospherique (visible) (K/jour)
125	! cool-----output-R- refroidissement dans l'IR (K/jour)
126	! albpla---output-R- albedo planetaire (entre 0 et 1)
127	! topsw----output-R- flux solaire net au sommet de l'atm.
128	! toplw----output-R- ray. IR montant au sommet de l'atmosphere
129	! solsw----output-R- flux solaire net a la surface
130	! solswfdiff----output-R- fraction de rayonnement diffus pour le flux solaire descendant a la surface
131	! sollw----output-R- ray. IR montant a la surface
132	! solswad---output-R- ray. solaire net absorbe a la surface (aerosol dir)
133	! topswad---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol dir)
134	! solswai---output-R- ray. solaire net absorbe a la surface (aerosol ind)
135	! topswai---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol ind)
136	!
137	! heat_volc-----output-R- echauffement atmospherique du au forcage volcanique (visible) (K/s)
138	! cool_volc-----output-R- refroidissement dans l'IR du au forcage volcanique (K/s)
139	!
140	! ATTENTION: swai and swad have to be interpreted in the following manner:
141	! ---------
142	! ok_ade=F & ok_aie=F -both are zero
143	! ok_ade=T & ok_aie=F -aerosol direct forcing is F_{AD} = topsw-topswad
144	! indirect is zero
145	! ok_ade=F & ok_aie=T -aerosol indirect forcing is F_{AI} = topsw-topswai
146	! direct is zero
147	! ok_ade=T & ok_aie=T -aerosol indirect forcing is F_{AI} = topsw-topswai
148	! aerosol direct forcing is F_{AD} = topswai-topswad
149	!
150	! --------- RRTM: output RECMWFL
151	! ZEMTD (KPROMA,KLEV+1) ; TOTAL DOWNWARD LONGWAVE EMISSIVITY
152	! ZEMTU (KPROMA,KLEV+1) ; TOTAL UPWARD LONGWAVE EMISSIVITY
153	! ZTRSO (KPROMA,KLEV+1) ; TOTAL SHORTWAVE TRANSMISSIVITY
154	! ZTH (KPROMA,KLEV+1) ; HALF LEVEL TEMPERATURE
155	! ZCTRSO(KPROMA,2) ; CLEAR-SKY SHORTWAVE TRANSMISSIVITY
156	! ZCEMTR(KPROMA,2) ; CLEAR-SKY NET LONGWAVE EMISSIVITY
157	! ZTRSOD(KPROMA) ; TOTAL-SKY SURFACE SW TRANSMISSITY
158	! ZLWFC (KPROMA,2) ; CLEAR-SKY LONGWAVE FLUXES
159	! ZLWFT (KPROMA,KLEV+1) ; TOTAL-SKY LONGWAVE FLUXES
160	! ZLWFT0(KPROMA,KLEV+1) ; CLEAR-SKY LONGWAVE FLUXES ! added by MPL 090109
161	! ZSWFC (KPROMA,2) ; CLEAR-SKY SHORTWAVE FLUXES
162	! ZSWFT (KPROMA,KLEV+1) ; TOTAL-SKY SHORTWAVE FLUXES
163	! ZSWFT0(KPROMA,KLEV+1) ; CLEAR-SKY SHORTWAVE FLUXES ! added by MPL 090109
164	! ZFLUX (KLON,2,KLEV+1) ; TOTAL LW FLUXES 1=up, 2=DWN ! added by MPL 080411
165	! ZFLUC (KLON,2,KLEV+1) ; CLEAR SKY LW FLUXES ! added by MPL 080411
166	! ZFSDWN(klon,KLEV+1) ; TOTAL SW DWN FLUXES ! added by MPL 080411
167	! ZFCDWN(klon,KLEV+1) ; CLEAR SKY SW DWN FLUXES ! added by MPL 080411
168	! ZFCCDWN(klon,KLEV+1) ; CLEAR SKY CLEAN (NO AEROSOL) SW DWN FLUXES ! added by OB 211117
169	! ZFSUP (klon,KLEV+1) ; TOTAL SW UP FLUXES ! added by MPL 080411
170	! ZFCUP (klon,KLEV+1) ; CLEAR SKY SW UP FLUXES ! added by MPL 080411
171	! ZFCCUP (klon,KLEV+1) ; CLEAR SKY CLEAN (NO AEROSOL) SW UP FLUXES ! added by OB 211117
172	! ZFLCCDWN(klon,KLEV+1) ; CLEAR SKY CLEAN (NO AEROSOL) LW DWN FLUXES ! added by OB 211117
173	! ZFLCCUP (klon,KLEV+1) ; CLEAR SKY CLEAN (NO AEROSOL) LW UP FLUXES ! added by OB 211117
174
175	!======================================================================
176
177	! ====================================================================
178	! Adapte au modele de chimie INCA par Celine Deandreis & Anne Cozic -- 2009
179	! 1 = ZERO
180	! 2 = AER total
181	! 3 = NAT
182	! 4 = BC
183	! 5 = SO4
184	! 6 = POM
185	! 7 = DUST
186	! 8 = SS
187	! 9 = NO3
188	!
189	! ====================================================================
190
191	! ==============
192	! DECLARATIONS
193	! ==============
194	include "YOETHF.h"
195	include "YOMCST.h"
196	include "clesphys.h"
197
198	! Input arguments
199	REAL, INTENT(in) :: dist
200	REAL, INTENT(in) :: rmu0(KLON), fract(KLON)
201	REAL, INTENT(in) :: paprs(KLON,KLEV+1), pplay(KLON,KLEV)
202	!albedo SB >>>
203	! REAL, INTENT(in) :: alb1(KLON), alb2(KLON), tsol(KLON)
204	REAL, INTENT(in) :: tsol(KLON)
205	REAL, INTENT(in) :: alb_dir(KLON,NSW),alb_dif(KLON,NSW)
206	REAL, INTENT(in) :: SFRWL(6)
207	!albedo SB <<<
208	REAL, INTENT(in) :: t(KLON,KLEV), q(KLON,KLEV)
209
210	REAL, INTENT(in):: wo(:, :, :) ! dimension(KLON,KLEV, 1 or 2)
211	! column-density of ozone in a layer, in kilo-Dobsons
212	! "wo(:, :, 1)" is for the average day-night field,
213	! "wo(:, :, 2)" is for daylight time.
214
215	LOGICAL, INTENT(in) :: ok_ade, ok_aie ! switches whether to use aerosol direct (indirect) effects or not
216	LOGICAL, INTENT(in) :: ok_volcan ! produce volcanic diags (SW/LW heat flux and rate)
217	INTEGER, INTENT(in) :: flag_volc_surfstrat ! allow to impose volcanic cooling rate at surf or heating in strato
218	LOGICAL :: lldebug=.false.
219	INTEGER, INTENT(in) :: flag_aerosol ! takes value 0 (no aerosol) or 1 to 6 (aerosols)
220	INTEGER, INTENT(in) :: flag_aerosol_strat ! use stratospheric aerosols
221	LOGICAL, INTENT(in) :: flag_aer_feedback ! activate aerosol radiative feedback
222	REAL, INTENT(in) :: cldfra(KLON,KLEV), cldemi(KLON,KLEV), cldtaupd(KLON,KLEV)
223	REAL, INTENT(in) :: tau_aero(KLON,KLEV,naero_grp,2) ! aerosol optical properties (see aeropt.F)
224	REAL, INTENT(in) :: piz_aero(KLON,KLEV,naero_grp,2) ! aerosol optical properties (see aeropt.F)
225	REAL, INTENT(in) :: cg_aero(KLON,KLEV,naero_grp,2) ! aerosol optical properties (see aeropt.F)
226	!--OB
227	REAL, INTENT(in) :: tau_aero_sw_rrtm(KLON,KLEV,2,NSW) ! aerosol optical properties RRTM
228	REAL, INTENT(in) :: piz_aero_sw_rrtm(KLON,KLEV,2,NSW) ! aerosol optical properties RRTM
229	REAL, INTENT(in) :: cg_aero_sw_rrtm(KLON,KLEV,2,NSW) ! aerosol optical properties RRTM
230	! AI
231	!--OB fin
232
233	!--C. Kleinschmitt
234	#ifdef CPP_RRTM
235	REAL, INTENT(in) :: tau_aero_lw_rrtm(KLON,KLEV,2,NLW) ! LW aerosol optical properties RRTM
236	#else
237	REAL, INTENT(in) :: tau_aero_lw_rrtm(KLON,KLEV,2,nbands_lw_rrtm)
238	#endif
239	!--C. Kleinschmitt end
240
241	REAL, INTENT(in) :: cldtaupi(KLON,KLEV) ! cloud optical thickness for pre-industrial aerosol concentrations
242	REAL, INTENT(in) :: qsat(klon,klev) ! Variable pour iflag_rrtm=1
243	REAL, INTENT(in) :: flwc(klon,klev) ! Variable pour iflag_rrtm=1
244	REAL, INTENT(in) :: fiwc(klon,klev) ! Variable pour iflag_rrtm=1
245	REAL, INTENT(in) :: ref_liq(klon,klev) ! cloud droplet radius present-day from newmicro
246	REAL, INTENT(in) :: ref_ice(klon,klev) ! ice crystal radius present-day from newmicro
247	REAL, INTENT(in) :: ref_liq_pi(klon,klev) ! cloud droplet radius pre-industrial from newmicro
248	REAL, INTENT(in) :: ref_ice_pi(klon,klev) ! ice crystal radius pre-industrial from newmicro
249
250	CHARACTER(len=512), INTENT(in) :: namelist_ecrad_file
251
252	! Output arguments
253	REAL, INTENT(out) :: heat(KLON,KLEV), cool(KLON,KLEV)
254	REAL, INTENT(out) :: heat0(KLON,KLEV), cool0(KLON,KLEV)
255	REAL, INTENT(out) :: heat_volc(KLON,KLEV), cool_volc(KLON,KLEV) !NL
256	REAL, INTENT(out) :: topsw(KLON), toplw(KLON)
257	REAL, INTENT(out) :: solsw(KLON), sollw(KLON), albpla(KLON), solswfdiff(KLON)
258	REAL, INTENT(out) :: topsw0(KLON), toplw0(KLON), solsw0(KLON), sollw0(KLON)
259	REAL, INTENT(out) :: sollwdown(KLON)
260	REAL, INTENT(out) :: swdn(KLON,kflev+1),swdn0(KLON,kflev+1), swdnc0(KLON,kflev+1)
261	REAL, INTENT(out) :: swup(KLON,kflev+1),swup0(KLON,kflev+1), swupc0(KLON,kflev+1)
262	REAL, INTENT(out) :: lwdn(KLON,kflev+1),lwdn0(KLON,kflev+1), lwdnc0(KLON,kflev+1)
263	REAL, INTENT(out) :: lwup(KLON,kflev+1),lwup0(KLON,kflev+1), lwupc0(KLON,kflev+1)
264	REAL, INTENT(out) :: topswad_aero(KLON), solswad_aero(KLON) ! output: aerosol direct forcing at TOA and surface
265	REAL, INTENT(out) :: topswai_aero(KLON), solswai_aero(KLON) ! output: aerosol indirect forcing atTOA and surface
266	REAL, INTENT(out) :: toplwad_aero(KLON), sollwad_aero(KLON) ! output: LW aerosol direct forcing at TOA and surface
267	REAL, INTENT(out) :: toplwai_aero(KLON), sollwai_aero(KLON) ! output: LW aerosol indirect forcing atTOA and surface
268	REAL, DIMENSION(klon), INTENT(out) :: topswad0_aero
269	REAL, DIMENSION(klon), INTENT(out) :: solswad0_aero
270	REAL, DIMENSION(klon), INTENT(out) :: toplwad0_aero
271	REAL, DIMENSION(klon), INTENT(out) :: sollwad0_aero
272	REAL, DIMENSION(kdlon,9), INTENT(out) :: topsw_aero
273	REAL, DIMENSION(kdlon,9), INTENT(out) :: topsw0_aero
274	REAL, DIMENSION(kdlon,9), INTENT(out) :: solsw_aero
275	REAL, DIMENSION(kdlon,9), INTENT(out) :: solsw0_aero
276	REAL, DIMENSION(kdlon,3), INTENT(out) :: topswcf_aero
277	REAL, DIMENSION(kdlon,3), INTENT(out) :: solswcf_aero
278	REAL, DIMENSION(kdlon,kflev+1), INTENT(out) :: ZSWFT0_i
279	REAL, DIMENSION(kdlon,kflev+1), INTENT(out) :: ZLWFT0_i
280
281	! Local variables
282	REAL(KIND=8) ZFSUP(KDLON,KFLEV+1)
283	REAL(KIND=8) ZFSDN(KDLON,KFLEV+1)
284	REAL(KIND=8) ZFSUP0(KDLON,KFLEV+1)
285	REAL(KIND=8) ZFSDN0(KDLON,KFLEV+1)
286	REAL(KIND=8) ZFSUPC0(KDLON,KFLEV+1)
287	REAL(KIND=8) ZFSDNC0(KDLON,KFLEV+1)
288	REAL(KIND=8) ZFLUP(KDLON,KFLEV+1)
289	REAL(KIND=8) ZFLDN(KDLON,KFLEV+1)
290	REAL(KIND=8) ZFLUP0(KDLON,KFLEV+1)
291	REAL(KIND=8) ZFLDN0(KDLON,KFLEV+1)
292	REAL(KIND=8) ZFLUPC0(KDLON,KFLEV+1)
293	REAL(KIND=8) ZFLDNC0(KDLON,KFLEV+1)
294	REAL(KIND=8) zx_alpha1, zx_alpha2
295	INTEGER k, kk, i, j, iof, nb_gr
296	INTEGER ist,iend,ktdia,kmode
297	REAL(KIND=8) PSCT
298	REAL(KIND=8) PALBD(kdlon,2), PALBP(kdlon,2)
299	! MPL 06.01.09: pour RRTM, creation de PALBD_NEW et PALBP_NEW
300	! avec NSW en deuxieme dimension
301	REAL(KIND=8) PALBD_NEW(kdlon,NSW), PALBP_NEW(kdlon,NSW)
302	REAL(KIND=8) PEMIS(kdlon), PDT0(kdlon), PVIEW(kdlon)
303	REAL(KIND=8) PPSOL(kdlon), PDP(kdlon,KLEV)
304	REAL(KIND=8) PTL(kdlon,kflev+1), PPMB(kdlon,kflev+1)
305	REAL(KIND=8) PTAVE(kdlon,kflev)
306	REAL(KIND=8) PWV(kdlon,kflev), PQS(kdlon,kflev)
307
308	REAL(KIND=8) cloud_cover_sw(klon)
309
310	!!!!!!! Declarations specifiques pour ECRAD !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
311	! AI 02.2021
312	#ifdef CPP_ECRAD
313	! ATTENTION les dimensions klon, kdlon ???
314	! INPUTS
315	REAL, DIMENSION(kdlon,kflev+1) :: ZSWFT0_ii, ZLWFT0_ii
316	REAL(KIND=8) ZEMISW(klon), & ! LW emissivity inside the window region
317	ZEMIS(klon) ! LW emissivity outside the window region
318	REAL(KIND=8) ZGELAM(klon), & ! longitudes en rad
319	ZGEMU(klon) ! sin(latitude)
320	REAL(KIND=8) ZCO2, & ! CO2 mass mixing ratios on full levels
321	ZCH4, & ! CH4 mass mixing ratios on full levels
322	ZN2O, & ! N2O mass mixing ratios on full levels
323	ZNO2, & ! NO2 mass mixing ratios on full levels
324	ZCFC11, & ! CFC11
325	ZCFC12, & ! CFC12
326	ZHCFC22, & ! HCFC22
327	ZCCL4, & ! CCL4
328	ZO2 ! O2
329
330	REAL(KIND=8) ZQ_RAIN(klon,klev), & ! Rain cloud mass mixing ratio (kg/kg) ?
331	ZQ_SNOW(klon,klev) ! Snow cloud mass mixing ratio (kg/kg) ?
332	REAL(KIND=8) ZAEROSOL_OLD(KLON,6,KLEV), & !
333	ZAEROSOL(KLON,KLEV,naero_grp) !
334	! OUTPUTS
335	REAL(KIND=8) ZFLUX_DIR(klon), & ! Direct compt of surf flux into horizontal plane
336	ZFLUX_DIR_CLEAR(klon), & ! CS Direct
337	ZFLUX_DIR_INTO_SUN(klon), & !
338	ZFLUX_UV(klon), & ! UV flux
339	ZFLUX_PAR(klon), & ! photosynthetically active radiation similarly
340	ZFLUX_PAR_CLEAR(klon), & ! CS photosynthetically
341	ZFLUX_SW_DN_TOA(klon), & ! DN SW flux at TOA
342	ZEMIS_OUT(klon) ! effective broadband emissivity
343
344	REAL(KIND=8) ZLWDERIVATIVE(klon,klev+1) ! LW derivatives
345	REAL(KIND=8) ZSWDIFFUSEBAND(klon,NSW), & ! SW DN flux in diffuse albedo band
346	ZSWDIRECTBAND(klon,NSW) ! SW DN flux in direct albedo band
347	REAL(KIND=8) SOLARIRAD
348	REAL(KIND=8) seuilmach
349	! AI 10 mars 22 : Pour les tests Offline
350	logical :: lldebug_for_offline = .false.
351	REAL(KIND=8) solaire_off(klon), &
352	ZCO2_off(klon,klev), &
353	ZCH4_off(klon,klev), & ! CH4 mass mixing ratios on full levels
354	ZN2O_off(klon,klev), & ! N2O mass mixing ratios on full levels
355	ZNO2_off(klon,klev), & ! NO2 mass mixing ratios on full levels
356	ZCFC11_off(klon,klev), & ! CFC11
357	ZCFC12_off(klon,klev), & ! CFC12
358	ZHCFC22_off(klon,klev), & ! HCFC22
359	ZCCL4_off(klon,klev), & ! CCL4
360	ZO2_off(klon,klev) ! O2#endif
361	#endif
362	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
363
364	REAL(kind=8) POZON(kdlon, kflev, size(wo, 3)) ! mass fraction of ozone
365	! "POZON(:, :, 1)" is for the average day-night field,
366	! "POZON(:, :, 2)" is for daylight time.
367	!!!!! Modif MPL 6.01.09 avec RRTM, on passe de 5 a 6
368	REAL(KIND=8) PAER(kdlon,kflev,6)
369	REAL(KIND=8) PCLDLD(kdlon,kflev)
370	REAL(KIND=8) PCLDLU(kdlon,kflev)
371	REAL(KIND=8) PCLDSW(kdlon,kflev)
372	REAL(KIND=8) PTAU(kdlon,2,kflev)
373	REAL(KIND=8) POMEGA(kdlon,2,kflev)
374	REAL(KIND=8) PCG(kdlon,2,kflev)
375	REAL(KIND=8) zfract(kdlon), zrmu0(kdlon), zdist
376	REAL(KIND=8) zheat(kdlon,kflev), zcool(kdlon,kflev)
377	REAL(KIND=8) zheat0(kdlon,kflev), zcool0(kdlon,kflev)
378	REAL(KIND=8) zheat_volc(kdlon,kflev), zcool_volc(kdlon,kflev) !NL
379	REAL(KIND=8) ztopsw(kdlon), ztoplw(kdlon)
380	REAL(KIND=8) zsolsw(kdlon), zsollw(kdlon), zalbpla(kdlon), zsolswfdiff(kdlon)
381	REAL(KIND=8) zsollwdown(kdlon)
382	REAL(KIND=8) ztopsw0(kdlon), ztoplw0(kdlon)
383	REAL(KIND=8) zsolsw0(kdlon), zsollw0(kdlon)
384	REAL(KIND=8) zznormcp
385	REAL(KIND=8) tauaero(kdlon,kflev,naero_grp,2) ! aer opt properties
386	REAL(KIND=8) pizaero(kdlon,kflev,naero_grp,2)
387	REAL(KIND=8) cgaero(kdlon,kflev,naero_grp,2)
388	REAL(KIND=8) PTAUA(kdlon,2,kflev) ! present-day value of cloud opt thickness (PTAU is pre-industrial value), local use
389	REAL(KIND=8) POMEGAA(kdlon,2,kflev) ! dito for single scatt albedo
390	REAL(KIND=8) ztopswadaero(kdlon), zsolswadaero(kdlon) ! Aerosol direct forcing at TOAand surface
391	REAL(KIND=8) ztopswad0aero(kdlon), zsolswad0aero(kdlon) ! Aerosol direct forcing at TOAand surface
392	REAL(KIND=8) ztopswaiaero(kdlon), zsolswaiaero(kdlon) ! dito, indirect
393	!--NL
394	REAL(KIND=8) zswadaero(kdlon,kflev+1) ! SW Aerosol direct forcing
395	REAL(KIND=8) zlwadaero(kdlon,kflev+1) ! LW Aerosol direct forcing
396	REAL(KIND=8) volmip_solsw(kdlon) ! SW clear sky in the case of VOLMIP
397	!-LW by CK
398	REAL(KIND=8) ztoplwadaero(kdlon), zsollwadaero(kdlon) ! LW Aerosol direct forcing at TOAand surface
399	REAL(KIND=8) ztoplwad0aero(kdlon), zsollwad0aero(kdlon) ! LW Aerosol direct forcing at TOAand surface
400	REAL(KIND=8) ztoplwaiaero(kdlon), zsollwaiaero(kdlon) ! dito, indirect
401	!-end
402	REAL(KIND=8) ztopsw_aero(kdlon,9), ztopsw0_aero(kdlon,9)
403	REAL(KIND=8) zsolsw_aero(kdlon,9), zsolsw0_aero(kdlon,9)
404	REAL(KIND=8) ztopswcf_aero(kdlon,3), zsolswcf_aero(kdlon,3)
405	! real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2 deje declare dans physiq.F MPL 20130618
406	!MPL input supplementaires pour RECMWFL
407	! flwc, fiwc = Liquid Water Content & Ice Water Content (kg/kg)
408	REAL(KIND=8) GEMU(klon)
409	!MPL input RECMWFL:
410	! Tableaux aux niveaux inverses pour respecter convention Arpege
411	REAL(KIND=8) ref_liq_i(klon,klev) ! cloud droplet radius present-day from newmicro (inverted)
412	REAL(KIND=8) ref_ice_i(klon,klev) ! ice crystal radius present-day from newmicro (inverted)
413	!--OB
414	REAL(KIND=8) ref_liq_pi_i(klon,klev) ! cloud droplet radius pre-industrial from newmicro (inverted)
415	REAL(KIND=8) ref_ice_pi_i(klon,klev) ! ice crystal radius pre-industrial from newmicro (inverted)
416	!--end OB
417	REAL(KIND=8) paprs_i(klon,klev+1)
418	REAL(KIND=8) pplay_i(klon,klev)
419	REAL(KIND=8) cldfra_i(klon,klev)
420	REAL(KIND=8) POZON_i(kdlon,kflev, size(wo, 3)) ! mass fraction of ozone
421	! "POZON(:, :, 1)" is for the average day-night field,
422	! "POZON(:, :, 2)" is for daylight time.
423	!!!!! Modif MPL 6.01.09 avec RRTM, on passe de 5 a 6
424	REAL(KIND=8) PAER_i(kdlon,kflev,6)
425	REAL(KIND=8) PDP_i(klon,klev)
426	REAL(KIND=8) t_i(klon,klev),q_i(klon,klev),qsat_i(klon,klev)
427	REAL(KIND=8) flwc_i(klon,klev),fiwc_i(klon,klev)
428	!MPL output RECMWFL:
429	REAL(KIND=8) ZEMTD (klon,klev+1),ZEMTD_i (klon,klev+1)
430	REAL(KIND=8) ZEMTU (klon,klev+1),ZEMTU_i (klon,klev+1)
431	REAL(KIND=8) ZTRSO (klon,klev+1),ZTRSO_i (klon,klev+1)
432	REAL(KIND=8) ZTH (klon,klev+1),ZTH_i (klon,klev+1)
433	REAL(KIND=8) ZCTRSO(klon,2)
434	REAL(KIND=8) ZCEMTR(klon,2)
435	REAL(KIND=8) ZTRSOD(klon)
436	REAL(KIND=8) ZLWFC (klon,2)
437	REAL(KIND=8) ZLWFT (klon,klev+1),ZLWFT_i (klon,klev+1)
438	REAL(KIND=8) ZSWFC (klon,2)
439	REAL(KIND=8) ZSWFT (klon,klev+1),ZSWFT_i (klon,klev+1)
440	REAL(KIND=8) ZFLUCDWN_i(klon,klev+1),ZFLUCUP_i(klon,klev+1)
441	REAL(KIND=8) PPIZA_TOT(klon,klev,NSW)
442	REAL(KIND=8) PCGA_TOT(klon,klev,NSW)
443	REAL(KIND=8) PTAU_TOT(klon,klev,NSW)
444	REAL(KIND=8) PPIZA_NAT(klon,klev,NSW)
445	REAL(KIND=8) PCGA_NAT(klon,klev,NSW)
446	REAL(KIND=8) PTAU_NAT(klon,klev,NSW)
447	#ifdef CPP_RRTM
448	REAL(KIND=8) PTAU_LW_TOT(klon,klev,NLW)
449	REAL(KIND=8) PTAU_LW_NAT(klon,klev,NLW)
450	#endif
451	REAL(KIND=8) PSFSWDIR(klon,NSW)
452	REAL(KIND=8) PSFSWDIF(klon,NSW)
453	REAL(KIND=8) PFSDNN(klon)
454	REAL(KIND=8) PFSDNV(klon)
455	!MPL On ne redefinit pas les tableaux ZFLUX,ZFLUC,
456	!MPL ZFSDWN,ZFCDWN,ZFSUP,ZFCUP car ils existent deja
457	!MPL sous les noms de ZFLDN,ZFLDN0,ZFLUP,ZFLUP0,
458	!MPL ZFSDN,ZFSDN0,ZFSUP,ZFSUP0
459	REAL(KIND=8) ZFLUX_i (klon,2,klev+1)
460	REAL(KIND=8) ZFLUC_i (klon,2,klev+1)
461	REAL(KIND=8) ZFSDWN_i (klon,klev+1)
462	REAL(KIND=8) ZFCDWN_i (klon,klev+1)
463	REAL(KIND=8) ZFCCDWN_i (klon,klev+1)
464	REAL(KIND=8) ZFSUP_i (klon,klev+1)
465	REAL(KIND=8) ZFCUP_i (klon,klev+1)
466	REAL(KIND=8) ZFCCUP_i (klon,klev+1)
467	REAL(KIND=8) ZFLCCDWN_i (klon,klev+1)
468	REAL(KIND=8) ZFLCCUP_i (klon,klev+1)
469	! 3 lignes suivantes a activer pour CCMVAL (MPL 20100412)
470	! REAL(KIND=8) RSUN(3,2)
471	! REAL(KIND=8) SUN(3)
472	! REAL(KIND=8) SUN_FRACT(2)
473	REAL, PARAMETER:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2
474	CHARACTER (LEN=80) :: abort_message
475	CHARACTER (LEN=80) :: modname='radlwsw_m'
476
477	REAL zdir, zdif
478
479	! ========= INITIALISATIONS ==============================================
480	IF (lldebug) THEN
481	print*,'Entree dans radlwsw '
482	print,'********** INITIALISATIONS ***************************'
483	print*,'klon, kdlon, klev, kflev =',klon, kdlon, klev, kflev
484	ENDIF
485
486	CALL assert(size(wo, 1) == klon, size(wo, 2) == klev, "radlwsw wo")
487
488	ist=1
489	iend=klon
490	ktdia=1
491	kmode=ist
492	! Aeros
493	tauaero(:,:,:,:)=0.
494	pizaero(:,:,:,:)=0.
495	cgaero(:,:,:,:)=0.
496	! lldebug=.FALSE.
497
498	ztopsw_aero(:,:) = 0. !ym missing init : warning : not initialized in SW_AEROAR4
499	ztopsw0_aero(:,:) = 0. !ym missing init : warning : not initialized in SW_AEROAR4
500	zsolsw_aero(:,:) = 0. !ym missing init : warning : not initialized in SW_AEROAR4
501	zsolsw0_aero(:,:) = 0. !ym missing init : warning : not initialized in SW_AEROAR4
502
503	ZTOPSWADAERO(:) = 0. !ym missing init
504	ZSOLSWADAERO(:) = 0. !ym missing init
505	ZTOPSWAD0AERO(:) = 0. !ym missing init
506	ZSOLSWAD0AERO(:) = 0. !ym missing init
507	ZTOPSWAIAERO(:) = 0. !ym missing init
508	ZSOLSWAIAERO(:) = 0. !ym missing init
509	ZTOPSWCF_AERO(:,:)= 0.!ym missing init
510	ZSOLSWCF_AERO(:,:) =0. !ym missing init
511
512	!
513	! AI 02.2021
514	#ifdef CPP_ECRAD
515	ZEMIS = 1.0
516	ZEMISW = 1.0
517	ZGELAM = longitude
518	ZGEMU = sin(latitude)
519	ZCO2 = RCO2
520	ZCH4 = RCH4
521	ZN2O = RN2O
522	ZNO2 = 0.0
523	ZCFC11 = RCFC11
524	ZCFC12 = RCFC12
525	ZHCFC22 = 0.0
526	ZO2 = 0.0
527	ZCCL4 = 0.0
528	ZQ_RAIN = 0.0
529	ZQ_SNOW = 0.0
530	ZAEROSOL_OLD = 0.0
531	ZAEROSOL = 0.0
532	seuilmach=tiny(seuilmach)
533	#endif
534
535	!-------------------------------------------
536	nb_gr = KLON / kdlon
537	IF (nb_gr*kdlon .NE. KLON) THEN
538	PRINT*, "kdlon mauvais:", KLON, kdlon, nb_gr
539	call abort_physic("radlwsw", "", 1)
540	ENDIF
541	IF (kflev .NE. KLEV) THEN
542	PRINT*, "kflev differe de KLEV, kflev, KLEV"
543	call abort_physic("radlwsw", "", 1)
544	ENDIF
545	!-------------------------------------------
546	DO k = 1, KLEV
547	DO i = 1, KLON
548	heat(i,k)=0.
549	cool(i,k)=0.
550	heat_volc(i,k)=0. !NL
551	cool_volc(i,k)=0. !NL
552	heat0(i,k)=0.
553	cool0(i,k)=0.
554	ENDDO
555	ENDDO
556	!
557	zdist = dist
558	!
559	PSCT = solaire/zdist/zdist
560
561	IF (type_trac == 'repr') THEN
562	#ifdef REPROBUS
563	IF (iflag_rrtm==0) THEN
564	IF (ok_SUNTIME) PSCT = solaireTIME/zdist/zdist
565	print,'Constante solaire: ',PSCTzdist*zdist
566	ENDIF
567	#endif
568	ENDIF
569
570	IF (lldebug) THEN
571	print,'************* Debut boucle de 1 a ', nb_gr
572	ENDIF
573
574	DO j = 1, nb_gr
575	iof = kdlon*(j-1)
576	DO i = 1, kdlon
577	zfract(i) = fract(iof+i)
578	zrmu0(i) = rmu0(iof+i)
579
580
581	IF (iflag_rrtm==0) THEN
582	! Albedo
583	PALBD(i,1)=alb_dif(iof+i,1)
584	PALBD(i,2)=alb_dif(iof+i,2)
585	PALBP(i,1)=alb_dir(iof+i,1)
586	PALBP(i,2)=alb_dir(iof+i,2)
587	! AI 02.2021 cas iflag_rrtm=1 et 2
588	ELSEIF (iflag_rrtm==1.OR.iflag_rrtm==2) THEN
589	DO kk=1,NSW
590	PALBD_NEW(i,kk)=alb_dif(iof+i,kk)
591	PALBP_NEW(i,kk)=alb_dir(iof+i,kk)
592	ENDDO
593	!
594	ENDIF
595	!albedo SB <<<
596
597	PEMIS(i) = 1.0 !!!!! A REVOIR (MPL)
598	PVIEW(i) = 1.66
599	PPSOL(i) = paprs(iof+i,1)
600	zx_alpha1 = (paprs(iof+i,1)-pplay(iof+i,2))/(pplay(iof+i,1)-pplay(iof+i,2))
601	zx_alpha2 = 1.0 - zx_alpha1
602	PTL(i,1) = t(iof+i,1) * zx_alpha1 + t(iof+i,2) * zx_alpha2
603	PTL(i,KLEV+1) = t(iof+i,KLEV)
604	PDT0(i) = tsol(iof+i) - PTL(i,1)
605	ENDDO
606	DO k = 2, kflev
607	DO i = 1, kdlon
608	PTL(i,k) = (t(iof+i,k)+t(iof+i,k-1))*0.5
609	ENDDO
610	ENDDO
611	DO k = 1, kflev
612	DO i = 1, kdlon
613	PDP(i,k) = paprs(iof+i,k)-paprs(iof+i,k+1)
614	PTAVE(i,k) = t(iof+i,k)
615	PWV(i,k) = MAX (q(iof+i,k), 1.0e-12)
616	PQS(i,k) = PWV(i,k)
617	! Confert from column density of ozone in a cell, in kDU, to a mass fraction
618	POZON(i,k, :) = wo(iof+i, k, :) * RG * dobson_u * 1e3 &
619	/ (paprs(iof+i, k) - paprs(iof+i, k+1))
620	! A activer pour CCMVAL on prend l'ozone impose (MPL 07042010)
621	! POZON(i,k,:) = wo(i,k,:)
622	! print *,'RADLWSW: POZON',k, POZON(i,k,1)
623	PCLDLD(i,k) = cldfra(iof+i,k)*cldemi(iof+i,k)
624	PCLDLU(i,k) = cldfra(iof+i,k)*cldemi(iof+i,k)
625	PCLDSW(i,k) = cldfra(iof+i,k)
626	PTAU(i,1,k) = MAX(cldtaupi(iof+i,k), 1.0e-05)! 1e-12 serait instable
627	PTAU(i,2,k) = MAX(cldtaupi(iof+i,k), 1.0e-05)! pour 32-bit machines
628	POMEGA(i,1,k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAU(i,1,k))
629	POMEGA(i,2,k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAU(i,2,k))
630	PCG(i,1,k) = 0.865
631	PCG(i,2,k) = 0.910
632	!-
633	! Introduced for aerosol indirect forcings.
634	! The following values use the cloud optical thickness calculated from
635	! present-day aerosol concentrations whereas the quantities without the
636	! "A" at the end are for pre-industial (natural-only) aerosol concentrations
637	!
638	PTAUA(i,1,k) = MAX(cldtaupd(iof+i,k), 1.0e-05)! 1e-12 serait instable
639	PTAUA(i,2,k) = MAX(cldtaupd(iof+i,k), 1.0e-05)! pour 32-bit machines
640	POMEGAA(i,1,k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAUA(i,1,k))
641	POMEGAA(i,2,k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAUA(i,2,k))
642	ENDDO
643	ENDDO
644
645	IF (type_trac == 'repr') THEN
646	#ifdef REPROBUS
647	ndimozon = size(wo, 3)
648	CALL RAD_INTERACTIF(POZON,iof)
649	#endif
650	ENDIF
651	!
652	DO k = 1, kflev+1
653	DO i = 1, kdlon
654	PPMB(i,k) = paprs(iof+i,k)/100.0
655	ENDDO
656	ENDDO
657	!
658	!!!!! Modif MPL 6.01.09 avec RRTM, on passe de 5 a 6
659	DO kk = 1, 6
660	DO k = 1, kflev
661	DO i = 1, kdlon
662	PAER(i,k,kk) = 1.0E-15 !!!!! A REVOIR (MPL)
663	ENDDO
664	ENDDO
665	ENDDO
666	DO k = 1, kflev
667	DO i = 1, kdlon
668	tauaero(i,k,:,1)=tau_aero(iof+i,k,:,1)
669	pizaero(i,k,:,1)=piz_aero(iof+i,k,:,1)
670	cgaero(i,k,:,1) =cg_aero(iof+i,k,:,1)
671	tauaero(i,k,:,2)=tau_aero(iof+i,k,:,2)
672	pizaero(i,k,:,2)=piz_aero(iof+i,k,:,2)
673	cgaero(i,k,:,2) =cg_aero(iof+i,k,:,2)
674	ENDDO
675	ENDDO
676	!
677	!===== iflag_rrtm ================================================
678	!
679	IF (iflag_rrtm == 0) THEN !!!! remettre 0 juste pour tester l'ancien rayt via rrtm
680	!
681	!--- Mise a zero des tableaux output du rayonnement LW-AR4 ----------
682	DO k = 1, kflev+1
683	DO i = 1, kdlon
684	! print *,'RADLWSW: boucle mise a zero i k',i,k
685	ZFLUP(i,k)=0.
686	ZFLDN(i,k)=0.
687	ZFLUP0(i,k)=0.
688	ZFLDN0(i,k)=0.
689	ZLWFT0_i(i,k)=0.
690	ZFLUCUP_i(i,k)=0.
691	ZFLUCDWN_i(i,k)=0.
692	ENDDO
693	ENDDO
694	DO k = 1, kflev
695	DO i = 1, kdlon
696	zcool(i,k)=0.
697	zcool_volc(i,k)=0. !NL
698	zcool0(i,k)=0.
699	ENDDO
700	ENDDO
701	DO i = 1, kdlon
702	ztoplw(i)=0.
703	zsollw(i)=0.
704	ztoplw0(i)=0.
705	zsollw0(i)=0.
706	zsollwdown(i)=0.
707	ENDDO
708	! Old radiation scheme, used for AR4 runs
709	! average day-night ozone for longwave
710	CALL LW_LMDAR4(&
711	PPMB, PDP,&
712	PPSOL,PDT0,PEMIS,&
713	PTL, PTAVE, PWV, POZON(:, :, 1), PAER,&
714	PCLDLD,PCLDLU,&
715	PVIEW,&
716	zcool, zcool0,&
717	ztoplw,zsollw,ztoplw0,zsollw0,&
718	zsollwdown,&
719	ZFLUP, ZFLDN, ZFLUP0,ZFLDN0)
720	!----- Mise a zero des tableaux output du rayonnement SW-AR4
721	DO k = 1, kflev+1
722	DO i = 1, kdlon
723	ZFSUP(i,k)=0.
724	ZFSDN(i,k)=0.
725	ZFSUP0(i,k)=0.
726	ZFSDN0(i,k)=0.
727	ZFSUPC0(i,k)=0.
728	ZFSDNC0(i,k)=0.
729	ZFLUPC0(i,k)=0.
730	ZFLDNC0(i,k)=0.
731	ZSWFT0_i(i,k)=0.
732	ZFCUP_i(i,k)=0.
733	ZFCDWN_i(i,k)=0.
734	ZFCCUP_i(i,k)=0.
735	ZFCCDWN_i(i,k)=0.
736	ZFLCCUP_i(i,k)=0.
737	ZFLCCDWN_i(i,k)=0.
738	zswadaero(i,k)=0. !--NL
739	ENDDO
740	ENDDO
741	DO k = 1, kflev
742	DO i = 1, kdlon
743	zheat(i,k)=0.
744	zheat_volc(i,k)=0.
745	zheat0(i,k)=0.
746	ENDDO
747	ENDDO
748	DO i = 1, kdlon
749	zalbpla(i)=0.
750	ztopsw(i)=0.
751	zsolsw(i)=0.
752	ztopsw0(i)=0.
753	zsolsw0(i)=0.
754	ztopswadaero(i)=0.
755	zsolswadaero(i)=0.
756	ztopswaiaero(i)=0.
757	zsolswaiaero(i)=0.
758	ENDDO
759
760	!--fraction of diffuse radiation in surface SW downward radiation
761	!--not computed with old radiation scheme
762	zsolswfdiff(:) = -999.999
763
764	! print *,'Avant SW_LMDAR4: PSCT zrmu0 zfract',PSCT, zrmu0, zfract
765	! daylight ozone, if we have it, for short wave
766	CALL SW_AEROAR4(PSCT, zrmu0, zfract,&
767	PPMB, PDP,&
768	PPSOL, PALBD, PALBP,&
769	PTAVE, PWV, PQS, POZON(:, :, size(wo, 3)), PAER,&
770	PCLDSW, PTAU, POMEGA, PCG,&
771	zheat, zheat0,&
772	zalbpla,ztopsw,zsolsw,ztopsw0,zsolsw0,&
773	ZFSUP,ZFSDN,ZFSUP0,ZFSDN0,&
774	tauaero, pizaero, cgaero, &
775	PTAUA, POMEGAA,&
776	ztopswadaero,zsolswadaero,&
777	ztopswad0aero,zsolswad0aero,&
778	ztopswaiaero,zsolswaiaero, &
779	ztopsw_aero,ztopsw0_aero,&
780	zsolsw_aero,zsolsw0_aero,&
781	ztopswcf_aero,zsolswcf_aero, &
782	ok_ade, ok_aie, flag_aerosol,flag_aerosol_strat)
783
784	ZSWFT0_i(:,:) = ZFSDN0(:,:)-ZFSUP0(:,:)
785	ZLWFT0_i(:,:) =-ZFLDN0(:,:)-ZFLUP0(:,:)
786
787	DO i=1,kdlon
788	DO k=1,kflev+1
789	lwdn0 ( iof+i,k) = ZFLDN0 ( i,k)
790	lwdn ( iof+i,k) = ZFLDN ( i,k)
791	lwup0 ( iof+i,k) = ZFLUP0 ( i,k)
792	lwup ( iof+i,k) = ZFLUP ( i,k)
793	swdn0 ( iof+i,k) = ZFSDN0 ( i,k)
794	swdn ( iof+i,k) = ZFSDN ( i,k)
795	swup0 ( iof+i,k) = ZFSUP0 ( i,k)
796	swup ( iof+i,k) = ZFSUP ( i,k)
797	ENDDO
798	ENDDO
799	!
800	ELSE IF (iflag_rrtm == 1) then
801	#ifdef CPP_RRTM
802	! if (prt_level.gt.10)write(lunout,*)'CPP_RRTM=.T.'
803	!===== iflag_rrtm=1, on passe dans SW via RECMWFL ===============
804
805	DO k = 1, kflev+1
806	DO i = 1, kdlon
807	ZEMTD_i(i,k)=0.
808	ZEMTU_i(i,k)=0.
809	ZTRSO_i(i,k)=0.
810	ZTH_i(i,k)=0.
811	ZLWFT_i(i,k)=0.
812	ZSWFT_i(i,k)=0.
813	ZFLUX_i(i,1,k)=0.
814	ZFLUX_i(i,2,k)=0.
815	ZFLUC_i(i,1,k)=0.
816	ZFLUC_i(i,2,k)=0.
817	ZFSDWN_i(i,k)=0.
818	ZFCDWN_i(i,k)=0.
819	ZFCCDWN_i(i,k)=0.
820	ZFSUP_i(i,k)=0.
821	ZFCUP_i(i,k)=0.
822	ZFCCUP_i(i,k)=0.
823	ZFLCCDWN_i(i,k)=0.
824	ZFLCCUP_i(i,k)=0.
825	ENDDO
826	ENDDO
827	!
828	!--OB
829	!--aerosol TOT - anthropogenic+natural - index 2
830	!--aerosol NAT - natural only - index 1
831	!
832	DO i = 1, kdlon
833	DO k = 1, kflev
834	DO kk=1, NSW
835	!
836	PTAU_TOT(i,kflev+1-k,kk)=tau_aero_sw_rrtm(i,k,2,kk)
837	PPIZA_TOT(i,kflev+1-k,kk)=piz_aero_sw_rrtm(i,k,2,kk)
838	PCGA_TOT(i,kflev+1-k,kk)=cg_aero_sw_rrtm(i,k,2,kk)
839	!
840	PTAU_NAT(i,kflev+1-k,kk)=tau_aero_sw_rrtm(i,k,1,kk)
841	PPIZA_NAT(i,kflev+1-k,kk)=piz_aero_sw_rrtm(i,k,1,kk)
842	PCGA_NAT(i,kflev+1-k,kk)=cg_aero_sw_rrtm(i,k,1,kk)
843	!
844	ENDDO
845	ENDDO
846	ENDDO
847	!-end OB
848	!
849	!--C. Kleinschmitt
850	!--aerosol TOT - anthropogenic+natural - index 2
851	!--aerosol NAT - natural only - index 1
852	!
853	DO i = 1, kdlon
854	DO k = 1, kflev
855	DO kk=1, NLW
856	!
857	PTAU_LW_TOT(i,kflev+1-k,kk)=tau_aero_lw_rrtm(i,k,2,kk)
858	PTAU_LW_NAT(i,kflev+1-k,kk)=tau_aero_lw_rrtm(i,k,1,kk)
859	!
860	ENDDO
861	ENDDO
862	ENDDO
863	!-end C. Kleinschmitt
864	!
865	DO i = 1, kdlon
866	ZCTRSO(i,1)=0.
867	ZCTRSO(i,2)=0.
868	ZCEMTR(i,1)=0.
869	ZCEMTR(i,2)=0.
870	ZTRSOD(i)=0.
871	ZLWFC(i,1)=0.
872	ZLWFC(i,2)=0.
873	ZSWFC(i,1)=0.
874	ZSWFC(i,2)=0.
875	PFSDNN(i)=0.
876	PFSDNV(i)=0.
877	DO kk = 1, NSW
878	PSFSWDIR(i,kk)=0.
879	PSFSWDIF(i,kk)=0.
880	ENDDO
881	ENDDO
882	!----- Fin des mises a zero des tableaux output de RECMWF -------------------
883	! GEMU(1:klon)=sin(rlatd(1:klon))
884	! On met les donnees dans l'ordre des niveaux arpege
885	paprs_i(:,1)=paprs(:,klev+1)
886	DO k=1,klev
887	paprs_i(1:klon,k+1) =paprs(1:klon,klev+1-k)
888	pplay_i(1:klon,k) =pplay(1:klon,klev+1-k)
889	cldfra_i(1:klon,k) =cldfra(1:klon,klev+1-k)
890	PDP_i(1:klon,k) =PDP(1:klon,klev+1-k)
891	t_i(1:klon,k) =t(1:klon,klev+1-k)
892	q_i(1:klon,k) =q(1:klon,klev+1-k)
893	qsat_i(1:klon,k) =qsat(1:klon,klev+1-k)
894	flwc_i(1:klon,k) =flwc(1:klon,klev+1-k)
895	fiwc_i(1:klon,k) =fiwc(1:klon,klev+1-k)
896	ref_liq_i(1:klon,k) =ref_liq(1:klon,klev+1-k)
897	ref_ice_i(1:klon,k) =ref_ice(1:klon,klev+1-k)
898	!-OB
899	ref_liq_pi_i(1:klon,k) =ref_liq_pi(1:klon,klev+1-k)
900	ref_ice_pi_i(1:klon,k) =ref_ice_pi(1:klon,klev+1-k)
901	ENDDO
902	DO k=1,kflev
903	POZON_i(1:klon,k,:)=POZON(1:klon,kflev+1-k,:)
904	!!! POZON_i(1:klon,k)=POZON(1:klon,k) !!! on laisse 1=sol et klev=top
905	! print *,'Juste avant RECMWFL: k tsol temp',k,tsol,t(1,k)
906	!!!!!!! Modif MPL 6.01.09 avec RRTM, on passe de 5 a 6
907	DO i=1,6
908	PAER_i(1:klon,k,i)=PAER(1:klon,kflev+1-k,i)
909	ENDDO
910	ENDDO
911
912	! print *,'RADLWSW: avant RECMWFL, RI0,rmu0=',solaire,rmu0
913
914	! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
915	! La version ARPEGE1D utilise differentes valeurs de la constante
916	! solaire suivant le rayonnement utilise.
917	! A controler ...
918	! SOLAR FLUX AT THE TOP (/YOMPHY3/)
919	! introduce season correction
920	!--------------------------------------
921	! RII0 = RIP0
922	! IF(LRAYFM)
923	! RII0 = RIP0M ! =rip0m if Morcrette non-each time step call.
924	! IF(LRAYFM15)
925	! RII0 = RIP0M15 ! =rip0m if Morcrette non-each time step call.
926	RII0=solaire/zdist/zdist
927	! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
928	! Ancien appel a RECMWF (celui du cy25)
929	! CALL RECMWF (ist , iend, klon , ktdia , klev , kmode ,
930	! s PALBD , PALBP , paprs_i , pplay_i , RCO2 , cldfra_i,
931	! s POZON_i , PAER_i , PDP_i , PEMIS , GEMU , rmu0,
932	! s q_i , qsat_i , fiwc_i , flwc_i , zmasq , t_i ,tsol,
933	! s ZEMTD_i , ZEMTU_i , ZTRSO_i ,
934	! s ZTH_i , ZCTRSO , ZCEMTR , ZTRSOD ,
935	! s ZLWFC , ZLWFT_i , ZSWFC , ZSWFT_i ,
936	! s ZFLUX_i , ZFLUC_i , ZFSDWN_i, ZFSUP_i , ZFCDWN_i,ZFCUP_i)
937	! s 'RECMWF ')
938	!
939	IF (lldebug) THEN
940	CALL writefield_phy('paprs_i',paprs_i,klev+1)
941	CALL writefield_phy('pplay_i',pplay_i,klev)
942	CALL writefield_phy('cldfra_i',cldfra_i,klev)
943	CALL writefield_phy('pozon_i',POZON_i,klev)
944	CALL writefield_phy('paer_i',PAER_i,klev)
945	CALL writefield_phy('pdp_i',PDP_i,klev)
946	CALL writefield_phy('q_i',q_i,klev)
947	CALL writefield_phy('qsat_i',qsat_i,klev)
948	CALL writefield_phy('fiwc_i',fiwc_i,klev)
949	CALL writefield_phy('flwc_i',flwc_i,klev)
950	CALL writefield_phy('t_i',t_i,klev)
951	CALL writefield_phy('palbd_new',PALBD_NEW,NSW)
952	CALL writefield_phy('palbp_new',PALBP_NEW,NSW)
953	ENDIF
954
955	! Nouvel appel a RECMWF (celui du cy32t0)
956	CALL RECMWF_AERO (ist , iend, klon , ktdia , klev , kmode ,&
957	PALBD_NEW,PALBP_NEW, paprs_i , pplay_i , RCO2 , cldfra_i,&
958	POZON_i , PAER_i , PDP_i , PEMIS , rmu0 ,&
959	q_i , qsat_i , fiwc_i , flwc_i , zmasq , t_i ,tsol,&
960	ref_liq_i, ref_ice_i, &
961	ref_liq_pi_i, ref_ice_pi_i, & ! rajoute par OB pour diagnostiquer effet indirect
962	ZEMTD_i , ZEMTU_i , ZTRSO_i ,&
963	ZTH_i , ZCTRSO , ZCEMTR , ZTRSOD ,&
964	ZLWFC , ZLWFT_i , ZSWFC , ZSWFT_i ,&
965	PSFSWDIR , PSFSWDIF, PFSDNN , PFSDNV ,&
966	PPIZA_TOT, PCGA_TOT,PTAU_TOT,&
967	PPIZA_NAT, PCGA_NAT,PTAU_NAT, & ! rajoute par OB pour diagnostiquer effet direct
968	PTAU_LW_TOT, PTAU_LW_NAT, & ! rajoute par C. Kleinschmitt
969	ZFLUX_i , ZFLUC_i ,&
970	ZFSDWN_i , ZFSUP_i , ZFCDWN_i, ZFCUP_i, ZFCCDWN_i, ZFCCUP_i, ZFLCCDWN_i, ZFLCCUP_i, &
971	ZTOPSWADAERO,ZSOLSWADAERO,& ! rajoute par OB pour diagnostics
972	ZTOPSWAD0AERO,ZSOLSWAD0AERO,&
973	ZTOPSWAIAERO,ZSOLSWAIAERO, &
974	ZTOPSWCF_AERO,ZSOLSWCF_AERO, &
975	ZSWADAERO, & !--NL
976	ZTOPLWADAERO,ZSOLLWADAERO,& ! rajoute par C. Kleinscmitt pour LW diagnostics
977	ZTOPLWAD0AERO,ZSOLLWAD0AERO,&
978	ZTOPLWAIAERO,ZSOLLWAIAERO, &
979	ZLWADAERO, & !--NL
980	volmip_solsw, flag_volc_surfstrat, & !--VOLMIP
981	ok_ade, ok_aie, ok_volcan, flag_aerosol,flag_aerosol_strat, flag_aer_feedback) ! flags aerosols
982
983	!--OB diagnostics
984	! & PTOPSWAIAERO,PSOLSWAIAERO,&
985	! & PTOPSWCFAERO,PSOLSWCFAERO,&
986	! & PSWADAERO,& !--NL
987	!!--LW diagnostics CK
988	! & PTOPLWADAERO,PSOLLWADAERO,&
989	! & PTOPLWAD0AERO,PSOLLWAD0AERO,&
990	! & PTOPLWAIAERO,PSOLLWAIAERO,&
991	! & PLWADAERO,& !--NL
992	!!..end
993	! & ok_ade, ok_aie, ok_volcan, flag_aerosol,flag_aerosol_strat,&
994	! & flag_aer_feedback)
995
996
997	! print *,'RADLWSW: apres RECMWF'
998	IF (lldebug) THEN
999	CALL writefield_phy('zemtd_i',ZEMTD_i,klev+1)
1000	CALL writefield_phy('zemtu_i',ZEMTU_i,klev+1)
1001	CALL writefield_phy('ztrso_i',ZTRSO_i,klev+1)
1002	CALL writefield_phy('zth_i',ZTH_i,klev+1)
1003	CALL writefield_phy('zctrso',ZCTRSO,2)
1004	CALL writefield_phy('zcemtr',ZCEMTR,2)
1005	CALL writefield_phy('ztrsod',ZTRSOD,1)
1006	CALL writefield_phy('zlwfc',ZLWFC,2)
1007	CALL writefield_phy('zlwft_i',ZLWFT_i,klev+1)
1008	CALL writefield_phy('zswfc',ZSWFC,2)
1009	CALL writefield_phy('zswft_i',ZSWFT_i,klev+1)
1010	CALL writefield_phy('psfswdir',PSFSWDIR,6)
1011	CALL writefield_phy('psfswdif',PSFSWDIF,6)
1012	CALL writefield_phy('pfsdnn',PFSDNN,1)
1013	CALL writefield_phy('pfsdnv',PFSDNV,1)
1014	CALL writefield_phy('ppiza_dst',PPIZA_TOT,klev)
1015	CALL writefield_phy('pcga_dst',PCGA_TOT,klev)
1016	CALL writefield_phy('ptaurel_dst',PTAU_TOT,klev)
1017	CALL writefield_phy('zflux_i',ZFLUX_i,klev+1)
1018	CALL writefield_phy('zfluc_i',ZFLUC_i,klev+1)
1019	CALL writefield_phy('zfsdwn_i',ZFSDWN_i,klev+1)
1020	CALL writefield_phy('zfsup_i',ZFSUP_i,klev+1)
1021	CALL writefield_phy('zfcdwn_i',ZFCDWN_i,klev+1)
1022	CALL writefield_phy('zfcup_i',ZFCUP_i,klev+1)
1023	ENDIF
1024
1025	! ---------
1026	! ---------
1027	! On retablit l'ordre des niveaux lmd pour les tableaux de sortie
1028	! D autre part, on multiplie les resultats SW par fract pour etre coherent
1029	! avec l ancien rayonnement AR4. Si nuit, fract=0 donc pas de
1030	! rayonnement SW. (MPL 260609)
1031	DO k=0,klev
1032	DO i=1,klon
1033	ZEMTD(i,k+1) = ZEMTD_i(i,k+1)
1034	ZEMTU(i,k+1) = ZEMTU_i(i,k+1)
1035	ZTRSO(i,k+1) = ZTRSO_i(i,k+1)
1036	ZTH(i,k+1) = ZTH_i(i,k+1)
1037	! ZLWFT(i,k+1) = ZLWFT_i(i,klev+1-k)
1038	! ZSWFT(i,k+1) = ZSWFT_i(i,klev+1-k)
1039	ZFLUP(i,k+1) = ZFLUX_i(i,1,k+1)
1040	ZFLDN(i,k+1) = ZFLUX_i(i,2,k+1)
1041	ZFLUP0(i,k+1) = ZFLUC_i(i,1,k+1)
1042	ZFLDN0(i,k+1) = ZFLUC_i(i,2,k+1)
1043	ZFSDN(i,k+1) = ZFSDWN_i(i,k+1)*fract(i)
1044	ZFSDN0(i,k+1) = ZFCDWN_i(i,k+1)*fract(i)
1045	ZFSDNC0(i,k+1)= ZFCCDWN_i(i,k+1)*fract(i)
1046	ZFSUP (i,k+1) = ZFSUP_i(i,k+1)*fract(i)
1047	ZFSUP0(i,k+1) = ZFCUP_i(i,k+1)*fract(i)
1048	ZFSUPC0(i,k+1)= ZFCCUP_i(i,k+1)*fract(i)
1049	ZFLDNC0(i,k+1)= ZFLCCDWN_i(i,k+1)
1050	ZFLUPC0(i,k+1)= ZFLCCUP_i(i,k+1)
1051	IF (ok_volcan) THEN
1052	ZSWADAERO(i,k+1)=ZSWADAERO(i,k+1)*fract(i) !--NL
1053	ENDIF
1054
1055	! Nouveau calcul car visiblement ZSWFT et ZSWFC sont nuls dans RRTM cy32
1056	! en sortie de radlsw.F90 - MPL 7.01.09
1057	ZSWFT(i,k+1) = (ZFSDWN_i(i,k+1)-ZFSUP_i(i,k+1))*fract(i)
1058	ZSWFT0_i(i,k+1) = (ZFCDWN_i(i,k+1)-ZFCUP_i(i,k+1))*fract(i)
1059	! WRITE(*,'("FSDN FSUP FCDN FCUP: ",4E12.5)') ZFSDWN_i(i,k+1),&
1060	! ZFSUP_i(i,k+1),ZFCDWN_i(i,k+1),ZFCUP_i(i,k+1)
1061	ZLWFT(i,k+1) =-ZFLUX_i(i,2,k+1)-ZFLUX_i(i,1,k+1)
1062	ZLWFT0_i(i,k+1)=-ZFLUC_i(i,2,k+1)-ZFLUC_i(i,1,k+1)
1063	! print *,'FLUX2 FLUX1 FLUC2 FLUC1',ZFLUX_i(i,2,k+1),&
1064	! & ZFLUX_i(i,1,k+1),ZFLUC_i(i,2,k+1),ZFLUC_i(i,1,k+1)
1065	ENDDO
1066	ENDDO
1067
1068	!--ajout OB
1069	ZTOPSWADAERO(:) =ZTOPSWADAERO(:) *fract(:)
1070	ZSOLSWADAERO(:) =ZSOLSWADAERO(:) *fract(:)
1071	ZTOPSWAD0AERO(:)=ZTOPSWAD0AERO(:)*fract(:)
1072	ZSOLSWAD0AERO(:)=ZSOLSWAD0AERO(:)*fract(:)
1073	ZTOPSWAIAERO(:) =ZTOPSWAIAERO(:) *fract(:)
1074	ZSOLSWAIAERO(:) =ZSOLSWAIAERO(:) *fract(:)
1075	ZTOPSWCF_AERO(:,1)=ZTOPSWCF_AERO(:,1)*fract(:)
1076	ZTOPSWCF_AERO(:,2)=ZTOPSWCF_AERO(:,2)*fract(:)
1077	ZTOPSWCF_AERO(:,3)=ZTOPSWCF_AERO(:,3)*fract(:)
1078	ZSOLSWCF_AERO(:,1)=ZSOLSWCF_AERO(:,1)*fract(:)
1079	ZSOLSWCF_AERO(:,2)=ZSOLSWCF_AERO(:,2)*fract(:)
1080	ZSOLSWCF_AERO(:,3)=ZSOLSWCF_AERO(:,3)*fract(:)
1081
1082	! ---------
1083	! ---------
1084	! On renseigne les champs LMDz, pour avoir la meme chose qu'en sortie de
1085	! LW_LMDAR4 et SW_LMDAR4
1086
1087	!--fraction of diffuse radiation in surface SW downward radiation
1088	DO i = 1, kdlon
1089	IF (fract(i).GT.0.0) THEN
1090	zdir=SUM(PSFSWDIR(i,:))
1091	zdif=SUM(PSFSWDIF(i,:))
1092	zsolswfdiff(i) = zdif/(zdir+zdif)
1093	ELSE !--night
1094	zsolswfdiff(i) = 1.0
1095	ENDIF
1096	ENDDO
1097	!
1098	DO i = 1, kdlon
1099	zsolsw(i) = ZSWFT(i,1)
1100	zsolsw0(i) = ZSWFT0_i(i,1)
1101	! zsolsw0(i) = ZFSDN0(i,1) -ZFSUP0(i,1)
1102	ztopsw(i) = ZSWFT(i,klev+1)
1103	ztopsw0(i) = ZSWFT0_i(i,klev+1)
1104	! ztopsw0(i) = ZFSDN0(i,klev+1)-ZFSUP0(i,klev+1)
1105	!
1106	! zsollw(i) = ZFLDN(i,1) -ZFLUP(i,1)
1107	! zsollw0(i) = ZFLDN0(i,1) -ZFLUP0(i,1)
1108	! ztoplw(i) = ZFLDN(i,klev+1) -ZFLUP(i,klev+1)
1109	! ztoplw0(i) = ZFLDN0(i,klev+1)-ZFLUP0(i,klev+1)
1110	zsollw(i) = ZLWFT(i,1)
1111	zsollw0(i) = ZLWFT0_i(i,1)
1112	ztoplw(i) = ZLWFT(i,klev+1)*(-1)
1113	ztoplw0(i) = ZLWFT0_i(i,klev+1)*(-1)
1114	!
1115	IF (fract(i) == 0.) THEN
1116	!!!!! A REVOIR MPL (20090630) ca n a pas de sens quand fract=0
1117	! pas plus que dans le sw_AR4
1118	zalbpla(i) = 1.0e+39
1119	ELSE
1120	zalbpla(i) = ZFSUP(i,klev+1)/ZFSDN(i,klev+1)
1121	ENDIF
1122	!!! 5 juin 2015
1123	!!! Correction MP bug RRTM
1124	zsollwdown(i)= -1.*ZFLDN(i,1)
1125	ENDDO
1126	! print*,'OK2'
1127
1128	!--add VOLMIP (surf cool or strat heat activate)
1129	IF (flag_volc_surfstrat > 0) THEN
1130	DO i = 1, kdlon
1131	zsolsw(i) = volmip_solsw(i)*fract(i)
1132	ENDDO
1133	ENDIF
1134
1135	! extrait de SW_AR4
1136	! DO k = 1, KFLEV
1137	! kpl1 = k+1
1138	! DO i = 1, KDLON
1139	! PHEAT(i,k) = -(ZFSUP(i,kpl1)-ZFSUP(i,k)) -(ZFSDN(i,k)-ZFSDN(i,kpl1))
1140	! PHEAT(i,k) = PHEAT(i,k) * RDAY*RG/RCPD / PDP(i,k)
1141	! ZLWFT(klon,k),ZSWFT
1142
1143	DO k=1,kflev
1144	DO i=1,kdlon
1145	zheat(i,k)=(ZSWFT(i,k+1)-ZSWFT(i,k))RDAYRG/RCPD/PDP(i,k)
1146	zheat0(i,k)=(ZSWFT0_i(i,k+1)-ZSWFT0_i(i,k))RDAYRG/RCPD/PDP(i,k)
1147	zcool(i,k)=(ZLWFT(i,k)-ZLWFT(i,k+1))RDAYRG/RCPD/PDP(i,k)
1148	zcool0(i,k)=(ZLWFT0_i(i,k)-ZLWFT0_i(i,k+1))RDAYRG/RCPD/PDP(i,k)
1149	IF (ok_volcan) THEN
1150	zheat_volc(i,k)=(ZSWADAERO(i,k+1)-ZSWADAERO(i,k))*RG/RCPD/PDP(i,k) !NL
1151	zcool_volc(i,k)=(ZLWADAERO(i,k)-ZLWADAERO(i,k+1))*RG/RCPD/PDP(i,k) !NL
1152	ENDIF
1153	! print *,'heat cool heat0 cool0 ',zheat(i,k),zcool(i,k),zheat0(i,k),zcool0(i,k)
1154	! ZFLUCUP_i(i,k)=ZFLUC_i(i,1,k)
1155	! ZFLUCDWN_i(i,k)=ZFLUC_i(i,2,k)
1156	ENDDO
1157	ENDDO
1158	#else
1159	abort_message="You should compile with -rrtm if running with iflag_rrtm=1"
1160	call abort_physic(modname, abort_message, 1)
1161	#endif
1162	!======================================================================
1163	! AI fev 2021
1164	ELSE IF(iflag_rrtm == 2) THEN
1165	print*,'Traitement cas iflag_rrtm = ',iflag_rrtm
1166	! print*,'Mise a zero des flux '
1167	#ifdef CPP_ECRAD
1168	DO k = 1, kflev+1
1169	DO i = 1, kdlon
1170	ZEMTD_i(i,k)=0.
1171	ZEMTU_i(i,k)=0.
1172	ZTRSO_i(i,k)=0.
1173	ZTH_i(i,k)=0.
1174	ZLWFT_i(i,k)=0.
1175	ZSWFT_i(i,k)=0.
1176	ZFLUX_i(i,1,k)=0.
1177	ZFLUX_i(i,2,k)=0.
1178	ZFLUC_i(i,1,k)=0.
1179	ZFLUC_i(i,2,k)=0.
1180	ZFSDWN_i(i,k)=0.
1181	ZFCDWN_i(i,k)=0.
1182	ZFCCDWN_i(i,k)=0.
1183	ZFSUP_i(i,k)=0.
1184	ZFCUP_i(i,k)=0.
1185	ZFCCUP_i(i,k)=0.
1186	ZFLCCDWN_i(i,k)=0.
1187	ZFLCCUP_i(i,k)=0.
1188	ENDDO
1189	ENDDO
1190	!
1191	! AI ATTENTION Aerosols A REVOIR
1192	DO i = 1, kdlon
1193	DO k = 1, kflev
1194	DO kk= 1, naero_grp
1195	! DO kk=1, NSW
1196	!
1197	! PTAU_TOT(i,kflev+1-k,kk)=tau_aero_sw_rrtm(i,k,2,kk)
1198	! PPIZA_TOT(i,kflev+1-k,kk)=piz_aero_sw_rrtm(i,k,2,kk)
1199	! PCGA_TOT(i,kflev+1-k,kk)=cg_aero_sw_rrtm(i,k,2,kk)
1200	!
1201	! PTAU_NAT(i,kflev+1-k,kk)=tau_aero_sw_rrtm(i,k,1,kk)
1202	! PPIZA_NAT(i,kflev+1-k,kk)=piz_aero_sw_rrtm(i,k,1,kk)
1203	! PCGA_NAT(i,kflev+1-k,kk)=cg_aero_sw_rrtm(i,k,1,kk)
1204	! ZAEROSOL(i,kflev+1-k,kk)=m_allaer(i,k,kk)
1205	ZAEROSOL(i,kflev+1-k,kk)=m_allaer(i,k,kk)
1206	!
1207	ENDDO
1208	ENDDO
1209	ENDDO
1210	!-end OB
1211	!
1212	! DO i = 1, kdlon
1213	! DO k = 1, kflev
1214	! DO kk=1, NLW
1215	!
1216	! PTAU_LW_TOT(i,kflev+1-k,kk)=tau_aero_lw_rrtm(i,k,2,kk)
1217	! PTAU_LW_NAT(i,kflev+1-k,kk)=tau_aero_lw_rrtm(i,k,1,kk)
1218	!
1219	! ENDDO
1220	! ENDDO
1221	! ENDDO
1222	!-end C. Kleinschmitt
1223	!
1224	DO i = 1, kdlon
1225	ZCTRSO(i,1)=0.
1226	ZCTRSO(i,2)=0.
1227	ZCEMTR(i,1)=0.
1228	ZCEMTR(i,2)=0.
1229	ZTRSOD(i)=0.
1230	ZLWFC(i,1)=0.
1231	ZLWFC(i,2)=0.
1232	ZSWFC(i,1)=0.
1233	ZSWFC(i,2)=0.
1234	PFSDNN(i)=0.
1235	PFSDNV(i)=0.
1236	DO kk = 1, NSW
1237	PSFSWDIR(i,kk)=0.
1238	PSFSWDIF(i,kk)=0.
1239	ENDDO
1240	ENDDO
1241	!----- Fin des mises a zero des tableaux output -------------------
1242
1243	! On met les donnees dans l'ordre des niveaux ecrad
1244	! print*,'On inverse sur la verticale '
1245	paprs_i(:,1)=paprs(:,klev+1)
1246	DO k=1,klev
1247	paprs_i(1:klon,k+1) =paprs(1:klon,klev+1-k)
1248	pplay_i(1:klon,k) =pplay(1:klon,klev+1-k)
1249	cldfra_i(1:klon,k) =cldfra(1:klon,klev+1-k)
1250	PDP_i(1:klon,k) =PDP(1:klon,klev+1-k)
1251	t_i(1:klon,k) =t(1:klon,klev+1-k)
1252	q_i(1:klon,k) =q(1:klon,klev+1-k)
1253	qsat_i(1:klon,k) =qsat(1:klon,klev+1-k)
1254	flwc_i(1:klon,k) =flwc(1:klon,klev+1-k)
1255	fiwc_i(1:klon,k) =fiwc(1:klon,klev+1-k)
1256	ref_liq_i(1:klon,k) =ref_liq(1:klon,klev+1-k)*1.0e-6
1257	ref_ice_i(1:klon,k) =ref_ice(1:klon,klev+1-k)*1.0e-6
1258	!-OB
1259	ref_liq_pi_i(1:klon,k) =ref_liq_pi(1:klon,klev+1-k)
1260	ref_ice_pi_i(1:klon,k) =ref_ice_pi(1:klon,klev+1-k)
1261	ENDDO
1262	DO k=1,kflev
1263	POZON_i(1:klon,k,:)=POZON(1:klon,kflev+1-k,:)
1264	! ZO3_DP_i(1:klon,k)=ZO3_DP(1:klon,kflev+1-k)
1265	! DO i=1,6
1266	PAER_i(1:klon,k,:)=PAER(1:klon,kflev+1-k,:)
1267	! ENDDO
1268	ENDDO
1269
1270	! AI 11.2021
1271	! Calcul de ZTH_i (temp aux interfaces 1:klev+1)
1272	! IFS currently sets the half-level temperature at the surface to be
1273	! equal to the skin temperature. The radiation scheme takes as input
1274	! only the half-level temperatures and assumes the Planck function to
1275	! vary linearly in optical depth between half levels. In the lowest
1276	! atmospheric layer, where the atmospheric temperature can be much
1277	! cooler than the skin temperature, this can lead to significant
1278	! differences between the effective temperature of this lowest layer
1279	! and the true value in the model.
1280	! We may approximate the temperature profile in the lowest model level
1281	! as piecewise linear between the top of the layer T[k-1/2], the
1282	! centre of the layer T[k] and the base of the layer Tskin. The mean
1283	! temperature of the layer is then 0.25T[k-1/2] + 0.5T[k] +
1284	! 0.25*Tskin, which can be achieved by setting the atmospheric
1285	! temperature at the half-level corresponding to the surface as
1286	! follows:
1287	! AI ATTENTION fais dans interface radlw
1288	!thermodynamics%temperature_hl(KIDIA:KFDIA,KLEV+1) &
1289	! & = PTEMPERATURE(KIDIA:KFDIA,KLEV) &
1290	! & + 0.5_JPRB * (PTEMPERATURE_H(KIDIA:KFDIA,KLEV+1) &
1291	! & -PTEMPERATURE_H(KIDIA:KFDIA,KLEV))
1292
1293	DO K=2,KLEV
1294	DO i = 1, kdlon
1295	ZTH_i(i,K)=&
1296	& (t_i(i,K-1)pplay_i(i,K-1)(pplay_i(i,K)-paprs_i(i,K))&
1297	& +t_i(i,K)pplay_i(i,K)(paprs_i(i,K)-pplay_i(i,K-1)))&
1298	& (1.0/(paprs_i(i,K)(pplay_i(i,K)-pplay_i(i,K-1))))
1299	ENDDO
1300	ENDDO
1301	DO i = 1, kdlon
1302	! Sommet
1303	ZTH_i(i,1)=t_i(i,1)-pplay_i(i,1)*(t_i(i,1)-ZTH_i(i,2))&
1304	& /(pplay_i(i,1)-paprs_i(i,2))
1305	! Vers le sol
1306	ZTH_i(i,KLEV+1)=t_i(i,KLEV) + 0.5 * &
1307	(tsol(i) - ZTH_i(i,KLEV))
1308	ENDDO
1309
1310
1311	print *,'RADLWSW: avant RADIATION_SCHEME '
1312
1313	! AI mars 2022
1314	SOLARIRAD = solaire/zdist/zdist
1315	!! diagnos pour la comparaison a la version offline
1316	!!! - Gas en VMR pour offline et MMR pour online
1317	!!! - on utilise pour solarirrad une valeur constante
1318	if (lldebug_for_offline) then
1319	SOLARIRAD = 1366.0896
1320	ZCH4_off = CH4_ppb*1e-9
1321	ZN2O_off = N2O_ppb*1e-9
1322	ZNO2_off = 0.0
1323	ZCFC11_off = CFC11_ppt*1e-12
1324	ZCFC12_off = CFC12_ppt*1e-12
1325	ZHCFC22_off = 0.0
1326	ZCCL4_off = 0.0
1327	ZO2_off = 0.0
1328	ZCO2_off = co2_ppm*1e-6
1329
1330	CALL writefield_phy('rmu0',rmu0,1)
1331	CALL writefield_phy('tsol',tsol,1)
1332	CALL writefield_phy('emissiv_out',ZEMIS,1)
1333	CALL writefield_phy('paprs_i',paprs_i,klev+1)
1334	CALL writefield_phy('ZTH_i',ZTH_i,klev+1)
1335	CALL writefield_phy('cldfra_i',cldfra_i,klev)
1336	CALL writefield_phy('q_i',q_i,klev)
1337	CALL writefield_phy('fiwc_i',fiwc_i,klev)
1338	CALL writefield_phy('flwc_i',flwc_i,klev)
1339	CALL writefield_phy('palbd_new',PALBD_NEW,NSW)
1340	CALL writefield_phy('palbp_new',PALBP_NEW,NSW)
1341	CALL writefield_phy('POZON',POZON_i(:,:,1),klev)
1342	CALL writefield_phy('ZCO2',ZCO2_off,klev)
1343	CALL writefield_phy('ZCH4',ZCH4_off,klev)
1344	CALL writefield_phy('ZN2O',ZN2O_off,klev)
1345	CALL writefield_phy('ZO2',ZO2_off,klev)
1346	CALL writefield_phy('ZNO2',ZNO2_off,klev)
1347	CALL writefield_phy('ZCFC11',ZCFC11_off,klev)
1348	CALL writefield_phy('ZCFC12',ZCFC12_off,klev)
1349	CALL writefield_phy('ZHCFC22',ZHCFC22_off,klev)
1350	CALL writefield_phy('ZCCL4',ZCCL4_off,klev)
1351	CALL writefield_phy('ref_liq_i',ref_liq_i,klev)
1352	CALL writefield_phy('ref_ice_i',ref_ice_i,klev)
1353	endif
1354	! lldebug_for_offline
1355
1356	CALL RADIATION_SCHEME &
1357	& (ist, iend, klon, klev, naero_grp, NSW, &
1358	& namelist_ecrad_file, ok_3Deffect, &
1359	& day_cur, current_time, &
1360	! Cste solaire/(d_Terre-Soleil)**2
1361	& SOLARIRAD, &
1362	! Cos(angle zin), temp sol
1363	& rmu0, tsol, &
1364	! Albedo diffuse et directe
1365	& PALBD_NEW,PALBP_NEW, &
1366	! Emessivite : PEMIS_WINDOW (???), &
1367	& ZEMIS, ZEMISW, &
1368	! longitude(rad), sin(latitude), PMASQ_ ???
1369	& ZGELAM, ZGEMU, &
1370	! Temp et pres aux interf, vapeur eau, Satur spec humid
1371	& paprs_i, ZTH_i, q_i, qsat_i, &
1372	! Gas
1373	& ZCO2, ZCH4, ZN2O, ZNO2, ZCFC11, ZCFC12, ZHCFC22, &
1374	& ZCCL4, POZON_i(:,:,1), ZO2, &
1375	! nuages :
1376	& cldfra_i, flwc_i, fiwc_i, ZQ_SNOW, &
1377	! rayons effectifs des gouttelettes
1378	& ref_liq_i, ref_ice_i, &
1379	! aerosols
1380	& ZAEROSOL_OLD, ZAEROSOL, &
1381	! Outputs
1382	! Net flux :
1383	& ZSWFT_i, ZLWFT_i, ZSWFT0_ii, ZLWFT0_ii, &
1384	! DWN flux :
1385	& ZFSDWN_i, ZFLUX_i(:,2,:), ZFCDWN_i, ZFLUC_i(:,2,:), &
1386	! UP flux :
1387	& ZFSUP_i, ZFLUX_i(:,1,:), ZFCUP_i, ZFLUC_i(:,1,:), &
1388	! Surf Direct flux : ATTENTION
1389	& ZFLUX_DIR, ZFLUX_DIR_CLEAR, ZFLUX_DIR_INTO_SUN, &
1390	! UV and para flux
1391	& ZFLUX_UV, ZFLUX_PAR, ZFLUX_PAR_CLEAR, &
1392	! & ZFLUX_SW_DN_TOA,
1393	& ZEMIS_OUT, ZLWDERIVATIVE, &
1394	& PSFSWDIF, PSFSWDIR, &
1395	& cloud_cover_sw)
1396
1397	print *,'========= RADLWSW: apres RADIATION_SCHEME ==================== '
1398
1399	if (lldebug_for_offline) then
1400	CALL writefield_phy('FLUX_LW',ZLWFT_i,klev+1)
1401	CALL writefield_phy('FLUX_LW_CLEAR',ZLWFT0_ii,klev+1)
1402	CALL writefield_phy('FLUX_SW',ZSWFT_i,klev+1)
1403	CALL writefield_phy('FLUX_SW_CLEAR',ZSWFT0_ii,klev+1)
1404	CALL writefield_phy('FLUX_DN_SW',ZFSDWN_i,klev+1)
1405	CALL writefield_phy('FLUX_DN_LW',ZFLUX_i(:,2,:),klev+1)
1406	CALL writefield_phy('FLUX_DN_SW_CLEAR',ZFCDWN_i,klev+1)
1407	CALL writefield_phy('FLUX_DN_LW_CLEAR',ZFLUC_i(:,2,:),klev+1)
1408	CALL writefield_phy('PSFSWDIR',PSFSWDIR,6)
1409	CALL writefield_phy('PSFSWDIF',PSFSWDIF,6)
1410	CALL writefield_phy('FLUX_UP_LW',ZFLUX_i(:,1,:),klev+1)
1411	CALL writefield_phy('FLUX_UP_LW_CLEAR',ZFLUC_i(:,1,:),klev+1)
1412	CALL writefield_phy('FLUX_UP_SW',ZFSUP_i,klev+1)
1413	CALL writefield_phy('FLUX_UP_SW_CLEAR',ZFCUP_i,klev+1)
1414	endif
1415
1416	! ---------
1417	! On retablit l'ordre des niveaux lmd pour les tableaux de sortie
1418	! D autre part, on multiplie les resultats SW par fract pour etre coherent
1419	! avec l ancien rayonnement AR4. Si nuit, fract=0 donc pas de
1420	! rayonnement SW. (MPL 260609)
1421	print*,'On retablit l ordre des niveaux verticaux pour LMDZ'
1422	print*,'On multiplie les flux SW par fract et LW dwn par -1'
1423	DO k=0,klev
1424	DO i=1,klon
1425	ZEMTD(i,k+1) = ZEMTD_i(i,klev+1-k)
1426	ZEMTU(i,k+1) = ZEMTU_i(i,klev+1-k)
1427	ZTRSO(i,k+1) = ZTRSO_i(i,klev+1-k)
1428	! ZTH(i,k+1) = ZTH_i(i,klev+1-k)
1429	! AI ATTENTION
1430	ZLWFT(i,k+1) = ZLWFT_i(i,klev+1-k)
1431	ZSWFT(i,k+1) = ZSWFT_i(i,klev+1-k)*fract(i)
1432	ZSWFT0_i(i,k+1) = ZSWFT0_ii(i,klev+1-k)*fract(i)
1433	ZLWFT0_i(i,k+1) = ZLWFT0_ii(i,klev+1-k)
1434	!
1435	ZFLUP(i,k+1) = ZFLUX_i(i,1,klev+1-k)
1436	ZFLDN(i,k+1) = -1.*ZFLUX_i(i,2,klev+1-k)
1437	ZFLUP0(i,k+1) = ZFLUC_i(i,1,klev+1-k)
1438	ZFLDN0(i,k+1) = -1.*ZFLUC_i(i,2,klev+1-k)
1439	ZFSDN(i,k+1) = ZFSDWN_i(i,klev+1-k)*fract(i)
1440	ZFSDN0(i,k+1) = ZFCDWN_i(i,klev+1-k)*fract(i)
1441	ZFSDNC0(i,k+1)= ZFCCDWN_i(i,klev+1-k)*fract(i)
1442	ZFSUP (i,k+1) = ZFSUP_i(i,klev+1-k)*fract(i)
1443	ZFSUP0(i,k+1) = ZFCUP_i(i,klev+1-k)*fract(i)
1444	ZFSUPC0(i,k+1)= ZFCCUP_i(i,klev+1-k)*fract(i)
1445	ZFLDNC0(i,k+1)= -1.*ZFLCCDWN_i(i,klev+1-k)
1446	ZFLUPC0(i,k+1)= ZFLCCUP_i(i,klev+1-k)
1447	IF (ok_volcan) THEN
1448	ZSWADAERO(i,k+1)=ZSWADAERO(i,klev+1-k)*fract(i) !--NL
1449	ENDIF
1450
1451	! Nouveau calcul car visiblement ZSWFT et ZSWFC sont nuls dans RRTM cy32
1452	! en sortie de radlsw.F90 - MPL 7.01.09
1453	! AI ATTENTION
1454	! ZSWFT(i,k+1) = (ZFSDWN_i(i,k+1)-ZFSUP_i(i,k+1))*fract(i)
1455	! ZSWFT0_i(i,k+1) = (ZFCDWN_i(i,k+1)-ZFCUP_i(i,k+1))*fract(i)
1456	! ZLWFT(i,k+1) =-ZFLUX_i(i,2,k+1)-ZFLUX_i(i,1,k+1)
1457	! ZLWFT0_i(i,k+1)=-ZFLUC_i(i,2,k+1)-ZFLUC_i(i,1,k+1)
1458	ENDDO
1459	ENDDO
1460
1461	!--ajout OB
1462	ZTOPSWADAERO(:) =ZTOPSWADAERO(:) *fract(:)
1463	ZSOLSWADAERO(:) =ZSOLSWADAERO(:) *fract(:)
1464	ZTOPSWAD0AERO(:)=ZTOPSWAD0AERO(:)*fract(:)
1465	ZSOLSWAD0AERO(:)=ZSOLSWAD0AERO(:)*fract(:)
1466	ZTOPSWAIAERO(:) =ZTOPSWAIAERO(:) *fract(:)
1467	ZSOLSWAIAERO(:) =ZSOLSWAIAERO(:) *fract(:)
1468	ZTOPSWCF_AERO(:,1)=ZTOPSWCF_AERO(:,1)*fract(:)
1469	ZTOPSWCF_AERO(:,2)=ZTOPSWCF_AERO(:,2)*fract(:)
1470	ZTOPSWCF_AERO(:,3)=ZTOPSWCF_AERO(:,3)*fract(:)
1471	ZSOLSWCF_AERO(:,1)=ZSOLSWCF_AERO(:,1)*fract(:)
1472	ZSOLSWCF_AERO(:,2)=ZSOLSWCF_AERO(:,2)*fract(:)
1473	ZSOLSWCF_AERO(:,3)=ZSOLSWCF_AERO(:,3)*fract(:)
1474
1475	! ---------
1476	! On renseigne les champs LMDz, pour avoir la meme chose qu'en sortie de
1477	! LW_LMDAR4 et SW_LMDAR4
1478
1479	!--fraction of diffuse radiation in surface SW downward radiation
1480	DO i = 1, kdlon
1481	zdir=SUM(PSFSWDIR(i,:))
1482	zdif=SUM(PSFSWDIF(i,:))
1483	IF (fract(i).GT.0.0.and.(zdir+zdif).gt.seuilmach) THEN
1484	zsolswfdiff(i) = zdif/(zdir+zdif)
1485	ELSE !--night
1486	zsolswfdiff(i) = 1.0
1487	ENDIF
1488	ENDDO
1489	!
1490	DO i = 1, kdlon
1491	zsolsw(i) = ZSWFT(i,1)
1492	zsolsw0(i) = ZSWFT0_i(i,1)
1493	ztopsw(i) = ZSWFT(i,klev+1)
1494	ztopsw0(i) = ZSWFT0_i(i,klev+1)
1495	zsollw(i) = ZLWFT(i,1)
1496	zsollw0(i) = ZLWFT0_i(i,1)
1497	ztoplw(i) = ZLWFT(i,klev+1)*(-1)
1498	ztoplw0(i) = ZLWFT0_i(i,klev+1)*(-1)
1499	!
1500	zsollwdown(i)= -1.*ZFLDN(i,1)
1501	ENDDO
1502
1503	DO k=1,kflev
1504	DO i=1,kdlon
1505	zheat(i,k)=(ZSWFT(i,k+1)-ZSWFT(i,k))RDAYRG/RCPD/PDP(i,k)
1506	zheat0(i,k)=(ZSWFT0_i(i,k+1)-ZSWFT0_i(i,k))RDAYRG/RCPD/PDP(i,k)
1507	zcool(i,k)=(ZLWFT(i,k)-ZLWFT(i,k+1))RDAYRG/RCPD/PDP(i,k)
1508	zcool0(i,k)=(ZLWFT0_i(i,k)-ZLWFT0_i(i,k+1))RDAYRG/RCPD/PDP(i,k)
1509	IF (ok_volcan) THEN
1510	zheat_volc(i,k)=(ZSWADAERO(i,k+1)-ZSWADAERO(i,k))*RG/RCPD/PDP(i,k) !NL
1511	zcool_volc(i,k)=(ZLWADAERO(i,k)-ZLWADAERO(i,k+1))*RG/RCPD/PDP(i,k) !NL
1512	ENDIF
1513	ENDDO
1514	ENDDO
1515	#endif
1516	print*,'Fin traitement ECRAD'
1517	! Fin ECRAD
1518	ENDIF ! iflag_rrtm
1519	! ecrad
1520	!======================================================================
1521
1522	DO i = 1, kdlon
1523	topsw(iof+i) = ztopsw(i)
1524	toplw(iof+i) = ztoplw(i)
1525	solsw(iof+i) = zsolsw(i)
1526	solswfdiff(iof+i) = zsolswfdiff(i)
1527	sollw(iof+i) = zsollw(i)
1528	sollwdown(iof+i) = zsollwdown(i)
1529	DO k = 1, kflev+1
1530	lwdn0 ( iof+i,k) = ZFLDN0 ( i,k)
1531	lwdn ( iof+i,k) = ZFLDN ( i,k)
1532	lwup0 ( iof+i,k) = ZFLUP0 ( i,k)
1533	lwup ( iof+i,k) = ZFLUP ( i,k)
1534	ENDDO
1535	topsw0(iof+i) = ztopsw0(i)
1536	toplw0(iof+i) = ztoplw0(i)
1537	solsw0(iof+i) = zsolsw0(i)
1538	sollw0(iof+i) = zsollw0(i)
1539	albpla(iof+i) = zalbpla(i)
1540
1541	DO k = 1, kflev+1
1542	swdnc0( iof+i,k) = ZFSDNC0( i,k)
1543	swdn0 ( iof+i,k) = ZFSDN0 ( i,k)
1544	swdn ( iof+i,k) = ZFSDN ( i,k)
1545	swupc0( iof+i,k) = ZFSUPC0( i,k)
1546	swup0 ( iof+i,k) = ZFSUP0 ( i,k)
1547	swup ( iof+i,k) = ZFSUP ( i,k)
1548	lwdnc0( iof+i,k) = ZFLDNC0( i,k)
1549	lwupc0( iof+i,k) = ZFLUPC0( i,k)
1550	ENDDO
1551	ENDDO
1552	!-transform the aerosol forcings, if they have
1553	! to be calculated
1554	IF (ok_ade) THEN
1555	DO i = 1, kdlon
1556	topswad_aero(iof+i) = ztopswadaero(i)
1557	topswad0_aero(iof+i) = ztopswad0aero(i)
1558	solswad_aero(iof+i) = zsolswadaero(i)
1559	solswad0_aero(iof+i) = zsolswad0aero(i)
1560	topsw_aero(iof+i,:) = ztopsw_aero(i,:)
1561	topsw0_aero(iof+i,:) = ztopsw0_aero(i,:)
1562	solsw_aero(iof+i,:) = zsolsw_aero(i,:)
1563	solsw0_aero(iof+i,:) = zsolsw0_aero(i,:)
1564	topswcf_aero(iof+i,:) = ztopswcf_aero(i,:)
1565	solswcf_aero(iof+i,:) = zsolswcf_aero(i,:)
1566	!-LW
1567	toplwad_aero(iof+i) = ztoplwadaero(i)
1568	toplwad0_aero(iof+i) = ztoplwad0aero(i)
1569	sollwad_aero(iof+i) = zsollwadaero(i)
1570	sollwad0_aero(iof+i) = zsollwad0aero(i)
1571	ENDDO
1572	ELSE
1573	DO i = 1, kdlon
1574	topswad_aero(iof+i) = 0.0
1575	solswad_aero(iof+i) = 0.0
1576	topswad0_aero(iof+i) = 0.0
1577	solswad0_aero(iof+i) = 0.0
1578	topsw_aero(iof+i,:) = 0.
1579	topsw0_aero(iof+i,:) =0.
1580	solsw_aero(iof+i,:) = 0.
1581	solsw0_aero(iof+i,:) = 0.
1582	!-LW
1583	toplwad_aero(iof+i) = 0.0
1584	sollwad_aero(iof+i) = 0.0
1585	toplwad0_aero(iof+i) = 0.0
1586	sollwad0_aero(iof+i) = 0.0
1587	ENDDO
1588	ENDIF
1589	IF (ok_aie) THEN
1590	DO i = 1, kdlon
1591	topswai_aero(iof+i) = ztopswaiaero(i)
1592	solswai_aero(iof+i) = zsolswaiaero(i)
1593	!-LW
1594	toplwai_aero(iof+i) = ztoplwaiaero(i)
1595	sollwai_aero(iof+i) = zsollwaiaero(i)
1596	ENDDO
1597	ELSE
1598	DO i = 1, kdlon
1599	topswai_aero(iof+i) = 0.0
1600	solswai_aero(iof+i) = 0.0
1601	!-LW
1602	toplwai_aero(iof+i) = 0.0
1603	sollwai_aero(iof+i) = 0.0
1604	ENDDO
1605	ENDIF
1606	DO k = 1, kflev
1607	DO i = 1, kdlon
1608	! scale factor to take into account the difference between
1609	! dry air and watter vapour scpecifi! heat capacity
1610	zznormcp=1.0+RVTMP2*PWV(i,k)
1611	heat(iof+i,k) = zheat(i,k)/zznormcp
1612	cool(iof+i,k) = zcool(i,k)/zznormcp
1613	heat0(iof+i,k) = zheat0(i,k)/zznormcp
1614	cool0(iof+i,k) = zcool0(i,k)/zznormcp
1615	IF(ok_volcan) THEN !NL
1616	heat_volc(iof+i,k) = zheat_volc(i,k)/zznormcp
1617	cool_volc(iof+i,k) = zcool_volc(i,k)/zznormcp
1618	ENDIF
1619	ENDDO
1620	ENDDO
1621
1622	ENDDO ! j = 1, nb_gr
1623
1624	IF (lldebug) THEN
1625	if (0.eq.1) then
1626	! Verifs dans le cas 1D
1627	print*,'================== Sortie de radlw ================='
1628	print,'***** LW LW LW *****************'
1629	print*,'ZLWFT =',ZLWFT
1630	print*,'ZLWFT0_i =',ZLWFT0_i
1631	print*,'ZFLUP0 =',ZFLUP0
1632	print*,'ZFLDN0 =',ZFLDN0
1633	print*,'ZFLDNC0 =',ZFLDNC0
1634	print*,'ZFLUPC0 =',ZFLUPC0
1635
1636	print,'***** SW SW SW *****************'
1637	print*,'ZSWFT =',ZSWFT
1638	print*,'ZSWFT0_i =',ZSWFT0_i
1639	print*,'ZFSDN =',ZFSDN
1640	print*,'ZFSDN0 =',ZFSDN0
1641	print*,'ZFSDNC0 =',ZFSDNC0
1642	print*,'ZFSUP =',ZFSUP
1643	print*,'ZFSUP0 =',ZFSUP0
1644	print*,'ZFSUPC0 =',ZFSUPC0
1645
1646	print,'***** LMDZ *****************'
1647	print*,'cool = ', cool
1648	print*,'heat = ', heat
1649	print*,'topsw = ', topsw
1650	print*,'toplw = ', toplw
1651	print*,'sollw = ', sollw
1652	print*,'solsw = ', solsw
1653	print*,'lwdn = ', lwdn
1654	print*,'lwup = ', lwup
1655	print*,'swdn = ', swdn
1656	print*,'swup =', swup
1657	endif
1658	ENDIF
1659
1660	END SUBROUTINE radlwsw
1661
1662	end module radlwsw_m

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

Download in other formats:

Original Format