Context Navigation

guide_loc_mod.f90 @ 5272

Last change on this file since 5272 was 5272, checked in by abarral, 3 months ago
Turn paramet.h into a module
Property copyright set to Name of program: LMDZ Creation date: 1984 Version: LMDZ5 License: CeCILL version 2 Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539 See the license file in the root directory
File size: 86.1 KB

Rev	Line
[5084]	1	!
	2	! $Id$
	3	!
[1632]	4	MODULE guide_loc_mod
	5
	6	!=======================================================================
	7	! Auteur: F.Hourdin
	8	! F. Codron 01/09
	9	!=======================================================================
	10
[3995]	11	USE getparam, only: ini_getparam, fin_getparam, getpar
[1632]	12	USE Write_Field_loc
[5270]	13	USE netcdf, ONLY: nf90_nowrite, nf90_open, nf90_inq_varid, nf90_close, &
	14	nf90_inq_dimid, nf90_inquire_dimension, nf90_inq_dimid, &
	15	nf90_inquire_dimension, nf90_enddef, nf90_def_dim, nf90_put_var, nf90_noerr, nf90_close, nf90_inq_varid, &
	16	nf90_redef, nf90_write, nf90_unlimited, nf90_float, nf90_clobber, nf90_64bit_offset, nf90_float, &
	17	nf90_create, nf90_def_var, nf90_open
[1823]	18	USE parallel_lmdz
[3995]	19	USE pres2lev_mod, only: pres2lev
[1632]	20
[5272]	21	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	22	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
	23	IMPLICIT NONE
[1632]	24
	25	! ---------------------------------------------
[5249]	26	! Declarations des cles logiques et parametres
[1632]	27	! ---------------------------------------------
	28	INTEGER, PRIVATE, SAVE :: iguide_read,iguide_int,iguide_sav
	29	INTEGER, PRIVATE, SAVE :: nlevnc, guide_plevs
	30	LOGICAL, PRIVATE, SAVE :: guide_u,guide_v,guide_T,guide_Q,guide_P
[5249]	31	LOGICAL, PRIVATE, SAVE :: guide_hr,guide_teta
	32	LOGICAL, PRIVATE, SAVE :: guide_BL,guide_reg,guide_add,gamma4,guide_zon
[1632]	33	LOGICAL, PRIVATE, SAVE :: invert_p,invert_y,ini_anal
	34	LOGICAL, PRIVATE, SAVE :: guide_2D,guide_sav,guide_modele
[3795]	35	!FC
	36	LOGICAL, PRIVATE, SAVE :: convert_Pa
[5249]	37
[1632]	38	REAL, PRIVATE, SAVE :: tau_min_u,tau_max_u
	39	REAL, PRIVATE, SAVE :: tau_min_v,tau_max_v
	40	REAL, PRIVATE, SAVE :: tau_min_T,tau_max_T
	41	REAL, PRIVATE, SAVE :: tau_min_Q,tau_max_Q
	42	REAL, PRIVATE, SAVE :: tau_min_P,tau_max_P
	43
	44	REAL, PRIVATE, SAVE :: lat_min_g,lat_max_g
	45	REAL, PRIVATE, SAVE :: lon_min_g,lon_max_g
	46	REAL, PRIVATE, SAVE :: tau_lon,tau_lat
	47
[4246]	48	REAL, PRIVATE, SAVE :: plim_guide_BL
	49
[5249]	50	REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE :: alpha_u,alpha_v
	51	REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE :: alpha_T,alpha_Q
[1632]	52	REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE :: alpha_P,alpha_pcor
[5249]	53
[1632]	54	! ---------------------------------------------
	55	! Variables de guidage
	56	! ---------------------------------------------
	57	! Variables des fichiers de guidage
	58	REAL, ALLOCATABLE, DIMENSION(:,:,:), PRIVATE, SAVE :: unat1,unat2
	59	REAL, ALLOCATABLE, DIMENSION(:,:,:), PRIVATE, SAVE :: vnat1,vnat2
	60	REAL, ALLOCATABLE, DIMENSION(:,:,:), PRIVATE, SAVE :: tnat1,tnat2
	61	REAL, ALLOCATABLE, DIMENSION(:,:,:), PRIVATE, SAVE :: qnat1,qnat2
	62	REAL, ALLOCATABLE, DIMENSION(:,:,:), PRIVATE, SAVE :: pnat1,pnat2
	63	REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: psnat1,psnat2
	64	REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE :: apnc,bpnc
	65	! Variables aux dimensions du modele
	66	REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: ugui1,ugui2
	67	REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: vgui1,vgui2
	68	REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: tgui1,tgui2
	69	REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: qgui1,qgui2
	70	REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE :: psgui1,psgui2
[5249]	71
[3995]	72	INTEGER,SAVE,PRIVATE :: ijbu,ijbv,ijeu,ijev !,ijnu,ijnv
[1632]	73	INTEGER,SAVE,PRIVATE :: jjbu,jjbv,jjeu,jjev,jjnu,jjnv
	74
	75
	76	CONTAINS
	77	!=======================================================================
	78
	79	SUBROUTINE guide_init
[2134]	80
[2598]	81	USE control_mod, ONLY: day_step
	82	USE serre_mod, ONLY: grossismx
[2134]	83
[5271]	84	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
	85	IMPLICIT NONE
[5249]	86
[5271]	87
[1632]	88
[5272]	89
[1632]	90	INTEGER :: error,ncidpl,rid,rcod
	91	CHARACTER (len = 80) :: abort_message
	92	CHARACTER (len = 20) :: modname = 'guide_init'
[3984]	93	CHARACTER (len = 20) :: namedim
[1632]	94
	95	! ---------------------------------------------
[5249]	96	! Lecture des parametres:
[1632]	97	! ---------------------------------------------
[2263]	98	call ini_getparam("nudging_parameters_out.txt")
[1632]	99	! Variables guidees
	100	CALL getpar('guide_u',.true.,guide_u,'guidage de u')
	101	CALL getpar('guide_v',.true.,guide_v,'guidage de v')
	102	CALL getpar('guide_T',.true.,guide_T,'guidage de T')
	103	CALL getpar('guide_P',.true.,guide_P,'guidage de P')
	104	CALL getpar('guide_Q',.true.,guide_Q,'guidage de Q')
	105	CALL getpar('guide_hr',.true.,guide_hr,'guidage de Q par H.R')
	106	CALL getpar('guide_teta',.false.,guide_teta,'guidage de T par Teta')
	107
	108	CALL getpar('guide_add',.false.,guide_add,'for�age constant?')
	109	CALL getpar('guide_zon',.false.,guide_zon,'guidage moy zonale')
[2134]	110	if (guide_zon .and. abs(grossismx - 1.) > 0.01) &
	111	call abort_gcm("guide_init", &
	112	"zonal nudging requires grid regular in longitude", 1)
[1632]	113
	114	! Constantes de rappel. Unite : fraction de jour
	115	CALL getpar('tau_min_u',0.02,tau_min_u,'Cste de rappel min, u')
	116	CALL getpar('tau_max_u', 10.,tau_max_u,'Cste de rappel max, u')
	117	CALL getpar('tau_min_v',0.02,tau_min_v,'Cste de rappel min, v')
	118	CALL getpar('tau_max_v', 10.,tau_max_v,'Cste de rappel max, v')
	119	CALL getpar('tau_min_T',0.02,tau_min_T,'Cste de rappel min, T')
	120	CALL getpar('tau_max_T', 10.,tau_max_T,'Cste de rappel max, T')
	121	CALL getpar('tau_min_Q',0.02,tau_min_Q,'Cste de rappel min, Q')
	122	CALL getpar('tau_max_Q', 10.,tau_max_Q,'Cste de rappel max, Q')
	123	CALL getpar('tau_min_P',0.02,tau_min_P,'Cste de rappel min, P')
	124	CALL getpar('tau_max_P', 10.,tau_max_P,'Cste de rappel max, P')
	125	CALL getpar('gamma4',.false.,gamma4,'Zone sans rappel elargie')
	126	CALL getpar('guide_BL',.true.,guide_BL,'guidage dans C.Lim')
[4246]	127	CALL getpar('plim_guide_BL',85000.,plim_guide_BL,'BL top presnivs value')
	128
[1632]	129	! Sauvegarde du for�age
	130	CALL getpar('guide_sav',.false.,guide_sav,'sauvegarde guidage')
	131	CALL getpar('iguide_sav',4,iguide_sav,'freq. sauvegarde guidage')
	132	! frequences f>0: fx/jour; f<0: tous les f jours; f=0: 1 seule fois.
[5084]	133	IF (iguide_sav.GT.0) THEN
[2124]	134	iguide_sav=day_step/iguide_sav
	135	ELSE if (iguide_sav == 0) then
	136	iguide_sav = huge(0)
[1632]	137	ELSE
[2124]	138	iguide_sav=day_step*iguide_sav
[1632]	139	ENDIF
	140
	141	! Guidage regional seulement (sinon constant ou suivant le zoom)
	142	CALL getpar('guide_reg',.false.,guide_reg,'guidage regional')
	143	CALL getpar('lat_min_g',-90.,lat_min_g,'Latitude mini guidage ')
	144	CALL getpar('lat_max_g', 90.,lat_max_g,'Latitude maxi guidage ')
	145	CALL getpar('lon_min_g',-180.,lon_min_g,'longitude mini guidage ')
	146	CALL getpar('lon_max_g', 180.,lon_max_g,'longitude maxi guidage ')
	147	CALL getpar('tau_lat', 5.,tau_lat,'raideur lat guide regional ')
	148	CALL getpar('tau_lon', 5.,tau_lon,'raideur lon guide regional ')
	149
	150	! Parametres pour lecture des fichiers
	151	CALL getpar('iguide_read',4,iguide_read,'freq. lecture guidage')
	152	CALL getpar('iguide_int',4,iguide_int,'freq. interpolation vert')
[5084]	153	IF (iguide_int.EQ.0) THEN
[1632]	154	iguide_int=1
[5084]	155	ELSEIF (iguide_int.GT.0) THEN
[1632]	156	iguide_int=day_step/iguide_int
	157	ELSE
	158	iguide_int=day_step*iguide_int
	159	ENDIF
	160	CALL getpar('guide_plevs',0,guide_plevs,'niveaux pression fichiers guidage')
	161	! Pour compatibilite avec ancienne version avec guide_modele
	162	CALL getpar('guide_modele',.false.,guide_modele,'niveaux pression ap+bp*psol')
	163	IF (guide_modele) THEN
	164	guide_plevs=1
	165	ENDIF
[3795]	166	!FC
	167	CALL getpar('convert_Pa',.true.,convert_Pa,'Convert Pressure levels in Pa')
[1632]	168	! Fin raccord
	169	CALL getpar('ini_anal',.false.,ini_anal,'Etat initial = analyse')
	170	CALL getpar('guide_invertp',.true.,invert_p,'niveaux p inverses')
	171	CALL getpar('guide_inverty',.true.,invert_y,'inversion N-S')
	172	CALL getpar('guide_2D',.false.,guide_2D,'fichier guidage lat-P')
	173
[2263]	174	call fin_getparam
[5249]	175
[1632]	176	! ---------------------------------------------
	177	! Determination du nombre de niveaux verticaux
	178	! des fichiers guidage
	179	! ---------------------------------------------
	180	ncidpl=-99
[5084]	181	if (guide_plevs.EQ.1) then
	182	if (ncidpl.eq.-99) then
[1902]	183	rcod=nf90_open('apbp.nc',Nf90_NOWRITe, ncidpl)
[5270]	184	if (rcod.NE.nf90_noerr) THEN
[3995]	185	abort_message=' Nudging error -> no file apbp.nc'
[1902]	186	CALL abort_gcm(modname,abort_message,1)
	187	endif
	188	endif
[5084]	189	elseif (guide_plevs.EQ.2) then
	190	if (ncidpl.EQ.-99) then
[1902]	191	rcod=nf90_open('P.nc',Nf90_NOWRITe,ncidpl)
[5270]	192	if (rcod.NE.nf90_noerr) THEN
[3995]	193	abort_message=' Nudging error -> no file P.nc'
[1902]	194	CALL abort_gcm(modname,abort_message,1)
	195	endif
	196	endif
[3984]	197
[1632]	198	elseif (guide_u) then
[5084]	199	if (ncidpl.eq.-99) then
[1902]	200	rcod=nf90_open('u.nc',Nf90_NOWRITe,ncidpl)
[5270]	201	if (rcod.NE.nf90_noerr) THEN
[3995]	202	abort_message=' Nudging error -> no file u.nc'
[1902]	203	CALL abort_gcm(modname,abort_message,1)
	204	endif
[5249]	205
[1902]	206	endif
[3984]	207
	208
[1632]	209	elseif (guide_v) then
[5084]	210	if (ncidpl.eq.-99) then
[1902]	211	rcod=nf90_open('v.nc',nf90_nowrite,ncidpl)
[5270]	212	if (rcod.NE.nf90_noerr) THEN
[3995]	213	abort_message=' Nudging error -> no file v.nc'
[1902]	214	CALL abort_gcm(modname,abort_message,1)
	215	endif
	216	endif
[3984]	217
[5249]	218
[1632]	219	elseif (guide_T) then
[5084]	220	if (ncidpl.eq.-99) then
[1902]	221	rcod=nf90_open('T.nc',nf90_nowrite,ncidpl)
[5270]	222	if (rcod.NE.nf90_noerr) THEN
[3995]	223	abort_message=' Nudging error -> no file T.nc'
[1902]	224	CALL abort_gcm(modname,abort_message,1)
	225	endif
	226	endif
[3984]	227
	228
	229
[1632]	230	elseif (guide_Q) then
[5084]	231	if (ncidpl.eq.-99) then
[1902]	232	rcod=nf90_open('hur.nc',nf90_nowrite, ncidpl)
[5270]	233	if (rcod.NE.nf90_noerr) THEN
[3995]	234	abort_message=' Nudging error -> no file hur.nc'
[1902]	235	CALL abort_gcm(modname,abort_message,1)
	236	endif
	237	endif
[3984]	238
	239
[5249]	240	endif
[5270]	241	error=nf90_inq_dimid(ncidpl,'LEVEL',rid)
	242	IF (error.NE.nf90_noerr) error=nf90_inq_dimid(ncidpl,'PRESSURE',rid)
	243	IF (error.NE.nf90_noerr) THEN
[3995]	244	abort_message='Nudging: error reading pressure levels'
[1632]	245	CALL abort_gcm(modname,abort_message,1)
	246	ENDIF
[5270]	247	error=nf90_inquire_dimension(ncidpl,rid,len=nlevnc)
[3995]	248	write(,)trim(modname)//' : number of vertical levels nlevnc', nlevnc
[1632]	249	rcod = nf90_close(ncidpl)
	250
	251	! ---------------------------------------------
	252	! Allocation des variables
	253	! ---------------------------------------------
[3995]	254	abort_message='nudging allocation error'
[1632]	255
	256	ALLOCATE(apnc(nlevnc), stat = error)
	257	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	258	ALLOCATE(bpnc(nlevnc), stat = error)
	259	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	260	apnc=0.;bpnc=0.
	261
	262	ALLOCATE(alpha_pcor(llm), stat = error)
	263	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	264	ALLOCATE(alpha_u(ijb_u:ije_u), stat = error)
	265	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	266	ALLOCATE(alpha_v(ijb_v:ije_v), stat = error)
	267	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	268	ALLOCATE(alpha_T(ijb_u:ije_u), stat = error)
	269	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	270	ALLOCATE(alpha_Q(ijb_u:ije_u), stat = error)
	271	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	272	ALLOCATE(alpha_P(ijb_u:ije_u), stat = error)
	273	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	274	alpha_u=0.;alpha_v=0;alpha_T=0;alpha_Q=0;alpha_P=0
