Context Navigation

guide_loc_mod.F90 @ 3836

Last change on this file since 3836 was 3803, checked in by Ehouarn Millour, 4 years ago
Make dynamics restart.nc file be in "64-bit offset" mode to enable generating large files. While at it also enforce this for limit.nc, guide_ins.nc and grilles_gcm.nc generated by the model. EM
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: 82.6 KB

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

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

Download in other formats:

Original Format