Context Navigation

radlwsw_m.F90 @ 4703

Last change on this file since 4703 was 4677, checked in by idelkadi, 15 months ago

Implementation in the LMDZ code of the double call of the ECRAD radiative transfer code to estimate the 3D radiative effect of clouds.

This double call of Ecrad is controlled by the ok_3Deffect logic key.
If this key is enabled, 2 files of parameter configuration "namelists" for ECRAD are required at runtime: namelist_ecrad and namelist_ecrad_s2.
If this key is deactivated, the configuration and initialization part (reading namelist and netcdf files) is performed only once during simulation (1st call to ECRAD). Otherwise, configuration and initialization are performed each time Ecrad is called.

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.6 KB

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

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

Download in other formats:

Original Format