Context Navigation

radlwsw_m.F90 @ 5306

Last change on this file since 5306 was 3756, checked in by oboucher, 4 years ago
The fraction of SW down radiation that is diffuse is extracted from radlwsw.F90 and passed on to pbl_surface for inclusion in fields_out for transfer to ORCHIDEE if it requested in input file coupling_fields.def. Hence there is nothing automatic here. The fraction of SW down radiation is also added as a diagnostic. To satisfy the case when this new quantity is used in ORCHIDEE, then it is added in the restart (and read in case it is there).
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: 47.1 KB

Line
1	!
2	! $Id: radlwsw_m.F90 3756 2020-07-10 21:50:17Z fairhead $
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_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 pour 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	heat,heat0,cool,cool0,albpla,&
27	heat_volc, cool_volc,&
28	topsw,toplw,solsw,solswfdiff,sollw,&
29	sollwdown,&
30	topsw0,toplw0,solsw0,sollw0,&
31	lwdnc0, lwdn0, lwdn, lwupc0, lwup0, lwup,&
32	swdnc0, swdn0, swdn, swupc0, swup0, swup,&
33	topswad_aero, solswad_aero,&
34	topswai_aero, solswai_aero, &
35	topswad0_aero, solswad0_aero,&
36	topsw_aero, topsw0_aero,&
37	solsw_aero, solsw0_aero, &
38	topswcf_aero, solswcf_aero,&
39	!-C. Kleinschmitt for LW diagnostics
40	toplwad_aero, sollwad_aero,&
41	toplwai_aero, sollwai_aero, &
42	toplwad0_aero, sollwad0_aero,&
43	!-end
44	ZLWFT0_i, ZFLDN0, ZFLUP0,&
45	ZSWFT0_i, ZFSDN0, ZFSUP0)
46
47
48
49	USE DIMPHY
50	USE assert_m, ONLY : assert
51	USE infotrac_phy, ONLY : type_trac
52	USE write_field_phy
53	#ifdef REPROBUS
54	USE CHEM_REP, ONLY : solaireTIME, ok_SUNTIME, ndimozon
55	#endif
56	#ifdef CPP_RRTM
57	! modules necessaires au rayonnement
58	! -----------------------------------------
59	! USE YOMCST , ONLY : RG ,RD ,RTT ,RPI
60	! USE YOERAD , ONLY : NSW ,LRRTM ,LINHOM , LCCNL,LCCNO,
61	! USE YOERAD , ONLY : NSW ,LRRTM ,LCCNL ,LCCNO ,&
62	! NSW mis dans .def MPL 20140211
63	! NLW ajoute par OB
64	USE YOERAD , ONLY : NLW, LRRTM ,LCCNL ,LCCNO ,&
65	NRADIP , NRADLP , NICEOPT, NLIQOPT ,RCCNLND , RCCNSEA
66	USE YOELW , ONLY : NSIL ,NTRA ,NUA ,TSTAND ,XP
67	USE YOESW , ONLY : RYFWCA ,RYFWCB ,RYFWCC ,RYFWCD,&
68	RYFWCE ,RYFWCF ,REBCUA ,REBCUB ,REBCUC,&
69	REBCUD ,REBCUE ,REBCUF ,REBCUI ,REBCUJ,&
70	REBCUG ,REBCUH ,RHSAVI ,RFULIO ,RFLAA0,&
71	RFLAA1 ,RFLBB0 ,RFLBB1 ,RFLBB2 ,RFLBB3,&
72	RFLCC0 ,RFLCC1 ,RFLCC2 ,RFLCC3 ,RFLDD0,&
73	RFLDD1 ,RFLDD2 ,RFLDD3 ,RFUETA ,RASWCA,&
74	RASWCB ,RASWCC ,RASWCD ,RASWCE ,RASWCF
75	! & RASWCB ,RASWCC ,RASWCD ,RASWCE ,RASWCF, RLINLI
76	USE YOERDU , ONLY : NUAER ,NTRAER ,REPLOG ,REPSC ,REPSCW ,DIFF
77	! USE YOETHF , ONLY : RTICE
78	USE YOERRTWN , ONLY : DELWAVE ,TOTPLNK
79	USE YOMPHY3 , ONLY : RII0
80	#endif
81	USE aero_mod
82
83	!======================================================================
84	! Auteur(s): Z.X. Li (LMD/CNRS) date: 19960719
85	! Objet: interface entre le modele et les rayonnements
86	! Arguments:
87	! dist-----input-R- distance astronomique terre-soleil
88	! rmu0-----input-R- cosinus de l'angle zenithal
89	! fract----input-R- duree d'ensoleillement normalisee
90	! co2_ppm--input-R- concentration du gaz carbonique (en ppm)
91	! paprs----input-R- pression a inter-couche (Pa)
92	! pplay----input-R- pression au milieu de couche (Pa)
93	! tsol-----input-R- temperature du sol (en K)
94	! alb1-----input-R- albedo du sol(entre 0 et 1) dans l'interval visible
95	! alb2-----input-R- albedo du sol(entre 0 et 1) dans l'interval proche infra-rouge
96	! t--------input-R- temperature (K)
97	! q--------input-R- vapeur d'eau (en kg/kg)
98	! cldfra---input-R- fraction nuageuse (entre 0 et 1)
99	! cldtaupd---input-R- epaisseur optique des nuages dans le visible (present-day value)
100	! cldemi---input-R- emissivite des nuages dans l'IR (entre 0 et 1)
101	! ok_ade---input-L- apply the Aerosol Direct Effect or not?
102	! ok_aie---input-L- apply the Aerosol Indirect Effect or not?
103	! ok_volcan-input-L- activate volcanic diags (SW heat & LW cool rate, SW & LW flux)
104	! flag_aerosol-input-I- aerosol flag from 0 to 6
105	! flag_aerosol_strat-input-I- use stratospheric aerosols flag (0, 1, 2)
106	! flag_aer_feedback-input-I- activate aerosol radiative feedback (T, F)
107	! tau_ae, piz_ae, cg_ae-input-R- aerosol optical properties (calculated in aeropt.F)
108	! cldtaupi-input-R- epaisseur optique des nuages dans le visible
109	! calculated for pre-industrial (pi) aerosol concentrations, i.e. with smaller
110	! droplet concentration, thus larger droplets, thus generally cdltaupi cldtaupd
111	! it is needed for the diagnostics of the aerosol indirect radiative forcing
112	!
113	! heat-----output-R- echauffement atmospherique (visible) (K/jour)
114	! cool-----output-R- refroidissement dans l'IR (K/jour)
115	! albpla---output-R- albedo planetaire (entre 0 et 1)
116	! topsw----output-R- flux solaire net au sommet de l'atm.
117	! toplw----output-R- ray. IR montant au sommet de l'atmosphere
118	! solsw----output-R- flux solaire net a la surface
119	! solswfdiff----output-R- fraction de rayonnement diffus pour le flux solaire descendant a la surface
120	! sollw----output-R- ray. IR montant a la surface
121	! solswad---output-R- ray. solaire net absorbe a la surface (aerosol dir)
122	! topswad---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol dir)
123	! solswai---output-R- ray. solaire net absorbe a la surface (aerosol ind)
124	! topswai---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol ind)
125	!
126	! heat_volc-----output-R- echauffement atmospherique du au forcage volcanique (visible) (K/s)
127	! cool_volc-----output-R- refroidissement dans l'IR du au forcage volcanique (K/s)
128	!
129	! ATTENTION: swai and swad have to be interpreted in the following manner:
130	! ---------
131	! ok_ade=F & ok_aie=F -both are zero
132	! ok_ade=T & ok_aie=F -aerosol direct forcing is F_{AD} = topsw-topswad
133	! indirect is zero
134	! ok_ade=F & ok_aie=T -aerosol indirect forcing is F_{AI} = topsw-topswai
135	! direct is zero
136	! ok_ade=T & ok_aie=T -aerosol indirect forcing is F_{AI} = topsw-topswai
137	! aerosol direct forcing is F_{AD} = topswai-topswad
138	!
139	! --------- RRTM: output RECMWFL
140	! ZEMTD (KPROMA,KLEV+1) ; TOTAL DOWNWARD LONGWAVE EMISSIVITY
141	! ZEMTU (KPROMA,KLEV+1) ; TOTAL UPWARD LONGWAVE EMISSIVITY
142	! ZTRSO (KPROMA,KLEV+1) ; TOTAL SHORTWAVE TRANSMISSIVITY
143	! ZTH (KPROMA,KLEV+1) ; HALF LEVEL TEMPERATURE
144	! ZCTRSO(KPROMA,2) ; CLEAR-SKY SHORTWAVE TRANSMISSIVITY
145	! ZCEMTR(KPROMA,2) ; CLEAR-SKY NET LONGWAVE EMISSIVITY
146	! ZTRSOD(KPROMA) ; TOTAL-SKY SURFACE SW TRANSMISSITY
147	! ZLWFC (KPROMA,2) ; CLEAR-SKY LONGWAVE FLUXES
148	! ZLWFT (KPROMA,KLEV+1) ; TOTAL-SKY LONGWAVE FLUXES
149	! ZLWFT0(KPROMA,KLEV+1) ; CLEAR-SKY LONGWAVE FLUXES ! added by MPL 090109
150	! ZSWFC (KPROMA,2) ; CLEAR-SKY SHORTWAVE FLUXES
151	! ZSWFT (KPROMA,KLEV+1) ; TOTAL-SKY SHORTWAVE FLUXES
152	! ZSWFT0(KPROMA,KLEV+1) ; CLEAR-SKY SHORTWAVE FLUXES ! added by MPL 090109
153	! ZFLUX (KLON,2,KLEV+1) ; TOTAL LW FLUXES 1=up, 2=DWN ! added by MPL 080411
154	! ZFLUC (KLON,2,KLEV+1) ; CLEAR SKY LW FLUXES ! added by MPL 080411
155	! ZFSDWN(klon,KLEV+1) ; TOTAL SW DWN FLUXES ! added by MPL 080411
156	! ZFCDWN(klon,KLEV+1) ; CLEAR SKY SW DWN FLUXES ! added by MPL 080411
157	! ZFCCDWN(klon,KLEV+1) ; CLEAR SKY CLEAN (NO AEROSOL) SW DWN FLUXES ! added by OB 211117
158	! ZFSUP (klon,KLEV+1) ; TOTAL SW UP FLUXES ! added by MPL 080411
159	! ZFCUP (klon,KLEV+1) ; CLEAR SKY SW UP FLUXES ! added by MPL 080411
160	! ZFCCUP (klon,KLEV+1) ; CLEAR SKY CLEAN (NO AEROSOL) SW UP FLUXES ! added by OB 211117
161	! ZFLCCDWN(klon,KLEV+1) ; CLEAR SKY CLEAN (NO AEROSOL) LW DWN FLUXES ! added by OB 211117
162	! ZFLCCUP (klon,KLEV+1) ; CLEAR SKY CLEAN (NO AEROSOL) LW UP FLUXES ! added by OB 211117
163
164	!======================================================================
165
166	! ====================================================================
167	! Adapte au modele de chimie INCA par Celine Deandreis & Anne Cozic -- 2009
168	! 1 = ZERO
169	! 2 = AER total
170	! 3 = NAT
171	! 4 = BC
172	! 5 = SO4
173	! 6 = POM
174	! 7 = DUST
175	! 8 = SS
176	! 9 = NO3
177	!
178	! ====================================================================
179	include "YOETHF.h"
180	include "YOMCST.h"
181	include "clesphys.h"
182
183	! Input arguments
184	REAL, INTENT(in) :: dist
185	REAL, INTENT(in) :: rmu0(KLON), fract(KLON)
186	REAL, INTENT(in) :: paprs(KLON,KLEV+1), pplay(KLON,KLEV)
187	!albedo SB >>>
188	! REAL, INTENT(in) :: alb1(KLON), alb2(KLON), tsol(KLON)
189	REAL, INTENT(in) :: tsol(KLON)
190	REAL, INTENT(in) :: alb_dir(KLON,NSW),alb_dif(KLON,NSW)
191	REAL, INTENT(in) :: SFRWL(6)
192	!albedo SB <<<
193	REAL, INTENT(in) :: t(KLON,KLEV), q(KLON,KLEV)
194
195	REAL, INTENT(in):: wo(:, :, :) ! dimension(KLON,KLEV, 1 or 2)
196	! column-density of ozone in a layer, in kilo-Dobsons
197	! "wo(:, :, 1)" is for the average day-night field,
198	! "wo(:, :, 2)" is for daylight time.
199
200	LOGICAL, INTENT(in) :: ok_ade, ok_aie ! switches whether to use aerosol direct (indirect) effects or not
201	LOGICAL, INTENT(in) :: ok_volcan ! produce volcanic diags (SW/LW heat flux and rate)
202	LOGICAL :: lldebug
203	INTEGER, INTENT(in) :: flag_aerosol ! takes value 0 (no aerosol) or 1 to 6 (aerosols)
204	INTEGER, INTENT(in) :: flag_aerosol_strat ! use stratospheric aerosols
205	LOGICAL, INTENT(in) :: flag_aer_feedback ! activate aerosol radiative feedback
206	REAL, INTENT(in) :: cldfra(KLON,KLEV), cldemi(KLON,KLEV), cldtaupd(KLON,KLEV)
207	REAL, INTENT(in) :: tau_aero(KLON,KLEV,naero_grp,2) ! aerosol optical properties (see aeropt.F)
208	REAL, INTENT(in) :: piz_aero(KLON,KLEV,naero_grp,2) ! aerosol optical properties (see aeropt.F)
209	REAL, INTENT(in) :: cg_aero(KLON,KLEV,naero_grp,2) ! aerosol optical properties (see aeropt.F)
210	!--OB
211	REAL, INTENT(in) :: tau_aero_sw_rrtm(KLON,KLEV,2,NSW) ! aerosol optical properties RRTM
212	REAL, INTENT(in) :: piz_aero_sw_rrtm(KLON,KLEV,2,NSW) ! aerosol optical properties RRTM
213	REAL, INTENT(in) :: cg_aero_sw_rrtm(KLON,KLEV,2,NSW) ! aerosol optical properties RRTM
214	!--OB fin
215
216	!--C. Kleinschmitt
217	#ifdef CPP_RRTM
218	REAL, INTENT(in) :: tau_aero_lw_rrtm(KLON,KLEV,2,NLW) ! LW aerosol optical properties RRTM
219	#else
220	REAL, INTENT(in) :: tau_aero_lw_rrtm(KLON,KLEV,2,nbands_lw_rrtm)
221	#endif
222	!--C. Kleinschmitt end
223
224	REAL, INTENT(in) :: cldtaupi(KLON,KLEV) ! cloud optical thickness for pre-industrial aerosol concentrations
225	REAL, INTENT(in) :: qsat(klon,klev) ! Variable pour iflag_rrtm=1
226	REAL, INTENT(in) :: flwc(klon,klev) ! Variable pour iflag_rrtm=1
227	REAL, INTENT(in) :: fiwc(klon,klev) ! Variable pour iflag_rrtm=1
228	REAL, INTENT(in) :: ref_liq(klon,klev) ! cloud droplet radius present-day from newmicro
229	REAL, INTENT(in) :: ref_ice(klon,klev) ! ice crystal radius present-day from newmicro
230	REAL, INTENT(in) :: ref_liq_pi(klon,klev) ! cloud droplet radius pre-industrial from newmicro
231	REAL, INTENT(in) :: ref_ice_pi(klon,klev) ! ice crystal radius pre-industrial from newmicro
232
233	! Output arguments
234	REAL, INTENT(out) :: heat(KLON,KLEV), cool(KLON,KLEV)
235	REAL, INTENT(out) :: heat0(KLON,KLEV), cool0(KLON,KLEV)
236	REAL, INTENT(out) :: heat_volc(KLON,KLEV), cool_volc(KLON,KLEV) !NL
237	REAL, INTENT(out) :: topsw(KLON), toplw(KLON)
238	REAL, INTENT(out) :: solsw(KLON), sollw(KLON), albpla(KLON), solswfdiff(KLON)
239	REAL, INTENT(out) :: topsw0(KLON), toplw0(KLON), solsw0(KLON), sollw0(KLON)
240	REAL, INTENT(out) :: sollwdown(KLON)
241	REAL, INTENT(out) :: swdn(KLON,kflev+1),swdn0(KLON,kflev+1), swdnc0(KLON,kflev+1)
242	REAL, INTENT(out) :: swup(KLON,kflev+1),swup0(KLON,kflev+1), swupc0(KLON,kflev+1)
243	REAL, INTENT(out) :: lwdn(KLON,kflev+1),lwdn0(KLON,kflev+1), lwdnc0(KLON,kflev+1)
244	REAL, INTENT(out) :: lwup(KLON,kflev+1),lwup0(KLON,kflev+1), lwupc0(KLON,kflev+1)
245	REAL, INTENT(out) :: topswad_aero(KLON), solswad_aero(KLON) ! output: aerosol direct forcing at TOA and surface
246	REAL, INTENT(out) :: topswai_aero(KLON), solswai_aero(KLON) ! output: aerosol indirect forcing atTOA and surface
247	REAL, INTENT(out) :: toplwad_aero(KLON), sollwad_aero(KLON) ! output: LW aerosol direct forcing at TOA and surface
248	REAL, INTENT(out) :: toplwai_aero(KLON), sollwai_aero(KLON) ! output: LW aerosol indirect forcing atTOA and surface
249	REAL, DIMENSION(klon), INTENT(out) :: topswad0_aero
250	REAL, DIMENSION(klon), INTENT(out) :: solswad0_aero
251	REAL, DIMENSION(klon), INTENT(out) :: toplwad0_aero
252	REAL, DIMENSION(klon), INTENT(out) :: sollwad0_aero
253	REAL, DIMENSION(kdlon,9), INTENT(out) :: topsw_aero
254	REAL, DIMENSION(kdlon,9), INTENT(out) :: topsw0_aero
255	REAL, DIMENSION(kdlon,9), INTENT(out) :: solsw_aero
256	REAL, DIMENSION(kdlon,9), INTENT(out) :: solsw0_aero
257	REAL, DIMENSION(kdlon,3), INTENT(out) :: topswcf_aero
258	REAL, DIMENSION(kdlon,3), INTENT(out) :: solswcf_aero
259	REAL, DIMENSION(kdlon,kflev+1), INTENT(out) :: ZSWFT0_i
260	REAL, DIMENSION(kdlon,kflev+1), INTENT(out) :: ZLWFT0_i
261
262	! Local variables
263	REAL(KIND=8) ZFSUP(KDLON,KFLEV+1)
264	REAL(KIND=8) ZFSDN(KDLON,KFLEV+1)
265	REAL(KIND=8) ZFSUP0(KDLON,KFLEV+1)
266	REAL(KIND=8) ZFSDN0(KDLON,KFLEV+1)
267	REAL(KIND=8) ZFSUPC0(KDLON,KFLEV+1)
268	REAL(KIND=8) ZFSDNC0(KDLON,KFLEV+1)
269	REAL(KIND=8) ZFLUP(KDLON,KFLEV+1)
270	REAL(KIND=8) ZFLDN(KDLON,KFLEV+1)
271	REAL(KIND=8) ZFLUP0(KDLON,KFLEV+1)
272	REAL(KIND=8) ZFLDN0(KDLON,KFLEV+1)
273	REAL(KIND=8) ZFLUPC0(KDLON,KFLEV+1)
274	REAL(KIND=8) ZFLDNC0(KDLON,KFLEV+1)
275	REAL(KIND=8) zx_alpha1, zx_alpha2
276	INTEGER k, kk, i, j, iof, nb_gr
277	INTEGER ist,iend,ktdia,kmode
278	REAL(KIND=8) PSCT
279	REAL(KIND=8) PALBD(kdlon,2), PALBP(kdlon,2)
280	! MPL 06.01.09: pour RRTM, creation de PALBD_NEW et PALBP_NEW
281	! avec NSW en deuxieme dimension
282	REAL(KIND=8) PALBD_NEW(kdlon,NSW), PALBP_NEW(kdlon,NSW)
283	REAL(KIND=8) PEMIS(kdlon), PDT0(kdlon), PVIEW(kdlon)
284	REAL(KIND=8) PPSOL(kdlon), PDP(kdlon,KLEV)
285	REAL(KIND=8) PTL(kdlon,kflev+1), PPMB(kdlon,kflev+1)
286	REAL(KIND=8) PTAVE(kdlon,kflev)
287	REAL(KIND=8) PWV(kdlon,kflev), PQS(kdlon,kflev)
288
289	REAL(kind=8) POZON(kdlon, kflev, size(wo, 3)) ! mass fraction of ozone
290	! "POZON(:, :, 1)" is for the average day-night field,
291	! "POZON(:, :, 2)" is for daylight time.
292	!!!!! Modif MPL 6.01.09 avec RRTM, on passe de 5 a 6
293	REAL(KIND=8) PAER(kdlon,kflev,6)
294	REAL(KIND=8) PCLDLD(kdlon,kflev)
295	REAL(KIND=8) PCLDLU(kdlon,kflev)
296	REAL(KIND=8) PCLDSW(kdlon,kflev)
297	REAL(KIND=8) PTAU(kdlon,2,kflev)
298	REAL(KIND=8) POMEGA(kdlon,2,kflev)
299	REAL(KIND=8) PCG(kdlon,2,kflev)
300	REAL(KIND=8) zfract(kdlon), zrmu0(kdlon), zdist
301	REAL(KIND=8) zheat(kdlon,kflev), zcool(kdlon,kflev)
302	REAL(KIND=8) zheat0(kdlon,kflev), zcool0(kdlon,kflev)
303	REAL(KIND=8) zheat_volc(kdlon,kflev), zcool_volc(kdlon,kflev) !NL
304	REAL(KIND=8) ztopsw(kdlon), ztoplw(kdlon)
305	REAL(KIND=8) zsolsw(kdlon), zsollw(kdlon), zalbpla(kdlon), zsolswfdiff(kdlon)
306	REAL(KIND=8) zsollwdown(kdlon)
307	REAL(KIND=8) ztopsw0(kdlon), ztoplw0(kdlon)
308	REAL(KIND=8) zsolsw0(kdlon), zsollw0(kdlon)
309	REAL(KIND=8) zznormcp
310	REAL(KIND=8) tauaero(kdlon,kflev,naero_grp,2) ! aer opt properties
311	REAL(KIND=8) pizaero(kdlon,kflev,naero_grp,2)
312	REAL(KIND=8) cgaero(kdlon,kflev,naero_grp,2)
313	REAL(KIND=8) PTAUA(kdlon,2,kflev) ! present-day value of cloud opt thickness (PTAU is pre-industrial value), local use
314	REAL(KIND=8) POMEGAA(kdlon,2,kflev) ! dito for single scatt albedo
315	REAL(KIND=8) ztopswadaero(kdlon), zsolswadaero(kdlon) ! Aerosol direct forcing at TOAand surface
316	REAL(KIND=8) ztopswad0aero(kdlon), zsolswad0aero(kdlon) ! Aerosol direct forcing at TOAand surface
317	REAL(KIND=8) ztopswaiaero(kdlon), zsolswaiaero(kdlon) ! dito, indirect
318	!--NL
319	REAL(KIND=8) zswadaero(kdlon,kflev+1) ! SW Aerosol direct forcing
320	REAL(KIND=8) zlwadaero(kdlon,kflev+1) ! LW Aerosol direct forcing
321	!-LW by CK
322	REAL(KIND=8) ztoplwadaero(kdlon), zsollwadaero(kdlon) ! LW Aerosol direct forcing at TOAand surface
323	REAL(KIND=8) ztoplwad0aero(kdlon), zsollwad0aero(kdlon) ! LW Aerosol direct forcing at TOAand surface
324	REAL(KIND=8) ztoplwaiaero(kdlon), zsollwaiaero(kdlon) ! dito, indirect
325	!-end
326	REAL(KIND=8) ztopsw_aero(kdlon,9), ztopsw0_aero(kdlon,9)
327	REAL(KIND=8) zsolsw_aero(kdlon,9), zsolsw0_aero(kdlon,9)
328	REAL(KIND=8) ztopswcf_aero(kdlon,3), zsolswcf_aero(kdlon,3)
329	! real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2 deje declare dans physiq.F MPL 20130618
330	!MPL input supplementaires pour RECMWFL
331	! flwc, fiwc = Liquid Water Content & Ice Water Content (kg/kg)
332	REAL(KIND=8) GEMU(klon)
333	!MPL input RECMWFL:
334	! Tableaux aux niveaux inverses pour respecter convention Arpege
335	REAL(KIND=8) ref_liq_i(klon,klev) ! cloud droplet radius present-day from newmicro (inverted)
336	REAL(KIND=8) ref_ice_i(klon,klev) ! ice crystal radius present-day from newmicro (inverted)
337	!--OB
338	REAL(KIND=8) ref_liq_pi_i(klon,klev) ! cloud droplet radius pre-industrial from newmicro (inverted)
339	REAL(KIND=8) ref_ice_pi_i(klon,klev) ! ice crystal radius pre-industrial from newmicro (inverted)
340	!--end OB
341	REAL(KIND=8) paprs_i(klon,klev+1)
342	REAL(KIND=8) pplay_i(klon,klev)
343	REAL(KIND=8) cldfra_i(klon,klev)
344	REAL(KIND=8) POZON_i(kdlon,kflev, size(wo, 3)) ! mass fraction of ozone
345	! "POZON(:, :, 1)" is for the average day-night field,
346	! "POZON(:, :, 2)" is for daylight time.
347	!!!!! Modif MPL 6.01.09 avec RRTM, on passe de 5 a 6
348	REAL(KIND=8) PAER_i(kdlon,kflev,6)
349	REAL(KIND=8) PDP_i(klon,klev)
350	REAL(KIND=8) t_i(klon,klev),q_i(klon,klev),qsat_i(klon,klev)
351	REAL(KIND=8) flwc_i(klon,klev),fiwc_i(klon,klev)
352	!MPL output RECMWFL:
353	REAL(KIND=8) ZEMTD (klon,klev+1),ZEMTD_i (klon,klev+1)
354	REAL(KIND=8) ZEMTU (klon,klev+1),ZEMTU_i (klon,klev+1)
355	REAL(KIND=8) ZTRSO (klon,klev+1),ZTRSO_i (klon,klev+1)
356	REAL(KIND=8) ZTH (klon,klev+1),ZTH_i (klon,klev+1)
357	REAL(KIND=8) ZCTRSO(klon,2)
358	REAL(KIND=8) ZCEMTR(klon,2)
359	REAL(KIND=8) ZTRSOD(klon)
360	REAL(KIND=8) ZLWFC (klon,2)
361	REAL(KIND=8) ZLWFT (klon,klev+1),ZLWFT_i (klon,klev+1)
362	REAL(KIND=8) ZSWFC (klon,2)
363	REAL(KIND=8) ZSWFT (klon,klev+1),ZSWFT_i (klon,klev+1)
364	REAL(KIND=8) ZFLUCDWN_i(klon,klev+1),ZFLUCUP_i(klon,klev+1)
365	REAL(KIND=8) PPIZA_TOT(klon,klev,NSW)
366	REAL(KIND=8) PCGA_TOT(klon,klev,NSW)
367	REAL(KIND=8) PTAU_TOT(klon,klev,NSW)
368	REAL(KIND=8) PPIZA_NAT(klon,klev,NSW)
369	REAL(KIND=8) PCGA_NAT(klon,klev,NSW)
370	REAL(KIND=8) PTAU_NAT(klon,klev,NSW)
371	#ifdef CPP_RRTM
372	REAL(KIND=8) PTAU_LW_TOT(klon,klev,NLW)
373	REAL(KIND=8) PTAU_LW_NAT(klon,klev,NLW)
374	#endif
375	REAL(KIND=8) PSFSWDIR(klon,NSW)
376	REAL(KIND=8) PSFSWDIF(klon,NSW)
377	REAL(KIND=8) PFSDNN(klon)
378	REAL(KIND=8) PFSDNV(klon)
379	!MPL On ne redefinit pas les tableaux ZFLUX,ZFLUC,
380	!MPL ZFSDWN,ZFCDWN,ZFSUP,ZFCUP car ils existent deja
381	!MPL sous les noms de ZFLDN,ZFLDN0,ZFLUP,ZFLUP0,
382	!MPL ZFSDN,ZFSDN0,ZFSUP,ZFSUP0
383	REAL(KIND=8) ZFLUX_i (klon,2,klev+1)
384	REAL(KIND=8) ZFLUC_i (klon,2,klev+1)
385	REAL(KIND=8) ZFSDWN_i (klon,klev+1)
386	REAL(KIND=8) ZFCDWN_i (klon,klev+1)
387	REAL(KIND=8) ZFCCDWN_i (klon,klev+1)
388	REAL(KIND=8) ZFSUP_i (klon,klev+1)
389	REAL(KIND=8) ZFCUP_i (klon,klev+1)
390	REAL(KIND=8) ZFCCUP_i (klon,klev+1)
391	REAL(KIND=8) ZFLCCDWN_i (klon,klev+1)
392	REAL(KIND=8) ZFLCCUP_i (klon,klev+1)
393	! 3 lignes suivantes a activer pour CCMVAL (MPL 20100412)
394	! REAL(KIND=8) RSUN(3,2)
395	! REAL(KIND=8) SUN(3)
396	! REAL(KIND=8) SUN_FRACT(2)
397	REAL, PARAMETER:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2
398	CHARACTER (LEN=80) :: abort_message
399	CHARACTER (LEN=80) :: modname='radlwsw_m'
400
401	REAL zdir, zdif
402
403	CALL assert(size(wo, 1) == klon, size(wo, 2) == klev, "radlwsw wo")
404	! initialisation
405	ist=1
406	iend=klon
407	ktdia=1
408	kmode=ist
409	tauaero(:,:,:,:)=0.
410	pizaero(:,:,:,:)=0.
411	cgaero(:,:,:,:)=0.
412	lldebug=.FALSE.
413
414	ztopsw_aero(:,:) = 0. !ym missing init : warning : not initialized in SW_AEROAR4
415	ztopsw0_aero(:,:) = 0. !ym missing init : warning : not initialized in SW_AEROAR4
416	zsolsw_aero(:,:) = 0. !ym missing init : warning : not initialized in SW_AEROAR4
417	zsolsw0_aero(:,:) = 0. !ym missing init : warning : not initialized in SW_AEROAR4
418
419	ZTOPSWADAERO(:) = 0. !ym missing init
420	ZSOLSWADAERO(:) = 0. !ym missing init
421	ZTOPSWAD0AERO(:) = 0. !ym missing init
422	ZSOLSWAD0AERO(:) = 0. !ym missing init
423	ZTOPSWAIAERO(:) = 0. !ym missing init
424	ZSOLSWAIAERO(:) = 0. !ym missing init
425	ZTOPSWCF_AERO(:,:)= 0.!ym missing init
426	ZSOLSWCF_AERO(:,:) =0. !ym missing init
427
428	!
429	!-------------------------------------------
430	nb_gr = KLON / kdlon
431	IF (nb_gr*kdlon .NE. KLON) THEN
432	PRINT*, "kdlon mauvais:", KLON, kdlon, nb_gr
433	call abort_physic("radlwsw", "", 1)
434	ENDIF
435	IF (kflev .NE. KLEV) THEN
436	PRINT*, "kflev differe de KLEV, kflev, KLEV"
437	call abort_physic("radlwsw", "", 1)
438	ENDIF
439	!-------------------------------------------
440	DO k = 1, KLEV
441	DO i = 1, KLON
442	heat(i,k)=0.
443	cool(i,k)=0.
444	heat_volc(i,k)=0. !NL
445	cool_volc(i,k)=0. !NL
446	heat0(i,k)=0.
447	cool0(i,k)=0.
448	ENDDO
449	ENDDO
450	!
451	zdist = dist
452	!
453	PSCT = solaire/zdist/zdist
454
455	IF (type_trac == 'repr') THEN
456	#ifdef REPROBUS
457	IF (iflag_rrtm==0) THEN
458	IF (ok_SUNTIME) PSCT = solaireTIME/zdist/zdist
459	print,'Constante solaire: ',PSCTzdist*zdist
460	ENDIF
461	#endif
462	ENDIF
463
464	DO j = 1, nb_gr
465	iof = kdlon*(j-1)
466	DO i = 1, kdlon
467	zfract(i) = fract(iof+i)
468	! zfract(i) = 1. !!!!!! essai MPL 19052010
469	zrmu0(i) = rmu0(iof+i)
470
471
472	!albedo SB >>>
473	!
474	IF (iflag_rrtm==0) THEN
475	!
476	PALBD(i,1)=alb_dif(iof+i,1)
477	PALBD(i,2)=alb_dif(iof+i,2)
478	PALBP(i,1)=alb_dir(iof+i,1)
479	PALBP(i,2)=alb_dir(iof+i,2)
480	!
481	ELSEIF (iflag_rrtm==1) THEn
482	!
483	DO kk=1,NSW
484	PALBD_NEW(i,kk)=alb_dif(iof+i,kk)
485	PALBP_NEW(i,kk)=alb_dir(iof+i,kk)
486	ENDDO
487	!
488	ENDIF
489	!albedo SB <<<
490
491
492	PEMIS(i) = 1.0 !!!!! A REVOIR (MPL)
493	PVIEW(i) = 1.66
494	PPSOL(i) = paprs(iof+i,1)
495	zx_alpha1 = (paprs(iof+i,1)-pplay(iof+i,2))/(pplay(iof+i,1)-pplay(iof+i,2))
496	zx_alpha2 = 1.0 - zx_alpha1
497	PTL(i,1) = t(iof+i,1) * zx_alpha1 + t(iof+i,2) * zx_alpha2
498	PTL(i,KLEV+1) = t(iof+i,KLEV)
499	PDT0(i) = tsol(iof+i) - PTL(i,1)
500	ENDDO
501	DO k = 2, kflev
502	DO i = 1, kdlon
503	PTL(i,k) = (t(iof+i,k)+t(iof+i,k-1))*0.5
504	ENDDO
505	ENDDO
506	DO k = 1, kflev
507	DO i = 1, kdlon
508	PDP(i,k) = paprs(iof+i,k)-paprs(iof+i,k+1)
509	PTAVE(i,k) = t(iof+i,k)
510	PWV(i,k) = MAX (q(iof+i,k), 1.0e-12)
511	PQS(i,k) = PWV(i,k)
512	! Confert from column density of ozone in a cell, in kDU, to a mass fraction
513	POZON(i,k, :) = wo(iof+i, k, :) * RG * dobson_u * 1e3 &
514	/ (paprs(iof+i, k) - paprs(iof+i, k+1))
515	! A activer pour CCMVAL on prend l'ozone impose (MPL 07042010)
516	! POZON(i,k,:) = wo(i,k,:)
517	! print *,'RADLWSW: POZON',k, POZON(i,k,1)
518	PCLDLD(i,k) = cldfra(iof+i,k)*cldemi(iof+i,k)
519	PCLDLU(i,k) = cldfra(iof+i,k)*cldemi(iof+i,k)
520	PCLDSW(i,k) = cldfra(iof+i,k)
521	PTAU(i,1,k) = MAX(cldtaupi(iof+i,k), 1.0e-05)! 1e-12 serait instable
522	PTAU(i,2,k) = MAX(cldtaupi(iof+i,k), 1.0e-05)! pour 32-bit machines
523	POMEGA(i,1,k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAU(i,1,k))
524	POMEGA(i,2,k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAU(i,2,k))
525	PCG(i,1,k) = 0.865
526	PCG(i,2,k) = 0.910
527	!-
528	! Introduced for aerosol indirect forcings.
529	! The following values use the cloud optical thickness calculated from
530	! present-day aerosol concentrations whereas the quantities without the
531	! "A" at the end are for pre-industial (natural-only) aerosol concentrations
532	!
533	PTAUA(i,1,k) = MAX(cldtaupd(iof+i,k), 1.0e-05)! 1e-12 serait instable
534	PTAUA(i,2,k) = MAX(cldtaupd(iof+i,k), 1.0e-05)! pour 32-bit machines
535	POMEGAA(i,1,k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAUA(i,1,k))
536	POMEGAA(i,2,k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAUA(i,2,k))
537	ENDDO
538	ENDDO
539
540	IF (type_trac == 'repr') THEN
541	#ifdef REPROBUS
542	ndimozon = size(wo, 3)
543	CALL RAD_INTERACTIF(POZON,iof)
544	#endif
545	ENDIF
546	!
547	DO k = 1, kflev+1
548	DO i = 1, kdlon
549	PPMB(i,k) = paprs(iof+i,k)/100.0
550	ENDDO
551	ENDDO
552	!
553	!!!!! Modif MPL 6.01.09 avec RRTM, on passe de 5 a 6
554	DO kk = 1, 6
555	DO k = 1, kflev
556	DO i = 1, kdlon
557	PAER(i,k,kk) = 1.0E-15 !!!!! A REVOIR (MPL)
558	ENDDO
559	ENDDO
560	ENDDO
561	DO k = 1, kflev
562	DO i = 1, kdlon
563	tauaero(i,k,:,1)=tau_aero(iof+i,k,:,1)
564	pizaero(i,k,:,1)=piz_aero(iof+i,k,:,1)
565	cgaero(i,k,:,1) =cg_aero(iof+i,k,:,1)
566	tauaero(i,k,:,2)=tau_aero(iof+i,k,:,2)
567	pizaero(i,k,:,2)=piz_aero(iof+i,k,:,2)
568	cgaero(i,k,:,2) =cg_aero(iof+i,k,:,2)
569	ENDDO
570	ENDDO
571	!
572	!===== iflag_rrtm ================================================
573	!
574	IF (iflag_rrtm == 0) THEN !!!! remettre 0 juste pour tester l'ancien rayt via rrtm
575	!
576	!--- Mise a zero des tableaux output du rayonnement LW-AR4 ----------
577	DO k = 1, kflev+1
578	DO i = 1, kdlon
579	! print *,'RADLWSW: boucle mise a zero i k',i,k
580	ZFLUP(i,k)=0.
581	ZFLDN(i,k)=0.
582	ZFLUP0(i,k)=0.
583	ZFLDN0(i,k)=0.
584	ZLWFT0_i(i,k)=0.
585	ZFLUCUP_i(i,k)=0.
586	ZFLUCDWN_i(i,k)=0.
587	ENDDO
588	ENDDO
589	DO k = 1, kflev
590	DO i = 1, kdlon
591	zcool(i,k)=0.
592	zcool_volc(i,k)=0. !NL
593	zcool0(i,k)=0.
594	ENDDO
595	ENDDO
596	DO i = 1, kdlon
597	ztoplw(i)=0.
598	zsollw(i)=0.
599	ztoplw0(i)=0.
600	zsollw0(i)=0.
601	zsollwdown(i)=0.
602	ENDDO
603	! Old radiation scheme, used for AR4 runs
604	! average day-night ozone for longwave
605	CALL LW_LMDAR4(&
606	PPMB, PDP,&
607	PPSOL,PDT0,PEMIS,&
608	PTL, PTAVE, PWV, POZON(:, :, 1), PAER,&
609	PCLDLD,PCLDLU,&
610	PVIEW,&
611	zcool, zcool0,&
612	ztoplw,zsollw,ztoplw0,zsollw0,&
613	zsollwdown,&
614	ZFLUP, ZFLDN, ZFLUP0,ZFLDN0)
615	!----- Mise a zero des tableaux output du rayonnement SW-AR4
616	DO k = 1, kflev+1
617	DO i = 1, kdlon
618	ZFSUP(i,k)=0.
619	ZFSDN(i,k)=0.
620	ZFSUP0(i,k)=0.
621	ZFSDN0(i,k)=0.
622	ZFSUPC0(i,k)=0.
623	ZFSDNC0(i,k)=0.
624	ZFLUPC0(i,k)=0.
625	ZFLDNC0(i,k)=0.
626	ZSWFT0_i(i,k)=0.
627	ZFCUP_i(i,k)=0.
628	ZFCDWN_i(i,k)=0.
629	ZFCCUP_i(i,k)=0.
630	ZFCCDWN_i(i,k)=0.
631	ZFLCCUP_i(i,k)=0.
632	ZFLCCDWN_i(i,k)=0.
633	zswadaero(i,k)=0. !--NL
634	ENDDO
635	ENDDO
636	DO k = 1, kflev
637	DO i = 1, kdlon
638	zheat(i,k)=0.
639	zheat_volc(i,k)=0.
640	zheat0(i,k)=0.
641	ENDDO
642	ENDDO
643	DO i = 1, kdlon
644	zalbpla(i)=0.
645	ztopsw(i)=0.
646	zsolsw(i)=0.
647	ztopsw0(i)=0.
648	zsolsw0(i)=0.
649	ztopswadaero(i)=0.
650	zsolswadaero(i)=0.
651	ztopswaiaero(i)=0.
652	zsolswaiaero(i)=0.
653	ENDDO
654
655	!--fraction of diffuse radiation in surface SW downward radiation
656	!--not computed with old radiation scheme
657	zsolswfdiff(:) = -999.999
658
659	! print *,'Avant SW_LMDAR4: PSCT zrmu0 zfract',PSCT, zrmu0, zfract
660	! daylight ozone, if we have it, for short wave
661	CALL SW_AEROAR4(PSCT, zrmu0, zfract,&
662	PPMB, PDP,&
663	PPSOL, PALBD, PALBP,&
664	PTAVE, PWV, PQS, POZON(:, :, size(wo, 3)), PAER,&
665	PCLDSW, PTAU, POMEGA, PCG,&
666	zheat, zheat0,&
667	zalbpla,ztopsw,zsolsw,ztopsw0,zsolsw0,&
668	ZFSUP,ZFSDN,ZFSUP0,ZFSDN0,&
669	tauaero, pizaero, cgaero, &
670	PTAUA, POMEGAA,&
671	ztopswadaero,zsolswadaero,&
672	ztopswad0aero,zsolswad0aero,&
673	ztopswaiaero,zsolswaiaero, &
674	ztopsw_aero,ztopsw0_aero,&
675	zsolsw_aero,zsolsw0_aero,&
676	ztopswcf_aero,zsolswcf_aero, &
677	ok_ade, ok_aie, flag_aerosol,flag_aerosol_strat)
678
679	ZSWFT0_i(:,:) = ZFSDN0(:,:)-ZFSUP0(:,:)
680	ZLWFT0_i(:,:) =-ZFLDN0(:,:)-ZFLUP0(:,:)
681
682	DO i=1,kdlon
683	DO k=1,kflev+1
684	lwdn0 ( iof+i,k) = ZFLDN0 ( i,k)
685	lwdn ( iof+i,k) = ZFLDN ( i,k)
686	lwup0 ( iof+i,k) = ZFLUP0 ( i,k)
687	lwup ( iof+i,k) = ZFLUP ( i,k)
688	swdn0 ( iof+i,k) = ZFSDN0 ( i,k)
689	swdn ( iof+i,k) = ZFSDN ( i,k)
690	swup0 ( iof+i,k) = ZFSUP0 ( i,k)
691	swup ( iof+i,k) = ZFSUP ( i,k)
692	ENDDO
693	ENDDO
694	!
695	ELSE
696	#ifdef CPP_RRTM
697	! if (prt_level.gt.10)write(lunout,*)'CPP_RRTM=.T.'
698	!===== iflag_rrtm=1, on passe dans SW via RECMWFL ===============
699
700	DO k = 1, kflev+1
701	DO i = 1, kdlon
702	ZEMTD_i(i,k)=0.
703	ZEMTU_i(i,k)=0.
704	ZTRSO_i(i,k)=0.
705	ZTH_i(i,k)=0.
706	ZLWFT_i(i,k)=0.
707	ZSWFT_i(i,k)=0.
708	ZFLUX_i(i,1,k)=0.
709	ZFLUX_i(i,2,k)=0.
710	ZFLUC_i(i,1,k)=0.
711	ZFLUC_i(i,2,k)=0.
712	ZFSDWN_i(i,k)=0.
713	ZFCDWN_i(i,k)=0.
714	ZFCCDWN_i(i,k)=0.
715	ZFSUP_i(i,k)=0.
716	ZFCUP_i(i,k)=0.
717	ZFCCUP_i(i,k)=0.
718	ZFLCCDWN_i(i,k)=0.
719	ZFLCCUP_i(i,k)=0.
720	ENDDO
721	ENDDO
722	!
723	!--OB
724	!--aerosol TOT - anthropogenic+natural - index 2
725	!--aerosol NAT - natural only - index 1
726	!
727	DO i = 1, kdlon
728	DO k = 1, kflev
729	DO kk=1, NSW
730	!
731	PTAU_TOT(i,kflev+1-k,kk)=tau_aero_sw_rrtm(i,k,2,kk)
732	PPIZA_TOT(i,kflev+1-k,kk)=piz_aero_sw_rrtm(i,k,2,kk)
733	PCGA_TOT(i,kflev+1-k,kk)=cg_aero_sw_rrtm(i,k,2,kk)
734	!
735	PTAU_NAT(i,kflev+1-k,kk)=tau_aero_sw_rrtm(i,k,1,kk)
736	PPIZA_NAT(i,kflev+1-k,kk)=piz_aero_sw_rrtm(i,k,1,kk)
737	PCGA_NAT(i,kflev+1-k,kk)=cg_aero_sw_rrtm(i,k,1,kk)
738	!
739	ENDDO
740	ENDDO
741	ENDDO
742	!-end OB
743	!
744	!--C. Kleinschmitt
745	!--aerosol TOT - anthropogenic+natural - index 2
746	!--aerosol NAT - natural only - index 1
747	!
748	DO i = 1, kdlon
749	DO k = 1, kflev
750	DO kk=1, NLW
751	!
752	PTAU_LW_TOT(i,kflev+1-k,kk)=tau_aero_lw_rrtm(i,k,2,kk)
753	PTAU_LW_NAT(i,kflev+1-k,kk)=tau_aero_lw_rrtm(i,k,1,kk)
754	!
755	ENDDO
756	ENDDO
757	ENDDO
758	!-end C. Kleinschmitt
759	!
760	DO i = 1, kdlon
761	ZCTRSO(i,1)=0.
762	ZCTRSO(i,2)=0.
763	ZCEMTR(i,1)=0.
764	ZCEMTR(i,2)=0.
765	ZTRSOD(i)=0.
766	ZLWFC(i,1)=0.
767	ZLWFC(i,2)=0.
768	ZSWFC(i,1)=0.
769	ZSWFC(i,2)=0.
770	PFSDNN(i)=0.
771	PFSDNV(i)=0.
772	DO kk = 1, NSW
773	PSFSWDIR(i,kk)=0.
774	PSFSWDIF(i,kk)=0.
775	ENDDO
776	ENDDO
777	!----- Fin des mises a zero des tableaux output de RECMWF -------------------
778	! GEMU(1:klon)=sin(rlatd(1:klon))
779	! On met les donnees dans l'ordre des niveaux arpege
780	paprs_i(:,1)=paprs(:,klev+1)
781	DO k=1,klev
782	paprs_i(1:klon,k+1) =paprs(1:klon,klev+1-k)
783	pplay_i(1:klon,k) =pplay(1:klon,klev+1-k)
784	cldfra_i(1:klon,k) =cldfra(1:klon,klev+1-k)
785	PDP_i(1:klon,k) =PDP(1:klon,klev+1-k)
786	t_i(1:klon,k) =t(1:klon,klev+1-k)
787	q_i(1:klon,k) =q(1:klon,klev+1-k)
788	qsat_i(1:klon,k) =qsat(1:klon,klev+1-k)
789	flwc_i(1:klon,k) =flwc(1:klon,klev+1-k)
790	fiwc_i(1:klon,k) =fiwc(1:klon,klev+1-k)
791	ref_liq_i(1:klon,k) =ref_liq(1:klon,klev+1-k)
792	ref_ice_i(1:klon,k) =ref_ice(1:klon,klev+1-k)
793	!-OB
794	ref_liq_pi_i(1:klon,k) =ref_liq_pi(1:klon,klev+1-k)
795	ref_ice_pi_i(1:klon,k) =ref_ice_pi(1:klon,klev+1-k)
796	ENDDO
797	DO k=1,kflev
798	POZON_i(1:klon,k,:)=POZON(1:klon,kflev+1-k,:)
799	!!! POZON_i(1:klon,k)=POZON(1:klon,k) !!! on laisse 1=sol et klev=top
800	! print *,'Juste avant RECMWFL: k tsol temp',k,tsol,t(1,k)
801	!!!!!!! Modif MPL 6.01.09 avec RRTM, on passe de 5 a 6
802	DO i=1,6
803	PAER_i(1:klon,k,i)=PAER(1:klon,kflev+1-k,i)
804	ENDDO
805	ENDDO
806	! print *,'RADLWSW: avant RECMWFL, RI0,rmu0=',solaire,rmu0
807
808	! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
809	! La version ARPEGE1D utilise differentes valeurs de la constante
810	! solaire suivant le rayonnement utilise.
811	! A controler ...
812	! SOLAR FLUX AT THE TOP (/YOMPHY3/)
813	! introduce season correction
814	!--------------------------------------
815	! RII0 = RIP0
816	! IF(LRAYFM)
817	! RII0 = RIP0M ! =rip0m if Morcrette non-each time step call.
818	! IF(LRAYFM15)
819	! RII0 = RIP0M15 ! =rip0m if Morcrette non-each time step call.
820	RII0=solaire/zdist/zdist
821	!print*,'+++ radlwsw: solaire ,RII0',solaire,RII0
822	! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
823	! Ancien appel a RECMWF (celui du cy25)
824	! CALL RECMWF (ist , iend, klon , ktdia , klev , kmode ,
825	! s PALBD , PALBP , paprs_i , pplay_i , RCO2 , cldfra_i,
826	! s POZON_i , PAER_i , PDP_i , PEMIS , GEMU , rmu0,
827	! s q_i , qsat_i , fiwc_i , flwc_i , zmasq , t_i ,tsol,
828	! s ZEMTD_i , ZEMTU_i , ZTRSO_i ,
829	! s ZTH_i , ZCTRSO , ZCEMTR , ZTRSOD ,
830	! s ZLWFC , ZLWFT_i , ZSWFC , ZSWFT_i ,
831	! s ZFLUX_i , ZFLUC_i , ZFSDWN_i, ZFSUP_i , ZFCDWN_i,ZFCUP_i)
832	! s 'RECMWF ')
833	!
834	IF (lldebug) THEN
835	CALL writefield_phy('paprs_i',paprs_i,klev+1)
836	CALL writefield_phy('pplay_i',pplay_i,klev)
837	CALL writefield_phy('cldfra_i',cldfra_i,klev)
838	CALL writefield_phy('pozon_i',POZON_i,klev)
839	CALL writefield_phy('paer_i',PAER_i,klev)
840	CALL writefield_phy('pdp_i',PDP_i,klev)
841	CALL writefield_phy('q_i',q_i,klev)
842	CALL writefield_phy('qsat_i',qsat_i,klev)
843	CALL writefield_phy('fiwc_i',fiwc_i,klev)
844	CALL writefield_phy('flwc_i',flwc_i,klev)
845	CALL writefield_phy('t_i',t_i,klev)
846	CALL writefield_phy('palbd_new',PALBD_NEW,NSW)
847	CALL writefield_phy('palbp_new',PALBP_NEW,NSW)
848	ENDIF
849
850	! Nouvel appel a RECMWF (celui du cy32t0)
851	CALL RECMWF_AERO (ist , iend, klon , ktdia , klev , kmode ,&
852	PALBD_NEW,PALBP_NEW, paprs_i , pplay_i , RCO2 , cldfra_i,&
853	POZON_i , PAER_i , PDP_i , PEMIS , rmu0 ,&
854	q_i , qsat_i , fiwc_i , flwc_i , zmasq , t_i ,tsol,&
855	ref_liq_i, ref_ice_i, &
856	ref_liq_pi_i, ref_ice_pi_i, & ! rajoute par OB pour diagnostiquer effet indirect
857	ZEMTD_i , ZEMTU_i , ZTRSO_i ,&
858	ZTH_i , ZCTRSO , ZCEMTR , ZTRSOD ,&
859	ZLWFC , ZLWFT_i , ZSWFC , ZSWFT_i ,&
860	PSFSWDIR , PSFSWDIF, PFSDNN , PFSDNV ,&
861	PPIZA_TOT, PCGA_TOT,PTAU_TOT,&
862	PPIZA_NAT, PCGA_NAT,PTAU_NAT, & ! rajoute par OB pour diagnostiquer effet direct
863	PTAU_LW_TOT, PTAU_LW_NAT, & ! rajoute par C. Kleinschmitt
864	ZFLUX_i , ZFLUC_i ,&
865	ZFSDWN_i , ZFSUP_i , ZFCDWN_i, ZFCUP_i, ZFCCDWN_i, ZFCCUP_i, ZFLCCDWN_i, ZFLCCUP_i, &
866	ZTOPSWADAERO,ZSOLSWADAERO,& ! rajoute par OB pour diagnostics
867	ZTOPSWAD0AERO,ZSOLSWAD0AERO,&
868	ZTOPSWAIAERO,ZSOLSWAIAERO, &
869	ZTOPSWCF_AERO,ZSOLSWCF_AERO, &
870	ZSWADAERO, & !--NL
871	ZTOPLWADAERO,ZSOLLWADAERO,& ! rajoute par C. Kleinscmitt pour LW diagnostics
872	ZTOPLWAD0AERO,ZSOLLWAD0AERO,&
873	ZTOPLWAIAERO,ZSOLLWAIAERO, &
874	ZLWADAERO, & !--NL
875	ok_ade, ok_aie, ok_volcan, flag_aerosol,flag_aerosol_strat, flag_aer_feedback) ! flags aerosols
876
877	! print *,'RADLWSW: apres RECMWF'
878	IF (lldebug) THEN
879	CALL writefield_phy('zemtd_i',ZEMTD_i,klev+1)
880	CALL writefield_phy('zemtu_i',ZEMTU_i,klev+1)
881	CALL writefield_phy('ztrso_i',ZTRSO_i,klev+1)
882	CALL writefield_phy('zth_i',ZTH_i,klev+1)
883	CALL writefield_phy('zctrso',ZCTRSO,2)
884	CALL writefield_phy('zcemtr',ZCEMTR,2)
885	CALL writefield_phy('ztrsod',ZTRSOD,1)
886	CALL writefield_phy('zlwfc',ZLWFC,2)
887	CALL writefield_phy('zlwft_i',ZLWFT_i,klev+1)
888	CALL writefield_phy('zswfc',ZSWFC,2)
889	CALL writefield_phy('zswft_i',ZSWFT_i,klev+1)
890	CALL writefield_phy('psfswdir',PSFSWDIR,6)
891	CALL writefield_phy('psfswdif',PSFSWDIF,6)
892	CALL writefield_phy('pfsdnn',PFSDNN,1)
893	CALL writefield_phy('pfsdnv',PFSDNV,1)
894	CALL writefield_phy('ppiza_dst',PPIZA_TOT,klev)
895	CALL writefield_phy('pcga_dst',PCGA_TOT,klev)
896	CALL writefield_phy('ptaurel_dst',PTAU_TOT,klev)
897	CALL writefield_phy('zflux_i',ZFLUX_i,klev+1)
898	CALL writefield_phy('zfluc_i',ZFLUC_i,klev+1)
899	CALL writefield_phy('zfsdwn_i',ZFSDWN_i,klev+1)
900	CALL writefield_phy('zfsup_i',ZFSUP_i,klev+1)
901	CALL writefield_phy('zfcdwn_i',ZFCDWN_i,klev+1)
902	CALL writefield_phy('zfcup_i',ZFCUP_i,klev+1)
903	ENDIF
904	! --------- output RECMWFL
905	! ZEMTD (KPROMA,KLEV+1) ; TOTAL DOWNWARD LONGWAVE EMISSIVITY
906	! ZEMTU (KPROMA,KLEV+1) ; TOTAL UPWARD LONGWAVE EMISSIVITY
907	! ZTRSO (KPROMA,KLEV+1) ; TOTAL SHORTWAVE TRANSMISSIVITY
908	! ZTH (KPROMA,KLEV+1) ; HALF LEVEL TEMPERATURE
909	! ZCTRSO (KPROMA,2) ; CLEAR-SKY SHORTWAVE TRANSMISSIVITY
910	! ZCEMTR (KPROMA,2) ; CLEAR-SKY NET LONGWAVE EMISSIVITY
911	! ZTRSOD (KPROMA) ; TOTAL-SKY SURFACE SW TRANSMISSITY
912	! ZLWFC (KPROMA,2) ; CLEAR-SKY LONGWAVE FLUXES
913	! ZLWFT (KPROMA,KLEV+1) ; TOTAL-SKY LONGWAVE FLUXES
914	! ZSWFC (KPROMA,2) ; CLEAR-SKY SHORTWAVE FLUXES
915	! ZSWFT (KPROMA,KLEV+1) ; TOTAL-SKY SHORTWAVE FLUXES
916	! PPIZA_TOT (KPROMA,KLEV,NSW); Single scattering albedo of total aerosols
917	! PCGA_TOT (KPROMA,KLEV,NSW); Assymetry factor for total aerosols
918	! PTAU_TOT (KPROMA,KLEV,NSW); Optical depth of total aerosols
919	! PPIZA_NAT (KPROMA,KLEV,NSW); Single scattering albedo of natural aerosols
920	! PCGA_NAT (KPROMA,KLEV,NSW); Assymetry factor for natural aerosols
921	! PTAU_NAT (KPROMA,KLEV,NSW); Optical depth of natiral aerosols
922	! PTAU_LW_TOT (KPROMA,KLEV,NLW); LW Optical depth of total aerosols
923	! PTAU_LW_NAT (KPROMA,KLEV,NLW); LW Optical depth of natural aerosols
924	! PSFSWDIR (KPROMA,NSW) ;
925	! PSFSWDIF (KPROMA,NSW) ;
926	! PFSDNN (KPROMA) ;
927	! PFSDNV (KPROMA) ;
928	! ---------
929	! ---------
930	! On retablit l'ordre des niveaux lmd pour les tableaux de sortie
931	! D autre part, on multiplie les resultats SW par fract pour etre coherent
932	! avec l ancien rayonnement AR4. Si nuit, fract=0 donc pas de
933	! rayonnement SW. (MPL 260609)
934	DO k=0,klev
935	DO i=1,klon
936	ZEMTD(i,k+1) = ZEMTD_i(i,k+1)
937	ZEMTU(i,k+1) = ZEMTU_i(i,k+1)
938	ZTRSO(i,k+1) = ZTRSO_i(i,k+1)
939	ZTH(i,k+1) = ZTH_i(i,k+1)
940	! ZLWFT(i,k+1) = ZLWFT_i(i,klev+1-k)
941	! ZSWFT(i,k+1) = ZSWFT_i(i,klev+1-k)
942	ZFLUP(i,k+1) = ZFLUX_i(i,1,k+1)
943	ZFLDN(i,k+1) = ZFLUX_i(i,2,k+1)
944	ZFLUP0(i,k+1) = ZFLUC_i(i,1,k+1)
945	ZFLDN0(i,k+1) = ZFLUC_i(i,2,k+1)
946	ZFSDN(i,k+1) = ZFSDWN_i(i,k+1)*fract(i)
947	ZFSDN0(i,k+1) = ZFCDWN_i(i,k+1)*fract(i)
948	ZFSDNC0(i,k+1)= ZFCCDWN_i(i,k+1)*fract(i)
949	ZFSUP (i,k+1) = ZFSUP_i(i,k+1)*fract(i)
950	ZFSUP0(i,k+1) = ZFCUP_i(i,k+1)*fract(i)
951	ZFSUPC0(i,k+1)= ZFCCUP_i(i,k+1)*fract(i)
952	ZFLDNC0(i,k+1)= ZFLCCDWN_i(i,k+1)
953	ZFLUPC0(i,k+1)= ZFLCCUP_i(i,k+1)
954	IF (ok_volcan) THEN
955	ZSWADAERO(i,k+1)=ZSWADAERO(i,k+1)*fract(i) !--NL
956	ENDIF
957
958	! Nouveau calcul car visiblement ZSWFT et ZSWFC sont nuls dans RRTM cy32
959	! en sortie de radlsw.F90 - MPL 7.01.09
960	ZSWFT(i,k+1) = (ZFSDWN_i(i,k+1)-ZFSUP_i(i,k+1))*fract(i)
961	ZSWFT0_i(i,k+1) = (ZFCDWN_i(i,k+1)-ZFCUP_i(i,k+1))*fract(i)
962	! WRITE(*,'("FSDN FSUP FCDN FCUP: ",4E12.5)') ZFSDWN_i(i,k+1),&
963	! ZFSUP_i(i,k+1),ZFCDWN_i(i,k+1),ZFCUP_i(i,k+1)
964	ZLWFT(i,k+1) =-ZFLUX_i(i,2,k+1)-ZFLUX_i(i,1,k+1)
965	ZLWFT0_i(i,k+1)=-ZFLUC_i(i,2,k+1)-ZFLUC_i(i,1,k+1)
966	! print *,'FLUX2 FLUX1 FLUC2 FLUC1',ZFLUX_i(i,2,k+1),&
967	! & ZFLUX_i(i,1,k+1),ZFLUC_i(i,2,k+1),ZFLUC_i(i,1,k+1)
968	ENDDO
969	ENDDO
970
971	!--ajout OB
972	ZTOPSWADAERO(:) =ZTOPSWADAERO(:) *fract(:)
973	ZSOLSWADAERO(:) =ZSOLSWADAERO(:) *fract(:)
974	ZTOPSWAD0AERO(:)=ZTOPSWAD0AERO(:)*fract(:)
975	ZSOLSWAD0AERO(:)=ZSOLSWAD0AERO(:)*fract(:)
976	ZTOPSWAIAERO(:) =ZTOPSWAIAERO(:) *fract(:)
977	ZSOLSWAIAERO(:) =ZSOLSWAIAERO(:) *fract(:)
978	ZTOPSWCF_AERO(:,1)=ZTOPSWCF_AERO(:,1)*fract(:)
979	ZTOPSWCF_AERO(:,2)=ZTOPSWCF_AERO(:,2)*fract(:)
980	ZTOPSWCF_AERO(:,3)=ZTOPSWCF_AERO(:,3)*fract(:)
981	ZSOLSWCF_AERO(:,1)=ZSOLSWCF_AERO(:,1)*fract(:)
982	ZSOLSWCF_AERO(:,2)=ZSOLSWCF_AERO(:,2)*fract(:)
983	ZSOLSWCF_AERO(:,3)=ZSOLSWCF_AERO(:,3)*fract(:)
984
985	! print*,'SW_RRTM ZFSDN0 1 , klev:',ZFSDN0(1:klon,1),ZFSDN0(1:klon,klev)
986	! print*,'SW_RRTM ZFSUP0 1 , klev:',ZFSUP0(1:klon,1),ZFSUP0(1:klon,klev)
987	! print*,'SW_RRTM ZFSDN 1 , klev:',ZFSDN(1:klon,1),ZFSDN(1:klon,klev)
988	! print*,'SW_RRTM ZFSUP 1 , klev:',ZFSUP(1:klon,1),ZFSUP(1:klon,klev)
989	! print*,'OK1'
990	! ---------
991	! ---------
992	! On renseigne les champs LMDz, pour avoir la meme chose qu'en sortie de
993	! LW_LMDAR4 et SW_LMDAR4
994
995	!--fraction of diffuse radiation in surface SW downward radiation
996	DO i = 1, kdlon
997	IF (fract(i).GT.0.0) THEN
998	zdir=SUM(PSFSWDIR(i,:))
999	zdif=SUM(PSFSWDIF(i,:))
1000	zsolswfdiff(i) = zdif/(zdir+zdif)
1001	ELSE !--night
1002	zsolswfdiff(i) = 1.0
1003	ENDIF
1004	ENDDO
1005	!
1006	DO i = 1, kdlon
1007	zsolsw(i) = ZSWFT(i,1)
1008	zsolsw0(i) = ZSWFT0_i(i,1)
1009	! zsolsw0(i) = ZFSDN0(i,1) -ZFSUP0(i,1)
1010	ztopsw(i) = ZSWFT(i,klev+1)
1011	ztopsw0(i) = ZSWFT0_i(i,klev+1)
1012	! ztopsw0(i) = ZFSDN0(i,klev+1)-ZFSUP0(i,klev+1)
1013	!
1014	! zsollw(i) = ZFLDN(i,1) -ZFLUP(i,1)
1015	! zsollw0(i) = ZFLDN0(i,1) -ZFLUP0(i,1)
1016	! ztoplw(i) = ZFLDN(i,klev+1) -ZFLUP(i,klev+1)
1017	! ztoplw0(i) = ZFLDN0(i,klev+1)-ZFLUP0(i,klev+1)
1018	zsollw(i) = ZLWFT(i,1)
1019	zsollw0(i) = ZLWFT0_i(i,1)
1020	ztoplw(i) = ZLWFT(i,klev+1)*(-1)
1021	ztoplw0(i) = ZLWFT0_i(i,klev+1)*(-1)
1022	!
1023	IF (fract(i) == 0.) THEN
1024	!!!!! A REVOIR MPL (20090630) ca n a pas de sens quand fract=0
1025	! pas plus que dans le sw_AR4
1026	zalbpla(i) = 1.0e+39
1027	ELSE
1028	zalbpla(i) = ZFSUP(i,klev+1)/ZFSDN(i,klev+1)
1029	ENDIF
1030	!!! 5 juin 2015
1031	!!! Correction MP bug RRTM
1032	zsollwdown(i)= -1.*ZFLDN(i,1)
1033	ENDDO
1034	! print*,'OK2'
1035
1036	! extrait de SW_AR4
1037	! DO k = 1, KFLEV
1038	! kpl1 = k+1
1039	! DO i = 1, KDLON
1040	! PHEAT(i,k) = -(ZFSUP(i,kpl1)-ZFSUP(i,k)) -(ZFSDN(i,k)-ZFSDN(i,kpl1))
1041	! PHEAT(i,k) = PHEAT(i,k) * RDAY*RG/RCPD / PDP(i,k)
1042	! ZLWFT(klon,k),ZSWFT
1043
1044	DO k=1,kflev
1045	DO i=1,kdlon
1046	zheat(i,k)=(ZSWFT(i,k+1)-ZSWFT(i,k))RDAYRG/RCPD/PDP(i,k)
1047	zheat0(i,k)=(ZSWFT0_i(i,k+1)-ZSWFT0_i(i,k))RDAYRG/RCPD/PDP(i,k)
1048	zcool(i,k)=(ZLWFT(i,k)-ZLWFT(i,k+1))RDAYRG/RCPD/PDP(i,k)
1049	zcool0(i,k)=(ZLWFT0_i(i,k)-ZLWFT0_i(i,k+1))RDAYRG/RCPD/PDP(i,k)
1050	IF (ok_volcan) THEN
1051	zheat_volc(i,k)=(ZSWADAERO(i,k+1)-ZSWADAERO(i,k))*RG/RCPD/PDP(i,k) !NL
1052	zcool_volc(i,k)=(ZLWADAERO(i,k)-ZLWADAERO(i,k+1))*RG/RCPD/PDP(i,k) !NL
1053	ENDIF
1054	! print *,'heat cool heat0 cool0 ',zheat(i,k),zcool(i,k),zheat0(i,k),zcool0(i,k)
1055	! ZFLUCUP_i(i,k)=ZFLUC_i(i,1,k)
1056	! ZFLUCDWN_i(i,k)=ZFLUC_i(i,2,k)
1057	ENDDO
1058	ENDDO
1059	#else
1060	abort_message="You should compile with -rrtm if running with iflag_rrtm=1"
1061	call abort_physic(modname, abort_message, 1)
1062	#endif
1063	ENDIF ! iflag_rrtm
1064	!======================================================================
1065
1066	DO i = 1, kdlon
1067	topsw(iof+i) = ztopsw(i)
1068	toplw(iof+i) = ztoplw(i)
1069	solsw(iof+i) = zsolsw(i)
1070	solswfdiff(iof+i) = zsolswfdiff(i)
1071	sollw(iof+i) = zsollw(i)
1072	sollwdown(iof+i) = zsollwdown(i)
1073	DO k = 1, kflev+1
1074	lwdn0 ( iof+i,k) = ZFLDN0 ( i,k)
1075	lwdn ( iof+i,k) = ZFLDN ( i,k)
1076	lwup0 ( iof+i,k) = ZFLUP0 ( i,k)
1077	lwup ( iof+i,k) = ZFLUP ( i,k)
1078	ENDDO
1079	topsw0(iof+i) = ztopsw0(i)
1080	toplw0(iof+i) = ztoplw0(i)
1081	solsw0(iof+i) = zsolsw0(i)
1082	sollw0(iof+i) = zsollw0(i)
1083	albpla(iof+i) = zalbpla(i)
1084
1085	DO k = 1, kflev+1
1086	swdnc0( iof+i,k) = ZFSDNC0( i,k)
1087	swdn0 ( iof+i,k) = ZFSDN0 ( i,k)
1088	swdn ( iof+i,k) = ZFSDN ( i,k)
1089	swupc0( iof+i,k) = ZFSUPC0( i,k)
1090	swup0 ( iof+i,k) = ZFSUP0 ( i,k)
1091	swup ( iof+i,k) = ZFSUP ( i,k)
1092	lwdnc0( iof+i,k) = ZFLDNC0( i,k)
1093	lwupc0( iof+i,k) = ZFLUPC0( i,k)
1094	ENDDO
1095	ENDDO
1096	!-transform the aerosol forcings, if they have
1097	! to be calculated
1098	IF (ok_ade) THEN
1099	DO i = 1, kdlon
1100	topswad_aero(iof+i) = ztopswadaero(i)
1101	topswad0_aero(iof+i) = ztopswad0aero(i)
1102	solswad_aero(iof+i) = zsolswadaero(i)
1103	solswad0_aero(iof+i) = zsolswad0aero(i)
1104	! MS the following lines seem to be wrong, why is iof on right hand side???
1105	! topsw_aero(iof+i,:) = ztopsw_aero(iof+i,:)
1106	! topsw0_aero(iof+i,:) = ztopsw0_aero(iof+i,:)
1107	! solsw_aero(iof+i,:) = zsolsw_aero(iof+i,:)
1108	! solsw0_aero(iof+i,:) = zsolsw0_aero(iof+i,:)
1109	topsw_aero(iof+i,:) = ztopsw_aero(i,:)
1110	topsw0_aero(iof+i,:) = ztopsw0_aero(i,:)
1111	solsw_aero(iof+i,:) = zsolsw_aero(i,:)
1112	solsw0_aero(iof+i,:) = zsolsw0_aero(i,:)
1113	topswcf_aero(iof+i,:) = ztopswcf_aero(i,:)
1114	solswcf_aero(iof+i,:) = zsolswcf_aero(i,:)
1115	!-LW
1116	toplwad_aero(iof+i) = ztoplwadaero(i)
1117	toplwad0_aero(iof+i) = ztoplwad0aero(i)
1118	sollwad_aero(iof+i) = zsollwadaero(i)
1119	sollwad0_aero(iof+i) = zsollwad0aero(i)
1120	ENDDO
1121	ELSE
1122	DO i = 1, kdlon
1123	topswad_aero(iof+i) = 0.0
1124	solswad_aero(iof+i) = 0.0
1125	topswad0_aero(iof+i) = 0.0
1126	solswad0_aero(iof+i) = 0.0
1127	topsw_aero(iof+i,:) = 0.
1128	topsw0_aero(iof+i,:) =0.
1129	solsw_aero(iof+i,:) = 0.
1130	solsw0_aero(iof+i,:) = 0.
1131	!-LW
1132	toplwad_aero(iof+i) = 0.0
1133	sollwad_aero(iof+i) = 0.0
1134	toplwad0_aero(iof+i) = 0.0
1135	sollwad0_aero(iof+i) = 0.0
1136	ENDDO
1137	ENDIF
1138	IF (ok_aie) THEN
1139	DO i = 1, kdlon
1140	topswai_aero(iof+i) = ztopswaiaero(i)
1141	solswai_aero(iof+i) = zsolswaiaero(i)
1142	!-LW
1143	toplwai_aero(iof+i) = ztoplwaiaero(i)
1144	sollwai_aero(iof+i) = zsollwaiaero(i)
1145	ENDDO
1146	ELSE
1147	DO i = 1, kdlon
1148	topswai_aero(iof+i) = 0.0
1149	solswai_aero(iof+i) = 0.0
1150	!-LW
1151	toplwai_aero(iof+i) = 0.0
1152	sollwai_aero(iof+i) = 0.0
1153	ENDDO
1154	ENDIF
1155	DO k = 1, kflev
1156	DO i = 1, kdlon
1157	! scale factor to take into account the difference between
1158	! dry air and watter vapour scpecifi! heat capacity
1159	zznormcp=1.0+RVTMP2*PWV(i,k)
1160	heat(iof+i,k) = zheat(i,k)/zznormcp
1161	cool(iof+i,k) = zcool(i,k)/zznormcp
1162	heat0(iof+i,k) = zheat0(i,k)/zznormcp
1163	cool0(iof+i,k) = zcool0(i,k)/zznormcp
1164	IF(ok_volcan) THEN !NL
1165	heat_volc(iof+i,k) = zheat_volc(i,k)/zznormcp
1166	cool_volc(iof+i,k) = zcool_volc(i,k)/zznormcp
1167	ENDIF
1168	ENDDO
1169	ENDDO
1170
1171	ENDDO ! j = 1, nb_gr
1172
1173	END SUBROUTINE radlwsw
1174
1175	end module radlwsw_m

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

Download in other formats:

Original Format