[5249]	275
[1632]	276	IF (guide_u) THEN
	277	ALLOCATE(unat1(iip1,jjb_u:jje_u,nlevnc), stat = error)
	278	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	279	ALLOCATE(ugui1(ijb_u:ije_u,llm), stat = error)
	280	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	281	ALLOCATE(unat2(iip1,jjb_u:jje_u,nlevnc), stat = error)
	282	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	283	ALLOCATE(ugui2(ijb_u:ije_u,llm), stat = error)
	284	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	285	unat1=0.;unat2=0.;ugui1=0.;ugui2=0.
	286	ENDIF
	287
	288	IF (guide_T) THEN
	289	ALLOCATE(tnat1(iip1,jjb_u:jje_u,nlevnc), stat = error)
	290	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	291	ALLOCATE(tgui1(ijb_u:ije_u,llm), stat = error)
	292	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	293	ALLOCATE(tnat2(iip1,jjb_u:jje_u,nlevnc), stat = error)
	294	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	295	ALLOCATE(tgui2(ijb_u:ije_u,llm), stat = error)
	296	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	297	tnat1=0.;tnat2=0.;tgui1=0.;tgui2=0.
	298	ENDIF
[5249]	299
[1632]	300	IF (guide_Q) THEN
	301	ALLOCATE(qnat1(iip1,jjb_u:jje_u,nlevnc), stat = error)
	302	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	303	ALLOCATE(qgui1(ijb_u:ije_u,llm), stat = error)
	304	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	305	ALLOCATE(qnat2(iip1,jjb_u:jje_u,nlevnc), stat = error)
	306	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	307	ALLOCATE(qgui2(ijb_u:ije_u,llm), stat = error)
	308	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	309	qnat1=0.;qnat2=0.;qgui1=0.;qgui2=0.
	310	ENDIF
	311
	312	IF (guide_v) THEN
	313	ALLOCATE(vnat1(iip1,jjb_v:jje_v,nlevnc), stat = error)
	314	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	315	ALLOCATE(vgui1(ijb_v:ije_v,llm), stat = error)
	316	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	317	ALLOCATE(vnat2(iip1,jjb_v:jje_v,nlevnc), stat = error)
	318	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	319	ALLOCATE(vgui2(ijb_v:ije_v,llm), stat = error)
	320	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	321	vnat1=0.;vnat2=0.;vgui1=0.;vgui2=0.
	322	ENDIF
	323
[5084]	324	IF (guide_plevs.EQ.2) THEN
[1632]	325	ALLOCATE(pnat1(iip1,jjb_u:jje_u,nlevnc), stat = error)
	326	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	327	ALLOCATE(pnat2(iip1,jjb_u:jje_u,nlevnc), stat = error)
	328	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	329	pnat1=0.;pnat2=0.;
	330	ENDIF
	331
[5084]	332	IF (guide_P.OR.guide_plevs.EQ.1) THEN
[1632]	333	ALLOCATE(psnat1(iip1,jjb_u:jje_u), stat = error)
	334	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	335	ALLOCATE(psnat2(iip1,jjb_u:jje_u), stat = error)
	336	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	337	psnat1=0.;psnat2=0.;
	338	ENDIF
	339	IF (guide_P) THEN
	340	ALLOCATE(psgui2(ijb_u:ije_u), stat = error)
	341	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	342	ALLOCATE(psgui1(ijb_u:ije_u), stat = error)
	343	IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
	344	psgui1=0.;psgui2=0.
	345	ENDIF
	346
	347	! ---------------------------------------------
	348	! Lecture du premier etat de guidage.
	349	! ---------------------------------------------
	350	IF (guide_2D) THEN
	351	CALL guide_read2D(1)
	352	ELSE
	353	CALL guide_read(1)
	354	ENDIF
	355	IF (guide_v) vnat1=vnat2
	356	IF (guide_u) unat1=unat2
	357	IF (guide_T) tnat1=tnat2
	358	IF (guide_Q) qnat1=qnat2
[5084]	359	IF (guide_plevs.EQ.2) pnat1=pnat2
	360	IF (guide_P.OR.guide_plevs.EQ.1) psnat1=psnat2
[1632]	361
	362	END SUBROUTINE guide_init
	363
	364	!=======================================================================
	365	SUBROUTINE guide_main(itau,ucov,vcov,teta,q,masse,ps)
[2021]	366	use exner_hyb_loc_m, only: exner_hyb_loc
	367	use exner_milieu_loc_m, only: exner_milieu_loc
[1823]	368	USE parallel_lmdz
[1632]	369	USE control_mod
[1806]	370	USE write_field_loc
[2597]	371	USE comconst_mod, ONLY: cpp, daysec, dtvr, kappa
[2600]	372	USE comvert_mod, ONLY: ap, bp, preff, presnivs, pressure_exner
[5249]	373
[5271]	374	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
[5272]	375	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	376	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
[5271]	377	IMPLICIT NONE
[5249]	378
[5271]	379
[1632]	380
[5272]	381
[1632]	382	! Variables entree
	383	INTEGER, INTENT(IN) :: itau !pas de temps
	384	REAL, DIMENSION (ijb_u:ije_u,llm), INTENT(INOUT) :: ucov,teta,q,masse
	385	REAL, DIMENSION (ijb_v:ije_v,llm), INTENT(INOUT) :: vcov
	386	REAL, DIMENSION (ijb_u:ije_u), INTENT(INOUT) :: ps
	387
	388	! Variables locales
	389	LOGICAL, SAVE :: first=.TRUE.
	390	!$OMP THREADPRIVATE(first)
	391	LOGICAL :: f_out ! sortie guidage
[1806]	392	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: f_addu ! var aux: champ de guidage
	393	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: f_addv ! var aux: champ de guidage
[1632]	394	! Variables pour fonction Exner (P milieu couche)
[2021]	395	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:,:) :: pk
[5249]	396	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: pks
[1632]	397	REAL :: unskap
[1806]	398	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: p ! besoin si guide_P
[1632]	399	! Compteurs temps:
	400	INTEGER, SAVE :: step_rea,count_no_rea,itau_test ! lecture guidage
	401	!$OMP THREADPRIVATE(step_rea,count_no_rea,itau_test)
	402	REAL :: ditau, dday_step
	403	REAL :: tau,reste ! position entre 2 etats de guidage
	404	REAL, SAVE :: factt ! pas de temps en fraction de jour
	405	!$OMP THREADPRIVATE(factt)
[5249]	406
[1632]	407	INTEGER :: i,j,l
[3995]	408	CHARACTER(LEN=20) :: modname="guide_main"
[5249]	409
	410	!$OMP MASTER
[3995]	411	ijbu=ij_begin ; ijeu=ij_end
[5249]	412	jjbu=jj_begin ; jjeu=jj_end ; jjnu=jjeu-jjbu+1
[3995]	413	ijbv=ij_begin ; ijev=ij_end
[5249]	414	jjbv=jj_begin ; jjev=jj_end ; jjnv=jjev-jjbv+1
[1632]	415	IF (pole_sud) THEN
[3995]	416	ijeu=ij_end-iip1
[1632]	417	ijev=ij_end-iip1
	418	jjev=jj_end-1
[5249]	419	jjnv=jjev-jjbv+1
[1632]	420	ENDIF
[3995]	421	IF (pole_nord) THEN
	422	ijbu=ij_begin+iip1
	423	ijbv=ij_begin
	424	ENDIF
[1632]	425	!$OMP END MASTER
	426	!$OMP BARRIER
[5249]	427
[2036]	428	! PRINT *,'---> on rentre dans guide_main'
[1632]	429	! CALL AllGather_Field(ucov,ip1jmp1,llm)
	430	! CALL AllGather_Field(vcov,ip1jm,llm)
	431	! CALL AllGather_Field(teta,ip1jmp1,llm)
	432	! CALL AllGather_Field(ps,ip1jmp1,1)
	433	! CALL AllGather_Field(q,ip1jmp1,llm)
[5249]	434
[1632]	435	!-----------------------------------------------------------------------
	436	! Initialisations au premier passage
	437	!-----------------------------------------------------------------------
	438
	439	IF (first) THEN
	440	first=.FALSE.
	441	!$OMP MASTER
[5249]	442	ALLOCATE(f_addu(ijb_u:ije_u,llm) )
	443	ALLOCATE(f_addv(ijb_v:ije_v,llm) )
	444	ALLOCATE(pk(iip1,jjb_u:jje_u,llm) )
	445	ALLOCATE(pks(iip1,jjb_u:jje_u) )
	446	ALLOCATE(p(ijb_u:ije_u,llmp1) )
	447	CALL guide_init
[1632]	448	!$OMP END MASTER
	449	!$OMP BARRIER
	450	itau_test=1001
	451	step_rea=1
	452	count_no_rea=0
	453	! Calcul des constantes de rappel
[5249]	454	factt=dtvr*iperiod/daysec
[1632]	455	!$OMP MASTER
[1806]	456	call tau2alpha(3, iip1, jjb_v, jje_v, factt, tau_min_v, tau_max_v, alpha_v)
	457	call tau2alpha(2, iip1, jjb_u, jje_u, factt, tau_min_u, tau_max_u, alpha_u)
	458	call tau2alpha(1, iip1, jjb_u, jje_u, factt, tau_min_T, tau_max_T, alpha_T)
	459	call tau2alpha(1, iip1, jjb_u, jje_u, factt, tau_min_P, tau_max_P, alpha_P)
	460	call tau2alpha(1, iip1, jjb_u, jje_u, factt, tau_min_Q, tau_max_Q, alpha_Q)
[1632]	461	! correction de rappel dans couche limite
	462	if (guide_BL) then
	463	alpha_pcor(:)=1.
	464	else
	465	do l=1,llm
[4246]	466	alpha_pcor(l)=(1.+tanh(((plim_guide_BL-presnivs(l))/preff)/0.05))/2.
[1632]	467	enddo
	468	endif
	469	!$OMP END MASTER
[1806]	470	!$OMP BARRIER
[5249]	471	! ini_anal: etat initial egal au guidage
[1632]	472	IF (ini_anal) THEN
	473	CALL guide_interp(ps,teta)
[5249]	474	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	475	DO l=1,llm
	476	IF (guide_u) ucov(ijbu:ijeu,l)=ugui2(ijbu:ijeu,l)
	477	IF (guide_v) vcov(ijbv:ijev,l)=ugui2(ijbv:ijev,l)
	478	IF (guide_T) teta(ijbu:ijeu,l)=tgui2(ijbu:ijeu,l)
	479	IF (guide_Q) q(ijbu:ijeu,l)=qgui2(ijbu:ijeu,l)
	480	ENDDO
[5249]	481
[1632]	482	IF (guide_P) THEN
[1806]	483	!$OMP MASTER
[1632]	484	ps(ijbu:ijeu)=psgui2(ijbu:ijeu)
[1806]	485	!$OMP END MASTER
	486	!$OMP BARRIER
[1632]	487	CALL pression_loc(ijnb_u,ap,bp,ps,p)
	488	CALL massdair_loc(p,masse)
[1806]	489	!$OMP BARRIER
[1632]	490	ENDIF
	491	RETURN
	492	ENDIF
	493
	494	ENDIF !first
	495
	496	!-----------------------------------------------------------------------
	497	! Lecture des fichiers de guidage ?
	498	!-----------------------------------------------------------------------
[5084]	499	IF (iguide_read.NE.0) THEN
[1632]	500	ditau=real(itau)
	501	dday_step=real(day_step)
[5084]	502	IF (iguide_read.LT.0) THEN
[1632]	503	tau=ditau/dday_step/REAL(iguide_read)
	504	ELSE
	505	tau=REAL(iguide_read)*ditau/dday_step
	506	ENDIF
	507	reste=tau-AINT(tau)
[5084]	508	IF (reste.EQ.0.) THEN
	509	IF (itau_test.EQ.itau) THEN
[3995]	510	write(,)trim(modname)//' second pass in advreel at itau=',&
	511	itau
[4469]	512	CALL abort_gcm("guide_loc_lod","stopped",1)
[1632]	513	ELSE
[1806]	514	!$OMP MASTER
[1632]	515	IF (guide_v) vnat1(:,jjbv:jjev,:)=vnat2(:,jjbv:jjev,:)
	516	IF (guide_u) unat1(:,jjbu:jjeu,:)=unat2(:,jjbu:jjeu,:)
	517	IF (guide_T) tnat1(:,jjbu:jjeu,:)=tnat2(:,jjbu:jjeu,:)
	518	IF (guide_Q) qnat1(:,jjbu:jjeu,:)=qnat2(:,jjbu:jjeu,:)
[5084]	519	IF (guide_plevs.EQ.2) pnat1(:,jjbu:jjeu,:)=pnat2(:,jjbu:jjeu,:)
	520	IF (guide_P.OR.guide_plevs.EQ.1) psnat1(:,jjbu:jjeu)=psnat2(:,jjbu:jjeu)
[1806]	521	!$OMP END MASTER
	522	!$OMP BARRIER
[1632]	523	step_rea=step_rea+1
	524	itau_test=itau
[3995]	525	if (is_master) then
	526	write(,)trim(modname)//' Reading nudging files, step ',&
	527	step_rea,'after ',count_no_rea,' skips'
	528	endif
[1632]	529	IF (guide_2D) THEN
[1806]	530	!$OMP MASTER
[1632]	531	CALL guide_read2D(step_rea)
[1806]	532	!$OMP END MASTER
	533	!$OMP BARRIER
[1632]	534	ELSE
[1806]	535	!$OMP MASTER
[1632]	536	CALL guide_read(step_rea)
[1806]	537	!$OMP END MASTER
	538	!$OMP BARRIER
[1632]	539	ENDIF
	540	count_no_rea=0
	541	ENDIF
	542	ELSE
	543	count_no_rea=count_no_rea+1
	544
	545	ENDIF
	546	ENDIF !iguide_read=0
	547
	548	!-----------------------------------------------------------------------
	549	! Interpolation et conversion des champs de guidage
	550	!-----------------------------------------------------------------------
[5084]	551	IF (MOD(itau,iguide_int).EQ.0) THEN
[1632]	552	CALL guide_interp(ps,teta)
	553	ENDIF
	554	! Repartition entre 2 etats de guidage
[5084]	555	IF (iguide_read.NE.0) THEN
[1632]	556	tau=reste
	557	ELSE
	558	tau=1.
	559	ENDIF
	560
[1806]	561	! CALL WriteField_u('ucov_guide',ucov)
	562	! CALL WriteField_v('vcov_guide',vcov)
	563	! CALL WriteField_u('teta_guide',teta)
	564	! CALL WriteField_u('masse_guide',masse)
[5249]	565
	566
[3995]	567	!-----------------------------------------------------------------------
[5249]	568	! Ajout des champs de guidage
[1632]	569	!-----------------------------------------------------------------------
	570	! Sauvegarde du guidage?
[5249]	571	f_out=((MOD(itau,iguide_sav).EQ.0).AND.guide_sav)
[1632]	572	IF (f_out) THEN
[2036]	573
	574	!$OMP BARRIER
	575	CALL pression_loc(ijnb_u,ap,bp,ps,p)
	576
	577	!$OMP BARRIER
	578	if (pressure_exner) then
	579	CALL exner_hyb_loc( ijnb_u, ps, p, pks, pk)
	580	else
	581	CALL exner_milieu_loc( ijnb_u, ps, p, pks, pk )
	582	endif
	583
	584	!$OMP BARRIER
	585
[1632]	586	unskap=1./kappa
[2036]	587	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[2028]	588	DO l = 1, llm
	589	DO j=jjbu,jjeu
	590	DO i =1, iip1
	591	p(i+(j-1)iip1,l) = preff ( pk(i,j,l)/cpp) ** unskap
	592	ENDDO
	593	ENDDO
	594	ENDDO
[2036]	595
	596	CALL guide_out("SP",jjp1,llm,p(ijb_u:ije_u,1:llm),1.)
[1632]	597	ENDIF
[5249]	598
[1632]	599	if (guide_u) then
	600	if (guide_add) then
[2036]	601	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	602	DO l=1,llm
	603	f_addu(ijbu:ijeu,l)=(1.-tau)ugui1(ijbu:ijeu,l)+tauugui2(ijbu:ijeu,l)
	604	ENDDO
[1632]	605	else
[2036]	606	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	607	DO l=1,llm
	608	f_addu(ijbu:ijeu,l)=(1.-tau)ugui1(ijbu:ijeu,l)+tauugui2(ijbu:ijeu,l)-ucov(ijbu:ijeu,l)
	609	ENDDO
