Context Navigation

startvar.F @ 1237

Last change on this file since 1237 was 1222, checked in by Ehouarn Millour, 15 years ago

Changes and cleanups to enable compiling without physics
and without ioipsl.

IOIPSL related cleanups:

bibio/writehist.F encapsulate the routine (which needs IOIPSL to function)

with #ifdef IOIPSL flag.

dyn3d/abort_gcm.F, dyn3dpar/abort_gcm.F and dyn3dpar/getparam.F90: use ioipsl_getincom module when not compiling with IOIPSL library, in order to always be able to use getin() routine.
removed unused "use IOIPSL" in dyn3dpar/guide_p_mod.F90
calendar related issue: Initialize day_ref and annee_ref in iniacademic.F (i.e. when they are not read from start.nc file)

Earth-specific programs/routines/modules:
create_etat0.F, fluxstokenc.F, limit_netcdf.F, startvar.F
(versions in dyn3d and dyn3dpar)
These routines and modules, which by design and porpose are made to function with
Earth physics are encapsulated with #CPP_EARTH cpp flag.

Earth-specific instructions:

calls to qminimum (specific treatment of first 2 tracers, i.e. water) in dyn3d/caladvtrac.F, dyn3d/integrd.F, dyn3dpar/caladvtrac_p.F, dyn3dpar/integrd_p.F only if (planet_type == 'earth')

Interaction with parallel physics:

routine dyn3dpar/parallel.F90 uses "surface_data" module (which is in the physics ...) to know value of "type_ocean" . Encapsulated that with #ifdef CPP_EARTH and set to a default type_ocean="dummy" otherwise.
So far, only Earth physics are parallelized, so all the interaction between parallel dynamics and parallel physics are encapsulated with #ifdef CCP_EARTH (this way we can run parallel without any physics). The (dyn3dpar) routines which contains such interaction are: bands.F90, gr_dyn_fi_p.F, gr_fi_dyn_p.F, mod_interface_dyn_phys.F90 This should later (when improving dyn/phys interface) be encapsulated with a more general and appropriate #ifdef CPP_PHYS cpp flag.

I checked that these changes do not alter results (on the simple
32x24x11 bench) on Ciclad (seq & mpi), Brodie (seq, mpi & omp) and
Vargas (seq, mpi & omp).

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 39.9 KB

Rev	Line
[630]	1	!
[1222]	2	! $Id: startvar.F 1222 2009-08-07 11:48:33Z fairhead $
[630]	3	!
	4	MODULE startvar
[1222]	5	#ifdef CPP_EARTH
	6	! This module is designed to work for Earth (and with ioipsl)
