Context Navigation

ts2IPCC.F90 @ 568

Last change on this file since 568 was 561, checked in by lmdzadmin, 20 years ago
Rajout d'une colonne indiquant dans quelle direction (up/down) les variables sont positives (pour les variables pour lesquelles cela a un sens) IM/LF
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 14.3 KB

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

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

Download in other formats:

Original Format