[5249]	610	endif
	611
[1806]	612	! CALL WriteField_u('f_addu',f_addu)
[1632]	613
[1806]	614	if (guide_zon) CALL guide_zonave_u(1,llm,f_addu)
	615	CALL guide_addfield_u(llm,f_addu,alpha_u)
[2028]	616	IF (f_out) CALL guide_out("ua",jjp1,llm,(1.-tau)ugui1(ijb_u:ije_u,:)+tauugui2(ijb_u:ije_u,:),factt)
	617	IF (f_out) CALL guide_out("u",jjp1,llm,ucov(ijb_u:ije_u,:),factt)
[3995]	618	IF (f_out) THEN
	619	! Ehouarn: fill the gaps adequately...
	620	IF (ijbu>ijb_u) f_addu(ijb_u:ijbu-1,:)=0
	621	IF (ijeu<ije_u) f_addu(ijeu+1:ije_u,:)=0
	622	CALL guide_out("ucov",jjp1,llm,f_addu(ijb_u:ije_u,:)/factt,factt)
	623	ENDIF
[2036]	624	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	625	DO l=1,llm
	626	ucov(ijbu:ijeu,l)=ucov(ijbu:ijeu,l)+f_addu(ijbu:ijeu,l)
	627	ENDDO
	628
[1632]	629	endif
	630
	631	if (guide_T) then
	632	if (guide_add) then
[2036]	633	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	634	DO l=1,llm
	635	f_addu(ijbu:ijeu,l)=(1.-tau)tgui1(ijbu:ijeu,l)+tautgui2(ijbu:ijeu,l)
	636	ENDDO
[1632]	637	else
[2036]	638	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	639	DO l=1,llm
	640	f_addu(ijbu:ijeu,l)=(1.-tau)tgui1(ijbu:ijeu,l)+tautgui2(ijbu:ijeu,l)-teta(ijbu:ijeu,l)
	641	ENDDO
[5249]	642	endif
[1806]	643	if (guide_zon) CALL guide_zonave_u(2,llm,f_addu)
	644	CALL guide_addfield_u(llm,f_addu,alpha_T)
[2028]	645	IF (f_out) CALL guide_out("teta",jjp1,llm,f_addu(:,:)/factt,factt)
[2036]	646	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	647	DO l=1,llm
	648	teta(ijbu:ijeu,l)=teta(ijbu:ijeu,l)+f_addu(ijbu:ijeu,l)
	649	ENDDO
[1632]	650	endif
	651
	652	if (guide_P) then
	653	if (guide_add) then
[1806]	654	!$OMP MASTER
	655	f_addu(ijbu:ijeu,1)=(1.-tau)psgui1(ijbu:ijeu)+taupsgui2(ijbu:ijeu)
	656	!$OMP END MASTER
	657	!$OMP BARRIER
[1632]	658	else
[1806]	659	!$OMP MASTER
	660	f_addu(ijbu:ijeu,1)=(1.-tau)psgui1(ijbu:ijeu)+taupsgui2(ijbu:ijeu)-ps(ijbu:ijeu)
	661	!$OMP END MASTER
	662	!$OMP BARRIER
[5249]	663	endif
[1806]	664	if (guide_zon) CALL guide_zonave_u(2,1,f_addu(ijb_u:ije_u,1))
	665	CALL guide_addfield_u(1,f_addu(ijb_u:ije_u,1),alpha_P)
[2036]	666	! IF (f_out) CALL guide_out("ps",jjp1,1,f_addu(ijb_u:ije_u,1)/factt,factt)
[1806]	667	!$OMP MASTER
	668	ps(ijbu:ijeu)=ps(ijbu:ijeu)+f_addu(ijbu:ijeu,1)
	669	!$OMP END MASTER
	670	!$OMP BARRIER
[1632]	671	CALL pression_loc(ijnb_u,ap,bp,ps,p)
	672	CALL massdair_loc(p,masse)
[1806]	673	!$OMP BARRIER
[1632]	674	endif
	675
	676	if (guide_Q) then
	677	if (guide_add) then
[2036]	678	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	679	DO l=1,llm
	680	f_addu(ijbu:ijeu,l)=(1.-tau)qgui1(ijbu:ijeu,l)+tauqgui2(ijbu:ijeu,l)
	681	ENDDO
[1632]	682	else
[2036]	683	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	684	DO l=1,llm
	685	f_addu(ijbu:ijeu,l)=(1.-tau)qgui1(ijbu:ijeu,l)+tauqgui2(ijbu:ijeu,l)-q(ijbu:ijeu,l)
	686	ENDDO
[5249]	687	endif
[1806]	688	if (guide_zon) CALL guide_zonave_u(2,llm,f_addu)
	689	CALL guide_addfield_u(llm,f_addu,alpha_Q)
[2028]	690	IF (f_out) CALL guide_out("q",jjp1,llm,f_addu(:,:)/factt,factt)
[1806]	691
[2036]	692	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	693	DO l=1,llm
	694	q(ijbu:ijeu,l)=q(ijbu:ijeu,l)+f_addu(ijbu:ijeu,l)
	695	ENDDO
[1632]	696	endif
	697
	698	if (guide_v) then
	699	if (guide_add) then
[2036]	700	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	701	DO l=1,llm
	702	f_addv(ijbv:ijev,l)=(1.-tau)vgui1(ijbv:ijev,l)+tauvgui2(ijbv:ijev,l)
	703	ENDDO
	704
[1632]	705	else
[2036]	706	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	707	DO l=1,llm
	708	f_addv(ijbv:ijev,l)=(1.-tau)vgui1(ijbv:ijev,l)+tauvgui2(ijbv:ijev,l)-vcov(ijbv:ijev,l)
	709	ENDDO
	710
[5249]	711	endif
	712
[1806]	713	if (guide_zon) CALL guide_zonave_v(2,jjm,llm,f_addv(ijb_v:ije_v,:))
[5249]	714
[1806]	715	CALL guide_addfield_v(llm,f_addv(ijb_v:ije_v,:),alpha_v)
[2028]	716	IF (f_out) CALL guide_out("v",jjm,llm,vcov(ijb_v:ije_v,:),factt)
	717	IF (f_out) CALL guide_out("va",jjm,llm,(1.-tau)vgui1(ijb_v:ije_v,:)+tauvgui2(ijb_v:ije_v,:),factt)
[3995]	718	IF (f_out) THEN
	719	! Ehouarn: Fill in the gaps adequately
[5249]	720	IF (ijbv>ijb_v) f_addv(ijb_v:ijbv-1,:)=0
[3995]	721	IF (ijev<ije_v) f_addv(ijev+1:ije_v,:)=0
	722	CALL guide_out("vcov",jjm,llm,f_addv(ijb_v:ije_v,:)/factt,factt)
	723	ENDIF
[1806]	724
[2036]	725	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	726	DO l=1,llm
	727	vcov(ijbv:ijev,l)=vcov(ijbv:ijev,l)+f_addv(ijbv:ijev,l)
	728	ENDDO
[1632]	729	endif
	730
	731	END SUBROUTINE guide_main
	732
	733
	734	SUBROUTINE guide_addfield_u(vsize,field,alpha)
	735	! field1=afield1+alphafield2
	736
[5271]	737	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
[5272]	738	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	739	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
[5271]	740	IMPLICIT NONE
	741
[1632]	742
[5272]	743
[1632]	744	! input variables
	745	INTEGER, INTENT(IN) :: vsize
[5249]	746	REAL, DIMENSION(ijb_u:ije_u), INTENT(IN) :: alpha
[1632]	747	REAL, DIMENSION(ijb_u:ije_u,vsize), INTENT(INOUT) :: field
	748
	749	! Local variables
	750	INTEGER :: l
	751
[2036]	752	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1632]	753	DO l=1,vsize
	754	field(ijbu:ijeu,l)=alpha(ijbu:ijeu)field(ijbu:ijeu,l)alpha_pcor(l)
	755	ENDDO
	756
	757	END SUBROUTINE guide_addfield_u
	758
	759
	760	SUBROUTINE guide_addfield_v(vsize,field,alpha)
	761	! field1=afield1+alphafield2
	762
[5271]	763	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
[5272]	764	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	765	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
[5271]	766	IMPLICIT NONE
	767
[1632]	768
[5272]	769
[1632]	770	! input variables
	771	INTEGER, INTENT(IN) :: vsize
[5249]	772	REAL, DIMENSION(ijb_v:ije_v), INTENT(IN) :: alpha
[1632]	773	REAL, DIMENSION(ijb_v:ije_v,vsize), INTENT(INOUT) :: field
	774
	775	! Local variables
	776	INTEGER :: l
	777
[2036]	778	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1632]	779	DO l=1,vsize
	780	field(ijbv:ijev,l)=alpha(ijbv:ijev)field(ijbv:ijev,l)alpha_pcor(l)
	781	ENDDO
	782
	783	END SUBROUTINE guide_addfield_v
[5249]	784
[1632]	785	!=======================================================================
	786
	787	SUBROUTINE guide_zonave_u(typ,vsize,field)
	788
[2597]	789	USE comconst_mod, ONLY: pi
[5249]	790
[5271]	791	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
[5272]	792	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	793	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
[5271]	794	IMPLICIT NONE
[1632]	795
[5271]	796
[5272]	797
[1632]	798	INCLUDE "comgeom.h"
[5249]	799
[1632]	800	! input/output variables
	801	INTEGER, INTENT(IN) :: typ
	802	INTEGER, INTENT(IN) :: vsize
	803	REAL, DIMENSION(ijb_u:ije_u,vsize), INTENT(INOUT) :: field
	804
	805	! Local variables
	806	LOGICAL, SAVE :: first=.TRUE.
[1806]	807	!$OMP THREADPRIVATE(first)
	808
[1632]	809	INTEGER, DIMENSION (2), SAVE :: imin, imax ! averaging domain
[5249]	810	!$OMP THREADPRIVATE(imin,imax)
[1632]	811	INTEGER :: i,j,l,ij
	812	REAL, DIMENSION (iip1) :: lond ! longitude in Deg.
	813	REAL, DIMENSION (jjb_u:jje_u,vsize):: fieldm ! zon-averaged field
	814
	815	IF (first) THEN
	816	first=.FALSE.
	817	!Compute domain for averaging
	818	lond=rlonu*180./pi
	819	imin(1)=1;imax(1)=iip1;
	820	imin(2)=1;imax(2)=iip1;
	821	IF (guide_reg) THEN
	822	DO i=1,iim
[5084]	823	IF (lond(i).LT.lon_min_g) imin(1)=i
	824	IF (lond(i).LE.lon_max_g) imax(1)=i
[1632]	825	ENDDO
	826	lond=rlonv*180./pi
	827	DO i=1,iim
[5084]	828	IF (lond(i).LT.lon_min_g) imin(2)=i
	829	IF (lond(i).LE.lon_max_g) imax(2)=i
[1632]	830	ENDDO
	831	ENDIF
	832	ENDIF
	833
[5249]	834
[2036]	835	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1632]	836	DO l=1,vsize
[1806]	837	fieldm(:,l)=0.
[1632]	838	! Compute zonal average
	839
	840	!correction bug ici
	841	! ---> a verifier
	842	! ym DO j=jjbv,jjev
	843	DO j=jjbu,jjeu
	844	DO i=imin(typ),imax(typ)
	845	ij=(j-1)*iip1+i
	846	fieldm(j,l)=fieldm(j,l)+field(ij,l)
	847	ENDDO
[5249]	848	ENDDO
[1632]	849	fieldm(:,l)=fieldm(:,l)/REAL(imax(typ)-imin(typ)+1)
	850	! Compute forcing
	851	DO j=jjbu,jjeu
	852	DO i=1,iip1
	853	ij=(j-1)*iip1+i
	854	field(ij,l)=fieldm(j,l)
	855	ENDDO
	856	ENDDO
	857	ENDDO
	858
	859	END SUBROUTINE guide_zonave_u
	860
	861
	862	SUBROUTINE guide_zonave_v(typ,hsize,vsize,field)
	863
[2597]	864	USE comconst_mod, ONLY: pi
[5249]	865
[5271]	866	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
[5272]	867	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	868	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
[5271]	869	IMPLICIT NONE
[1632]	870
[5271]	871
[5272]	872
[1632]	873	INCLUDE "comgeom.h"
[5249]	874
[1632]	875	! input/output variables
	876	INTEGER, INTENT(IN) :: typ
	877	INTEGER, INTENT(IN) :: vsize
	878	INTEGER, INTENT(IN) :: hsize
	879	REAL, DIMENSION(ijb_v:ije_v,vsize), INTENT(INOUT) :: field
	880
	881	! Local variables
	882	LOGICAL, SAVE :: first=.TRUE.
[1806]	883	!$OMP THREADPRIVATE(first)
[1632]	884	INTEGER, DIMENSION (2), SAVE :: imin, imax ! averaging domain
[1806]	885	!$OMP THREADPRIVATE(imin, imax)
[1632]	886	INTEGER :: i,j,l,ij
	887	REAL, DIMENSION (iip1) :: lond ! longitude in Deg.
	888	REAL, DIMENSION (jjb_v:jjev,vsize):: fieldm ! zon-averaged field
	889
	890	IF (first) THEN
	891	first=.FALSE.
	892	!Compute domain for averaging
	893	lond=rlonu*180./pi
	894	imin(1)=1;imax(1)=iip1;
	895	imin(2)=1;imax(2)=iip1;
	896	IF (guide_reg) THEN
	897	DO i=1,iim
[5084]	898	IF (lond(i).LT.lon_min_g) imin(1)=i
	899	IF (lond(i).LE.lon_max_g) imax(1)=i
[1632]	900	ENDDO
	901	lond=rlonv*180./pi
	902	DO i=1,iim
[5084]	903	IF (lond(i).LT.lon_min_g) imin(2)=i
	904	IF (lond(i).LE.lon_max_g) imax(2)=i
[1632]	905	ENDDO
	906	ENDIF
	907	ENDIF
	908
[2036]	909	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1632]	910	DO l=1,vsize
	911	! Compute zonal average
[1806]	912	fieldm(:,l)=0.
[1632]	913	DO j=jjbv,jjev
	914	DO i=imin(typ),imax(typ)
	915	ij=(j-1)*iip1+i
	916	fieldm(j,l)=fieldm(j,l)+field(ij,l)
	917	ENDDO
[5249]	918	ENDDO
[1632]	919	fieldm(:,l)=fieldm(:,l)/REAL(imax(typ)-imin(typ)+1)
	920	! Compute forcing
	921	DO j=jjbv,jjev
	922	DO i=1,iip1
	923	ij=(j-1)*iip1+i
	924	field(ij,l)=fieldm(j,l)
	925	ENDDO
	926	ENDDO
	927	ENDDO
	928
	929
	930	END SUBROUTINE guide_zonave_v
[5249]	931
[1632]	932	!=======================================================================
	933	SUBROUTINE guide_interp(psi,teta)
[2021]	934	use exner_hyb_loc_m, only: exner_hyb_loc
	935	use exner_milieu_loc_m, only: exner_milieu_loc
[1823]	936	USE parallel_lmdz
[1632]	937	USE mod_hallo
	938	USE Bands
[2597]	939	USE comconst_mod, ONLY: cpp, kappa
[2600]	940	USE comvert_mod, ONLY: preff, pressure_exner, bp, ap, disvert_type
[5271]	941	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
[5272]	942	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	943	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
[5271]	944	IMPLICIT NONE
[1632]	945
[5271]	946
[5272]	947
[1632]	948	include "comgeom2.h"
	949
	950	REAL, DIMENSION (iip1,jjb_u:jje_u), INTENT(IN) :: psi ! Psol gcm
	951	REAL, DIMENSION (iip1,jjb_u:jje_u,llm), INTENT(IN) :: teta ! Temp. Pot. gcm
	952
	953	LOGICAL, SAVE :: first=.TRUE.
