[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]

hdf5impex.hxx VIGRA

1 /************************************************************************/
2 /* */
3 /* Copyright 2009 by Michael Hanselmann and Ullrich Koethe */
4 /* */
5 /* This file is part of the VIGRA computer vision library. */
6 /* The VIGRA Website is */
7 /* http://hci.iwr.uni-heidelberg.de/vigra/ */
8 /* Please direct questions, bug reports, and contributions to */
9 /* ullrich.koethe@iwr.uni-heidelberg.de or */
10 /* vigra@informatik.uni-hamburg.de */
11 /* */
12 /* Permission is hereby granted, free of charge, to any person */
13 /* obtaining a copy of this software and associated documentation */
14 /* files (the "Software"), to deal in the Software without */
15 /* restriction, including without limitation the rights to use, */
16 /* copy, modify, merge, publish, distribute, sublicense, and/or */
17 /* sell copies of the Software, and to permit persons to whom the */
18 /* Software is furnished to do so, subject to the following */
19 /* conditions: */
20 /* */
21 /* The above copyright notice and this permission notice shall be */
22 /* included in all copies or substantial portions of the */
23 /* Software. */
24 /* */
25 /* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND */
26 /* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES */
27 /* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND */
28 /* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT */
29 /* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, */
30 /* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING */
31 /* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR */
32 /* OTHER DEALINGS IN THE SOFTWARE. */
33 /* */
34 /************************************************************************/
35 
36 #ifndef VIGRA_HDF5IMPEX_HXX
37 #define VIGRA_HDF5IMPEX_HXX
38 
39 #include <string>
40 #include <algorithm>
41 #include <utility>
42 
43 #define H5Gcreate_vers 2
44 #define H5Gopen_vers 2
45 #define H5Dopen_vers 2
46 #define H5Dcreate_vers 2
47 #define H5Acreate_vers 2
48 #define H5Eset_auto_vers 2
49 #define H5Eget_auto_vers 2
50 
51 #include <hdf5.h>
52 
53 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6)
54 # ifndef H5Gopen
55 # define H5Gopen(a, b, c) H5Gopen(a, b)
56 # endif
57 # ifndef H5Gcreate
58 # define H5Gcreate(a, b, c, d, e) H5Gcreate(a, b, 1)
59 # endif
60 # ifndef H5Dopen
61 # define H5Dopen(a, b, c) H5Dopen(a, b)
62 # endif
63 # ifndef H5Dcreate
64 # define H5Dcreate(a, b, c, d, e, f, g) H5Dcreate(a, b, c, d, f)
65 # endif
66 # ifndef H5Acreate
67 # define H5Acreate(a, b, c, d, e, f) H5Acreate(a, b, c, d, e)
68 # endif
69 # ifndef H5Pset_obj_track_times
70 # define H5Pset_obj_track_times(a, b) do {} while (0)
71 # endif
72 # include <H5LT.h>
73 #else
74 # include <hdf5_hl.h>
75 #endif
76 
77 #include "impex.hxx"
78 #include "multi_array.hxx"
79 #include "multi_iterator_coupled.hxx"
80 #include "multi_impex.hxx"
81 #include "utilities.hxx"
82 #include "error.hxx"
83 
84 #if defined(_MSC_VER)
85 # include <io.h>
86 #else
87 # include <unistd.h>
88 #endif
89 
90 namespace vigra {
91 
92 /** \addtogroup VigraHDF5Impex Import/Export of Images and Arrays in HDF5 Format
93 
94  Supports arrays with arbitrary element types and arbitrary many dimensions.
95  See the <a href="http://www.hdfgroup.org/HDF5/">HDF5 Website</a> for more
96  information on the HDF5 file format.
97 */
98 //@{
99 
100  /** \brief Check if given filename refers to a HDF5 file.
101  */
102 inline bool isHDF5(char const * filename)
103 {
104 #ifdef _MSC_VER
105  return _access(filename, 0) != -1 && H5Fis_hdf5(filename);
106 #else
107  return access(filename, F_OK) == 0 && H5Fis_hdf5(filename);
108 #endif
109 }
110 
111  /** \brief Temporarily disable HDF5's native error output.
112 
113  This should be used when you want to call an HDF5 function
114  that is known to fail (e.g. during testing), or when you want
115  to use an alternative error reporting mechanism (e.g. exceptions).
116 
117  <b>Usage:</b>
118 
119  <b>\#include</b> <vigra/hdf5impex.hxx><br>
120  Namespace: vigra
121  \code
122  {
123  HDF5DisableErrorOutput hdf5DisableErrorOutput;
124 
125  ... // call some HDF5 function
126 
127  } // restore the original error reporting in the destructor of HDF5DisableErrorOutput
128  \endcode
129  */
131 {
132  H5E_auto2_t old_func_;
133  void *old_client_data_;
134 
136  HDF5DisableErrorOutput & operator=(HDF5DisableErrorOutput const &);
137 
138  public:
140  : old_func_(0)
141  , old_client_data_(0)
142  {
143  H5Eget_auto2(H5E_DEFAULT, &old_func_, &old_client_data_);
144  H5Eset_auto2(H5E_DEFAULT, NULL, NULL);
145  }
146 
148  {
149  H5Eset_auto2(H5E_DEFAULT, old_func_, old_client_data_);
150  }
151 };
152 
153  /** \brief Wrapper for unique hid_t objects.
154 
155  This class offers the functionality of <tt>std::unique_ptr</tt> for HDF5 handles
156  (type <tt>hid_t</tt>). Unfortunately, <tt>std::unique_ptr</tt> cannot be used directly
157  for this purpose because it only works with pointers, whereas <tt>hid_t</tt> is an integer type.
158 
159  Newly created or opened HDF5 handles are stored as objects of type <tt>hid_t</tt>. When the handle
160  is no longer needed, the appropriate close function must be called. However, if a function is
161  aborted by an exception, this is difficult to ensure. Class HDF5Handle is a smart pointer that
162  solves this problem by calling the close function in the destructor (This is analogous to how
163  <tt>std::unique_ptr</tt> calls 'delete' on the contained pointer). A pointer to the close function
164  must be passed to the constructor, along with an error message that is raised when
165  creation/opening fails.
166 
167  When a <tt>HDF5Handle</tt> is created or assigned from another one, ownership passes on to the
168  left-hand-side handle object, and the right-hand-side objects is resest to a NULL handle.
169 
170  Since <tt>HDF5Handle</tt> objects are convertible to <tt>hid_t</tt>, they can be used in the code
171  in place of the latter.
172 
173  <b>Usage:</b>
174 
175  \code
176  HDF5Handle file_id(H5Fopen(filename, H5F_ACC_RDWR, H5P_DEFAULT),
177  &H5Fclose,
178  "Error message when H5Fopen() fails.");
179 
180  ... // use file_id in the same way as a plain hid_t object
181 
182  // pass ownership to a new handle object
183  HDF5Handle new_handle(file_id);
184 
185  assert(file_id.get() == 0);
186  \endcode
187 
188  <b>\#include</b> <vigra/hdf5impex.hxx><br>
189  Namespace: vigra
190  */
192 {
193 public:
194  typedef herr_t (*Destructor)(hid_t);
195 
196 private:
197  hid_t handle_;
198  Destructor destructor_;
199 
200 public:
201 
202  /** \brief Default constructor.
203  Creates a NULL handle.
204  **/
206  : handle_( 0 ),
207  destructor_(0)
208  {}
209 
210  /** \brief Create a wrapper for a hid_t object.
211 
212  The hid_t object \a h is assumed to be the return value of an open or create function.
213  It will be closed with the given close function \a destructor as soon as this
214  HDF5Handle is destructed, except when \a destructor is a NULL pointer (in which
215  case nothing happens at destruction time). If \a h has a value that indicates
216  failed opening or creation (by HDF5 convention, this means that \a h is negative),
217  an exception is raised by calling <tt>vigra_fail(error_message)</tt>.
218 
219  <b>Usage:</b>
220 
221  \code
222  HDF5Handle file_id(H5Fopen(filename, H5F_ACC_RDWR, H5P_DEFAULT),
223  &H5Fclose,
224  "Error message.");
225 
226  ... // use file_id in the same way
227  \endcode
228  */
229  HDF5Handle(hid_t h, Destructor destructor, const char * error_message)
230  : handle_( h ),
231  destructor_(destructor)
232  {
233  if(handle_ < 0)
234  vigra_fail(error_message);
235  }
236 
237  /** \brief Copy constructor.
238 
239  Hands over ownership of the RHS handle (analogous to <tt>std::unique_pt</tt>).
240  */
242  : handle_( h.handle_ ),
243  destructor_(h.destructor_)
244  {
245  const_cast<HDF5Handle &>(h).handle_ = 0;
246  }
247 
248  /** \brief Assignment.
249  Calls close() for the LHS handle and hands over ownership of the
250  RHS handle (analogous to <tt>std::unique_pt</tt>).
251  */
253  {
254  if(h.handle_ != handle_)
255  {
256  close();
257  handle_ = h.handle_;
258  destructor_ = h.destructor_;
259  const_cast<HDF5Handle &>(h).handle_ = 0;
260  }
261  return *this;
262  }
263 
264  /** \brief Destructor.
265  Calls close() for the contained handle.
266  */
268  {
269  close();
270  }
271 
272  /** \brief Explicitly call the stored destructor (if one has been registered in the
273  constructor) for the contained handle and set the wrapper to NULL. Returns
274  a negative value when the destructor call for the handle fails, and
275  a non-negative value otherwise.
276  */
277  herr_t close()
278  {
279  herr_t res = 1;
280  if(handle_ && destructor_)
281  res = (*destructor_)(handle_);
282  handle_ = 0;
283  destructor_ = 0;
284  return res;
285  }
286 
287  /** \brief Return the contained handle and set the wrapper to NULL
288  without calling <tt>close()</tt>.
289  */
290  hid_t release()
291  {
292  hid_t res = handle_;
293  handle_ = 0;
294  destructor_ = 0;
295  return res;
296  }
297 
298  /** \brief Reset the wrapper to a new handle.
299 
300  Equivalent to <tt>handle = HDF5Handle(h, destructor, error_message)</tt>.
301  */
302  void reset(hid_t h, Destructor destructor, const char * error_message)
303  {
304  if(h < 0)
305  vigra_fail(error_message);
306  if(h != handle_)
307  {
308  close();
309  handle_ = h;
310  destructor_ = destructor;
311  }
312  }
313 
314  /** \brief Swap the contents of two handle wrappers.
315 
316  Also available as <tt>std::swap(handle1, handle2)</tt>.
317  */
318  void swap(HDF5Handle & h)
319  {
320  std::swap(handle_, h.handle_);
321  std::swap(destructor_, h.destructor_);
322  }
323 
324  /** \brief Get a temporary hid_t object for the contained handle.
325  Do not call a close function on the return value - a crash will be likely
326  otherwise.
327  */
328  hid_t get() const
329  {
330  return handle_;
331  }
332 
333  /** \brief Convert to a plain hid_t object.
334 
335  This function ensures that hid_t objects can be transparently replaced with
336  HDF5Handle objects in user code. Do not call a close function on the return
337  value - a crash will be likely otherwise.
338  */
339  operator hid_t() const
340  {
341  return handle_;
342  }
343 
344  /** \brief Equality comparison of the contained handle.
345  */
346  bool operator==(HDF5Handle const & h) const
347  {
348  return handle_ == h.handle_;
349  }
350 
351  /** \brief Equality comparison of the contained handle.
352  */
353  bool operator==(hid_t h) const
354  {
355  return handle_ == h;
356  }
357 
358  /** \brief Inequality comparison of the contained handle.
359  */
360  bool operator!=(HDF5Handle const & h) const
361  {
362  return handle_ != h.handle_;
363  }
364 
365  /** \brief Inequality comparison of the contained handle.
366  */
367  bool operator!=(hid_t h) const
368  {
369  return handle_ != h;
370  }
371 };
372 
373 
374  /** \brief Wrapper for shared hid_t objects.
375 
376  This class offers the functionality of <tt>std::shared_ptr</tt> for HDF5 handles
377  (type <tt>hid_t</tt>). Unfortunately, <tt>std::shared_ptr</tt> cannot be used directly
378  for this purpose because it only works with pointers, whereas <tt>hid_t</tt> is an integer type.
379 
380  Newly created or opened HDF5 handles are stored as objects of type <tt>hid_t</tt>. When the handle
381  is no longer needed, the appropriate close function must be called. However, if a function is
382  aborted by an exception, this is difficult to ensure. Class HDF5HandleShared is a smart pointer
383  that solves this problem by calling the close function in the destructor of the handle's last
384  owner (This is analogous to how <tt>std::shared_ptr</tt> calls 'delete' on the contained
385  pointer). A pointer to the close function must be passed to the constructor, along with an error
386  message that is raised when creation/opening fails.
387 
388  When a <tt>HDF5HandleShared</tt> is created or assigned from another one, ownership is shared
389  between the two handles, and the value returned by <tt>use_count()</tt> increases by one.
390 
391  Since <tt>HDF5HandleShared</tt> objects are convertible to <tt>hid_t</tt>, they can be used in the code
392  in place of the latter.
393 
394  <b>Usage:</b>
395 
396  \code
397  HDF5HandleShared file_id(H5Fopen(filename, H5F_ACC_RDWR, H5P_DEFAULT),
398  &H5Fclose,
399  "Error message when H5Fopen() fails.");
400 
401  ... // use file_id in the same way as a plain hid_t object
402 
403  // share ownership between same_id and file_id
404  HDF5HandleShared same_id(file_id);
405  assert(same_id.use_count() == 2);
406  assert(same_id.get() == file_id.get());
407  \endcode
408 
409  <b>\#include</b> <vigra/hdf5impex.hxx><br>
410  Namespace: vigra
411  */
413 {
414 public:
415  typedef herr_t (*Destructor)(hid_t);
416 
417 private:
418  hid_t handle_;
419  Destructor destructor_;
420  size_t * refcount_;
421 
422 public:
423 
424  /** \brief Default constructor.
425  Creates a NULL handle.
426  **/
428  : handle_( 0 ),
429  destructor_(0),
430  refcount_(0)
431  {}
432 
433  /** \brief Create a shared wrapper for a plain hid_t object.
434 
435  The hid_t object \a h is assumed to be the return value of an open or create function.
436  It will be closed with the given close function \a destructor as soon as this
437  HDF5HandleShared is destructed, except when \a destructor is a NULL pointer (in which
438  case nothing happens at destruction time). If \a h has a value that indicates
439  failed opening or creation (by HDF5 convention, this means that \a h is negative),
440  an exception is raised by calling <tt>vigra_fail(error_message)</tt>.
441 
442  <b>Usage:</b>
443 
444  \code
445  HDF5HandleShared file_id(H5Fopen(filename, H5F_ACC_RDWR, H5P_DEFAULT),
446  &H5Fclose,
447  "Error message.");
448 
449  ... // use file_id in the same way
450  \endcode
451  */
452  HDF5HandleShared(hid_t h, Destructor destructor, const char * error_message)
453  : handle_( h ),
454  destructor_(destructor),
455  refcount_(0)
456  {
457  if(handle_ < 0)
458  vigra_fail(error_message);
459  if(handle_ > 0)
460  refcount_ = new size_t(1);
461  }
462 
463  /** \brief Copy constructor.
464  Shares ownership with the RHS handle (analogous to <tt>std::shared_ptr</tt>).
465  */
467  : handle_( h.handle_ ),
468  destructor_(h.destructor_),
469  refcount_(h.refcount_)
470  {
471  if(refcount_)
472  ++(*refcount_);
473  }
474 
475  /** \brief Assignment.
476  Call close() for the present LHS handle and share ownership with the
477  RHS handle (analogous to <tt>std::shared_ptr</tt>).
478  */
480  {
481  if(h.handle_ != handle_)
482  {
483  close();
484  handle_ = h.handle_;
485  destructor_ = h.destructor_;
486  refcount_ = h.refcount_;
487  if(refcount_)
488  ++(*refcount_);
489  }
490  return *this;
491  }
492 
493  /** \brief Destructor (calls close())
494  */
496  {
497  close();
498  }
499 
500  /** \brief Close the handle if this is the unique (i.e. last) owner.
501 
502  Decrements the reference counter and calls the destructor function of
503  the handle (if one has been registered in the constructor) when the counter
504  reaches zero. Sets this wrapper to NULL in any case. Returns
505  a negative value when the destructor call for the handle fails, and
506  a non-negative value otherwise.
507  */
508  herr_t close()
509  {
510  herr_t res = 1;
511  if(refcount_)
512  {
513  --(*refcount_);
514  if(*refcount_ == 0)
515  {
516  if(destructor_)
517  res = (*destructor_)(handle_);
518  delete refcount_;
519  }
520  }
521  handle_ = 0;
522  destructor_ = 0;
523  refcount_ = 0;
524  return res;
525  }
526 
527  /** \brief Reset the handle to a new value.
528 
529  Equivalent to <tt>handle = HDF5HandleShared(h, destructor, error_message)</tt>.
530  */
531  void reset(hid_t h, Destructor destructor, const char * error_message)
532  {
533  if(h < 0)
534  vigra_fail(error_message);
535  if(h != handle_)
536  {
537  close();
538  handle_ = h;
539  destructor_ = destructor;
540  if(handle_ > 0)
541  refcount_ = new size_t(1);
542  }
543  }
544 
545  /** \brief Get the number of owners of the contained handle.
546  */
547  size_t use_count() const
548  {
549  return refcount_
550  ? *refcount_
551  : 0;
552  }
553 
554  /** \brief Check if this is the unique owner of the conatined handle.
555 
556  Equivalent to <tt>handle.use_count() == 1</tt>.
557  */
558  bool unique() const
559  {
560  return use_count() == 1;
561  }
562 
563  /** \brief Swap the contents of two handle wrappers.
564 
565  Also available as <tt>std::swap(handle1, handle2)</tt>.
566  */
568  {
569  std::swap(handle_, h.handle_);
570  std::swap(destructor_, h.destructor_);
571  std::swap(refcount_, h.refcount_);
572  }
573 
574  /** \brief Get a temporary hid_t object for the contained handle.
575  Do not call a close function on the return value - a crash will be likely
576  otherwise.
577  */
578  hid_t get() const
579  {
580  return handle_;
581  }
582 
583  /** \brief Convert to a plain hid_t object.
584 
585  This function ensures that hid_t objects can be transparently replaced with
586  HDF5HandleShared objects in user code. Do not call a close function on the return
587  value - a crash will be likely otherwise.
588  */
589  operator hid_t() const
590  {
591  return handle_;
592  }
593 
594  /** \brief Equality comparison of the contained handle.
595  */
596  bool operator==(HDF5HandleShared const & h) const
597  {
598  return handle_ == h.handle_;
599  }
600 
601  /** \brief Equality comparison of the contained handle.
602  */
603  bool operator==(hid_t h) const
604  {
605  return handle_ == h;
606  }
607 
608  /** \brief Inequality comparison of the contained handle.
609  */
610  bool operator!=(HDF5HandleShared const & h) const
611  {
612  return handle_ != h.handle_;
613  }
614 
615  /** \brief Inequality comparison of the contained handle.
616  */
617  bool operator!=(hid_t h) const
618  {
619  return handle_ != h;
620  }
621 };
622 
623 //@}
624 
625 } // namespace vigra
626 
627 namespace std {
628 
629 inline void swap(vigra::HDF5Handle & l, vigra::HDF5Handle & r)
630 {
631  l.swap(r);
632 }
633 
634 inline void swap(vigra::HDF5HandleShared & l, vigra::HDF5HandleShared & r)
635 {
636  l.swap(r);
637 }
638 
639 } // namespace std
640 
641 namespace vigra {
642 
643 /** \addtogroup VigraHDF5Impex
644 */
645 //@{
646 
647 
648 /********************************************************/
649 /* */
650 /* HDF5ImportInfo */
651 /* */
652 /********************************************************/
653 
654 /** \brief Argument object for the function readHDF5().
655 
656 See \ref readHDF5() for a usage example. This object must be
657 used to read an image or array from an HDF5 file
658 and enquire about its properties.
659 
660 <b>\#include</b> <vigra/hdf5impex.hxx><br>
661 Namespace: vigra
662 */
664 {
665  public:
666  enum PixelType { UINT8, UINT16, UINT32, UINT64,
667  INT8, INT16, INT32, INT64,
668  FLOAT, DOUBLE };
669 
670  /** Construct HDF5ImportInfo object.
671 
672  The dataset \a pathInFile in the HDF5 file \a filename is accessed to
673  read its properties. \a pathInFile may contain '/'-separated group
674  names, but must end with the name of the desired dataset:
675 
676  \code
677  HDF5ImportInfo info(filename, "/group1/group2/my_dataset");
678  \endcode
679  */
680  VIGRA_EXPORT HDF5ImportInfo( const char* filePath, const char* pathInFile );
681 
682  VIGRA_EXPORT ~HDF5ImportInfo();
683 
684  /** Get the filename of this HDF5 object.
685  */
686  VIGRA_EXPORT const std::string& getFilePath() const;
687 
688  /** Get the dataset's full name in the HDF5 file.
689  */
690  VIGRA_EXPORT const std::string& getPathInFile() const;
691 
692  /** Get a handle to the file represented by this info object.
693  */
694  VIGRA_EXPORT hid_t getH5FileHandle() const;
695 
696  /** Get a handle to the dataset represented by this info object.
697  */
698  VIGRA_EXPORT hid_t getDatasetHandle() const;
699 
700  /** Get the number of dimensions of the dataset represented by this info object.
701  */
702  VIGRA_EXPORT MultiArrayIndex numDimensions() const;
703 
704  /** Get the shape of the dataset represented by this info object.
705 
706  Note that the memory order between VIGRA and HDF5 files differs: VIGRA uses
707  Fortran-order, while HDF5 uses C-order. This function therefore reverses the axis
708  order relative to the file contents. That is, when the axes in the file are
709  ordered as 'z', 'y', 'x', this function will return the shape in the order
710  'x', 'y', 'z'.
711  */
712  VIGRA_EXPORT ArrayVector<hsize_t> const & shape() const
713  {
714  return m_dims;
715  }
716 
717  /** Get the shape (length) of the dataset along dimension \a dim.
718 
719  Note that the memory order between VIGRA and HDF5 files differs: VIGRA uses
720  Fortran-order, while HDF5 uses C-order. This function therefore reverses the axis
721  order relative to the file contents. That is, when the axes in the file are
722  ordered as 'z', 'y', 'x', this function will return the shape in the order
723  'x', 'y', 'z'.
724  */
725  VIGRA_EXPORT MultiArrayIndex shapeOfDimension(const int dim) const;
726 
727  /** Query the pixel type of the dataset.
728 
729  Possible values are:
730  <DL>
731  <DT>"INT8"<DD> 8-bit signed integer (unsigned char)
732  <DT>"UINT8"<DD> 8-bit unsigned integer (unsigned char)
733  <DT>"INT16"<DD> 16-bit signed integer (short)
734  <DT>"UINT16"<DD> 16-bit unsigned integer (unsigned short)
735  <DT>"INT32"<DD> 32-bit signed integer (long)
736  <DT>"UINT32"<DD> 32-bit unsigned integer (unsigned long)
737  <DT>"INT64"<DD> 64-bit signed integer (long long)
738  <DT>"UINT64"<DD> 64-bit unsigned integer (unsigned long long)
739  <DT>"FLOAT"<DD> 32-bit floating point (float)
740  <DT>"DOUBLE"<DD> 64-bit floating point (double)
741  </DL>
742  */
743  VIGRA_EXPORT const char * getPixelType() const;
744 
745  /** Query the pixel type of the dataset.
746 
747  Same as getPixelType(), but the result is returned as a
748  ImageImportInfo::PixelType enum. This is useful to implement
749  a switch() on the pixel type.
750 
751  Possible values are:
752  <DL>
753  <DT>UINT8<DD> 8-bit unsigned integer (unsigned char)
754  <DT>INT16<DD> 16-bit signed integer (short)
755  <DT>UINT16<DD> 16-bit unsigned integer (unsigned short)
756  <DT>INT32<DD> 32-bit signed integer (long)
757  <DT>UINT32<DD> 32-bit unsigned integer (unsigned long)
758  <DT>FLOAT<DD> 32-bit floating point (float)
759  <DT>DOUBLE<DD> 64-bit floating point (double)
760  </DL>
761  */
762  VIGRA_EXPORT PixelType pixelType() const;
763 
764  private:
765  HDF5HandleShared m_file_handle, m_dataset_handle;
766  std::string m_filename, m_path, m_pixeltype;
767  hssize_t m_dimensions;
768  ArrayVector<hsize_t> m_dims;
769 };
770 
771 
772 namespace detail {
773 
774 template <class T>
775 struct HDF5TypeTraits
776 {
777  static hid_t getH5DataType()
778  {
779  std::runtime_error("getH5DataType(): invalid type");
780  return 0;
781  }
782 
783  static int numberOfBands()
784  {
785  std::runtime_error("numberOfBands(): invalid type");
786  return 0;
787  }
788 };
789 
790 template<class T>
791 inline hid_t getH5DataType()
792 {
793  return HDF5TypeTraits<T>::getH5DataType();
794 }
795 
796 #define VIGRA_H5_DATATYPE(type, h5type) \
797 template <> \
798 struct HDF5TypeTraits<type> \
799 { \
800  static hid_t getH5DataType() \
801  { \
802  return h5type; \
803  } \
804  static int numberOfBands() \
805  { \
806  return 1; \
807  } \
808  typedef type value_type; \
809 }; \
810 template <int M> \
811 struct HDF5TypeTraits<TinyVector<type, M> > \
812 { \
813  static hid_t getH5DataType() \
814  { \
815  return h5type; \
816  } \
817  static int numberOfBands() \
818  { \
819  return M; \
820  } \
821  typedef type value_type; \
822 }; \
823 template <> \
824 struct HDF5TypeTraits<RGBValue<type> > \
825 { \
826  static hid_t getH5DataType() \
827  { \
828  return h5type; \
829  } \
830  static int numberOfBands() \
831  { \
832  return 3; \
833  } \
834  typedef type value_type; \
835 }; \
836 
837 VIGRA_H5_DATATYPE(char, H5T_NATIVE_CHAR)
838 VIGRA_H5_DATATYPE(signed char, H5T_NATIVE_SCHAR)
839 VIGRA_H5_DATATYPE(unsigned char, H5T_NATIVE_UCHAR)
840 VIGRA_H5_DATATYPE(signed short, H5T_NATIVE_SHORT)
841 VIGRA_H5_DATATYPE(unsigned short, H5T_NATIVE_USHORT)
842 VIGRA_H5_DATATYPE(signed int, H5T_NATIVE_INT)
843 VIGRA_H5_DATATYPE(unsigned int, H5T_NATIVE_UINT)
844 VIGRA_H5_DATATYPE(signed long, H5T_NATIVE_LONG)
845 VIGRA_H5_DATATYPE(unsigned long, H5T_NATIVE_ULONG)
846 VIGRA_H5_DATATYPE(signed long long, H5T_NATIVE_LLONG)
847 VIGRA_H5_DATATYPE(unsigned long long, H5T_NATIVE_ULLONG)
848 VIGRA_H5_DATATYPE(float, H5T_NATIVE_FLOAT)
849 VIGRA_H5_DATATYPE(double, H5T_NATIVE_DOUBLE)
850 VIGRA_H5_DATATYPE(long double, H5T_NATIVE_LDOUBLE)
851 
852 // char arrays with flexible length require 'handcrafted' H5 datatype
853 template <>
854 struct HDF5TypeTraits<char*>
855 {
856  static hid_t getH5DataType()
857  {
858  hid_t stringtype = H5Tcopy (H5T_C_S1);
859  H5Tset_size(stringtype, H5T_VARIABLE);
860  return stringtype;
861  }
862 
863  static int numberOfBands()
864  {
865  return 1;
866  }
867 };
868 
869 template <>
870 struct HDF5TypeTraits<const char*>
871 {
872  static hid_t getH5DataType()
873  {
874  hid_t stringtype = H5Tcopy (H5T_C_S1);
875  H5Tset_size(stringtype, H5T_VARIABLE);
876  return stringtype;
877  }
878 
879  static int numberOfBands()
880  {
881  return 1;
882  }
883 };
884 
885 #undef VIGRA_H5_DATATYPE
886 
887 #if 0
888 template<>
889 inline hid_t getH5DataType<FFTWComplex<float> >()
890 {
891  hid_t complex_id = H5Tcreate (H5T_COMPOUND, sizeof (FFTWComplex<float>));
892  H5Tinsert (complex_id, "real", 0, H5T_NATIVE_FLOAT);
893  H5Tinsert (complex_id, "imaginary", sizeof(float), H5T_NATIVE_FLOAT);
894  return complex_id;
895 }
896 
897 template<>
898 inline hid_t getH5DataType<FFTWComplex<double> >()
899 {
900  hid_t complex_id = H5Tcreate (H5T_COMPOUND, sizeof (FFTWComplex<double>));
901  H5Tinsert (complex_id, "real", 0, H5T_NATIVE_DOUBLE);
902  H5Tinsert (complex_id, "imaginary", sizeof(double), H5T_NATIVE_DOUBLE);
903  return complex_id;
904 }
905 #endif
906 
907 
908 } // namespace detail
909 
910 // helper friend function for callback HDF5_ls_inserter_callback()
911 void HDF5_ls_insert(void*, const std::string &);
912 // callback function for ls(), called via HDF5File::H5Literate()
913 // see http://www.parashift.com/c++-faq-lite/pointers-to-members.html#faq-33.2
914 // for as to why.
915 
916 VIGRA_EXPORT H5O_type_t HDF5_get_type(hid_t, const char*);
917 extern "C" VIGRA_EXPORT herr_t HDF5_ls_inserter_callback(hid_t, const char*, const H5L_info_t*, void*);
918 extern "C" VIGRA_EXPORT herr_t HDF5_listAttributes_inserter_callback(hid_t, const char*, const H5A_info_t*, void*);
919 
920 /********************************************************/
921 /* */
922 /* HDF5File */
923 /* */
924 /********************************************************/
925 
926 
927 /** \brief Access to HDF5 files
928 
929 HDF5File provides a convenient way of accessing data in HDF5 files. vigra::MultiArray
930 structures of any dimension can be stored to / loaded from HDF5 files. Typical
931 HDF5 features like subvolume access, chunks and data compression are available,
932 string attributes can be attached to any dataset or group. Group- or dataset-handles
933 are encapsulated in the class and managed automatically. The internal file-system like
934 structure can be accessed by functions like "cd()" or "mkdir()".
935 
936 Note that the memory order between VIGRA and HDF5 files differs: VIGRA uses
937 Fortran-order, while HDF5 uses C-order. This means that a VIGRA MultiArray,
938 whose indices represent the 'x'-, 'y'-, and 'z'-axis in that order, is reversed
939 upon writing to an HDF5 file, i.e. in the file the axis order is 'z', 'y', 'x'.
940 Likewise, the order is reversed upon reading.
941 
942 <b>Example:</b>
943 Write the MultiArray out_multi_array to file. Change the current directory to
944 "/group" and read in the same MultiArray as in_multi_array.
945 \code
946 HDF5File file("/path/to/file",HDF5File::New);
947 file.mkdir("group");
948 file.write("/group/dataset", out_multi_array);
949 
950 file.cd("/group");
951 file.read("dataset", in_multi_array);
952 
953 \endcode
954 
955 <b>\#include</b> <vigra/hdf5impex.hxx><br>
956 Namespace: vigra
957 */
958 class HDF5File
959 {
960  protected:
961  HDF5HandleShared fileHandle_;
962 
963  // current group handle
964  HDF5Handle cGroupHandle_;
965 
966  private:
967  // time tagging of datasets, turned off (= 0) by default.
968  int track_time;
969 
970  bool read_only_;
971 
972  // helper classes for ls() and listAttributes()
973  struct ls_closure
974  {
975  virtual void insert(const std::string &) = 0;
976  virtual ~ls_closure() {}
977  };
978 
979  // datastructure to hold a std::vector<std::string>
980  struct lsOpData : public ls_closure
981  {
982  std::vector<std::string> & objects;
983  lsOpData(std::vector<std::string> & o) : objects(o) {}
984  void insert(const std::string & x)
985  {
986  objects.push_back(x);
987  }
988  };
989 
990  // datastructure to hold an associative container
991  template<class Container>
992  struct ls_container_data : public ls_closure
993  {
994  Container & objects;
995  ls_container_data(Container & o) : objects(o) {}
996  void insert(const std::string & x)
997  {
998  objects.insert(std::string(x));
999  }
1000  };
1001 
1002  public:
1003 
1004  // helper for callback HDF5_ls_inserter_callback(), used by ls()
1005  friend void HDF5_ls_insert(void*, const std::string &);
1006 
1007  /** \brief Set how a file is opened.
1008 
1009  OpenMode::New creates a new file. If the file already exists, it is overwritten.
1010 
1011  OpenMode::ReadWrite opens a file for reading/writing. The file will be created if it doesn't exist.
1012 
1013  OpenMode::ReadOnly opens a file for reading. The file as well as any dataset to be accessed must already exist.
1014  */
1015  enum OpenMode {
1016  New, // Create new empty file (existing file will be deleted).
1017  Open, // Open file. Create if not existing.
1018  ReadWrite = Open, // Alias for Open.
1019  OpenReadOnly, // Open file in read-only mode.
1020  ReadOnly = OpenReadOnly, // Alias for OpenReadOnly
1021  Replace, // for ChunkedArrayHDF5: replace dataset if it exists, create otherwise
1022  Default // for ChunkedArrayHDF5: use New if file doesn't exist,
1023  // ReadOnly if file and dataset exist
1024  // Open otherwise
1025  };
1026 
1027  /** \brief Default constructor.
1028 
1029  A file can later be opened via the open() function. Time tagging of datasets is disabled.
1030  */
1032  : track_time(0)
1033  {}
1034 
1035  /** \brief Construct with time tagging of datasets enabled.
1036 
1037  If \a track_creation_times is 'true', time tagging of datasets will be enabled.
1038  */
1039  explicit HDF5File(bool track_creation_times)
1040  : track_time(track_creation_times ? 1 : 0),
1041  read_only_(true)
1042  {}
1043 
1044  /** \brief Open or create an HDF5File object.
1045 
1046  Creates or opens HDF5 file with given filename.
1047  The current group is set to "/".
1048 
1049  Note that the HDF5File class is not copyable (the copy constructor is
1050  private to enforce this).
1051  */
1052  HDF5File(std::string filePath, OpenMode mode, bool track_creation_times = false)
1053  : track_time(track_creation_times ? 1 : 0)
1054  {
1055  open(filePath, mode);
1056  }
1057 
1058  /** \brief Initialize an HDF5File object from HDF5 file handle
1059 
1060  Initializes an HDF5File object corresponding to the HDF5 file
1061  opened elsewhere. If \a fileHandle is constructed with a
1062  <tt>NULL</tt> destructor, ownership is not transferred
1063  to the new HDF5File object, and you must ensure that the file is
1064  not closed while the new HDF5File object is in use. Otherwise,
1065  ownership will be shared.
1066 
1067  The current group is set to the specified \a pathname. If
1068  \a read_only is 'true', you cannot create new datasets or
1069  overwrite data.
1070 
1071  \warning In case the underlying HDF5 library used by Vigra is not
1072  exactly the same library used to open the file with the given id,
1073  this method will lead to crashes.
1074  */
1075  explicit HDF5File(HDF5HandleShared const & fileHandle,
1076  const std::string & pathname = "",
1077  bool read_only = false)
1078  : fileHandle_(fileHandle),
1079  read_only_(read_only)
1080 
1081  {
1082  // get group handle for given pathname
1083  // calling openCreateGroup_ without setting a valid cGroupHandle does
1084  // not work. Starting from root() is a safe bet.
1085  root();
1086  cGroupHandle_ = HDF5Handle(openCreateGroup_(pathname), &H5Gclose,
1087  "HDF5File(fileHandle, pathname): Failed to open group");
1088 
1089  // extract track_time attribute
1090  hbool_t track_times_tmp;
1091  HDF5Handle plist_id(H5Fget_create_plist(fileHandle_), &H5Pclose,
1092  "HDF5File(fileHandle, pathname): Failed to open file creation property list");
1093  herr_t status = H5Pget_obj_track_times(plist_id, &track_times_tmp );
1094  vigra_postcondition(status >= 0,
1095  "HDF5File(fileHandle, pathname): cannot access track time attribute");
1096  track_time = track_times_tmp;
1097  }
1098 
1099  /** \brief Copy a HDF5File object.
1100 
1101  The new object will refer to the same file and group as \a other.
1102  */
1103  HDF5File(HDF5File const & other)
1104  : fileHandle_(other.fileHandle_),
1105  track_time(other.track_time),
1106  read_only_(other.read_only_)
1107  {
1108  cGroupHandle_ = HDF5Handle(openCreateGroup_(other.currentGroupName_()), &H5Gclose,
1109  "HDF5File(HDF5File const &): Failed to open group.");
1110  }
1111 
1112  /** \brief The destructor flushes and closes the file.
1113  */
1115  {
1116  // The members fileHandle_ and cGroupHandle_ are automatically closed
1117  // as they are of type HDF5Handle and are properly initialised.
1118  // The closing of fileHandle_ implies flushing the file to
1119  // the operating system, see
1120  // http://www.hdfgroup.org/HDF5/doc/RM/RM_H5F.html#File-Close .
1121  }
1122 
1123  /** \brief Assign a HDF5File object.
1124 
1125  Calls close() on the present file and The new object will refer to the same file and group as \a other.
1126  */
1127  HDF5File & operator=(HDF5File const & other)
1128  {
1129  if(this != &other)
1130  {
1131  close();
1132  fileHandle_ = other.fileHandle_;
1133  cGroupHandle_ = HDF5Handle(openCreateGroup_(other.currentGroupName_()), &H5Gclose,
1134  "HDF5File::operator=(): Failed to open group.");
1135  track_time = other.track_time;
1136  read_only_ = other.read_only_;
1137  }
1138  return *this;
1139  }
1140 
1141  int file_use_count() const
1142  {
1143  return fileHandle_.use_count();
1144  }
1145 
1146  bool isOpen() const
1147  {
1148  return fileHandle_ != 0;
1149  }
1150 
1151  bool isReadOnly() const
1152  {
1153  return read_only_;
1154  }
1155 
1156  void setReadOnly(bool stat=true)
1157  {
1158  read_only_ = stat;
1159  }
1160 
1161  /** \brief Open or create the given file in the given mode and set the group to "/".
1162  If another file is currently open, it is first closed.
1163  */
1164  void open(std::string filePath, OpenMode mode)
1165  {
1166  close();
1167 
1168  std::string errorMessage = "HDF5File.open(): Could not open or create file '" + filePath + "'.";
1169  fileHandle_ = HDF5HandleShared(createFile_(filePath, mode), &H5Fclose, errorMessage.c_str());
1170  cGroupHandle_ = HDF5Handle(openCreateGroup_("/"), &H5Gclose, "HDF5File.open(): Failed to open root group.");
1171  setReadOnly(mode == OpenReadOnly);
1172  }
1173 
1174  /** \brief Close the current file.
1175  */
1176  void close()
1177  {
1178  bool success = cGroupHandle_.close() >= 0 && fileHandle_.close() >= 0;
1179  vigra_postcondition(success, "HDF5File.close() failed.");
1180  }
1181 
1182  /** \brief Change current group to "/".
1183  */
1184  inline void root()
1185  {
1186  std::string message = "HDF5File::root(): Could not open group '/'.";
1187  cGroupHandle_ = HDF5Handle(H5Gopen(fileHandle_, "/", H5P_DEFAULT),&H5Gclose,message.c_str());
1188  }
1189 
1190  /** \brief Change the current group.
1191  Both absolute and relative group names are allowed.
1192  */
1193  inline void cd(std::string groupName)
1194  {
1195  cGroupHandle_ = getGroupHandle(groupName, "HDF5File::cd()");
1196  }
1197 
1198  /** \brief Change the current group to its parent group.
1199  Returns true if successful, false otherwise. If unsuccessful,
1200  the group will not change.
1201  */
1202  inline bool cd_up()
1203  {
1204  std::string groupName = currentGroupName_();
1205 
1206  //do not try to move up if we already in "/"
1207  if(groupName == "/"){
1208  return false;
1209  }
1210 
1211  size_t lastSlash = groupName.find_last_of('/');
1212 
1213  std::string parentGroup (groupName.begin(), groupName.begin()+lastSlash+1);
1214 
1215  cd(parentGroup);
1216 
1217  return true;
1218  }
1219 
1220  /** \brief Change the current group to its parent group.
1221  Returns true if successful, false otherwise. If unsuccessful,
1222  the group will not change.
1223  */
1224  inline bool cd_up(int levels)
1225  {
1226  std::string groupName = currentGroupName_();
1227 
1228  for(int i = 0; i<levels; i++)
1229  {
1230  if(!cd_up())
1231  {
1232  // restore old group if neccessary
1233  if(groupName != currentGroupName_())
1234  cd(groupName);
1235  return false;
1236  }
1237  }
1238  return true;
1239  }
1240 
1241  /** \brief Create a new group.
1242  If the first character is a "/", the path will be interpreted as absolute path,
1243  otherwise it will be interpreted as path relative to the current group.
1244  */
1245  inline void mkdir(std::string groupName)
1246  {
1247  vigra_precondition(!isReadOnly(),
1248  "HDF5File::mkdir(): file is read-only.");
1249 
1250  std::string message = "HDF5File::mkdir(): Could not create group '" + groupName + "'.\n";
1251 
1252  // make groupName clean
1253  groupName = get_absolute_path(groupName);
1254 
1255  HDF5Handle(openCreateGroup_(groupName.c_str()),&H5Gclose,message.c_str());
1256  }
1257 
1258  /** \brief Change the current group; create it if necessary.
1259  If the first character is a "/", the path will be interpreted as absolute path,
1260  otherwise it will be interpreted as path relative to the current group.
1261  */
1262  inline void cd_mk(std::string groupName)
1263  {
1264  vigra_precondition(!isReadOnly(),
1265  "HDF5File::cd_mk(): file is read-only.");
1266 
1267  std::string message = "HDF5File::cd_mk(): Could not create group '" + groupName + "'.";
1268 
1269  // make groupName clean
1270  groupName = get_absolute_path(groupName);
1271 
1272  cGroupHandle_ = HDF5Handle(openCreateGroup_(groupName.c_str()),&H5Gclose,message.c_str());
1273  }
1274 
1275  // helper function for the various ls() variants.
1276  void ls_H5Literate(ls_closure & data) const
1277  {
1278  H5Literate(cGroupHandle_, H5_INDEX_NAME, H5_ITER_NATIVE, NULL,
1279  HDF5_ls_inserter_callback, static_cast<void*>(&data));
1280  }
1281 
1282  /** \brief List the contents of the current group.
1283  The function returns a vector of strings holding the entries of the
1284  current group. Only datasets and groups are listed, other objects
1285  (e.g. datatypes) are ignored. Group names always have a trailing "/".
1286  */
1287  inline std::vector<std::string> ls() const
1288  {
1289  std::vector<std::string> list;
1290  lsOpData data(list);
1291  ls_H5Literate(data);
1292  return list;
1293  }
1294 
1295  /** \brief List the contents of the current group into a container-like
1296  object via insert().
1297 
1298  Only datasets and groups are inserted, other objects (e.g., datatypes) are ignored.
1299  Group names always have a trailing "/".
1300 
1301  The argument cont is presumably an associative container, however,
1302  only its member function <tt>cont.insert(std::string)</tt> will be
1303  called.
1304  \param cont reference to a container supplying a member function
1305  <tt>insert(const i_type &)</tt>, where <tt>i_type</tt>
1306  is convertible to <tt>std::string</tt>.
1307  */
1308  template<class Container>
1309  void ls(Container & cont) const
1310  {
1311  ls_container_data<Container> data(cont);
1312  ls_H5Literate(data);
1313  }
1314 
1315  /** \brief Get the path of the current group.
1316  */
1317  inline std::string pwd() const
1318  {
1319  return currentGroupName_();
1320  }
1321 
1322  /** \brief Get the name of the associated file.
1323  */
1324  inline std::string filename() const
1325  {
1326  return fileName_();
1327  }
1328 
1329  /** \brief Check if given datasetName exists.
1330  */
1331  inline bool existsDataset(std::string datasetName) const
1332  {
1333  // make datasetName clean
1334  datasetName = get_absolute_path(datasetName);
1335  return (H5Lexists(fileHandle_, datasetName.c_str(), H5P_DEFAULT) > 0);
1336  }
1337 
1338  /** \brief Get the number of dimensions of a certain dataset
1339  If the first character is a "/", the path will be interpreted as absolute path,
1340  otherwise it will be interpreted as path relative to the current group.
1341  */
1342  hssize_t getDatasetDimensions(std::string datasetName) const
1343  {
1344  HDF5Handle datasetHandle = getDatasetHandle(datasetName);
1345 
1346  return getDatasetDimensions_(datasetHandle);
1347  }
1348 
1349  hssize_t getDatasetDimensions_(hid_t dataset) const
1350  {
1351  std::string errorMessage = "HDF5File::getDatasetDimensions(): Unable to access dataspace.";
1352  HDF5Handle dataspaceHandle(H5Dget_space(dataset), &H5Sclose, errorMessage.c_str());
1353 
1354  //return dimension information
1355  return H5Sget_simple_extent_ndims(dataspaceHandle);
1356  }
1357 
1358  /** \brief Get the shape of each dimension of a certain dataset.
1359 
1360  Normally, this function is called after determining the dimension of the
1361  dataset using \ref getDatasetDimensions().
1362  If the first character is a "/", the path will be interpreted as absolute path,
1363  otherwise it will be interpreted as path relative to the current group.
1364 
1365  Note that the memory order between VIGRA and HDF5 files differs: VIGRA uses
1366  Fortran-order, while HDF5 uses C-order. This function therefore reverses the axis
1367  order relative to the file contents. That is, when the axes in the file are
1368  ordered as 'z', 'y', 'x', this function will return the shape in the order
1369  'x', 'y', 'z'.
1370  */
1371  ArrayVector<hsize_t> getDatasetShape(std::string datasetName) const
1372  {
1373  // make datasetName clean
1374  datasetName = get_absolute_path(datasetName);
1375 
1376  //Open dataset and dataspace
1377  std::string errorMessage = "HDF5File::getDatasetShape(): Unable to open dataset '" + datasetName + "'.";
1378  HDF5Handle datasetHandle = HDF5Handle(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
1379 
1380  errorMessage = "HDF5File::getDatasetShape(): Unable to access dataspace.";
1381  HDF5Handle dataspaceHandle(H5Dget_space(datasetHandle), &H5Sclose, errorMessage.c_str());
1382 
1383  //get dimension information
1384  ArrayVector<hsize_t>::size_type dimensions = H5Sget_simple_extent_ndims(dataspaceHandle);
1385 
1386  ArrayVector<hsize_t> shape(dimensions);
1387  ArrayVector<hsize_t> maxdims(dimensions);
1388  H5Sget_simple_extent_dims(dataspaceHandle, shape.data(), maxdims.data());
1389 
1390  // invert the dimensions to guarantee VIGRA-compatible order.
1391  std::reverse(shape.begin(), shape.end());
1392  return shape;
1393  }
1394 
1395  /** Query the pixel type of the dataset.
1396 
1397  Possible values are:
1398  <DL>
1399  <DT>"INT8"<DD> 8-bit signed integer (unsigned char)
1400  <DT>"UINT8"<DD> 8-bit unsigned integer (unsigned char)
1401  <DT>"INT16"<DD> 16-bit signed integer (short)
1402  <DT>"UINT16"<DD> 16-bit unsigned integer (unsigned short)
1403  <DT>"INT32"<DD> 32-bit signed integer (long)
1404  <DT>"UINT32"<DD> 32-bit unsigned integer (unsigned long)
1405  <DT>"INT64"<DD> 64-bit signed integer (long long)
1406  <DT>"UINT64"<DD> 64-bit unsigned integer (unsigned long long)
1407  <DT>"FLOAT"<DD> 32-bit floating point (float)
1408  <DT>"DOUBLE"<DD> 64-bit floating point (double)
1409  <DT>"UNKNOWN"<DD> any other type
1410  </DL>
1411  */
1412  std::string getDatasetType(std::string const & datasetName) const
1413  {
1414  HDF5Handle datasetHandle = getDatasetHandle(datasetName);
1415 
1416  hid_t datatype = H5Dget_type(datasetHandle);
1417  H5T_class_t dataclass = H5Tget_class(datatype);
1418  size_t datasize = H5Tget_size(datatype);
1419  H5T_sign_t datasign = H5Tget_sign(datatype);
1420 
1421  if(dataclass == H5T_FLOAT)
1422  {
1423  if(datasize == 4)
1424  return "FLOAT";
1425  else if(datasize == 8)
1426  return "DOUBLE";
1427  }
1428  else if(dataclass == H5T_INTEGER)
1429  {
1430  if(datasign == H5T_SGN_NONE)
1431  {
1432  if(datasize == 1)
1433  return "UINT8";
1434  else if(datasize == 2)
1435  return "UINT16";
1436  else if(datasize == 4)
1437  return "UINT32";
1438  else if(datasize == 8)
1439  return "UINT64";
1440  }
1441  else
1442  {
1443  if(datasize == 1)
1444  return "INT8";
1445  else if(datasize == 2)
1446  return "INT16";
1447  else if(datasize == 4)
1448  return "INT32";
1449  else if(datasize == 8)
1450  return "INT64";
1451  }
1452  }
1453  return "UNKNOWN";
1454  }
1455 
1456  /** \brief Obtain the HDF5 handle of a dataset.
1457  */
1458  HDF5Handle getDatasetHandle(std::string const & datasetName) const
1459  {
1460  std::string errorMessage = "HDF5File::getDatasetHandle(): Unable to open dataset '" + datasetName + "'.";
1461  return HDF5Handle(getDatasetHandle_(get_absolute_path(datasetName)), &H5Dclose, errorMessage.c_str());
1462  }
1463 
1464  /** \brief Obtain a shared HDF5 handle of a dataset.
1465  */
1466  HDF5HandleShared getDatasetHandleShared(std::string const & datasetName) const
1467  {
1468  std::string errorMessage = "HDF5File::getDatasetHandle(): Unable to open dataset '" + datasetName + "'.";
1469  return HDF5HandleShared(getDatasetHandle_(get_absolute_path(datasetName)), &H5Dclose, errorMessage.c_str());
1470  }
1471 
1472  /** \brief Obtain the HDF5 handle of a group (create the group if it doesn't exist).
1473  */
1474  HDF5Handle getGroupHandle(std::string group_name,
1475  std::string function_name = "HDF5File::getGroupHandle()")
1476  {
1477  std::string errorMessage = function_name + ": Group '" + group_name + "' not found.";
1478 
1479  // make group_name clean
1480  group_name = get_absolute_path(group_name);
1481 
1482  // group must exist
1483  vigra_precondition(group_name == "/" || H5Lexists(fileHandle_, group_name.c_str(), H5P_DEFAULT) != 0,
1484  errorMessage.c_str());
1485 
1486  // open group and return group handle
1487  return HDF5Handle(openCreateGroup_(group_name), &H5Gclose, "Internal error");
1488  }
1489 
1490  // helper function for the various listAttributes() variants.
1491  void ls_H5Aiterate(std::string const & group_or_dataset, ls_closure & data) const
1492  {
1493  H5O_type_t h5_type = get_object_type_(group_or_dataset);
1494  vigra_precondition(h5_type == H5O_TYPE_GROUP || h5_type == H5O_TYPE_DATASET,
1495  "HDF5File::listAttributes(): object \"" + group_or_dataset + "\" is neither a group nor a dataset.");
1496  // get object handle
1497  HDF5Handle object_handle(h5_type == H5O_TYPE_GROUP
1498  ? const_cast<HDF5File*>(this)->openCreateGroup_(group_or_dataset)
1499  : getDatasetHandle_(group_or_dataset),
1500  h5_type == H5O_TYPE_GROUP
1501  ? &H5Gclose
1502  : &H5Dclose,
1503  "HDF5File::listAttributes(): unable to open object.");
1504  hsize_t n = 0;
1505  H5Aiterate2(object_handle, H5_INDEX_NAME, H5_ITER_NATIVE, &n,
1506  HDF5_listAttributes_inserter_callback, static_cast<void*>(&data));
1507  }
1508 
1509  /** \brief List the attribute names of the given group or dataset.
1510 
1511  If \a group_or_dataset is empty or <tt>"."</tt> (a dot), the command
1512  refers to the current group of this file object.
1513  */
1514  inline std::vector<std::string> listAttributes(std::string const & group_or_dataset) const
1515  {
1516  std::vector<std::string> list;
1517  lsOpData data(list);
1518  ls_H5Aiterate(group_or_dataset, data);
1519  return list;
1520  }
1521 
1522  /** \brief Insert the attribute names of the given group or dataset into the given
1523  \a container by calling <tt>container.insert(std::string)</tt>.
1524 
1525  If \a group_or_dataset is empty or <tt>"."</tt> (a dot), the command
1526  refers to the current group of this file object.
1527  */
1528  template<class Container>
1529  void listAttributes(std::string const & group_or_dataset, Container & container) const
1530  {
1531  ls_container_data<Container> data(container);
1532  ls_H5Aiterate(group_or_dataset, data);
1533  }
1534 
1535  /** \brief Obtain the HDF5 handle of a attribute.
1536  */
1537  HDF5Handle getAttributeHandle(std::string dataset_name, std::string attribute_name) const
1538  {
1539  std::string message = "HDF5File::getAttributeHandle(): Attribute '" + attribute_name + "' not found.";
1540  return HDF5Handle(H5Aopen(getDatasetHandle(dataset_name), attribute_name.c_str(), H5P_DEFAULT),
1541  &H5Aclose, message.c_str());
1542  }
1543 
1544  /* Writing Attributes */
1545 
1546  /** \brief Write MultiArray Attributes.
1547  * In contrast to datasets, subarray access, chunks and compression are not available.
1548  */
1549  template<unsigned int N, class T, class Stride>
1550  inline void writeAttribute(std::string object_name,
1551  std::string attribute_name,
1552  const MultiArrayView<N, T, Stride> & array)
1553  {
1554  // make object_name clean
1555  object_name = get_absolute_path(object_name);
1556 
1557  write_attribute_(object_name, attribute_name, array, detail::getH5DataType<T>(), 1);
1558  }
1559 
1560  template<unsigned int N, class T, int SIZE, class Stride>
1561  inline void writeAttribute(std::string datasetName,
1562  std::string attributeName,
1563  const MultiArrayView<N, TinyVector<T, SIZE>, Stride> & array)
1564  {
1565  // make datasetName clean
1566  datasetName = get_absolute_path(datasetName);
1567 
1568  write_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), SIZE);
1569  }
1570 
1571  template<unsigned int N, class T, class Stride>
1572  inline void writeAttribute(std::string datasetName,
1573  std::string attributeName,
1574  const MultiArrayView<N, RGBValue<T>, Stride> & array)
1575  {
1576  // make datasetName clean
1577  datasetName = get_absolute_path(datasetName);
1578 
1579  write_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), 3);
1580  }
1581 
1582  /** \brief Write a single value.
1583  Specialization of the write function for simple datatypes
1584  */
1585  inline void writeAttribute(std::string object_name, std::string attribute_name, char data)
1586  { writeAtomicAttribute(object_name,attribute_name,data); }
1587  inline void writeAttribute(std::string datasetName, std::string attributeName, signed char data)
1588  { writeAtomicAttribute(datasetName,attributeName,data); }
1589  inline void writeAttribute(std::string datasetName, std::string attributeName, signed short data)
1590  { writeAtomicAttribute(datasetName,attributeName,data); }
1591  inline void writeAttribute(std::string datasetName, std::string attributeName, signed int data)
1592  { writeAtomicAttribute(datasetName,attributeName,data); }
1593  inline void writeAttribute(std::string datasetName, std::string attributeName, signed long data)
1594  { writeAtomicAttribute(datasetName,attributeName,data); }
1595  inline void writeAttribute(std::string datasetName, std::string attributeName, signed long long data)
1596  { writeAtomicAttribute(datasetName,attributeName,data); }
1597  inline void writeAttribute(std::string datasetName, std::string attributeName, unsigned char data)
1598  { writeAtomicAttribute(datasetName,attributeName,data); }
1599  inline void writeAttribute(std::string datasetName, std::string attributeName, unsigned short data)
1600  { writeAtomicAttribute(datasetName,attributeName,data); }
1601  inline void writeAttribute(std::string datasetName, std::string attributeName, unsigned int data)
1602  { writeAtomicAttribute(datasetName,attributeName,data); }
1603  inline void writeAttribute(std::string datasetName, std::string attributeName, unsigned long data)
1604  { writeAtomicAttribute(datasetName,attributeName,data); }
1605  inline void writeAttribute(std::string datasetName, std::string attributeName, unsigned long long data)
1606  { writeAtomicAttribute(datasetName,attributeName,data); }
1607  inline void writeAttribute(std::string datasetName, std::string attributeName, float data)
1608  { writeAtomicAttribute(datasetName,attributeName,data); }
1609  inline void writeAttribute(std::string datasetName, std::string attributeName, double data)
1610  { writeAtomicAttribute(datasetName,attributeName,data); }
1611  inline void writeAttribute(std::string datasetName, std::string attributeName, long double data)
1612  { writeAtomicAttribute(datasetName,attributeName,data); }
1613  inline void writeAttribute(std::string datasetName, std::string attributeName, const char* data)
1614  { writeAtomicAttribute(datasetName,attributeName,data); }
1615  inline void writeAttribute(std::string datasetName, std::string attributeName, std::string const & data)
1616  { writeAtomicAttribute(datasetName,attributeName,data.c_str()); }
1617 
1618  /** \brief Test if attribute exists.
1619  */
1620  bool existsAttribute(std::string object_name, std::string attribute_name)
1621  {
1622  std::string obj_path = get_absolute_path(object_name);
1623  htri_t exists = H5Aexists_by_name(fileHandle_, obj_path.c_str(),
1624  attribute_name.c_str(), H5P_DEFAULT);
1625  vigra_precondition(exists >= 0, "HDF5File::existsAttribute(): "
1626  "object '" + object_name + "' "
1627  "not found.");
1628  return exists != 0;
1629  }
1630 
1631  // Reading Attributes
1632 
1633  /** \brief Read MultiArray Attributes.
1634  * In contrast to datasets, subarray access is not available.
1635  */
1636  template<unsigned int N, class T, class Stride>
1637  inline void readAttribute(std::string object_name,
1638  std::string attribute_name,
1640  {
1641  // make object_name clean
1642  object_name = get_absolute_path(object_name);
1643 
1644  read_attribute_(object_name, attribute_name, array, detail::getH5DataType<T>(), 1);
1645  }
1646 
1647  template<unsigned int N, class T, int SIZE, class Stride>
1648  inline void readAttribute(std::string datasetName,
1649  std::string attributeName,
1650  MultiArrayView<N, TinyVector<T, SIZE>, Stride> array)
1651  {
1652  // make datasetName clean
1653  datasetName = get_absolute_path(datasetName);
1654 
1655  read_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), SIZE);
1656  }
1657 
1658  template<unsigned int N, class T, class Stride>
1659  inline void readAttribute(std::string datasetName,
1660  std::string attributeName,
1661  MultiArrayView<N, RGBValue<T>, Stride> array)
1662  {
1663  // make datasetName clean
1664  datasetName = get_absolute_path(datasetName);
1665 
1666  read_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), 3);
1667  }
1668 
1669  /** \brief Read a single value.
1670  Specialization of the read function for simple datatypes
1671  */
1672  inline void readAttribute(std::string object_name, std::string attribute_name, char &data)
1673  { readAtomicAttribute(object_name,attribute_name,data); }
1674  inline void readAttribute(std::string datasetName, std::string attributeName, signed char &data)
1675  { readAtomicAttribute(datasetName,attributeName,data); }
1676  inline void readAttribute(std::string datasetName, std::string attributeName, signed short &data)
1677  { readAtomicAttribute(datasetName,attributeName,data); }
1678  inline void readAttribute(std::string datasetName, std::string attributeName, signed int &data)
1679  { readAtomicAttribute(datasetName,attributeName,data); }
1680  inline void readAttribute(std::string datasetName, std::string attributeName, signed long &data)
1681  { readAtomicAttribute(datasetName,attributeName,data); }
1682  inline void readAttribute(std::string datasetName, std::string attributeName, signed long long &data)
1683  { readAtomicAttribute(datasetName,attributeName,data); }
1684  inline void readAttribute(std::string datasetName, std::string attributeName, unsigned char &data)
1685  { readAtomicAttribute(datasetName,attributeName,data); }
1686  inline void readAttribute(std::string datasetName, std::string attributeName, unsigned short &data)
1687  { readAtomicAttribute(datasetName,attributeName,data); }
1688  inline void readAttribute(std::string datasetName, std::string attributeName, unsigned int &data)
1689  { readAtomicAttribute(datasetName,attributeName,data); }
1690  inline void readAttribute(std::string datasetName, std::string attributeName, unsigned long &data)
1691  { readAtomicAttribute(datasetName,attributeName,data); }
1692  inline void readAttribute(std::string datasetName, std::string attributeName, unsigned long long &data)
1693  { readAtomicAttribute(datasetName,attributeName,data); }
1694  inline void readAttribute(std::string datasetName, std::string attributeName, float &data)
1695  { readAtomicAttribute(datasetName,attributeName,data); }
1696  inline void readAttribute(std::string datasetName, std::string attributeName, double &data)
1697  { readAtomicAttribute(datasetName,attributeName,data); }
1698  inline void readAttribute(std::string datasetName, std::string attributeName, long double &data)
1699  { readAtomicAttribute(datasetName,attributeName,data); }
1700  inline void readAttribute(std::string datasetName, std::string attributeName, std::string &data)
1701  { readAtomicAttribute(datasetName,attributeName,data); }
1702 
1703  // Writing data
1704 
1705  /** \brief Write multi arrays.
1706 
1707  Chunks can be activated by setting
1708  \code iChunkSize = size; //size > 0
1709  \endcode .
1710  The chunks will be hypercubes with edge length size. When <tt>iChunkSize == 0</tt>
1711  (default), the behavior depends on the <tt>compression</tt> setting: If no
1712  compression is requested, the data is written without chunking. Otherwise,
1713  chuning is required, and the chunk size is automatically selected such that
1714  each chunk contains about 300k pixels.
1715 
1716  Compression can be activated by setting
1717  \code compression = parameter; // 0 < parameter <= 9
1718  \endcode
1719  where 0 stands for no compression and 9 for maximum compression.
1720 
1721  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
1722  otherwise it will be interpreted as path relative to the current group.
1723 
1724  Note that the memory order between VIGRA and HDF5 files differs: VIGRA uses
1725  Fortran-order, while HDF5 uses C-order. This means that a VIGRA MultiArray,
1726  whose indices represent the 'x'-, 'y'-, and 'z'-axis in that order, is reversed
1727  upon writing to an HDF5 file, i.e. in the file the axis order is 'z', 'y', 'x'.
1728  */
1729  template<unsigned int N, class T, class Stride>
1730  inline void write(std::string datasetName,
1731  const MultiArrayView<N, T, Stride> & array,
1732  int iChunkSize = 0, int compression = 0)
1733  {
1734  // make datasetName clean
1735  datasetName = get_absolute_path(datasetName);
1736 
1737  typename MultiArrayShape<N>::type chunkSize;
1738  for(unsigned int i = 0; i < N; i++){
1739  chunkSize[i] = iChunkSize;
1740  }
1741  write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize, compression);
1742  }
1743 
1744  /** \brief Write multi arrays.
1745  Chunks can be activated by providing a MultiArrayShape as chunkSize.
1746  chunkSize must have equal dimension as array.
1747 
1748  Compression can be activated by setting
1749  \code compression = parameter; // 0 < parameter <= 9
1750  \endcode
1751  where 0 stands for no compression and 9 for maximum compression.
1752 
1753  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
1754  otherwise it will be interpreted as path relative to the current group.
1755 
1756  Note that the memory order between VIGRA and HDF5 files differs: VIGRA uses
1757  Fortran-order, while HDF5 uses C-order. This means that a VIGRA MultiArray,
1758  whose indices represent the 'x'-, 'y'-, and 'z'-axis in that order, is reversed
1759  upon writing to an HDF5 file, i.e. in the file the axis order is 'z', 'y', 'x'.
1760  */
1761  template<unsigned int N, class T, class Stride>
1762  inline void write(std::string datasetName,
1763  const MultiArrayView<N, T, Stride> & array,
1764  typename MultiArrayShape<N>::type chunkSize, int compression = 0)
1765  {
1766  // make datasetName clean
1767  datasetName = get_absolute_path(datasetName);
1768 
1769  write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize, compression);
1770  }
1771 
1772  /** \brief Write a multi array into a larger volume.
1773  blockOffset determines the position, where array is written.
1774 
1775  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
1776  otherwise it will be interpreted as path relative to the current group.
1777 
1778  Note that the memory order between VIGRA and HDF5 files differs: VIGRA uses
1779  Fortran-order, while HDF5 uses C-order. This means that a VIGRA MultiArray,
1780  whose indices represent the 'x'-, 'y'-, and 'z'-axis in that order, is reversed
1781  upon writing to an HDF5 file, i.e. in the file the axis order is 'z', 'y', 'x'.
1782  */
1783  template<unsigned int N, class T, class Stride>
1784  inline void writeBlock(std::string datasetName,
1785  typename MultiArrayShape<N>::type blockOffset,
1786  const MultiArrayView<N, T, Stride> & array)
1787  {
1788  // make datasetName clean
1789  datasetName = get_absolute_path(datasetName);
1790  typedef detail::HDF5TypeTraits<T> TypeTraits;
1791  writeBlock_(datasetName, blockOffset, array,
1792  TypeTraits::getH5DataType(), TypeTraits::numberOfBands());
1793  }
1794 
1795  template<unsigned int N, class T, class Stride>
1796  inline herr_t writeBlock(HDF5HandleShared dataset,
1797  typename MultiArrayShape<N>::type blockOffset,
1798  const MultiArrayView<N, T, Stride> & array)
1799  {
1800  typedef detail::HDF5TypeTraits<T> TypeTraits;
1801  return writeBlock_(dataset, blockOffset, array,
1802  TypeTraits::getH5DataType(), TypeTraits::numberOfBands());
1803  }
1804 
1805  // non-scalar (TinyVector) and unstrided multi arrays
1806  template<unsigned int N, class T, int SIZE, class Stride>
1807  inline void write(std::string datasetName,
1808  const MultiArrayView<N, TinyVector<T, SIZE>, Stride> & array,
1809  int iChunkSize = 0, int compression = 0)
1810  {
1811  // make datasetName clean
1812  datasetName = get_absolute_path(datasetName);
1813 
1814  typename MultiArrayShape<N>::type chunkSize;
1815  for(unsigned i = 0; i < N; i++){
1816  chunkSize[i] = iChunkSize;
1817  }
1818  write_(datasetName, array, detail::getH5DataType<T>(), SIZE, chunkSize, compression);
1819  }
1820 
1821  template<unsigned int N, class T, int SIZE, class Stride>
1822  inline void write(std::string datasetName,
1823  const MultiArrayView<N, TinyVector<T, SIZE>, Stride> & array,
1824  typename MultiArrayShape<N>::type chunkSize, int compression = 0)
1825  {
1826  // make datasetName clean
1827  datasetName = get_absolute_path(datasetName);
1828 
1829  write_(datasetName, array, detail::getH5DataType<T>(), SIZE, chunkSize, compression);
1830  }
1831 
1832  /** \brief Write array vectors.
1833 
1834  Compression can be activated by setting
1835  \code compression = parameter; // 0 < parameter <= 9
1836  \endcode
1837  where 0 stands for no compression and 9 for maximum compression.
1838 
1839  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
1840  otherwise it will be interpreted as path relative to the current group.
1841  */
1842  template<class T>
1843  void write(const std::string & datasetName,
1844  const ArrayVectorView<T> & array,
1845  int compression = 0)
1846  {
1847  // convert to a (trivial) MultiArrayView and forward.
1848  MultiArrayShape<1>::type shape(static_cast<MultiArrayIndex>(array.size()));
1849  const MultiArrayView<1, T> m_array(shape, const_cast<T*>(array.data()));
1850  write(datasetName, m_array, compression);
1851  }
1852 
1853  // non-scalar (RGBValue) and unstrided multi arrays
1854  template<unsigned int N, class T, class Stride>
1855  inline void write(std::string datasetName,
1856  const MultiArrayView<N, RGBValue<T>, Stride> & array,
1857  int iChunkSize = 0, int compression = 0)
1858  {
1859  // make datasetName clean
1860  datasetName = get_absolute_path(datasetName);
1861 
1862  typename MultiArrayShape<N>::type chunkSize;
1863  for(unsigned i = 0; i < N; i++){
1864  chunkSize[i] = iChunkSize;
1865  }
1866  write_(datasetName, array, detail::getH5DataType<T>(), 3, chunkSize, compression);
1867  }
1868 
1869  template<unsigned int N, class T, class Stride>
1870  inline void write(std::string datasetName,
1871  const MultiArrayView<N, RGBValue<T>, Stride> & array,
1872  typename MultiArrayShape<N>::type chunkSize, int compression = 0)
1873  {
1874  // make datasetName clean
1875  datasetName = get_absolute_path(datasetName);
1876 
1877  write_(datasetName, array, detail::getH5DataType<T>(), 3, chunkSize, compression);
1878  }
1879 
1880  /** \brief Write a single value.
1881  Specialization of the write function for simple datatypes
1882  */
1883  inline void write(std::string datasetName, char data) { writeAtomic(datasetName,data); }
1884  inline void write(std::string datasetName, signed char data) { writeAtomic(datasetName,data); }
1885  inline void write(std::string datasetName, signed short data) { writeAtomic(datasetName,data); }
1886  inline void write(std::string datasetName, signed int data) { writeAtomic(datasetName,data); }
1887  inline void write(std::string datasetName, signed long data) { writeAtomic(datasetName,data); }
1888  inline void write(std::string datasetName, signed long long data) { writeAtomic(datasetName,data); }
1889  inline void write(std::string datasetName, unsigned char data) { writeAtomic(datasetName,data); }
1890  inline void write(std::string datasetName, unsigned short data) { writeAtomic(datasetName,data); }
1891  inline void write(std::string datasetName, unsigned int data) { writeAtomic(datasetName,data); }
1892  inline void write(std::string datasetName, unsigned long data) { writeAtomic(datasetName,data); }
1893  inline void write(std::string datasetName, unsigned long long data) { writeAtomic(datasetName,data); }
1894  inline void write(std::string datasetName, float data) { writeAtomic(datasetName,data); }
1895  inline void write(std::string datasetName, double data) { writeAtomic(datasetName,data); }
1896  inline void write(std::string datasetName, long double data) { writeAtomic(datasetName,data); }
1897  inline void write(std::string datasetName, const char* data) { writeAtomic(datasetName,data); }
1898  inline void write(std::string datasetName, std::string const & data) { writeAtomic(datasetName,data.c_str()); }
1899 
1900  // Reading data
1901 
1902  /** \brief Read data into a multi array.
1903  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
1904  otherwise it will be interpreted as path relative to the current group.
1905 
1906  Note that the memory order between VIGRA and HDF5 files differs: VIGRA uses
1907  Fortran-order, while HDF5 uses C-order. This means that a HDF5 dataset,
1908  whose indices represent the 'z'-, 'y'-, and 'x'-axis in that order, is reversed
1909  upon reading into a MultiArrayView, i.e. in the array axis order must be 'x', 'y', 'z'.
1910  */
1911  template<unsigned int N, class T, class Stride>
1912  inline void read(std::string datasetName, MultiArrayView<N, T, Stride> array)
1913  {
1914  // make datasetName clean
1915  datasetName = get_absolute_path(datasetName);
1916 
1917  read_(datasetName, array, detail::getH5DataType<T>(), 1);
1918  }
1919 
1920  /** \brief Read data into a MultiArray. Resize MultiArray to the correct size.
1921  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
1922  otherwise it will be interpreted as path relative to the current group.
1923 
1924  Note that the memory order between VIGRA and HDF5 files differs: VIGRA uses
1925  Fortran-order, while HDF5 uses C-order. This means that a HDF5 dataset,
1926  whose indices represent the 'z'-, 'y'-, and 'x'-axis in that order, is reversed
1927  upon reading into a MultiArray, i.e. in the array axis order will be 'x', 'y', 'z'.
1928  */
1929  template<unsigned int N, class T, class Alloc>
1930  inline void readAndResize(std::string datasetName, MultiArray<N, T, Alloc> & array)
1931  {
1932  // make datasetName clean
1933  datasetName = get_absolute_path(datasetName);
1934 
1935  // get dataset dimension
1936  ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName);
1937 
1938  // check if dimensions are correct
1939  vigra_precondition(N == MultiArrayIndex(dimshape.size()), // the object in the HDF5 file may have one additional dimension which we then interpret as the pixel type bands
1940  "HDF5File::readAndResize(): Array dimension disagrees with dataset dimension.");
1941 
1942  // reshape target MultiArray
1943  typename MultiArrayShape<N>::type shape;
1944  for(int k=0; k < static_cast<int>(dimshape.size()); ++k)
1945  shape[k] = static_cast<MultiArrayIndex>(dimshape[k]);
1946  array.reshape(shape);
1947 
1948  read_(datasetName, array, detail::getH5DataType<T>(), 1);
1949  }
1950 
1951  /** \brief Read data into an array vector.
1952  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
1953  otherwise it will be interpreted as path relative to the current group.
1954  */
1955  template<class T>
1956  inline void read(const std::string & datasetName, ArrayVectorView<T> array)
1957  {
1958  // convert to a (trivial) MultiArrayView and forward.
1959  MultiArrayShape<1>::type shape(array.size());
1960  MultiArrayView<1, T> m_array(shape, (array.data()));
1961  read(datasetName, m_array);
1962  }
1963 
1964  /** \brief Read data into an array vector. Resize the array vector to the correct size.
1965  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
1966  otherwise it will be interpreted as path relative to the current group.
1967  */
1968  template<class T>
1969  inline void readAndResize(std::string datasetName,
1970  ArrayVector<T> & array)
1971  {
1972  // make dataset name clean
1973  datasetName = get_absolute_path(datasetName);
1974 
1975  // get dataset dimension
1976  ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName);
1977 
1978  // check if dimensions are correct
1979  vigra_precondition(1 == MultiArrayIndex(dimshape.size()),
1980  "HDF5File::readAndResize(): Array dimension disagrees with Dataset dimension must equal one for vigra::ArrayVector.");
1981 
1982  // resize target array vector
1983  array.resize((typename ArrayVector<T>::size_type)dimshape[0]);
1984  // convert to a (trivial) MultiArrayView and forward.
1985  MultiArrayShape<1>::type shape(static_cast<MultiArrayIndex>(array.size()));
1986  MultiArrayView<1, T> m_array(shape, (array.data()));
1987 
1988  read_(datasetName, m_array, detail::getH5DataType<T>(), 1);
1989  }
1990 
1991  /** \brief Read a block of data into a multi array.
1992  This function allows to read a small block out of a larger volume stored
1993  in an HDF5 dataset.
1994 
1995  blockOffset determines the position of the block.
1996  blockSize determines the size in each dimension of the block.
1997 
1998  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
1999  otherwise it will be interpreted as path relative to the current group.
2000 
2001  Note that the memory order between VIGRA and HDF5 files differs: VIGRA uses
2002  Fortran-order, while HDF5 uses C-order. This means that a HDF5 dataset,
2003  whose indices represent the 'z'-, 'y'-, and 'x'-axis in that order, is reversed
2004  upon reading into a MultiArray, i.e. in the array axis order will be 'x', 'y', 'z'.
2005  */
2006  template<unsigned int N, class T, class Stride>
2007  inline void readBlock(std::string datasetName,
2008  typename MultiArrayShape<N>::type blockOffset,
2009  typename MultiArrayShape<N>::type blockShape,
2011  {
2012  // make datasetName clean
2013  datasetName = get_absolute_path(datasetName);
2014  typedef detail::HDF5TypeTraits<T> TypeTraits;
2015  readBlock_(datasetName, blockOffset, blockShape, array,
2016  TypeTraits::getH5DataType(), TypeTraits::numberOfBands());
2017  }
2018 
2019  template<unsigned int N, class T, class Stride>
2020  inline herr_t readBlock(HDF5HandleShared dataset,
2021  typename MultiArrayShape<N>::type blockOffset,
2022  typename MultiArrayShape<N>::type blockShape,
2024  {
2025  typedef detail::HDF5TypeTraits<T> TypeTraits;
2026  return readBlock_(dataset, blockOffset, blockShape, array,
2027  TypeTraits::getH5DataType(), TypeTraits::numberOfBands());
2028  }
2029 
2030  // non-scalar (TinyVector) and unstrided target MultiArrayView
2031  template<unsigned int N, class T, int SIZE, class Stride>
2032  inline void read(std::string datasetName, MultiArrayView<N, TinyVector<T, SIZE>, Stride> array)
2033  {
2034  // make datasetName clean
2035  datasetName = get_absolute_path(datasetName);
2036 
2037  read_(datasetName, array, detail::getH5DataType<T>(), SIZE);
2038  }
2039 
2040  // non-scalar (TinyVector) MultiArray
2041  template<unsigned int N, class T, int SIZE, class Alloc>
2042  inline void readAndResize(std::string datasetName, MultiArray<N, TinyVector<T, SIZE>, Alloc> & array)
2043  {
2044  // make datasetName clean
2045  datasetName = get_absolute_path(datasetName);
2046 
2047  // get dataset dimension
2048  ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName);
2049 
2050  // check if dimensions are correct
2051  vigra_precondition((N+1) == MultiArrayIndex(dimshape.size()) &&
2052  SIZE == dimshape[0], // the object in the HDF5 file must have one additional dimension which we interpret as the pixel type bands
2053  "HDF5File::readAndResize(): Array dimension disagrees with dataset dimension.");
2054 
2055  // reshape target MultiArray
2056  typename MultiArrayShape<N>::type shape;
2057  for(int k=1; k < static_cast<int>(dimshape.size()); ++k)
2058  shape[k-1] = static_cast<MultiArrayIndex>(dimshape[k]);
2059  array.reshape(shape);
2060 
2061  read_(datasetName, array, detail::getH5DataType<T>(), SIZE);
2062  }
2063 
2064  // non-scalar (RGBValue) and unstrided target MultiArrayView
2065  template<unsigned int N, class T, class Stride>
2066  inline void read(std::string datasetName, MultiArrayView<N, RGBValue<T>, Stride> array)
2067  {
2068  // make datasetName clean
2069  datasetName = get_absolute_path(datasetName);
2070 
2071  read_(datasetName, array, detail::getH5DataType<T>(), 3);
2072  }
2073 
2074  // non-scalar (RGBValue) MultiArray
2075  template<unsigned int N, class T, class Alloc>
2076  inline void readAndResize(std::string datasetName, MultiArray<N, RGBValue<T>, Alloc> & array)
2077  {
2078  // make datasetName clean
2079  datasetName = get_absolute_path(datasetName);
2080 
2081  // get dataset dimension
2082  ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName);
2083 
2084  // check if dimensions are correct
2085  vigra_precondition((N+1) == MultiArrayIndex(dimshape.size()) &&
2086  3 == dimshape[0], // the object in the HDF5 file must have one additional dimension which we interpret as the pixel type bands
2087  "HDF5File::readAndResize(): Array dimension disagrees with dataset dimension.");
2088 
2089  // reshape target MultiArray
2090  typename MultiArrayShape<N>::type shape;
2091  for(int k=1; k < static_cast<int>(dimshape.size()); ++k)
2092  shape[k-1] = static_cast<MultiArrayIndex>(dimshape[k]);
2093  array.reshape(shape);
2094 
2095  read_(datasetName, array, detail::getH5DataType<T>(), 3);
2096  }
2097 
2098  /** \brief Read a single value.
2099  Specialization of the read function for simple datatypes
2100  */
2101  inline void read(std::string datasetName, char &data) { readAtomic(datasetName,data); }
2102  inline void read(std::string datasetName, signed char &data) { readAtomic(datasetName,data); }
2103  inline void read(std::string datasetName, signed short &data) { readAtomic(datasetName,data); }
2104  inline void read(std::string datasetName, signed int &data) { readAtomic(datasetName,data); }
2105  inline void read(std::string datasetName, signed long &data) { readAtomic(datasetName,data); }
2106  inline void read(std::string datasetName, signed long long &data) { readAtomic(datasetName,data); }
2107  inline void read(std::string datasetName, unsigned char &data) { readAtomic(datasetName,data); }
2108  inline void read(std::string datasetName, unsigned short &data) { readAtomic(datasetName,data); }
2109  inline void read(std::string datasetName, unsigned int &data) { readAtomic(datasetName,data); }
2110  inline void read(std::string datasetName, unsigned long &data) { readAtomic(datasetName,data); }
2111  inline void read(std::string datasetName, unsigned long long &data) { readAtomic(datasetName,data); }
2112  inline void read(std::string datasetName, float &data) { readAtomic(datasetName,data); }
2113  inline void read(std::string datasetName, double &data) { readAtomic(datasetName,data); }
2114  inline void read(std::string datasetName, long double &data) { readAtomic(datasetName,data); }
2115  inline void read(std::string datasetName, std::string &data) { readAtomic(datasetName,data); }
2116 
2117  /** \brief Create a new dataset.
2118  This function can be used to create a dataset filled with a default value \a init,
2119  for example before writing data into it using \ref writeBlock().
2120 
2121  shape determines the dimension and the size of the dataset.
2122 
2123  Chunks can be activated by providing a MultiArrayShape as chunkSize.
2124  chunkSize must have equal dimension as array.
2125 
2126  Compression can be activated by setting
2127  \code compression = parameter; // 0 < parameter <= 9
2128  \endcode
2129  where 0 stands for no compression and 9 for maximum compression. If
2130  a non-zero compression level is specified, but the chunk size is zero,
2131  a default chunk size will be chosen (compression always requires chunks).
2132 
2133  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
2134  otherwise it will be interpreted as path relative to the current group.
2135 
2136  Note that the memory order between VIGRA and HDF5 files differs: VIGRA uses
2137  Fortran-order, while HDF5 uses C-order. This means that a VIGRA MultiArray,
2138  whose indices represent the 'x'-, 'y'-, and 'z'-axis in that order, is reversed
2139  upon writing to an HDF5 file, i.e. in the file the axis order is 'z', 'y', 'x'.
2140  */
2141  template<int N, class T>
2142  HDF5HandleShared
2143  createDataset(std::string datasetName,
2144  TinyVector<MultiArrayIndex, N> const & shape,
2145  typename detail::HDF5TypeTraits<T>::value_type init =
2146  typename detail::HDF5TypeTraits<T>::value_type(),
2147 #ifdef _MSC_VER
2148  TinyVector<MultiArrayIndex, N> const & chunkSize = TinyVector<MultiArrayIndex, N>(),
2149 #else
2150  TinyVector<MultiArrayIndex, N> const & chunkSize = (TinyVector<MultiArrayIndex, N>()),
2151 #endif
2152  int compressionParameter = 0);
2153 
2154  // for backwards compatibility
2155  template<int N, class T>
2156  HDF5HandleShared
2157  createDataset(std::string datasetName,
2158  TinyVector<MultiArrayIndex, N> const & shape,
2159  T init,
2160  int iChunkSize,
2161  int compressionParameter = 0)
2162  {
2163  typename MultiArrayShape<N>::type chunkSize;
2164  for(int i = 0; i < N; i++){
2165  chunkSize[i] = iChunkSize;
2166  }
2167  return this->template createDataset<N, T>(datasetName, shape, init,
2168  chunkSize, compressionParameter);
2169  }
2170 
2171  /** \brief Immediately write all data to disk
2172  */
2173  inline void flushToDisk()
2174  {
2175  if(fileHandle_)
2176  H5Fflush(fileHandle_, H5F_SCOPE_GLOBAL);
2177  }
2178 
2179  private:
2180 
2181  /* Simple extension of std::string for splitting into two parts
2182  *
2183  * Strings (in particular: file/dataset paths) will be split into two
2184  * parts. The split is made at the last occurrence of the delimiter.
2185  *
2186  * For example, "/path/to/some/file" will be split (delimiter = "/") into
2187  * first() = "/path/to/some" and last() = "file".
2188  */
2189  class SplitString: public std::string {
2190  public:
2191  SplitString(std::string &sstring): std::string(sstring) {};
2192 
2193  // return the part of the string before the delimiter
2194  std::string first(char delimiter = '/')
2195  {
2196  size_t lastPos = find_last_of(delimiter);
2197  if(lastPos == std::string::npos) // delimiter not found --> no first
2198  return "";
2199 
2200  return std::string(begin(), begin()+lastPos+1);
2201  }
2202 
2203  // return the part of the string after the delimiter
2204  std::string last(char delimiter = '/')
2205  {
2206  size_t lastPos = find_last_of(delimiter);
2207  if(lastPos == std::string::npos) // delimiter not found --> only last
2208  return std::string(*this);
2209  return std::string(begin()+lastPos+1, end());
2210  }
2211  };
2212 
2213  template <class Shape>
2214  ArrayVector<hsize_t>
2215  defineChunks(Shape chunks, Shape const & shape, int numBands, int compression = 0)
2216  {
2217  if(prod(chunks) > 0)
2218  {
2219  ArrayVector<hsize_t> res(chunks.begin(), chunks.end());
2220  if(numBands > 1)
2221  res.insert(res.begin(), static_cast<hsize_t>(numBands));
2222  return res;
2223  }
2224  else if(compression > 0)
2225  {
2226  // set default chunks to enable compression
2227  chunks = min(detail::ChunkShape<Shape::static_size>::defaultShape(), shape);
2228  ArrayVector<hsize_t> res(chunks.begin(), chunks.end());
2229  if(numBands > 1)
2230  res.insert(res.begin(), static_cast<hsize_t>(numBands));
2231  return res;
2232  }
2233  else
2234  {
2235  return ArrayVector<hsize_t>();
2236  }
2237  }
2238 
2239  public:
2240 
2241  /** \brief takes any path and converts it into an absolute path
2242  in the current file.
2243 
2244  Elements like "." and ".." are treated as expected.
2245  Links are not supported or resolved.
2246  */
2247  inline std::string get_absolute_path(std::string path) const {
2248  // check for empty input or "." and return the current folder
2249  if(path.length() == 0 || path == "."){
2250  return currentGroupName_();
2251  }
2252 
2253  std::string str;
2254  // convert to absolute path
2255  if(relativePath_(path)){
2256  std::string cname = currentGroupName_();
2257  if (cname == "/")
2258  str = currentGroupName_()+path;
2259  else
2260  str = currentGroupName_()+"/"+path;
2261  }else{
2262  str = path;
2263  }
2264 
2265  // cut out "./"
2266  std::string::size_type startpos = 0;
2267  while(str.find(std::string("./"), startpos) != std::string::npos){
2268  std::string::size_type pos = str.find(std::string("./"), startpos);
2269  startpos = pos+1;
2270  // only cut if "./" is not part of "../" (see below)
2271  if(str.substr(pos-1,3) != "../"){
2272  // cut out part of the string
2273  str = str.substr(0,pos) + str.substr(pos+2,str.length()-pos-2);
2274  startpos = pos;
2275  }
2276  }
2277 
2278  // cut out pairs of "bla/../"
2279  while(str.find(std::string("..")) != std::string::npos){
2280  std::string::size_type pos = str.find(std::string(".."));
2281 
2282  // find first slash after ".."
2283  std::string::size_type end = str.find("/",pos);
2284  if(end != std::string::npos){
2285  // also include slash
2286  end++;
2287  }else{
2288  // no "/" after ".." --> this is a group, add a "/"
2289  str = str + "/";
2290  end = str.length();
2291  }
2292 
2293  // find first slash before ".."
2294  std::string::size_type prev_slash = str.rfind("/",pos);
2295  // if the root slash is the first before ".." --> Error
2296  vigra_invariant(prev_slash != 0 && prev_slash != std::string::npos,
2297  "Error parsing path: "+str);
2298  // find second slash before ".."
2299  std::string::size_type begin = str.rfind("/",prev_slash-1);
2300 
2301  // cut out part of the string
2302  str = str.substr(0,begin+1) + str.substr(end,str.length()-end);
2303  }
2304 
2305  return str;
2306  }
2307 
2308  protected:
2309 
2310  /* checks if the given path is a relative path.
2311  */
2312  inline bool relativePath_(std::string & path) const
2313  {
2314  std::string::size_type pos = path.find('/') ;
2315  if(pos == 0)
2316  return false;
2317 
2318  return true;
2319  }
2320 
2321  /* return the name of the current group
2322  */
2323  inline std::string currentGroupName_() const
2324  {
2325  int len = H5Iget_name(cGroupHandle_,NULL,1000);
2326  ArrayVector<char> name (len+1,0);
2327  H5Iget_name(cGroupHandle_,name.begin(),len+1);
2328 
2329  return std::string(name.begin());
2330  }
2331 
2332  /* return the name of the current file
2333  */
2334  inline std::string fileName_() const
2335  {
2336  int len = H5Fget_name(fileHandle_,NULL,1000);
2337  ArrayVector<char> name (len+1,0);
2338  H5Fget_name(fileHandle_,name.begin(),len+1);
2339 
2340  return std::string(name.begin());
2341  }
2342 
2343  /* create an empty file or open an existing one
2344  */
2345  inline hid_t createFile_(std::string filePath, OpenMode mode = Open)
2346  {
2347  // try to open file
2348  FILE * pFile;
2349  pFile = fopen ( filePath.c_str(), "r" );
2350  hid_t fileId;
2351 
2352  // check if opening was successful (= file exists)
2353  if ( pFile == NULL )
2354  {
2355  vigra_precondition(mode != OpenReadOnly,
2356  "HDF5File::open(): cannot open non-existing file in read-only mode.");
2357  fileId = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
2358  }
2359  else
2360  {
2361  fclose(pFile);
2362  if(mode == OpenReadOnly)
2363  {
2364  fileId = H5Fopen(filePath.c_str(), H5F_ACC_RDONLY, H5P_DEFAULT);
2365  }
2366  else if(mode == New)
2367  {
2368  std::remove(filePath.c_str());
2369  fileId = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
2370  }
2371  else
2372  {
2373  fileId = H5Fopen(filePath.c_str(), H5F_ACC_RDWR, H5P_DEFAULT);
2374  }
2375  }
2376  return fileId;
2377  }
2378 
2379  /* \brief Open a group.
2380 
2381  A negative value is returned when the group does not exist or when opening
2382  fails for other reasons.
2383  */
2384  hid_t openGroup_(std::string groupName) const
2385  {
2386  return const_cast<HDF5File *>(this)->openCreateGroup_(groupName, false);
2387  }
2388 
2389  /* \brief Open or create a group.
2390 
2391  If \a create is <tt>true</tt> and the group does not exist, it will be created,
2392  including all necessary parent groups. If group creation fails, a negative
2393  value is returned. Likewise, a negative value is returned when \a create
2394  is <tt>false</tt> and the group does not exist or when opening of the group
2395  fails for other reasons.
2396  */
2397  hid_t openCreateGroup_(std::string groupName, bool create = true)
2398  {
2399  // make groupName clean
2400  groupName = get_absolute_path(groupName);
2401 
2402  // open root group
2403  hid_t parent = H5Gopen(fileHandle_, "/", H5P_DEFAULT);
2404  if(groupName == "/")
2405  {
2406  return parent;
2407  }
2408 
2409  // remove leading /
2410  groupName = std::string(groupName.begin()+1, groupName.end());
2411 
2412  // check if the groupName has finishing slash
2413  if( groupName.size() != 0 && *groupName.rbegin() != '/')
2414  {
2415  groupName = groupName + '/';
2416  }
2417 
2418  // We determine if the group exists by checking the return value of H5Gopen.
2419  // To do so, we must temporarily disable error reporting.
2420  // Alternatively, we could use H5LTfind_dataset(), but this is much slower.
2421 
2422  // Save current error handling.
2423  H5E_auto2_t error_handler_ori;
2424  void *error_data_ori;
2425  H5Eget_auto(H5E_DEFAULT, &error_handler_ori, &error_data_ori);
2426 
2427  // Turn off error handling and register function to restore it upon return.
2428  H5Eset_auto(H5E_DEFAULT, NULL, NULL);
2429  VIGRA_FINALLY(
2430  H5Eset_auto(H5E_DEFAULT, error_handler_ori, error_data_ori));
2431 
2432  // Open or create subgroups one by one
2433  std::string::size_type begin = 0, end = groupName.find('/');
2434  while (end != std::string::npos)
2435  {
2436  std::string group(groupName.begin()+begin, groupName.begin()+end);
2437 
2438  hid_t prevParent = parent;
2439  parent = H5Gopen(prevParent, group.c_str(), H5P_DEFAULT);
2440  if(parent < 0 && create) // group doesn't exist, but we are supposed to create it
2441  parent = H5Gcreate(prevParent, group.c_str(), H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
2442  H5Gclose(prevParent);
2443 
2444  if(parent < 0)
2445  break;
2446 
2447  begin = end + 1;
2448  end = groupName.find('/', begin);
2449  }
2450 
2451  return parent;
2452  }
2453 
2454  /* delete a dataset by unlinking it from the file structure. This does not
2455  delete the data!
2456  */
2457  inline void deleteDataset_(hid_t parent, std::string datasetName)
2458  {
2459  // delete existing data and create new dataset
2460  if(H5LTfind_dataset(parent, datasetName.c_str()))
2461  {
2462 
2463  #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6)
2464  if(H5Gunlink(parent, datasetName.c_str()) < 0)
2465  {
2466  vigra_postcondition(false, "HDF5File::deleteDataset_(): Unable to delete existing data.");
2467  }
2468  #else
2469  if(H5Ldelete(parent, datasetName.c_str(), H5P_DEFAULT ) < 0)
2470  {
2471  vigra_postcondition(false, "HDF5File::deleteDataset_(): Unable to delete existing data.");
2472  }
2473  #endif
2474  }
2475  }
2476 
2477  /* get the handle of a dataset specified by a string
2478  */
2479  hid_t getDatasetHandle_(std::string datasetName) const
2480  {
2481  // make datasetName clean
2482  datasetName = get_absolute_path(datasetName);
2483 
2484  std::string groupname = SplitString(datasetName).first();
2485  std::string setname = SplitString(datasetName).last();
2486 
2487  if(H5Lexists(fileHandle_, datasetName.c_str(), H5P_DEFAULT) <= 0)
2488  {
2489  std::cerr << "HDF5File::getDatasetHandle_(): Dataset '" << datasetName << "' does not exist.\n";
2490  return -1;
2491  }
2492 
2493  // Open parent group
2494  HDF5Handle groupHandle(openGroup_(groupname), &H5Gclose, "HDF5File::getDatasetHandle_(): Internal error");
2495 
2496  return H5Dopen(groupHandle, setname.c_str(), H5P_DEFAULT);
2497  }
2498 
2499  /* get the type of an object specified by a string
2500  */
2501  H5O_type_t get_object_type_(std::string name) const
2502  {
2503  name = get_absolute_path(name);
2504  std::string group_name = SplitString(name).first();
2505  std::string object_name = SplitString(name).last();
2506  if (!object_name.size())
2507  return H5O_TYPE_GROUP;
2508 
2509  htri_t exists = H5Lexists(fileHandle_, name.c_str(), H5P_DEFAULT);
2510  vigra_precondition(exists > 0, "HDF5File::get_object_type_(): "
2511  "object \"" + name + "\" "
2512  "not found.");
2513  // open parent group
2514  HDF5Handle group_handle(openGroup_(group_name), &H5Gclose, "Internal error");
2515  return HDF5_get_type(group_handle, name.c_str());
2516  }
2517 
2518  /* low-level write function to write vigra MultiArray data as an attribute
2519  */
2520  template<unsigned int N, class T, class Stride>
2521  void write_attribute_(std::string name,
2522  const std::string & attribute_name,
2523  const MultiArrayView<N, T, Stride> & array,
2524  const hid_t datatype,
2525  const int numBandsOfType);
2526 
2527  /* Write single value attribute
2528  This function allows to write data of atomic datatypes (int, long, double)
2529  as an attribute in the HDF5 file. So it is not necessary to create a MultiArray
2530  of size 1 to write a single number.
2531  */
2532  template<class T>
2533  inline void writeAtomicAttribute(std::string datasetName, std::string attributeName, const T data)
2534  {
2535  // make datasetName clean
2536  datasetName = get_absolute_path(datasetName);
2537 
2538  typename MultiArrayShape<1>::type chunkSize;
2539  chunkSize[0] = 0;
2540  MultiArray<1,T> array(MultiArrayShape<1>::type(1));
2541  array[0] = data;
2542  write_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), 1);
2543  }
2544 
2545  /* low-level read function to read vigra MultiArray data from attributes
2546  */
2547  template<unsigned int N, class T, class Stride>
2548  void read_attribute_(std::string datasetName,
2549  std::string attributeName,
2550  MultiArrayView<N, T, Stride> array,
2551  const hid_t datatype, const int numBandsOfType);
2552 
2553  /* Read a single value attribute.
2554  This functions allows to read a single value attribute of atomic datatype (int, long, double)
2555  from the HDF5 file. So it is not necessary to create a MultiArray
2556  of size 1 to read a single number.
2557  */
2558  template<class T>
2559  inline void readAtomicAttribute(std::string datasetName, std::string attributeName, T & data)
2560  {
2561  // make datasetName clean
2562  datasetName = get_absolute_path(datasetName);
2563 
2564  MultiArray<1,T> array(MultiArrayShape<1>::type(1));
2565  read_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), 1);
2566  data = array[0];
2567  }
2568 
2569  inline void readAtomicAttribute(std::string datasetName, std::string attributeName, std::string & data)
2570  {
2571  // make datasetName clean
2572  datasetName = get_absolute_path(datasetName);
2573 
2574  MultiArray<1,const char *> array(MultiArrayShape<1>::type(1));
2575  read_attribute_(datasetName, attributeName, array, detail::getH5DataType<const char *>(), 1);
2576  data = std::string(array[0]);
2577  }
2578 
2579  /* low-level write function to write vigra unstrided MultiArray data
2580  */
2581  template<unsigned int N, class T, class Stride>
2582  void write_(std::string &datasetName,
2583  const MultiArrayView<N, T, Stride> & array,
2584  const hid_t datatype,
2585  const int numBandsOfType,
2586  typename MultiArrayShape<N>::type &chunkSize,
2587  int compressionParameter = 0);
2588 
2589  /* Write single value as dataset.
2590  This functions allows to write data of atomic datatypes (int, long, double)
2591  as a dataset in the HDF5 file. So it is not necessary to create a MultiArray
2592  of size 1 to write a single number.
2593 
2594  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
2595  otherwise it will be interpreted as path relative to the current group.
2596  */
2597  template<class T>
2598  inline void writeAtomic(std::string datasetName, const T data)
2599  {
2600  // make datasetName clean
2601  datasetName = get_absolute_path(datasetName);
2602 
2603  typename MultiArrayShape<1>::type chunkSize;
2604  chunkSize[0] = 0;
2605  MultiArray<1,T> array(MultiArrayShape<1>::type(1));
2606  array[0] = data;
2607  write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize,0);
2608  }
2609 
2610  /* low-level read function to read vigra unstrided MultiArray data
2611  */
2612  template<unsigned int N, class T, class Stride>
2613  void read_(std::string datasetName,
2614  MultiArrayView<N, T, Stride> array,
2615  const hid_t datatype, const int numBandsOfType);
2616 
2617  /* Read a single value.
2618  This functions allows to read a single datum of atomic datatype (int, long, double)
2619  from the HDF5 file. So it is not necessary to create a MultiArray
2620  of size 1 to read a single number.
2621 
2622  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
2623  otherwise it will be interpreted as path relative to the current group.
2624  */
2625  template<class T>
2626  inline void readAtomic(std::string datasetName, T & data)
2627  {
2628  // make datasetName clean
2629  datasetName = get_absolute_path(datasetName);
2630 
2631  MultiArray<1,T> array(MultiArrayShape<1>::type(1));
2632  read_(datasetName, array, detail::getH5DataType<T>(), 1);
2633  data = array[0];
2634  }
2635 
2636  inline void readAtomic(std::string datasetName, std::string & data)
2637  {
2638  // make datasetName clean
2639  datasetName = get_absolute_path(datasetName);
2640 
2641  MultiArray<1,const char *> array(MultiArrayShape<1>::type(1));
2642  read_(datasetName, array, detail::getH5DataType<const char *>(), 1);
2643  data = std::string(array[0]);
2644  }
2645 
2646  /* low-level write function to write vigra unstrided MultiArray data into a sub-block of a dataset
2647  */
2648  template<unsigned int N, class T, class Stride>
2649  void writeBlock_(std::string datasetName,
2650  typename MultiArrayShape<N>::type &blockOffset,
2651  const MultiArrayView<N, T, Stride> & array,
2652  const hid_t datatype,
2653  const int numBandsOfType)
2654  {
2655  // open dataset if it exists
2656  std::string errorMessage = "HDF5File::writeBlock(): Error opening dataset '" + datasetName + "'.";
2657  HDF5HandleShared dataset(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
2658  herr_t status = writeBlock_(dataset, blockOffset, array, datatype, numBandsOfType);
2659  vigra_postcondition(status >= 0,
2660  "HDF5File::writeBlock(): write to dataset '" + datasetName + "' via H5Dwrite() failed.");
2661  }
2662 
2663  /* low-level write function to write vigra unstrided MultiArray data into a
2664  sub-block of a dataset. Returns the result of the internal call
2665  to <tt>H5Dwrite()</tt>.
2666  */
2667  template<unsigned int N, class T, class Stride>
2668  herr_t writeBlock_(HDF5HandleShared dataset,
2669  typename MultiArrayShape<N>::type &blockOffset,
2670  const MultiArrayView<N, T, Stride> & array,
2671  const hid_t datatype,
2672  const int numBandsOfType);
2673 
2674  /* low-level read function to read vigra unstrided MultiArray data from a sub-block of a dataset.
2675 
2676  The array must have the same shape as the block.
2677  */
2678  template<unsigned int N, class T, class Stride>
2679  void readBlock_(std::string datasetName,
2680  typename MultiArrayShape<N>::type &blockOffset,
2681  typename MultiArrayShape<N>::type &blockShape,
2682  MultiArrayView<N, T, Stride> array,
2683  const hid_t datatype, const int numBandsOfType)
2684  {
2685  std::string errorMessage ("HDF5File::readBlock(): Unable to open dataset '" + datasetName + "'.");
2686  HDF5HandleShared dataset(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
2687  herr_t status = readBlock_(dataset, blockOffset, blockShape, array, datatype, numBandsOfType);
2688  vigra_postcondition(status >= 0,
2689  "HDF5File::readBlock(): read from dataset '" + datasetName + "' via H5Dread() failed.");
2690  }
2691 
2692  /* low-level read function to read vigra unstrided MultiArray data from a sub-block of a dataset.
2693 
2694  The array must have the same shape as the block. Returns the result of the internal call
2695  to <tt>H5Dread()</tt>.
2696  */
2697  template<unsigned int N, class T, class Stride>
2698  herr_t readBlock_(HDF5HandleShared dataset,
2699  typename MultiArrayShape<N>::type &blockOffset,
2700  typename MultiArrayShape<N>::type &blockShape,
2701  MultiArrayView<N, T, Stride> array,
2702  const hid_t datatype, const int numBandsOfType);
2703 }; /* class HDF5File */
2704 
2705 /********************************************************************/
2706 
2707 template<int N, class T>
2708 HDF5HandleShared
2709 HDF5File::createDataset(std::string datasetName,
2710  TinyVector<MultiArrayIndex, N> const & shape,
2711  typename detail::HDF5TypeTraits<T>::value_type init,
2712  TinyVector<MultiArrayIndex, N> const & chunkSize,
2713  int compressionParameter)
2714 {
2715  vigra_precondition(!isReadOnly(),
2716  "HDF5File::createDataset(): file is read-only.");
2717 
2718  // make datasetName clean
2719  datasetName = get_absolute_path(datasetName);
2720 
2721  std::string groupname = SplitString(datasetName).first();
2722  std::string setname = SplitString(datasetName).last();
2723 
2724  hid_t parent = openCreateGroup_(groupname);
2725 
2726  // delete the dataset if it already exists
2727  deleteDataset_(parent, setname);
2728 
2729  // invert dimensions to guarantee c-order
2730  // add an extra dimension in case that the data is non-scalar
2731  typedef detail::HDF5TypeTraits<T> TypeTraits;
2732  ArrayVector<hsize_t> shape_inv;
2733  if(TypeTraits::numberOfBands() > 1)
2734  {
2735  shape_inv.resize(N+1);
2736  shape_inv[N] = TypeTraits::numberOfBands();
2737  }
2738  else
2739  {
2740  shape_inv.resize(N);
2741  }
2742  for(int k=0; k<N; ++k)
2743  shape_inv[N-1-k] = shape[k];
2744 
2745  // create dataspace
2746  HDF5Handle
2747  dataspaceHandle = HDF5Handle(H5Screate_simple(shape_inv.size(), shape_inv.data(), NULL),
2748  &H5Sclose, "HDF5File::createDataset(): unable to create dataspace for scalar data.");
2749 
2750  // set fill value
2751  HDF5Handle plist ( H5Pcreate(H5P_DATASET_CREATE), &H5Pclose, "HDF5File::createDataset(): unable to create property list." );
2752  H5Pset_fill_value(plist, TypeTraits::getH5DataType(), &init);
2753 
2754  // turn off time tagging of datasets by default.
2755  H5Pset_obj_track_times(plist, track_time);
2756 
2757  // enable chunks
2758  ArrayVector<hsize_t> chunks(defineChunks(chunkSize, shape, TypeTraits::numberOfBands(), compressionParameter));
2759  if(chunks.size() > 0)
2760  {
2761  std::reverse(chunks.begin(), chunks.end());
2762  H5Pset_chunk (plist, chunks.size(), chunks.begin());
2763  }
2764 
2765  // enable compression
2766  if(compressionParameter > 0)
2767  {
2768  H5Pset_deflate(plist, compressionParameter);
2769  }
2770 
2771  //create the dataset.
2772  HDF5HandleShared datasetHandle(H5Dcreate(parent, setname.c_str(),
2773  TypeTraits::getH5DataType(),
2774  dataspaceHandle, H5P_DEFAULT, plist, H5P_DEFAULT),
2775  &H5Dclose,
2776  "HDF5File::createDataset(): unable to create dataset.");
2777  if(parent != cGroupHandle_)
2778  H5Gclose(parent);
2779 
2780  return datasetHandle;
2781 }
2782 
2783 /********************************************************************/
2784 
2785 template<unsigned int N, class T, class Stride>
2786 void HDF5File::write_(std::string &datasetName,
2787  const MultiArrayView<N, T, Stride> & array,
2788  const hid_t datatype,
2789  const int numBandsOfType,
2790  typename MultiArrayShape<N>::type &chunkSize,
2791  int compressionParameter)
2792 {
2793  vigra_precondition(!isReadOnly(),
2794  "HDF5File::write(): file is read-only.");
2795 
2796  std::string groupname = SplitString(datasetName).first();
2797  std::string setname = SplitString(datasetName).last();
2798 
2799  // shape of the array. Add one dimension, if array contains non-scalars.
2800  ArrayVector<hsize_t> shape(array.shape().begin(), array.shape().end());
2801  std::reverse(shape.begin(), shape.end());
2802 
2803  if(numBandsOfType > 1)
2804  shape.push_back(numBandsOfType);
2805 
2806  HDF5Handle dataspace(H5Screate_simple(shape.size(), shape.begin(), NULL), &H5Sclose,
2807  "HDF5File::write(): Can not create dataspace.");
2808 
2809  // create and open group:
2810  std::string errorMessage ("HDF5File::write(): can not create group '" + groupname + "'.");
2811  HDF5Handle groupHandle(openCreateGroup_(groupname), &H5Gclose, errorMessage.c_str());
2812 
2813  // delete dataset, if it already exists
2814  deleteDataset_(groupHandle, setname.c_str());
2815 
2816  // set up properties list
2817  HDF5Handle plist(H5Pcreate(H5P_DATASET_CREATE), &H5Pclose,
2818  "HDF5File::write(): unable to create property list." );
2819 
2820  // turn off time tagging of datasets by default.
2821  H5Pset_obj_track_times(plist, track_time);
2822 
2823  // enable chunks
2824  ArrayVector<hsize_t> chunks(defineChunks(chunkSize, array.shape(), numBandsOfType, compressionParameter));
2825  if(chunks.size() > 0)
2826  {
2827  std::reverse(chunks.begin(), chunks.end());
2828  H5Pset_chunk (plist, chunks.size(), chunks.begin());
2829  }
2830 
2831  // enable compression
2832  if(compressionParameter > 0)
2833  {
2834  H5Pset_deflate(plist, compressionParameter);
2835  }
2836 
2837  // create dataset
2838  HDF5Handle datasetHandle(H5Dcreate(groupHandle, setname.c_str(), datatype, dataspace,H5P_DEFAULT, plist, H5P_DEFAULT),
2839  &H5Dclose, "HDF5File::write(): Can not create dataset.");
2840 
2841  herr_t status = 0;
2842  if(array.isUnstrided())
2843  {
2844  // Write the data directly from the array data buffer
2845  status = H5Dwrite(datasetHandle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data());
2846  }
2847  else
2848  {
2849  // otherwise, we need an intermediate buffer
2850  // FIXME: right now, the buffer has the same size as the array to be read
2851  // incomplete code for better solutions is below
2852  // MultiArray<N, T> buffer(array);
2853  // status = H5Dwrite(datasetHandle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, buffer.data());
2854 
2855  int offset = numBandsOfType > 1 ? 1 : 0;
2856  std::reverse(shape.begin(), shape.end());
2857  if(chunks.size() > 0)
2858  {
2859  // if the file is chunked, we use a buffer that matches the chunk size.
2860  std::reverse(chunks.begin(), chunks.end());
2861  }
2862  else
2863  {
2864  // otherwise, we compute a suitable chunk size.
2865  ArrayVector<hsize_t>(shape.size(), 1).swap(chunks);
2866  chunks[0] = numBandsOfType;
2867  MultiArrayIndex prod = 1;
2868  for(unsigned int k=0; k<N; ++k)
2869  {
2870  chunks[k+offset] = array.shape(k);
2871  prod *= array.shape(k);
2872  if(prod > 300000)
2873  break;
2874  }
2875  }
2876 
2877  ArrayVector<hsize_t> null(shape.size(), 0),
2878  start(shape.size(), 0),
2879  count(shape.size(), 1);
2880 
2881  count[N-1-offset] = numBandsOfType;
2882 
2883  typedef typename MultiArrayShape<N>::type Shape;
2884  Shape chunkCount, chunkMaxShape;
2885  for(unsigned int k=offset; k<chunks.size(); ++k)
2886  {
2887  chunkMaxShape[k-offset] = chunks[k];
2888  chunkCount[k-offset] = static_cast<MultiArrayIndex>(std::ceil(double(shape[k]) / chunks[k]));
2889  }
2890 
2891  typename CoupledIteratorType<N>::type chunkIter = createCoupledIterator(chunkCount),
2892  chunkEnd = chunkIter.getEndIterator();
2893  for(; chunkIter != chunkEnd; ++chunkIter)
2894  {
2895  Shape chunkStart(chunkIter.point() * chunkMaxShape),
2896  chunkStop(min(chunkStart + chunkMaxShape, array.shape()));
2897  MultiArray<N, T> buffer(array.subarray(chunkStart, chunkStop));
2898 
2899  for(unsigned int k=0; k<N; ++k)
2900  {
2901  start[N-1-k] = chunkStart[k];
2902  count[N-1-k] = buffer.shape(k);
2903  }
2904  if(offset == 1)
2905  {
2906  start[N] = 0;
2907  count[N] = numBandsOfType;
2908  }
2909  HDF5Handle filespace(H5Dget_space(datasetHandle),
2910  &H5Sclose, "HDF5File::write(): unable to create hyperslabs.");
2911  status = H5Sselect_hyperslab(filespace, H5S_SELECT_SET, start.data(), NULL, count.data(), NULL);
2912  if(status < 0)
2913  break;
2914 
2915  HDF5Handle dataspace2(H5Screate_simple(count.size(), count.data(), NULL),
2916  &H5Sclose, "HDF5File::write(): unable to create hyperslabs.");
2917  status = H5Sselect_hyperslab(dataspace2, H5S_SELECT_SET, null.data(), NULL, count.data(), NULL);
2918  if(status < 0)
2919  break;
2920 
2921  status = H5Dwrite(datasetHandle, datatype, dataspace2, filespace, H5P_DEFAULT, buffer.data());
2922  if(status < 0)
2923  break;
2924  }
2925  }
2926  vigra_postcondition(status >= 0,
2927  "HDF5File::write(): write to dataset '" + datasetName + "' via H5Dwrite() failed.");
2928 }
2929 
2930 /********************************************************************/
2931 
2932 template<unsigned int N, class T, class Stride>
2933 herr_t HDF5File::writeBlock_(HDF5HandleShared datasetHandle,
2934  typename MultiArrayShape<N>::type &blockOffset,
2935  const MultiArrayView<N, T, Stride> & array,
2936  const hid_t datatype,
2937  const int numBandsOfType)
2938 {
2939  vigra_precondition(!isReadOnly(),
2940  "HDF5File::writeBlock(): file is read-only.");
2941 
2942  ArrayVector<hsize_t> boffset, bshape, bones(N+1, 1);
2943  hssize_t dimensions = getDatasetDimensions_(datasetHandle);
2944  if(numBandsOfType > 1)
2945  {
2946  vigra_precondition(N+1 == dimensions,
2947  "HDF5File::readBlock(): Array dimension disagrees with data dimension.");
2948  bshape.resize(N+1);
2949  boffset.resize(N+1);
2950  bshape[N] = numBandsOfType;
2951  boffset[N] = 0;
2952  }
2953  else
2954  {
2955  vigra_precondition(N == dimensions,
2956  "HDF5File::readBlock(): Array dimension disagrees with data dimension.");
2957  bshape.resize(N);
2958  boffset.resize(N);
2959  }
2960 
2961  for(unsigned i = 0; i < N; ++i)
2962  {
2963  // vigra and hdf5 use different indexing
2964  bshape[N-1-i] = array.shape(i);
2965  boffset[N-1-i] = blockOffset[i];
2966  }
2967 
2968  // create a target dataspace in memory with the shape of the desired block
2969  HDF5Handle memspace_handle (H5Screate_simple(bshape.size(), bshape.data(), NULL),
2970  &H5Sclose,
2971  "Unable to get origin dataspace");
2972 
2973  // get file dataspace and select the desired block
2974  HDF5Handle dataspaceHandle (H5Dget_space(datasetHandle),&H5Sclose,"Unable to create target dataspace");
2975  H5Sselect_hyperslab(dataspaceHandle, H5S_SELECT_SET,
2976  boffset.data(), bones.data(), bones.data(), bshape.data());
2977 
2978  herr_t status = 0;
2979  if(array.isUnstrided())
2980  {
2981  // when the array is unstrided, we can read the data directly from the array buffer
2982  status = H5Dwrite( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, array.data());
2983  }
2984  else
2985  {
2986  // otherwise, we must copy the data into an unstrided extra buffer
2987  MultiArray<N, T> buffer(array);
2988  status = H5Dwrite( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, buffer.data());
2989  }
2990  return status;
2991 }
2992 
2993 /********************************************************************/
2994 
2995 template<unsigned int N, class T, class Stride>
2996 void HDF5File::write_attribute_(std::string name,
2997  const std::string & attribute_name,
2998  const MultiArrayView<N, T, Stride> & array,
2999  const hid_t datatype,
3000  const int numBandsOfType)
3001 {
3002  vigra_precondition(!isReadOnly(),
3003  "HDF5File::writeAttribute(): file is read-only.");
3004 
3005  // shape of the array. Add one dimension, if array contains non-scalars.
3006  ArrayVector<hsize_t> shape(array.shape().begin(), array.shape().end());
3007  std::reverse(shape.begin(), shape.end());
3008  if(numBandsOfType > 1)
3009  shape.push_back(numBandsOfType);
3010 
3011  HDF5Handle dataspace(H5Screate_simple(shape.size(),
3012  shape.begin(), NULL),
3013  &H5Sclose, "HDF5File::writeAttribute(): Can not"
3014  " create dataspace.");
3015 
3016  std::string errorMessage ("HDF5File::writeAttribute(): can not find "
3017  "object '" + name + "'.");
3018 
3019  H5O_type_t h5_type = get_object_type_(name);
3020  bool is_group = h5_type == H5O_TYPE_GROUP;
3021  if (!is_group && h5_type != H5O_TYPE_DATASET)
3022  vigra_precondition(0, "HDF5File::writeAttribute(): object \""
3023  + name + "\" is neither a group nor a "
3024  "dataset.");
3025  // get parent object handle
3026  HDF5Handle object_handle(is_group
3027  ? openCreateGroup_(name)
3028  : getDatasetHandle_(name),
3029  is_group
3030  ? &H5Gclose
3031  : &H5Dclose,
3032  errorMessage.c_str());
3033  // create / open attribute
3034  bool exists = existsAttribute(name, attribute_name);
3035  HDF5Handle attributeHandle(exists
3036  ? H5Aopen(object_handle,
3037  attribute_name.c_str(),
3038  H5P_DEFAULT)
3039  : H5Acreate(object_handle,
3040  attribute_name.c_str(), datatype,
3041  dataspace, H5P_DEFAULT,
3042  H5P_DEFAULT),
3043  &H5Aclose,
3044  "HDF5File::writeAttribute(): Can not create"
3045  " attribute.");
3046  herr_t status = 0;
3047  if(array.isUnstrided())
3048  {
3049  // write the data directly from the array data buffer
3050  status = H5Awrite(attributeHandle, datatype, array.data());
3051  }
3052  else
3053  {
3054  // write the data via an unstrided copy
3055  // (we assume that attributes are small arrays, so that the wasted memory is uncritical)
3056  MultiArray<N, T> buffer(array);
3057  status = H5Awrite(attributeHandle, datatype, buffer.data());
3058  }
3059  vigra_postcondition(status >= 0,
3060  "HDF5File::writeAttribute(): write to attribute '" + attribute_name + "' via H5Awrite() failed.");
3061 }
3062 
3063 /********************************************************************/
3064 
3065 template<unsigned int N, class T, class Stride>
3066 void HDF5File::read_(std::string datasetName,
3067  MultiArrayView<N, T, Stride> array,
3068  const hid_t datatype, const int numBandsOfType)
3069 {
3070  //Prepare to read without using HDF5ImportInfo
3071  ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName);
3072 
3073  std::string errorMessage ("HDF5File::read(): Unable to open dataset '" + datasetName + "'.");
3074  HDF5Handle datasetHandle(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
3075 
3076  // the object in the HDF5 file may have one additional dimension which we
3077  // interprete as the pixel type's bands
3078  int offset = (numBandsOfType > 1)
3079  ? 1
3080  : 0;
3081 
3082  vigra_precondition(MultiArrayIndex(N + offset) == MultiArrayIndex(dimshape.size()),
3083  "HDF5File::read(): Array dimension disagrees with dataset dimension.");
3084 
3085  typename MultiArrayShape<N>::type shape;
3086  for(int k=offset; k < (int)dimshape.size(); ++k)
3087  shape[k-offset] = (MultiArrayIndex)dimshape[k];
3088 
3089  vigra_precondition(shape == array.shape(),
3090  "HDF5File::read(): Array shape disagrees with dataset shape.");
3091  if (offset)
3092  vigra_precondition(dimshape[0] == static_cast<hsize_t>(numBandsOfType),
3093  "HDF5File::read(): Band count doesn't match destination array compound type.");
3094 
3095  herr_t status = 0;
3096  if(array.isUnstrided())
3097  {
3098  // when the array is unstrided, we can read the data directly into the array buffer
3099  status = H5Dread(datasetHandle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data() );
3100  }
3101  else
3102  {
3103  // otherwise, we need an intermediate buffer
3104 
3105  ArrayVector<hsize_t> null(dimshape.size(), 0),
3106  chunks(dimshape.size(), 1),
3107  start(dimshape.size(), 0),
3108  count(dimshape.size(), 1);
3109 
3110  HDF5Handle properties(H5Dget_create_plist(datasetHandle),
3111  &H5Pclose, "HDF5File::read(): failed to get property list");
3112  if(H5D_CHUNKED == H5Pget_layout(properties))
3113  {
3114  // if the file is chunked, we use a buffer that matches the chunk size.
3115  H5Pget_chunk(properties, static_cast<int>(chunks.size()), chunks.data());
3116  std::reverse(chunks.begin(), chunks.end());
3117  }
3118  else
3119  {
3120  // otherwise, we compute a suitable chunk size.
3121  chunks[0] = numBandsOfType;
3122  MultiArrayIndex prod = 1;
3123  for(unsigned int k=0; k<N; ++k)
3124  {
3125  chunks[k+offset] = array.shape(k);
3126  prod *= array.shape(k);
3127  if(prod > 300000)
3128  break;
3129  }
3130  }
3131 
3132  count[N-1-offset] = static_cast<hsize_t>(numBandsOfType);
3133 
3134  typedef typename MultiArrayShape<N>::type Shape;
3135  Shape chunkCount, chunkMaxShape;
3136  for(unsigned int k=offset; k<chunks.size(); ++k)
3137  {
3138  chunkMaxShape[k-offset] = chunks[k];
3139  chunkCount[k-offset] = (MultiArrayIndex)std::ceil(double(dimshape[k]) / chunks[k]);
3140  }
3141 
3142  typename CoupledIteratorType<N>::type chunkIter = createCoupledIterator(chunkCount),
3143  chunkEnd = chunkIter.getEndIterator();
3144  for(; chunkIter != chunkEnd; ++chunkIter)
3145  {
3146  Shape chunkStart(chunkIter.point() * chunkMaxShape),
3147  chunkStop(min(chunkStart + chunkMaxShape, array.shape()));
3148  MultiArray<N, T> buffer(chunkStop - chunkStart);
3149 
3150  for(unsigned int k=0; k<N; ++k)
3151  {
3152  start[N-1-k] = chunkStart[k];
3153  count[N-1-k] = buffer.shape(k);
3154  }
3155  if(offset == 1)
3156  {
3157  start[N] = 0;
3158  count[N] = numBandsOfType;
3159  }
3160  HDF5Handle filespace(H5Dget_space(datasetHandle),
3161  &H5Sclose, "HDF5File::read(): unable to create hyperslabs.");
3162  status = H5Sselect_hyperslab(filespace, H5S_SELECT_SET, start.data(), NULL, count.data(), NULL);
3163  if(status < 0)
3164  break;
3165 
3166  HDF5Handle dataspace(H5Screate_simple(count.size(), count.data(), NULL),
3167  &H5Sclose, "HDF5File::read(): unable to create hyperslabs.");
3168  status = H5Sselect_hyperslab(dataspace, H5S_SELECT_SET, null.data(), NULL, count.data(), NULL);
3169  if(status < 0)
3170  break;
3171 
3172  status = H5Dread(datasetHandle, datatype, dataspace, filespace, H5P_DEFAULT, buffer.data());
3173  if(status < 0)
3174  break;
3175 
3176  array.subarray(chunkStart, chunkStop) = buffer;
3177  }
3178  }
3179  vigra_postcondition(status >= 0,
3180  "HDF5File::read(): read from dataset '" + datasetName + "' via H5Dread() failed.");
3181 }
3182 
3183 /********************************************************************/
3184 
3185 template<unsigned int N, class T, class Stride>
3186 herr_t HDF5File::readBlock_(HDF5HandleShared datasetHandle,
3187  typename MultiArrayShape<N>::type &blockOffset,
3188  typename MultiArrayShape<N>::type &blockShape,
3189  MultiArrayView<N, T, Stride> array,
3190  const hid_t datatype, const int numBandsOfType)
3191 {
3192  vigra_precondition(blockShape == array.shape(),
3193  "HDF5File::readBlock(): Array shape disagrees with block size.");
3194 
3195  ArrayVector<hsize_t> boffset, bshape, bones(N+1, 1);
3196  hssize_t dimensions = getDatasetDimensions_(datasetHandle);
3197  if(numBandsOfType > 1)
3198  {
3199  vigra_precondition(N+1 == dimensions,
3200  "HDF5File::readBlock(): Array dimension disagrees with data dimension.");
3201  bshape.resize(N+1);
3202  boffset.resize(N+1);
3203  bshape[N] = numBandsOfType;
3204  boffset[N] = 0;
3205  }
3206  else
3207  {
3208  vigra_precondition(N == dimensions,
3209  "HDF5File::readBlock(): Array dimension disagrees with data dimension.");
3210  bshape.resize(N);
3211  boffset.resize(N);
3212  }
3213 
3214  for(unsigned i = 0; i < N; ++i)
3215  {
3216  // vigra and hdf5 use different indexing
3217  bshape[N-1-i] = blockShape[i];
3218  boffset[N-1-i] = blockOffset[i];
3219  }
3220 
3221  // create a target dataspace in memory with the shape of the desired block
3222  HDF5Handle memspace_handle(H5Screate_simple(bshape.size(), bshape.data(), NULL),
3223  &H5Sclose,
3224  "Unable to create target dataspace");
3225 
3226  // get file dataspace and select the desired block
3227  HDF5Handle dataspaceHandle(H5Dget_space(datasetHandle), &H5Sclose,
3228  "Unable to get dataspace");
3229  H5Sselect_hyperslab(dataspaceHandle, H5S_SELECT_SET,
3230  boffset.data(), bones.data(), bones.data(), bshape.data());
3231 
3232  herr_t status = 0;
3233  if(array.isUnstrided())
3234  {
3235  // when the array is unstrided, we can read the data directly into the array buffer
3236  status = H5Dread( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, array.data());
3237  }
3238  else
3239  {
3240  // otherwise, we need an unstrided extra buffer ...
3241  MultiArray<N, T> buffer(array.shape());
3242  status = H5Dread( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, buffer.data());
3243  // ... and must copy the values
3244  if(status >= 0)
3245  array = buffer;
3246  }
3247  return status;
3248 }
3249 
3250 /********************************************************************/
3251 
3252 template<unsigned int N, class T, class Stride>
3253 void HDF5File::read_attribute_(std::string datasetName,
3254  std::string attributeName,
3255  MultiArrayView<N, T, Stride> array,
3256  const hid_t datatype, const int numBandsOfType)
3257 {
3258  std::string dataset_path = get_absolute_path(datasetName);
3259  // open Attribute handle
3260  std::string message = "HDF5File::readAttribute(): could not get handle for attribute '"+attributeName+"'' of object '"+dataset_path+"'.";
3261  HDF5Handle attr_handle (H5Aopen_by_name(fileHandle_,dataset_path.c_str(),attributeName.c_str(),H5P_DEFAULT,H5P_DEFAULT),&H5Aclose, message.c_str());
3262 
3263  // get Attribute dataspace
3264  message = "HDF5File::readAttribute(): could not get dataspace for attribute '"+attributeName+"'' of object '"+dataset_path+"'.";
3265  HDF5Handle attr_dataspace_handle (H5Aget_space(attr_handle),&H5Sclose,message.c_str());
3266 
3267  // obtain Attribute shape
3268  int raw_dims = H5Sget_simple_extent_ndims(attr_dataspace_handle);
3269  int dims = std::max(raw_dims, 1); // scalar attributes may be stored with raw_dims==0
3270  ArrayVector<hsize_t> dimshape(dims);
3271  if(raw_dims > 0)
3272  H5Sget_simple_extent_dims(attr_dataspace_handle, dimshape.data(), NULL);
3273  else
3274  dimshape[0] = 1;
3275 
3276  // invert the dimensions to guarantee VIGRA-compatible order
3277  std::reverse(dimshape.begin(), dimshape.end());
3278 
3279  int offset = (numBandsOfType > 1)
3280  ? 1
3281  : 0;
3282  message = "HDF5File::readAttribute(): Array dimension disagrees with dataset dimension.";
3283  // the object in the HDF5 file may have one additional dimension which we then interpret as the pixel type bands
3284  vigra_precondition(MultiArrayIndex(N + offset) == MultiArrayIndex(dims), message);
3285 
3286  for(int k=offset; k < (int)dimshape.size(); ++k)
3287  vigra_precondition(array.shape()[k-offset] == (MultiArrayIndex)dimshape[k],
3288  "HDF5File::readAttribute(): Array shape disagrees with dataset shape");
3289 
3290  herr_t status = 0;
3291  if(array.isUnstrided())
3292  {
3293  // when the array is unstrided, we can read the data directly into the array buffer
3294  status = H5Aread( attr_handle, datatype, array.data());
3295  }
3296  else
3297  {
3298  // otherwise, we need an unstrided extra buffer ...
3299  // (we assume that attributes are small arrays, so that the wasted memory is uncritical)
3300  MultiArray<N, T> buffer(array.shape());
3301  status = H5Aread( attr_handle, datatype, buffer.data() );
3302  // ... and must copy the values
3303  if(status >= 0)
3304  array = buffer;
3305  }
3306  vigra_postcondition(status >= 0,
3307  "HDF5File::readAttribute(): read from attribute '" + attributeName + "' via H5Aread() failed.");
3308 }
3309 
3310 /********************************************************************/
3311 
3312 /** \brief Read the data specified by the given \ref vigra::HDF5ImportInfo object
3313  and write the into the given 'array'.
3314 
3315  The array must have the correct number of dimensions and shape for the dataset
3316  represented by 'info'. When the element type of 'array' differs from the stored element
3317  type, HDF5 will convert the type on the fly (except when the HDF5 version is 1.6 or below,
3318  in which case an error will result). Multi-channel element types (i.e. \ref vigra::RGBValue,
3319  \ref vigra::TinyVector, and \ref vigra::FFTWComplex) are recognized and handled correctly.
3320 
3321  <b> Declaration:</b>
3322 
3323  \code
3324  namespace vigra {
3325  template<unsigned int N, class T, class StrideTag>
3326  void
3327  readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, StrideTag> array);
3328  }
3329  \endcode
3330 
3331  <b> Usage:</b>
3332 
3333  <b>\#include</b> <vigra/hdf5impex.hxx><br>
3334  Namespace: vigra
3335 
3336  \code
3337 
3338  HDF5ImportInfo info(filename, dataset_name);
3339  vigra_precondition(info.numDimensions() == 3, "Dataset must be 3-dimensional.");
3340 
3341  MultiArrayShape<3>::type shape(info.shape().begin());
3342  MultiArray<3, int> array(shape);
3343 
3344  readHDF5(info, array);
3345  \endcode
3346 */
3347 doxygen_overloaded_function(template <...> void readHDF5)
3348 
3349 template<unsigned int N, class T, class StrideTag>
3350 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, StrideTag> array)
3351 {
3352  HDF5File file(info.getFilePath(), HDF5File::OpenReadOnly);
3353  file.read(info.getPathInFile(), array);
3354 }
3355 
3356 inline hid_t openGroup(hid_t parent, std::string group_name)
3357 {
3358  //std::cout << group_name << std::endl;
3359  size_t last_slash = group_name.find_last_of('/');
3360  if (last_slash == std::string::npos || last_slash != group_name.size() - 1)
3361  group_name = group_name + '/';
3362  std::string::size_type begin = 0, end = group_name.find('/');
3363  int ii = 0;
3364  while (end != std::string::npos)
3365  {
3366  std::string group(group_name.begin()+begin, group_name.begin()+end);
3367  hid_t prev_parent = parent;
3368  parent = H5Gopen(prev_parent, group.c_str(), H5P_DEFAULT);
3369 
3370  if(ii != 0) H5Gclose(prev_parent);
3371  if(parent < 0) return parent;
3372  ++ii;
3373  begin = end + 1;
3374  end = group_name.find('/', begin);
3375  }
3376  return parent;
3377 }
3378 
3379 inline hid_t createGroup(hid_t parent, std::string group_name)
3380 {
3381  if(group_name.size() == 0 ||*group_name.rbegin() != '/')
3382  group_name = group_name + '/';
3383  if(group_name == "/")
3384  return H5Gopen(parent, group_name.c_str(), H5P_DEFAULT);
3385 
3386  std::string::size_type begin = 0, end = group_name.find('/');
3387  int ii = 0;
3388  while (end != std::string::npos)
3389  {
3390  std::string group(group_name.begin()+begin, group_name.begin()+end);
3391  hid_t prev_parent = parent;
3392 
3393  if(H5LTfind_dataset(parent, group.c_str()) == 0)
3394  {
3395  parent = H5Gcreate(prev_parent, group.c_str(), H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
3396  } else {
3397  parent = H5Gopen(prev_parent, group.c_str(), H5P_DEFAULT);
3398  }
3399 
3400  if(ii != 0) H5Gclose(prev_parent);
3401  if(parent < 0) return parent;
3402  ++ii;
3403  begin = end + 1;
3404  end = group_name.find('/', begin);
3405  }
3406  return parent;
3407 }
3408 
3409 inline void deleteDataset(hid_t parent, std::string dataset_name)
3410 {
3411  // delete existing data and create new dataset
3412  if(H5LTfind_dataset(parent, dataset_name.c_str()))
3413  {
3414  //std::cout << "dataset already exists" << std::endl;
3415 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6)
3416  if(H5Gunlink(parent, dataset_name.c_str()) < 0)
3417  {
3418  vigra_postcondition(false, "writeToHDF5File(): Unable to delete existing data.");
3419  }
3420 #else
3421  if(H5Ldelete(parent, dataset_name.c_str(), H5P_DEFAULT ) < 0)
3422  {
3423  vigra_postcondition(false, "createDataset(): Unable to delete existing data.");
3424  }
3425 #endif
3426  }
3427 }
3428 
3429 inline hid_t createFile(std::string filePath, bool append_ = true)
3430 {
3431  FILE * pFile;
3432  pFile = fopen ( filePath.c_str(), "r" );
3433  hid_t file_id;
3434  if ( pFile == NULL )
3435  {
3436  file_id = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
3437  }
3438  else if(append_)
3439  {
3440  fclose( pFile );
3441  file_id = H5Fopen(filePath.c_str(), H5F_ACC_RDWR, H5P_DEFAULT);
3442  }
3443  else
3444  {
3445  fclose(pFile);
3446  std::remove(filePath.c_str());
3447  file_id = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
3448  }
3449  return file_id;
3450 }
3451 
3452 template<unsigned int N, class T, class Tag>
3453 void createDataset(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, Tag> & array, const hid_t datatype, const int numBandsOfType, HDF5Handle & file_handle, HDF5Handle & dataset_handle)
3454 {
3455  std::string path_name(pathInFile), group_name, data_set_name, message;
3456  std::string::size_type delimiter = path_name.rfind('/');
3457 
3458  //create or open file
3459  file_handle = HDF5Handle(createFile(filePath), &H5Fclose,
3460  "createDataset(): unable to open output file.");
3461 
3462  // get the groupname and the filename
3463  if(delimiter == std::string::npos)
3464  {
3465  group_name = "/";
3466  data_set_name = path_name;
3467  }
3468  else
3469  {
3470  group_name = std::string(path_name.begin(), path_name.begin()+delimiter);
3471  data_set_name = std::string(path_name.begin()+delimiter+1, path_name.end());
3472  }
3473 
3474  // create all groups
3475  HDF5Handle group(createGroup(file_handle, group_name), &H5Gclose,
3476  "createDataset(): Unable to create and open group. generic v");
3477 
3478  // delete the dataset if it already exists
3479  deleteDataset(group, data_set_name);
3480 
3481  // create dataspace
3482  // add an extra dimension in case that the data is non-scalar
3483  HDF5Handle dataspace_handle;
3484  if(numBandsOfType > 1) {
3485  // invert dimensions to guarantee c-order
3486  hsize_t shape_inv[N+1]; // one additional dimension for pixel type channel(s)
3487  for(unsigned int k=0; k<N; ++k) {
3488  shape_inv[N-1-k] = array.shape(k); // the channels (eg of an RGB image) are represented by the first dimension (before inversion)
3489  //std::cout << shape_inv[N-k] << " (" << N << ")";
3490  }
3491  shape_inv[N] = numBandsOfType;
3492 
3493  // create dataspace
3494  dataspace_handle = HDF5Handle(H5Screate_simple(N+1, shape_inv, NULL),
3495  &H5Sclose, "createDataset(): unable to create dataspace for non-scalar data.");
3496  } else {
3497  // invert dimensions to guarantee c-order
3498  hsize_t shape_inv[N];
3499  for(unsigned int k=0; k<N; ++k)
3500  shape_inv[N-1-k] = array.shape(k);
3501 
3502  // create dataspace
3503  dataspace_handle = HDF5Handle(H5Screate_simple(N, shape_inv, NULL),
3504  &H5Sclose, "createDataset(): unable to create dataspace for scalar data.");
3505  }
3506 
3507  //alloc memory for dataset.
3508  dataset_handle = HDF5Handle(H5Dcreate(group,
3509  data_set_name.c_str(),
3510  datatype,
3511  dataspace_handle,
3512  H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT),
3513  &H5Dclose, "createDataset(): unable to create dataset.");
3514 }
3515 
3516 
3517 
3518 
3519 /** \brief Store array data in an HDF5 file.
3520 
3521  The number of dimensions, shape and element type of the stored dataset is automatically
3522  determined from the properties of the given \a array. Strided arrays are stored in an
3523  unstrided way, i.e. in contiguous scan-order. Multi-channel element types
3524  (i.e. \ref vigra::RGBValue, \ref vigra::TinyVector and \ref vigra::FFTWComplex)
3525  are recognized and handled correctly
3526  (in particular, the will form the innermost dimension of the stored dataset).
3527  \a pathInFile may contain '/'-separated group names, but must end with the name
3528  of the dataset to be created.
3529 
3530  <b> Declaration:</b>
3531 
3532  \code
3533  namespace vigra {
3534  template<unsigned int N, class T, class StrideTag>
3535  void
3536  writeHDF5(const char* filePath, const char* pathInFile,
3537  MultiArrayView<N, T, StrideTag>const & array);
3538  }
3539  \endcode
3540 
3541  <b> Usage:</b>
3542 
3543  <b>\#include</b> <vigra/hdf5impex.hxx><br>
3544  Namespace: vigra
3545 
3546  \code
3547  MultiArrayShape<3>::type shape(100, 200, 20);
3548  MultiArray<3, int> array(shape);
3549  ... // fill array with data
3550 
3551  writeHDF5("mydata.h5", "/group1/my_dataset", array);
3552  \endcode
3553 */
3554 doxygen_overloaded_function(template <...> void writeHDF5)
3555 
3556 template<unsigned int N, class T, class StrideTag>
3557 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, StrideTag> & array)
3558 {
3559  HDF5File file(filePath, HDF5File::Open);
3560  file.write(pathInFile, array);
3561 }
3562 
3563 namespace detail
3564 {
3565 struct MaxSizeFnc
3566 {
3567  size_t size;
3568 
3569  MaxSizeFnc()
3570  : size(0)
3571  {}
3572 
3573  void operator()(std::string const & in)
3574  {
3575  size = in.size() > size ?
3576  in.size() :
3577  size;
3578  }
3579 };
3580 }
3581 
3582 
3583 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR == 8) || DOXYGEN
3584 /** Write a numeric MultiArray as an attribute with name \a name
3585  of the dataset specified by the handle \a loc.
3586 
3587  <b>\#include</b> <vigra/hdf5impex.hxx><br>
3588  Namespace: vigra
3589 */
3590 template<size_t N, class T, class C>
3591 void writeHDF5Attr(hid_t loc,
3592  const char* name,
3593  MultiArrayView<N, T, C> const & array)
3594 {
3595  if(H5Aexists(loc, name) > 0)
3596  H5Adelete(loc, name);
3597 
3598  ArrayVector<hsize_t> shape(array.shape().begin(),
3599  array.shape().end());
3600  HDF5Handle
3601  dataspace_handle(H5Screate_simple(N, shape.data(), NULL),
3602  &H5Sclose,
3603  "writeToHDF5File(): unable to create dataspace.");
3604 
3605  HDF5Handle attr(H5Acreate(loc,
3606  name,
3607  detail::getH5DataType<T>(),
3608  dataspace_handle,
3609  H5P_DEFAULT ,H5P_DEFAULT ),
3610  &H5Aclose,
3611  "writeHDF5Attr: unable to create Attribute");
3612 
3613  //copy data - since attributes are small - who cares!
3614  ArrayVector<T> buffer;
3615  for(int ii = 0; ii < array.size(); ++ii)
3616  buffer.push_back(array[ii]);
3617  H5Awrite(attr, detail::getH5DataType<T>(), buffer.data());
3618 }
3619 
3620 /** Write a string MultiArray as an attribute with name \a name
3621  of the dataset specified by the handle \a loc.
3622 
3623  <b>\#include</b> <vigra/hdf5impex.hxx><br>
3624  Namespace: vigra
3625 */
3626 template<size_t N, class C>
3627 void writeHDF5Attr(hid_t loc,
3628  const char* name,
3629  MultiArrayView<N, std::string, C> const & array)
3630 {
3631  if(H5Aexists(loc, name) > 0)
3632  H5Adelete(loc, name);
3633 
3634  ArrayVector<hsize_t> shape(array.shape().begin(),
3635  array.shape().end());
3636  HDF5Handle
3637  dataspace_handle(H5Screate_simple(N, shape.data(), NULL),
3638  &H5Sclose,
3639  "writeToHDF5File(): unable to create dataspace.");
3640 
3641  HDF5Handle atype(H5Tcopy (H5T_C_S1),
3642  &H5Tclose,
3643  "writeToHDF5File(): unable to create type.");
3644 
3645  detail::MaxSizeFnc max_size;
3646  max_size = std::for_each(array.data(),array.data()+ array.size(), max_size);
3647  H5Tset_size (atype, max_size.size);
3648 
3649  HDF5Handle attr(H5Acreate(loc,
3650  name,
3651  atype,
3652  dataspace_handle,
3653  H5P_DEFAULT ,H5P_DEFAULT ),
3654  &H5Aclose,
3655  "writeHDF5Attr: unable to create Attribute");
3656 
3657  std::string buf ="";
3658  for(int ii = 0; ii < array.size(); ++ii)
3659  {
3660  buf = buf + array[ii]
3661  + std::string(max_size.size - array[ii].size(), ' ');
3662  }
3663  H5Awrite(attr, atype, buf.c_str());
3664 }
3665 
3666 /** Write a numeric ArrayVectorView as an attribute with name \a name
3667  of the dataset specified by the handle \a loc.
3668 
3669  <b>\#include</b> <vigra/hdf5impex.hxx><br>
3670  Namespace: vigra
3671 */
3672 template<class T>
3673 inline void writeHDF5Attr( hid_t loc,
3674  const char* name,
3675  ArrayVectorView<T> & array)
3676 {
3677  writeHDF5Attr(loc, name,
3679  array.data()));
3680 }
3681 
3682 /** write an Attribute given a file and a path in the file.
3683  the path in the file should have the format
3684  [attribute] or /[subgroups/]dataset.attribute or
3685  /[subgroups/]group.attribute.
3686  The attribute is written to the root group, a dataset or a subgroup
3687  respectively
3688 */
3689 template<class Arr>
3690 inline void writeHDF5Attr( std::string filePath,
3691  std::string pathInFile,
3692  Arr & ar)
3693 {
3694  std::string path_name(pathInFile), group_name, data_set_name, message, attr_name;
3695  std::string::size_type delimiter = path_name.rfind('/');
3696 
3697  //create or open file
3698  HDF5Handle file_id(createFile(filePath), &H5Fclose,
3699  "writeToHDF5File(): unable to open output file.");
3700 
3701  // get the groupname and the filename
3702  if(delimiter == std::string::npos)
3703  {
3704  group_name = "/";
3705  data_set_name = path_name;
3706  }
3707 
3708  else
3709  {
3710  group_name = std::string(path_name.begin(), path_name.begin()+delimiter);
3711  data_set_name = std::string(path_name.begin()+delimiter+1, path_name.end());
3712  }
3713  delimiter = data_set_name.rfind('.');
3714  if(delimiter == std::string::npos)
3715  {
3716  attr_name = path_name;
3717  data_set_name = "/";
3718  }
3719  else
3720  {
3721  attr_name = std::string(data_set_name.begin()+delimiter+1, data_set_name.end());
3722  data_set_name = std::string(data_set_name.begin(), data_set_name.begin()+delimiter);
3723  }
3724 
3725  HDF5Handle group(openGroup(file_id, group_name), &H5Gclose,
3726  "writeToHDF5File(): Unable to create and open group. attr ver");
3727 
3728  if(data_set_name != "/")
3729  {
3730  HDF5Handle dset(H5Dopen(group, data_set_name.c_str(), H5P_DEFAULT), &H5Dclose,
3731  "writeHDF5Attr():unable to open dataset");
3732  writeHDF5Attr(hid_t(dset), attr_name.c_str(), ar);
3733  }
3734  else
3735  {
3736  writeHDF5Attr(hid_t(group), attr_name.c_str(), ar);
3737  }
3738 
3739 }
3740 #endif
3741 
3742 //@}
3743 
3744 } // namespace vigra
3745 
3746 #endif // VIGRA_HDF5IMPEX_HXX
HDF5File & operator=(HDF5File const &other)
Assign a HDF5File object.
Definition: hdf5impex.hxx:1127
void read(const std::string &datasetName, ArrayVectorView< T > array)
Read data into an array vector. If the first character of datasetName is a "/", the path will be inte...
Definition: hdf5impex.hxx:1956
HDF5Handle()
Default constructor. Creates a NULL handle.
Definition: hdf5impex.hxx:205
void write(const std::string &datasetName, const ArrayVectorView< T > &array, int compression=0)
Write array vectors.
Definition: hdf5impex.hxx:1843
bool cd_up(int levels)
Change the current group to its parent group. Returns true if successful, false otherwise. If unsuccessful, the group will not change.
Definition: hdf5impex.hxx:1224
std::vector< std::string > listAttributes(std::string const &group_or_dataset) const
List the attribute names of the given group or dataset.
Definition: hdf5impex.hxx:1514
bool operator==(HDF5Handle const &h) const
Equality comparison of the contained handle.
Definition: hdf5impex.hxx:346
Temporarily disable HDF5's native error output.
Definition: hdf5impex.hxx:130
ArrayVector< hsize_t > getDatasetShape(std::string datasetName) const
Get the shape of each dimension of a certain dataset.
Definition: hdf5impex.hxx:1371
Wrapper for unique hid_t objects.
Definition: hdf5impex.hxx:191
HDF5File(std::string filePath, OpenMode mode, bool track_creation_times=false)
Open or create an HDF5File object.
Definition: hdf5impex.hxx:1052
MultiArrayIndex numDimensions() const
bool operator!=(hid_t h) const
Inequality comparison of the contained handle.
Definition: hdf5impex.hxx:617
HDF5HandleShared()
Default constructor. Creates a NULL handle.
Definition: hdf5impex.hxx:427
const difference_type & shape() const
Definition: multi_array.hxx:1594
void readAttribute(std::string object_name, std::string attribute_name, MultiArrayView< N, T, Stride > array)
Read MultiArray Attributes. In contrast to datasets, subarray access is not available.
Definition: hdf5impex.hxx:1637
void cd_mk(std::string groupName)
Change the current group; create it if necessary. If the first character is a "/", the path will be interpreted as absolute path, otherwise it will be interpreted as path relative to the current group.
Definition: hdf5impex.hxx:1262
hssize_t getDatasetDimensions(std::string datasetName) const
Get the number of dimensions of a certain dataset If the first character is a "/", the path will be interpreted as absolute path, otherwise it will be interpreted as path relative to the current group.
Definition: hdf5impex.hxx:1342
Definition: array_vector.hxx:76
const_iterator begin() const
Definition: array_vector.hxx:223
pointer data() const
Definition: multi_array.hxx:1844
HDF5ImportInfo(const char *filePath, const char *pathInFile)
HDF5Handle(hid_t h, Destructor destructor, const char *error_message)
Create a wrapper for a hid_t object.
Definition: hdf5impex.hxx:229
std::string get_absolute_path(std::string path) const
takes any path and converts it into an absolute path in the current file.
Definition: hdf5impex.hxx:2247
void swap(HDF5HandleShared &h)
Swap the contents of two handle wrappers.
Definition: hdf5impex.hxx:567
void reshape(const difference_type &shape)
Definition: multi_array.hxx:2807
HDF5HandleShared getDatasetHandleShared(std::string const &datasetName) const
Obtain a shared HDF5 handle of a dataset.
Definition: hdf5impex.hxx:1466
void mkdir(std::string groupName)
Create a new group. If the first character is a "/", the path will be interpreted as absolute path...
Definition: hdf5impex.hxx:1245
bool existsAttribute(std::string object_name, std::string attribute_name)
Test if attribute exists.
Definition: hdf5impex.hxx:1620
const std::string & getFilePath() const
std::ptrdiff_t MultiArrayIndex
Definition: multi_fwd.hxx:60
void writeHDF5Attr(hid_t loc, const char *name, MultiArrayView< N, T, C > const &array)
Definition: hdf5impex.hxx:3591
HDF5HandleShared createDataset(std::string datasetName, TinyVector< MultiArrayIndex, N > const &shape, typename detail::HDF5TypeTraits< T >::value_type init=typename detail::HDF5TypeTraits< T >::value_type(), TinyVector< MultiArrayIndex, N > const &chunkSize=(TinyVector< MultiArrayIndex, N >()), int compressionParameter=0)
Create a new dataset. This function can be used to create a dataset filled with a default value init...
Definition: hdf5impex.hxx:2709
bool operator==(hid_t h) const
Equality comparison of the contained handle.
Definition: hdf5impex.hxx:353
hid_t release()
Return the contained handle and set the wrapper to NULL without calling close().
Definition: hdf5impex.hxx:290
const std::string & getPathInFile() const
herr_t close()
Close the handle if this is the unique (i.e. last) owner.
Definition: hdf5impex.hxx:508
#define VIGRA_FINALLY(destructor)
Definition: utilities.hxx:206
void ls(Container &cont) const
List the contents of the current group into a container-like object via insert(). ...
Definition: hdf5impex.hxx:1309
void write(std::string datasetName, const MultiArrayView< N, T, Stride > &array, typename MultiArrayShape< N >::type chunkSize, int compression=0)
Write multi arrays. Chunks can be activated by providing a MultiArrayShape as chunkSize. chunkSize must have equal dimension as array.
Definition: hdf5impex.hxx:1762
Argument object for the function readHDF5().
Definition: hdf5impex.hxx:663
bool isHDF5(char const *filename)
Check if given filename refers to a HDF5 file.
Definition: hdf5impex.hxx:102
void reset(hid_t h, Destructor destructor, const char *error_message)
Reset the handle to a new value.
Definition: hdf5impex.hxx:531
HDF5File(HDF5File const &other)
Copy a HDF5File object.
Definition: hdf5impex.hxx:1103
HDF5HandleShared & operator=(HDF5HandleShared const &h)
Assignment. Call close() for the present LHS handle and share ownership with the RHS handle (analogou...
Definition: hdf5impex.hxx:479
HDF5HandleShared(HDF5HandleShared const &h)
Copy constructor. Shares ownership with the RHS handle (analogous to std::shared_ptr).
Definition: hdf5impex.hxx:466
bool operator!=(hid_t h) const
Inequality comparison of the contained handle.
Definition: hdf5impex.hxx:367
bool operator!=(HDF5Handle const &h) const
Inequality comparison of the contained handle.
Definition: hdf5impex.hxx:360
difference_type_1 size() const
Definition: multi_array.hxx:1587
Wrapper for shared hid_t objects.
Definition: hdf5impex.hxx:412
std::string filename() const
Get the name of the associated file.
Definition: hdf5impex.hxx:1324
HDF5File()
Default constructor.
Definition: hdf5impex.hxx:1031
HDF5Handle & operator=(HDF5Handle const &h)
Assignment. Calls close() for the LHS handle and hands over ownership of the RHS handle (analogous to...
Definition: hdf5impex.hxx:252
~HDF5File()
The destructor flushes and closes the file.
Definition: hdf5impex.hxx:1114
NumericTraits< V >::Promote prod(TinyVectorBase< V, SIZE, D1, D2 > const &l)
product of the vector's elements
Definition: tinyvector.hxx:2097
bool operator==(HDF5HandleShared const &h) const
Equality comparison of the contained handle.
Definition: hdf5impex.hxx:596
HDF5File(HDF5HandleShared const &fileHandle, const std::string &pathname="", bool read_only=false)
Initialize an HDF5File object from HDF5 file handle.
Definition: hdf5impex.hxx:1075
HDF5Handle getAttributeHandle(std::string dataset_name, std::string attribute_name) const
Obtain the HDF5 handle of a attribute.
Definition: hdf5impex.hxx:1537
void readAndResize(std::string datasetName, MultiArray< N, T, Alloc > &array)
Read data into a MultiArray. Resize MultiArray to the correct size. If the first character of dataset...
Definition: hdf5impex.hxx:1930
bool operator==(hid_t h) const
Equality comparison of the contained handle.
Definition: hdf5impex.hxx:603
void readHDF5(...)
Read the data specified by the given vigra::HDF5ImportInfo object and write the into the given 'array...
void readAndResize(std::string datasetName, ArrayVector< T > &array)
Read data into an array vector. Resize the array vector to the correct size. If the first character o...
Definition: hdf5impex.hxx:1969
~HDF5Handle()
Destructor. Calls close() for the contained handle.
Definition: hdf5impex.hxx:267
bool operator!=(HDF5HandleShared const &h) const
Inequality comparison of the contained handle.
Definition: hdf5impex.hxx:610
void read(std::string datasetName, char &data)
Read a single value. Specialization of the read function for simple datatypes.
Definition: hdf5impex.hxx:2101
doxygen_overloaded_function(template<...> void separableConvolveBlockwise) template< unsigned int N
Separated convolution on ChunkedArrays.
void open(std::string filePath, OpenMode mode)
Open or create the given file in the given mode and set the group to "/". If another file is currentl...
Definition: hdf5impex.hxx:1164
PixelType pixelType() const
OpenMode
Set how a file is opened.
Definition: hdf5impex.hxx:1015
void readAttribute(std::string object_name, std::string attribute_name, char &data)
Read a single value. Specialization of the read function for simple datatypes.
Definition: hdf5impex.hxx:1672
void writeBlock(std::string datasetName, typename MultiArrayShape< N >::type blockOffset, const MultiArrayView< N, T, Stride > &array)
Write a multi array into a larger volume. blockOffset determines the position, where array is written...
Definition: hdf5impex.hxx:1784
void flushToDisk()
Immediately write all data to disk.
Definition: hdf5impex.hxx:2173
HDF5Handle(HDF5Handle const &h)
Copy constructor.
Definition: hdf5impex.hxx:241
void reset(hid_t h, Destructor destructor, const char *error_message)
Reset the wrapper to a new handle.
Definition: hdf5impex.hxx:302
Class for fixed size vectors.This class contains an array of size SIZE of the specified VALUETYPE...
Definition: accessor.hxx:940
void writeHDF5(...)
Store array data in an HDF5 file.
MultiArrayIndex shapeOfDimension(const int dim) const
HDF5Handle getGroupHandle(std::string group_name, std::string function_name="HDF5File::getGroupHandle()")
Obtain the HDF5 handle of a group (create the group if it doesn't exist).
Definition: hdf5impex.hxx:1474
std::string pwd() const
Get the path of the current group.
Definition: hdf5impex.hxx:1317
void read(std::string datasetName, MultiArrayView< N, T, Stride > array)
Read data into a multi array. If the first character of datasetName is a "/", the path will be interp...
Definition: hdf5impex.hxx:1912
void listAttributes(std::string const &group_or_dataset, Container &container) const
Insert the attribute names of the given group or dataset into the given container by calling containe...
Definition: hdf5impex.hxx:1529
bool unique() const
Check if this is the unique owner of the conatined handle.
Definition: hdf5impex.hxx:558
void write(std::string datasetName, char data)
Write a single value. Specialization of the write function for simple datatypes.
Definition: hdf5impex.hxx:1883
const char * getPixelType() const
void writeAttribute(std::string object_name, std::string attribute_name, char data)
Write a single value. Specialization of the write function for simple datatypes.
Definition: hdf5impex.hxx:1585
void writeAttribute(std::string object_name, std::string attribute_name, const MultiArrayView< N, T, Stride > &array)
Write MultiArray Attributes. In contrast to datasets, subarray access, chunks and compression are not...
Definition: hdf5impex.hxx:1550
herr_t close()
Explicitly call the stored destructor (if one has been registered in the constructor) for the contain...
Definition: hdf5impex.hxx:277
void write(std::string datasetName, const MultiArrayView< N, T, Stride > &array, int iChunkSize=0, int compression=0)
Write multi arrays.
Definition: hdf5impex.hxx:1730
image import and export functions
bool isUnstrided(unsigned int dimension=N-1) const
Definition: multi_array.hxx:1232
void close()
Close the current file.
Definition: hdf5impex.hxx:1176
bool cd_up()
Change the current group to its parent group. Returns true if successful, false otherwise. If unsuccessful, the group will not change.
Definition: hdf5impex.hxx:1202
Base class for, and view to, vigra::MultiArray.
Definition: multi_array.hxx:650
std::string getDatasetType(std::string const &datasetName) const
Definition: hdf5impex.hxx:1412
bool existsDataset(std::string datasetName) const
Check if given datasetName exists.
Definition: hdf5impex.hxx:1331
const_iterator end() const
Definition: array_vector.hxx:237
const_pointer data() const
Definition: array_vector.hxx:209
ArrayVector< hsize_t > const & shape() const
Definition: hdf5impex.hxx:712
size_type size() const
Definition: array_vector.hxx:358
MultiArrayView subarray(difference_type p, difference_type q) const
Definition: multi_array.hxx:1474
void readBlock(std::string datasetName, typename MultiArrayShape< N >::type blockOffset, typename MultiArrayShape< N >::type blockShape, MultiArrayView< N, T, Stride > array)
Read a block of data into a multi array. This function allows to read a small block out of a larger v...
Definition: hdf5impex.hxx:2007
Class for a single RGB value.
Definition: accessor.hxx:938
int ceil(FixedPoint< IntBits, FracBits > v)
rounding up.
Definition: fixedpoint.hxx:675
hid_t getDatasetHandle() const
HDF5HandleShared(hid_t h, Destructor destructor, const char *error_message)
Create a shared wrapper for a plain hid_t object.
Definition: hdf5impex.hxx:452
hid_t getH5FileHandle() const
~HDF5HandleShared()
Destructor (calls close())
Definition: hdf5impex.hxx:495
HDF5Handle getDatasetHandle(std::string const &datasetName) const
Obtain the HDF5 handle of a dataset.
Definition: hdf5impex.hxx:1458
void root()
Change current group to "/".
Definition: hdf5impex.hxx:1184
HDF5File(bool track_creation_times)
Construct with time tagging of datasets enabled.
Definition: hdf5impex.hxx:1039
void cd(std::string groupName)
Change the current group. Both absolute and relative group names are allowed.
Definition: hdf5impex.hxx:1193
size_t use_count() const
Get the number of owners of the contained handle.
Definition: hdf5impex.hxx:547
void swap(HDF5Handle &h)
Swap the contents of two handle wrappers.
Definition: hdf5impex.hxx:318
Access to HDF5 files.
Definition: hdf5impex.hxx:958
std::vector< std::string > ls() const
List the contents of the current group. The function returns a vector of strings holding the entries ...
Definition: hdf5impex.hxx:1287

© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de)
Heidelberg Collaboratory for Image Processing, University of Heidelberg, Germany

html generated using doxygen and Python
vigra 1.11.0 (Thu Mar 17 2016)