[630]	7	!
	8	!
	9	! There are three ways to access data from the database of atmospheric data which
	10	! can be used to initialize the model. This depends on the type of field which needs
	11	! to be extracted. In any case the call should come after a restget and should be of the type :
	12	! CALL startget(...)
	13	!
	14	! We will details the possible arguments to startget here :
	15	!
	16	! - A 2D variable on the dynamical grid :
	17	! CALL startget(varname, iml, jml, lon_in, lat_in, champ, val_ex, jml2, lon_in2, lat_in2, interbar )
	18	!
	19	! - A 1D variable on the physical grid :
	20	! CALL startget(varname, iml, jml, lon_in, lat_in, nbindex, champ, val_exp, jml2, lon_in2, lat_in2, interbar )
	21	!
	22	!
	23	! - A 3D variable on the dynamical grid :
	24	! CALL startget(varname, iml, jml, lon_in, lat_in, lml, pls, workvar, champ, val_exp, jml2, lon_in2, lat_in2, interbar )
	25	!
	26	!
	27	! There is special constraint on the atmospheric data base except that the
	28	! the data needs to be in netCDF and the variables should have the the following
	29	! names in the file :
	30	!
	31	! 'RELIEF' : High resolution orography
	32	! 'ST' : Surface temperature
	33	! 'CDSW' : Soil moisture
	34	! 'Z' : Surface geopotential
	35	! 'SP' : Surface pressure
	36	! 'U' : East ward wind
	37	! 'V' : Northward wind
	38	! 'TEMP' : Temperature
	39	! 'R' : Relative humidity
	40	!
	41	USE ioipsl
	42	!
	43	!
	44	IMPLICIT NONE
	45	!
	46	!
	47	PRIVATE
	48	PUBLIC startget
	49	!
	50	!
	51	INTERFACE startget
	52	MODULE PROCEDURE startget_phys2d, startget_phys1d, startget_dyn
	53	END INTERFACE
	54	!
	55	INTEGER, SAVE :: fid_phys, fid_dyn
	56	INTEGER, SAVE :: iml_phys, iml_rel, iml_dyn
	57	INTEGER, SAVE :: jml_phys, jml_rel, jml_dyn
	58	INTEGER, SAVE :: llm_dyn, ttm_dyn
	59	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: lon_phys, lon_rug,
	60	. lon_alb, lon_rel, lon_dyn
	61	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: lat_phys, lat_rug,
	62	. lat_alb, lat_rel, lat_dyn
	63	REAL, ALLOCATABLE, SAVE, DIMENSION (:) :: levdyn_ini
	64	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: relief, zstd, zsig,
	65	. zgam, zthe, zpic, zval
	66	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: rugo, masque, phis
	67	!
	68	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: tsol, qsol, psol_dyn
	69	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:,:) :: var_ana3d
	70	!
	71	CONTAINS
	72	!
	73	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	74	!
	75	SUBROUTINE startget_phys2d(varname, iml, jml, lon_in, lat_in,
	76	. champ, val_exp, jml2, lon_in2, lat_in2 , interbar, masque_lu )
	77	!
	78	! There is a big mess with the size in logitude, should it be iml or iml+1.
	79	! I have chosen to use the iml+1 as an argument to this routine and we declare
	80	! internaly smaler fields when needed. This needs to be cleared once and for all in LMDZ.
	81	! A convention is required.
	82	!
	83	!
	84	CHARACTER(), INTENT(in) :: varname
	85	INTEGER, INTENT(in) :: iml, jml ,jml2
	86	REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
	87	REAL, INTENT(in) :: lon_in2(iml), lat_in2(jml2)
	88	REAL, INTENT(inout) :: champ(iml,jml)
	89	REAL, INTENT(in) :: val_exp
	90	REAL, INTENT(in), optional :: masque_lu(iml,jml)
	91	LOGICAL interbar
	92	!
	93	! This routine only works if the variable does not exist or is constant
	94	!
	95	IF ( MINVAL(champ(:,:)).EQ.MAXVAL(champ(:,:)) .AND.
	96	.MINVAL(champ(:,:)).EQ.val_exp ) THEN
	97	!
	98	SELECTCASE(varname)
	99	!
	100	CASE ('relief')
	101	!
	102	! If we do not have the orography we need to get it
	103	!
	104	IF ( .NOT.ALLOCATED(relief)) THEN
	105	!
	106	if (present(masque_lu)) then
	107	CALL start_init_orog( iml, jml, lon_in, lat_in,
	108	. jml2,lon_in2,lat_in2, interbar, masque_lu )
	109	else
	110	CALL start_init_orog( iml, jml, lon_in, lat_in,
	111	. jml2,lon_in2,lat_in2, interbar)
	112	endif
	113	!
	114	ENDIF
	115	!
	116	IF (SIZE(relief) .NE. SIZE(champ)) THEN
	117	!
	118	WRITE(,) 'STARTVAR module has been',
	119	.' initialized to the wrong size'
	120	STOP
	121	!
	122	ENDIF
	123	!
	124	champ(:,:) = relief(:,:)
	125	!
	126	CASE ('rugosite')
	127	!
	128	! If we do not have the orography we need to get it
	129	!
	130	IF ( .NOT.ALLOCATED(rugo)) THEN
	131	!
	132	CALL start_init_orog( iml, jml, lon_in, lat_in,
	133	. jml2,lon_in2,lat_in2 , interbar )
	134	!
	135	ENDIF
	136	!
	137	IF (SIZE(rugo) .NE. SIZE(champ)) THEN
	138	!
	139	WRITE(,)
	140	. 'STARTVAR module has been initialized to the wrong size'
	141	STOP
	142	!
	143	ENDIF
	144	!
	145	champ(:,:) = rugo(:,:)
	146	!
	147	CASE ('masque')
	148	!
	149	! If we do not have the orography we need to get it
	150	!
	151	IF ( .NOT.ALLOCATED(masque)) THEN
	152	!
	153	CALL start_init_orog( iml, jml, lon_in, lat_in,
	154	. jml2,lon_in2,lat_in2 , interbar )
	155	!
	156	ENDIF
	157	!
	158	IF (SIZE(masque) .NE. SIZE(champ)) THEN
	159	!
	160	WRITE(,)
	161	. 'STARTVAR module has been initialized to the wrong size'
	162	STOP
	163	!
	164	ENDIF
	165	!
	166	champ(:,:) = masque(:,:)
	167	!
	168	CASE ('surfgeo')
	169	!
	170	! If we do not have the orography we need to get it
	171	!
	172	IF ( .NOT.ALLOCATED(phis)) THEN
	173	!
	174	CALL start_init_orog( iml, jml, lon_in, lat_in,
	175	. jml2,lon_in2, lat_in2 , interbar )
	176	!
	177	ENDIF
	178	!
	179	IF (SIZE(phis) .NE. SIZE(champ)) THEN
	180	!
	181	WRITE(,)
	182	. 'STARTVAR module has been initialized to the wrong size'
	183	STOP
	184	!
	185	ENDIF
	186	!
	187	champ(:,:) = phis(:,:)
	188	!
	189	CASE ('psol')
	190	!
	191	! If we do not have the orography we need to get it
	192	!
	193	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	194	!
	195	CALL start_init_dyn( iml, jml, lon_in, lat_in,
	196	. jml2,lon_in2, lat_in2 , interbar )
	197	!
	198	ENDIF
	199	!
	200	IF (SIZE(psol_dyn) .NE. SIZE(champ)) THEN
	201	!
	202	WRITE(,)
	203	. 'STARTVAR module has been initialized to the wrong size'
	204	STOP
	205	!
	206	ENDIF
	207	!
	208	champ(:,:) = psol_dyn(:,:)
	209	!
	210	CASE DEFAULT
	211	!
	212	WRITE(,) 'startget_phys2d'
	213	WRITE(,) 'No rule is present to extract variable',
	214	. varname(:LEN_TRIM(varname)),' from any data set'
	215	STOP
	216	!
	217	END SELECT
	218	!
	219	ELSE
	220	!
	221	! There are a few fields we might need if we need to interpolate 3D filed. Thus if they come through here we
	222	! will catch them
	223	!
	224	SELECTCASE(varname)
	225	!
	226	CASE ('surfgeo')
	227	!
	228	IF ( .NOT.ALLOCATED(phis)) THEN
	229	ALLOCATE(phis(iml,jml))
	230	ENDIF
	231	!
	232	IF (SIZE(phis) .NE. SIZE(champ)) THEN
	233	!
	234	WRITE(,)
	235	. 'STARTVAR module has been initialized to the wrong size'
	236	STOP
	237	!
	238	ENDIF
	239	!
	240	phis(:,:) = champ(:,:)
	241	!
	242	END SELECT
	243	!
	244	ENDIF
	245	!
	246	END SUBROUTINE startget_phys2d
	247	!
	248	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	249	!
	250	SUBROUTINE start_init_orog ( iml,jml,lon_in, lat_in,jml2,lon_in2 ,
	251	, lat_in2 , interbar, masque_lu )
	252	!
	253	INTEGER, INTENT(in) :: iml, jml, jml2
	254	REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
	255	REAL, INTENT(in) :: lon_in2(iml), lat_in2(jml2)
	256	REAL, intent(in), optional :: masque_lu(iml,jml)
	257	LOGICAL interbar
	258	!
	259	! LOCAL
	260	!
	261	LOGICAL interbar2
	262	REAL :: lev(1), date, dt,chmin,chmax
	263	INTEGER :: itau(1), fid
	264	INTEGER :: llm_tmp, ttm_tmp
	265	INTEGER :: i, j
	266	INTEGER :: iret
	267	CHARACTER*25 title
	268	REAL, ALLOCATABLE :: relief_hi(:,:)
	269	REAL, ALLOCATABLE :: lon_rad(:), lat_rad(:)
	270	REAL, ALLOCATABLE :: lon_ini(:), lat_ini(:)
	271	REAL, ALLOCATABLE :: tmp_var(:,:)
	272	INTEGER, ALLOCATABLE :: tmp_int(:,:)
	273	!
	274	CHARACTER*120 :: orogfname
	275	LOGICAL :: check=.TRUE.
	276	!
	277	!
	278	orogfname = 'Relief.nc'
	279	!
	280	IF ( check ) WRITE(,) 'Reading the high resolution orography'
	281	!
	282	CALL flininfo(orogfname,iml_rel, jml_rel, llm_tmp, ttm_tmp, fid)
	283	!
	284	ALLOCATE (lat_rel(iml_rel,jml_rel), stat=iret)
	285	ALLOCATE (lon_rel(iml_rel,jml_rel), stat=iret)
	286	ALLOCATE (relief_hi(iml_rel,jml_rel), stat=iret)
	287	!
	288	CALL flinopen(orogfname, .FALSE., iml_rel, jml_rel,
	289	.llm_tmp, lon_rel, lat_rel, lev, ttm_tmp,
	290	. itau, date, dt, fid)
	291	!
	292	CALL flinget(fid, 'RELIEF', iml_rel, jml_rel, llm_tmp,
	293	. ttm_tmp, 1, 1, relief_hi)
	294	!
	295	CALL flinclo(fid)
	296	!
	297	! In case we have a file which is in degrees we do the transformation
	298	!
	299	ALLOCATE(lon_rad(iml_rel))
	300	ALLOCATE(lon_ini(iml_rel))
	301
	302	IF ( MAXVAL(lon_rel(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	303	lon_ini(:) = lon_rel(:,1) * 2.0 * ASIN(1.0) / 180.0
	304	ELSE
	305	lon_ini(:) = lon_rel(:,1)
	306	ENDIF
	307
	308	ALLOCATE(lat_rad(jml_rel))
	309	ALLOCATE(lat_ini(jml_rel))
	310
	311	IF ( MAXVAL(lat_rel(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	312	lat_ini(:) = lat_rel(1,:) * 2.0 * ASIN(1.0) / 180.0
	313	ELSE
	314	lat_ini(:) = lat_rel(1,:)
	315	ENDIF
	316	!
	317	!
	318
	319	title='RELIEF'
	320
	321	interbar2 = .FALSE.
	322	CALL conf_dat2d(title,iml_rel, jml_rel, lon_ini, lat_ini,
	323	. lon_rad, lat_rad, relief_hi , interbar2 )
	324
	325	IF ( check ) WRITE(,) 'Computes all the parameters needed',
	326	.' for the gravity wave drag code'
	327	!
	328	! Allocate the data we need to put in the interpolated fields
	329	!
	330	! RELIEF: orographie moyenne
	331	ALLOCATE(relief(iml,jml))
	332	! zphi : orographie moyenne
	333	ALLOCATE(phis(iml,jml))
	334	! zstd: deviation standard de l'orographie sous-maille
	335	ALLOCATE(zstd(iml,jml))
	336	! zsig: pente de l'orographie sous-maille
	337	ALLOCATE(zsig(iml,jml))
	338	! zgam: anisotropy de l'orographie sous maille
	339	ALLOCATE(zgam(iml,jml))
	340	! zthe: orientation de l'axe oriente dans la direction
	341	! de plus grande pente de l'orographie sous maille
	342	ALLOCATE(zthe(iml,jml))
	343	! zpic: hauteur pics de la SSO
	344	ALLOCATE(zpic(iml,jml))
	345	! zval: hauteur vallees de la SSO
	346	ALLOCATE(zval(iml,jml))
	347	! masque : Masque terre ocean
	348	ALLOCATE(tmp_int(iml,jml))
	349	ALLOCATE(masque(iml,jml))
	350
	351	masque = -99999.
	352	if (present(masque_lu)) then
	353	masque = masque_lu
	354	endif
	355	!
	356	CALL grid_noro(iml_rel, jml_rel, lon_rad, lat_rad, relief_hi,
	357	. iml-1, jml, lon_in, lat_in,
	358	. phis, relief, zstd, zsig, zgam, zthe, zpic, zval, masque)
	359	phis = phis * 9.81
	360	!
	361	! masque(:,:) = FLOAT(tmp_int(:,:))
	362	!
	363	! Compute surface roughness
	364	!
	365	IF ( check ) WRITE(,)
	366	.'Compute surface roughness induced by the orography'
	367	!
	368	ALLOCATE(rugo(iml,jml))
	369	ALLOCATE(tmp_var(iml-1,jml))
	370	!
	371	CALL rugsoro(iml_rel, jml_rel, lon_rad, lat_rad, relief_hi,
	372	. iml-1, jml, lon_in, lat_in, tmp_var)
	373	!
	374	DO j = 1, jml
	375	DO i = 1, iml-1
	376	rugo(i,j) = tmp_var(i,j)
	377	ENDDO
	378	rugo(iml,j) = tmp_var(1,j)
	379	ENDDO
	380	c
	381	cc * rugo n'est pas utilise pour l'instant ****
	382	!
	383	! Build land-sea mask
	384	!
	385	!
	386	RETURN
	387	!
	388	END SUBROUTINE start_init_orog
	389	!
	390	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	391	!
	392	SUBROUTINE startget_phys1d(varname, iml, jml, lon_in,
	393	.lat_in, nbindex, champ, val_exp ,jml2, lon_in2, lat_in2,interbar)
	394	!
	395	CHARACTER(), INTENT(in) :: varname
	396	INTEGER, INTENT(in) :: iml, jml, nbindex, jml2
	397	REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
	398	REAL, INTENT(in) :: lon_in2(iml), lat_in2(jml2)
	399	REAL, INTENT(inout) :: champ(nbindex)
	400	REAL, INTENT(in) :: val_exp
	401	LOGICAL interbar
	402	!
	403	!
	404	! This routine only works if the variable does not exist or is constant
	405	!
	406	IF ( MINVAL(champ(:)).EQ.MAXVAL(champ(:)) .AND.
	407	.MINVAL(champ(:)).EQ.val_exp ) THEN
	408	SELECTCASE(varname)
	409	CASE ('tsol')
	410	IF ( .NOT.ALLOCATED(tsol)) THEN
	411	CALL start_init_phys( iml, jml, lon_in, lat_in,
	412	. jml2, lon_in2, lat_in2, interbar )
	413	ENDIF
	414	IF ( SIZE(tsol) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	415	WRITE(,)
	416	. 'STARTVAR module has been initialized to the wrong size'
	417	STOP
	418	ENDIF
	419	CALL gr_dyn_fi(1, iml, jml, nbindex, tsol, champ)
	420	CASE ('qsol')
	421	IF ( .NOT.ALLOCATED(qsol)) THEN
	422	CALL start_init_phys( iml, jml, lon_in, lat_in,
	423	. jml2, lon_in2,lat_in2 , interbar )
	424	ENDIF
	425	IF ( SIZE(qsol) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	426	WRITE(,)
	427	. 'STARTVAR module has been initialized to the wrong size'
	428	STOP
	429	ENDIF
	430	CALL gr_dyn_fi(1, iml, jml, nbindex, qsol, champ)
	431	CASE ('psol')
	432	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	433	CALL start_init_dyn( iml, jml, lon_in, lat_in,
	434	. jml2, lon_in2,lat_in2 , interbar )
	435	ENDIF
	436	IF (SIZE(psol_dyn) .NE. SIZE(lon_in)*SIZE(lat_in)) THEN
	437	WRITE(,)
	438	. 'STARTVAR module has been initialized to the wrong size'
	439	STOP
	440	ENDIF
	441	CALL gr_dyn_fi(1, iml, jml, nbindex, psol_dyn, champ)
	442	CASE ('zmea')
	443	IF ( .NOT.ALLOCATED(relief)) THEN
	444	CALL start_init_orog( iml, jml, lon_in, lat_in,
	445	. jml2, lon_in2,lat_in2 , interbar )
	446	ENDIF
	447	IF ( SIZE(relief) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	448	WRITE(,)
	449	. 'STARTVAR module has been initialized to the wrong size'
	450	STOP
	451	ENDIF
	452	CALL gr_dyn_fi(1, iml, jml, nbindex, relief, champ)
	453	CASE ('zstd')
	454	IF ( .NOT.ALLOCATED(zstd)) THEN
	455	CALL start_init_orog( iml, jml, lon_in, lat_in,
	456	. jml2, lon_in2,lat_in2 , interbar )
	457	ENDIF
	458	IF ( SIZE(zstd) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	459	WRITE(,)
	460	. 'STARTVAR module has been initialized to the wrong size'
	461	STOP
	462	ENDIF
	463	CALL gr_dyn_fi(1, iml, jml, nbindex,zstd, champ)
	464	CASE ('zsig')
	465	IF ( .NOT.ALLOCATED(zsig)) THEN
	466	CALL start_init_orog( iml, jml, lon_in, lat_in,
	467	. jml2, lon_in2,lat_in2 , interbar )
	468	ENDIF
	469	IF ( SIZE(zsig) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	470	WRITE(,)
	471	. 'STARTVAR module has been initialized to the wrong size'
	472	STOP
	473	ENDIF
	474	CALL gr_dyn_fi(1, iml, jml, nbindex,zsig, champ)
	475	CASE ('zgam')
	476	IF ( .NOT.ALLOCATED(zgam)) THEN
	477	CALL start_init_orog( iml, jml, lon_in, lat_in,
	478	. jml2, lon_in2,lat_in2 , interbar )
	479	ENDIF
	480	IF ( SIZE(zgam) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	481	WRITE(,)
	482	. 'STARTVAR module has been initialized to the wrong size'
	483	STOP
	484	ENDIF
	485	CALL gr_dyn_fi(1, iml, jml, nbindex,zgam, champ)
	486	CASE ('zthe')
	487	IF ( .NOT.ALLOCATED(zthe)) THEN
	488	CALL start_init_orog( iml, jml, lon_in, lat_in,
	489	. jml2, lon_in2,lat_in2 , interbar )
	490	ENDIF
	491	IF ( SIZE(zthe) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	492	WRITE(,)
	493	. 'STARTVAR module has been initialized to the wrong size'
	494	STOP
	495	ENDIF
	496	CALL gr_dyn_fi(1, iml, jml, nbindex,zthe, champ)
	497	CASE ('zpic')
	498	IF ( .NOT.ALLOCATED(zpic)) THEN
	499	CALL start_init_orog( iml, jml, lon_in, lat_in,
	500	. jml2, lon_in2,lat_in2 , interbar )
	501	ENDIF
	502	IF ( SIZE(zpic) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	503	WRITE(,)
	504	. 'STARTVAR module has been initialized to the wrong size'
	505	STOP
	506	ENDIF
	507	CALL gr_dyn_fi(1, iml, jml, nbindex,zpic, champ)
	508	CASE ('zval')
	509	IF ( .NOT.ALLOCATED(zval)) THEN
	510	CALL start_init_orog( iml, jml, lon_in, lat_in,
	511	. jml2, lon_in2,lat_in2 , interbar )
	512	ENDIF
	513	IF ( SIZE(zval) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	514	WRITE(,)
	515	. 'STARTVAR module has been initialized to the wrong size'
	516	STOP
	517	ENDIF
	518	CALL gr_dyn_fi(1, iml, jml, nbindex,zval, champ)
	519	CASE ('rads')
	520	champ(:) = 0.0
	521	CASE ('snow')
	522	champ(:) = 0.0
[1000]	523	cIM "slab" ocean
	524	CASE ('tslab')
	525	champ(:) = 0.0
	526	CASE ('seaice')
[630]	527	champ(:) = 0.0
	528	CASE ('rugmer')
	529	champ(:) = 0.001
	530	CASE ('agsno')
	531	champ(:) = 50.0
	532	CASE DEFAULT
	533	WRITE(,) 'startget_phys1d'
	534	WRITE(,) 'No rule is present to extract variable ',
	535	. varname(:LEN_TRIM(varname)),' from any data set'
	536	STOP
	537	END SELECT
	538	ELSE
	539	!
	540	! If we see tsol we catch it as we may need it for a 3D interpolation
	541	!
	542	SELECTCASE(varname)
	543	CASE ('tsol')
	544	IF ( .NOT.ALLOCATED(tsol)) THEN
	545	ALLOCATE(tsol(SIZE(lon_in),SIZE(lat_in) ))
	546	ENDIF
	547	CALL gr_fi_dyn(1, iml, jml, nbindex, champ, tsol)
	548	END SELECT
	549	ENDIF
	550	END SUBROUTINE startget_phys1d
	551	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	552	!
	553	SUBROUTINE start_init_phys( iml, jml, lon_in, lat_in, jml2,
	554	. lon_in2, lat_in2 , interbar )
	555	!
	556	INTEGER, INTENT(in) :: iml, jml ,jml2
	557	REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
	558	REAL, INTENT(in) :: lon_in2(iml), lat_in2(jml2)
	559	LOGICAL interbar
	560	!
	561	! LOCAL
	562	!
	563	!ac REAL :: lev(1), date, dt
	564	REAL :: date, dt
	565	REAL, DIMENSION(:), ALLOCATABLE :: levphys_ini
	566	!ac
	567	INTEGER :: itau(1)
	568	INTEGER :: llm_tmp, ttm_tmp
	569	INTEGER :: i, j
	570	!
	571	CHARACTER*25 title
	572	CHARACTER*120 :: physfname
	573	LOGICAL :: check=.TRUE.
	574	!
	575	REAL, ALLOCATABLE :: lon_rad(:), lat_rad(:)
	576	REAL, ALLOCATABLE :: lon_ini(:), lat_ini(:)
	577	REAL, ALLOCATABLE :: var_ana(:,:), tmp_var(:,:)
	578	!
	579	physfname = 'ECPHY.nc'
	580	!
	581	IF ( check ) WRITE(,) 'Opening the surface analysis'
	582	!
	583	CALL flininfo(physfname, iml_phys, jml_phys, llm_tmp,
	584	. ttm_tmp, fid_phys)
	585	!
	586	ALLOCATE (lat_phys(iml_phys,jml_phys))
	587	ALLOCATE (lon_phys(iml_phys,jml_phys))
	588	!ac
	589	ALLOCATE (levphys_ini(llm_tmp))
	590	!
	591	! CALL flinopen(physfname, .FALSE., iml_phys, jml_phys,
	592	! . llm_tmp, lon_phys, lat_phys, lev, ttm_tmp,
	593	! . itau, date, dt, fid_phys)
	594	!
	595	CALL flinopen(physfname, .FALSE., iml_phys, jml_phys,
	596	. llm_tmp, lon_phys, lat_phys, levphys_ini, ttm_tmp,
	597	. itau, date, dt, fid_phys)
	598	!
	599	DEALLOCATE (levphys_ini)
	600	!ac
	601	!
	602	! Allocate the space we will need to get the data out of this file
	603	!
	604	ALLOCATE(var_ana(iml_phys, jml_phys))
	605	!
	606	! In case we have a file which is in degrees we do the transformation
	607	!
	608	ALLOCATE(lon_rad(iml_phys))
	609	ALLOCATE(lon_ini(iml_phys))
	610
	611	IF ( MAXVAL(lon_phys(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	612	lon_ini(:) = lon_phys(:,1) * 2.0 * ASIN(1.0) / 180.0
	613	ELSE
	614	lon_ini(:) = lon_phys(:,1)
	615	ENDIF
	616
	617	ALLOCATE(lat_rad(jml_phys))
	618	ALLOCATE(lat_ini(jml_phys))
	619
	620	IF ( MAXVAL(lat_phys(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	621	lat_ini(:) = lat_phys(1,:) * 2.0 * ASIN(1.0) / 180.0
	622	ELSE
	623	lat_ini(:) = lat_phys(1,:)
	624	ENDIF
	625
	626
	627	!
	628	! We get the two standard varibales
	629	! Surface temperature
	630	!
	631	ALLOCATE(tsol(iml,jml))
	632	ALLOCATE(tmp_var(iml-1,jml))
	633	!
	634	!
	635
	636	CALL flinget(fid_phys, 'ST', iml_phys, jml_phys,
	637	.llm_tmp, ttm_tmp, 1, 1, var_ana)
	638
	639	title='ST'
	640	CALL conf_dat2d(title,iml_phys, jml_phys, lon_ini, lat_ini,
	641	. lon_rad, lat_rad, var_ana , interbar )
	642
	643	IF ( interbar ) THEN
	644	WRITE(6,*) '-------------------------------------------------',
	645	,'--------------'
	646	WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
	647	, ' pour ST $$$ '
	648	WRITE(6,*) '-------------------------------------------------',
	649	,'--------------'
	650	CALL inter_barxy ( iml_phys,jml_phys -1,lon_rad,lat_rad ,
	651	, var_ana, iml-1, jml-1, lon_in2, lat_in2, jml, tmp_var )
	652	ELSE
	653	CALL grille_m(iml_phys, jml_phys, lon_rad, lat_rad,
	654	. var_ana, iml-1, jml, lon_in, lat_in, tmp_var )
	655	ENDIF
	656
	657	CALL gr_int_dyn(tmp_var, tsol, iml-1, jml)
	658	!
	659	! Soil moisture
	660	!
	661	ALLOCATE(qsol(iml,jml))
	662	CALL flinget(fid_phys, 'CDSW', iml_phys, jml_phys,
	663	. llm_tmp, ttm_tmp, 1, 1, var_ana)
	664
	665	title='CDSW'
	666	CALL conf_dat2d(title,iml_phys, jml_phys, lon_ini, lat_ini,
	667	. lon_rad, lat_rad, var_ana, interbar )
	668
	669	IF ( interbar ) THEN
	670	WRITE(6,*) '-------------------------------------------------',
	671	,'--------------'
	672	WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
	673	, ' pour CDSW $$$ '
	674	WRITE(6,*) '-------------------------------------------------',
	675	,'--------------'
	676	CALL inter_barxy ( iml_phys,jml_phys -1,lon_rad,lat_rad ,
	677	, var_ana, iml-1, jml-1, lon_in2, lat_in2, jml, tmp_var )
	678	ELSE
	679	CALL grille_m(iml_phys, jml_phys, lon_rad, lat_rad,
	680	. var_ana, iml-1, jml, lon_in, lat_in, tmp_var )
	681	ENDIF
	682	c
	683	CALL gr_int_dyn(tmp_var, qsol, iml-1, jml)
	684	!
	685	CALL flinclo(fid_phys)
	686	!
	687	END SUBROUTINE start_init_phys
	688	!
	689	!
	690	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	691	!
	692	!
	693	SUBROUTINE startget_dyn(varname, iml, jml, lon_in, lat_in,
	694	. lml, pls, workvar, champ, val_exp,jml2, lon_in2, lat_in2 ,
	695	, interbar )
	696	!
	697	! ARGUMENTS
	698	!
	699	CHARACTER(), INTENT(in) :: varname
	700	INTEGER, INTENT(in) :: iml, jml, lml, jml2
	701	REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
	702	REAL, INTENT(in) :: lon_in2(iml), lat_in2(jml2)
	703	REAL, INTENT(in) :: pls(iml, jml, lml)
	704	REAL, INTENT(in) :: workvar(iml, jml, lml)
	705	REAL, INTENT(inout) :: champ(iml, jml, lml)
	706	REAL, INTENT(in) :: val_exp
	707	LOGICAL interbar
	708	!
	709	! LOCAL
	710	!
	711	INTEGER :: il, ij, ii
	712	REAL :: xppn, xpps
	713	!
	714	#include "dimensions.h"
	715	#include "paramet.h"
	716	#include "comgeom2.h"
	717	#include "comconst.h"
	718	!
	719	! This routine only works if the variable does not exist or is constant
	720	!
	721	IF ( MINVAL(champ(:,:,:)).EQ.MAXVAL(champ(:,:,:)) .AND.
	722	. MINVAL(champ(:,:,:)).EQ.val_exp ) THEN
	723	!
	724	SELECTCASE(varname)
	725	CASE ('u')
	726	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	727	CALL start_init_dyn( iml, jml, lon_in, lat_in, jml2 ,
	728	. lon_in2,lat_in2 , interbar )
	729	ENDIF
	730	CALL start_inter_3d('U', iml, jml, lml, lon_in,
	731	. lat_in, jml2, lon_in2, lat_in2, pls, champ,interbar )
	732	DO il=1,lml
	733	DO ij=1,jml
	734	DO ii=1,iml-1
	735	champ(ii,ij,il) = champ(ii,ij,il) * cu(ii,ij)
	736	ENDDO
	737	champ(iml,ij, il) = champ(1,ij, il)
	738	ENDDO
	739	ENDDO
	740	CASE ('v')
	741	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	742	CALL start_init_dyn( iml, jml, lon_in, lat_in , jml2,
	743	. lon_in2, lat_in2 , interbar )
	744	ENDIF
	745	CALL start_inter_3d('V', iml, jml, lml, lon_in,
	746	. lat_in, jml2, lon_in2, lat_in2, pls, champ, interbar )
	747	DO il=1,lml
	748	DO ij=1,jml
	749	DO ii=1,iml-1
	750	champ(ii,ij,il) = champ(ii,ij,il) * cv(ii,ij)
	751	ENDDO
	752	champ(iml,ij, il) = champ(1,ij, il)
	753	ENDDO
	754	ENDDO
	755	CASE ('t')
	756	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	757	CALL start_init_dyn( iml, jml, lon_in, lat_in, jml2 ,
	758	. lon_in2, lat_in2 ,interbar )
	759	ENDIF
	760	CALL start_inter_3d('TEMP', iml, jml, lml, lon_in,
	761	. lat_in, jml2, lon_in2, lat_in2, pls, champ, interbar )
	762
	763	CASE ('tpot')
	764	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	765	CALL start_init_dyn( iml, jml, lon_in, lat_in , jml2 ,
	766	. lon_in2, lat_in2 , interbar )
	767	ENDIF
	768	CALL start_inter_3d('TEMP', iml, jml, lml, lon_in,
	769	. lat_in, jml2, lon_in2, lat_in2, pls, champ, interbar )
	770	IF ( MINVAL(workvar(:,:,:)) .NE. MAXVAL(workvar(:,:,:)) )
	771	. THEN
	772	DO il=1,lml
	773	DO ij=1,jml
	774	DO ii=1,iml-1
	775	champ(ii,ij,il) = champ(ii,ij,il) * cpp
	776	. / workvar(ii,ij,il)
	777	ENDDO
	778	champ(iml,ij,il) = champ(1,ij,il)
	779	ENDDO
	780	ENDDO
	781	DO il=1,lml
	782	xppn = SUM(aire(:,1)*champ(:,1,il))/apoln
	783	xpps = SUM(aire(:,jml)*champ(:,jml,il))/apols
	784	champ(:,1,il) = xppn
	785	champ(:,jml,il) = xpps
	786	ENDDO
	787	ELSE
	788	WRITE(,)'Could not compute potential temperature as the'
	789	WRITE(,)'Exner function is missing or constant.'
	790	STOP
	791	ENDIF
	792	CASE ('q')
	793	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	794	CALL start_init_dyn( iml, jml, lon_in, lat_in, jml2 ,
	795	. lon_in2, lat_in2 , interbar )
	796	ENDIF
	797	CALL start_inter_3d('R', iml, jml, lml, lon_in, lat_in,
	798	. jml2, lon_in2, lat_in2, pls, champ, interbar )
	799	IF ( MINVAL(workvar(:,:,:)) .NE. MAXVAL(workvar(:,:,:)) )
	800	. THEN
	801	DO il=1,lml
	802	DO ij=1,jml
	803	DO ii=1,iml-1
	804	champ(ii,ij,il) = 0.01 * champ(ii,ij,il) *
	805	. workvar(ii,ij,il)
	806	ENDDO
	807	champ(iml,ij,il) = champ(1,ij,il)
	808	ENDDO
	809	ENDDO
	810	WHERE ( champ .LT. 0.) champ = 1.0E-10
	811	DO il=1,lml
	812	xppn = SUM(aire(:,1)*champ(:,1,il))/apoln
	813	xpps = SUM(aire(:,jml)*champ(:,jml,il))/apols
	814	champ(:,1,il) = xppn
	815	champ(:,jml,il) = xpps
	816	ENDDO
	817	ELSE
	818	WRITE(,)'Could not compute specific humidity as the'
	819	WRITE(,)'saturated humidity is missing or constant.'
	820	STOP
	821	ENDIF
	822	CASE DEFAULT
	823	WRITE(,) 'startget_dyn'
	824	WRITE(,) 'No rule is present to extract variable ',
	825	. varname(:LEN_TRIM(varname)),' from any data set'
	826	STOP
	827	END SELECT
	828	ENDIF
	829	END SUBROUTINE startget_dyn
	830	!
	831	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	832	!
	833	SUBROUTINE start_init_dyn( iml, jml, lon_in, lat_in,jml2,lon_in2 ,
	834	, lat_in2 , interbar )
	835	!
	836	INTEGER, INTENT(in) :: iml, jml, jml2
	837	REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
	838	REAL, INTENT(in) :: lon_in2(iml), lat_in2(jml2)
	839	LOGICAL interbar
	840	!
	841	! LOCAL
	842	!
	843	REAL :: lev(1), date, dt
	844	INTEGER :: itau(1)
	845	INTEGER :: i, j
	846	integer :: iret
	847	!
	848	CHARACTER*120 :: physfname
	849	LOGICAL :: check=.TRUE.
	850	!
	851	REAL, ALLOCATABLE :: lon_rad(:), lat_rad(:)
	852	REAL, ALLOCATABLE :: lon_ini(:), lat_ini(:)
	853	REAL, ALLOCATABLE :: var_ana(:,:), tmp_var(:,:), z(:,:)
	854	REAL, ALLOCATABLE :: xppn(:), xpps(:)
	855	LOGICAL :: allo
	856	!
	857	!
	858	#include "dimensions.h"
	859	#include "paramet.h"
	860	#include "comgeom2.h"
	861
	862	CHARACTER*25 title
	863
	864	!
	865	physfname = 'ECDYN.nc'
	866	!
	867	IF ( check ) WRITE(,) 'Opening the surface analysis'
	868	!
	869	CALL flininfo(physfname, iml_dyn, jml_dyn, llm_dyn,
	870	. ttm_dyn, fid_dyn)
	871	IF ( check ) WRITE(,) 'Values read: ', iml_dyn, jml_dyn,
	872	. llm_dyn, ttm_dyn
	873	!
	874	ALLOCATE (lat_dyn(iml_dyn,jml_dyn), stat=iret)
	875	ALLOCATE (lon_dyn(iml_dyn,jml_dyn), stat=iret)
	876	ALLOCATE (levdyn_ini(llm_dyn), stat=iret)
	877	!
	878	CALL flinopen(physfname, .FALSE., iml_dyn, jml_dyn, llm_dyn,
	879	. lon_dyn, lat_dyn, levdyn_ini, ttm_dyn,
	880	. itau, date, dt, fid_dyn)
	881	!
	882
	883	allo = allocated (var_ana)
	884	if (allo) then
	885	DEALLOCATE(var_ana, stat=iret)
	886	endif
	887	ALLOCATE(var_ana(iml_dyn, jml_dyn), stat=iret)
	888
	889	allo = allocated (lon_rad)
	890	if (allo) then
	891	DEALLOCATE(lon_rad, stat=iret)
	892	endif
	893
	894	ALLOCATE(lon_rad(iml_dyn), stat=iret)
	895	ALLOCATE(lon_ini(iml_dyn))
	896
	897	IF ( MAXVAL(lon_dyn(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	898	lon_ini(:) = lon_dyn(:,1) * 2.0 * ASIN(1.0) / 180.0
	899	ELSE
	900	lon_ini(:) = lon_dyn(:,1)
	901	ENDIF
	902
	903	ALLOCATE(lat_rad(jml_dyn))
	904	ALLOCATE(lat_ini(jml_dyn))
	905
	906	IF ( MAXVAL(lat_dyn(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	907	lat_ini(:) = lat_dyn(1,:) * 2.0 * ASIN(1.0) / 180.0
	908	ELSE
	909	lat_ini(:) = lat_dyn(1,:)
	910	ENDIF
	911	!
	912
	913
	914	ALLOCATE(z(iml, jml))
	915	ALLOCATE(tmp_var(iml-1,jml))
	916	!
	917	CALL flinget(fid_dyn, 'Z', iml_dyn, jml_dyn, 0, ttm_dyn,
	918	. 1, 1, var_ana)
	919	c
	920	title='Z'
	921	CALL conf_dat2d( title,iml_dyn, jml_dyn,lon_ini, lat_ini,
	922	. lon_rad, lat_rad, var_ana, interbar )
	923	c
	924	IF ( interbar ) THEN
	925	WRITE(6,*) '-------------------------------------------------',
	926	,'--------------'
	927	WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
	928	, ' pour Z $$$ '
	929	WRITE(6,*) '-------------------------------------------------',
	930	,'--------------'
	931	CALL inter_barxy ( iml_dyn,jml_dyn -1,lon_rad,lat_rad ,
	932	, var_ana, iml-1, jml-1, lon_in2, lat_in2, jml, tmp_var)
	933	ELSE
	934	CALL grille_m(iml_dyn, jml_dyn , lon_rad, lat_rad, var_ana,
	935	. iml-1, jml, lon_in, lat_in, tmp_var)
	936	ENDIF
	937
	938	CALL gr_int_dyn(tmp_var, z, iml-1, jml)
	939	!
	940	ALLOCATE(psol_dyn(iml, jml))
	941	!
	942	CALL flinget(fid_dyn, 'SP', iml_dyn, jml_dyn, 0, ttm_dyn,
	943	. 1, 1, var_ana)
	944
	945	title='SP'
	946	CALL conf_dat2d( title,iml_dyn, jml_dyn,lon_ini, lat_ini,
	947	. lon_rad, lat_rad, var_ana, interbar )
	948
	949	IF ( interbar ) THEN
	950	WRITE(6,*) '-------------------------------------------------',
	951	,'--------------'
	952	WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
	953	, ' pour SP $$$ '
	954	WRITE(6,*) '-------------------------------------------------',
	955	,'--------------'
	956	CALL inter_barxy ( iml_dyn,jml_dyn -1,lon_rad,lat_rad ,
	957	, var_ana, iml-1, jml-1, lon_in2, lat_in2, jml, tmp_var)
	958	ELSE
	959	CALL grille_m(iml_dyn, jml_dyn , lon_rad, lat_rad, var_ana,
	960	. iml-1, jml, lon_in, lat_in, tmp_var )
	961	ENDIF
	962
	963	CALL gr_int_dyn(tmp_var, psol_dyn, iml-1, jml)
	964	!
	965	IF ( .NOT.ALLOCATED(tsol)) THEN
	966	! These variables may have been allocated by the need to
	967	! create a start field for them or by the varibale
	968	! coming out of the restart file. In case we dor have it we will initialize it.
	969	!
	970	CALL start_init_phys( iml, jml, lon_in, lat_in,jml2,lon_in2,
	971	. lat_in2 , interbar )
	972	ELSE
	973	IF ( SIZE(tsol) .NE. SIZE(psol_dyn) ) THEN
	974	WRITE(,) 'start_init_dyn :'
	975	WRITE(,) 'The temperature field we have does not ',
	976	. 'have the right size'
	977	STOP
	978	ENDIF
	979	ENDIF
	980	IF ( .NOT.ALLOCATED(phis)) THEN
	981	!
	982	! These variables may have been allocated by the need to create a start field for them or by the varibale
	983	! coming out of the restart file. In case we dor have it we will initialize it.
	984	!
	985	CALL start_init_orog( iml, jml, lon_in, lat_in, jml2, lon_in2 ,
	986	. lat_in2 , interbar )
	987	!
	988	ELSE
	989	!
	990	IF (SIZE(phis) .NE. SIZE(psol_dyn)) THEN
	991	!
	992	WRITE(,) 'start_init_dyn :'
	993	WRITE(,) 'The orography field we have does not ',
	994	. ' have the right size'
	995	STOP
	996	ENDIF
	997	!
	998	ENDIF
	999	!
	1000	! PSOL is computed in Pascals
	1001	!
	1002	!
	1003	DO j = 1, jml
	1004	DO i = 1, iml-1
	1005	psol_dyn(i,j) = psol_dyn(i,j)*(1.0+(z(i,j)-phis(i,j))
	1006	. /287.0/tsol(i,j))
	1007	ENDDO
	1008	psol_dyn(iml,j) = psol_dyn(1,j)
	1009	ENDDO
	1010	!
	1011	!
	1012	ALLOCATE(xppn(iml-1))
	1013	ALLOCATE(xpps(iml-1))
	1014	!
	1015	DO i = 1, iml-1
	1016	xppn(i) = aire( i,1) * psol_dyn( i,1)
	1017	xpps(i) = aire( i,jml) * psol_dyn( i,jml)
	1018	ENDDO
	1019	!
	1020	DO i = 1, iml
	1021	psol_dyn(i,1 ) = SUM(xppn)/apoln
	1022	psol_dyn(i,jml) = SUM(xpps)/apols
	1023	ENDDO
	1024	!
	1025	RETURN
	1026	!
	1027	END SUBROUTINE start_init_dyn
	1028	!
	1029	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	1030	!
	1031	SUBROUTINE start_inter_3d(varname, iml, jml, lml, lon_in,
	1032	. lat_in, jml2, lon_in2, lat_in2, pls_in, var3d, interbar )
	1033	!
	1034	! This subroutine gets a variables from a 3D file and does the interpolations needed
	1035	!
	1036	!
	1037	! ARGUMENTS
	1038	!
	1039	CHARACTER() :: varname
	1040	INTEGER :: iml, jml, lml, jml2
	1041	REAL :: lon_in(iml), lat_in(jml), pls_in(iml, jml, lml)
	1042	REAL :: lon_in2(iml) , lat_in2(jml2)
	1043	REAL :: var3d(iml, jml, lml)
	1044	LOGICAL interbar
	1045	real chmin,chmax
	1046	!
	1047	! LOCAL
	1048	!
	1049	CHARACTER*25 title
	1050	INTEGER :: ii, ij, il, jsort,i,j,l
	1051	REAL :: bx, by
	1052	REAL, ALLOCATABLE :: lon_rad(:), lat_rad(:)
	1053	REAL, ALLOCATABLE :: lon_ini(:), lat_ini(:) , lev_dyn(:)
	1054	REAL, ALLOCATABLE :: var_tmp2d(:,:), var_tmp3d(:,:,:)
	1055	REAL, ALLOCATABLE :: ax(:), ay(:), yder(:)
[1000]	1056	! REAL, ALLOCATABLE :: varrr(:,:,:)
[630]	1057	INTEGER, ALLOCATABLE :: lind(:)
	1058	!
	1059	LOGICAL :: check = .TRUE.
	1060	!
	1061	IF ( .NOT. ALLOCATED(var_ana3d)) THEN
	1062	ALLOCATE(var_ana3d(iml_dyn, jml_dyn, llm_dyn))
	1063	ENDIF
[1000]	1064	! ALLOCATE(varrr(iml_dyn, jml_dyn, llm_dyn))
[630]	1065	!
	1066	!
	1067	IF ( check) WRITE(,) 'Going into flinget to extract the 3D ',
	1068	. ' field.', fid_dyn
	1069	IF ( check) WRITE(,) fid_dyn, varname, iml_dyn, jml_dyn,
	1070	. llm_dyn,ttm_dyn
	1071	!
	1072	CALL flinget(fid_dyn, varname, iml_dyn, jml_dyn, llm_dyn,
	1073	. ttm_dyn, 1, 1, var_ana3d)
	1074	!
	1075	IF ( check) WRITE(,) 'Allocating space for the interpolation',
	1076	. iml, jml, llm_dyn
	1077	!
	1078	ALLOCATE(lon_rad(iml_dyn))
	1079	ALLOCATE(lon_ini(iml_dyn))
	1080
	1081	IF ( MAXVAL(lon_dyn(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	1082	lon_ini(:) = lon_dyn(:,1) * 2.0 * ASIN(1.0) / 180.0
	1083	ELSE
	1084	lon_ini(:) = lon_dyn(:,1)
	1085	ENDIF
	1086
	1087	ALLOCATE(lat_rad(jml_dyn))
	1088	ALLOCATE(lat_ini(jml_dyn))
	1089
	1090	ALLOCATE(lev_dyn(llm_dyn))
	1091
	1092	IF ( MAXVAL(lat_dyn(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	1093	lat_ini(:) = lat_dyn(1,:) * 2.0 * ASIN(1.0) / 180.0
	1094	ELSE
	1095	lat_ini(:) = lat_dyn(1,:)
	1096	ENDIF
	1097	!
	1098
	1099	CALL conf_dat3d ( varname,iml_dyn, jml_dyn, llm_dyn, lon_ini,
	1100	. lat_ini, levdyn_ini, lon_rad, lat_rad, lev_dyn, var_ana3d ,
	1101	, interbar )
	1102
	1103	ALLOCATE(var_tmp2d(iml-1, jml))
	1104	ALLOCATE(var_tmp3d(iml, jml, llm_dyn))
	1105	ALLOCATE(ax(llm_dyn))
	1106	ALLOCATE(ay(llm_dyn))
	1107	ALLOCATE(yder(llm_dyn))
	1108	ALLOCATE(lind(llm_dyn))
	1109	!
	1110
	1111	DO il=1,llm_dyn
	1112	!
	1113	IF( interbar ) THEN
	1114	IF( il.EQ.1 ) THEN
	1115	WRITE(6,*) '-------------------------------------------------',
	1116	,'--------------'
	1117	WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
	1118	, ' pour ', varname
	1119	WRITE(6,*) '-------------------------------------------------',
	1120	,'--------------'
	1121	ENDIF
	1122	CALL inter_barxy ( iml_dyn, jml_dyn -1,lon_rad, lat_rad,
	1123	, var_ana3d(:,:,il),iml-1, jml2, lon_in2, lat_in2,jml,var_tmp2d )
	1124	ELSE
	1125	CALL grille_m(iml_dyn, jml_dyn, lon_rad, lat_rad,
	1126	. var_ana3d(:,:,il), iml-1, jml, lon_in, lat_in, var_tmp2d )
	1127	ENDIF
	1128	!
	1129	CALL gr_int_dyn(var_tmp2d, var_tmp3d(:,:,il), iml-1, jml)
	1130	!
	1131	ENDDO
	1132	!
	1133	DO il=1,llm_dyn
	1134	lind(il) = llm_dyn-il+1
	1135	ENDDO
	1136	!
	1137	c
	1138	c ... Pour l'interpolation verticale ,on interpole du haut de l'atmosphere
	1139	c vers le sol ...
	1140	c
	1141	DO ij=1,jml
	1142	DO ii=1,iml-1
	1143	!
	1144	ax(:) = lev_dyn(lind(:))
	1145	ay(:) = var_tmp3d(ii, ij, lind(:))
	1146	!
	1147
	1148	CALL SPLINE(ax, ay, llm_dyn, 1.e30, 1.e30, yder)
	1149	!
	1150	DO il=1,lml
	1151	bx = pls_in(ii, ij, il)
	1152	CALL SPLINT(ax, ay, yder, llm_dyn, bx, by)
	1153	var3d(ii, ij, il) = by
	1154	ENDDO
	1155	!
	1156	ENDDO
	1157	var3d(iml, ij, :) = var3d(1, ij, :)
	1158	ENDDO
	1159
	1160	do il=1,lml
	1161	call minmax(iml*jml,var3d(1,1,il),chmin,chmax)
	1162	SELECTCASE(varname)
	1163	CASE('U')
	1164	WRITE(,) ' U min max l ',il,chmin,chmax
	1165	CASE('V')
	1166	WRITE(,) ' V min max l ',il,chmin,chmax
	1167	CASE('TEMP')
	1168	WRITE(,) ' TEMP min max l ',il,chmin,chmax
	1169	CASE('R')
	1170	WRITE(,) ' R min max l ',il,chmin,chmax
	1171	END SELECT
	1172	enddo
	1173
	1174	DEALLOCATE(lon_rad)
[1000]	1175	DEALLOCATE(lon_ini)
[630]	1176	DEALLOCATE(lat_rad)
[1000]	1177	DEALLOCATE(lat_ini)
	1178	DEALLOCATE(lev_dyn)
[630]	1179	DEALLOCATE(var_tmp2d)
	1180	DEALLOCATE(var_tmp3d)
	1181	DEALLOCATE(ax)
	1182	DEALLOCATE(ay)
	1183	DEALLOCATE(yder)
	1184	DEALLOCATE(lind)
	1185
	1186	!
	1187	RETURN
	1188	!
	1189	END SUBROUTINE start_inter_3d
	1190	!
[1222]	1191	#endif
	1192	! of #ifdef CPP_EARTH
[630]	1193	END MODULE startvar

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

Download in other formats:

Original Format