[1806]	954	!$OMP THREADPRIVATE(first)
[1632]	955	! Variables pour niveaux pression:
[1806]	956	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:,:) :: plnc1,plnc2 !niveaux pression guidage
	957	REAL, ALLOCATABLE, SAVE,DIMENSION (:,:,:) :: plunc,plsnc !niveaux pression modele
	958	REAL, ALLOCATABLE, SAVE,DIMENSION (:,:,:) :: plvnc !niveaux pression modele
[5249]	959	REAL, ALLOCATABLE, SAVE,DIMENSION (:,:,:) :: p ! pression intercouches
[1806]	960	REAL, ALLOCATABLE, SAVE,DIMENSION (:,:,:) :: pls, pext ! var intermediaire
[5249]	961	REAL, ALLOCATABLE, SAVE,DIMENSION (:,:,:) :: pbarx
	962	REAL, ALLOCATABLE, SAVE,DIMENSION (:,:,:) :: pbary
[1632]	963	! Variables pour fonction Exner (P milieu couche)
[2021]	964	REAL, ALLOCATABLE, SAVE,DIMENSION (:,:,:) :: pk
[5249]	965	REAL ,ALLOCATABLE, SAVE, DIMENSION (:,:) :: pks
[1632]	966	REAL :: unskap
	967	! Pression de vapeur saturante
[1806]	968	REAL, ALLOCATABLE, SAVE,DIMENSION (:,:) :: qsat
[1632]	969	!Variables intermediaires interpolation
[5249]	970	REAL, ALLOCATABLE, SAVE,DIMENSION (:,:,:) :: zu1,zu2
[1806]	971	REAL, ALLOCATABLE, SAVE,DIMENSION (:,:,:) :: zv1,zv2
[5249]	972
[1632]	973	INTEGER :: i,j,l,ij
[3995]	974	CHARACTER(LEN=20),PARAMETER :: modname="guide_interp"
[5249]	975	TYPE(Request),SAVE :: Req
[1848]	976	!$OMP THREADPRIVATE(Req)
[5249]	977
[3995]	978	if (is_master) write(,)trim(modname)//': interpolate nudging variables'
[1632]	979	! -----------------------------------------------------------------
	980	! Calcul des niveaux de pression champs guidage (pour T et Q)
	981	! -----------------------------------------------------------------
[1806]	982	IF (first) THEN
	983	!$OMP MASTER
[5249]	984	ALLOCATE(plnc1(iip1,jjb_u:jje_u,nlevnc) )
	985	ALLOCATE(plnc2(iip1,jjb_u:jje_u,nlevnc) )
	986	ALLOCATE(plunc(iip1,jjb_u:jje_u,llm) )
	987	ALLOCATE(plsnc(iip1,jjb_u:jje_u,llm) )
	988	ALLOCATE(plvnc(iip1,jjb_v:jje_v,llm) )
	989	ALLOCATE(p(iip1,jjb_u:jje_u,llmp1) )
	990	ALLOCATE(pls(iip1,jjb_u:jje_u,llm) )
	991	ALLOCATE(pext(iip1,jjb_u:jje_u,llm) )
	992	ALLOCATE(pbarx(iip1,jjb_u:jje_u,llm) )
	993	ALLOCATE(pbary(iip1,jjb_v:jje_v,llm) )
	994	ALLOCATE(pk(iip1,jjb_u:jje_u,llm) )
	995	ALLOCATE(pks (iip1,jjb_u:jje_u) )
	996	ALLOCATE(qsat(ijb_u:ije_u,llm) )
	997	ALLOCATE(zu1(iip1,jjb_u:jje_u,llm) )
	998	ALLOCATE(zu2(iip1,jjb_u:jje_u,llm) )
	999	ALLOCATE(zv1(iip1,jjb_v:jje_v,llm) )
[1806]	1000	ALLOCATE(zv2(iip1,jjb_v:jje_v,llm) )
	1001	!$OMP END MASTER
	1002	!$OMP BARRIER
[5249]	1003	ENDIF
[1806]	1004
[5249]	1005
	1006
	1007
[5084]	1008	IF (guide_plevs.EQ.0) THEN
[1806]	1009	!$OMP DO
[1632]	1010	DO l=1,nlevnc
	1011	DO j=jjbu,jjeu
	1012	DO i=1,iip1
	1013	plnc2(i,j,l)=apnc(l)
	1014	plnc1(i,j,l)=apnc(l)
	1015	ENDDO
	1016	ENDDO
	1017	ENDDO
[5249]	1018	ENDIF
[1632]	1019
	1020	if (first) then
	1021	first=.FALSE.
[1806]	1022	!$OMP MASTER
[3995]	1023	write(,)trim(modname)//' : check vertical level order'
	1024	write(,)trim(modname)//' LMDZ :'
[1632]	1025	do l=1,llm
[3995]	1026	write(,)trim(modname)//' PL(',l,')=',(ap(l)+ap(l+1))/2. &
[1632]	1027	+psi(1,jjeu)*(bp(l)+bp(l+1))/2.
	1028	enddo
[3995]	1029	write(,)trim(modname)//' nudging file :'
[1632]	1030	SELECT CASE (guide_plevs)
[5249]	1031	CASE (0)
[1632]	1032	do l=1,nlevnc
[3995]	1033	write(,)trim(modname)//' PL(',l,')=',plnc2(1,jjbu,l)
[1632]	1034	enddo
	1035	CASE (1)
	1036	DO l=1,nlevnc
[3995]	1037	write(,)trim(modname)//' PL(',l,')=',&
[4041]	1038	apnc(l)+bpnc(l)*psnat2(1,jjbu)
[3995]	1039	ENDDO
[1632]	1040	CASE (2)
	1041	do l=1,nlevnc
[3995]	1042	write(,)trim(modname)//' PL(',l,')=',pnat2(1,jjbu,l)
[1632]	1043	enddo
	1044	END SELECT
[3995]	1045	write(,)trim(modname)//' invert ordering: invert_p=',invert_p
[1632]	1046	if (guide_u) then
	1047	do l=1,nlevnc
[3995]	1048	write(,)trim(modname)//' U(',l,')=',unat2(1,jjbu,l)
[1632]	1049	enddo
	1050	endif
	1051	if (guide_T) then
	1052	do l=1,nlevnc
[3995]	1053	write(,)trim(modname)//' T(',l,')=',tnat2(1,jjbu,l)
[1632]	1054	enddo
	1055	endif
[1806]	1056	!$OMP END MASTER
[3995]	1057	endif ! of if (first)
[4042]	1058
[4054]	1059	if (guide_plevs /= 1 .or. guide_t .and. .not. guide_teta &
	1060	.or. guide_q .and. guide_hr) then
	1061	CALL pression_loc( ijnb_u, ap, bp, psi, p )
	1062	if (disvert_type==1) then
	1063	CALL exner_hyb_loc(ijnb_u,psi,p,pks,pk)
	1064	else ! we assume that we are in the disvert_type==2 case
	1065	CALL exner_milieu_loc(ijnb_u,psi,p,pks,pk)
	1066	endif
	1067	end if
	1068
[1632]	1069	! -----------------------------------------------------------------
[5249]	1070	! Calcul niveaux pression modele
[1632]	1071	! -----------------------------------------------------------------
	1072
	1073	! .... Calcul de pls , pression au milieu des couches ,en Pascals
[5084]	1074	IF (guide_plevs.EQ.1) THEN
[2036]	1075	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1632]	1076	DO l=1,llm
[2597]	1077	DO j=jjbu,jjeu
	1078	DO i =1, iip1
[1632]	1079	pls(i,j,l)=(ap(l)+ap(l+1))/2.+psi(i,j)*(bp(l)+bp(l+1))/2.
[2028]	1080	ENDDO
	1081	ENDDO
[1632]	1082	ENDDO
	1083	ELSE
[2028]	1084	unskap=1./kappa
[1806]	1085	!$OMP BARRIER
[2036]	1086	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[2028]	1087	DO l = 1, llm
	1088	DO j=jjbu,jjeu
[2597]	1089	DO i =1, iip1
	1090	pls(i,j,l) = preff * ( pk(i,j,l)/cpp) ** unskap
	1091	ENDDO
[2028]	1092	ENDDO
	1093	ENDDO
[1632]	1094	ENDIF
	1095
	1096	! calcul des pressions pour les grilles u et v
[2036]	1097	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1632]	1098	do l=1,llm
	1099	do j=jjbu,jjeu
	1100	do i=1,iip1
	1101	pext(i,j,l)=pls(i,j,l)*aire(i,j)
	1102	enddo
	1103	enddo
	1104	enddo
	1105
	1106	CALL Register_Hallo_u(pext,llm,1,2,2,1,Req)
	1107	CALL SendRequest(Req)
[1806]	1108	!$OMP BARRIER
[1632]	1109	CALL WaitRequest(Req)
[1806]	1110	!$OMP BARRIER
[1632]	1111
	1112	call massbar_loc(pext, pbarx, pbary )
[1806]	1113	!$OMP BARRIER
[2036]	1114	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1632]	1115	do l=1,llm
	1116	do j=jjbu,jjeu
	1117	do i=1,iip1
	1118	plunc(i,j,l)=pbarx(i,j,l)/aireu(i,j)
	1119	plsnc(i,j,l)=pls(i,j,l)
	1120	enddo
	1121	enddo
	1122	enddo
[2036]	1123	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1632]	1124	do l=1,llm
	1125	do j=jjbv,jjev
	1126	do i=1,iip1
	1127	plvnc(i,j,l)=pbary(i,j,l)/airev(i,j)
	1128	enddo
	1129	enddo
	1130	enddo
	1131
	1132	! -----------------------------------------------------------------
	1133	! Interpolation verticale champs guidage sur niveaux modele
	1134	! Conversion en variables gcm (ucov, vcov...)
	1135	! -----------------------------------------------------------------
	1136	if (guide_P) then
[1806]	1137	!$OMP MASTER
[1632]	1138	do j=jjbu,jjeu
	1139	do i=1,iim
	1140	ij=(j-1)*iip1+i
	1141	psgui1(ij)=psnat1(i,j)
	1142	psgui2(ij)=psnat2(i,j)
	1143	enddo
	1144	psgui1(iip1*j)=psnat1(1,j)
	1145	psgui2(iip1*j)=psnat2(1,j)
	1146	enddo
[1806]	1147	!$OMP END MASTER
	1148	!$OMP BARRIER
[1632]	1149	endif
	1150
	1151	IF (guide_T) THEN
	1152	! Calcul des nouvelles valeurs des niveaux de pression du guidage
[5084]	1153	IF (guide_plevs.EQ.1) THEN
[1806]	1154	!$OMP DO
[1632]	1155	DO l=1,nlevnc
	1156	DO j=jjbu,jjeu
	1157	DO i=1,iip1
	1158	plnc2(i,j,l)=apnc(l)+bpnc(l)*psnat2(i,j)
	1159	plnc1(i,j,l)=apnc(l)+bpnc(l)*psnat1(i,j)
	1160	ENDDO
	1161	ENDDO
	1162	ENDDO
[5084]	1163	ELSE IF (guide_plevs.EQ.2) THEN
[1806]	1164	!$OMP DO
[1632]	1165	DO l=1,nlevnc
	1166	DO j=jjbu,jjeu
	1167	DO i=1,iip1
	1168	plnc2(i,j,l)=pnat2(i,j,l)
	1169	plnc1(i,j,l)=pnat1(i,j,l)
	1170	ENDDO
	1171	ENDDO
	1172	ENDDO
	1173	ENDIF
	1174
	1175	! Interpolation verticale
[1806]	1176	!$OMP MASTER
[1632]	1177	CALL pres2lev(tnat1(:,jjbu:jjeu,:),zu1(:,jjbu:jjeu,:),nlevnc,llm, &
	1178	plnc1(:,jjbu:jjeu,:),plsnc(:,jjbu:jjeu,:),iip1,jjnu,invert_p)
	1179	CALL pres2lev(tnat2(:,jjbu:jjeu,:),zu2(:,jjbu:jjeu,:),nlevnc,llm, &
	1180	plnc2(:,jjbu:jjeu,:),plsnc(:,jjbu:jjeu,:),iip1,jjnu,invert_p)
[1806]	1181	!$OMP END MASTER
	1182	!$OMP BARRIER
[1632]	1183	! Conversion en variables GCM
[2036]	1184	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1632]	1185	do l=1,llm
	1186	do j=jjbu,jjeu
	1187	IF (guide_teta) THEN
	1188	do i=1,iim
	1189	ij=(j-1)*iip1+i
	1190	tgui1(ij,l)=zu1(i,j,l)
	1191	tgui2(ij,l)=zu2(i,j,l)
	1192	enddo
	1193	ELSE
	1194	do i=1,iim
	1195	ij=(j-1)*iip1+i
	1196	tgui1(ij,l)=zu1(i,j,l)*cpp/pk(i,j,l)
	1197	tgui2(ij,l)=zu2(i,j,l)*cpp/pk(i,j,l)
	1198	enddo
	1199	ENDIF
[5249]	1200	tgui1(jiip1,l)=tgui1((j-1)iip1+1,l)
	1201	tgui2(jiip1,l)=tgui2((j-1)iip1+1,l)
[1632]	1202	enddo
	1203	if (pole_nord) then
	1204	do i=1,iip1
	1205	tgui1(i,l)=tgui1(1,l)
	1206	tgui2(i,l)=tgui2(1,l)
	1207	enddo
	1208	endif
	1209	if (pole_sud) then
	1210	do i=1,iip1
[5249]	1211	tgui1(ip1jm+i,l)=tgui1(ip1jm+1,l)
	1212	tgui2(ip1jm+i,l)=tgui2(ip1jm+1,l)
[1632]	1213	enddo
	1214	endif
	1215	enddo
	1216	ENDIF
	1217
	1218	IF (guide_Q) THEN
	1219	! Calcul des nouvelles valeurs des niveaux de pression du guidage
[5084]	1220	IF (guide_plevs.EQ.1) THEN
[1806]	1221	!$OMP DO
[1632]	1222	DO l=1,nlevnc
	1223	DO j=jjbu,jjeu
	1224	DO i=1,iip1
	1225	plnc2(i,j,l)=apnc(l)+bpnc(l)*psnat2(i,j)
	1226	plnc1(i,j,l)=apnc(l)+bpnc(l)*psnat1(i,j)
	1227	ENDDO
	1228	ENDDO
	1229	ENDDO
[5084]	1230	ELSE IF (guide_plevs.EQ.2) THEN
[1806]	1231	!$OMP DO
[1632]	1232	DO l=1,nlevnc
	1233	DO j=jjbu,jjeu
	1234	DO i=1,iip1
	1235	plnc2(i,j,l)=pnat2(i,j,l)
	1236	plnc1(i,j,l)=pnat1(i,j,l)
	1237	ENDDO
	1238	ENDDO
	1239	ENDDO
	1240	ENDIF
	1241
	1242	! Interpolation verticale
[1806]	1243	!$OMP MASTER
[1632]	1244	CALL pres2lev(qnat1(:,jjbu:jjeu,:),zu1(:,jjbu:jjeu,:),nlevnc,llm, &
	1245	plnc1(:,jjbu:jjeu,:),plsnc(:,jjbu:jjeu,:),iip1,jjnu,invert_p)
	1246	CALL pres2lev(qnat2(:,jjbu:jjeu,:),zu2(:,jjbu:jjeu,:),nlevnc,llm, &
	1247	plnc2(:,jjbu:jjeu,:),plsnc(:,jjbu:jjeu,:),iip1,jjnu,invert_p)
[1806]	1248	!$OMP END MASTER
	1249	!$OMP BARRIER
[1632]	1250
	1251	! Conversion en variables GCM
	1252	! On suppose qu'on a la bonne variable dans le fichier de guidage:
	1253	! Hum.Rel si guide_hr, Hum.Spec. sinon.
[2036]	1254	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1632]	1255	do l=1,llm
	1256	do j=jjbu,jjeu
	1257	do i=1,iim
	1258	ij=(j-1)*iip1+i
	1259	qgui1(ij,l)=zu1(i,j,l)
	1260	qgui2(ij,l)=zu2(i,j,l)
	1261	enddo
[5249]	1262	qgui1(jiip1,l)=qgui1((j-1)iip1+1,l)
	1263	qgui2(jiip1,l)=qgui2((j-1)iip1+1,l)
