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

© 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.10.0 (Thu Jan 8 2015)