Context Navigation

limit_netcdf.F @ 507

Last change on this file since 507 was 425, checked in by lmdz, 22 years ago

Dans le cas "oldice" la routine "sea_ice" utilise seulement l'interpol.
grille_m Levan

LF

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

File size: 28.9 KB

Line
1	C
2	C $Header$
3	C
4	SUBROUTINE limit_netcdf ( interbar, extrap, oldice )
5	c
6	IMPLICIT none
7	c
8	c-------------------------------------------------------------
9	C Author : L. Fairhead
10	C Date : 27/01/94
11	C Objet : Construction des fichiers de conditions aux limites
12	C pour le nouveau
13	C modele a partir de fichiers de climatologie. Les deux
14	C grilles doivent etre regulieres
15	c
16	c Modifie par z.x.li (le23mars1994)
17	c Modifie par L. Fairhead (fairhead@lmd.jussieu.fr) septembre 1999
18	c pour lecture netcdf dans LMDZ.3.3
19	c Modifie par P;Le Van , juillet 2001
20	c-------------------------------------------------------------
21	c
22	#include "dimensions.h"
23	#include "paramet.h"
24	#include "control.h"
25	#include "logic.h"
26	#include "netcdf.inc"
27	#include "comvert.h"
28	#include "comgeom2.h"
29	#include "comconst.h"
30	c
31	c-----------------------------------------------------------------------
32	LOGICAL interbar, extrap, oldice
33
34	INTEGER KIDIA, KFDIA, KLON, KLEV
35	PARAMETER (KIDIA=1,KFDIA=iim*(jjm-1)+2,
36	. KLON=KFDIA-KIDIA+1,KLEV=llm)
37	c-----------------------------------------------------------------------
38	REAL phy_nat(klon,360), phy_nat0(klon)
39	REAL phy_alb(klon,360)
40	REAL phy_sst(klon,360)
41	REAL phy_bil(klon,360)
42	REAL phy_rug(klon,360)
43	REAL phy_ice(klon,360)
44	c
45	REAL masque(iip1,jjp1)
46	REAL mask(iim,jjp1)
47
48	C Declarations pour le champ de depart
49	INTEGER imdep, jmdep,lmdep
50	INTEGER tbid
51	PARAMETER ( tbid = 60 ) ! >52 semaines
52	REAL timecoord(tbid)
53	c
54	REAL , ALLOCATABLE :: dlon_msk(:), dlat_msk(:)
55	REAL , ALLOCATABLE :: lonmsk_ini(:), latmsk_ini(:)
56	REAL , ALLOCATABLE :: dlon(:), dlat(:)
57	REAL , ALLOCATABLE :: dlon_ini(:), dlat_ini(:)
58	REAL , ALLOCATABLE :: champ_msk(:), champ(:)
59	REAL , ALLOCATABLE :: work(:,:)
60
61	CHARACTER*25 title
62
63	C Declarations pour le champ interpole 2D
64	REAL champint(iim,jjp1)
65	real chmin,chmax
66
67	C Declarations pour le champ interpole 3D
68	REAL champtime(iim,jjp1,tbid)
69	REAL timeyear(tbid)
70	REAL champan(iip1,jjp1,366)
71
72	C Declarations pour l'inteprolation verticale
73	REAL ax(tbid), ay(tbid)
74	REAL by
75	REAL yder(tbid)
76
77
78	INTEGER ierr
79	INTEGER dimfirst(3)
80	INTEGER dimlast(3)
81	c
82	INTEGER nid, ndim, ntim
83	INTEGER dims(2), debut(2), epais(2)
84	INTEGER id_tim
85	INTEGER id_NAT, id_SST, id_BILS, id_RUG, id_ALB
86
87	INTEGER i, j, k, l
88	c Diverses variables locales
89	REAL time
90
91	INTEGER longcles
92	PARAMETER ( longcles = 20 )
93	REAL clesphy0 ( longcles )
94	#include "serre.h"
95	INTEGER ncid,varid,ndimid(4),dimid
96	character*30 namedim
97
98	pi = 4. * ATAN(1.)
99	rad = 6 371 229.
100	omeg = 4.* ASIN(1.)/(24.*3600.)
101	g = 9.8
102	daysec = 86400.
103	kappa = 0.2857143
104	cpp = 1004.70885
105	dtvr = daysec/FLOAT(day_step)
106	c
107	C Traitement du relief au sol
108	c
109	write(,) 'Traitement du relief au sol pour fabriquer masque'
110	ierr = NF_OPEN('Relief.nc', NF_NOWRITE, ncid)
111
112	if (ierr.ne.0) then
113	print *, NF_STRERROR(ierr)
114	STOP
115	ENDIF
116
117	ierr = NF_INQ_VARID(ncid,'RELIEF',varid)
118	if (ierr.ne.0) then
119	print *, NF_STRERROR(ierr)
120	STOP
121	ENDIF
122	ierr = NF_INQ_VARDIMID (ncid,varid,ndimid)
123	if (ierr.ne.0) then
124	print *, NF_STRERROR(ierr)
125	STOP
126	ENDIF
127	ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep)
128	if (ierr.ne.0) then
129	print *, NF_STRERROR(ierr)
130	STOP
131	ENDIF
132	print*,'variable ', namedim, 'dimension ', imdep
133	ierr = NF_INQ_VARID(ncid,namedim,dimid)
134	if (ierr.ne.0) then
135	print *, NF_STRERROR(ierr)
136	STOP
137	ENDIF
138
139	ALLOCATE( lonmsk_ini(imdep) )
140	ALLOCATE( dlon_msk(imdep) )
141
142	#ifdef NC_DOUBLE
143	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,lonmsk_ini)
144	#else
145	ierr = NF_GET_VAR_REAL(ncid,dimid,lonmsk_ini)
146	#endif
147
148	c
149	if (ierr.ne.0) then
150	print *, NF_STRERROR(ierr)
151	STOP
152	ENDIF
153	ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep)
154	if (ierr.ne.0) then
155	print *, NF_STRERROR(ierr)
156	STOP
157	ENDIF
158	print*,'variable ', namedim, 'dimension ', jmdep
159	ierr = NF_INQ_VARID(ncid,namedim,dimid)
160	if (ierr.ne.0) then
161	print *, NF_STRERROR(ierr)
162	STOP
163	ENDIF
164
165	ALLOCATE( latmsk_ini(jmdep) )
166	ALLOCATE( dlat_msk(jmdep) )
167	ALLOCATE( champ_msk(imdep*jmdep) )
168
169	#ifdef NC_DOUBLE
170	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,latmsk_ini)
171	#else
172	ierr = NF_GET_VAR_REAL(ncid,dimid,latmsk_ini)
173	#endif
174	c
175	if (ierr.ne.0) then
176	print *, NF_STRERROR(ierr)
177	STOP
178	ENDIF
179	#ifdef NC_DOUBLE
180	ierr = NF_GET_VAR_DOUBLE(ncid,varid,champ_msk)
181	#else
182	ierr = NF_GET_VAR_REAL(ncid,varid,champ_msk)
183	#endif
184	c
185	if (ierr.ne.0) then
186	print *, NF_STRERROR(ierr)
187	STOP
188	ENDIF
189	c
190	title='RELIEF'
191
192	CALL conf_dat2d(title,imdep, jmdep, lonmsk_ini, latmsk_ini,
193	. dlon_msk, dlat_msk, champ_msk, interbar )
194
195	CALL mask_c_o(imdep, jmdep, dlon_msk, dlat_msk,champ_msk,
196	. iim, jjp1, rlonv, rlatu, champint)
197	CALL gr_int_dyn(champint, masque, iim, jjp1)
198	DO i = 1, iim
199	masque(i,1) = FLOAT(NINT(masque(i,1)))
200	masque(i,jjp1) = FLOAT(NINT(masque(i,jjp1)))
201	ENDDO
202	DO i = 1, iim
203	DO j = 1, jjp1
204	mask(i,j) = champint(i,j)
205	ENDDO
206	ENDDO
207	CALL gr_dyn_fi(1, iip1, jjp1, klon, masque, phy_nat0)
208	ierr = NF_CLOSE(ncid)
209	c
210	c
211	C Traitement de la rugosite
212	c
213	PRINT*, 'Traitement de la rugosite'
214	ierr = NF_OPEN('Rugos.nc', NF_NOWRITE, ncid)
215	if (ierr.ne.0) then
216	print *, NF_STRERROR(ierr)
217	STOP
218	ENDIF
219
220	ierr = NF_INQ_VARID(ncid,'RUGOS',varid)
221	if (ierr.ne.0) then
222	print *, NF_STRERROR(ierr)
223	STOP
224	ENDIF
225	ierr = NF_INQ_VARDIMID (ncid,varid,ndimid)
226	if (ierr.ne.0) then
227	print *, NF_STRERROR(ierr)
228	STOP
229	ENDIF
230	ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep)
231	if (ierr.ne.0) then
232	print *, NF_STRERROR(ierr)
233	STOP
234	ENDIF
235	print*,'variable ', namedim, 'dimension ', imdep
236	ierr = NF_INQ_VARID(ncid,namedim,dimid)
237	if (ierr.ne.0) then
238	print *, NF_STRERROR(ierr)
239	STOP
240	ENDIF
241
242	ALLOCATE( dlon_ini(imdep) )
243	ALLOCATE( dlon(imdep) )
244
245	#ifdef NC_DOUBLE
246	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlon_ini)
247	#else
248	ierr = NF_GET_VAR_REAL(ncid,dimid,dlon_ini)
249	#endif
250	if (ierr.ne.0) then
251	print *, NF_STRERROR(ierr)
252	STOP
253	ENDIF
254	ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep)
255	if (ierr.ne.0) then
256	print *, NF_STRERROR(ierr)
257	STOP
258	ENDIF
259	print*,'variable ', namedim, 'dimension ', jmdep
260	ierr = NF_INQ_VARID(ncid,namedim,dimid)
261	if (ierr.ne.0) then
262	print *, NF_STRERROR(ierr)
263	STOP
264	ENDIF
265
266	ALLOCATE( dlat_ini(jmdep) )
267	ALLOCATE( dlat(jmdep) )
268
269	#ifdef NC_DOUBLE
270	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlat_ini)
271	#else
272	ierr = NF_GET_VAR_REAL(ncid,dimid,dlat_ini)
273	#endif
274	if (ierr.ne.0) then
275	print *, NF_STRERROR(ierr)
276	STOP
277	ENDIF
278	ierr = NF_INQ_DIM(ncid,ndimid(3), namedim, lmdep)
279	if (ierr.ne.0) then
280	print *, NF_STRERROR(ierr)
281	STOP
282	ENDIF
283	print*,'variable ', namedim, 'dimension ', lmdep
284	ierr = NF_INQ_VARID(ncid,namedim,dimid)
285	if (ierr.ne.0) then
286	print *, NF_STRERROR(ierr)
287	STOP
288	ENDIF
289	#ifdef NC_DOUBLE
290	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,timecoord)
291	#else
292	ierr = NF_GET_VAR_REAL(ncid,dimid,timecoord)
293	#endif
294	if (ierr.ne.0) then
295	print *, NF_STRERROR(ierr)
296	STOP
297	ENDIF
298	c
299	ALLOCATE( champ(imdep*jmdep) )
300
301	DO 200 l = 1, lmdep
302	dimfirst(1) = 1
303	dimfirst(2) = 1
304	dimfirst(3) = l
305	c
306	dimlast(1) = imdep
307	dimlast(2) = jmdep
308	dimlast(3) = 1
309	c
310	PRINT*,'Lecture temporelle et int. horizontale ',l,timecoord(l)
311	print*,dimfirst,dimlast
312	#ifdef NC_DOUBLE
313	ierr = NF_GET_VARA_DOUBLE(ncid,varid,dimfirst,dimlast,champ)
314	#else
315	ierr = NF_GET_VARA_REAL(ncid,varid,dimfirst,dimlast,champ)
316	#endif
317	if (ierr.ne.0) then
318	print *, NF_STRERROR(ierr)
319	STOP
320	ENDIF
321
322	title = 'Rugosite Amip '
323	c
324	CALL conf_dat2d(title,imdep, jmdep, dlon_ini, dlat_ini,
325	. dlon, dlat, champ, interbar )
326
327	IF ( interbar ) THEN
328	DO j = 1, imdep * jmdep
329	champ(j) = LOG(champ(j))
330	ENDDO
331
332	IF( l.EQ.1 ) THEN
333	WRITE(6,*) '-------------------------------------------------',
334	,'------------------------'
335	WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
336	, ' pour la rugosite $$$ '
337	WRITE(6,*) '-------------------------------------------------',
338	,'------------------------'
339	ENDIF
340	CALL inter_barxy ( imdep,jmdep -1,dlon,dlat,champ ,
341	, iim,jjm,rlonu,rlatv, jjp1,champint )
342	DO j=1,jjp1
343	DO i=1,iim
344	champint(i,j)=EXP(champint(i,j))
345	ENDDO
346	ENDDO
347
348	DO j = 1, jjp1
349	DO i = 1, iim
350	IF(NINT(mask(i,j)).NE.1) THEN
351	champint( i,j ) = 0.001
352	ENDIF
353	ENDDO
354	ENDDO
355	ELSE
356	CALL rugosite(imdep, jmdep, dlon, dlat, champ,
357	. iim, jjp1, rlonv, rlatu, champint, mask)
358	ENDIF
359	DO j = 1,jjp1
360	DO i = 1, iim
361	champtime (i,j,l) = champint(i,j)
362	ENDDO
363	ENDDO
364	200 CONTINUE
365	c
366	DO l = 1, lmdep
367	timeyear(l) = timecoord(l)
368	ENDDO
369
370	PRINT 222, timeyear
371	222 FORMAT(2x,' Time year ',10f6.1)
372	c
373
374	PRINT*, 'Interpolation temporelle dans l annee'
375
376	DO j = 1, jjp1
377	DO i = 1, iim
378	DO l = 1, lmdep
379	ax(l) = timeyear(l)
380	ay(l) = champtime (i,j,l)
381	ENDDO
382	CALL SPLINE(ax,ay,lmdep,1.e30,1.e30,yder)
383	DO k = 1, 360
384	time = FLOAT(k-1)
385	CALL SPLINT(ax,ay,yder,lmdep,time,by)
386	champan(i,j,k) = by
387	ENDDO
388	ENDDO
389	ENDDO
390	DO k = 1, 360
391	DO j = 1, jjp1
392	champan(iip1,j,k) = champan(1,j,k)
393	ENDDO
394	IF ( k.EQ.10 ) THEN
395	DO j = 1, jjp1
396	CALL minmax( iip1,champan(1,j,10),chmin,chmax )
397	PRINT *,' Rugosite au temps 10 ', chmin,chmax,j
398	ENDDO
399	ENDIF
400	ENDDO
401	c
402	DO k = 1, 360
403	CALL gr_dyn_fi(1,iip1,jjp1,klon,champan(1,1,k), phy_rug(1,k))
404	ENDDO
405	c
406	ierr = NF_CLOSE(ncid)
407
408	DEALLOCATE( dlon )
409	DEALLOCATE( dlon_ini )
410	DEALLOCATE( dlat )
411	DEALLOCATE( dlat_ini )
412	DEALLOCATE( champ )
413	c
414	c
415	C Traitement de la glace oceanique
416	c
417	PRINT*, 'Traitement de la glace oceanique'
418	ierr = NF_OPEN('AMIP.nc', NF_NOWRITE, ncid)
419	if (ierr.ne.0) then
420	print *, NF_STRERROR(ierr)
421	STOP
422	ENDIF
423
424	ierr = NF_INQ_VARID(ncid,'SEA_ICE',varid)
425	if (ierr.ne.0) then
426	print *, NF_STRERROR(ierr)
427	STOP
428	ENDIF
429	ierr = NF_INQ_VARDIMID (ncid,varid,ndimid)
430	if (ierr.ne.0) then
431	print *, NF_STRERROR(ierr)
432	STOP
433	ENDIF
434	ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep)
435	if (ierr.ne.0) then
436	print *, NF_STRERROR(ierr)
437	STOP
438	ENDIF
439	print*,'variable ', namedim, 'dimension ', imdep
440	ierr = NF_INQ_VARID(ncid,namedim,dimid)
441	if (ierr.ne.0) then
442	print *, NF_STRERROR(ierr)
443	STOP
444	ENDIF
445
446	ALLOCATE ( dlon_ini(imdep) )
447	ALLOCATE ( dlon(imdep) )
448
449	#ifdef NC_DOUBLE
450	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlon_ini)
451	#else
452	ierr = NF_GET_VAR_REAL(ncid,dimid,dlon_ini)
453	#endif
454	if (ierr.ne.0) then
455	print *, NF_STRERROR(ierr)
456	STOP
457	ENDIF
458	ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep)
459	if (ierr.ne.0) then
460	print *, NF_STRERROR(ierr)
461	STOP
462	ENDIF
463	print*,'variable ', namedim, jmdep
464	ierr = NF_INQ_VARID(ncid,namedim,dimid)
465	if (ierr.ne.0) then
466	print *, NF_STRERROR(ierr)
467	STOP
468	ENDIF
469
470	ALLOCATE ( dlat_ini(jmdep) )
471	ALLOCATE ( dlat(jmdep) )
472
473	#ifdef NC_DOUBLE
474	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlat_ini)
475	#else
476	ierr = NF_GET_VAR_REAL(ncid,dimid,dlat_ini)
477	#endif
478	if (ierr.ne.0) then
479	print *, NF_STRERROR(ierr)
480	STOP
481	ENDIF
482	ierr = NF_INQ_DIM(ncid,ndimid(3), namedim, lmdep)
483	if (ierr.ne.0) then
484	print *, NF_STRERROR(ierr)
485	STOP
486	ENDIF
487	print*,'variable ', namedim, lmdep
488	ierr = NF_INQ_VARID(ncid,namedim,dimid)
489	if (ierr.ne.0) then
490	print *, NF_STRERROR(ierr)
491	STOP
492	ENDIF
493	#ifdef NC_DOUBLE
494	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,timecoord)
495	#else
496	ierr = NF_GET_VAR_REAL(ncid,dimid,timecoord)
497	#endif
498	if (ierr.ne.0) then
499	print *, NF_STRERROR(ierr)
500	STOP
501	ENDIF
502	c
503	ALLOCATE ( champ(imdep*jmdep) )
504
505	DO l = 1, lmdep
506	dimfirst(1) = 1
507	dimfirst(2) = 1
508	dimfirst(3) = l
509	c
510	dimlast(1) = imdep
511	dimlast(2) = jmdep
512	dimlast(3) = 1
513	c
514	PRINT*,'Lecture temporelle et int. horizontale ',l,timecoord(l)
515	#ifdef NC_DOUBLE
516	ierr = NF_GET_VARA_DOUBLE(ncid,varid,dimfirst,dimlast,champ)
517	#else
518	ierr = NF_GET_VARA_REAL(ncid,varid,dimfirst,dimlast,champ)
519	#endif
520	if (ierr.ne.0) then
521	print *, NF_STRERROR(ierr)
522	STOP
523	ENDIF
524
525	title = 'Sea-ice Amip '
526	c
527	IF( oldice ) THEN
528	CALL conf_dat2d(title,imdep, jmdep, dlon_ini, dlat_ini,
529	. dlon, dlat, champ, .FALSE. )
530	ELSE
531	CALL conf_dat2d(title,imdep, jmdep, dlon_ini, dlat_ini,
532	. dlon, dlat, champ, interbar )
533	ENDIF
534	c
535	IF( oldice ) THEN
536	CALL sea_ice(imdep, jmdep, dlon, dlat, champ,
537	. iim, jjp1, rlonv, rlatu, champint )
538	ELSEIF ( interbar ) THEN
539	IF( l.EQ.1 ) THEN
540	WRITE(6,*) '-------------------------------------------------',
541	,'------------------------'
542	WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
543	, ' pour Sea-ice Amip $$$ '
544	WRITE(6,*) '-------------------------------------------------',
545	,'------------------------'
546	ENDIF
547
548	CALL inter_barxy ( imdep,jmdep -1,dlon, dlat ,
549	, champ, iim, jjm, rlonu, rlatv, jjp1, champint )
550	ELSE
551	CALL sea_ice(imdep, jmdep, dlon, dlat, champ,
552	. iim, jjp1, rlonv, rlatu, champint )
553	ENDIF
554	DO j = 1,jjp1
555	DO i = 1, iim
556	champtime (i,j,l) = champint(i,j)
557	ENDDO
558	ENDDO
559	ENDDO
560	c
561	DO l = 1, lmdep
562	timeyear(l) = timecoord(l)
563	ENDDO
564	PRINT 222, timeyear
565	c
566	PRINT*, 'Interpolation temporelle'
567	DO j = 1, jjp1
568	DO i = 1, iim
569	DO l = 1, lmdep
570	ax(l) = timeyear(l)
571	ay(l) = champtime (i,j,l)
572	ENDDO
573	CALL SPLINE(ax,ay,lmdep,1.e30,1.e30,yder)
574	DO k = 1, 360
575	time = FLOAT(k-1)
576	CALL SPLINT(ax,ay,yder,lmdep,time,by)
577	champan(i,j,k) = by
578	ENDDO
579	ENDDO
580	ENDDO
581	DO k = 1, 360
582	DO j = 1, jjp1
583	champan(iip1, j, k) = champan(1, j, k)
584	ENDDO
585	IF ( k.EQ.10 ) THEN
586	DO j = 1, jjp1
587	CALL minmax( iip1,champan(1,j,10),chmin,chmax )
588	PRINT *,' Sea ice au temps 10 ', chmin,chmax,j
589	ENDDO
590	ENDIF
591	ENDDO
592	c
593
594	DO k = 1, 360
595	CALL gr_dyn_fi(1, iip1, jjp1, klon,
596	. champan(1,1,k), phy_ice(1,k))
597	DO i = 1, klon
598	phy_nat(i,k) = phy_nat0(i)
599	IF ( (phy_ice(i,k) - 0.5).GE.1.e-5 ) THEN
600	IF (NINT(phy_nat0(i)).EQ.0) THEN
601	phy_nat(i,k) = 3.0
602	ELSE
603	phy_nat(i,k) = 2.0
604	ENDIF
605	ENDIF
606	IF( NINT(phy_nat(i,k)).EQ.0 ) THEN
607	IF ( phy_rug(i,k).NE.0.001 ) phy_rug(i,k) = 0.001
608	ENDIF
609	ENDDO
610
611	ENDDO
612	c
613
614	ierr = NF_CLOSE(ncid)
615	c
616	DEALLOCATE( dlon )
617	DEALLOCATE( dlon_ini )
618	DEALLOCATE( dlat )
619	DEALLOCATE( dlat_ini )
620	DEALLOCATE( champ )
621
622	477 continue
623	c
624	C Traitement de la sst
625	c
626	PRINT*, 'Traitement de la sst'
627	ierr = NF_OPEN('AMIP.nc', NF_NOWRITE, ncid)
628	if (ierr.ne.0) then
629	print *, NF_STRERROR(ierr)
630	STOP
631	ENDIF
632
633	ierr = NF_INQ_VARID(ncid,'SST',varid)
634	if (ierr.ne.0) then
635	print *, NF_STRERROR(ierr)
636	STOP
637	ENDIF
638	ierr = NF_INQ_VARDIMID (ncid,varid,ndimid)
639	if (ierr.ne.0) then
640	print *, NF_STRERROR(ierr)
641	STOP
642	ENDIF
643	ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep)
644	if (ierr.ne.0) then
645	print *, NF_STRERROR(ierr)
646	STOP
647	ENDIF
648	print*,'variable SST ', namedim,'dimension ', imdep
649	ierr = NF_INQ_VARID(ncid,namedim,dimid)
650	if (ierr.ne.0) then
651	print *, NF_STRERROR(ierr)
652	STOP
653	ENDIF
654
655	ALLOCATE( dlon_ini(imdep) )
656	ALLOCATE( dlon(imdep) )
657
658	#ifdef NC_DOUBLE
659	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlon_ini)
660	#else
661	ierr = NF_GET_VAR_REAL(ncid,dimid,dlon_ini)
662	#endif
663
664	if (ierr.ne.0) then
665	print *, NF_STRERROR(ierr)
666	STOP
667	ENDIF
668	ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep)
669	if (ierr.ne.0) then
670	print *, NF_STRERROR(ierr)
671	STOP
672	ENDIF
673	print*,'variable SST ', namedim, 'dimension ', jmdep
674	ierr = NF_INQ_VARID(ncid,namedim,dimid)
675	if (ierr.ne.0) then
676	print *, NF_STRERROR(ierr)
677	STOP
678	ENDIF
679
680	ALLOCATE( dlat_ini(jmdep) )
681	ALLOCATE( dlat(jmdep) )
682
683	#ifdef NC_DOUBLE
684	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlat_ini)
685	#else
686	ierr = NF_GET_VAR_REAL(ncid,dimid,dlat_ini)
687	#endif
688	if (ierr.ne.0) then
689	print *, NF_STRERROR(ierr)
690	STOP
691	ENDIF
692	ierr = NF_INQ_DIM(ncid,ndimid(3), namedim, lmdep)
693	if (ierr.ne.0) then
694	print *, NF_STRERROR(ierr)
695	STOP
696	ENDIF
697	print*,'variable ', namedim, 'dimension ', lmdep
698	ierr = NF_INQ_VARID(ncid,namedim,dimid)
699	if (ierr.ne.0) then
700	print *, NF_STRERROR(ierr)
701	STOP
702	ENDIF
703	#ifdef NC_DOUBLE
704	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,timecoord)
705	#else
706	ierr = NF_GET_VAR_REAL(ncid,dimid,timecoord)
707	#endif
708	if (ierr.ne.0) then
709	print *, NF_STRERROR(ierr)
710	STOP
711	ENDIF
712
713	ALLOCATE( champ(imdep*jmdep) )
714	IF( extrap ) THEN
715	ALLOCATE ( work(imdep,jmdep) )
716	ENDIF
717	c
718	DO l = 1, lmdep
719	dimfirst(1) = 1
720	dimfirst(2) = 1
721	dimfirst(3) = l
722	c
723	dimlast(1) = imdep
724	dimlast(2) = jmdep
725	dimlast(3) = 1
726	c
727	PRINT*,'Lecture temporelle et int. horizontale ',l,timecoord(l)
728	#ifdef NC_DOUBLE
729	ierr = NF_GET_VARA_DOUBLE(ncid,varid,dimfirst,dimlast,champ)
730	#else
731	ierr = NF_GET_VARA_REAL(ncid,varid,dimfirst,dimlast,champ)
732	#endif
733	if (ierr.ne.0) then
734	print *, NF_STRERROR(ierr)
735	STOP
736	ENDIF
737
738	title='Sst Amip'
739	c
740	CALL conf_dat2d(title,imdep, jmdep, dlon_ini, dlat_ini,
741	. dlon, dlat, champ, interbar )
742	IF ( extrap ) THEN
743	CALL extrapol(champ, imdep, jmdep, 999999.,.TRUE.,.TRUE.,2,work)
744	ENDIF
745	c
746
747	IF ( interbar ) THEN
748	IF( l.EQ.1 ) THEN
749	WRITE(6,*) '-------------------------------------------------',
750	,'------------------------'
751	WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
752	, ' pour la Sst Amip $$$ '
753	WRITE(6,*) '-------------------------------------------------',
754	,'------------------------'
755	ENDIF
756	CALL inter_barxy ( imdep,jmdep -1,dlon, dlat ,
757	, champ, iim, jjm, rlonu, rlatv, jjp1, champint )
758	ELSE
759	CALL grille_m(imdep, jmdep, dlon, dlat, champ,
760	. iim, jjp1, rlonv, rlatu, champint )
761	ENDIF
762
763	DO j = 1,jjp1
764	DO i = 1, iim
765	champtime (i,j,l) = champint(i,j)
766	ENDDO
767	ENDDO
768	ENDDO
769	c
770	DO l = 1, lmdep
771	timeyear(l) = timecoord(l)
772	ENDDO
773	print 222, timeyear
774	c
775	C interpolation temporelle
776	DO j = 1, jjp1
777	DO i = 1, iim
778	DO l = 1, lmdep
779	ax(l) = timeyear(l)
780	ay(l) = champtime (i,j,l)
781	ENDDO
782	CALL SPLINE(ax,ay,lmdep,1.e30,1.e30,yder)
783	DO k = 1, 360
784	time = FLOAT(k-1)
785	CALL SPLINT(ax,ay,yder,lmdep,time,by)
786	champan(i,j,k) = by
787	ENDDO
788	ENDDO
789	ENDDO
790	DO k = 1, 360
791	DO j = 1, jjp1
792	champan(iip1,j,k) = champan(1,j,k)
793	ENDDO
794	IF ( k.EQ.10 ) THEN
795	DO j = 1, jjp1
796	CALL minmax( iip1,champan(1,j,10),chmin,chmax )
797	PRINT *,' SST au temps 10 ', chmin,chmax,j
798	ENDDO
799	ENDIF
800	ENDDO
801	c
802	DO k = 1, 360
803	CALL gr_dyn_fi(1, iip1, jjp1, klon,
804	. champan(1,1,k), phy_sst(1,k))
805	ENDDO
806	c
807	ierr = NF_CLOSE(ncid)
808	c
809	DEALLOCATE( dlon )
810	DEALLOCATE( dlon_ini )
811	DEALLOCATE( dlat )
812	DEALLOCATE( dlat_ini )
813	DEALLOCATE( champ )
814	c
815	C Traitement de l'albedo
816	c
817	PRINT*, 'Traitement de l albedo'
818	ierr = NF_OPEN('Albedo.nc', NF_NOWRITE, ncid)
819	if (ierr.ne.0) then
820	print *, NF_STRERROR(ierr)
821	STOP
822	ENDIF
823	ierr = NF_INQ_VARID(ncid,'ALBEDO',varid)
824	if (ierr.ne.0) then
825	print *, NF_STRERROR(ierr)
826	STOP
827	ENDIF
828	ierr = NF_INQ_VARDIMID (ncid,varid,ndimid)
829	if (ierr.ne.0) then
830	print *, NF_STRERROR(ierr)
831	STOP
832	ENDIF
833	ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep)
834	if (ierr.ne.0) then
835	print *, NF_STRERROR(ierr)
836	STOP
837	ENDIF
838	print*,'variable ', namedim, 'dimension ', imdep
839	ierr = NF_INQ_VARID(ncid,namedim,dimid)
840	if (ierr.ne.0) then
841	print *, NF_STRERROR(ierr)
842	STOP
843	ENDIF
844
845	ALLOCATE ( dlon_ini(imdep) )
846	ALLOCATE ( dlon(imdep) )
847
848	#ifdef NC_DOUBLE
849	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlon_ini)
850	#else
851	ierr = NF_GET_VAR_REAL(ncid,dimid,dlon_ini)
852	#endif
853	if (ierr.ne.0) then
854	print *, NF_STRERROR(ierr)
855	STOP
856	ENDIF
857	ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep)
858	if (ierr.ne.0) then
859	print *, NF_STRERROR(ierr)
860	STOP
861	ENDIF
862	print*,'variable ', namedim, 'dimension ', jmdep
863	ierr = NF_INQ_VARID(ncid,namedim,dimid)
864	if (ierr.ne.0) then
865	print *, NF_STRERROR(ierr)
866	STOP
867	ENDIF
868
869	ALLOCATE ( dlat_ini(jmdep) )
870	ALLOCATE ( dlat(jmdep) )
871
872	#ifdef NC_DOUBLE
873	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlat_ini)
874	#else
875	ierr = NF_GET_VAR_REAL(ncid,dimid,dlat_ini)
876	#endif
877	if (ierr.ne.0) then
878	print *, NF_STRERROR(ierr)
879	STOP
880	ENDIF
881	ierr = NF_INQ_DIM(ncid,ndimid(3), namedim, lmdep)
882	if (ierr.ne.0) then
883	print *, NF_STRERROR(ierr)
884	STOP
885	ENDIF
886	print*,'variable ', namedim, 'dimension ', lmdep
887	ierr = NF_INQ_VARID(ncid,namedim,dimid)
888	if (ierr.ne.0) then
889	print *, NF_STRERROR(ierr)
890	STOP
891	ENDIF
892	#ifdef NC_DOUBLE
893	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,timecoord)
894	#else
895	ierr = NF_GET_VAR_REAL(ncid,dimid,timecoord)
896	#endif
897	if (ierr.ne.0) then
898	print *, NF_STRERROR(ierr)
899	STOP
900	ENDIF
901	c
902	ALLOCATE ( champ(imdep*jmdep) )
903
904	DO l = 1, lmdep
905	dimfirst(1) = 1
906	dimfirst(2) = 1
907	dimfirst(3) = l
908	c
909	dimlast(1) = imdep
910	dimlast(2) = jmdep
911	dimlast(3) = 1
912	c
913	PRINT*,'Lecture temporelle et int. horizontale ',l,timecoord(l)
914	#ifdef NC_DOUBLE
915	ierr = NF_GET_VARA_DOUBLE(ncid,varid,dimfirst,dimlast,champ)
916	#else
917	ierr = NF_GET_VARA_REAL(ncid,varid,dimfirst,dimlast,champ)
918	#endif
919	if (ierr.ne.0) then
920	print *, NF_STRERROR(ierr)
921	STOP
922	ENDIF
923
924	title='Albedo Amip'
925	c
926	CALL conf_dat2d(title,imdep, jmdep, dlon_ini, dlat_ini,
927	. dlon, dlat, champ, interbar )
928	c
929	c
930	IF ( interbar ) THEN
931	IF( l.EQ.1 ) THEN
932	WRITE(6,*) '-------------------------------------------------',
933	,'------------------------'
934	WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
935	, ' pour l Albedo Amip $$$ '
936	WRITE(6,*) '-------------------------------------------------',
937	,'------------------------'
938	ENDIF
939
940	CALL inter_barxy ( imdep,jmdep -1,dlon, dlat ,
941	, champ, iim, jjm, rlonu, rlatv, jjp1, champint )
942	ELSE
943	CALL grille_m(imdep, jmdep, dlon, dlat, champ,
944	. iim, jjp1, rlonv, rlatu, champint )
945	ENDIF
946	c
947	DO j = 1,jjp1
948	DO i = 1, iim
949	champtime (i, j, l) = champint(i, j)
950	ENDDO
951	ENDDO
952	ENDDO
953	c
954	DO l = 1, lmdep
955	timeyear(l) = timecoord(l)
956	ENDDO
957	print 222, timeyear
958	c
959	C interpolation temporelle
960	DO j = 1, jjp1
961	DO i = 1, iim
962	DO l = 1, lmdep
963	ax(l) = timeyear(l)
964	ay(l) = champtime (i, j, l)
965	ENDDO
966	CALL SPLINE(ax,ay,lmdep,1.e30,1.e30,yder)
967	DO k = 1, 360
968	time = FLOAT(k-1)
969	CALL SPLINT(ax,ay,yder,lmdep,time,by)
970	champan(i,j,k) = by
971	ENDDO
972	ENDDO
973	ENDDO
974	DO k = 1, 360
975	DO j = 1, jjp1
976	champan(iip1, j, k) = champan(1, j, k)
977	ENDDO
978	IF ( k.EQ.10 ) THEN
979	DO j = 1, jjp1
980	CALL minmax( iip1,champan(1,j,10),chmin,chmax )
981	PRINT *,' Albedo au temps 10 ', chmin,chmax,j
982	ENDDO
983	ENDIF
984	ENDDO
985	c
986	DO k = 1, 360
987	CALL gr_dyn_fi(1, iip1, jjp1, klon,
988	. champan(1,1,k), phy_alb(1,k))
989	ENDDO
990	c
991	ierr = NF_CLOSE(ncid)
992	c
993	c
994	DO k = 1, 360
995	DO i = 1, klon
996	phy_bil(i,k) = 0.0
997	ENDDO
998	ENDDO
999	c
1000	PRINT*, 'Ecriture du fichier limit'
1001	c
1002	ierr = NF_CREATE ("limit.nc", NF_CLOBBER, nid)
1003	c
1004	ierr = NF_PUT_ATT_TEXT (nid, NF_GLOBAL, "title", 30,
1005	. "Fichier conditions aux limites")
1006	ierr = NF_DEF_DIM (nid, "points_physiques", klon, ndim)
1007	ierr = NF_DEF_DIM (nid, "time", NF_UNLIMITED, ntim)
1008	c
1009	dims(1) = ndim
1010	dims(2) = ntim
1011	c
1012	ccc ierr = NF_DEF_VAR (nid, "TEMPS", NF_DOUBLE, 1,ntim, id_tim)
1013	ierr = NF_DEF_VAR (nid, "TEMPS", NF_FLOAT, 1,ntim, id_tim)
1014	ierr = NF_PUT_ATT_TEXT (nid, id_tim, "title", 17,
1015	. "Jour dans l annee")
1016	ccc ierr = NF_DEF_VAR (nid, "NAT", NF_DOUBLE, 2,dims, id_NAT)
1017	ierr = NF_DEF_VAR (nid, "NAT", NF_FLOAT, 2,dims, id_NAT)
1018	ierr = NF_PUT_ATT_TEXT (nid, id_NAT, "title", 23,
1019	. "Nature du sol (0,1,2,3)")
1020	ccc ierr = NF_DEF_VAR (nid, "SST", NF_DOUBLE, 2,dims, id_SST)
1021	ierr = NF_DEF_VAR (nid, "SST", NF_FLOAT, 2,dims, id_SST)
1022	ierr = NF_PUT_ATT_TEXT (nid, id_SST, "title", 35,
1023	. "Temperature superficielle de la mer")
1024	ccc ierr = NF_DEF_VAR (nid, "BILS", NF_DOUBLE, 2,dims, id_BILS)
1025	ierr = NF_DEF_VAR (nid, "BILS", NF_FLOAT, 2,dims, id_BILS)
1026	ierr = NF_PUT_ATT_TEXT (nid, id_BILS, "title", 32,
1027	. "Reference flux de chaleur au sol")
1028	ccc ierr = NF_DEF_VAR (nid, "ALB", NF_DOUBLE, 2,dims, id_ALB)
1029	ierr = NF_DEF_VAR (nid, "ALB", NF_FLOAT, 2,dims, id_ALB)
1030	ierr = NF_PUT_ATT_TEXT (nid, id_ALB, "title", 19,
1031	. "Albedo a la surface")
1032	ccc ierr = NF_DEF_VAR (nid, "RUG", NF_DOUBLE, 2,dims, id_RUG)
1033	ierr = NF_DEF_VAR (nid, "RUG", NF_FLOAT, 2,dims, id_RUG)
1034	ierr = NF_PUT_ATT_TEXT (nid, id_RUG, "title", 8,
1035	. "Rugosite")
1036	c
1037	ierr = NF_ENDDEF(nid)
1038	c
1039	DO k = 1, 360
1040	c
1041	debut(1) = 1
1042	debut(2) = k
1043	epais(1) = klon
1044	epais(2) = 1
1045	c
1046	#ifdef NC_DOUBLE
1047	ierr = NF_PUT_VAR1_DOUBLE (nid,id_tim,k,DBLE(k))
1048	ierr = NF_PUT_VARA_DOUBLE (nid,id_NAT,debut,epais,phy_nat(1,k))
1049	ierr = NF_PUT_VARA_DOUBLE (nid,id_SST,debut,epais,phy_sst(1,k))
1050	ierr = NF_PUT_VARA_DOUBLE (nid,id_BILS,debut,epais,phy_bil(1,k))
1051	ierr = NF_PUT_VARA_DOUBLE (nid,id_ALB,debut,epais,phy_alb(1,k))
1052	ierr = NF_PUT_VARA_DOUBLE (nid,id_RUG,debut,epais,phy_rug(1,k))
1053	#else
1054	ierr = NF_PUT_VAR1_REAL (nid,id_tim,k,FLOAT(k))
1055	ierr = NF_PUT_VARA_REAL (nid,id_NAT,debut,epais,phy_nat(1,k))
1056	ierr = NF_PUT_VARA_REAL (nid,id_SST,debut,epais,phy_sst(1,k))
1057	ierr = NF_PUT_VARA_REAL (nid,id_BILS,debut,epais,phy_bil(1,k))
1058	ierr = NF_PUT_VARA_REAL (nid,id_ALB,debut,epais,phy_alb(1,k))
1059	ierr = NF_PUT_VARA_REAL (nid,id_RUG,debut,epais,phy_rug(1,k))
1060	#endif
1061	c
1062	ENDDO
1063	c
1064	ierr = NF_CLOSE(nid)
1065	c
1066	STOP
1067	END

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

Download in other formats:

Original Format