[1632]	1264	enddo
	1265	if (pole_nord) then
	1266	do i=1,iip1
	1267	qgui1(i,l)=qgui1(1,l)
	1268	qgui2(i,l)=qgui2(1,l)
	1269	enddo
	1270	endif
[3547]	1271	if (pole_sud) then
[1632]	1272	do i=1,iip1
[5249]	1273	qgui1(ip1jm+i,l)=qgui1(ip1jm+1,l)
	1274	qgui2(ip1jm+i,l)=qgui2(ip1jm+1,l)
[1632]	1275	enddo
	1276	endif
	1277	enddo
	1278	IF (guide_hr) THEN
[2036]	1279	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1806]	1280	do l=1,llm
	1281	CALL q_sat(iip1jjnu,teta(:,jjbu:jjeu,l)pk(:,jjbu:jjeu,l)/cpp, &
	1282	plsnc(:,jjbu:jjeu,l),qsat(ijbu:ijeu,l))
	1283	qgui1(ijbu:ijeu,l)=qgui1(ijbu:ijeu,l)qsat(ijbu:ijeu,l)0.01 !hum. rel. en %
[5249]	1284	qgui2(ijbu:ijeu,l)=qgui2(ijbu:ijeu,l)qsat(ijbu:ijeu,l)0.01
[1806]	1285	enddo
	1286
[1632]	1287	ENDIF
	1288	ENDIF
	1289
	1290	IF (guide_u) THEN
	1291	! Calcul des nouvelles valeurs des niveaux de pression du guidage
[5084]	1292	IF (guide_plevs.EQ.1) THEN
[1806]	1293	!$OMP DO
[1632]	1294	DO l=1,nlevnc
	1295	DO j=jjbu,jjeu
	1296	DO i=1,iim
	1297	plnc2(i,j,l)=apnc(l)+bpnc(l)(psnat2(i,j)aire(i,j)*alpha1p2(i,j) &
	1298	& +psnat2(i+1,j)aire(i+1,j)alpha3p4(i+1,j))/aireu(i,j)
	1299	plnc1(i,j,l)=apnc(l)+bpnc(l)(psnat1(i,j)aire(i,j)*alpha1p2(i,j) &
	1300	& +psnat1(i+1,j)aire(i+1,j)alpha3p4(i+1,j))/aireu(i,j)
	1301	ENDDO
	1302	plnc2(iip1,j,l)=plnc2(1,j,l)
	1303	plnc1(iip1,j,l)=plnc1(1,j,l)
	1304	ENDDO
	1305	ENDDO
[5084]	1306	ELSE IF (guide_plevs.EQ.2) THEN
[1806]	1307	!$OMP DO
[1632]	1308	DO l=1,nlevnc
	1309	DO j=jjbu,jjeu
	1310	DO i=1,iim
	1311	plnc2(i,j,l)=(pnat2(i,j,l)aire(i,j)alpha1p2(i,j) &
	1312	& +pnat2(i+1,j,l)aire(i,j)alpha3p4(i+1,j))/aireu(i,j)
	1313	plnc1(i,j,l)=(pnat1(i,j,l)aire(i,j)alpha1p2(i,j) &
	1314	& +pnat1(i+1,j,l)aire(i,j)alpha3p4(i+1,j))/aireu(i,j)
	1315	ENDDO
	1316	plnc2(iip1,j,l)=plnc2(1,j,l)
	1317	plnc1(iip1,j,l)=plnc1(1,j,l)
	1318	ENDDO
	1319	ENDDO
	1320	ENDIF
[5249]	1321
[1632]	1322	! Interpolation verticale
[1806]	1323	!$OMP MASTER
[1632]	1324	CALL pres2lev(unat1(:,jjbu:jjeu,:),zu1(:,jjbu:jjeu,:),nlevnc,llm, &
	1325	plnc1(:,jjbu:jjeu,:),plunc(:,jjbu:jjeu,:),iip1,jjnu,invert_p)
	1326	CALL pres2lev(unat2(:,jjbu:jjeu,:),zu2(:,jjbu:jjeu,:),nlevnc,llm, &
	1327	plnc2(:,jjbu:jjeu,:),plunc(:,jjbu:jjeu,:),iip1,jjnu,invert_p)
[1806]	1328	!$OMP END MASTER
	1329	!$OMP BARRIER
[1632]	1330
	1331	! Conversion en variables GCM
[2036]	1332	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1632]	1333	do l=1,llm
	1334	do j=jjbu,jjeu
	1335	do i=1,iim
	1336	ij=(j-1)*iip1+i
	1337	ugui1(ij,l)=zu1(i,j,l)*cu(i,j)
	1338	ugui2(ij,l)=zu2(i,j,l)*cu(i,j)
	1339	enddo
[5249]	1340	ugui1(jiip1,l)=ugui1((j-1)iip1+1,l)
	1341	ugui2(jiip1,l)=ugui2((j-1)iip1+1,l)
[1632]	1342	enddo
	1343	if (pole_nord) then
	1344	do i=1,iip1
	1345	ugui1(i,l)=0.
	1346	ugui2(i,l)=0.
	1347	enddo
	1348	endif
	1349	if (pole_sud) then
	1350	do i=1,iip1
	1351	ugui1(ip1jm+i,l)=0.
	1352	ugui2(ip1jm+i,l)=0.
	1353	enddo
	1354	endif
	1355	enddo
	1356	ENDIF
[5249]	1357
[1632]	1358	IF (guide_v) THEN
	1359	! Calcul des nouvelles valeurs des niveaux de pression du guidage
[5084]	1360	IF (guide_plevs.EQ.1) THEN
[2597]	1361	CALL Register_Hallo_u(psnat1,1,1,2,2,1,Req)
	1362	CALL Register_Hallo_u(psnat2,1,1,2,2,1,Req)
	1363	CALL SendRequest(Req)
[1806]	1364	!$OMP BARRIER
[2597]	1365	CALL WaitRequest(Req)
[1806]	1366	!$OMP BARRIER
	1367	!$OMP DO
[1632]	1368	DO l=1,nlevnc
	1369	DO j=jjbv,jjev
	1370	DO i=1,iip1
	1371	plnc2(i,j,l)=apnc(l)+bpnc(l)(psnat2(i,j)aire(i,j)*alpha2p3(i,j) &
	1372	& +psnat2(i,j+1)aire(i,j+1)alpha1p4(i,j+1))/airev(i,j)
	1373	plnc1(i,j,l)=apnc(l)+bpnc(l)(psnat1(i,j)aire(i,j)*alpha2p3(i,j) &
	1374	& +psnat1(i,j+1)aire(i,j+1)alpha1p4(i,j+1))/airev(i,j)
	1375	ENDDO
	1376	ENDDO
	1377	ENDDO
[5084]	1378	ELSE IF (guide_plevs.EQ.2) THEN
[2597]	1379	CALL Register_Hallo_u(pnat1,llm,1,2,2,1,Req)
	1380	CALL Register_Hallo_u(pnat2,llm,1,2,2,1,Req)
	1381	CALL SendRequest(Req)
[1806]	1382	!$OMP BARRIER
[2597]	1383	CALL WaitRequest(Req)
[1806]	1384	!$OMP BARRIER
	1385	!$OMP DO
[1632]	1386	DO l=1,nlevnc
	1387	DO j=jjbv,jjev
	1388	DO i=1,iip1
	1389	plnc2(i,j,l)=(pnat2(i,j,l)aire(i,j)alpha2p3(i,j) &
	1390	& +pnat2(i,j+1,l)aire(i,j)alpha1p4(i,j+1))/airev(i,j)
	1391	plnc1(i,j,l)=(pnat1(i,j,l)aire(i,j)alpha2p3(i,j) &
	1392	& +pnat1(i,j+1,l)aire(i,j)alpha1p4(i,j+1))/airev(i,j)
	1393	ENDDO
	1394	ENDDO
	1395	ENDDO
	1396	ENDIF
	1397	! Interpolation verticale
[1806]	1398
	1399	!$OMP MASTER
[1632]	1400	CALL pres2lev(vnat1(:,jjbv:jjev,:),zv1(:,jjbv:jjev,:),nlevnc,llm, &
	1401	plnc1(:,jjbv:jjev,:),plvnc(:,jjbv:jjev,:),iip1,jjnv,invert_p)
	1402	CALL pres2lev(vnat2(:,jjbv:jjev,:),zv2(:,jjbv:jjev,:),nlevnc,llm, &
	1403	plnc2(:,jjbv:jjev,:),plvnc(:,jjbv:jjev,:),iip1,jjnv,invert_p)
[1806]	1404	!$OMP END MASTER
	1405	!$OMP BARRIER
[1632]	1406	! Conversion en variables GCM
[2036]	1407	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[1632]	1408	do l=1,llm
	1409	do j=jjbv,jjev
	1410	do i=1,iim
	1411	ij=(j-1)*iip1+i
	1412	vgui1(ij,l)=zv1(i,j,l)*cv(i,j)
	1413	vgui2(ij,l)=zv2(i,j,l)*cv(i,j)
	1414	enddo
[5249]	1415	vgui1(jiip1,l)=vgui1((j-1)iip1+1,l)
	1416	vgui2(jiip1,l)=vgui2((j-1)iip1+1,l)
[1632]	1417	enddo
	1418	enddo
	1419	ENDIF
	1420
[5249]	1421
[1632]	1422	END SUBROUTINE guide_interp
	1423
	1424	!=======================================================================
[1806]	1425	SUBROUTINE tau2alpha(typ,pim,jjb,jje,factt,taumin,taumax,alpha)
[1632]	1426
	1427	! Calcul des constantes de rappel alpha (=1/tau)
	1428
[2597]	1429	use comconst_mod, only: pi
[2598]	1430	use serre_mod, only: clat, clon, grossismx, grossismy
[5249]	1431
[5271]	1432	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
[5272]	1433	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	1434	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
[5271]	1435	implicit none
[1632]	1436
[5271]	1437
[5272]	1438
[1632]	1439	include "comgeom2.h"
	1440
	1441	! input arguments :
	1442	INTEGER, INTENT(IN) :: typ ! u(2),v(3), ou scalaire(1)
[1806]	1443	INTEGER, INTENT(IN) :: pim ! dimensions en lon
	1444	INTEGER, INTENT(IN) :: jjb,jje ! dimensions en lat
[1632]	1445	REAL, INTENT(IN) :: factt ! pas de temps en fraction de jour
	1446	REAL, INTENT(IN) :: taumin,taumax
	1447	! output arguments:
[5249]	1448	REAL, DIMENSION(pim,jjb:jje), INTENT(OUT) :: alpha
	1449
[1632]	1450	! local variables:
	1451	LOGICAL, SAVE :: first=.TRUE.
	1452	REAL, SAVE :: gamma,dxdy_min,dxdy_max
	1453	REAL, DIMENSION (iip1,jjp1) :: zdx,zdy
	1454	REAL, DIMENSION (iip1,jjp1) :: dxdys,dxdyu
	1455	REAL, DIMENSION (iip1,jjm) :: dxdyv
	1456	real dxdy_
	1457	real zlat,zlon
	1458	real alphamin,alphamax,xi
	1459	integer i,j,ilon,ilat
[3995]	1460	character(len=20),parameter :: modname="tau2alpha"
[1632]	1461
	1462
	1463	alphamin=factt/taumax
	1464	alphamax=factt/taumin
	1465	IF (guide_reg.OR.guide_add) THEN
	1466	alpha=alphamax
	1467	!-----------------------------------------------------------------------
	1468	! guide_reg: alpha=alpha_min dans region, 0. sinon.
	1469	!-----------------------------------------------------------------------
	1470	IF (guide_reg) THEN
[1806]	1471	do j=jjb,jje
[1632]	1472	do i=1,pim
[5084]	1473	if (typ.eq.2) then
[1632]	1474	zlat=rlatu(j)*180./pi
	1475	zlon=rlonu(i)*180./pi
[5084]	1476	elseif (typ.eq.1) then
[1632]	1477	zlat=rlatu(j)*180./pi
	1478	zlon=rlonv(i)*180./pi
[5084]	1479	elseif (typ.eq.3) then
[1632]	1480	zlat=rlatv(j)*180./pi
	1481	zlon=rlonv(i)*180./pi
	1482	endif
	1483	alpha(i,j)=alphamax/16.* &
	1484	(1.+tanh((zlat-lat_min_g)/tau_lat))* &
	1485	(1.+tanh((lat_max_g-zlat)/tau_lat))* &
	1486	(1.+tanh((zlon-lon_min_g)/tau_lon))* &
	1487	(1.+tanh((lon_max_g-zlon)/tau_lon))
	1488	enddo
	1489	enddo
	1490	ENDIF
	1491	ELSE
	1492	!-----------------------------------------------------------------------
	1493	! Sinon, alpha varie entre alpha_min et alpha_max suivant le zoom.
	1494	!-----------------------------------------------------------------------
	1495	!Calcul de l'aire des mailles
	1496	do j=2,jjm
	1497	do i=2,iip1
	1498	zdx(i,j)=0.5*(cu(i-1,j)+cu(i,j))/cos(rlatu(j))
	1499	enddo
	1500	zdx(1,j)=zdx(iip1,j)
	1501	enddo
	1502	do j=2,jjm
	1503	do i=1,iip1
	1504	zdy(i,j)=0.5*(cv(i,j-1)+cv(i,j))
	1505	enddo
	1506	enddo
	1507	do i=1,iip1
	1508	zdx(i,1)=zdx(i,2)
	1509	zdx(i,jjp1)=zdx(i,jjm)
	1510	zdy(i,1)=zdy(i,2)
	1511	zdy(i,jjp1)=zdy(i,jjm)
	1512	enddo
	1513	do j=1,jjp1
	1514	do i=1,iip1
	1515	dxdys(i,j)=sqrt(zdx(i,j)zdx(i,j)+zdy(i,j)zdy(i,j))
	1516	enddo
	1517	enddo
[5084]	1518	IF (typ.EQ.2) THEN
[1632]	1519	do j=1,jjp1
	1520	do i=1,iim
	1521	dxdyu(i,j)=0.5*(dxdys(i,j)+dxdys(i+1,j))
	1522	enddo
	1523	dxdyu(iip1,j)=dxdyu(1,j)
	1524	enddo
	1525	ENDIF
[5084]	1526	IF (typ.EQ.3) THEN
[1632]	1527	do j=1,jjm
	1528	do i=1,iip1
	1529	dxdyv(i,j)=0.5*(dxdys(i,j)+dxdys(i,j+1))
	1530	enddo
	1531	enddo
	1532	ENDIF
	1533	! Premier appel: calcul des aires min et max et de gamma.
[5249]	1534	IF (first) THEN
[1632]	1535	first=.FALSE.
	1536	! coordonnees du centre du zoom
[5249]	1537	CALL coordij(clon,clat,ilon,ilat)
[1632]	1538	! aire de la maille au centre du zoom
	1539	dxdy_min=dxdys(ilon,ilat)
	1540	! dxdy maximale de la maille
	1541	dxdy_max=0.
	1542	do j=1,jjp1
	1543	do i=1,iip1
	1544	dxdy_max=max(dxdy_max,dxdys(i,j))
	1545	enddo
	1546	enddo
	1547	! Calcul de gamma
[5084]	1548	if (abs(grossismx-1.).lt.0.1.or.abs(grossismy-1.).lt.0.1) then
[3995]	1549	write(,)trim(modname)//' ATTENTION modele peu zoome'
	1550	write(,)trim(modname)//' ATTENTION on prend une constante de guidage cste'
	1551	gamma=0.
[1632]	1552	else
[3995]	1553	gamma=(dxdy_max-2.*dxdy_min)/(dxdy_max-dxdy_min)
	1554	write(,)trim(modname)//' gamma=',gamma
[5084]	1555	if (gamma.lt.1.e-5) then
[3995]	1556	write(,)trim(modname)//' gamma =',gamma,'<1e-5'
[4469]	1557	CALL abort_gcm("guide_loc_mod","stopped",1)
[3995]	1558	endif
	1559	gamma=log(0.5)/log(gamma)
