Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

lint_datasets.F90 @ 3607

Last change on this file since 3607 was 3560, checked in by debatzbr, 7 weeks ago
Addition of the microphysics model in moments.
File size: 71.5 KB

Line
1	! Copyright (c) Université de Reims Champagnne-Ardenne (2010-2015)
2	! Contributor: Jérémie Burgalat (jeremie.burgalat@univ-reims.fr)
3	!
4	! This software is a computer program whose purpose is to compute multi-variate
5	! linear interpolation.
6	!
7	! This software is governed by the CeCILL-B license under French law and
8	! abiding by the rules of distribution of free software. You can use,
9	! modify and/ or redistribute the software under the terms of the CeCILL-B
10	! license as circulated by CEA, CNRS and INRIA at the following URL
11	! "http://www.cecill.info".
12	!
13	! As a counterpart to the access to the source code and rights to copy,
14	! modify and redistribute granted by the license, users are provided only
15	! with a limited warranty and the software's author, the holder of the
16	! economic rights, and the successive licensors have only limited
17	! liability.
18	!
19	! In this respect, the user's attention is drawn to the risks associated
20	! with loading, using, modifying and/or developing or reproducing the
21	! software by the user in light of its specific status of free software,
22	! that may mean that it is complicated to manipulate, and that also
23	! therefore means that it is reserved for developers and experienced
24	! professionals having in-depth computer knowledge. Users are therefore
25	! encouraged to load and test the software's suitability as regards their
26	! requirements in conditions enabling the security of their systems and/or
27	! data to be ensured and, more generally, to use and operate it in the
28	! same conditions as regards security.
29	!
30	! The fact that you are presently reading this means that you have had
31	! knowledge of the CeCILL-B license and that you accept its terms.
32
33	!! File: lint_datasets.F90
34	!! Summary: Dataset module definition file
35	!! Author: J. Burgalat
36	!! Date: 2014-2015
37
38	MODULE LINT_DATASETS
39	!! dataset definitions module
40	!!
41	!! This module defines simple derived types that encapsulate data set variables. For a N-dimension
42	!! dataset, N vectors of \(n_{i}\) elements and a single N-dimensional array of
43	!! \(\prod_{i}^{N} n_{i}\) elements are defined.
44	!!
45	!! If the data set is to be used within [[lint_dset(module)]] module then each coordinate values of the
46	!! dataset must be sorted either in ascending or descending order with no duplicates. The module does
47	!! not provide functionnality to sort and to check such requirements.
48	!!
49	!! The module also provides two interfaces to initialize data sets from input data file which can be
50	!! a NetCDF file (if NetCDF support is enabled at compilation time) or an ASCII file. In the latter
51	!! case, the file must contain a header which is formatted as follows:
52	!!
53	!! - The first line must contain only one value which is the number of coordinates (__N__)
54	!! - The second line must contain all the dimension sizes (that is __N__ values)
55	!! - Each other lines should contain __N+1__ columns with, for the __N__ first columns, the
56	!! coordinates values and finally the point value in the last column.
57	!!
58	!! @note
59	!! Note that for ASCII file, data must be ordered so first dimensions vary first. Same requirement
60	!! is needed for NetCDF file but in most cases, it is implicitly done (if dimensions are ordered).
61	USE NETCDF
62	IMPLICIT NONE
63
64	PRIVATE
65
66	PUBLIC :: read_dset, write_dset, clear_dset, is_in, has_data, debug
67
68	LOGICAL :: debug = .false. !! A control flag to enable verbose mode
69
70	!> Initialize a data set from either an ASCII or a NetCDF file
71	!!
72	!! Whatever the kind of file, this interface assumes that you know the dimensionnality of the
73	!! extracted variable. If file content is not compatible with the kind of data set given as
74	!! output argument, or if some I/O error occured, the dataset is cleared.
75	!! @note
76	!! Netcdf reader interface is available only if the library has been compiled with NetCDF support.
77	INTERFACE read_dset
78	MODULE PROCEDURE ncdf_rd_1d,ncdf_rd_2d,ncdf_rd_3d,ncdf_rd_4d,ncdf_rd_5d
79	#if HAVE_NC4_FTN
80	MODULE PROCEDURE ncdf4_rd_1d,ncdf4_rd_2d,ncdf4_rd_3d,ncdf4_rd_4d,ncdf4_rd_5d
81	#endif
82	MODULE PROCEDURE ascii_rd_1d,ascii_rd_2d,ascii_rd_3d,ascii_rd_4d,ascii_rd_5d
83	END INTERFACE
84
85	!> Write a dataset in a netcdf (classic) file.
86	INTERFACE write_dset
87	MODULE PROCEDURE nc_wr_1d,nc_wr_2d,nc_wr_3d,nc_wr_4d,nc_wr_5d
88	END INTERFACE
89
90	!> Clear the given data set
91	INTERFACE clear_dset
92	MODULE PROCEDURE clr_1d_set,clr_2d_set,clr_3d_set,clr_4d_set,clr_5d_set
93	END INTERFACE
94
95	!> Check if given point is within the data set
96	INTERFACE is_in
97	MODULE PROCEDURE is_in_1d,is_in_2d,is_in_3d,is_in_4d,is_in_5d
98	END INTERFACE
99
100	!> Private interface to netcdf informations getters
101	INTERFACE get_nc_info
102	MODULE PROCEDURE get_nc3_info
103	#if HAVE_NC4_FTN
104	MODULE PROCEDURE get_nc4_info
105	#endif
106	END INTERFACE
107
108	INTERFACE has_data
109	MODULE PROCEDURE has_d_1d,has_d_2d,has_d_3d,has_d_4d,has_d_5d
110	END INTERFACE
111
112
113	TYPE, PUBLIC :: DSET1D
114	!! A 1D data set
115	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: x !! X coordinate tabulated values
116	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: data !! Tabulated function's value at each coordinate
117	CHARACTER(len=NF90_MAX_NAME) :: xname = "X" !! Name of the X coordinate
118	CHARACTER(len=NF90_MAX_NAME) :: dname = "data" !! Name of the data block
119	END TYPE DSET1D
120
121	TYPE, PUBLIC :: DSET2D
122	!! A 2D data set
123	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: x !! X coordinate tabulated values
124	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: y !! Y coordinate tabulated values
125	REAL(kind=8), DIMENSION(:,:), ALLOCATABLE :: data !! Tabulated function's value at each coordinate
126	CHARACTER(len=NF90_MAX_NAME) :: xname = "X" !! Name of the X coordinate
127	CHARACTER(len=NF90_MAX_NAME) :: yname = "Y" !! Name of the Y coordinate
128	CHARACTER(len=NF90_MAX_NAME) :: dname = "data" !! Name of the data block
129	END TYPE DSET2D
130
131	TYPE, PUBLIC :: DSET3D
132	!! A 3D data set
133	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: x !! X coordinate tabulated values
134	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: y !! Y coordinate tabulated values
135	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: z !! Z coordinate tabulated values
136	REAL(kind=8), DIMENSION(:,:,:), ALLOCATABLE :: data !! Tabulated function's value at each coordinate
137	CHARACTER(len=NF90_MAX_NAME) :: xname = "X" !! Name of the X coordinate
138	CHARACTER(len=NF90_MAX_NAME) :: yname = "Y" !! Name of the Y coordinate
139	CHARACTER(len=NF90_MAX_NAME) :: zname = "Z" !! Name of the Z coordinate
140	CHARACTER(len=NF90_MAX_NAME) :: dname = "data" !! Name of the data block
141	END TYPE DSET3D
142
143	TYPE, PUBLIC :: DSET4D
144	!! A 4D data set
145	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: x !! X coordinate tabulated values
146	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: y !! Y coordinate tabulated values
147	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: z !! Z coordinate tabulated values
148	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: t !! T coordinate tabulated values
149	REAL(kind=8), DIMENSION(:,:,:,:), ALLOCATABLE :: data !! Tabulated function's value at each coordinate
150	CHARACTER(len=NF90_MAX_NAME) :: xname = "X" !! Name of the X coordinate
151	CHARACTER(len=NF90_MAX_NAME) :: yname = "Y" !! Name of the Y coordinate
152	CHARACTER(len=NF90_MAX_NAME) :: zname = "Z" !! Name of the Z coordinate
153	CHARACTER(len=NF90_MAX_NAME) :: tname = "T" !! Name of the T coordinate
154	CHARACTER(len=NF90_MAX_NAME) :: dname = "data" !! Name of the data block
155	END TYPE DSET4D
156
157	TYPE, PUBLIC :: DSET5D
158	!! A 5D data set
159	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: x !! X coordinate tabulated values
160	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: y !! Y coordinate tabulated values
161	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: z !! Z coordinate tabulated values
162	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: t !! T coordinate tabulated values
163	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: w !! W coordinate tabulated values
164	REAL(kind=8), DIMENSION(:,:,:,:,:), ALLOCATABLE :: data !! Tabulated function's value at each coordinate
165	CHARACTER(len=NF90_MAX_NAME) :: xname = "X" !! Name of the X coordinate
166	CHARACTER(len=NF90_MAX_NAME) :: yname = "Y" !! Name of the Y coordinate
167	CHARACTER(len=NF90_MAX_NAME) :: zname = "Z" !! Name of the Z coordinate
168	CHARACTER(len=NF90_MAX_NAME) :: tname = "T" !! Name of the T coordinate
169	CHARACTER(len=NF90_MAX_NAME) :: wname = "W" !! Name of the W coordinate
170	CHARACTER(len=NF90_MAX_NAME) :: dname = "data" !! Name of the data block
171	END TYPE DSET5D
172
173	CONTAINS
174
175
176	FUNCTION is_in_1d(set,x) RESULT(ret)
177	!! Check if point is in the 1D data set
178	TYPE(DSET1D), INTENT(in) :: set !! Dataset object to search in
179	REAL(kind=8), INTENT(in) :: x !! coordinate of the point to check
180	LOGICAL :: ret !! .true. if the point is in the data set, .false. otherwise
181	REAL(kind=8) :: l,u
182	ret=.true.
183	l = set%x(1) ; u= set%x(size(set%x))
184	IF ((x>=l.EQV.x<=u).OR.(x<=l.EQV.x>=u)) ret=.false.
185	RETURN
186	END FUNCTION is_in_1d
187
188	FUNCTION is_in_2d(set,x,y) RESULT(ret)
189	!! Check if point is in the 2D data set
190	TYPE(DSET2D), INTENT(in) :: set !! Dataset object to search in
191	REAL(kind=8), INTENT(in) :: x !! X coordinate of the point to check
192	REAL(kind=8), INTENT(in) :: y !! Y coordinate of the point to check
193	LOGICAL :: ret !! .true. if the point is in the data set, .false. otherwise
194	REAL(kind=8) :: l,u
195	ret=.false.
196	l = set%x(1) ; u= set%x(size(set%x))
197	IF ((x>l.AND.x>u).OR.(x<l.AND.x<u)) RETURN
198	IF ((x>l.AND.x>u).OR.(x<l.AND.x<u)) RETURN
199	l = set%y(1) ; u= set%y(size(set%y))
200	IF ((y>l.AND.y>u).OR.(y<l.AND.y<u)) RETURN
201	ret=.true.
202	RETURN
203	END FUNCTION is_in_2d
204
205	FUNCTION is_in_3d(set,x,y,z) RESULT(ret)
206	!! Check if point is in the 3D data set
207	TYPE(DSET3D), INTENT(in) :: set !! Dataset object to search in
208	REAL(kind=8), INTENT(in) :: x !! X coordinate of the point to check
209	REAL(kind=8), INTENT(in) :: y !! Y coordinate of the point to check
210	REAL(kind=8), INTENT(in) :: z !! Z coordinate of the point to check
211	LOGICAL :: ret !! .true. if the point is in the data set, .false. otherwise
212	REAL(kind=8) :: l,u
213	ret=.false.
214	l = set%x(1) ; u= set%x(size(set%x))
215	IF ((x>l.AND.x>u).OR.(x<l.AND.x<u)) RETURN
216	l = set%y(1) ; u= set%y(size(set%y))
217	IF ((y>l.AND.y>u).OR.(y<l.AND.y<u)) RETURN
218	l = set%z(1) ; u= set%z(size(set%z))
219	IF ((z>l.AND.z>u).OR.(z<l.AND.z<u)) RETURN
220	ret=.true.
221	RETURN
222	END FUNCTION is_in_3d
223
224	FUNCTION is_in_4d(set,x,y,z,t) RESULT(ret)
225	!! Check if point is in the 4D data set
226	TYPE(DSET4D), INTENT(in) :: set !! Dataset object to search in
227	REAL(kind=8), INTENT(in) :: x !! X coordinate of the point to check
228	REAL(kind=8), INTENT(in) :: y !! Y coordinate of the point to check
229	REAL(kind=8), INTENT(in) :: z !! Z coordinate of the point to check
230	REAL(kind=8), INTENT(in) :: t !! T coordinate of the point to check
231	LOGICAL :: ret !! .true. if the point is in the data set, .false. otherwise
232	REAL(kind=8) :: l,u
233	ret=.false.
234	l = set%x(1) ; u= set%x(size(set%x))
235	IF ((x>l.AND.x>u).OR.(x<l.AND.x<u)) RETURN
236	l = set%y(1) ; u= set%y(size(set%y))
237	IF ((y>l.AND.y>u).OR.(y<l.AND.y<u)) RETURN
238	l = set%z(1) ; u= set%z(size(set%z))
239	IF ((z>l.AND.z>u).OR.(z<l.AND.z<u)) RETURN
240	l = set%t(1) ; u= set%t(size(set%t))
241	IF ((t>l.AND.t>u).OR.(t<l.AND.t<u)) RETURN
242	ret=.true.
243	RETURN
244	END FUNCTION is_in_4d
245
246	FUNCTION is_in_5d(set,x,y,z,t,w) RESULT(ret)
247	!! Check if point is in the 4D data set
248	TYPE(DSET5D), INTENT(in) :: set !! Dataset object to search in
249	REAL(kind=8), INTENT(in) :: x !! X coordinate of the point to check
250	REAL(kind=8), INTENT(in) :: y !! Y coordinate of the point to check
251	REAL(kind=8), INTENT(in) :: z !! Z coordinate of the point to check
252	REAL(kind=8), INTENT(in) :: t !! T coordinate of the point to check
253	REAL(kind=8), INTENT(in) :: w !! W coordinate of the point to check
254	LOGICAL :: ret !! .true. if the point is in the data set, .false. otherwise
255	REAL(kind=8) :: l,u
256	ret=.false.
257	l = set%x(1) ; u= set%x(size(set%x))
258	IF ((x>l.AND.x>u).OR.(x<l.AND.x<u)) RETURN
259	l = set%y(1) ; u= set%y(size(set%y))
260	IF ((y>l.AND.y>u).OR.(y<l.AND.y<u)) RETURN
261	l = set%z(1) ; u= set%z(size(set%z))
262	IF ((z>l.AND.z>u).OR.(z<l.AND.z<u)) RETURN
263	l = set%t(1) ; u= set%t(size(set%t))
264	IF ((t>l.AND.t>u).OR.(t<l.AND.t<u)) RETURN
265	l = set%w(1) ; u= set%w(size(set%w))
266	IF ((w>l.AND.w>u).OR.(w<l.AND.w<u)) RETURN
267	ret=.true.
268	RETURN
269	END FUNCTION is_in_5d
270
271	FUNCTION ascii_header(path,ds,dp) RESULT(ret)
272	!! Read ASCII file header
273	!!
274	!! The method assumes the header is on two lines :
275	!!
276	!! - the first line must contain a single value which is the number of
277	!! dimensions.
278	!! - the second must contain N values with the size of each dimensions (N
279	!! referring to the first line number).
280	!!
281	!! The file remains open in the logical unit 666 unless some error occured.
282	CHARACTER(len=*), INTENT(in) :: path !! Path of the ASCII file to read
283	INTEGER, INTENT(out), DIMENSION(:) :: ds !! Size of each dimensions
284	INTEGER, INTENT(out), DIMENSION(:) :: dp !! Product of each lower dimension size (1 for 1st)
285	LOGICAL :: ret !! .true. if no error occured, .false. otherwise
286	INTEGER :: nd,i,e
287	CHARACTER(len=15) :: i2s
288	INQUIRE(666,OPENED=ret)
289	IF (ret) THEN
290	WRITE(,) 'ERROR: LUN 666 already used...'
291	ret = .false. ; RETURN
292	ENDIF
293	ret = .false.
294	OPEN(666,file=TRIM(path))
295	READ(666,*) nd
296	IF (nd /= SIZE(ds)) THEN
297	WRITE(i2s,*) nd ; i2s=ADJUSTL(i2s)
298	WRITE(*,'(a)') "ERROR: Incompatible size, DSET should be "//TRIM(i2s)//"D"
299	RETURN
300	ENDIF
301	READ(666,*,iostat=e) ds
302	IF (e /= 0) THEN
303	WRITE(*,'(a)') 'ERROR: Cannot get dimensions size'
304	CLOSE(666) ; RETURN
305	ENDIF
306	dp(1) = 1
307	DO i=2,nd
308	dp(i)=PRODUCT(ds(:i-1))
309	ENDDO
310	ret = .true.
311	END FUNCTION ascii_header
312
313	FUNCTION nc_get_dim_by_name(fid,dn,values,ds,did) RESULT(ret)
314	!! Get informations about a dimension.
315	!!
316	!! The method gets the values and size of a given dimension.
317	!!
318	!! Errors occur if:
319	!!
320	!! - the given name is not a dimension name.
321	!! - the size of the dimension is less than 2 (the method assumes dimension has an afferent variable with values).
322	!! - the method can not retrieve the values of the afferent variable.
323	INTEGER, INTENT(in) :: fid
324	!! Id of the NetCDF file
325	CHARACTER(len=NF90_MAX_NAME), INTENT(in) :: dn
326	!! Name of the dimension.
327	REAL(kind=8), INTENT(out), DIMENSION(:), ALLOCATABLE :: values
328	!! Values of the dimension
329	INTEGER, INTENT(out) :: ds
330	!! Size of __values__
331	INTEGER, INTENT(out) :: did
332	!! Id of the NetCDF dimension
333	LOGICAL :: ret
334	!! .true. if no error(s) occured, .false. otherwise
335	INTEGER :: vid,err
336	CHARACTER(len=15) :: i2s
337	ret = .false.
338	! --- Get dimension informations
339	err = NF90_INQ_DIMID(fid,dn,did)
340	IF (err /= NF90_NOERR) RETURN
341	err = NF90_INQUIRE_DIMENSION(fid,did,len=ds)
342	IF (err /= NF90_NOERR) RETURN
343	IF (ds < 2) RETURN
344	err = NF90_INQ_VARID(fid,TRIM(dn),vid)
345	IF (err /= NF90_NOERR) RETURN
346	ALLOCATE(values(ds))
347	err = NF90_GET_VAR(fid,vid,values)
348	IF (err /= NF90_NOERR) RETURN
349	ret = .true.
350	END FUNCTION nc_get_dim_by_name
351
352	FUNCTION nc_get_dim_by_id(fid,did,values,ds,dn) RESULT(ret)
353	!! Get informations about a dimension.
354	!!
355	!! The method gets the values and size of a given dimension.
356	!!
357	!! Errors occur if:
358	!!
359	!! - the given id is not a dimension identifier.
360	!! - the size of the dimension is less than 2 (the method assumes dimension has an afferent variable with values).
361	!! - the method can not retrieve the values of the afferent variable.
362	INTEGER, INTENT(in) :: fid
363	!! Id of the NetCDF file
364	INTEGER, INTENT(in) :: did
365	!! Id of the NetCDF dimension
366	REAL(kind=8), INTENT(out), DIMENSION(:), ALLOCATABLE :: values
367	!! Values of the dimension
368	INTEGER, INTENT(out) :: ds
369	!! Size of __values__
370	CHARACTER(len=NF90_MAX_NAME), INTENT(out) :: dn
371	!! Name of the dimension.
372	LOGICAL :: ret
373	!! .true. if no error(s) occured, .false. otherwise
374	INTEGER :: vid,err
375	CHARACTER(len=15) :: i2s
376	ret = .false.
377	! --- Get dimension informations
378	IF (NF90_INQUIRE_DIMENSION(fid,did,dn,ds) /= NF90_NOERR) RETURN
379	IF (ds < 2) RETURN
380	IF (NF90_INQ_VARID(fid,TRIM(dn),vid) /= NF90_NOERR) RETURN
381	ALLOCATE(values(ds))
382	IF (NF90_GET_VAR(fid,vid,values) /= NF90_NOERR) RETURN
383	ret = .true.
384	END FUNCTION nc_get_dim_by_id
385
386	FUNCTION get_nc3_info(path,variable,ncid,vid,dimids,verbose) RESULT(ret)
387	!! Get variable informations from NetCDF file
388	!!
389	!! The method attempts to read the given NetCDF file header and retrieves
390	!! variable's informations from it:
391	!!
392	!! - parent file ID
393	!! - variable ID
394	!! - variable's dimensions ID
395	!! - variable's dimensions sizes
396	!!
397	!! The method always opens and closes the file. If the file cannot be opened,
398	!! __ncid__ is set to -1.
399	CHARACTER(len=*), INTENT(in) :: path
400	!! Path of the NetCDF file
401	CHARACTER(len=*), INTENT(in) :: variable
402	!! Name of the NetCDF variable to extract
403	INTEGER, INTENT(out) :: ncid
404	!! NetCDF file ID
405	INTEGER, INTENT(out) :: vid
406	!! NetCDF Variable ID
407	INTEGER, INTENT(out), DIMENSION(:), ALLOCATABLE :: dimids
408	!! Id of the variable's dimensions. __dimids__ is not allocated on error.
409	LOGICAL, INTENT(in), OPTIONAL :: verbose
410	!! True to print out message on error (default to False).
411	LOGICAL :: ret
412	!! .true. if no errors occured, .false. otherwise
413	INTEGER :: fid,ty,nc,err
414	LOGICAL :: zlog
415	zlog = .false. ; IF (PRESENT(verbose)) zlog = verbose
416	ret = .false.
417	ncid = -1
418	! Opens file
419	IF (NF90_OPEN(TRIM(path),NF90_NOWRITE,fid) /= NF90_NOERR) THEN
420	IF (zlog) WRITE(*,'(a)') 'ERROR: Cannot open '//trim(path)
421	RETURN
422	ENDIF
423	ncid = fid
424	! Searches for variable
425	IF (NF90_INQ_VARID(ncid,TRIM(variable),vid) /= NF90_NOERR) THEN
426	IF (zlog) WRITE(*,'(a)') 'Cannot find '//TRIM(variable)//' in '//trim(path)
427	nc = NF90_CLOSE(fid)
428	RETURN
429	ENDIF
430	! Get variable infos
431	! 1st call to get type and number of dimensions)
432	IF (NF90_INQUIRE_VARIABLE(ncid,vid,xtype=ty,ndims=nc) /= NF90_NOERR) THEN
433	IF (zlog) WRITE(*,'(a)') 'Cannot access to '//TRIM(variable)//' informations'
434	nc = NF90_CLOSE(fid)
435	RETURN
436	ELSE
437	! Checks type
438	IF (ty == NF90_CHAR) THEN
439	IF (zlog) WRITE(*,'(a)') 'Inconsistent variable type (should be numeric)'
440	nc = NF90_CLOSE(fid)
441	RETURN
442	ENDIF
443	ALLOCATE(dimids(nc))
444	ENDIF
445	! Gets variable's dimensions informations
446	! first get dimensions id
447	IF (NF90_INQUIRE_VARIABLE(ncid,vid,dimids=dimids) /= NF90_NOERR) THEN
448	IF (zlog) WRITE(*,'(a)') 'Cannot access to '//TRIM(variable)//' informations'
449	nc = NF90_CLOSE(fid)
450	DEALLOCATE(dimids)
451	RETURN
452	ENDIF
453	nc = NF90_CLOSE(fid)
454	ret = .true.
455	END FUNCTION get_nc3_info
456
457	#if HAVE_NC4_FTN
458	FUNCTION get_nc4_info(path,variable,group,ncid,vid,dimids,verbose) RESULT(ret)
459	!! Get variable informations from NetCDF4 file
460	!!
461	!! The method attempts to read the given NetCDF file header and retrieves
462	!! variable's informations from it :
463	!!
464	!! - parent group/file ID
465	!! - variable ID
466	!! - variable's dimensions ID
467	!! - variable's dimensions sizes
468	!!
469	!! The method always opens and closes the file. If the file cannot be opened,
470	!! __ncid__ is set to -1.
471	CHARACTER(len=*), INTENT(in) :: path
472	!! Path of the NetCDF file
473	CHARACTER(len=*), INTENT(in) :: variable
474	!! Name of the NetCDF variable to extract
475	CHARACTER(len=*), INTENT(in) :: group
476	!! Fullname of the variable's parent group (should be set to empty string for root group)
477	INTEGER, INTENT(out) :: ncid
478	!! NetCDF group ID
479	INTEGER, INTENT(out) :: vid
480	!! NetCDF Variable ID
481	INTEGER, INTENT(out), DIMENSION(:), ALLOCATABLE :: dimids
482	!! Id of the variable's dimensions. On error, __dimids__ is not allocated.
483	LOGICAL, INTENT(in), OPTIONAL :: verbose
484	!! True to print out message on error (default to False).
485	LOGICAL :: ret
486	!! .true. if no errors occured, .false. otherwise
487	INTEGER :: fid,ty,nc,err
488	LOGICAL :: zlog
489	zlog = .false. ; IF (PRESENT(verbose)) zlog = verbose
490	ret = .false.
491	ncid = -1
492	! Opens file
493	IF (NF90_OPEN(TRIM(path),NF90_NOWRITE,fid) /= NF90_NOERR) THEN
494	IF (zlog) WRITE(*,'(a)') 'ERROR: Cannot open '//trim(path)
495	RETURN
496	ENDIF
497	! Searches for group based on its fullname
498	IF (LEN_TRIM(group) == 0) THEN
499	nc = NF90_NOERR ; ncid = fid
500	ELSE
501	nc = NF90_INQ_GRP_FULL_NCID(fid,TRIM(group),ncid)
502	ENDIF
503	SELECT CASE(nc)
504	! NF90_ENOGRP is missing from Netcdf-Fortran4.2 : its value=-125
505	CASE(-125)
506	IF (zlog) WRITE(*,'(a)') TRIM(group)//' does not exist in '//TRIM(path)
507	nc = NF90_CLOSE(fid) ; RETURN
508	CASE(NF90_ENOTNC4,NF90_ESTRICTNC3)
509	IF (zlog) WRITE(*,'(a)') TRIM(path)//' is not a NetCDF-4 file with "group" feature'
510	nc = NF90_CLOSE(fid) ; RETURN
511	CASE(NF90_EHDFERR)
512	IF (zlog) WRITE(*,'(a)') "Too bad, an HDF5 error has been reported..."
513	nc = NF90_CLOSE(fid) ; RETURN
514	END SELECT
515	! Searches for variable
516	IF (NF90_INQ_VARID(ncid,TRIM(variable),vid) /= NF90_NOERR) THEN
517	IF (zlog) WRITE(*,'(a)') 'Cannot find '//TRIM(variable)//' in '//trim(path)
518	nc = NF90_CLOSE(fid)
519	RETURN
520	ENDIF
521	! Get variable infos
522	! 1st call to get type and number of dimensions)
523	IF (NF90_INQUIRE_VARIABLE(ncid,vid,xtype=ty,ndims=nc) /= NF90_NOERR) THEN
524	IF (zlog) WRITE(*,'(a)') 'Cannot access to '//TRIM(variable)//' informations'
525	nc = NF90_CLOSE(fid)
526	RETURN
527	ELSE
528	! Checks type
529	IF (ty == NF90_CHAR) THEN
530	IF (zlog) WRITE(*,'(a)') 'Inconsistent variable type (should be numeric)'
531	nc = NF90_CLOSE(fid)
532	RETURN
533	ENDIF
534	ALLOCATE(dimids(nc))
535	ENDIF
536	! Gets variable's dimensions informations
537	! first get dimensions id
538	IF (NF90_INQUIRE_VARIABLE(ncid,vid,dimids=dimids) /= NF90_NOERR) THEN
539	IF (zlog) WRITE(*,'(a)') 'Cannot access to '//TRIM(variable)//' informations'
540	nc = NF90_CLOSE(fid)
541	DEALLOCATE(dimids)
542	RETURN
543	ENDIF
544	nc = NF90_CLOSE(fid)
545	ret = .true.
546	END FUNCTION get_nc4_info
547	#endif
548
549	!-------------------------
550	! NetCDF data file readers
551	!-------------------------
552
553	FUNCTION ncdf_rd_1d(path,variable,set) RESULT(ret)
554	!! Read a 1D data set from a NetCDF file
555	CHARACTER(len=*), INTENT(in) :: path !! Path of the NetCDF file
556	CHARACTER(len=*), INTENT(in) :: variable !! Name of the NetCDF variable to extract
557	TYPE(DSET1D), INTENT(out) :: set !! Output dataset object
558	LOGICAL :: ret !! .true. if no errors occured, .false. otherwise
559	INTEGER :: fi,vi,nd,iret
560	INTEGER, DIMENSION(:), ALLOCATABLE :: di,ds
561	CHARACTER(len=NF90_MAX_NAME) :: dn
562	CHARACTER(len=15) :: i2s
563	ret = .false.
564	! --- Reset the data set
565	CALL clear_dset(set)
566	! --- Check NetCDF file info
567	IF (.NOT.get_nc_info(path,variable,fi,vi,di)) RETURN
568	iret = NF90_OPEN(path,NF90_NOWRITE,fi)
569	! --- Check dimension size
570	nd = SIZE(di)
571	IF (nd /= 1) THEN
572	IF (debug) THEN
573	WRITE(i2s,*) nd ; i2s=ADJUSTL(i2s)
574	WRITE(*,'(a)') "ERROR: Incompatible size, DSET should be"//TRIM(i2s)//"D"
575	ENDIF
576	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
577	ENDIF
578	! --- Get coordinates values
579	ALLOCATE(ds(nd))
580	! ------ X coordinate
581	IF (.NOT.nc_get_dim_by_id(fi,di(1),set%x,ds(1),set%xname)) THEN
582	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
583	ENDIF
584	! --- Read data
585	ALLOCATE(set%data(ds(1)))
586	IF (NF90_GET_VAR(fi,vi,set%data) /= NF90_NOERR) THEN
587	IF (debug) WRITE(*,'(a)') "ERROR:"//TRIM(variable)//": Cannot get values"
588	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
589	ENDIF
590	nd = NF90_CLOSE(fi)
591	set%dname = variable
592	ret = .true.
593	RETURN
594	END FUNCTION ncdf_rd_1d
595
596	FUNCTION ncdf_rd_2d(path,variable,set) RESULT(ret)
597	!! Read a 2D data set from a NetCDF file
598	CHARACTER(len=*), INTENT(in) :: path !! Path of the NetCDF file
599	CHARACTER(len=*), INTENT(in) :: variable !! Name of the NetCDF variable to extract
600	TYPE(DSET2D), INTENT(out) :: set !! Output dataset object
601	LOGICAL :: ret !! .true. if no errors occured, .false. otherwise
602	INTEGER :: fi,vi,nd,iret
603	INTEGER, DIMENSION(:), ALLOCATABLE :: di,ds
604	CHARACTER(len=15) :: i2s
605	ret = .false.
606	! --- Reset the data set
607	CALL clear_dset(set)
608	! --- Check NetCDF file info
609	IF (.NOT.get_nc_info(path,variable,fi,vi,di)) RETURN
610	iret = NF90_OPEN(path,NF90_NOWRITE,fi)
611	! --- Check dimension size
612	nd = SIZE(di)
613	IF (nd /= 2) THEN
614	IF (debug) THEN
615	WRITE(i2s,*) nd ; i2s=ADJUSTL(i2s)
616	WRITE(*,'(a)') "ERROR: Incompatible size, DSET should be"//TRIM(i2s)//"D"
617	ENDIF
618	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
619	ENDIF
620	! --- Get dimensions informations
621	ALLOCATE(ds(nd))
622	! ------ X coordinate
623	IF (.NOT.nc_get_dim_by_id(fi,di(1),set%x,ds(1),set%xname)) THEN
624	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
625	ENDIF
626	! ------ Y coordinate
627	IF (.NOT.nc_get_dim_by_id(fi,di(2),set%y,ds(2),set%yname)) THEN
628	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
629	ENDIF
630	! --- Read data
631	ALLOCATE(set%data(ds(1),ds(2)))
632	IF (NF90_GET_VAR(fi,vi,set%data) /= NF90_NOERR) THEN
633	IF (debug) WRITE(*,'(a)') "ERROR:"//TRIM(variable)//": Cannot get values"
634	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
635	ENDIF
636	nd = NF90_CLOSE(fi)
637	set%dname = variable
638	ret = .true.
639	RETURN
640	END FUNCTION ncdf_rd_2d
641
642	FUNCTION ncdf_rd_3d(path,variable,set) RESULT(ret)
643	!! Read a 3D data set from a NetCDF file
644	CHARACTER(len=*), INTENT(in) :: path !! Path of the NetCDF file
645	CHARACTER(len=*), INTENT(in) :: variable !! Name of the NetCDF variable to extract
646	TYPE(DSET3D), INTENT(out) :: set !! Output dataset object
647	LOGICAL :: ret !! .true. if no errors occured, .false. otherwise
648	INTEGER :: fi,vi,nd,iret
649	INTEGER, DIMENSION(:), ALLOCATABLE :: di,ds
650	CHARACTER(len=15) :: i2s
651	ret = .false.
652	! --- Reset the data set
653	CALL clear_dset(set)
654	! --- Check NetCDF file info
655	IF (.NOT.get_nc_info(path,variable,fi,vi,di)) RETURN
656	iret = NF90_OPEN(path,NF90_NOWRITE,fi)
657	! --- Check dimension size
658	nd = SIZE(di)
659	IF (nd /= 3) THEN
660	IF (debug) THEN
661	WRITE(i2s,*) nd ; i2s=ADJUSTL(i2s)
662	WRITE(*,'(a)') "ERROR: Incompatible size, DSET should be"//TRIM(i2s)//"D"
663	ENDIF
664	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
665	ENDIF
666	! --- Get dimensions informations
667	ALLOCATE(ds(nd))
668	! ------ X coordinate
669	IF (.NOT.nc_get_dim_by_id(fi,di(1),set%x,ds(1),set%xname)) THEN
670	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
671	ENDIF
672	! ------ Y coordinate
673	IF (.NOT.nc_get_dim_by_id(fi,di(2),set%y,ds(2),set%yname)) THEN
674	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
675	ENDIF
676	! ------ Z coordinate
677	IF (.NOT.nc_get_dim_by_id(fi,di(3),set%z,ds(3),set%zname)) THEN
678	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
679	ENDIF
680	! --- Read data
681	ALLOCATE(set%data(ds(1),ds(2),ds(3)))
682	IF (NF90_GET_VAR(fi,vi,set%data) /= NF90_NOERR) THEN
683	IF (debug) WRITE(*,'(a)') "ERROR:"//TRIM(variable)//": Cannot get values"
684	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
685	ENDIF
686	nd = NF90_CLOSE(fi)
687	set%dname = variable
688	ret = .true.
689	RETURN
690	END FUNCTION ncdf_rd_3d
691
692	FUNCTION ncdf_rd_4d(path,variable,set) RESULT(ret)
693	!! Read a 4D data set from a NetCDF file
694	CHARACTER(len=*), INTENT(in) :: path !! Path of the NetCDF file
695	CHARACTER(len=*), INTENT(in) :: variable !! Name of the NetCDF variable to extract
696	TYPE(DSET4D), INTENT(out) :: set !! Output dataset object
697	LOGICAL :: ret !! .true. if no errors occured, .false. otherwise
698	INTEGER :: fi,vi,nd,iret
699	INTEGER, DIMENSION(:), ALLOCATABLE :: di,ds
700	CHARACTER(len=15) :: i2s
701	ret = .false.
702	! --- Reset the data set
703	CALL clear_dset(set)
704	! --- Check NetCDF file info
705	IF (.NOT.get_nc_info(path,variable,fi,vi,di)) RETURN
706	iret = NF90_OPEN(path,NF90_NOWRITE,fi)
707	! --- Check dimension size
708	nd = SIZE(di)
709	IF (nd /= 4) THEN
710	IF (debug) THEN
711	WRITE(i2s,*) nd ; i2s=ADJUSTL(i2s)
712	WRITE(*,'(a)') "ERROR: Incompatible size, DSET should be"//TRIM(i2s)//"D"
713	ENDIF
714	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
715	ENDIF
716	! --- Get dimensions informations
717	ALLOCATE(ds(nd))
718	! ------ X coordinate
719	IF (.NOT.nc_get_dim_by_id(fi,di(1),set%x,ds(1),set%xname)) THEN
720	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
721	ENDIF
722	! ------ Y coordinate
723	IF (.NOT.nc_get_dim_by_id(fi,di(2),set%y,ds(2),set%yname)) THEN
724	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
725	ENDIF
726	! ------ Z coordinate
727	IF (.NOT.nc_get_dim_by_id(fi,di(3),set%z,ds(3),set%zname)) THEN
728	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
729	ENDIF
730	! ------ T coordinate
731	IF (.NOT.nc_get_dim_by_id(fi,di(4),set%t,ds(4),set%tname)) THEN
732	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
733	ENDIF
734	! --- Read data
735	ALLOCATE(set%data(ds(1),ds(2),ds(3),ds(4)))
736	IF (NF90_GET_VAR(fi,vi,set%data) /= NF90_NOERR) THEN
737	IF (debug) WRITE(*,'(a)') "ERROR:"//TRIM(variable)//": Cannot get values"
738	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
739	ENDIF
740	nd = NF90_CLOSE(fi)
741	set%dname = variable
742	ret = .true.
743	RETURN
744	END FUNCTION ncdf_rd_4d
745
746	FUNCTION ncdf_rd_5d(path,variable,set) RESULT(ret)
747	!! Read a 5D data set from a NetCDF file
748	CHARACTER(len=*), INTENT(in) :: path !! Path of the NetCDF file
749	CHARACTER(len=*), INTENT(in) :: variable !! Name of the NetCDF variable to extract
750	TYPE(DSET5D), INTENT(out) :: set !! Output dataset object
751	LOGICAL :: ret !! .true. if no errors occured, .false. otherwise
752	INTEGER :: fi,vi,nd,iret
753	INTEGER, DIMENSION(:), ALLOCATABLE :: di,ds
754	CHARACTER(len=15) :: i2s
755	ret = .false.
756	! --- Reset the data set
757	CALL clear_dset(set)
758	! --- Check NetCDF file info
759	IF (.NOT.get_nc_info(path,variable,fi,vi,di)) RETURN
760	iret = NF90_OPEN(path,NF90_NOWRITE,fi)
761	! --- Check dimension size
762	nd = SIZE(di)
763	IF (nd /= 5) THEN
764	IF (debug) THEN
765	WRITE(i2s,*) nd ; i2s=ADJUSTL(i2s)
766	WRITE(*,'(a)') "ERROR: Incompatible size, DSET should be"//TRIM(i2s)//"D"
767	ENDIF
768	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
769	ENDIF
770	! --- Get dimensions informations
771	ALLOCATE(ds(nd))
772	! ------ X coordinate
773	IF (.NOT.nc_get_dim_by_id(fi,di(1),set%x,ds(1),set%xname)) THEN
774	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
775	ENDIF
776	! ------ Y coordinate
777	IF (.NOT.nc_get_dim_by_id(fi,di(2),set%y,ds(2),set%yname)) THEN
778	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
779	ENDIF
780	! ------ Z coordinate
781	IF (.NOT.nc_get_dim_by_id(fi,di(3),set%z,ds(3),set%zname)) THEN
782	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
783	ENDIF
784	! ------ T coordinate
785	IF (.NOT.nc_get_dim_by_id(fi,di(4),set%t,ds(4),set%tname)) THEN
786	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
787	ENDIF
788	! ------ W coordinate
789	IF (.NOT.nc_get_dim_by_id(fi,di(5),set%w,ds(5),set%wname)) THEN
790	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
791	ENDIF
792	! --- Read data
793	ALLOCATE(set%data(ds(1),ds(2),ds(3),ds(4),ds(5)))
794	IF (NF90_GET_VAR(fi,vi,set%data) /= NF90_NOERR) THEN
795	IF (debug) WRITE(*,'(a)') "ERROR:"//TRIM(variable)//": Cannot get values"
796	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
797	ENDIF
798	nd = NF90_CLOSE(fi)
799	set%dname = variable
800	ret = .true.
801	RETURN
802	END FUNCTION ncdf_rd_5d
803
804	! NetCDF4 with groups !
805	#if HAVE_NC4_FTN
806
807	FUNCTION ncdf4_rd_1d(path,group,variable,set) RESULT(ret)
808	!! Read a 1D data set from a NetCDF4 file
809	CHARACTER(len=*), INTENT(in) :: path !! Path of the NetCDF file
810	CHARACTER(len=*), INTENT(in) :: group !! Parent group of the variable (empty for root)
811	CHARACTER(len=*), INTENT(in) :: variable !! Name of the NetCDF variable to extract
812	TYPE(DSET1D), INTENT(out) :: set !! Output dataset
813	LOGICAL :: ret !! .true. if no errors occured, .false. otherwise
814	INTEGER :: fi,vi,nd,iret
815	INTEGER, DIMENSION(:), ALLOCATABLE :: di,ds
816	CHARACTER(len=NF90_MAX_NAME) :: dn
817	CHARACTER(len=15) :: i2s
818	ret = .false.
819	! --- Reset the data set
820	CALL clear_dset(set)
821	! --- Check NetCDF file info
822	IF (.NOT.get_nc_info(path,group,variable,fi,vi,di)) RETURN
823	iret = NF90_OPEN(path,NF90_NOWRITE,fi)
824	! --- Check dimension size
825	nd = SIZE(di)
826	IF (nd /= 1) THEN
827	IF (debug) THEN
828	WRITE(i2s,*) nd ; i2s=ADJUSTL(i2s)
829	WRITE(*,'(a)') "ERROR: Incompatible size, DSET should be"//TRIM(i2s)//"D"
830	ENDIF
831	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
832	ENDIF
833	! --- Get coordinates values
834	ALLOCATE(ds(nd))
835	! ------ X coordinate
836	IF (.NOT.nc_get_dim_by_id(fi,di(1),set%x,ds(1),set%xname)) THEN
837	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
838	ENDIF
839	! --- Read data
840	ALLOCATE(set%data(ds(1)))
841	IF (NF90_GET_VAR(fi,vi,set%data) /= NF90_NOERR) THEN
842	IF (debug) WRITE(*,'(a)') "ERROR:"//TRIM(variable)//": Cannot get values"
843	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
844	ENDIF
845	nd = NF90_CLOSE(fi)
846	set%dname = variable
847	ret = .true.
848	RETURN
849	END FUNCTION ncdf4_rd_1d
850
851	FUNCTION ncdf4_rd_2d(path,group,variable,set) RESULT(ret)
852	!! Read a 2D data set from a NetCDF4 file
853	CHARACTER(len=*), INTENT(in) :: path !! Path of the NetCDF file
854	CHARACTER(len=*), INTENT(in) :: group !! Parent group of the variable (empty for root)
855	CHARACTER(len=*), INTENT(in) :: variable !! Name of the NetCDF variable to extract
856	TYPE(DSET2D), INTENT(out) :: set !! Output dataset
857	LOGICAL :: ret !! .true. if no errors occured, .false. otherwise
858	INTEGER :: fi,vi,nd,iret
859	INTEGER, DIMENSION(:), ALLOCATABLE :: di,ds
860	CHARACTER(len=15) :: i2s
861	ret = .false.
862	! --- Reset the data set
863	CALL clear_dset(set)
864	! --- Check NetCDF file info
865	IF (.NOT.get_nc_info(path,group,variable,fi,vi,di)) rETURN
866	iret = NF90_OPEN(path,NF90_NOWRITE,fi)
867	! --- Check dimension size
868	nd = SIZE(di)
869	IF (nd /= 2) THEN
870	IF (debug) THEN
871	WRITE(i2s,*) nd ; i2s=ADJUSTL(i2s)
872	WRITE(*,'(a)') "ERROR: Incompatible size, DSET should be"//TRIM(i2s)//"D"
873	ENDIF
874	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
875	ENDIF
876	! --- Get dimensions informations
877	ALLOCATE(ds(nd))
878	! ------ X coordinate
879	IF (.NOT.nc_get_dim_by_id(fi,di(1),set%x,ds(1),set%xname)) THEN
880	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
881	ENDIF
882	! ------ Y coordinate
883	IF (.NOT.nc_get_dim_by_id(fi,di(2),set%y,ds(2),set%yname)) THEN
884	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
885	ENDIF
886	! --- Read data
887	ALLOCATE(set%data(ds(1),ds(2)))
888	IF (NF90_GET_VAR(fi,vi,set%data) /= NF90_NOERR) THEN
889	IF (debug) WRITE(*,'(a)') "ERROR:"//TRIM(variable)//": Cannot get values"
890	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
891	ENDIF
892	nd = NF90_CLOSE(fi)
893	set%dname = variable
894	ret = .true.
895	RETURN
896	END FUNCTION ncdf4_rd_2d
897
898	FUNCTION ncdf4_rd_3d(path,group,variable,set) RESULT(ret)
899	!! Read a 3D data set from a NetCDF4 file
900	CHARACTER(len=*), INTENT(in) :: path !! Path of the NetCDF file
901	CHARACTER(len=*), INTENT(in) :: group !! Parent group of the variable (empty for root)
902	CHARACTER(len=*), INTENT(in) :: variable !! Name of the NetCDF variable to extract
903	TYPE(DSET3D), INTENT(out) :: set !! Output dataset
904	LOGICAL :: ret !! .true. if no errors occured, .false. otherwise
905	INTEGER :: fi,vi,nd,iret
906	INTEGER, DIMENSION(:), ALLOCATABLE :: di,ds
907	CHARACTER(len=15) :: i2s
908	ret = .false.
909	! --- Reset the data set
910	CALL clear_dset(set)
911	! --- Check NetCDF file info
912	IF (.NOT.get_nc_info(path,group,variable,fi,vi,di)) RETURN
913	iret = NF90_OPEN(path,NF90_NOWRITE,fi)
914	! --- Check dimension size
915	nd = SIZE(di)
916	IF (nd /= 3) THEN
917	IF (debug) THEN
918	WRITE(i2s,*) nd ; i2s=ADJUSTL(i2s)
919	WRITE(*,'(a)') "ERROR: Incompatible size, DSET should be"//TRIM(i2s)//"D"
920	ENDIF
921	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
922	ENDIF
923	! --- Get dimensions informations
924	ALLOCATE(ds(nd))
925	! ------ X coordinate
926	IF (.NOT.nc_get_dim_by_id(fi,di(1),set%x,ds(1),set%xname)) THEN
927	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
928	ENDIF
929	! ------ Y coordinate
930	IF (.NOT.nc_get_dim_by_id(fi,di(2),set%y,ds(2),set%yname)) THEN
931	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
932	ENDIF
933	! ------ Z coordinate
934	IF (.NOT.nc_get_dim_by_id(fi,di(3),set%z,ds(3),set%zname)) THEN
935	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
936	ENDIF
937	! --- Read data
938	ALLOCATE(set%data(ds(1),ds(2),ds(3)))
939	IF (NF90_GET_VAR(fi,vi,set%data) /= NF90_NOERR) THEN
940	IF (debug) WRITE(*,'(a)') "ERROR:"//TRIM(variable)//": Cannot get values"
941	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
942	ENDIF
943	nd = NF90_CLOSE(fi)
944	set%dname = variable
945	ret = .true.
946	RETURN
947	END FUNCTION ncdf4_rd_3d
948
949	FUNCTION ncdf4_rd_4d(path,group,variable,set) RESULT(ret)
950	!! Read a 4D data set from a NetCDF4 file
951	CHARACTER(len=*), INTENT(in) :: path !! Path of the NetCDF file
952	CHARACTER(len=*), INTENT(in) :: group !! Parent group of the variable (empty for root)
953	CHARACTER(len=*), INTENT(in) :: variable !! Name of the NetCDF variable to extract
954	TYPE(DSET4D), INTENT(out) :: set !! Output dataset
955	LOGICAL :: ret !! .true. if no errors occured, .false. otherwise
956	INTEGER :: fi,vi,nd,iret
957	INTEGER, DIMENSION(:), ALLOCATABLE :: di,ds
958	CHARACTER(len=15) :: i2s
959	ret = .false.
960	! --- Reset the data set
961	CALL clear_dset(set)
962	! --- Check NetCDF file info
963	IF (.NOT.get_nc_info(path,group,variable,fi,vi,di)) RETURN
964	iret = NF90_OPEN(path,NF90_NOWRITE,fi)
965	! --- Check dimension size
966	nd = SIZE(di)
967	IF (nd /= 4) THEN
968	IF (debug) THEN
969	WRITE(i2s,*) nd ; i2s=ADJUSTL(i2s)
970	WRITE(*,'(a)') "ERROR: Incompatible size, DSET should be"//TRIM(i2s)//"D"
971	ENDIF
972	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
973	ENDIF
974	! --- Get dimensions informations
975	ALLOCATE(ds(nd))
976	! ------ X coordinate
977	IF (.NOT.nc_get_dim_by_id(fi,di(1),set%x,ds(1),set%xname)) THEN
978	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
979	ENDIF
980	! ------ Y coordinate
981	IF (.NOT.nc_get_dim_by_id(fi,di(2),set%y,ds(2),set%yname)) THEN
982	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
983	ENDIF
984	! ------ Z coordinate
985	IF (.NOT.nc_get_dim_by_id(fi,di(3),set%z,ds(3),set%zname)) THEN
986	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
987	ENDIF
988	! ------ T coordinate
989	IF (.NOT.nc_get_dim_by_id(fi,di(4),set%t,ds(4),set%tname)) THEN
990	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
991	ENDIF
992	! --- Read data
993	ALLOCATE(set%data(ds(1),ds(2),ds(3),ds(4)))
994	IF (NF90_GET_VAR(fi,vi,set%data) /= NF90_NOERR) THEN
995	IF (debug) WRITE(*,'(a)') "ERROR:"//TRIM(variable)//": Cannot get values"
996	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
997	ENDIF
998	nd = NF90_CLOSE(fi)
999	set%dname = variable
1000	ret = .true.
1001	RETURN
1002	END FUNCTION ncdf4_rd_4d
1003
1004	FUNCTION ncdf4_rd_5d(path,group,variable,set) RESULT(ret)
1005	!! Read a 5D data set from a NetCDF4 file
1006	CHARACTER(len=*), INTENT(in) :: path !! Path of the NetCDF file
1007	CHARACTER(len=*), INTENT(in) :: group !! Parent group of the variable (empty for root)
1008	CHARACTER(len=*), INTENT(in) :: variable !! Name of the NetCDF variable to extract
1009	TYPE(DSET5D), INTENT(out) :: set !! Output dataset
1010	LOGICAL :: ret !! .true. if no errors occured, .false. otherwise
1011	INTEGER :: fi,vi,nd,iret
1012	INTEGER, DIMENSION(:), ALLOCATABLE :: di,ds
1013	CHARACTER(len=15) :: i2s
1014	ret = .false.
1015	! --- Reset the data set
1016	CALL clear_dset(set)
1017	! --- Check NetCDF file info
1018	IF (.NOT.get_nc_info(path,group,variable,fi,vi,di)) RETURN
1019	iret = NF90_OPEN(path,NF90_NOWRITE,fi)
1020	! --- Check dimension size
1021	nd = SIZE(di)
1022	IF (nd /= 5) THEN
1023	IF (debug) THEN
1024	WRITE(i2s,*) nd ; i2s=ADJUSTL(i2s)
1025	WRITE(*,'(a)') "ERROR: Incompatible size, DSET should be"//TRIM(i2s)//"D"
1026	ENDIF
1027	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
1028	ENDIF
1029	! --- Get dimensions informations
1030	ALLOCATE(ds(nd))
1031	! ------ X coordinate
1032	IF (.NOT.nc_get_dim_by_id(fi,di(1),set%x,ds(1),set%xname)) THEN
1033	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
1034	ENDIF
1035	! ------ Y coordinate
1036	IF (.NOT.nc_get_dim_by_id(fi,di(2),set%y,ds(2),set%yname)) THEN
1037	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
1038	ENDIF
1039	! ------ Z coordinate
1040	IF (.NOT.nc_get_dim_by_id(fi,di(3),set%z,ds(3),set%zname)) THEN
1041	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
1042	ENDIF
1043	! ------ T coordinate
1044	IF (.NOT.nc_get_dim_by_id(fi,di(4),set%t,ds(4),set%tname)) THEN
1045	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
1046	ENDIF
1047	! ------ W coordinate
1048	IF (.NOT.nc_get_dim_by_id(fi,di(5),set%w,ds(5),set%wname)) THEN
1049	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
1050	CALL clear_dset(set) ; RETURN
1051	ENDIF
1052	! --- Read data
1053	ALLOCATE(set%data(ds(1),ds(2),ds(3),ds(4),ds(5)))
1054	IF (NF90_GET_VAR(fi,vi,set%data) /= NF90_NOERR) THEN
1055	IF (debug) WRITE(*,'(a)') "ERROR:"//TRIM(variable)//": Cannot get values"
1056	nd = NF90_CLOSE(fi) ; CALL clear_dset(set) ; RETURN
1057	ENDIF
1058	nd = NF90_CLOSE(fi)
1059	set%dname = variable
1060	ret = .true.
1061	RETURN
1062	END FUNCTION ncdf4_rd_5d
1063	#endif
1064
1065	!------------------------
1066	! ASCII data file readers
1067	!------------------------
1068
1069	FUNCTION ascii_rd_1d(path,set) RESULT(ret)
1070	!! Read a 1D data set from a ASCII file
1071	CHARACTER(len=*), INTENT(in) :: path !! Path of the ASCII file to read
1072	TYPE(dset1d), INTENT(out) :: set !! output 1D dataset
1073	LOGICAL :: ret !! .true. if no error occured, .false. otherwise
1074	INTEGER :: i,e
1075	REAL(kind=8), DIMENSION(2) :: vl
1076	INTEGER, DIMENSION(1) :: cc,ds,dp
1077	ret = .false.
1078	CALL clear_dset(set)
1079	IF (.NOT.ascii_header(path,ds,dp)) RETURN
1080	ALLOCATE(set%x(ds(1)), &
1081	set%data(ds(1)))
1082	cc(:) = 1
1083	DO i=1,ds(1)
1084	READ(666,*) vl(:)
1085	set%x(cc(1)) = vl(1)
1086	set%data(cc(1)) = vl(2)
1087	cc(1) = cc(1) + 1
1088	ENDDO
1089	READ(666,*,iostat=e) vl(1)
1090	IF (e == 0) THEN
1091	IF (debug) WRITE(,) 'ERROR: Extra value found...'
1092	ret = .false. ; call clear_dset(set) ; CLOSE(666) ; RETURN
1093	ENDIF
1094	ret = .true.
1095	CLOSE(666)
1096	RETURN
1097	END FUNCTION ascii_rd_1d
1098
1099	FUNCTION ascii_rd_2d(path,set) RESULT(ret)
1100	!! Read a 2D data set from a ASCII file
1101	CHARACTER(len=*), INTENT(in) :: path !! Path of the ASCII file to read
1102	TYPE(dset2d), INTENT(out) :: set !! output 2D dataset
1103	LOGICAL :: ret !! .true. if no error occured, .false. otherwise
1104	INTEGER :: i,e
1105	REAL(kind=8), DIMENSION(3) :: vl
1106	INTEGER, DIMENSION(2) :: cc,ds,dp
1107	ret = .false.
1108	CALL clear_dset(set)
1109	IF (.NOT.ascii_header(path,ds,dp)) RETURN
1110	ALLOCATE(set%x(ds(1)),set%y(ds(2)), &
1111	set%data(ds(1),ds(2)))
1112	cc(:) = 1
1113	DO i=1,PRODUCT(ds)
1114	READ(666,*,iostat=e) vl ! Reads the line
1115	IF (e /= 0) THEN
1116	IF (debug) WRITE(*,'(a)') 'ERROR: Malformed file, line: ',i+2
1117	ret = .false. ; call clear_dset(set) ; CLOSE(666) ; RETURN
1118	ENDIF
1119	! --- X coordinate
1120	set%x(cc(1)) = vl(1)
1121	! --- Y coordinate
1122	set%y(cc(2)) = vl(2)
1123	! --- Data
1124	set%data(cc(1),cc(2)) = vl(3)
1125	! - Update counters
1126	IF (mod(i,dp(1))==0) cc(1) = cc(1)+1 ; IF (cc(1) > ds(1)) cc(1)=1
1127	IF (mod(i,dp(2))==0) cc(2) = cc(2)+1 ; IF (cc(2) > ds(2)) cc(2)=1
1128	ENDDO
1129	READ(666,*,iostat=e) vl(1)
1130	IF (e == 0) THEN
1131	IF (debug) WRITE(*,'(a)') 'ERROR: Extra value found'
1132	ret = .false. ; call clear_dset(set) ; CLOSE(666) ; RETURN
1133	ENDIF
1134	ret = .true.
1135	CLOSE(666)
1136	RETURN
1137	END FUNCTION ascii_rd_2d
1138
1139	FUNCTION ascii_rd_3d(path,set) RESULT(ret)
1140	!! Read a 3D data set from a ASCII file
1141	CHARACTER(len=*), INTENT(in) :: path !! Path of the ASCII file to read
1142	TYPE(dset3d), INTENT(out) :: set !! output 3D dataset
1143	LOGICAL :: ret !! .true. if no error occured, .false. otherwise
1144	INTEGER :: i,e
1145	REAL(kind=8), DIMENSION(4) :: vl
1146	INTEGER, DIMENSION(3) :: cc,ds,dp
1147	ret = .false.
1148	CALL clear_dset(set)
1149	IF (.NOT.ascii_header(path,ds,dp)) RETURN
1150	ALLOCATE(set%x(ds(1)),set%y(ds(2)),set%z(ds(3)), &
1151	set%data(ds(1),ds(2),ds(3)))
1152	cc(:) = 1
1153	DO i=1,PRODUCT(ds)
1154	READ(666,*,iostat=e) vl ! Reads the line
1155	IF (e /= 0) THEN
1156	IF (debug) WRITE(*,'(a)') 'ERROR: Malformed file, line: ',i+2
1157	ret = .false. ; call clear_dset(set) ; CLOSE(666) ; RETURN
1158	ENDIF
1159	! --- X coordinate
1160	set%x(cc(1)) = vl(1)
1161	! --- Y coordinate
1162	set%y(cc(2)) = vl(2)
1163	! --- Z coordinate
1164	set%z(cc(3)) = vl(3)
1165	! --- Data
1166	set%data(cc(1),cc(2),cc(3)) = vl(4)
1167	! - Update counters
1168	IF (mod(i,dp(1))==0) cc(1) = cc(1)+1 ; IF (cc(1) > ds(1)) cc(1)=1
1169	IF (mod(i,dp(2))==0) cc(2) = cc(2)+1 ; IF (cc(2) > ds(2)) cc(2)=1
1170	IF (mod(i,dp(3))==0) cc(3) = cc(3)+1 ; IF (cc(3) > ds(3)) cc(3)=1
1171	ENDDO
1172	READ(666,*,iostat=e) vl(1)
1173	IF (e == 0) THEN
1174	IF (debug) WRITE(*,'(a)') 'ERROR: Extra value found'
1175	ret = .false. ; call clear_dset(set) ; CLOSE(666) ; RETURN
1176	ENDIF
1177	ret = .true.
1178	CLOSE(666)
1179	RETURN
1180	END FUNCTION ascii_rd_3d
1181
1182	FUNCTION ascii_rd_4d(path,set) RESULT(ret)
1183	!! Read a 4D data set from a ASCII file
1184	CHARACTER(len=*), INTENT(in) :: path !! Path of the ASCII file to read
1185	TYPE(dset4d), INTENT(out) :: set !! output 4D dataset
1186	LOGICAL :: ret !! .true. if no error occured, .false. otherwise
1187	INTEGER :: i,e
1188	REAL(kind=8), DIMENSION(5) :: vl
1189	INTEGER, DIMENSION(4) :: cc,ds,dp
1190	ret = .false.
1191	CALL clear_dset(set)
1192	IF (.NOT.ascii_header(path,ds,dp)) RETURN
1193	ALLOCATE(set%x(ds(1)),set%y(ds(2)),set%z(ds(3)),set%t(ds(4)), &
1194	set%data(ds(1),ds(2),ds(3),ds(4)))
1195	cc(:) = 1
1196	DO i=1,PRODUCT(ds)
1197	READ(666,*,iostat=e) vl ! Reads the line
1198	IF (e /= 0) THEN
1199	IF (debug) WRITE(*,'(a)') 'ERROR: Malformed file, line: ',i+2
1200	ret = .false. ; call clear_dset(set) ; CLOSE(666) ; RETURN
1201	ENDIF
1202	! --- X coordinate
1203	set%x(cc(1)) = vl(1)
1204	! --- Y coordinate
1205	set%y(cc(2)) = vl(2)
1206	! --- Z coordinate
1207	set%z(cc(3)) = vl(3)
1208	! --- T coordinate
1209	set%t(cc(4)) = vl(4)
1210	! --- Data
1211	set%data(cc(1),cc(2),cc(3),cc(4)) = vl(5)
1212	! - Update counters
1213	IF (mod(i,dp(1))==0) cc(1) = cc(1)+1 ; IF (cc(1) > ds(1)) cc(1)=1
1214	IF (mod(i,dp(2))==0) cc(2) = cc(2)+1 ; IF (cc(2) > ds(2)) cc(2)=1
1215	IF (mod(i,dp(3))==0) cc(3) = cc(3)+1 ; IF (cc(3) > ds(3)) cc(3)=1
1216	IF (mod(i,dp(4))==0) cc(4) = cc(4)+1 ; IF (cc(4) > ds(4)) cc(4)=1
1217	ENDDO
1218	READ(666,*,iostat=e) vl(1)
1219	IF (e == 0) THEN
1220	IF (debug) WRITE(*,'(a)') 'ERROR: Extra value found'
1221	ret = .false. ; call clear_dset(set) ; CLOSE(666) ; RETURN
1222	ENDIF
1223	CLOSE(666)
1224	ret = .true.
1225	RETURN
1226	END FUNCTION ascii_rd_4d
1227
1228	FUNCTION ascii_rd_5d(path,set) RESULT(ret)
1229	!! Read a 5D data set from a ASCII file
1230	CHARACTER(len=*), INTENT(in) :: path !! Path of the ASCII file to read
1231	TYPE(dset5d), INTENT(out) :: set !! output 5D dataset
1232	LOGICAL :: ret !! .true. if no error occured, .false. otherwise
1233	INTEGER :: i,e
1234	REAL(kind=8), DIMENSION(6) :: vl
1235	INTEGER, DIMENSION(5) :: cc,ds,dp
1236	ret = .false.
1237	CALL clear_dset(set)
1238	IF (.NOT.ascii_header(path,ds,dp)) RETURN
1239	ALLOCATE(set%x(ds(1)),set%y(ds(2)),set%z(ds(3)),set%t(ds(4)),set%w(ds(5)), &
1240	set%data(ds(1),ds(2),ds(3),ds(4),ds(5)))
1241	cc(:) = 1
1242	DO i=1,PRODUCT(ds)
1243	READ(666,*,iostat=e) vl ! Reads the line
1244	IF (e /= 0) THEN
1245	IF (debug) WRITE(*,'(a)') 'ERROR: Malformed file, line: ',i+2
1246	ret = .false. ; call clear_dset(set) ; CLOSE(666) ; RETURN
1247	ENDIF
1248	! --- X coordinate
1249	set%x(cc(1)) = vl(1)
1250	! --- Y coordinate
1251	set%y(cc(2)) = vl(2)
1252	! --- Z coordinate
1253	set%z(cc(3)) = vl(3)
1254	! --- T coordinate
1255	set%t(cc(4)) = vl(4)
1256	! --- W coordinate
1257	set%w(cc(5)) = vl(5)
1258	! --- Data
1259	set%data(cc(1),cc(2),cc(3),cc(4),cc(5)) = vl(6)
1260	! - Update counters
1261	IF (mod(i,dp(1)) == 0) cc(1)=cc(1)+1 ; IF (cc(1) > ds(1)) cc(1)=1
1262	IF (mod(i,dp(2)) == 0) cc(2)=cc(2)+1 ; IF (cc(2) > ds(2)) cc(2)=1
1263	IF (mod(i,dp(3)) == 0) cc(3)=cc(3)+1 ; IF (cc(3) > ds(3)) cc(3)=1
1264	IF (mod(i,dp(4)) == 0) cc(4)=cc(4)+1 ; IF (cc(4) > ds(4)) cc(4)=1
1265	IF (mod(i,dp(5)) == 0) cc(5)=cc(5)+1 ; IF (cc(5) > ds(5)) cc(5)=1
1266	ENDDO
1267	READ(666,*,iostat=e) vl(1)
1268	IF (e == 0) THEN
1269	IF (debug) WRITE(*,'(a)') 'ERROR: Extra value found'
1270	ret = .false. ; call clear_dset(set) ; CLOSE(666) ; RETURN
1271	ENDIF
1272	CLOSE(666)
1273	ret = .true.
1274	RETURN
1275	END FUNCTION ascii_rd_5d
1276
1277	SUBROUTINE clr_1d_set(set)
1278	!! Clear the given 1D data set
1279	TYPE(DSET1D), INTENT(inout) :: set !! dataset object to clear
1280	IF (ALLOCATED(set%data)) DEALLOCATE(set%data)
1281	IF (ALLOCATED(set%x)) DEALLOCATE(set%x)
1282	RETURN
1283	END SUBROUTINE clr_1d_set
1284
1285	SUBROUTINE clr_2d_set(set)
1286	!! Clear the given 2D data set
1287	TYPE(DSET2D), INTENT(inout) :: set !! dataset object to clear
1288	IF (ALLOCATED(set%data)) DEALLOCATE(set%data)
1289	IF (ALLOCATED(set%x)) DEALLOCATE(set%x)
1290	IF (ALLOCATED(set%y)) DEALLOCATE(set%y)
1291	RETURN
1292	END SUBROUTINE clr_2d_set
1293
1294	SUBROUTINE clr_3d_set(set)
1295	!! Clear the given 3D data set
1296	TYPE(DSET3D), INTENT(inout) :: set !! dataset object to clear
1297	IF (ALLOCATED(set%data)) DEALLOCATE(set%data)
1298	IF (ALLOCATED(set%x)) DEALLOCATE(set%x)
1299	IF (ALLOCATED(set%y)) DEALLOCATE(set%y)
1300	IF (ALLOCATED(set%z)) DEALLOCATE(set%z)
1301	RETURN
1302	END SUBROUTINE clr_3d_set
1303
1304	SUBROUTINE clr_4d_set(set)
1305	!! Clear the given 4D data set
1306	TYPE(DSET4D), INTENT(inout) :: set !! dataset object to clear
1307	IF (ALLOCATED(set%data)) DEALLOCATE(set%data)
1308	IF (ALLOCATED(set%x)) DEALLOCATE(set%x)
1309	IF (ALLOCATED(set%y)) DEALLOCATE(set%y)
1310	IF (ALLOCATED(set%z)) DEALLOCATE(set%z)
1311	IF (ALLOCATED(set%t)) DEALLOCATE(set%t)
1312	RETURN
1313	END SUBROUTINE clr_4d_set
1314
1315	SUBROUTINE clr_5d_set(set)
1316	!! Clear the given 5D data set
1317	TYPE(DSET5D), INTENT(inout) :: set !! dataset object to clear
1318	IF (ALLOCATED(set%data)) DEALLOCATE(set%data)
1319	IF (ALLOCATED(set%x)) DEALLOCATE(set%x)
1320	IF (ALLOCATED(set%y)) DEALLOCATE(set%y)
1321	IF (ALLOCATED(set%z)) DEALLOCATE(set%z)
1322	IF (ALLOCATED(set%t)) DEALLOCATE(set%t)
1323	IF (ALLOCATED(set%w)) DEALLOCATE(set%w)
1324	RETURN
1325	END SUBROUTINE clr_5d_set
1326
1327	FUNCTION has_d_1d(dset) RESULT(yes)
1328	!! Check whether or not the dataset has data.
1329	TYPE(DSET1D), INTENT(in) :: dset !! Dataset to check
1330	LOGICAL :: yes !! return status
1331	yes = (ALLOCATED(dset%data).AND. &
1332	ALLOCATED(dset%x))
1333	END FUNCTION has_d_1d
1334
1335	FUNCTION has_d_2d(dset) RESULT(yes)
1336	!! Check whether or not the dataset has data.
1337	TYPE(DSET2D), INTENT(in) :: dset !! Dataset to check
1338	LOGICAL :: yes !! return status
1339	yes = (ALLOCATED(dset%data).AND. &
1340	ALLOCATED(dset%x) .AND. &
1341	ALLOCATED(dset%y))
1342	END FUNCTION has_d_2d
1343
1344	FUNCTION has_d_3d(dset) RESULT(yes)
1345	!! Check whether or not the dataset has data.
1346	TYPE(DSET3D), INTENT(in) :: dset !! Dataset to check
1347	LOGICAL :: yes !! return status
1348	yes = (ALLOCATED(dset%data).AND. &
1349	ALLOCATED(dset%x) .AND. &
1350	ALLOCATED(dset%y) .AND. &
1351	ALLOCATED(dset%z))
1352	END FUNCTION has_d_3d
1353
1354	FUNCTION has_d_4d(dset) RESULT(yes)
1355	!! Check whether or not the dataset has data.
1356	TYPE(DSET4D), INTENT(in) :: dset !! Dataset to check
1357	LOGICAL :: yes !! return status
1358	yes = (ALLOCATED(dset%data).AND. &
1359	ALLOCATED(dset%x) .AND. &
1360	ALLOCATED(dset%y) .AND. &
1361	ALLOCATED(dset%z) .AND. &
1362	ALLOCATED(dset%t))
1363	END FUNCTION has_d_4d
1364
1365	FUNCTION has_d_5d(dset) RESULt(yes)
1366	!! Check whether or not the dataset has data.
1367	TYPE(DSET5D), INTENT(in) :: dset !! Dataset to check
1368	LOGICAL :: yes !! return status
1369	yes = (ALLOCATED(dset%data).AND. &
1370	ALLOCATED(dset%x) .AND. &
1371	ALLOCATED(dset%y) .AND. &
1372	ALLOCATED(dset%z) .AND. &
1373	ALLOCATED(dset%t) .AND. &
1374	ALLOCATED(dset%w))
1375	END FUNCTION has_d_5d
1376
1377	FUNCTION nc_wr_1d(path,dset) RESULT(ret)
1378	CHARACTER(len=*), INTENT(in) :: path !! Path of the netcdf file.
1379	TYPE(DSET1D), INTENT(in) :: dset !! Dataset to write
1380	LOGICAL :: ret !! Return status
1381	INTEGER :: fi,vi,nd,ds,iret
1382	INTEGER, DIMENSION(:), ALLOCATABLE :: di
1383	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: tv
1384	CHARACTER(len=15) :: i2s
1385	LOGICAL :: hxd,ok
1386	INTEGER :: xid,vid
1387	ret = has_data(dset)
1388	IF (.NOT.ret) RETURN
1389	ret = .false.
1390	INQUIRE(FILE=TRIM(path),EXIST=ok)
1391	IF (ok) THEN
1392	IF (NF90_OPEN(TRIM(path), NF90_WRITE, fi) /= NF90_NOERR) RETURN
1393	ELSE
1394	IF (NF90_CREATE(TRIM(path), NF90_NOCLOBBER, fi) /= NF90_NOERR) RETURN
1395	ENDIF
1396	iret = NF90_CLOSE(fi)
1397	! if variable already exist: get out of here !
1398	IF (get_nc_info(path,dset%dname,fi,vi,di)) RETURN
1399
1400	iret = NF90_OPEN(path,NF90_WRITE,fi)
1401
1402	hxd = nc_get_dim_by_name(fi,dset%xname,tv,ds,xid)
1403	IF (hxd) THEN
1404	IF(.NOT.compare(dset%x,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1405	ENDIF
1406
1407	IF (.NOT.hxd) THEN
1408	ret = nc_set_dim(fi,dset%xname,dset%x,xid)
1409	IF (.NOT.ret) THEN ; iret= NF90_CLOSE(fi) ; RETURN ; ENDIF
1410	ENDIF
1411	ret = .false.
1412	iret = nf90_redef(fi)
1413	iret = nf90_def_var(fi, TRIM(dset%dname), nf90_double,(/xid/),vid)
1414	IF (iret /= 0) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1415	iret = nf90_enddef(fi)
1416	iret = nf90_put_var(fi,vid,dset%data)
1417	ret = (iret == 0)
1418	iret=NF90_CLOSE(fi)
1419	ret = .true.
1420	RETURN
1421	END FUNCTION nc_wr_1d
1422
1423	FUNCTION nc_wr_2d(path,dset) RESULT(ret)
1424	CHARACTER(len=*), INTENT(in) :: path !! Path of the netcdf file.
1425	TYPE(DSET2D), INTENT(in) :: dset !! Dataset to write
1426	LOGICAL :: ret !! Return status
1427	INTEGER :: fi,vi,nd,ds,iret,nferr
1428	INTEGER, DIMENSION(:), ALLOCATABLE :: di
1429	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: tv
1430	CHARACTER(len=15) :: i2s
1431	LOGICAL :: hxd,hyd,ok
1432	INTEGER :: xid,yid,vid
1433	ret = has_data(dset)
1434	IF (.NOT.ret) RETURN
1435	ret = .false.
1436	INQUIRE(FILE=TRIM(path),EXIST=ok)
1437	IF (ok) THEN
1438	IF (NF90_OPEN(TRIM(path), NF90_WRITE, fi) /= NF90_NOERR) RETURN
1439	ELSE
1440	IF (NF90_CREATE(TRIM(path), NF90_NOCLOBBER, fi) /= NF90_NOERR) RETURN
1441	ENDIF
1442	iret = NF90_CLOSE(fi)
1443
1444	! if variable already exist: get out of here !
1445	IF (get_nc_info(path,dset%dname,fi,vi,di)) RETURN
1446
1447	iret = NF90_OPEN(path,NF90_WRITE,fi)
1448
1449	hxd = nc_get_dim_by_name(fi,dset%xname,tv,ds,xid)
1450	IF (hxd) THEN
1451	IF(.NOT.compare(dset%x,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1452	ENDIF
1453	hyd = nc_get_dim_by_name(fi,dset%yname,tv,ds,yid)
1454	IF (hyd) THEN
1455	IF(.NOT.compare(dset%y,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1456	ENDIF
1457
1458	IF (.NOT.hxd) THEN
1459	ret = nc_set_dim(fi,dset%xname,dset%x,xid)
1460	IF (.NOT.ret) THEN ; iret= NF90_CLOSE(fi) ; RETURN ; ENDIF
1461	ENDIF
1462	IF (.NOT.hyd) THEN
1463	ret = nc_set_dim(fi,dset%yname,dset%y,yid)
1464	IF (.NOT.ret) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1465	ENDIF
1466	ret = .false.
1467	iret = nf90_redef(fi)
1468	iret = nf90_def_var(fi, TRIM(dset%dname), nf90_double,(/xid,yid/),vid)
1469	IF (iret /= 0) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1470	iret = nf90_enddef(fi)
1471	iret = nf90_put_var(fi,vid,dset%data)
1472	ret = (iret == 0)
1473	iret=NF90_CLOSE(fi)
1474	ret = .true.
1475	RETURN
1476	END FUNCTION nc_wr_2d
1477
1478	FUNCTION nc_wr_3d(path,dset) RESULT(ret)
1479	CHARACTER(len=*), INTENT(in) :: path !! Path of the netcdf file.
1480	TYPE(DSET3D), INTENT(in) :: dset !! Dataset to write
1481	LOGICAL :: ret !! Return status
1482	INTEGER :: fi,vi,nd,ds,iret
1483	INTEGER, DIMENSION(:), ALLOCATABLE :: di
1484	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: tv
1485	CHARACTER(len=15) :: i2s
1486	LOGICAL :: hxd,hyd,hzd,ok
1487	INTEGER :: xid,yid,zid,vid
1488	ret = has_data(dset)
1489	IF (.NOT.ret) RETURN
1490	ret = .false.
1491	INQUIRE(FILE=TRIM(path),EXIST=ok)
1492	IF (ok) THEN
1493	IF (NF90_OPEN(TRIM(path), NF90_WRITE, fi) /= NF90_NOERR) RETURN
1494	ELSE
1495	IF (NF90_CREATE(TRIM(path), NF90_NOCLOBBER, fi) /= NF90_NOERR) RETURN
1496	ENDIF
1497	iret = NF90_CLOSE(fi)
1498
1499	! if variable already exist: get out of here !
1500	IF (get_nc_info(path,dset%dname,fi,vi,di)) RETURN
1501
1502	iret = NF90_OPEN(path,NF90_WRITE,fi)
1503
1504	hxd = nc_get_dim_by_name(fi,dset%xname,tv,ds,xid)
1505	IF (hxd) THEN
1506	IF(.NOT.compare(dset%x,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1507	ENDIF
1508	hyd = nc_get_dim_by_name(fi,dset%yname,tv,ds,yid)
1509	IF (hyd) THEN
1510	IF(.NOT.compare(dset%y,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1511	ENDIF
1512	hzd = nc_get_dim_by_name(fi,dset%zname,tv,ds,zid)
1513	IF (hzd) THEN
1514	IF(.NOT.compare(dset%z,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1515	ENDIF
1516
1517	IF (.NOT.hxd) THEN
1518	ret = nc_set_dim(fi,dset%xname,dset%x,xid)
1519	IF (.NOT.ret) THEN ; iret= NF90_CLOSE(fi) ; RETURN ; ENDIF
1520	ENDIF
1521	IF (.NOT.hyd) THEN
1522	ret = nc_set_dim(fi,dset%yname,dset%y,yid)
1523	IF (.NOT.ret) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1524	ENDIF
1525	IF (.NOT.hzd) THEN
1526	ret = nc_set_dim(fi,dset%zname,dset%z,zid)
1527	IF (.NOT.ret) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1528	ENDIF
1529	ret = .false.
1530	iret = nf90_redef(fi)
1531	iret = nf90_def_var(fi, TRIM(dset%dname), nf90_double,(/xid,yid,zid/),vid)
1532	IF (iret /= 0) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1533	iret = nf90_enddef(fi)
1534	iret = nf90_put_var(fi,vid,dset%data)
1535	ret = (iret == 0)
1536	iret=NF90_CLOSE(fi)
1537	ret = .true.
1538	RETURN
1539	END FUNCTION nc_wr_3d
1540
1541	FUNCTION nc_wr_4d(path,dset) RESULT(ret)
1542	CHARACTER(len=*), INTENT(in) :: path !! Path of the netcdf file.
1543	TYPE(DSET4D), INTENT(in) :: dset !! Dataset to write
1544	LOGICAL :: ret !! Return status
1545	INTEGER :: fi,vi,nd,ds,iret
1546	INTEGER, DIMENSION(:), ALLOCATABLE :: di
1547	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: tv
1548	CHARACTER(len=15) :: i2s
1549	LOGICAL :: hxd,hyd,hzd,htd,ok
1550	INTEGER :: xid,yid,zid,tid,vid
1551	ret = has_data(dset)
1552	IF (.NOT.ret) RETURN
1553	ret = .false.
1554	INQUIRE(FILE=TRIM(path),EXIST=ok)
1555	IF (ok) THEN
1556	IF (NF90_OPEN(TRIM(path), NF90_WRITE, fi) /= NF90_NOERR) RETURN
1557	ELSE
1558	IF (NF90_CREATE(TRIM(path), NF90_NOCLOBBER, fi) /= NF90_NOERR) RETURN
1559	ENDIF
1560	iret = NF90_CLOSE(fi)
1561
1562	! if variable already exist: get out of here !
1563	IF (get_nc_info(path,dset%dname,fi,vi,di)) RETURN
1564
1565	iret = NF90_OPEN(path,NF90_WRITE,fi)
1566
1567	hxd = nc_get_dim_by_name(fi,dset%xname,tv,ds,xid)
1568	IF (hxd) THEN
1569	IF(.NOT.compare(dset%x,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1570	ENDIF
1571	hyd = nc_get_dim_by_name(fi,dset%yname,tv,ds,yid)
1572	IF (hyd) THEN
1573	IF(.NOT.compare(dset%y,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1574	ENDIF
1575	hzd = nc_get_dim_by_name(fi,dset%zname,tv,ds,zid)
1576	IF (hzd) THEN
1577	IF(.NOT.compare(dset%z,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1578	ENDIF
1579	htd = nc_get_dim_by_name(fi,dset%tname,tv,ds,tid)
1580	IF (htd) THEN
1581	IF(.NOT.compare(dset%t,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1582	ENDIF
1583
1584	IF (.NOT.hxd) THEN
1585	ret = nc_set_dim(fi,dset%xname,dset%x,xid)
1586	IF (.NOT.ret) THEN ; iret= NF90_CLOSE(fi) ; RETURN ; ENDIF
1587	ENDIF
1588	IF (.NOT.hyd) THEN
1589	ret = nc_set_dim(fi,dset%yname,dset%y,yid)
1590	IF (.NOT.ret) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1591	ENDIF
1592	IF (.NOT.hzd) THEN
1593	ret = nc_set_dim(fi,dset%zname,dset%z,zid)
1594	IF (.NOT.ret) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1595	ENDIF
1596	IF (.NOT.htd) THEN
1597	ret = nc_set_dim(fi,dset%tname,dset%t,tid)
1598	IF (.NOT.ret) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1599	ENDIF
1600	ret = .false.
1601	iret = nf90_redef(fi)
1602	iret = nf90_def_var(fi, TRIM(dset%dname), nf90_double,(/xid,yid,zid,tid/),vid)
1603	IF (iret /= 0) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1604	iret = nf90_enddef(fi)
1605	iret = nf90_put_var(fi,vid,dset%data)
1606	ret = (iret == 0)
1607	iret=NF90_CLOSE(fi)
1608	ret = .true.
1609	RETURN
1610	END FUNCTION nc_wr_4d
1611
1612	FUNCTION nc_wr_5d(path,dset) RESULT(ret)
1613	CHARACTER(len=*), INTENT(in) :: path !! Path of the netcdf file.
1614	TYPE(DSET5D), INTENT(in) :: dset !! Dataset to write
1615	LOGICAL :: ret !! Return status
1616	INTEGER :: fi,vi,nd,ds,iret
1617	INTEGER, DIMENSION(:), ALLOCATABLE :: di
1618	REAL(kind=8), DIMENSION(:), ALLOCATABLE :: tv
1619	CHARACTER(len=15) :: i2s
1620	LOGICAL :: hxd,hyd,hzd,htd,hwd,ok
1621	INTEGER :: xid,yid,zid,tid,wid,vid
1622	ret = has_data(dset)
1623	IF (.NOT.ret) RETURN
1624	ret = .false.
1625	INQUIRE(FILE=TRIM(path),EXIST=ok)
1626	IF (ok) THEN
1627	IF (NF90_OPEN(TRIM(path), NF90_WRITE, fi) /= NF90_NOERR) RETURN
1628	ELSE
1629	IF (NF90_CREATE(TRIM(path), NF90_NOCLOBBER, fi) /= NF90_NOERR) RETURN
1630	ENDIF
1631	iret = NF90_CLOSE(fi)
1632
1633	! if variable already exist: get out of here !
1634	IF (get_nc_info(path,dset%dname,fi,vi,di)) RETURN
1635
1636	iret = NF90_OPEN(path,NF90_WRITE,fi)
1637
1638	hxd = nc_get_dim_by_name(fi,dset%xname,tv,ds,xid)
1639	IF (hxd) THEN
1640	IF(.NOT.compare(dset%x,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1641	ENDIF
1642	hyd = nc_get_dim_by_name(fi,dset%yname,tv,ds,yid)
1643	IF (hyd) THEN
1644	IF(.NOT.compare(dset%y,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1645	ENDIF
1646	hzd = nc_get_dim_by_name(fi,dset%zname,tv,ds,zid)
1647	IF (hzd) THEN
1648	IF(.NOT.compare(dset%z,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1649	ENDIF
1650	htd = nc_get_dim_by_name(fi,dset%tname,tv,ds,tid)
1651	IF (htd) THEN
1652	IF(.NOT.compare(dset%t,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1653	ENDIF
1654	hwd = nc_get_dim_by_name(fi,dset%wname,tv,ds,wid)
1655	IF (hwd) THEN
1656	IF (.NOT.compare(dset%w,tv)) THEN ; ret = .false. ; RETURN ; ENDIF
1657	ENDIF
1658
1659	IF (.NOT.hxd) THEN
1660	ret = nc_set_dim(fi,dset%xname,dset%x,xid)
1661	IF (.NOT.ret) THEN ; iret= NF90_CLOSE(fi) ; RETURN ; ENDIF
1662	ENDIF
1663	IF (.NOT.hyd) THEN
1664	ret = nc_set_dim(fi,dset%yname,dset%y,yid)
1665	IF (.NOT.ret) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1666	ENDIF
1667	IF (.NOT.hzd) THEN
1668	ret = nc_set_dim(fi,dset%zname,dset%z,zid)
1669	IF (.NOT.ret) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1670	ENDIF
1671	IF (.NOT.htd) THEN
1672	ret = nc_set_dim(fi,dset%tname,dset%t,tid)
1673	IF (.NOT.ret) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1674	ENDIF
1675	IF (.NOT.hwd) THEN
1676	ret = nc_set_dim(fi,dset%wname,dset%w,wid)
1677	IF (.NOT.ret) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1678	ENDIF
1679	ret = .false.
1680	iret = nf90_redef(fi)
1681	iret = nf90_def_var(fi, TRIM(dset%dname), nf90_double,(/xid,yid,zid,tid,wid/),vid)
1682	IF (iret /= 0) THEN ; iret = NF90_CLOSE(fi) ; RETURN ; ENDIF
1683	iret = nf90_enddef(fi)
1684	iret = nf90_put_var(fi,vid,dset%data)
1685	ret = (iret == 0)
1686	iret=NF90_CLOSE(fi)
1687	ret = .true.
1688	RETURN
1689	END FUNCTION nc_wr_5d
1690
1691	FUNCTION nc_set_dim(fid,name,values,dimid) RESULT(ret)
1692	!! Set a new dimension (and its associated values) in a NetCDF file.
1693	INTEGER, INTENT(in) :: fid !! Netcdf file id.
1694	CHARACTER(len=*), INTENT(in) :: name !! Name of the dimension to set.
1695	REAL(kind=8), DIMENSION(:), INTENT(in) :: values !! Associated values.
1696	INTEGER, INTENT(out) :: dimid !! Output dimension id.
1697	LOGICAL :: ret !! Return status
1698	INTEGER :: iret,vid
1699	ret = .false.
1700	iret = nf90_redef(fid)
1701	iret = nf90_def_dim(fid, TRIM(name), size(values), dimid)
1702	WRITE(*,'(a)') "nc_set_dim: "//TRIM(name)//" --> "//TRIM(nf90_strerror(iret))
1703	IF (iret /=0) RETURN
1704	iret = nf90_def_var(fid, TRIM(name), nf90_double,(/dimid/),vid)
1705	IF (iret /=0) RETURN
1706	iret = nf90_enddef(fid)
1707	iret = nf90_put_var(fid, vid,values)
1708	ret = .true.
1709	END FUNCTION nc_set_dim
1710
1711	FUNCTION compare(vec1,vec2,tol) RESULT(ret)
1712	!! Compare two vector of double.
1713	!!
1714	!! The test is checked against the difference (element wise) of the two vector compared to
1715	!! to a given tolerance.
1716	REAL(kind=8), DIMENSION(:), INTENT(in) :: vec1 !! First vector to compare
1717	REAL(kind=8), DIMENSION(:), INTENT(in) :: vec2 !! Second vector to compare
1718	REAL(kind=8), INTENT(in), OPTIONAL :: tol !! Tolerance to apply for the comparison
1719	LOGICAL :: ret !! .true. if both vectors are equal, .false. otherwise.
1720	REAL(kind=8) :: ztol
1721	INTEGER :: i
1722	ztol = 1d-10 ; IF (PRESENT(tol)) ztol = abs(tol)
1723	ret = .false.
1724	IF (size(vec1) /= size(vec2)) RETURN
1725
1726	DO i=1,size(vec1) ; IF (abs(vec1(i)-vec2(i)) > ztol) RETURN ; ENDDO
1727	ret = .true.
1728	END FUNCTION compare
1729
1730	END MODULE LINT_DATASETS

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

Download in other formats:

Original Format