Context Navigation

ts2IPCC.F90 @ 605

Last change on this file since 605 was 605, checked in by Laurent Fairhead, 19 years ago
Un probleme sur la boucle de lecture du fichier de correspondance variable IPSL - variable IPCC empechait l'execution sur PC LF
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 14.3 KB

Line
1	!
2	! $Header$
3	!
4	PROGRAM ts2IPCC
5
6	!
7	! Filtre permettant de transformer les fichiers de serie temporelle a une
8	! variable de l'IPSL en fichiers acceptables par le PCMDI/IPCC.
9	!
10	! Utilisation de la bibliothèque CMOR du PCMDI
11	!
12	! L. Fairhead 2004/08
13	!
14	! Ce programme est appelé avec un argument, le nom du fichier à traiter
15	! Il nécessite aussi un fichier config.def contenant diverses informations
16	! (voir plus bas) et l'accès à un tableau faisant la correspondance entre
17	! les noms de variables du modèle et les noms imposés par l'IPCC
18	! Pour l'instant on ne traite que les fichiers contenant la serie
19	! temporelle d'une seule variable
20
21	use cmor_users_functions
22	use netcdf
23
24	implicit none
25
26	#include "netcdf.inc"
27
28	CHARACTER (len=256) :: orig_file ! nom du fichier à traiter
29	character (len=512) :: line_read
30	CHARACTER (len=128) :: inpath, contact, hist_gen,repert,instit
31	CHARACTER (len=128) :: hist_var,expt_id,source,comment,refs
32	CHARACTER (len=20) :: action
33	character (len=20) :: first_part
34	character (len=1004) :: second_part
35	CHARACTER (len=20), DIMENSION(100) :: ipsl_name, ipsl_units, ipsl_pos
36	CHARACTER (len=20), DIMENSION(100) :: ipcc_name, ipcc_table
37	CHARACTER (len=80) :: varname, units, namedim
38
39	INTEGER :: orig_file_id, nvars, ndims
40	INTEGER :: verbos, exit_ctl, realis, indice,index_table
41	INTEGER :: iargc, iostat, ierr
42	INTEGER :: i,idim
43	INTEGER, ALLOCATABLE, DIMENSION(:) :: dimids,axis_ids,lendim
44	INTEGER :: latid, lonid, vertid, timeid
45	INTEGER :: varid, cmorvarid
46	INTEGER :: ilat, ilon, ivert, itime
47	INTEGER :: lunout ! device de sortie
48
49	logical :: found = .false.
50
51	REAL, ALLOCATABLE, DIMENSION(:) :: lon, lat, vert, time
52	REAL, ALLOCATABLE, DIMENSION(:) :: lon_bounds, lat_bounds
53	REAL, ALLOCATABLE, DIMENSION(:,:,:,:) :: donnees
54	DOUBLE PRECISION, DIMENSION(1) :: rdate
55	real :: missing_value
56
57	external iargc
58
59	!
60	! quelques initialisations
61	lunout = 6
62	varname = 'xxxxxxxx'
63
64	!
65	! On vérifie que l'appel au programme a bien un argument:
66	CALL getarg(1, orig_file)
67	IF (iargc() == 0 .OR. orig_file == '-h') then
68	WRITE(lunout,*)' '
69	WRITE(lunout,*)' Utilisation de ce programme: '
70	WRITE(lunout,*)' ./ts2IPCC nom_de_fichier [variable]'
71	WRITE(lunout,*)' ou nom_de_fichier est le nom du fichier a traiter'
72	WRITE(lunout,*)' et variable la variable a traiter [optionel]'
73	WRITE(lunout,*)' '
74	WRITE(lunout,*)' ./ts2IPCC -h sort ce message'
75	WRITE(lunout,*)' '
76	stop
77	ENDIF
78	if (iargc() == 2) then
79	CALL getarg(2, varname)
80	endif
81
82	!
83	! Lecture du fichier de configuration
84	OPEN (20, IOSTAT=iostat, file='config.def',form='formatted')
85	IF (iostat /= 0) then
86	WRITE(lunout,*)'Erreur ouverture du fichier config.def'
87	stop
88	endif
89
90	config: do
91	READ(20,'(A)',iostat=iostat)line_read
92	IF (iostat /= 0) exit
93	line_read = TRIM(line_read)
94	IF (INDEX(line_read, '#') /= 1) THEN
95	first_part = trim(line_read(1:INDEX(line_read, '=')-1))
96	second_part = trim(line_read(INDEX(line_read, '=')+1:))
97	selectcase(first_part)
98	case('inpath')
99	inpath = trim(second_part)
100	case('file_action')
101	action = trim(second_part)
102	case('verbosity')
103	READ(second_part,'(i)') verbos
104	case('exit_control')
105	READ(second_part,'(i)') exit_ctl
106	case('repertoire')
107	repert = trim(second_part)
108	case('experiment_ID')
109	expt_id = trim(second_part)
110	case('institut')
111	instit = trim(second_part)
112	case('source')
113	source = trim(second_part)
114	case('realisation')
115	READ(second_part,'(i)') realis
116	case('hist_gen')
117	hist_gen = trim(second_part)
118	case('comment')
119	comment = trim(second_part)
120	case('refs')
121	refs = trim(second_part)
122	case('hist_var')
123	hist_var = trim(second_part)
124	case('contact')
125	contact = trim(second_part)
126	end select
127	endif
128	enddo config
129	if (iostat > 0) then
130	WRITE(lunout,*)'Probleme de lecture du fichier config.def, iostat = ',iostat
131	stop
132	endif
133	close(20)
134
135	!
136	! Lecture du tableau de correspondance nom IPSL <=> nom IPCC
137	OPEN (20, IOSTAT=iostat, file='table.def',form='formatted')
138	IF (iostat /= 0) then
139	WRITE(lunout,*)'Erreur ouverture du fichier table.def'
140	stop
141	endif
142	indice = 0
143	table: do
144	READ(20,'(A)',iostat=iostat)line_read
145	IF (iostat /= 0) exit
146	line_read = TRIM(line_read)
147	IF (INDEX(line_read, '#') /= 1) THEN
148	indice = indice + 1
149	ipsl_name(indice) = trim(line_read(1:INDEX(line_read, '\|')-1))
150	line_read = trim(line_read(INDEX(line_read, '\|')+1:))
151	ipsl_units(indice) = trim(line_read(1:INDEX(line_read, '\|')-1))
152	line_read = trim(line_read(INDEX(line_read, '\|')+1:))
153	ipsl_pos(indice) = trim(line_read(1:INDEX(line_read, '\|')-1))
154	line_read = trim(line_read(INDEX(line_read, '\|')+1:))
155	ipcc_name(indice) = trim(line_read(1:INDEX(line_read, '\|')-1))
156	ipcc_table(indice) = trim(line_read(INDEX(line_read, '\|')+1:))
157	endif
158	enddo table
159	indice = indice - 1
160	close(20)
161	! DO i = 1, indice
162	! WRITE(lunout,*)ipsl_name(i),ipsl_units(i),ipcc_name(i),ipcc_table(i)
163	! enddo
164
165	!
166	! Ouverture du fichier a traiter
167	ierr = nf_open(orig_file, NF_NOWRITE, orig_file_id)
168	IF (ierr /= NF_NOERR) then
169	WRITE(lunout,*)NF_STRERROR(ierr)
170	stop
171	endif
172	!
173	! trouver la variable a traiter, c'est une variable a 3 ou 4 dimensions
174	ierr = nf_inq_nvars(orig_file_id, nvars)
175	IF (ierr /= NF_NOERR) then
176	WRITE(lunout,*)NF_STRERROR(ierr)
177	stop
178	endif
179
180	i = 0
181	if (varname == 'xxxxxxxx') then
182	DO while (.not.found)
183	i = i + 1
184	if (i > nvars) then
185	WRITE(lunout,*)' pas de variable 3d ou 4d trouvee'
186	stop
187	endif
188	ierr = nf_inq_varname(orig_file_id, i, varname)
189	IF (ierr /= NF_NOERR) then
190	WRITE(lunout,*)NF_STRERROR(ierr)
191	stop
192	endif
193	ierr = nf_inq_varndims(orig_file_id, i, ndims)
194	IF (ierr /= NF_NOERR) then
195	WRITE(lunout,*)NF_STRERROR(ierr)
196	stop
197	endif
198	if (ndims > 2) found = .true.
199	enddo
200	else
201	ierr = nf_inq_varid(orig_file_id, varname, varid)
202	IF (ierr /= NF_NOERR) then
203	WRITE(lunout,*)NF_STRERROR(ierr)
204	stop
205	endif
206	endif
207
208	!
209	! recherche de la correspondance nom IPSL <=> nom IPCC
210	found = .false.
211	i = 0
212	do while (.not. found)
213	i = i + 1
214	if (i > indice) then
215	WRITE(lunout,*)'La variable ',trim(varname),' n''est pas dans le tableau de correspondance table.def'
216	stop
217	endif
218	IF (varname == ipsl_name(i)) THEN
219	index_table = i
220	found = .true.
221	endif
222	enddo
223
224	WRITE(lunout,*)' found variable = ', trim(varname)
225	WRITE(lunout,*)' ipcc_name = ', trim(ipcc_name(index_table))
226
227
228	!
229	! Initialisation CMOR
230	ierr = cmor_setup(inpath=inpath, netcdf_file_action=action,set_verbosity=verbos,&
231	& exit_control=exit_ctl)
232	IF (ierr /= 0) then
233	WRITE(lunout,*)'Probleme dans cmor_setup, ierr = ', ierr
234	endif
235
236	!
237	! Initialisation dataset
238	ierr = cmor_dataset(outpath=repert, &
239	& experiment_id=expt_id, &
240	& institution=instit, &
241	& source=source, &
242	& calendar='360_day', &
243	& realization=realis, &
244	& contact=contact, &
245	& history=hist_gen, &
246	& comment=comment, &
247	& references=refs)
248	IF (ierr /= 0) then
249	WRITE(lunout,*)'Probleme dans cmor_dataset, ierr = ', ierr
250	endif
251
252	!
253	! Definition des axes
254
255	ierr = nf90_inq_varid(orig_file_id,TRIM(varname), varid)
256	if (ierr /= 0) call handle_err(ierr)
257	ierr = nf90_Inquire_Variable(orig_file_id, varid, ndims = ndims)
258	if (ierr /= 0) call handle_err(ierr)
259	allocate (dimids(ndims))
260	allocate (axis_ids(ndims))
261	allocate (lendim(ndims))
262	ierr = nf90_Inquire_Variable(orig_file_id, varid, dimids = dimids)
263	if (ierr /= 0) call handle_err(ierr)
264
265	do idim = 1, ndims
266	ierr = nf90_Inquire_Dimension(orig_file_id, dimids(idim), &
267	name = namedim, len = lendim(idim))
268	if (ierr /= 0) call handle_err(ierr)
269	units=' '
270	selectcase(trim(namedim))
271	case('lat')
272	! lecture de la latitude:
273	allocate(lat(lendim(idim)))
274	ierr = nf_inq_varid(orig_file_id, namedim, latid)
275	if (ierr /= 0) call handle_err(ierr)
276	ierr = nf_get_var_real(orig_file_id, latid, lat)
277	if (ierr /= 0) call handle_err(ierr)
278	ierr = nf_get_att_text(orig_file_id, latid, 'units', units)
279	if (ierr /= 0) call handle_err(ierr)
280	allocate(lat_bounds(lendim(idim)+1))
281	DO i = 2, lendim(idim)
282	lat_bounds(i) = lat(i-1) - (lat(i-1) - lat(i))/2
283	enddo
284	lat_bounds(1) = lat(1)
285	lat_bounds(lendim(idim)+1) = lat(lendim(idim))
286	! definition de la latitude
287	axis_ids(idim) = cmor_axis( &
288	table=trim(ipcc_table(index_table)), &
289	table_entry='latitude', &
290	units=trim(units), &
291	length=lendim(idim), &
292	coord_vals=lat, &
293	cell_bounds=lat_bounds)
294	!
295	!
296	case('lon')
297	! lecture de la longitude:
298	allocate(lon(lendim(idim)))
299	ierr = nf_inq_varid(orig_file_id, namedim, lonid)
300	if (ierr /= 0) call handle_err(ierr)
301	ierr = nf_get_var_real(orig_file_id, lonid, lon)
302	if (ierr /= 0) call handle_err(ierr)
303	ierr = nf_get_att_text(orig_file_id, lonid, 'units', units)
304	if (ierr /= 0) call handle_err(ierr)
305	ALLOCATE(lon_bounds(lendim(idim)+1))
306	DO i = 2, lendim(idim)
307	lon_bounds(i) = lon(i-1) - (lon(i-1) - lon(i))/2
308	enddo
309	lon_bounds(1) = lon(1) - (lon_bounds(3) -lon_bounds(2))/2.
310	lon_bounds(lendim(idim)+1) = lon(lendim(idim)) + (lon_bounds(lendim(idim))-lon_bounds(lendim(idim)-1))/2.
311
312	! definition de la longitude
313	axis_ids(idim) = cmor_axis( &
314	table=trim(ipcc_table(index_table)), &
315	table_entry='longitude', &
316	units=trim(units), &
317	length=lendim(idim), &
318	coord_vals=lon, &
319	cell_bounds=lon_bounds)
320	!
321	!
322	case('presnivs')
323	! lecture de la verticale:
324	allocate(vert(lendim(idim)))
325	ierr = nf_inq_varid(orig_file_id, namedim, vertid)
326	if (ierr /= 0) call handle_err(ierr)
327	ierr = nf_get_var_real(orig_file_id, vertid, vert)
328	if (ierr /= 0) call handle_err(ierr)
329	ierr = nf_get_att_text(orig_file_id, vertid, 'units', units)
330	if (ierr /= 0) call handle_err(ierr)
331	!
332	! definition de la verticale
333	if (units == 'mb') units='Pa'
334	axis_ids(idim) = cmor_axis( &
335	table=trim(ipcc_table(index_table)), &
336	table_entry='pressure', &
337	units=trim(units), &
338	length=lendim(idim), &
339	coord_vals=vert)
340	!
341	!
342	case('time_counter')
343	! definition du temps
344	if (idim /= ndims) then
345	write(lunout,*)'la dimension temps doit etre la derniere dimension'
346	stop
347	endif
348	allocate(time(lendim(idim)))
349	ierr = nf_inq_varid(orig_file_id,namedim,timeid)
350	if (ierr /=0) call handle_err(ierr)
351	ierr = nf_get_var_real(orig_file_id, timeid, time)
352	if (ierr /=0) call handle_err(ierr)
353	ierr = nf_get_att_text(orig_file_id,timeid, 'units', units)
354	if (ierr /=0) call handle_err(ierr)
355	axis_ids(idim) = cmor_axis( &
356	table=trim(ipcc_table(index_table)), &
357	table_entry='time', &
358	units=trim(units), &
359	length=lendim(idim), &
360	interval='30 minutes')
361	itime= lendim(idim)
362	case default
363	write(lunout,*)'Dimension: ', trim(namedim),' non reconnue'
364	stop
365	endselect
366	enddo
367
368	!
369	! Definition de la variable a ecrire
370	units=' '
371	ierr = nf_inq_varid(orig_file_id,TRIM(varname), varid)
372	if (ierr /= 0) call handle_err(ierr)
373	ierr = nf_get_att_text(orig_file_id, varid, 'units', units)
374	if (ierr /= 0) call handle_err(ierr)
375	ierr = nf_get_att_real(orig_file_id, varid, 'missing_value', missing_value)
376	if (ierr /= 0) call handle_err(ierr)
377	cmorvarid = cmor_variable( &
378	table=trim(ipcc_table(index_table)), &
379	table_entry=trim(ipcc_name(index_table)), &
380	units=trim(units), &
381	axis_ids= axis_ids, &
382	missing_value=real(missing_value), &
383	positive = trim(ipsl_pos(index_table)), &
384	original_name=trim(varname))
385	!
386	! Lecture de la variable
387	if (ndims == 3) then
388	ALLOCATE (donnees(lendim(1), lendim(2), 1, lendim(3) ))
389	else if (ndims ==4) then
390	allocate (donnees(lendim(1), lendim(2), lendim(3), lendim(4) ))
391	endif
392	ierr = nf_get_var_real(orig_file_id, varid, donnees)
393	!
394	! Ecriture de la variable
395
396	DO i = 1, itime
397	rdate(1) = dble(time(i))
398	ierr = cmor_write( &
399	var_id = cmorvarid, &
400	DATA = REAL(donnees(:,:,:,i)), &
401	ntimes_passed = 1, &
402	time_vals = rdate)
403	enddo
404	!
405	! Fin CMOR
406	ierr = cmor_close()
407	IF (ierr /= 0) then
408	WRITE(lunout,*)'Probleme dans cmor_close, ierr = ', ierr
409	endif
410
411	!
412	! fermeture fichier originel
413	ierr = nf_close(orig_file_id)
414	IF (ierr /= NF_NOERR) then
415	WRITE(lunout,*)NF_STRERROR(ierr)
416	stop
417	endif
418
419	contains
420
421	subroutine handle_err(status)
422	use netcdf
423
424	implicit none
425	integer, intent(in) :: status
426
427	write(lunout,*)nf90_strerror(status)
428	stop
429
430	end subroutine handle_err
431
432	END PROGRAM ts2IPCC
433

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

Download in other formats:

Original Format