[5249]	1560	if (gamma4) then
[3995]	1561	gamma=min(gamma,4.)
	1562	endif
	1563	write(,)trim(modname)//' gamma=',gamma
[1632]	1564	endif
	1565	ENDIF !first
	1566
[1806]	1567	do j=jjb,jje
[1632]	1568	do i=1,pim
[5084]	1569	if (typ.eq.1) then
[1632]	1570	dxdy_=dxdys(i,j)
	1571	zlat=rlatu(j)*180./pi
[5084]	1572	elseif (typ.eq.2) then
[1632]	1573	dxdy_=dxdyu(i,j)
	1574	zlat=rlatu(j)*180./pi
[5084]	1575	elseif (typ.eq.3) then
[1632]	1576	dxdy_=dxdyv(i,j)
	1577	zlat=rlatv(j)*180./pi
	1578	endif
[5084]	1579	if (abs(grossismx-1.).lt.0.1.or.abs(grossismy-1.).lt.0.1) then
[1632]	1580	! pour une grille reguliere, xi=xxx**0=1 -> alpha=alphamin
	1581	alpha(i,j)=alphamin
	1582	else
	1583	xi=((dxdy_max-dxdy_)/(dxdy_max-dxdy_min))**gamma
	1584	xi=min(xi,1.)
[5084]	1585	if(lat_min_g.le.zlat .and. zlat.le.lat_max_g) then
[1632]	1586	alpha(i,j)=xialphamin+(1.-xi)alphamax
	1587	else
	1588	alpha(i,j)=0.
	1589	endif
	1590	endif
	1591	enddo
	1592	enddo
	1593	ENDIF ! guide_reg
	1594
[2134]	1595	if (.not. guide_add) alpha = 1. - exp(- alpha)
	1596
[1632]	1597	END SUBROUTINE tau2alpha
	1598
	1599	!=======================================================================
	1600	SUBROUTINE guide_read(timestep)
[5249]	1601	USE netcdf, ONLY: nf90_put_var
[5271]	1602	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
[5272]	1603	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	1604	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
[5271]	1605	IMPLICIT NONE
	1606
[1632]	1607
[5272]	1608
[1632]	1609	INTEGER, INTENT(IN) :: timestep
	1610
	1611	LOGICAL, SAVE :: first=.TRUE.
	1612	! Identification fichiers et variables NetCDF:
	1613	INTEGER, SAVE :: ncidu,varidu,ncidv,varidv,ncidp,varidp
	1614	INTEGER, SAVE :: ncidQ,varidQ,ncidt,varidt,ncidps,varidps
[3984]	1615	INTEGER :: ncidpl,varidpl,varidap,varidbp,dimid,lendim
[1632]	1616	! Variables auxiliaires NetCDF:
	1617	INTEGER, DIMENSION(4) :: start,count
	1618	INTEGER :: status,rcode
	1619	CHARACTER (len = 80) :: abort_message
	1620	CHARACTER (len = 20) :: modname = 'guide_read'
[3984]	1621	CHARACTER (len = 20) :: namedim
[1632]	1622	abort_message='pb in guide_read'
	1623
	1624	! -----------------------------------------------------------------
	1625	! Premier appel: initialisation de la lecture des fichiers
	1626	! -----------------------------------------------------------------
	1627	if (first) then
	1628	ncidpl=-99
[3995]	1629	write(,),trim(modname)//': opening nudging files '
[1632]	1630	! Ap et Bp si Niveaux de pression hybrides
[5084]	1631	if (guide_plevs.EQ.1) then
[3995]	1632	write(,),trim(modname)//' Reading nudging on model levels'
[1632]	1633	rcode = nf90_open('apbp.nc', nf90_nowrite, ncidpl)
[5270]	1634	IF (rcode.NE.nf90_noerr) THEN
[3995]	1635	abort_message='Nudging: error -> no file apbp.nc'
[1902]	1636	CALL abort_gcm(modname,abort_message,1)
	1637	ENDIF
[1632]	1638	rcode = nf90_inq_varid(ncidpl, 'AP', varidap)
[5270]	1639	IF (rcode.NE.nf90_noerr) THEN
[3995]	1640	abort_message='Nudging: error -> no AP variable in file apbp.nc'
[1902]	1641	CALL abort_gcm(modname,abort_message,1)
	1642	ENDIF
[1632]	1643	rcode = nf90_inq_varid(ncidpl, 'BP', varidbp)
[5270]	1644	IF (rcode.NE.nf90_noerr) THEN
[3995]	1645	abort_message='Nudging: error -> no BP variable in file apbp.nc'
[1902]	1646	CALL abort_gcm(modname,abort_message,1)
	1647	ENDIF
[3995]	1648	write(,),trim(modname)//' ncidpl,varidap',ncidpl,varidap
[1632]	1649	endif
[5249]	1650
[1632]	1651	! Pression si guidage sur niveaux P variables
[5084]	1652	if (guide_plevs.EQ.2) then
[1632]	1653	rcode = nf90_open('P.nc', nf90_nowrite, ncidp)
[5270]	1654	IF (rcode.NE.nf90_noerr) THEN
[3995]	1655	abort_message='Nudging: error -> no file P.nc'
[1902]	1656	CALL abort_gcm(modname,abort_message,1)
	1657	ENDIF
[1632]	1658	rcode = nf90_inq_varid(ncidp, 'PRES', varidp)
[5270]	1659	IF (rcode.NE.nf90_noerr) THEN
[3995]	1660	abort_message='Nudging: error -> no PRES variable in file P.nc'
[1902]	1661	CALL abort_gcm(modname,abort_message,1)
	1662	ENDIF
[3995]	1663	write(,),trim(modname)//' ncidp,varidp',ncidp,varidp
[5084]	1664	if (ncidpl.eq.-99) ncidpl=ncidp
[1632]	1665	endif
[3984]	1666
[1632]	1667	! Vent zonal
	1668	if (guide_u) then
	1669	rcode = nf90_open('u.nc', nf90_nowrite, ncidu)
[5270]	1670	IF (rcode.NE.nf90_noerr) THEN
[3995]	1671	abort_message='Nudging: error -> no file u.nc'
[1902]	1672	CALL abort_gcm(modname,abort_message,1)
	1673	ENDIF
[1632]	1674	rcode = nf90_inq_varid(ncidu, 'UWND', varidu)
[5270]	1675	IF (rcode.NE.nf90_noerr) THEN
[3995]	1676	abort_message='Nudging: error -> no UWND variable in file u.nc'
[1902]	1677	CALL abort_gcm(modname,abort_message,1)
	1678	ENDIF
[3995]	1679	write(,),trim(modname)//' ncidu,varidu',ncidu,varidu
[5084]	1680	if (ncidpl.eq.-99) ncidpl=ncidu
[3984]	1681
[5249]	1682
[3984]	1683	status=NF90_INQ_DIMID(ncidu, "LONU", dimid)
	1684	status=NF90_INQUIRE_DIMENSION(ncidu,dimid,namedim,lendim)
[5084]	1685	IF (lendim .NE. iip1) THEN
[3995]	1686	abort_message='dimension LONU different from iip1 in u.nc'
[3984]	1687	CALL abort_gcm(modname,abort_message,1)
	1688	ENDIF
	1689
	1690	status=NF90_INQ_DIMID(ncidu, "LATU", dimid)
	1691	status=NF90_INQUIRE_DIMENSION(ncidu,dimid,namedim,lendim)
[5084]	1692	IF (lendim .NE. jjp1) THEN
[3995]	1693	abort_message='dimension LATU different from jjp1 in u.nc'
[3984]	1694	CALL abort_gcm(modname,abort_message,1)
	1695	ENDIF
[5249]	1696
[1632]	1697	endif
[3984]	1698
[1632]	1699	! Vent meridien
	1700	if (guide_v) then
	1701	rcode = nf90_open('v.nc', nf90_nowrite, ncidv)
[5270]	1702	IF (rcode.NE.nf90_noerr) THEN
[3995]	1703	abort_message='Nudging: error -> no file v.nc'
[1902]	1704	CALL abort_gcm(modname,abort_message,1)
	1705	ENDIF
[1632]	1706	rcode = nf90_inq_varid(ncidv, 'VWND', varidv)
[5270]	1707	IF (rcode.NE.nf90_noerr) THEN
[3995]	1708	abort_message='Nudging: error -> no VWND variable in file v.nc'
[1902]	1709	CALL abort_gcm(modname,abort_message,1)
	1710	ENDIF
[3995]	1711	write(,),trim(modname)//' ncidv,varidv',ncidv,varidv
[5084]	1712	if (ncidpl.eq.-99) ncidpl=ncidv
[5249]	1713
[3984]	1714	status=NF90_INQ_DIMID(ncidv, "LONV", dimid)
	1715	status=NF90_INQUIRE_DIMENSION(ncidv,dimid,namedim,lendim)
[5249]	1716
[5084]	1717	IF (lendim .NE. iip1) THEN
[3995]	1718	abort_message='dimension LONV different from iip1 in v.nc'
[3984]	1719	CALL abort_gcm(modname,abort_message,1)
	1720	ENDIF
	1721
	1722
	1723	status=NF90_INQ_DIMID(ncidv, "LATV", dimid)
	1724	status=NF90_INQUIRE_DIMENSION(ncidv,dimid,namedim,lendim)
[5084]	1725	IF (lendim .NE. jjm) THEN
[3995]	1726	abort_message='dimension LATV different from jjm in v.nc'
[3984]	1727	CALL abort_gcm(modname,abort_message,1)
	1728	ENDIF
[5249]	1729
[3984]	1730	endif
	1731
[1632]	1732	! Temperature
	1733	if (guide_T) then
	1734	rcode = nf90_open('T.nc', nf90_nowrite, ncidt)
[5270]	1735	IF (rcode.NE.nf90_noerr) THEN
[3995]	1736	abort_message='Nudging: error -> no file T.nc'
[1902]	1737	CALL abort_gcm(modname,abort_message,1)
	1738	ENDIF
[1632]	1739	rcode = nf90_inq_varid(ncidt, 'AIR', varidt)
[5270]	1740	IF (rcode.NE.nf90_noerr) THEN
[3995]	1741	abort_message='Nudging: error -> no AIR variable in file T.nc'
[1902]	1742	CALL abort_gcm(modname,abort_message,1)
	1743	ENDIF
[3995]	1744	write(,),trim(modname)//' ncidT,varidT',ncidt,varidt
[5084]	1745	if (ncidpl.eq.-99) ncidpl=ncidt
[3984]	1746
	1747	status=NF90_INQ_DIMID(ncidt, "LONV", dimid)
	1748	status=NF90_INQUIRE_DIMENSION(ncidt,dimid,namedim,lendim)
[5084]	1749	IF (lendim .NE. iip1) THEN
[3995]	1750	abort_message='dimension LONV different from iip1 in T.nc'
[3984]	1751	CALL abort_gcm(modname,abort_message,1)
	1752	ENDIF
	1753
	1754	status=NF90_INQ_DIMID(ncidt, "LATU", dimid)
	1755	status=NF90_INQUIRE_DIMENSION(ncidt,dimid,namedim,lendim)
[5084]	1756	IF (lendim .NE. jjp1) THEN
[3995]	1757	abort_message='dimension LATU different from jjp1 in T.nc'
[3984]	1758	CALL abort_gcm(modname,abort_message,1)
	1759	ENDIF
	1760
[1632]	1761	endif
[3984]	1762
[1632]	1763	! Humidite
	1764	if (guide_Q) then
	1765	rcode = nf90_open('hur.nc', nf90_nowrite, ncidQ)
[5270]	1766	IF (rcode.NE.nf90_noerr) THEN
[3995]	1767	abort_message='Nudging: error -> no file hur.nc'
[1902]	1768	CALL abort_gcm(modname,abort_message,1)
	1769	ENDIF
[1632]	1770	rcode = nf90_inq_varid(ncidQ, 'RH', varidQ)
[5270]	1771	IF (rcode.NE.nf90_noerr) THEN
[3995]	1772	abort_message='Nudging: error -> no RH variable in file hur.nc'
[1902]	1773	CALL abort_gcm(modname,abort_message,1)
	1774	ENDIF
[3995]	1775	write(,),trim(modname)//' ncidQ,varidQ',ncidQ,varidQ
[5084]	1776	if (ncidpl.eq.-99) ncidpl=ncidQ
[3984]	1777
	1778
	1779	status=NF90_INQ_DIMID(ncidQ, "LONV", dimid)
	1780	status=NF90_INQUIRE_DIMENSION(ncidQ,dimid,namedim,lendim)
[5084]	1781	IF (lendim .NE. iip1) THEN
[3995]	1782	abort_message='dimension LONV different from iip1 in hur.nc'
[3984]	1783	CALL abort_gcm(modname,abort_message,1)
	1784	ENDIF
	1785
	1786	status=NF90_INQ_DIMID(ncidQ, "LATU", dimid)
	1787	status=NF90_INQUIRE_DIMENSION(ncidQ,dimid,namedim,lendim)
[5084]	1788	IF (lendim .NE. jjp1) THEN
[3995]	1789	abort_message='dimension LATU different from jjp1 in hur.nc'
[3984]	1790	CALL abort_gcm(modname,abort_message,1)
	1791	ENDIF
	1792
	1793
[1632]	1794	endif
	1795	! Pression de surface
[5084]	1796	if ((guide_P).OR.(guide_plevs.EQ.1)) then
[1632]	1797	rcode = nf90_open('ps.nc', nf90_nowrite, ncidps)
[5270]	1798	IF (rcode.NE.nf90_noerr) THEN
[3995]	1799	abort_message='Nudging: error -> no file ps.nc'
[1902]	1800	CALL abort_gcm(modname,abort_message,1)
	1801	ENDIF
[1632]	1802	rcode = nf90_inq_varid(ncidps, 'SP', varidps)
[5270]	1803	IF (rcode.NE.nf90_noerr) THEN
[3995]	1804	abort_message='Nudging: error -> no SP variable in file ps.nc'
[1902]	1805	CALL abort_gcm(modname,abort_message,1)
	1806	ENDIF
[3995]	1807	write(,),trim(modname)//' ncidps,varidps',ncidps,varidps
[1632]	1808	endif
	1809	! Coordonnee verticale
[5084]	1810	if (guide_plevs.EQ.0) then
[1632]	1811	rcode = nf90_inq_varid(ncidpl, 'LEVEL', varidpl)
[5084]	1812	IF (rcode.NE.0) rcode = nf90_inq_varid(ncidpl, 'PRESSURE', varidpl)
[3995]	1813	write(,),trim(modname)//' ncidpl,varidpl',ncidpl,varidpl
[1632]	1814	endif
	1815	! Coefs ap, bp pour calcul de la pression aux differents niveaux
[5084]	1816	IF (guide_plevs.EQ.1) THEN
[5249]	1817	status = nf90_put_var(ncidpl, varidap, apnc, [1], [nlevnc])
	1818	status = nf90_put_var(ncidpl, varidbp, bpnc, [1], [nlevnc])
[5084]	1819	ELSEIF (guide_plevs.EQ.0) THEN
[5249]	1820	status = nf90_put_var(ncidpl, varidpl, apnc, [1], [nlevnc])
[3795]	1821	!FC Pour les corrections la pression est deja en Pascals on commente la ligne ci-dessous
	1822	IF(convert_Pa) apnc=apnc*100.! conversion en Pascals
[1632]	1823	bpnc(:)=0.
	1824	ENDIF
	1825	first=.FALSE.
	1826	ENDIF ! (first)
	1827
	1828	! -----------------------------------------------------------------
	1829	! lecture des champs u, v, T, Q, ps
	1830	! -----------------------------------------------------------------
	1831
	1832	! dimensions pour les champs scalaires et le vent zonal
	1833	start(1)=1
	1834	start(2)=jjb_u
	1835	start(3)=1
	1836	start(4)=timestep
	1837
	1838	count(1)=iip1
	1839	count(2)=jjnb_u
	1840	count(3)=nlevnc
	1841	count(4)=1
	1842
	1843	IF (invert_y) start(2)=jjp1-jje_u+1
