Context Navigation

limit_netcdf.F @ 415

Last change on this file since 415 was 259, checked in by lmdz, 23 years ago
Nouveaux programmes pour la creation des etats initiaux et des conditions aux limites. 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.7 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	CALL conf_dat2d(title,imdep, jmdep, dlon_ini, dlat_ini,
528	. dlon, dlat, champ, interbar )
529	c
530	IF( oldice ) THEN
531	CALL sea_ice(imdep, jmdep, dlon, dlat, champ,
532	. iim, jjp1, rlonv, rlatu, champint )
533	ELSEIF ( interbar ) THEN
534	IF( l.EQ.1 ) THEN
535	WRITE(6,*) '-------------------------------------------------',
536	,'------------------------'
537	WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
538	, ' pour Sea-ice Amip $$$ '
539	WRITE(6,*) '-------------------------------------------------',
540	,'------------------------'
541	ENDIF
542
543	CALL inter_barxy ( imdep,jmdep -1,dlon, dlat ,
544	, champ, iim, jjm, rlonu, rlatv, jjp1, champint )
545	ELSE
546	CALL sea_ice(imdep, jmdep, dlon, dlat, champ,
547	. iim, jjp1, rlonv, rlatu, champint )
548	ENDIF
549	DO j = 1,jjp1
550	DO i = 1, iim
551	champtime (i,j,l) = champint(i,j)
552	ENDDO
553	ENDDO
554	ENDDO
555	c
556	DO l = 1, lmdep
557	timeyear(l) = timecoord(l)
558	ENDDO
559	PRINT 222, timeyear
560	c
561	PRINT*, 'Interpolation temporelle'
562	DO j = 1, jjp1
563	DO i = 1, iim
564	DO l = 1, lmdep
565	ax(l) = timeyear(l)
566	ay(l) = champtime (i,j,l)
567	ENDDO
568	CALL SPLINE(ax,ay,lmdep,1.e30,1.e30,yder)
569	DO k = 1, 360
570	time = FLOAT(k-1)
571	CALL SPLINT(ax,ay,yder,lmdep,time,by)
572	champan(i,j,k) = by
573	ENDDO
574	ENDDO
575	ENDDO
576	DO k = 1, 360
577	DO j = 1, jjp1
578	champan(iip1, j, k) = champan(1, j, k)
579	ENDDO
580	IF ( k.EQ.10 ) THEN
581	DO j = 1, jjp1
582	CALL minmax( iip1,champan(1,j,10),chmin,chmax )
583	PRINT *,' Sea ice au temps 10 ', chmin,chmax,j
584	ENDDO
585	ENDIF
586	ENDDO
587	c
588
589	DO k = 1, 360
590	CALL gr_dyn_fi(1, iip1, jjp1, klon,
591	. champan(1,1,k), phy_ice(1,k))
592	DO i = 1, klon
593	phy_nat(i,k) = phy_nat0(i)
594	IF ( (phy_ice(i,k) - 0.5).GE.1.e-5 ) THEN
595	IF (NINT(phy_nat0(i)).EQ.0) THEN
596	phy_nat(i,k) = 3.0
597	ELSE
598	phy_nat(i,k) = 2.0
599	ENDIF
600	ENDIF
601	IF( NINT(phy_nat(i,k)).EQ.0 ) THEN
602	IF ( phy_rug(i,k).NE.0.001 ) phy_rug(i,k) = 0.001
603	ENDIF
604	ENDDO
605
606	ENDDO
607	c
608
609	ierr = NF_CLOSE(ncid)
610	c
611	DEALLOCATE( dlon )
612	DEALLOCATE( dlon_ini )
613	DEALLOCATE( dlat )
614	DEALLOCATE( dlat_ini )
615	DEALLOCATE( champ )
616
617	477 continue
618	c
619	C Traitement de la sst
620	c
621	PRINT*, 'Traitement de la sst'
622	ierr = NF_OPEN('AMIP.nc', NF_NOWRITE, ncid)
623	if (ierr.ne.0) then
624	print *, NF_STRERROR(ierr)
625	STOP
626	ENDIF
627
628	ierr = NF_INQ_VARID(ncid,'SST',varid)
629	if (ierr.ne.0) then
630	print *, NF_STRERROR(ierr)
631	STOP
632	ENDIF
633	ierr = NF_INQ_VARDIMID (ncid,varid,ndimid)
634	if (ierr.ne.0) then
635	print *, NF_STRERROR(ierr)
636	STOP
637	ENDIF
638	ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep)
639	if (ierr.ne.0) then
640	print *, NF_STRERROR(ierr)
641	STOP
642	ENDIF
643	print*,'variable SST ', namedim,'dimension ', imdep
644	ierr = NF_INQ_VARID(ncid,namedim,dimid)
645	if (ierr.ne.0) then
646	print *, NF_STRERROR(ierr)
647	STOP
648	ENDIF
649
650	ALLOCATE( dlon_ini(imdep) )
651	ALLOCATE( dlon(imdep) )
652
653	#ifdef NC_DOUBLE
654	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlon_ini)
655	#else
656	ierr = NF_GET_VAR_REAL(ncid,dimid,dlon_ini)
657	#endif
658
659	if (ierr.ne.0) then
660	print *, NF_STRERROR(ierr)
661	STOP
662	ENDIF
663	ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep)
664	if (ierr.ne.0) then
665	print *, NF_STRERROR(ierr)
666	STOP
667	ENDIF
668	print*,'variable SST ', namedim, 'dimension ', jmdep
669	ierr = NF_INQ_VARID(ncid,namedim,dimid)
670	if (ierr.ne.0) then
671	print *, NF_STRERROR(ierr)
672	STOP
673	ENDIF
674
675	ALLOCATE( dlat_ini(jmdep) )
676	ALLOCATE( dlat(jmdep) )
677
678	#ifdef NC_DOUBLE
679	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlat_ini)
680	#else
681	ierr = NF_GET_VAR_REAL(ncid,dimid,dlat_ini)
682	#endif
683	if (ierr.ne.0) then
684	print *, NF_STRERROR(ierr)
685	STOP
686	ENDIF
687	ierr = NF_INQ_DIM(ncid,ndimid(3), namedim, lmdep)
688	if (ierr.ne.0) then
689	print *, NF_STRERROR(ierr)
690	STOP
691	ENDIF
692	print*,'variable ', namedim, 'dimension ', lmdep
693	ierr = NF_INQ_VARID(ncid,namedim,dimid)
694	if (ierr.ne.0) then
695	print *, NF_STRERROR(ierr)
696	STOP
697	ENDIF
698	#ifdef NC_DOUBLE
699	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,timecoord)
700	#else
701	ierr = NF_GET_VAR_REAL(ncid,dimid,timecoord)
702	#endif
703	if (ierr.ne.0) then
704	print *, NF_STRERROR(ierr)
705	STOP
706	ENDIF
707
708	ALLOCATE( champ(imdep*jmdep) )
709	IF( extrap ) THEN
710	ALLOCATE ( work(imdep,jmdep) )
711	ENDIF
712	c
713	DO l = 1, lmdep
714	dimfirst(1) = 1
715	dimfirst(2) = 1
716	dimfirst(3) = l
717	c
718	dimlast(1) = imdep
719	dimlast(2) = jmdep
720	dimlast(3) = 1
721	c
722	PRINT*,'Lecture temporelle et int. horizontale ',l,timecoord(l)
723	#ifdef NC_DOUBLE
724	ierr = NF_GET_VARA_DOUBLE(ncid,varid,dimfirst,dimlast,champ)
725	#else
726	ierr = NF_GET_VARA_REAL(ncid,varid,dimfirst,dimlast,champ)
727	#endif
728	if (ierr.ne.0) then
729	print *, NF_STRERROR(ierr)
730	STOP
731	ENDIF
732
733	title='Sst Amip'
734	c
735	CALL conf_dat2d(title,imdep, jmdep, dlon_ini, dlat_ini,
736	. dlon, dlat, champ, interbar )
737	IF ( extrap ) THEN
738	CALL extrapol(champ, imdep, jmdep, 999999.,.TRUE.,.TRUE.,2,work)
739	ENDIF
740	c
741
742	IF ( interbar ) THEN
743	IF( l.EQ.1 ) THEN
744	WRITE(6,*) '-------------------------------------------------',
745	,'------------------------'
746	WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
747	, ' pour la Sst Amip $$$ '
748	WRITE(6,*) '-------------------------------------------------',
749	,'------------------------'
750	ENDIF
751	CALL inter_barxy ( imdep,jmdep -1,dlon, dlat ,
752	, champ, iim, jjm, rlonu, rlatv, jjp1, champint )
753	ELSE
754	CALL grille_m(imdep, jmdep, dlon, dlat, champ,
755	. iim, jjp1, rlonv, rlatu, champint )
756	ENDIF
757
758	DO j = 1,jjp1
759	DO i = 1, iim
760	champtime (i,j,l) = champint(i,j)
761	ENDDO
762	ENDDO
763	ENDDO
764	c
765	DO l = 1, lmdep
766	timeyear(l) = timecoord(l)
767	ENDDO
768	print 222, timeyear
769	c
770	C interpolation temporelle
771	DO j = 1, jjp1
772	DO i = 1, iim
773	DO l = 1, lmdep
774	ax(l) = timeyear(l)
775	ay(l) = champtime (i,j,l)
776	ENDDO
777	CALL SPLINE(ax,ay,lmdep,1.e30,1.e30,yder)
778	DO k = 1, 360
779	time = FLOAT(k-1)
780	CALL SPLINT(ax,ay,yder,lmdep,time,by)
781	champan(i,j,k) = by
782	ENDDO
783	ENDDO
784	ENDDO
785	DO k = 1, 360
786	DO j = 1, jjp1
787	champan(iip1,j,k) = champan(1,j,k)
788	ENDDO
789	IF ( k.EQ.10 ) THEN
790	DO j = 1, jjp1
791	CALL minmax( iip1,champan(1,j,10),chmin,chmax )
792	PRINT *,' SST au temps 10 ', chmin,chmax,j
793	ENDDO
794	ENDIF
795	ENDDO
796	c
797	DO k = 1, 360
798	CALL gr_dyn_fi(1, iip1, jjp1, klon,
799	. champan(1,1,k), phy_sst(1,k))
800	ENDDO
801	c
802	ierr = NF_CLOSE(ncid)
803	c
804	DEALLOCATE( dlon )
805	DEALLOCATE( dlon_ini )
806	DEALLOCATE( dlat )
807	DEALLOCATE( dlat_ini )
808	DEALLOCATE( champ )
809	c
810	C Traitement de l'albedo
811	c
812	PRINT*, 'Traitement de l albedo'
813	ierr = NF_OPEN('Albedo.nc', NF_NOWRITE, ncid)
814	if (ierr.ne.0) then
815	print *, NF_STRERROR(ierr)
816	STOP
817	ENDIF
818	ierr = NF_INQ_VARID(ncid,'ALBEDO',varid)
819	if (ierr.ne.0) then
820	print *, NF_STRERROR(ierr)
821	STOP
822	ENDIF
823	ierr = NF_INQ_VARDIMID (ncid,varid,ndimid)
824	if (ierr.ne.0) then
825	print *, NF_STRERROR(ierr)
826	STOP
827	ENDIF
828	ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep)
829	if (ierr.ne.0) then
830	print *, NF_STRERROR(ierr)
831	STOP
832	ENDIF
833	print*,'variable ', namedim, 'dimension ', imdep
834	ierr = NF_INQ_VARID(ncid,namedim,dimid)
835	if (ierr.ne.0) then
836	print *, NF_STRERROR(ierr)
837	STOP
838	ENDIF
839
840	ALLOCATE ( dlon_ini(imdep) )
841	ALLOCATE ( dlon(imdep) )
842
843	#ifdef NC_DOUBLE
844	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlon_ini)
845	#else
846	ierr = NF_GET_VAR_REAL(ncid,dimid,dlon_ini)
847	#endif
848	if (ierr.ne.0) then
849	print *, NF_STRERROR(ierr)
850	STOP
851	ENDIF
852	ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep)
853	if (ierr.ne.0) then
854	print *, NF_STRERROR(ierr)
855	STOP
856	ENDIF
857	print*,'variable ', namedim, 'dimension ', jmdep
858	ierr = NF_INQ_VARID(ncid,namedim,dimid)
859	if (ierr.ne.0) then
860	print *, NF_STRERROR(ierr)
861	STOP
862	ENDIF
863
864	ALLOCATE ( dlat_ini(jmdep) )
865	ALLOCATE ( dlat(jmdep) )
866
867	#ifdef NC_DOUBLE
868	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlat_ini)
869	#else
870	ierr = NF_GET_VAR_REAL(ncid,dimid,dlat_ini)
871	#endif
872	if (ierr.ne.0) then
873	print *, NF_STRERROR(ierr)
874	STOP
875	ENDIF
876	ierr = NF_INQ_DIM(ncid,ndimid(3), namedim, lmdep)
877	if (ierr.ne.0) then
878	print *, NF_STRERROR(ierr)
879	STOP
880	ENDIF
881	print*,'variable ', namedim, 'dimension ', lmdep
882	ierr = NF_INQ_VARID(ncid,namedim,dimid)
883	if (ierr.ne.0) then
884	print *, NF_STRERROR(ierr)
885	STOP
886	ENDIF
887	#ifdef NC_DOUBLE
888	ierr = NF_GET_VAR_DOUBLE(ncid,dimid,timecoord)
889	#else
890	ierr = NF_GET_VAR_REAL(ncid,dimid,timecoord)
891	#endif
892	if (ierr.ne.0) then
893	print *, NF_STRERROR(ierr)
894	STOP
895	ENDIF
896	c
897	ALLOCATE ( champ(imdep*jmdep) )
898
899	DO l = 1, lmdep
900	dimfirst(1) = 1
901	dimfirst(2) = 1
902	dimfirst(3) = l
903	c
904	dimlast(1) = imdep
905	dimlast(2) = jmdep
906	dimlast(3) = 1
907	c
908	PRINT*,'Lecture temporelle et int. horizontale ',l,timecoord(l)
909	#ifdef NC_DOUBLE
910	ierr = NF_GET_VARA_DOUBLE(ncid,varid,dimfirst,dimlast,champ)
911	#else
912	ierr = NF_GET_VARA_REAL(ncid,varid,dimfirst,dimlast,champ)
913	#endif
914	if (ierr.ne.0) then
915	print *, NF_STRERROR(ierr)
916	STOP
917	ENDIF
918
919	title='Albedo Amip'
920	c
921	CALL conf_dat2d(title,imdep, jmdep, dlon_ini, dlat_ini,
922	. dlon, dlat, champ, interbar )
923	c
924	c
925	IF ( interbar ) THEN
926	IF( l.EQ.1 ) THEN
927	WRITE(6,*) '-------------------------------------------------',
928	,'------------------------'
929	WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
930	, ' pour l Albedo Amip $$$ '
931	WRITE(6,*) '-------------------------------------------------',
932	,'------------------------'
933	ENDIF
934
935	CALL inter_barxy ( imdep,jmdep -1,dlon, dlat ,
936	, champ, iim, jjm, rlonu, rlatv, jjp1, champint )
937	ELSE
938	CALL grille_m(imdep, jmdep, dlon, dlat, champ,
939	. iim, jjp1, rlonv, rlatu, champint )
940	ENDIF
941	c
942	DO j = 1,jjp1
943	DO i = 1, iim
944	champtime (i, j, l) = champint(i, j)
945	ENDDO
946	ENDDO
947	ENDDO
948	c
949	DO l = 1, lmdep
950	timeyear(l) = timecoord(l)
951	ENDDO
952	print 222, timeyear
953	c
954	C interpolation temporelle
955	DO j = 1, jjp1
956	DO i = 1, iim
957	DO l = 1, lmdep
958	ax(l) = timeyear(l)
959	ay(l) = champtime (i, j, l)
960	ENDDO
961	CALL SPLINE(ax,ay,lmdep,1.e30,1.e30,yder)
962	DO k = 1, 360
963	time = FLOAT(k-1)
964	CALL SPLINT(ax,ay,yder,lmdep,time,by)
965	champan(i,j,k) = by
966	ENDDO
967	ENDDO
968	ENDDO
969	DO k = 1, 360
970	DO j = 1, jjp1
971	champan(iip1, j, k) = champan(1, j, k)
972	ENDDO
973	IF ( k.EQ.10 ) THEN
974	DO j = 1, jjp1
975	CALL minmax( iip1,champan(1,j,10),chmin,chmax )
976	PRINT *,' Albedo au temps 10 ', chmin,chmax,j
977	ENDDO
978	ENDIF
979	ENDDO
980	c
981	DO k = 1, 360
982	CALL gr_dyn_fi(1, iip1, jjp1, klon,
983	. champan(1,1,k), phy_alb(1,k))
984	ENDDO
985	c
986	ierr = NF_CLOSE(ncid)
987	c
988	c
989	DO k = 1, 360
990	DO i = 1, klon
991	phy_bil(i,k) = 0.0
992	ENDDO
993	ENDDO
994	c
995	PRINT*, 'Ecriture du fichier limit'
996	c
997	ierr = NF_CREATE ("limit.nc", NF_CLOBBER, nid)
998	c
999	ierr = NF_PUT_ATT_TEXT (nid, NF_GLOBAL, "title", 30,
1000	. "Fichier conditions aux limites")
1001	ierr = NF_DEF_DIM (nid, "points_physiques", klon, ndim)
1002	ierr = NF_DEF_DIM (nid, "time", NF_UNLIMITED, ntim)
1003	c
1004	dims(1) = ndim
1005	dims(2) = ntim
1006	c
1007	ccc ierr = NF_DEF_VAR (nid, "TEMPS", NF_DOUBLE, 1,ntim, id_tim)
1008	ierr = NF_DEF_VAR (nid, "TEMPS", NF_FLOAT, 1,ntim, id_tim)
1009	ierr = NF_PUT_ATT_TEXT (nid, id_tim, "title", 17,
1010	. "Jour dans l annee")
1011	ccc ierr = NF_DEF_VAR (nid, "NAT", NF_DOUBLE, 2,dims, id_NAT)
1012	ierr = NF_DEF_VAR (nid, "NAT", NF_FLOAT, 2,dims, id_NAT)
1013	ierr = NF_PUT_ATT_TEXT (nid, id_NAT, "title", 23,
1014	. "Nature du sol (0,1,2,3)")
1015	ccc ierr = NF_DEF_VAR (nid, "SST", NF_DOUBLE, 2,dims, id_SST)
1016	ierr = NF_DEF_VAR (nid, "SST", NF_FLOAT, 2,dims, id_SST)
1017	ierr = NF_PUT_ATT_TEXT (nid, id_SST, "title", 35,
1018	. "Temperature superficielle de la mer")
1019	ccc ierr = NF_DEF_VAR (nid, "BILS", NF_DOUBLE, 2,dims, id_BILS)
1020	ierr = NF_DEF_VAR (nid, "BILS", NF_FLOAT, 2,dims, id_BILS)
1021	ierr = NF_PUT_ATT_TEXT (nid, id_BILS, "title", 32,
1022	. "Reference flux de chaleur au sol")
1023	ccc ierr = NF_DEF_VAR (nid, "ALB", NF_DOUBLE, 2,dims, id_ALB)
1024	ierr = NF_DEF_VAR (nid, "ALB", NF_FLOAT, 2,dims, id_ALB)
1025	ierr = NF_PUT_ATT_TEXT (nid, id_ALB, "title", 19,
1026	. "Albedo a la surface")
1027	ccc ierr = NF_DEF_VAR (nid, "RUG", NF_DOUBLE, 2,dims, id_RUG)
1028	ierr = NF_DEF_VAR (nid, "RUG", NF_FLOAT, 2,dims, id_RUG)
1029	ierr = NF_PUT_ATT_TEXT (nid, id_RUG, "title", 8,
1030	. "Rugosite")
1031	c
1032	ierr = NF_ENDDEF(nid)
1033	c
1034	DO k = 1, 360
1035	c
1036	debut(1) = 1
1037	debut(2) = k
1038	epais(1) = klon
1039	epais(2) = 1
1040	c
1041	#ifdef NC_DOUBLE
1042	ierr = NF_PUT_VAR1_DOUBLE (nid,id_tim,k,DBLE(k))
1043	ierr = NF_PUT_VARA_DOUBLE (nid,id_NAT,debut,epais,phy_nat(1,k))
1044	ierr = NF_PUT_VARA_DOUBLE (nid,id_SST,debut,epais,phy_sst(1,k))
1045	ierr = NF_PUT_VARA_DOUBLE (nid,id_BILS,debut,epais,phy_bil(1,k))
1046	ierr = NF_PUT_VARA_DOUBLE (nid,id_ALB,debut,epais,phy_alb(1,k))
1047	ierr = NF_PUT_VARA_DOUBLE (nid,id_RUG,debut,epais,phy_rug(1,k))
1048	#else
1049	ierr = NF_PUT_VAR1_REAL (nid,id_tim,k,FLOAT(k))
1050	ierr = NF_PUT_VARA_REAL (nid,id_NAT,debut,epais,phy_nat(1,k))
1051	ierr = NF_PUT_VARA_REAL (nid,id_SST,debut,epais,phy_sst(1,k))
1052	ierr = NF_PUT_VARA_REAL (nid,id_BILS,debut,epais,phy_bil(1,k))
1053	ierr = NF_PUT_VARA_REAL (nid,id_ALB,debut,epais,phy_alb(1,k))
1054	ierr = NF_PUT_VARA_REAL (nid,id_RUG,debut,epais,phy_rug(1,k))
1055	#endif
1056	c
1057	ENDDO
1058	c
1059	ierr = NF_CLOSE(nid)
1060	c
1061	STOP
1062	END

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

Download in other formats:

Original Format