DISPATCH
h5_mod.f90
1 !===============================================================================
2 !> Module interface to the HDF5 library
3 !===============================================================================
4 
5 MODULE h5_mod
6  USE hdf5
7  USE io_unit_mod
8  USE omp_lock_mod
9  implicit none
10  private
11  type, public:: h5_t
12  integer:: err
13  integer:: verbose=0
14  integer(HID_T):: fid=0, gid=0
15  type(lock_t):: lock
16  contains
17  procedure:: open
18  procedure:: close
19  procedure:: group_open
20  procedure:: set_open
21  procedure:: set_write
22  procedure:: set_close
23  procedure:: att_open
24  procedure:: att_write
25  procedure:: att_close
26  procedure:: exists
27  procedure, private:: ints
28  procedure, private:: ints_1d
29  procedure, private:: ints_2d
30  procedure, private:: ints_3d
31  procedure, private:: ints_4d
32  procedure, private:: real
33  procedure, private:: real_1d
34  procedure, private:: real_2d
35  procedure, private:: real_3d
36  procedure, private:: real_4d
37  generic, public:: output => ints, ints_1d, ints_2d, ints_3d, ints_4d, &
38  real, real_1d, real_2d, real_3d, real_4d
39  end type
40  type(h5_t), public:: h5
41 CONTAINS
42 
43 !> -----------------------------------------------------------------------------
44 !> Interface FUNCTION to open or create an HDF5 file
45 !> -----------------------------------------------------------------------------
46 INTEGER(HID_T) FUNCTION open (self, name, new, verbose) result (fid)
47  class(h5_t):: self
48  character(len=*):: name
49  logical, optional:: new
50  integer, optional:: verbose
51  logical:: exists
52  !.............................................................................
53  call h5open_f (self%err)
54  inquire (file=name, exist=exists)
55  if (exists) then
56  call h5fopen_f (name, h5f_acc_rdwr_f, fid, self%err)
57  else
58  call h5fcreate_f (name, h5f_acc_excl_f, fid, self%err)
59  end if
60  self%fid = fid
61  if (present(new)) then
62  new = .not.exists
63  end if
64  if (present(verbose)) &
65  self%verbose = verbose
66  if (self%verbose > 1) &
67  print *, 'open ', trim(name), exists, self%err
68 END FUNCTION open
69 
70 !> -----------------------------------------------------------------------------
71 !> Close an HDF5 file
72 !> -----------------------------------------------------------------------------
73 SUBROUTINE close (self, id)
74  class(h5_t):: self
75  integer(HID_T):: id
76  !.............................................................................
77  call h5fclose_f (id, self%err)
78  call h5close_f (self%err)
79  if (self%verbose > 1) &
80  print *, 'close', self%err
81 END SUBROUTINE close
82 
83 !> -----------------------------------------------------------------------------
84 !> Open group
85 !> -----------------------------------------------------------------------------
86 INTEGER(HID_T) FUNCTION group_open (self, id, name, new) result (gid)
87  class(h5_t):: self
88  integer(HID_T):: id
89  character(len=*):: name
90  logical, optional:: new
91  integer:: type, links, corder
92  logical:: exists
93  !.............................................................................
94  call h5lexists_f (id, name, exists, self%err)
95  if (exists) then
96  call h5gopen_f (id, name, gid, self%err)
97  else
98  call h5gcreate_f (id, name, gid, self%err)
99  end if
100  self%gid = gid
101  if (present(new)) then
102  new = .not.exists
103  end if
104  if (self%verbose > 1) &
105  print *, 'group_open ', name, new, self%err
106 END FUNCTION group_open
107 
108 !> -----------------------------------------------------------------------------
109 !> Open or create dataset
110 !> -----------------------------------------------------------------------------
111 INTEGER(HID_T) FUNCTION set_open (self, id, name, a, new) result (did)
112  class(h5_t):: self
113  integer(HID_T):: id, sid
114  character(len=*):: name
115  real:: a(:,:,:)
116  integer:: type, links, corder
117  logical, optional:: new
118  logical:: exists
119  !.............................................................................
120  call h5lexists_f (id, name, exists, self%err)
121  if (exists) then
122  call h5dopen_f (id, name, did, self%err)
123  else
124  call h5screate_simple_f (SIZE(shape(a)), shape(a,kind=hsize_t), sid, self%err)
125  call h5dcreate_f (id, name, h5t_native_real, sid, did, self%err)
126  end if
127  if (present(new)) then
128  new = .not.exists
129  end if
130  if (self%verbose > 1) &
131  print *, 'set_open ', name, exists, self%err
132 END FUNCTION set_open
133 
134 !> -----------------------------------------------------------------------------
135 !> Write a data set
136 !> -----------------------------------------------------------------------------
137 SUBROUTINE set_write (self, id, a)
138  class(h5_t):: self
139  integer(HID_T):: id
140  real:: a(:,:,:)
141  !.............................................................................
142  call h5dwrite_f (id, h5t_native_real, a, shape(a,kind=hsize_t), self%err)
143  if (self%verbose > 1) &
144  print *, 'set_write', shape(a), self%err
145 END SUBROUTINE set_write
146 
147 !> -----------------------------------------------------------------------------
148 !> Close a data set
149 !> -----------------------------------------------------------------------------
150 SUBROUTINE set_close (self, id)
151  class(h5_t):: self
152  integer(HID_T):: id
153  !.............................................................................
154  call h5dclose_f (id, self%err)
155  if (self%verbose > 1) &
156  print *, 'set_close', self%err
157 END SUBROUTINE set_close
158 
159 !> -----------------------------------------------------------------------------
160 !> Open or create dataset
161 !> -----------------------------------------------------------------------------
162 INTEGER(HID_T) FUNCTION att_open (self, id, name, a, new) result (did)
163  class(h5_t):: self
164  integer(HID_T):: id, sid
165  character(len=*):: name
166  real:: a(:,:,:)
167  integer:: type, links, corder
168  logical, optional:: new
169  logical:: exists
170  !.............................................................................
171  call h5aexists_f (id, name, exists, self%err)
172  if (exists) then
173  call h5aopen_f (id, name, did, self%err)
174  else
175  call h5screate_simple_f (1, shape(a,kind=hsize_t), sid, self%err)
176  call h5acreate_f (id, name, h5t_native_real, sid, did, self%err)
177  end if
178  if (present(new)) then
179  new = .not.exists
180  end if
181  if (self%verbose > 1) &
182  print *, 'att_open ', name, exists, self%err
183 END FUNCTION att_open
184 
185 !> -----------------------------------------------------------------------------
186 !> Write a data set
187 !> -----------------------------------------------------------------------------
188 SUBROUTINE att_write (self, id, a)
189  class(h5_t):: self
190  integer(HID_T):: id
191  real:: a(:,:,:)
192  !.............................................................................
193  call h5awrite_f (id, h5t_native_real, a, shape(a,kind=hsize_t), self%err)
194  if (self%verbose > 1) &
195  print *, 'att_write', shape(a), self%err
196 END SUBROUTINE att_write
197 
198 !> -----------------------------------------------------------------------------
199 !> Close a data set
200 !> -----------------------------------------------------------------------------
201 SUBROUTINE att_close (self, id)
202  class(h5_t):: self
203  integer(HID_T):: id
204  !.............................................................................
205  call h5aclose_f (id, self%err)
206  if (self%verbose > 1) &
207  print *, 'att_close', self%err
208 END SUBROUTINE att_close
209 
210 !> -----------------------------------------------------------------------------
211 !> Check if /sequence/RRRRRR/variable
212 !> -----------------------------------------------------------------------------
213 LOGICAL FUNCTION exists (self, sequence, record, name, id, attribute)
214  class(h5_t) :: self
215  character(len=*) :: sequence, name
216  integer :: record
217  character(len=6) :: srecord
218  integer(HID_T) :: id
219  logical, optional:: attribute
220  !.............................................................................
221  if (self%fid == 0) then
222  self%fid = self%open (trim(io_unit%outputname)//'/hdf5.dat')
223  end if
224  write (srecord,'(i6.6)') record
225  id = self%group_open (self%group_open (self%fid, sequence), srecord)
226  if (present(attribute)) then
227  call h5aexists_f (id, name, exists, self%err)
228  else
229  call h5lexists_f (id, name, exists, self%err)
230  end if
231 END FUNCTION exists
232 
233 !> -----------------------------------------------------------------------------
234 SUBROUTINE ints (self, seq, record, name, a)
235  class(h5_t) :: self
236  character(len=*) :: seq, name
237  integer :: record
238  integer :: a
239  integer(HID_T) :: id, sid, did, rank(1)
240  !.............................................................................
241  if (self%exists(seq, record, name, id, attribute=.true.)) then
242  call h5aopen_f (id, name, did, self%err)
243  else
244  rank = 1
245  call h5screate_simple_f (size(rank), rank, sid, self%err)
246  call h5acreate_f (id, name, h5t_native_integer, sid, did, self%err)
247  end if
248  call h5awrite_f (did, h5t_native_integer, a, rank, self%err)
249  if (self%verbose > 1) &
250  print *, 'h5_t%ints: attribute =', name, a, self%err
251  call h5aclose_f (did, self%err)
252 END SUBROUTINE ints
253 !> -----------------------------------------------------------------------------
254 SUBROUTINE ints_1d (self, seq, record, name, a)
255  class(h5_t) :: self
256  character(len=*) :: seq, name
257  integer :: record
258  integer :: a(:)
259  integer(HID_T) :: id, sid, did
260  !.............................................................................
261  if (self%exists(seq, record, name, id)) then
262  call h5dopen_f (id, name, did, self%err)
263  else
264  call h5screate_simple_f (SIZE(shape(a)), shape(a,kind=hsize_t), sid, self%err)
265  call h5dcreate_f (id, name, h5t_native_integer, sid, did, self%err)
266  end if
267  call h5dwrite_f (did, h5t_native_integer, a, shape(a,kind=hsize_t), self%err)
268  call h5dclose_f (did, self%err)
269 END SUBROUTINE ints_1d
270 !> -----------------------------------------------------------------------------
271 SUBROUTINE ints_2d (self, seq, record, name, a)
272  class(h5_t) :: self
273  character(len=*) :: seq, name
274  integer :: record
275  integer :: a(:,:)
276  integer(HID_T) :: id, sid, did
277  !.............................................................................
278  if (self%exists(seq, record, name, id)) then
279  call h5dopen_f (id, name, did, self%err)
280  else
281  call h5screate_simple_f (SIZE(shape(a)), shape(a,kind=hsize_t), sid, self%err)
282  call h5dcreate_f (id, name, h5t_native_integer, sid, did, self%err)
283  end if
284  call h5dwrite_f (did, h5t_native_integer, a, shape(a,kind=hsize_t), self%err)
285  call h5dclose_f (did, self%err)
286 END SUBROUTINE ints_2d
287 !> -----------------------------------------------------------------------------
288 SUBROUTINE ints_3d (self, seq, record, name, a)
289  class(h5_t) :: self
290  character(len=*) :: seq, name
291  integer :: record
292  integer :: a(:,:,:)
293  integer(HID_T) :: id, sid, did
294  !.............................................................................
295  if (self%exists(seq, record, name, id)) then
296  call h5dopen_f (id, name, did, self%err)
297  else
298  call h5screate_simple_f (SIZE(shape(a)), shape(a,kind=hsize_t), sid, self%err)
299  call h5dcreate_f (id, name, h5t_native_integer, sid, did, self%err)
300  end if
301  call h5dwrite_f (did, h5t_native_integer, a, shape(a,kind=hsize_t), self%err)
302  call h5dclose_f (did, self%err)
303 END SUBROUTINE ints_3d
304 !> -----------------------------------------------------------------------------
305 SUBROUTINE ints_4d (self, seq, record, name, a)
306  class(h5_t) :: self
307  character(len=*) :: seq, name
308  integer :: record
309  integer :: a(:,:,:,:)
310  integer(HID_T) :: id, sid, did
311  !.............................................................................
312  if (self%exists(seq, record, name, id)) then
313  call h5dopen_f (id, name, did, self%err)
314  else
315  call h5screate_simple_f (SIZE(shape(a)), shape(a,kind=hsize_t), sid, self%err)
316  call h5dcreate_f (id, name, h5t_native_integer, sid, did, self%err)
317  end if
318  call h5dwrite_f (did, h5t_native_integer, a, shape(a,kind=hsize_t), self%err)
319  call h5dclose_f (did, self%err)
320 END SUBROUTINE ints_4d
321 
322 !> -----------------------------------------------------------------------------
323 SUBROUTINE real (self, seq, record, name, a)
324  class(h5_t) :: self
325  character(len=*) :: seq, name
326  integer :: record
327  real :: a
328  integer(HID_T) :: id, sid, did, rank(1)
329  !.............................................................................
330  if (self%exists(seq, record, name, id, attribute=.true.)) then
331  call h5aopen_f (id, name, did, self%err)
332  else
333  rank = 1
334  call h5screate_simple_f (size(rank), rank, sid, self%err)
335  call h5acreate_f (id, name, h5t_native_real, sid, did, self%err)
336  end if
337  call h5awrite_f (did, h5t_native_real, a, rank, self%err)
338  if (self%verbose > 1) &
339  print *, 'h5_t%real: attribute =', name, a, self%err
340  call h5aclose_f (did, self%err)
341 END SUBROUTINE real
342 !> -----------------------------------------------------------------------------
343 SUBROUTINE real_1d (self, seq, record, name, a)
344  class(h5_t) :: self
345  character(len=*) :: seq, name
346  integer :: record
347  real :: a(:)
348  integer(HID_T) :: id, sid, did
349  !.............................................................................
350  if (self%exists(seq, record, name, id)) then
351  call h5dopen_f (id, name, did, self%err)
352  else
353  call h5screate_simple_f (SIZE(shape(a)), shape(a,kind=hsize_t), sid, self%err)
354  call h5dcreate_f (id, name, h5t_native_real, sid, did, self%err)
355  end if
356  call h5dwrite_f (did, h5t_native_real, a, shape(a,kind=hsize_t), self%err)
357  call h5dclose_f (did, self%err)
358 END SUBROUTINE real_1d
359 !> -----------------------------------------------------------------------------
360 SUBROUTINE real_2d (self, seq, record, name, a)
361  class(h5_t) :: self
362  character(len=*) :: seq, name
363  integer :: record
364  real :: a(:,:)
365  integer(HID_T) :: id, sid, did
366  !.............................................................................
367  if (self%exists(seq, record, name, id)) then
368  call h5dopen_f (id, name, did, self%err)
369  else
370  call h5screate_simple_f (SIZE(shape(a)), shape(a,kind=hsize_t), sid, self%err)
371  call h5dcreate_f (id, name, h5t_native_real, sid, did, self%err)
372  end if
373  call h5dwrite_f (did, h5t_native_real, a, shape(a,kind=hsize_t), self%err)
374  call h5dclose_f (did, self%err)
375 END SUBROUTINE real_2d
376 !> -----------------------------------------------------------------------------
377 SUBROUTINE real_3d (self, seq, record, name, a)
378  class(h5_t) :: self
379  character(len=*) :: seq, name
380  integer :: record
381  real :: a(:,:,:)
382  integer(HID_T) :: id, sid, did
383  !.............................................................................
384  if (self%exists(seq, record, name, id)) then
385  call h5dopen_f (id, name, did, self%err)
386  else
387  call h5screate_simple_f (SIZE(shape(a)), shape(a,kind=hsize_t), sid, self%err)
388  call h5dcreate_f (id, name, h5t_native_real, sid, did, self%err)
389  end if
390  call h5dwrite_f (did, h5t_native_real, a, shape(a,kind=hsize_t), self%err)
391  call h5dclose_f (did, self%err)
392 END SUBROUTINE real_3d
393 !> -----------------------------------------------------------------------------
394 SUBROUTINE real_4d (self, seq, record, name, a)
395  class(h5_t) :: self
396  character(len=*) :: seq, name
397  integer :: record
398  real :: a(:,:,:,:)
399  integer(HID_T) :: id, sid, did
400  !.............................................................................
401  if (self%exists(seq, record, name, id)) then
402  call h5dopen_f (id, name, did, self%err)
403  else
404  call h5screate_simple_f (SIZE(shape(a)), shape(a,kind=hsize_t), sid, self%err)
405  call h5dcreate_f (id, name, h5t_native_real, sid, did, self%err)
406  end if
407  call h5dwrite_f (did, h5t_native_real, a, shape(a,kind=hsize_t), self%err)
408  call h5dclose_f (did, self%err)
409 END SUBROUTINE real_4d
410 
411 END MODULE h5_mod
Module interface to the HDF5 library.
Definition: h5_mod.f90:5
The lock module uses nested locks, to allow versatile use of locks, where a procedure may want to mak...