[5249]	1844	! Pression
[5084]	1845	if (guide_plevs.EQ.2) then
[5249]	1846	status = nf90_put_var(ncidp, varidp, pnat2, start, count)
[1632]	1847	IF (invert_y) THEN
	1848	! PRINT*,"Invertion impossible actuellement"
	1849	! CALL abort_gcm(modname,abort_message,1)
	1850	CALL invert_lat(iip1,jjnb_u,nlevnc,pnat2)
	1851	ENDIF
	1852	endif
	1853
	1854	! Vent zonal
	1855	if (guide_u) then
[5249]	1856	status = nf90_put_var(ncidu, varidu, unat2, start, count)
[1632]	1857	IF (invert_y) THEN
	1858	! PRINT*,"Invertion impossible actuellement"
	1859	! CALL abort_gcm(modname,abort_message,1)
	1860	CALL invert_lat(iip1,jjnb_u,nlevnc,unat2)
	1861	ENDIF
	1862
	1863	endif
	1864
[2028]	1865
[1632]	1866	! Temperature
	1867	if (guide_T) then
[5249]	1868	status = nf90_put_var(ncidt, varidt, tnat2, start, count)
[1632]	1869	IF (invert_y) THEN
	1870	! PRINT*,"Invertion impossible actuellement"
	1871	! CALL abort_gcm(modname,abort_message,1)
	1872	CALL invert_lat(iip1,jjnb_u,nlevnc,tnat2)
	1873	ENDIF
	1874	endif
	1875
	1876	! Humidite
	1877	if (guide_Q) then
[5249]	1878	status = nf90_put_var(ncidQ, varidQ, qnat2, start, count)
[1632]	1879	IF (invert_y) THEN
	1880	! PRINT*,"Invertion impossible actuellement"
	1881	! CALL abort_gcm(modname,abort_message,1)
	1882	CALL invert_lat(iip1,jjnb_u,nlevnc,qnat2)
	1883	ENDIF
	1884
	1885	endif
	1886
	1887	! Vent meridien
	1888	if (guide_v) then
	1889	start(2)=jjb_v
	1890	count(2)=jjnb_v
	1891	IF (invert_y) start(2)=jjm-jje_v+1
[5249]	1892	status = nf90_put_var(ncidv, varidv, vnat2, start, count)
[1632]	1893	IF (invert_y) THEN
	1894	! PRINT*,"Invertion impossible actuellement"
	1895	! CALL abort_gcm(modname,abort_message,1)
	1896	CALL invert_lat(iip1,jjnb_v,nlevnc,vnat2)
	1897	ENDIF
	1898	endif
	1899
	1900	! Pression de surface
[5084]	1901	if ((guide_P).OR.(guide_plevs.EQ.1)) then
[1632]	1902	start(2)=jjb_u
	1903	start(3)=timestep
	1904	start(4)=0
	1905	count(2)=jjnb_u
	1906	count(3)=1
	1907	count(4)=0
	1908	IF (invert_y) start(2)=jjp1-jje_u+1
[5249]	1909	status = nf90_put_var(ncidps, varidps, psnat2, start, count)
[1632]	1910	IF (invert_y) THEN
	1911	! PRINT*,"Invertion impossible actuellement"
	1912	! CALL abort_gcm(modname,abort_message,1)
	1913	CALL invert_lat(iip1,jjnb_u,1,psnat2)
	1914	ENDIF
	1915	endif
	1916
	1917	END SUBROUTINE guide_read
	1918
	1919	!=======================================================================
	1920	SUBROUTINE guide_read2D(timestep)
[5249]	1921	USE netcdf, ONLY: nf90_put_var
[5271]	1922	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
[5272]	1923	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	1924	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
[5271]	1925	IMPLICIT NONE
	1926
[1632]	1927
[5272]	1928
[1632]	1929	INTEGER, INTENT(IN) :: timestep
	1930
	1931	LOGICAL, SAVE :: first=.TRUE.
	1932	! Identification fichiers et variables NetCDF:
	1933	INTEGER, SAVE :: ncidu,varidu,ncidv,varidv,ncidp,varidp
	1934	INTEGER, SAVE :: ncidQ,varidQ,ncidt,varidt,ncidps,varidps
	1935	INTEGER :: ncidpl,varidpl,varidap,varidbp
	1936	! Variables auxiliaires NetCDF:
	1937	INTEGER, DIMENSION(4) :: start,count
	1938	INTEGER :: status,rcode
	1939	! Variables for 3D extension:
	1940	REAL, DIMENSION (jjb_u:jje_u,llm) :: zu
	1941	REAL, DIMENSION (jjb_v:jje_v,llm) :: zv
	1942	INTEGER :: i
	1943	CHARACTER (len = 80) :: abort_message
[1906]	1944	CHARACTER (len = 20) :: modname = 'guide_read2D'
[1632]	1945	abort_message='pb in guide_read2D'
	1946
	1947	! -----------------------------------------------------------------
	1948	! Premier appel: initialisation de la lecture des fichiers
	1949	! -----------------------------------------------------------------
	1950	if (first) then
	1951	ncidpl=-99
[3995]	1952	write(,)trim(modname)//' : opening nudging files '
[1632]	1953	! Ap et Bp si niveaux de pression hybrides
[5084]	1954	if (guide_plevs.EQ.1) then
[3995]	1955	write(,)trim(modname)//' Reading nudging on model levels'
	1956	rcode = nf90_open('apbp.nc', nf90_nowrite, ncidpl)
[5270]	1957	IF (rcode.NE.nf90_noerr) THEN
[3995]	1958	abort_message='Nudging: error -> no file apbp.nc'
	1959	CALL abort_gcm(modname,abort_message,1)
	1960	ENDIF
	1961	rcode = nf90_inq_varid(ncidpl, 'AP', varidap)
[5270]	1962	IF (rcode.NE.nf90_noerr) THEN
[3995]	1963	abort_message='Nudging: error -> no AP variable in file apbp.nc'
	1964	CALL abort_gcm(modname,abort_message,1)
	1965	ENDIF
	1966	rcode = nf90_inq_varid(ncidpl, 'BP', varidbp)
[5270]	1967	IF (rcode.NE.nf90_noerr) THEN
[3995]	1968	abort_message='Nudging: error -> no BP variable in file apbp.nc'
	1969	CALL abort_gcm(modname,abort_message,1)
	1970	ENDIF
	1971	write(,)trim(modname)//'ncidpl,varidap',ncidpl,varidap
[1632]	1972	endif
	1973	! Pression
[5084]	1974	if (guide_plevs.EQ.2) then
[3995]	1975	rcode = nf90_open('P.nc', nf90_nowrite, ncidp)
[5270]	1976	IF (rcode.NE.nf90_noerr) THEN
[3995]	1977	abort_message='Nudging: error -> no file P.nc'
	1978	CALL abort_gcm(modname,abort_message,1)
	1979	ENDIF
	1980	rcode = nf90_inq_varid(ncidp, 'PRES', varidp)
[5270]	1981	IF (rcode.NE.nf90_noerr) THEN
[3995]	1982	abort_message='Nudging: error -> no PRES variable in file P.nc'
	1983	CALL abort_gcm(modname,abort_message,1)
	1984	ENDIF
	1985	write(,)trim(modname)//' ncidp,varidp',ncidp,varidp
[5084]	1986	if (ncidpl.eq.-99) ncidpl=ncidp
[1632]	1987	endif
	1988	! Vent zonal
	1989	if (guide_u) then
[3995]	1990	rcode = nf90_open('u.nc', nf90_nowrite, ncidu)
[5270]	1991	IF (rcode.NE.nf90_noerr) THEN
[3995]	1992	abort_message='Nudging: error -> no file u.nc'
	1993	CALL abort_gcm(modname,abort_message,1)
	1994	ENDIF
	1995	rcode = nf90_inq_varid(ncidu, 'UWND', varidu)
[5270]	1996	IF (rcode.NE.nf90_noerr) THEN
[3995]	1997	abort_message='Nudging: error -> no UWND variable in file u.nc'
	1998	CALL abort_gcm(modname,abort_message,1)
	1999	ENDIF
	2000	write(,)trim(modname)//' ncidu,varidu',ncidu,varidu
[5084]	2001	if (ncidpl.eq.-99) ncidpl=ncidu
[1632]	2002	endif
[2028]	2003
[1632]	2004	! Vent meridien
	2005	if (guide_v) then
[3995]	2006	rcode = nf90_open('v.nc', nf90_nowrite, ncidv)
[5270]	2007	IF (rcode.NE.nf90_noerr) THEN
[3995]	2008	abort_message='Nudging: error -> no file v.nc'
	2009	CALL abort_gcm(modname,abort_message,1)
	2010	ENDIF
	2011	rcode = nf90_inq_varid(ncidv, 'VWND', varidv)
[5270]	2012	IF (rcode.NE.nf90_noerr) THEN
[3995]	2013	abort_message='Nudging: error -> no VWND variable in file v.nc'
	2014	CALL abort_gcm(modname,abort_message,1)
	2015	ENDIF
	2016	write(,)trim(modname)//' ncidv,varidv',ncidv,varidv
[5084]	2017	if (ncidpl.eq.-99) ncidpl=ncidv
[3995]	2018	endif
[1632]	2019	! Temperature
	2020	if (guide_T) then
[3995]	2021	rcode = nf90_open('T.nc', nf90_nowrite, ncidt)
[5270]	2022	IF (rcode.NE.nf90_noerr) THEN
[3995]	2023	abort_message='Nudging: error -> no file T.nc'
	2024	CALL abort_gcm(modname,abort_message,1)
	2025	ENDIF
	2026	rcode = nf90_inq_varid(ncidt, 'AIR', varidt)
[5270]	2027	IF (rcode.NE.nf90_noerr) THEN
[3995]	2028	abort_message='Nudging: error -> no AIR variable in file T.nc'
	2029	CALL abort_gcm(modname,abort_message,1)
	2030	ENDIF
	2031	write(,)trim(modname)//' ncidT,varidT',ncidt,varidt
[5084]	2032	if (ncidpl.eq.-99) ncidpl=ncidt
[1632]	2033	endif
	2034	! Humidite
	2035	if (guide_Q) then
[3995]	2036	rcode = nf90_open('hur.nc', nf90_nowrite, ncidQ)
[5270]	2037	IF (rcode.NE.nf90_noerr) THEN
[3995]	2038	abort_message='Nudging: error -> no file hur.nc'
	2039	CALL abort_gcm(modname,abort_message,1)
	2040	ENDIF
	2041	rcode = nf90_inq_varid(ncidQ, 'RH', varidQ)
[5270]	2042	IF (rcode.NE.nf90_noerr) THEN
[3995]	2043	abort_message='Nudging: error -> no RH,variable in file hur.nc'
	2044	CALL abort_gcm(modname,abort_message,1)
	2045	ENDIF
	2046	write(,)trim(modname)//' ncidQ,varidQ',ncidQ,varidQ
[5084]	2047	if (ncidpl.eq.-99) ncidpl=ncidQ
[1632]	2048	endif
	2049	! Pression de surface
[5084]	2050	if ((guide_P).OR.(guide_plevs.EQ.1)) then
[3995]	2051	rcode = nf90_open('ps.nc', nf90_nowrite, ncidps)
[5270]	2052	IF (rcode.NE.nf90_noerr) THEN
[3995]	2053	abort_message='Nudging: error -> no file ps.nc'
	2054	CALL abort_gcm(modname,abort_message,1)
	2055	ENDIF
	2056	rcode = nf90_inq_varid(ncidps, 'SP', varidps)
[5270]	2057	IF (rcode.NE.nf90_noerr) THEN
[3995]	2058	abort_message='Nudging: error -> no SP variable in file ps.nc'
	2059	CALL abort_gcm(modname,abort_message,1)
	2060	ENDIF
	2061	write(,)trim(modname)//' ncidps,varidps',ncidps,varidps
[1632]	2062	endif
	2063	! Coordonnee verticale
[5084]	2064	if (guide_plevs.EQ.0) then
[3995]	2065	rcode = nf90_inq_varid(ncidpl, 'LEVEL', varidpl)
[5084]	2066	IF (rcode.NE.0) rcode = nf90_inq_varid(ncidpl, 'PRESSURE', varidpl)
[3995]	2067	write(,)trim(modname)//' ncidpl,varidpl',ncidpl,varidpl
[1632]	2068	endif
	2069	! Coefs ap, bp pour calcul de la pression aux differents niveaux
[5084]	2070	if (guide_plevs.EQ.1) then
[5249]	2071	status = nf90_put_var(ncidpl, varidap, apnc, [1], [nlevnc])
	2072	status = nf90_put_var(ncidpl, varidbp, bpnc, [1], [nlevnc])
[5084]	2073	elseif (guide_plevs.EQ.0) THEN
[5249]	2074	status = nf90_put_var(ncidpl, varidpl, apnc, [1], [nlevnc])
[1632]	2075	apnc=apnc*100.! conversion en Pascals
	2076	bpnc(:)=0.
	2077	endif
	2078	first=.FALSE.
	2079	endif ! (first)
	2080
	2081	! -----------------------------------------------------------------
	2082	! lecture des champs u, v, T, Q, ps
	2083	! -----------------------------------------------------------------
	2084
	2085	! dimensions pour les champs scalaires et le vent zonal
	2086	start(1)=1
	2087	start(2)=jjb_u
	2088	start(3)=1
	2089	start(4)=timestep
	2090
	2091	count(1)=1
	2092	count(2)=jjnb_u
	2093	count(3)=nlevnc
	2094	count(4)=1
	2095
	2096	IF (invert_y) start(2)=jjp1-jje_u+1
	2097	! Pression
[5084]	2098	if (guide_plevs.EQ.2) then
[5249]	2099	status = nf90_put_var(ncidp, varidp, zu, start, count)
[1632]	2100	DO i=1,iip1
	2101	pnat2(i,:,:)=zu(:,:)
	2102	ENDDO
	2103
	2104	IF (invert_y) THEN
	2105	! PRINT*,"Invertion impossible actuellement"
	2106	! CALL abort_gcm(modname,abort_message,1)
	2107	CALL invert_lat(iip1,jjnb_u,nlevnc,pnat2)
	2108	ENDIF
	2109	endif
	2110	! Vent zonal
	2111	if (guide_u) then
[5249]	2112	status = nf90_put_var(ncidu, varidu, zu, start, count)
[1632]	2113	DO i=1,iip1
	2114	unat2(i,:,:)=zu(:,:)
	2115	ENDDO
	2116
	2117	IF (invert_y) THEN
	2118	! PRINT*,"Invertion impossible actuellement"
	2119	! CALL abort_gcm(modname,abort_message,1)
	2120	CALL invert_lat(iip1,jjnb_u,nlevnc,unat2)
	2121	ENDIF
	2122	endif
	2123
[2028]	2124
[1632]	2125	! Temperature
	2126	if (guide_T) then
[5249]	2127	status = nf90_put_var(ncidt, varidt, zu, start, count)
[1632]	2128	DO i=1,iip1
	2129	tnat2(i,:,:)=zu(:,:)
	2130	ENDDO
	2131
	2132	IF (invert_y) THEN
	2133	! PRINT*,"Invertion impossible actuellement"
	2134	! CALL abort_gcm(modname,abort_message,1)
	2135	CALL invert_lat(iip1,jjnb_u,nlevnc,tnat2)
	2136	ENDIF
	2137	endif
	2138
	2139	! Humidite
	2140	if (guide_Q) then
[5249]	2141	status = nf90_put_var(ncidQ, varidQ, zu, start, count)
[1632]	2142	DO i=1,iip1
	2143	qnat2(i,:,:)=zu(:,:)
	2144	ENDDO
[5249]	2145
[1632]	2146	IF (invert_y) THEN
	2147	! PRINT*,"Invertion impossible actuellement"
	2148	! CALL abort_gcm(modname,abort_message,1)
	2149	CALL invert_lat(iip1,jjnb_u,nlevnc,qnat2)
	2150	ENDIF
	2151	endif
	2152
	2153	! Vent meridien
	2154	if (guide_v) then
	2155	start(2)=jjb_v
	2156	count(2)=jjnb_v
	2157	IF (invert_y) start(2)=jjm-jje_v+1
