Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

startvar.F @ 575

Last change on this file since 575 was 575, checked in by (none), 19 years ago
This commit was manufactured by cvs2svn to create tag 'IPSL-CM4_CMIP'.
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 39.7 KB

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

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

Download in other formats:

Original Format