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

© 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 (Fri Feb 21 2014)