[5249]	2158	status = nf90_put_var(ncidv, varidv, zv, start, count)
[1632]	2159	DO i=1,iip1
	2160	vnat2(i,:,:)=zv(:,:)
	2161	ENDDO
	2162
	2163	IF (invert_y) THEN
[5249]	2164
[1632]	2165	! PRINT*,"Invertion impossible actuellement"
	2166	! CALL abort_gcm(modname,abort_message,1)
	2167	CALL invert_lat(iip1,jjnb_v,nlevnc,vnat2)
	2168	ENDIF
	2169	endif
	2170
	2171	! Pression de surface
[5084]	2172	if ((guide_P).OR.(guide_plevs.EQ.1)) then
[1632]	2173	start(2)=jjb_u
	2174	start(3)=timestep
	2175	start(4)=0
	2176	count(2)=jjnb_u
	2177	count(3)=1
	2178	count(4)=0
	2179	IF (invert_y) start(2)=jjp1-jje_u+1
[5249]	2180	status = nf90_put_var(ncidps, varidps, zu(:, 1), start, count)
[1632]	2181	DO i=1,iip1
	2182	psnat2(i,:)=zu(:,1)
	2183	ENDDO
	2184
	2185	IF (invert_y) THEN
	2186	! PRINT*,"Invertion impossible actuellement"
	2187	! CALL abort_gcm(modname,abort_message,1)
	2188	CALL invert_lat(iip1,jjnb_u,1,psnat2)
	2189	ENDIF
	2190	endif
	2191
	2192	END SUBROUTINE guide_read2D
[5249]	2193
[1632]	2194	!=======================================================================
[2028]	2195	SUBROUTINE guide_out(varname,hsize,vsize,field_loc,factt)
[1823]	2196	USE parallel_lmdz
[2028]	2197	USE mod_hallo, ONLY : gather_field_u, gather_field_v
[2597]	2198	USE comconst_mod, ONLY: pi
[2600]	2199	USE comvert_mod, ONLY: presnivs
[2740]	2200	use netcdf95, only: nf95_def_var, nf95_put_var
[5249]	2201	use netcdf, only: nf90_float, nf90_put_var
[2740]	2202
[5271]	2203	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
[5272]	2204	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	2205	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
[5271]	2206	IMPLICIT NONE
[1632]	2207
[5271]	2208
[5272]	2209
[1632]	2210	INCLUDE "comgeom2.h"
[5249]	2211
[1632]	2212	! Variables entree
[2028]	2213	CHARACTER(), INTENT(IN) :: varname
[1632]	2214	INTEGER, INTENT (IN) :: hsize,vsize
[2036]	2215	! REAL, DIMENSION (iip1,hsize,vsize), INTENT(IN) :: field_loc
	2216	REAL, DIMENSION (:,:), INTENT(IN) :: field_loc
	2217	REAL factt
[1632]	2218
	2219	! Variables locales
	2220	INTEGER, SAVE :: timestep=0
	2221	! Identites fichier netcdf
	2222	INTEGER :: nid, id_lonu, id_lonv, id_latu, id_latv, id_tim, id_lev
	2223	INTEGER :: vid_lonu,vid_lonv,vid_latu,vid_latv,vid_cu,vid_cv,vid_lev
[2740]	2224	INTEGER :: vid_au,vid_av, varid_alpha_t, varid_alpha_q
[1632]	2225	INTEGER, DIMENSION (3) :: dim3
	2226	INTEGER, DIMENSION (4) :: dim4,count,start
[2028]	2227	INTEGER :: ierr, varid,l
[2740]	2228	REAL zu(ip1jmp1),zv(ip1jm), zt(iip1, jjp1), zq(iip1, jjp1)
[2036]	2229	REAL, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: field_glo
[3995]	2230	CHARACTER(LEN=20),PARAMETER :: modname="guide_out"
[5249]	2231
[2028]	2232	!$OMP MASTER
[2036]	2233	ALLOCATE(field_glo(iip1,hsize,vsize))
	2234	!$OMP END MASTER
	2235	!$OMP BARRIER
[2028]	2236
[3995]	2237	! write(,)trim(modname)//' after allocation ',hsize,vsize
[2036]	2238
[2028]	2239	IF (hsize==jjp1) THEN
[2036]	2240	CALL gather_field_u(field_loc,field_glo,vsize)
[2028]	2241	ELSE IF (hsize==jjm) THEN
	2242	CALL gather_field_v(field_loc,field_glo, vsize)
	2243	ENDIF
	2244
[3995]	2245	! write(,)trim(modname)//' after gather '
[2028]	2246	CALL Gather_field_u(alpha_u,zu,1)
[2740]	2247	CALL Gather_field_u(alpha_t,zt,1)
	2248	CALL Gather_field_u(alpha_q,zq,1)
[2028]	2249	CALL Gather_field_v(alpha_v,zv,1)
[2036]	2250
	2251	IF (mpi_rank > 0) THEN
	2252	!$OMP MASTER
	2253	DEALLOCATE(field_glo)
[2034]	2254	!$OMP END MASTER
	2255	!$OMP BARRIER
[2028]	2256
[2036]	2257	RETURN
	2258	ENDIF
[5249]	2259
[2034]	2260	!$OMP MASTER
[5249]	2261	IF (timestep.EQ.0) THEN
[1632]	2262	! ----------------------------------------------
	2263	! initialisation fichier de sortie
	2264	! ----------------------------------------------
	2265	! Ouverture du fichier
[5270]	2266	ierr=nf90_create("guide_ins.nc",IOR(nf90_clobber,nf90_64bit_offset),nid)
[1632]	2267	! Definition des dimensions
[5270]	2268	ierr=nf90_def_dim(nid,"LONU",iip1,id_lonu)
	2269	ierr=nf90_def_dim(nid,"LONV",iip1,id_lonv)
	2270	ierr=nf90_def_dim(nid,"LATU",jjp1,id_latu)
	2271	ierr=nf90_def_dim(nid,"LATV",jjm,id_latv)
	2272	ierr=nf90_def_dim(nid,"LEVEL",llm,id_lev)
	2273	ierr=nf90_def_dim(nid,"TIME",nf90_unlimited,id_tim)
[1632]	2274
	2275	! Creation des variables dimensions
[5270]	2276	ierr=nf90_def_var(nid,"LONU",nf90_float,id_lonu,vid_lonu)
	2277	ierr=nf90_def_var(nid,"LONV",nf90_float,id_lonv,vid_lonv)
	2278	ierr=nf90_def_var(nid,"LATU",nf90_float,id_latu,vid_latu)
	2279	ierr=nf90_def_var(nid,"LATV",nf90_float,id_latv,vid_latv)
	2280	ierr=nf90_def_var(nid,"LEVEL",nf90_float,id_lev,vid_lev)
	2281	ierr=nf90_def_var(nid,"cu",nf90_float,(/id_lonu,id_latu/),vid_cu)
	2282	ierr=nf90_def_var(nid,"cv",nf90_float,(/id_lonv,id_latv/),vid_cv)
	2283	ierr=nf90_def_var(nid,"au",nf90_float,(/id_lonu,id_latu/),vid_au)
	2284	ierr=nf90_def_var(nid,"av",nf90_float,(/id_lonv,id_latv/),vid_av)
[2740]	2285	call nf95_def_var(nid, "alpha_T", nf90_float, (/id_lonv, id_latu/), &
	2286	varid_alpha_t)
	2287	call nf95_def_var(nid, "alpha_q", nf90_float, (/id_lonv, id_latu/), &
	2288	varid_alpha_q)
[5249]	2289
[5270]	2290	ierr=nf90_enddef(nid)
[1632]	2291
	2292	! Enregistrement des variables dimensions
[5249]	2293	ierr = nf90_put_var(nid, vid_lonu, rlonu * 180. / pi)
	2294	ierr = nf90_put_var(nid, vid_lonv, rlonv * 180. / pi)
	2295	ierr = nf90_put_var(nid, vid_latu, rlatu * 180. / pi)
	2296	ierr = nf90_put_var(nid, vid_latv, rlatv * 180. / pi)
	2297	ierr = nf90_put_var(nid, vid_lev, presnivs)
	2298	ierr = nf90_put_var(nid, vid_cu, cu)
	2299	ierr = nf90_put_var(nid, vid_cv, cv)
	2300	ierr = nf90_put_var(nid, vid_au, zu)
	2301	ierr = nf90_put_var(nid, vid_av, zv)
[2740]	2302	call nf95_put_var(nid, varid_alpha_t, zt)
	2303	call nf95_put_var(nid, varid_alpha_q, zq)
[1632]	2304	! --------------------------------------------------------------------
	2305	! Cr�ation des variables sauvegard�es
	2306	! --------------------------------------------------------------------
[5270]	2307	ierr = nf90_redef(nid)
[1632]	2308	! Pressure (GCM)
	2309	dim4=(/id_lonv,id_latu,id_lev,id_tim/)
[5270]	2310	ierr = nf90_def_var(nid,"SP",nf90_float,dim4,varid)
[1632]	2311	! Surface pressure (guidage)
	2312	IF (guide_P) THEN
	2313	dim3=(/id_lonv,id_latu,id_tim/)
[5270]	2314	ierr = nf90_def_var(nid,"ps",nf90_float,dim3,varid)
[1632]	2315	ENDIF
	2316	! Zonal wind
	2317	IF (guide_u) THEN
	2318	dim4=(/id_lonu,id_latu,id_lev,id_tim/)
[5270]	2319	ierr = nf90_def_var(nid,"u",nf90_float,dim4,varid)
	2320	ierr = nf90_def_var(nid,"ua",nf90_float,dim4,varid)
	2321	ierr = nf90_def_var(nid,"ucov",nf90_float,dim4,varid)
[1632]	2322	ENDIF
	2323	! Merid. wind
	2324	IF (guide_v) THEN
	2325	dim4=(/id_lonv,id_latv,id_lev,id_tim/)
[5270]	2326	ierr = nf90_def_var(nid,"v",nf90_float,dim4,varid)
	2327	ierr = nf90_def_var(nid,"va",nf90_float,dim4,varid)
	2328	ierr = nf90_def_var(nid,"vcov",nf90_float,dim4,varid)
[1632]	2329	ENDIF
	2330	! Pot. Temperature
	2331	IF (guide_T) THEN
	2332	dim4=(/id_lonv,id_latu,id_lev,id_tim/)
[5270]	2333	ierr = nf90_def_var(nid,"teta",nf90_float,dim4,varid)
[1632]	2334	ENDIF
	2335	! Specific Humidity
	2336	IF (guide_Q) THEN
	2337	dim4=(/id_lonv,id_latu,id_lev,id_tim/)
[5270]	2338	ierr = nf90_def_var(nid,"q",nf90_float,dim4,varid)
[1632]	2339	ENDIF
[5249]	2340
[5270]	2341	ierr = nf90_enddef(nid)
	2342	ierr = nf90_close(nid)
[1632]	2343	ENDIF ! timestep=0
	2344
	2345	! --------------------------------------------------------------------
	2346	! Enregistrement du champ
	2347	! --------------------------------------------------------------------
[5249]	2348
[5270]	2349	ierr=nf90_open("guide_ins.nc",nf90_write,nid)
[1632]	2350
[2028]	2351	IF (varname=="SP") timestep=timestep+1
	2352
[5270]	2353	ierr = nf90_inq_varid(nid,varname,varid)
[1632]	2354	SELECT CASE (varname)
[2028]	2355	CASE ("SP","ps")
[1632]	2356	start=(/1,1,1,timestep/)
	2357	count=(/iip1,jjp1,llm,1/)
[2028]	2358	CASE ("v","va","vcov")
[1632]	2359	start=(/1,1,1,timestep/)
	2360	count=(/iip1,jjm,llm,1/)
[2028]	2361	CASE DEFAULT
[1632]	2362	start=(/1,1,1,timestep/)
	2363	count=(/iip1,jjp1,llm,1/)
[2028]	2364	END SELECT
	2365
[2036]	2366	!$OMP END MASTER
	2367	!$OMP BARRIER
	2368
[2028]	2369	SELECT CASE (varname)
	2370
	2371	CASE("u","ua")
[2036]	2372	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[5249]	2373	DO l=1,llm
[2039]	2374	field_glo(:,2:jjm,l)=field_glo(:,2:jjm,l)/cu(:,2:jjm)
	2375	field_glo(:,1,l)=0. ; field_glo(:,jjp1,l)=0.
	2376	ENDDO
[2028]	2377	CASE("v","va")
[2036]	2378	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[2039]	2379	DO l=1,llm
	2380	field_glo(:,:,l)=field_glo(:,:,l)/cv(:,:)
	2381	ENDDO
[2028]	2382	END SELECT
	2383
	2384	! if (varname=="ua") then
	2385	! call dump2d(iip1,jjp1,field_glo,'ua gui1 1ere couche ')
	2386	! call dump2d(iip1,jjp1,field_glo(:,:,llm),'ua gui1 llm ')
	2387	! endif
[2036]	2388
	2389	!$OMP MASTER
	2390
[5249]	2391	ierr = nf90_put_var(nid, varid, field_glo, start, count)
[5270]	2392	ierr = nf90_close(nid)
[1632]	2393
[2036]	2394	DEALLOCATE(field_glo)
[2028]	2395	!$OMP END MASTER
	2396	!$OMP BARRIER
	2397
[1632]	2398	END SUBROUTINE guide_out
[5249]	2399
	2400
[1632]	2401	!===========================================================================
	2402	subroutine correctbid(iim,nl,x)
	2403	integer iim,nl
	2404	real x(iim+1,nl)
	2405	integer i,l
	2406	real zz
	2407
	2408	do l=1,nl
	2409	do i=2,iim-1
[5084]	2410	if(abs(x(i,l)).gt.1.e10) then
[1632]	2411	zz=0.5*(x(i-1,l)+x(i+1,l))
	2412	print*,'correction ',i,l,x(i,l),zz
	2413	x(i,l)=zz
	2414	endif
	2415	enddo
	2416	enddo
	2417	return
	2418	end subroutine correctbid
	2419
[2028]	2420
	2421	!====================================================================
	2422	! Ascii debug output. Could be reactivated
	2423	!====================================================================
	2424
	2425	subroutine dump2du(var,varname)
	2426	use parallel_lmdz
	2427	use mod_hallo
[5271]	2428	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
[5272]	2429	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	2430	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
[2028]	2431	implicit none
[5271]	2432
[2028]	2433
[5272]	2434
[2028]	2435	CHARACTER (len=*) :: varname
	2436
	2437
	2438	real, dimension(ijb_u:ije_u) :: var
	2439
	2440	real, dimension(ip1jmp1) :: var_glob
	2441
	2442	RETURN
	2443
	2444	call barrier
	2445	CALL Gather_field_u(var,var_glob,1)
	2446	call barrier
	2447
	2448	if (mpi_rank==0) then
	2449	call dump2d(iip1,jjp1,var_glob,varname)
	2450	endif
	2451
	2452	call barrier
	2453
	2454	return
	2455	end subroutine dump2du
	2456
	2457	!====================================================================
	2458	! Ascii debug output. Could be reactivated
	2459	!====================================================================
	2460	subroutine dumpall
[5271]	2461	USE dimensions_mod, ONLY: iim, jjm, llm, ndm
[5272]	2462	USE paramet_mod_h, ONLY: iip1, iip2, iip3, jjp1, llmp1, llmp2, llmm1, kftd, ip1jm, ip1jmp1, &
	2463	ip1jmi1, ijp1llm, ijmllm, mvar, jcfil, jcfllm
[5271]	2464	implicit none
	2465
[5272]	2466
[2028]	2467	include "comgeom.h"
	2468	call barrier
	2469	call dump2du(alpha_u(ijb_u:ije_u),' alpha_u couche 1')
	2470	call dump2du(unat2(:,jjbu:jjeu,nlevnc),' unat2 couche nlevnc')
	2471	call dump2du(ugui1(ijb_u:ije_u,1)*sqrt(unscu2(ijb_u:ije_u)),' ugui1 couche 1')
	2472	return
	2473	end subroutine dumpall
	2474
[1632]	2475	!===========================================================================
	2476	END MODULE guide_loc_mod

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

Download in other formats:

Original Format