Context Navigation

startvar.F @ 3400

Last change on this file since 3400 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

Line
1	!
2	! $Id: startvar.F 1222 2009-08-07 11:48:33Z fairhead $
3	!
4	MODULE startvar
5	#ifdef CPP_EARTH
6	! This module is designed to work for Earth (and with ioipsl)
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
523	cIM "slab" ocean
524	CASE ('tslab')
525	champ(:) = 0.0
526	CASE ('seaice')
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(:)
1056	! REAL, ALLOCATABLE :: varrr(:,:,:)
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
1064	! ALLOCATE(varrr(iml_dyn, jml_dyn, llm_dyn))
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)
1175	DEALLOCATE(lon_ini)
1176	DEALLOCATE(lat_rad)
1177	DEALLOCATE(lat_ini)
1178	DEALLOCATE(lev_dyn)
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	!
1191	#endif
1192	! of #ifdef CPP_EARTH
1193	END MODULE startvar

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

Download in other formats:

Original Format