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

hdf5impex.hxx
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 # include <H5LT.h>
66 #else
67 # include <hdf5_hl.h>
68 #endif
69 
70 #include "impex.hxx"
71 #include "multi_array.hxx"
72 #include "multi_impex.hxx"
73 
74 namespace vigra {
75 
76 /** \addtogroup VigraHDF5Impex Import/Export of Images and Arrays in HDF5 Format
77 
78  Supports arrays with arbitrary element types and arbitrary many dimensions.
79  See the <a href="http://www.hdfgroup.org/HDF5/">HDF5 Website</a> for more
80  information on the HDF5 file format.
81 **/
82 //@{
83 
84  /** \brief Wrapper for hid_t objects.
85 
86  Newly created or opened HDF5 handles are usually stored as objects of type 'hid_t'. When the handle
87  is no longer needed, the appropriate close function must be called. However, if a function is
88  aborted by an exception, this is difficult to ensure. Class HDF5Handle is a smart pointer that
89  solves this problem by calling the close function in the destructor (This is analogous to how
90  std::auto_ptr calls 'delete' on the contained pointer). A pointer to the close function must be
91  passed to the constructor, along with an error message that is raised when creation/opening fails.
92 
93  Since HDF5Handle objects are convertible to hid_t, they can be used in the code in place
94  of the latter.
95 
96  <b>Usage:</b>
97 
98  \code
99  HDF5Handle file_id(H5Fopen(filename, H5F_ACC_RDWR, H5P_DEFAULT),
100  &H5Fclose,
101  "Error message.");
102 
103  ... // use file_id in the same way as a plain hid_t object
104  \endcode
105 
106  <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br>
107  Namespace: vigra
108  **/
110 {
111 public:
112  typedef herr_t (*Destructor)(hid_t);
113 
114 private:
115  hid_t handle_;
116  Destructor destructor_;
117 
118 public:
119 
120  /** \brief Default constuctor.
121  Creates a NULL handle.
122  **/
124  : handle_( 0 ),
125  destructor_(0)
126  {}
127 
128  /** \brief Create a wrapper for a hid_t object.
129 
130  The hid_t object \a h is assumed to be the return value of an open or create function.
131  It will be closed with the given close function \a destructor as soon as this
132  HDF5Handle is destructed, except when \a destructor is a NULL pointer (in which
133  case nothing happens at destruction time). If \a h has a value that indicates
134  failed opening or creation (by HDF5 convention, this means if it is a negative number),
135  an exception is raised by calling <tt>vigra_fail(error_message)</tt>.
136 
137  <b>Usage:</b>
138 
139  \code
140  HDF5Handle file_id(H5Fopen(filename, H5F_ACC_RDWR, H5P_DEFAULT),
141  &H5Fclose,
142  "Error message.");
143 
144  ... // use file_id in the same way
145  \endcode
146  **/
147  HDF5Handle(hid_t h, Destructor destructor, const char * error_message)
148  : handle_( h ),
149  destructor_(destructor)
150  {
151  if(handle_ < 0)
152  vigra_fail(error_message);
153  }
154 
155  /** \brief Copy constructor.
156  Hands over ownership of the RHS handle (analogous to std::auto_ptr).
157  **/
159  : handle_( h.handle_ ),
160  destructor_(h.destructor_)
161  {
162  const_cast<HDF5Handle &>(h).handle_ = 0;
163  }
164 
165  /** \brief Assignment.
166  Calls close() for the LHS handle and hands over ownership of the
167  RHS handle (analogous to std::auto_ptr).
168  **/
170  {
171  if(h.handle_ != handle_)
172  {
173  close();
174  handle_ = h.handle_;
175  destructor_ = h.destructor_;
176  const_cast<HDF5Handle &>(h).handle_ = 0;
177  }
178  return *this;
179  }
180 
181  /** \brief Destreuctor.
182  Calls close() for the contained handle.
183  **/
185  {
186  close();
187  }
188 
189  /** \brief Explicitly call the stored function (if one has been stored within
190  this object) for the contained handle and set the handle to NULL.
191  **/
192  herr_t close()
193  {
194  herr_t res = 1;
195  if(handle_ && destructor_)
196  res = (*destructor_)(handle_);
197  handle_ = 0;
198  return res;
199  }
200 
201  /** \brief Get a temporary hid_t object for the contained handle.
202  Do not call a close function on the return value - a crash will be likely
203  otherwise.
204  **/
205  hid_t get() const
206  {
207  return handle_;
208  }
209 
210  /** \brief Convert to a plain hid_t object.
211 
212  This function ensures that hid_t objects can be transparently replaced with
213  HDF5Handle objects in user code. Do not call a close function on the return
214  value - a crash will be likely otherwise.
215  **/
216  operator hid_t() const
217  {
218  return handle_;
219  }
220 
221  /** \brief Equality comparison of the contained handle.
222  **/
223  bool operator==(HDF5Handle const & h) const
224  {
225  return handle_ == h.handle_;
226  }
227 
228  /** \brief Equality comparison of the contained handle.
229  **/
230  bool operator==(hid_t h) const
231  {
232  return handle_ == h;
233  }
234 
235  /** \brief Unequality comparison of the contained handle.
236  **/
237  bool operator!=(HDF5Handle const & h) const
238  {
239  return handle_ != h.handle_;
240  }
241 
242  /** \brief Unequality comparison of the contained handle.
243  **/
244  bool operator!=(hid_t h) const
245  {
246  return handle_ != h;
247  }
248 };
249 
250 
251 /********************************************************/
252 /* */
253 /* HDF5ImportInfo */
254 /* */
255 /********************************************************/
256 
257 /** \brief Argument object for the function readHDF5().
258 
259 See \ref readHDF5() for a usage example. This object must be
260 used to read an image or array from a HDF5 file
261 and enquire about its properties.
262 
263 <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br>
264 Namespace: vigra
265 **/
267 {
268  public:
269  enum PixelType { UINT8, UINT16, UINT32, UINT64,
270  INT8, INT16, INT32, INT64,
271  FLOAT, DOUBLE };
272 
273  /** Construct HDF5ImportInfo object.
274 
275  The dataset \a pathInFile in the HDF5 file \a filename is accessed to
276  read its properties. \a pathInFile may contain '/'-separated group
277  names, but must end with the name of the desired dataset:
278 
279  \code
280  HDF5ImportInfo info(filename, "/group1/group2/my_dataset");
281  \endcode
282  **/
283  VIGRA_EXPORT HDF5ImportInfo( const char* filePath, const char* pathInFile );
284 
285  VIGRA_EXPORT ~HDF5ImportInfo();
286 
287  /** Get the filename of this HDF5 object.
288  **/
289  VIGRA_EXPORT const std::string& getFilePath() const;
290 
291  /** Get the dataset's full name in the HDF5 file.
292  **/
293  VIGRA_EXPORT const std::string& getPathInFile() const;
294 
295  /** Get a handle to the file represented by this info object.
296  **/
297  VIGRA_EXPORT hid_t getH5FileHandle() const;
298 
299  /** Get a handle to the dataset represented by this info object.
300  **/
301  VIGRA_EXPORT hid_t getDatasetHandle() const;
302 
303  /** Get the number of dimensions of the dataset represented by this info object.
304  **/
305  VIGRA_EXPORT MultiArrayIndex numDimensions() const;
306 
307  /** Get the shape of the dataset represented by this info object.
308  **/
309  VIGRA_EXPORT ArrayVector<hsize_t> const & shape() const
310  {
311  return m_dims;
312  }
313 
314  /** Get the shape (length) of the dataset along dimension \a dim.
315  **/
316  VIGRA_EXPORT MultiArrayIndex shapeOfDimension(const int dim) const;
317 
318  /** Query the pixel type of the dataset.
319 
320  Possible values are:
321  <DL>
322  <DT>"UINT8"<DD> 8-bit unsigned integer (unsigned char)
323  <DT>"INT16"<DD> 16-bit signed integer (short)
324  <DT>"UINT16"<DD> 16-bit unsigned integer (unsigned short)
325  <DT>"INT32"<DD> 32-bit signed integer (long)
326  <DT>"UINT32"<DD> 32-bit unsigned integer (unsigned long)
327  <DT>"FLOAT"<DD> 32-bit floating point (float)
328  <DT>"DOUBLE"<DD> 64-bit floating point (double)
329  </DL>
330  **/
331  VIGRA_EXPORT const char * getPixelType() const;
332 
333  /** Query the pixel type of the dataset.
334 
335  Same as getPixelType(), but the result is returned as a
336  ImageImportInfo::PixelType enum. This is useful to implement
337  a switch() on the pixel type.
338 
339  Possible values are:
340  <DL>
341  <DT>UINT8<DD> 8-bit unsigned integer (unsigned char)
342  <DT>INT16<DD> 16-bit signed integer (short)
343  <DT>UINT16<DD> 16-bit unsigned integer (unsigned short)
344  <DT>INT32<DD> 32-bit signed integer (long)
345  <DT>UINT32<DD> 32-bit unsigned integer (unsigned long)
346  <DT>FLOAT<DD> 32-bit floating point (float)
347  <DT>DOUBLE<DD> 64-bit floating point (double)
348  </DL>
349  **/
350  VIGRA_EXPORT PixelType pixelType() const;
351 
352  private:
353  HDF5Handle m_file_handle, m_dataset_handle;
354  std::string m_filename, m_path, m_pixeltype;
355  hssize_t m_dimensions;
356  ArrayVector<hsize_t> m_dims;
357 };
358 
359 
360 namespace detail {
361 
362 template<class type>
363 inline hid_t getH5DataType()
364 {
365  std::runtime_error("getH5DataType(): invalid type");
366  return 0;
367 }
368 
369 #define VIGRA_H5_DATATYPE(type, h5type) \
370 template<> \
371 inline hid_t getH5DataType<type>() \
372 { return h5type;}
373 VIGRA_H5_DATATYPE(char, H5T_NATIVE_CHAR)
374 VIGRA_H5_DATATYPE(Int8, H5T_NATIVE_INT8)
375 VIGRA_H5_DATATYPE(Int16, H5T_NATIVE_INT16)
376 VIGRA_H5_DATATYPE(Int32, H5T_NATIVE_INT32)
377 VIGRA_H5_DATATYPE(Int64, H5T_NATIVE_INT64)
378 VIGRA_H5_DATATYPE(UInt8, H5T_NATIVE_UINT8)
379 VIGRA_H5_DATATYPE(UInt16, H5T_NATIVE_UINT16)
380 VIGRA_H5_DATATYPE(UInt32, H5T_NATIVE_UINT32)
381 VIGRA_H5_DATATYPE(UInt64, H5T_NATIVE_UINT64)
382 VIGRA_H5_DATATYPE(float, H5T_NATIVE_FLOAT)
383 VIGRA_H5_DATATYPE(double, H5T_NATIVE_DOUBLE)
384 VIGRA_H5_DATATYPE(long double, H5T_NATIVE_LDOUBLE)
385 
386 #undef VIGRA_H5_DATATYPE
387 
388 } // namespace detail
389 
390 
391 
392 
393 /********************************************************/
394 /* */
395 /* HDF5File */
396 /* */
397 /********************************************************/
398 
399 
400 /** \brief Access to HDF5 files
401 
402 HDF5File proviedes a convenient way of accessing data in HDF5 files. vigra::MultiArray
403 structures of any dimension can be stored to / loaded from HDF5 files. Typical
404 HDF5 features like subvolume access, chunks and data compression are available,
405 string attributes can be attached to any dataset or group. Group- or dataset-handles
406 are encapsulated in the class and managed automatically. The internal file-system like
407 structure can be accessed by functions like "cd()" or "mkdir()".
408 
409 
410 <b>Example:</b>
411 Write the MultiArray out_multi_array to file. Change the current directory to
412 "/group" and read in the same MultiArray as in_multi_array.
413 \code
414 HDF5File file("/path/to/file",HDF5File::New);
415 file.mkdir("group");
416 file.write("/group/dataset", out_multi_array);
417 
418 file.cd("/group");
419 file.read("dataset", in_multi_array);
420 
421 \endcode
422 
423 <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br>
424 Namespace: vigra
425 **/
426 class HDF5File
427 {
428  private:
429  HDF5Handle fileHandle_;
430 
431  // current group handle
432  HDF5Handle cGroupHandle_;
433 
434 
435  // datastructure to hold a list of dataset and group names
436  struct lsOpData
437  {
438  std::vector<std::string> objects;
439  };
440 
441  // operator function, used by H5Literate
442  static herr_t opFunc (hid_t loc_id, const char *name, const H5L_info_t *info, void *operator_data)
443  {
444  // get information about object
445  H5O_info_t infobuf;
446  H5Oget_info_by_name (loc_id, name, &infobuf, H5P_DEFAULT);
447 
448  // add name to list, if object is a dataset or a group
449  if(infobuf.type == H5O_TYPE_GROUP)
450  {
451  (*(struct lsOpData *) operator_data).objects.push_back(std::string(name)+"/");
452  }
453  if(infobuf.type == H5O_TYPE_DATASET)
454  {
455  (*(struct lsOpData *) operator_data).objects.push_back(std::string(name));
456  }
457 
458  return 0;
459  }
460 
461  public:
462  /** \brief Set how a file is opened.
463  OpenMode::New creates a new file. If the file already exists, overwrite it.
464 
465  OpenMode::Open opens a file for reading/writing. The file will be created,
466  if nescessary.
467  **/
468  enum OpenMode {
469  New, // Create new empty file (existing file will be deleted).
470  Open // Open file. Create if not existing.
471  };
472 
473 
474  /** \brief Create a HDF5File object.
475 
476  Creates or opens HDF5 file at position filename. The current group is set
477  to "/".
478  **/
479  HDF5File(std::string filename, OpenMode mode)
480  {
481  std::string errorMessage = "HDF5File: Could not create file '" + filename + "'.";
482  fileHandle_ = HDF5Handle(createFile_(filename, mode), &H5Fclose, errorMessage.c_str());
483  cGroupHandle_ = HDF5Handle(openCreateGroup_("/"), &H5Gclose, "HDF5File(): Failed to open root group.");
484  }
485 
486 
487  /** \brief Destructor to make sure that all data is flushed before closing the file.
488  */
490  {
491  //Write everything to disk before closing
492  H5Fflush(fileHandle_, H5F_SCOPE_GLOBAL);
493  }
494 
495 
496  /** \brief Change current group to "/".
497  */
498  void root()
499  {
500  std::string message = "HDF5File::root(): Could not open group '/'.";
501  cGroupHandle_ = HDF5Handle(H5Gopen(fileHandle_, "/", H5P_DEFAULT),&H5Gclose,message.c_str());
502  }
503 
504 
505  /** \brief Change the current group.
506  If the first character is a "/", the path will be interpreted as absolute path,
507  otherwise it will be interpreted as path relative to the current group.
508  */
509  void cd(std::string groupName)
510  {
511  std::string message = "HDF5File::cd(): Could not open group '" + groupName + "'.\n";
512  if(groupName == "/")
513  {
514  cGroupHandle_ = HDF5Handle(openCreateGroup_("/"),&H5Gclose,message.c_str());
515  return;
516  }
517  else if(groupName =="..")
518  {
519  cd_up();
520  }
521  else{
522  if(relativePath_(groupName))
523  {
524  if (H5Lexists(cGroupHandle_, groupName.c_str(), H5P_DEFAULT) == 0)
525  {
526  std::cerr << message;
527  return;
528  }
529  cGroupHandle_ = HDF5Handle(openCreateGroup_(groupName),&H5Gclose,message.c_str());
530  }
531  else
532  {
533  if (H5Lexists(fileHandle_, groupName.c_str(), H5P_DEFAULT) == 0)
534  {
535  std::cerr << message;
536  return;
537  }
538  cGroupHandle_ = HDF5Handle(openCreateGroup_(groupName),&H5Gclose,message.c_str());
539  }
540  }
541  }
542 
543  /** \brief Change the current group to its parent group.
544  returns true if successful, false otherwise.
545  */
546  bool cd_up()
547  {
548  std::string groupName = currentGroupName_();
549 
550  //do not try to move up if we already in "/"
551  if(groupName == "/"){
552  return false;
553  }
554 
555  size_t lastSlash = groupName.find_last_of('/');
556 
557  std::string parentGroup (groupName.begin(), groupName.begin()+lastSlash+1);
558 
559  cd(parentGroup);
560 
561  return true;
562  }
563  void cd_up(int levels)
564  {
565  for(int i = 0; i<levels; i++)
566  cd_up();
567  }
568 
569 
570  /** \brief Create a new group.
571  If the first character is a "/", the path will be interpreted as absolute path,
572  otherwise it will be interpreted as path relative to the current group.
573  */
574  void mkdir(std::string groupName)
575  {
576  hid_t handle = openCreateGroup_(groupName.c_str());
577  if (handle != cGroupHandle_){
578  H5Gclose(handle);
579  }
580  }
581 
582  /** \brief Change the current group; create it if nescessary.
583  If the first character is a "/", the path will be interpreted as absolute path,
584  otherwise it will be interpreted as path relative to the current group.
585  */
586  void cd_mk(std::string groupName)
587  {
588  std::string message = "HDF5File::cd_mk(): Could not create group '" + groupName + "'.";
589  cGroupHandle_ = HDF5Handle(openCreateGroup_(groupName.c_str()),&H5Gclose,message.c_str());
590  }
591 
592 
593 
594  /** \brief List the content of the current group.
595  The function returns a vector of strings holding the entries of the current
596  group. Only datasets and groups are listed, other objects (e.g. datatypes)
597  are ignored. Group names always have an ending "/".
598  */
599  std::vector<std::string> ls()
600  {
601  lsOpData data;
602  H5Literate(cGroupHandle_,H5_INDEX_NAME,H5_ITER_NATIVE,NULL, &opFunc, (void *) &data);
603 
604  return data.objects;
605  }
606 
607 
608  /** \brief Get the path of the current group.
609  */
610  std::string pwd()
611  {
612  return currentGroupName_();
613  }
614 
615  /** \brief Get the name of the associated file.
616  */
617  std::string filename()
618  {
619  return fileName_();
620  }
621 
622  /** \brief Get the number of dimensions of a certain dataset
623  If the first character is a "/", the path will be interpreted as absolute path,
624  otherwise it will be interpreted as path relative to the current group.
625  */
626  hssize_t getDatasetDimensions(std::string datasetName)
627  {
628  //Open dataset and dataspace
629  std::string errorMessage = "HDF5File::getDatasetDimensions(): Unable to open dataset '" + datasetName + "'.";
630  HDF5Handle datasetHandle = HDF5Handle(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
631 
632  errorMessage = "HDF5File::getDatasetDimensions(): Unable to access dataspace.";
633  HDF5Handle dataspaceHandle(H5Dget_space(datasetHandle), &H5Sclose, errorMessage.c_str());
634 
635  //return dimension information
636  return H5Sget_simple_extent_ndims(dataspaceHandle);
637  }
638 
639  /** \brief Get the shape of each dimension of a certain dataset.
640  Normally, this function is called after determining the dimension of the
641  dataset using \ref getDatasetDimensions().
642  If the first character is a "/", the path will be interpreted as absolute path,
643  otherwise it will be interpreted as path relative to the current group.
644  */
645  ArrayVector<hsize_t> getDatasetShape(std::string datasetName)
646  {
647  //Open dataset and dataspace
648  std::string errorMessage = "HDF5File::getDatasetShape(): Unable to open dataset '" + datasetName + "'.";
649  HDF5Handle datasetHandle = HDF5Handle(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
650 
651  errorMessage = "HDF5File::getDatasetShape(): Unable to access dataspace.";
652  HDF5Handle dataspaceHandle(H5Dget_space(datasetHandle), &H5Sclose, errorMessage.c_str());
653 
654  //get dimension information
655  ArrayVector<hsize_t>::size_type dimensions = H5Sget_simple_extent_ndims(dataspaceHandle);
656 
657  ArrayVector<hsize_t> shape(dimensions);
658  ArrayVector<hsize_t> maxdims(dimensions);
659  H5Sget_simple_extent_dims(dataspaceHandle, shape.data(), maxdims.data());
660 
661 
662  // invert the dimensions to guarantee c-order
663  ArrayVector<hsize_t> shape_inv(dimensions);
664  for(ArrayVector<hsize_t>::size_type i=0; i<shape.size(); i++) {
665  shape_inv[i] = shape[dimensions-1-i];
666  }
667 
668  return shape_inv;
669  }
670 
671  /** \brief Attach a string attribute to an existing object.
672  The attribute can be attached to datasets and groups. The string may have arbitrary length.
673  */
674  void setAttribute(std::string datasetName, std::string attributeName, std::string text)
675  {
676  std::string groupname = SplitString(datasetName).first();
677  std::string setname = SplitString(datasetName).last();
678 
679  std::string errorMessage ("HDF5File::setAttribute(): Unable to open group '" + groupname + "'.");
680  HDF5Handle groupHandle (openCreateGroup_(groupname), &H5Gclose, errorMessage.c_str());
681 
682  H5LTset_attribute_string(groupHandle,setname.c_str(), attributeName.c_str(),text.c_str());
683  }
684 
685 
686  /** \brief Get an attribute string of an object.
687  */
688  std::string getAttribute(std::string datasetName, std::string attributeName)
689  {
690  typedef ArrayVector<char>::size_type SizeType;
691  if (H5Lexists(fileHandle_, datasetName.c_str(), H5P_DEFAULT) == 0)
692  {
693  std::cerr << "HDF5File::getAttribute(): Dataset '" << datasetName << "' does not exist.\n";
694  return "error - dataset not found";
695  }
696 
697  std::string groupname = SplitString(datasetName).first();
698  std::string setname = SplitString(datasetName).last();
699 
700  std::string errorMessage ("HDF5File::setAttribute(): Unable to open group '" + groupname + "'.");
701  HDF5Handle groupHandle (openCreateGroup_(groupname), &H5Gclose, errorMessage.c_str());
702 
703  // get the size of the attribute
704  HDF5Handle AttrHandle (H5Aopen_by_name(groupHandle,setname.c_str(),attributeName.c_str(),H5P_DEFAULT, H5P_DEFAULT),&H5Aclose, "HDF5File::getAttribute(): Unable to open attribute.");
705  SizeType len = (SizeType)H5Aget_storage_size(AttrHandle);
706 
707  //read the attribute
708  ArrayVector<char> text (len+1,0);
709  H5LTget_attribute_string(groupHandle, setname.c_str(), attributeName.c_str(), text.begin());
710 
711  return std::string(text.begin());
712  }
713 
714 
715  // Writing data
716 
717  /** \brief Write multi arrays.
718  Chunks can be activated by setting \code iChunkSize = size; //size > 0 \endcode .
719  The chunks will be hypercubes with edge length size.
720 
721  Compression can be activated by setting \code compression = parameter; // 0 < parameter <= 9 \endcode
722  where 0 stands for no compression and 9 for maximum compression.
723 
724  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
725  otherwise it will be interpreted as path relative to the current group.
726  */
727  template<unsigned int N, class T>
728  inline void write(std::string datasetName, const MultiArrayView<N, T, UnstridedArrayTag> & array, int iChunkSize = 0, int compression = 0)
729  {
730  typename MultiArrayShape<N>::type chunkSize;
731  for(int i = 0; i < N; i++){
732  chunkSize[i] = iChunkSize;
733  }
734  write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize, compression);
735  }
736 
737  /** \brief Write multi arrays.
738  Chunks can be activated by providing a MultiArrayShape as chunkSize.
739  chunkSize must have equal dimension as array.
740 
741  Compression can be activated by setting \code compression = parameter; // 0 < parameter <= 9 \endcode
742  where 0 stands for no compression and 9 for maximum compression.
743 
744  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
745  otherwise it will be interpreted as path relative to the current group.
746  */
747  template<unsigned int N, class T>
748  inline void write(std::string datasetName, const MultiArrayView<N, T, UnstridedArrayTag> & array, typename MultiArrayShape<N>::type chunkSize, int compression = 0)
749  {
750  write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize, compression);
751  }
752 
753  /** \brief Write a multi array into a larger volume.
754  blockOffset determines the position, where array is written.
755 
756  Chunks can be activated by providing a MultiArrayShape as chunkSize.
757  chunkSize must have equal dimension as array.
758 
759  Compression can be activated by setting \code compression = parameter; // 0 < parameter <= 9 \endcode
760  where 0 stands for no compression and 9 for maximum compression.
761 
762  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
763  otherwise it will be interpreted as path relative to the current group.
764  */
765  template<unsigned int N, class T>
766  inline void writeBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, const MultiArrayView<N, T, UnstridedArrayTag> & array)
767  {
768  writeBlock_(datasetName, blockOffset, array, detail::getH5DataType<T>(), 1);
769  }
770 
771  // non-scalar (TinyVector) and unstrided multi arrays
772  template<unsigned int N, class T, int SIZE>
773  inline void write(std::string datasetName, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array, int iChunkSize = 0, int compression = 0)
774  {
775  typename MultiArrayShape<N>::type chunkSize;
776  for(int i = 0; i < N; i++){
777  chunkSize[i] = iChunkSize;
778  }
779  write_(datasetName, array, detail::getH5DataType<T>(), SIZE, chunkSize, compression);
780  }
781 
782  template<unsigned int N, class T, int SIZE>
783  inline void write(std::string datasetName, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array, typename MultiArrayShape<N>::type chunkSize, int compression = 0)
784  {
785  write_(datasetName, array, detail::getH5DataType<T>(), SIZE, chunkSize, compression);
786  }
787 
788  template<unsigned int N, class T, int SIZE>
789  inline void writeBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array)
790  {
791  writeBlock_(datasetName, blockOffset, array, detail::getH5DataType<T>(), SIZE);
792  }
793 
794  // non-scalar (RGBValue) and unstrided multi arrays
795  template<unsigned int N, class T>
796  inline void write(std::string datasetName, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array, int iChunkSize = 0, int compression = 0)
797  {
798  typename MultiArrayShape<N>::type chunkSize;
799  for(int i = 0; i < N; i++){
800  chunkSize[i] = iChunkSize;
801  }
802  write_(datasetName, array, detail::getH5DataType<T>(), 3, chunkSize, compression);
803  }
804 
805  template<unsigned int N, class T>
806  inline void write(std::string datasetName, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array, typename MultiArrayShape<N>::type chunkSize, int compression = 0)
807  {
808  write_(datasetName, array, detail::getH5DataType<T>(), 3, chunkSize, compression);
809  }
810 
811  template<unsigned int N, class T>
812  inline void writeBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array)
813  {
814  writeBlock_(datasetName, blockOffset, array, detail::getH5DataType<T>(), 3);
815  }
816 
817  /** \brief Write single value as dataset.
818  This functions allows to write data of atomic datatypes (int, long, double)
819  as a dataset in the HDF5 file. So it is not nescessary to create a MultiArray
820  of size 1 to write a single number.
821 
822  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
823  otherwise it will be interpreted as path relative to the current group.
824  */
825  template<class T>
826  inline void writeAtomic(std::string datasetName, const T data)
827  {
828  typename MultiArrayShape<1>::type chunkSize;
829  chunkSize[0] = 0;
831  array[0] = data;
832  write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize,0);
833  }
834 
835 
836  /* Reading data. */
837 
838  /** \brief Read data into a multi array.
839  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
840  otherwise it will be interpreted as path relative to the current group.
841  */
842  template<unsigned int N, class T>
843  inline void read(std::string datasetName, MultiArrayView<N, T, UnstridedArrayTag> array)
844  {
845  read_(datasetName, array, detail::getH5DataType<T>(), 1);
846  }
847 
848  /** \brief Read a block of data into s multi array.
849  This function allows to read a small block out of a larger volume stored
850  in a HDF5 dataset.
851 
852  blockOffset determines the position of the block.
853  blockSize determines the size in each dimension of the block.
854 
855  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
856  otherwise it will be interpreted as path relative to the current group.
857  */
858  template<unsigned int N, class T>
859  inline void readBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, typename MultiArrayShape<N>::type blockShape, MultiArrayView<N, T, UnstridedArrayTag> array)
860  {
861  readBlock_(datasetName, blockOffset, blockShape, array, detail::getH5DataType<T>(), 1);
862  }
863 
864  // non-scalar (TinyVector) and unstrided target multi array
865  template<unsigned int N, class T, int SIZE>
866  inline void read(std::string datasetName, MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> array)
867  {
868  read_(datasetName, array, detail::getH5DataType<T>(), SIZE);
869  }
870 
871  template<unsigned int N, class T, int SIZE>
872  inline void readBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, typename MultiArrayShape<N>::type blockShape, MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> array)
873  {
874  readBlock_(datasetName, blockOffset, blockShape, array, detail::getH5DataType<T>(), SIZE);
875  }
876 
877  // non-scalar (RGBValue) and unstrided target multi array
878  template<unsigned int N, class T>
879  inline void read(std::string datasetName, MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> array)
880  {
881  read_(datasetName, array, detail::getH5DataType<T>(), 3);
882  }
883 
884  template<unsigned int N, class T>
885  inline void readBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, typename MultiArrayShape<N>::type blockShape, MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> array)
886  {
887  readBlock_(datasetName, blockOffset, blockShape, array, detail::getH5DataType<T>(), 3);
888  }
889 
890  /** \brief Read a single value.
891  This functions allows to read a single datum of atomic datatype (int, long, double)
892  from the HDF5 file. So it is not nescessary to create a MultiArray
893  of size 1 to read a single number.
894 
895  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
896  otherwise it will be interpreted as path relative to the current group.
897  */
898  template<class T>
899  inline void readAtomic(std::string datasetName, T & data)
900  {
902  read_(datasetName, array, detail::getH5DataType<T>(), 1);
903  data = array[0];
904  }
905 
906 
907  /** \brief Create a new dataset.
908  This function can be used to create a dataset filled with a default value,
909  for example before writing data into it using \ref writeBlock().
910  Attention: only atomic datatypes are provided. For spectral data, add an
911  dimension (case RGB: add one dimension of size 3).
912 
913  shape determines the dimension and the size of the dataset.
914 
915  Chunks can be activated by providing a MultiArrayShape as chunkSize.
916  chunkSize must have equal dimension as array.
917 
918  Compression can be activated by setting \code compression = parameter; // 0 < parameter <= 9 \endcode
919  where 0 stands for no compression and 9 for maximum compression.
920 
921  If the first character of datasetName is a "/", the path will be interpreted as absolute path,
922  otherwise it will be interpreted as path relative to the current group.
923  */
924  template<unsigned int N, class T>
925  inline void createDataset(std::string datasetName, typename MultiArrayShape<N>::type shape, T init, int iChunkSize = 0, int compressionParameter = 0)
926  {
927  typename MultiArrayShape<N>::type chunkSize;
928  for(int i = 0; i < N; i++){
929  chunkSize[i] = iChunkSize;
930  }
931  createDataset<N,T>(datasetName, shape, init, chunkSize, compressionParameter);
932  }
933 
934  template<unsigned int N, class T>
935  inline void createDataset(std::string datasetName, typename MultiArrayShape<N>::type shape, T init, typename MultiArrayShape<N>::type chunkSize, int compressionParameter = 0)
936  {
937  std::string groupname = SplitString(datasetName).first();
938  std::string setname = SplitString(datasetName).last();
939 
940  hid_t parent = openCreateGroup_(groupname);
941 
942  // delete the dataset if it already exists
943  deleteDataset_(parent, setname);
944 
945  // create dataspace
946  // add an extra dimension in case that the data is non-scalar
947  HDF5Handle dataspaceHandle;
948 
949  // invert dimensions to guarantee c-order
950  hsize_t shape_inv[N];
951  for(unsigned int k=0; k<N; ++k)
952  shape_inv[N-1-k] = shape[k];
953 
954  // create dataspace
955  dataspaceHandle = HDF5Handle(H5Screate_simple(N, shape_inv, NULL),
956  &H5Sclose, "HDF5File::createDataset(): unable to create dataspace for scalar data.");
957 
958  // set fill value
959  HDF5Handle plist ( H5Pcreate(H5P_DATASET_CREATE), &H5Pclose, "HDF5File::createDataset(): unable to create property list." );
960  H5Pset_fill_value(plist,detail::getH5DataType<T>(), &init);
961 
962  // enable chunks
963  if(chunkSize[0] > 0)
964  {
965  hsize_t cSize [N];
966  for(int i = 0; i<N; i++)
967  {
968  cSize[i] = chunkSize[N-1-i];
969  }
970  H5Pset_chunk (plist, N, cSize);
971  }
972 
973  // enable compression
974  if(compressionParameter > 0)
975  {
976  H5Pset_deflate(plist, compressionParameter);
977  }
978 
979  //create the dataset.
980  HDF5Handle datasetHandle ( H5Dcreate(parent, setname.c_str(), detail::getH5DataType<T>(), dataspaceHandle, H5P_DEFAULT, plist, H5P_DEFAULT),
981  &H5Dclose, "HDF5File::createDataset(): unable to create dataset.");
982  if(parent != cGroupHandle_)
983  H5Gclose(parent);
984  }
985 
986  /** \brief Immediately write all data to disk
987  */
988  void flushToDisk()
989  {
990  H5Fflush(fileHandle_, H5F_SCOPE_GLOBAL);
991  }
992 
993 
994  private:
995 
996  /* Simple extension of std::string for splitting into two parts
997  *
998  * Strings (in particular: file/dataset paths) will be split into two
999  * parts. The split is made at the last occurance of the delimiter.
1000  *
1001  * For example, "/path/to/some/file" will be split (delimiter = "/") into
1002  * first() = "/path/to/some" and last() = "file".
1003  */
1004  class SplitString: public std::string {
1005  public:
1006  SplitString(std::string &sstring): std::string(sstring) {};
1007 
1008  // return the part of the string before the delimiter
1009  std::string first(char delimiter = '/')
1010  {
1011  size_t last = find_last_of(delimiter);
1012  if(last == std::string::npos) // delimiter not found --> no first
1013  return "";
1014 
1015  return std::string(begin(), begin()+last+1);
1016  }
1017 
1018  // return the part of the string after the delimiter
1019  std::string last(char delimiter = '/')
1020  {
1021  size_t last = find_last_of(delimiter);
1022  if(last == std::string::npos) // delimiter not found --> only last
1023  return std::string(*this);
1024  return std::string(begin()+last+1, end());
1025  }
1026  };
1027 
1028  inline bool relativePath_(std::string & path)
1029  {
1030  std::string::size_type pos = path.find('/') ;
1031  if(pos == 0)
1032  return false;
1033 
1034  return true;
1035  }
1036 
1037 
1038  std::string currentGroupName_()
1039  {
1040  int len = H5Iget_name(cGroupHandle_,NULL,1000);
1041  ArrayVector<char> name (len+1,0);
1042  H5Iget_name(cGroupHandle_,name.begin(),len+1);
1043 
1044  return std::string(name.begin());
1045  }
1046 
1047  std::string fileName_()
1048  {
1049  int len = H5Fget_name(fileHandle_,NULL,1000);
1050  ArrayVector<char> name (len+1,0);
1051  H5Fget_name(fileHandle_,name.begin(),len+1);
1052 
1053  return std::string(name.begin());
1054  }
1055 
1056 
1057  inline hid_t createFile_(std::string filePath, OpenMode mode = Open)
1058  {
1059  // try to open file
1060  FILE * pFile;
1061  pFile = fopen ( filePath.c_str(), "r" );
1062  hid_t fileId;
1063 
1064  // check if opening was successful (= file exists)
1065  if ( pFile == NULL )
1066  {
1067  fileId = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
1068  }
1069  else if(mode == Open)
1070  {
1071  fclose( pFile );
1072  fileId = H5Fopen(filePath.c_str(), H5F_ACC_RDWR, H5P_DEFAULT);
1073  }
1074  else
1075  {
1076  fclose(pFile);
1077  std::remove(filePath.c_str());
1078  fileId = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
1079  }
1080  return fileId;
1081  }
1082 
1083 
1084  // open a group and subgroups. Create if nescessary.
1085  inline hid_t openCreateGroup_(std::string groupName)
1086  {
1087  // check if groupName is absolute or relative
1088  hid_t parent;
1089  if(relativePath_(groupName))
1090  {
1091  parent = cGroupHandle_;
1092  }
1093  else
1094  {
1095  // open root group
1096  parent = H5Gopen(fileHandle_, "/", H5P_DEFAULT);
1097  if(groupName == "/")
1098  {
1099  return parent;
1100  }
1101 
1102  // remove leading /
1103  groupName = std::string(groupName.begin()+1, groupName.end());
1104  }
1105 
1106 
1107  // check if the groupName has the right format
1108  if( groupName.size() != 0 && *groupName.rbegin() != '/')
1109  {
1110  groupName = groupName + '/';
1111  }
1112 
1113 
1114  //create subgroups one by one
1115  std::string::size_type begin = 0, end = groupName.find('/');
1116  int ii = 0;
1117  while (end != std::string::npos)
1118  {
1119  std::string group(groupName.begin()+begin, groupName.begin()+end);
1120  hid_t prevParent = parent;
1121 
1122  if(H5LTfind_dataset(parent, group.c_str()) == 0)
1123  {
1124  parent = H5Gcreate(prevParent, group.c_str(), H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
1125  } else {
1126  parent = H5Gopen(prevParent, group.c_str(), H5P_DEFAULT);
1127  }
1128 
1129  if(ii != 0)
1130  {
1131  H5Gclose(prevParent);
1132  }
1133  if(parent < 0)
1134  {
1135  return parent;
1136  }
1137  ++ii;
1138  begin = end + 1;
1139  end = groupName.find('/', begin);
1140  }
1141 
1142 
1143  return parent;
1144  }
1145 
1146 
1147  inline void deleteDataset_(hid_t parent, std::string datasetName)
1148  {
1149  // delete existing data and create new dataset
1150  if(H5LTfind_dataset(parent, datasetName.c_str()))
1151  {
1152  //std::cout << "dataset already exists" << std::endl;
1153  #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6)
1154  if(H5Gunlink(parent, datasetName.c_str()) < 0)
1155  {
1156  vigra_postcondition(false, "HDF5File::deleteDataset_(): Unable to delete existing data.");
1157  }
1158  #else
1159  if(H5Ldelete(parent, datasetName.c_str(), H5P_DEFAULT ) < 0)
1160  {
1161  vigra_postcondition(false, "HDF5File::deleteDataset_(): Unable to delete existing data.");
1162  }
1163  #endif
1164  }
1165  }
1166 
1167  hid_t getDatasetHandle_(std::string datasetName)
1168  {
1169  std::string groupname = SplitString(datasetName).first();
1170  std::string setname = SplitString(datasetName).last();
1171 
1172  if(relativePath_(datasetName))
1173  {
1174  if (H5Lexists(cGroupHandle_, datasetName.c_str(), H5P_DEFAULT) != 1)
1175  {
1176  std::cerr << "HDF5File::getDatasetHandle_(): Dataset '" << datasetName << "' does not exist.\n";
1177  return -1;
1178  }
1179  //Open parent group
1180  hid_t groupHandle = openCreateGroup_(groupname);
1181 
1182  hid_t datasetHandle = H5Dopen(groupHandle, setname.c_str(), H5P_DEFAULT);
1183 
1184  if(groupHandle != cGroupHandle_){
1185  H5Gclose(groupHandle);
1186  }
1187 
1188  //return dataset handle
1189  return datasetHandle;
1190  }
1191  else
1192  {
1193  if (H5Lexists(fileHandle_, datasetName.c_str(), H5P_DEFAULT) != 1)
1194  {
1195  std::cerr << "HDF5File::getDatasetHandle_(): Dataset '" << datasetName << "' does not exist.\n";
1196  return -1;
1197  }
1198 
1199  //Open parent group
1200  hid_t groupHandle = openCreateGroup_(groupname);
1201 
1202  hid_t datasetHandle = H5Dopen(groupHandle, setname.c_str(), H5P_DEFAULT);
1203 
1204  if(groupHandle != cGroupHandle_){
1205  H5Gclose(groupHandle);
1206  }
1207 
1208  //return dataset handle
1209  return datasetHandle;
1210  }
1211  }
1212 
1213 
1214  // unstrided multi arrays
1215  template<unsigned int N, class T>
1216  void write_(std::string &datasetName, const MultiArrayView<N, T, UnstridedArrayTag> & array, const hid_t datatype, const int numBandsOfType, typename MultiArrayShape<N>::type &chunkSize, int compressionParameter = 0)
1217  {
1218  std::string groupname = SplitString(datasetName).first();
1219  std::string setname = SplitString(datasetName).last();
1220 
1221  // shape of the array. Add one dimension, if array contains non-scalars.
1222  ArrayVector<hsize_t> shape(N + (numBandsOfType > 1),0);
1223  for(int i = 0; i < N; i++){
1224  shape[N-1-i] = array.shape(i); // reverse order
1225  }
1226 
1227  if(numBandsOfType > 1)
1228  shape[N] = numBandsOfType;
1229 
1230  HDF5Handle dataspace ( H5Screate_simple(N + (numBandsOfType > 1), shape.begin(), NULL), &H5Sclose, "HDF5File::write(): Can not create dataspace.");
1231 
1232  // create and open group:
1233  std::string errorMessage ("HDF5File::write(): can not create group '" + groupname + "'.");
1234  hid_t groupHandle = openCreateGroup_(groupname);
1235  if(groupHandle <= 0)
1236  {
1237  std::cerr << errorMessage << "\n";
1238  }
1239 
1240  // delete dataset, if it already exists
1241  deleteDataset_(groupHandle, setname.c_str());
1242 
1243  // set up properties list
1244  HDF5Handle plist ( H5Pcreate(H5P_DATASET_CREATE), &H5Pclose, "HDF5File::write(): unable to create property list." );
1245 
1246  // enable chunks
1247  if(chunkSize[0] > 0)
1248  {
1249  ArrayVector<hsize_t> cSize(N + (numBandsOfType > 1),0);
1250  for(int i = 0; i<N; i++)
1251  {
1252  cSize[i] = chunkSize[N-1-i];
1253  }
1254  if(numBandsOfType > 1)
1255  cSize[N] = numBandsOfType;
1256 
1257  H5Pset_chunk (plist, N + (numBandsOfType > 1), cSize.begin());
1258  }
1259 
1260  // enable compression
1261  if(compressionParameter > 0)
1262  {
1263  H5Pset_deflate(plist, compressionParameter);
1264  }
1265 
1266  // create dataset
1267  HDF5Handle datasetHandle (H5Dcreate(groupHandle, setname.c_str(), datatype, dataspace,H5P_DEFAULT, plist, H5P_DEFAULT), &H5Dclose, "HDF5File::write(): Can not create dataset.");
1268 
1269  // Write the data to the HDF5 dataset as is
1270  H5Dwrite( datasetHandle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data());
1271 
1272  if(groupHandle != cGroupHandle_)
1273  {
1274  H5Gclose(groupHandle);
1275  }
1276  }
1277 
1278  // unstrided target multi array
1279  template<unsigned int N, class T>
1280  void read_(std::string datasetName, MultiArrayView<N, T, UnstridedArrayTag> array, const hid_t datatype, const int numBandsOfType)
1281  {
1282  //Prepare to read without using HDF5ImportInfo
1283  ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName);
1284  hssize_t dimensions = getDatasetDimensions(datasetName);
1285 
1286  std::string errorMessage ("HDF5File::read(): Unable to open dataset '" + datasetName + "'.");
1287  HDF5Handle datasetHandle (getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
1288 
1289  int offset = (numBandsOfType > 1);
1290 
1291  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
1292  "HDF5File::read(): Array dimension disagrees with dataset dimension.");
1293 
1294  typename MultiArrayShape<N>::type shape;
1295  for(int k=offset; k< MultiArrayIndex(dimensions); ++k) {
1296  shape[k-offset] = MultiArrayIndex(dimshape[k]);
1297  }
1298 
1299  vigra_precondition(shape == array.shape(),
1300  "HDF5File::read(): Array shape disagrees with dataset shape.");
1301 
1302  // simply read in the data as is
1303  H5Dread( datasetHandle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data() ); // .data() possible since void pointer!
1304  }
1305 
1306  template<unsigned int N, class T>
1307  void writeBlock_(std::string datasetName, typename MultiArrayShape<N>::type &blockOffset, const MultiArrayView<N, T, UnstridedArrayTag> & array, const hid_t datatype, const int numBandsOfType)
1308  {
1309  // open dataset if it exists
1310  std::string errorMessage = "HDF5File::writeBlock(): Error opening dataset '" + datasetName + "'.";
1311  HDF5Handle datasetHandle (getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
1312 
1313 
1314  // hyperslab parameters for position, size, ...
1315  hsize_t boffset [N];
1316  hsize_t bshape [N];
1317  hsize_t bones [N];
1318 
1319  for(int i = 0; i < N; i++){
1320  boffset[i] = blockOffset[N-1-i];
1321  bshape[i] = array.size(N-1-i);
1322  bones[i] = 1;
1323  }
1324 
1325  // create a target dataspace in memory with the shape of the desired block
1326  HDF5Handle memspace_handle (H5Screate_simple(N,bshape,NULL),&H5Sclose,"Unable to get origin dataspace");
1327 
1328  // get file dataspace and select the desired block
1329  HDF5Handle dataspaceHandle (H5Dget_space(datasetHandle),&H5Sclose,"Unable to create target dataspace");
1330  H5Sselect_hyperslab(dataspaceHandle, H5S_SELECT_SET, boffset, bones, bones, bshape);
1331 
1332  // Write the data to the HDF5 dataset as is
1333  H5Dwrite( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, array.data()); // .data() possible since void pointer!
1334  }
1335 
1336 
1337 
1338  // unstrided target multi array
1339  // read a block of a HDF5 dataset into a MultiArray
1340  template<unsigned int N, class T>
1341  void readBlock_(std::string datasetName, typename MultiArrayShape<N>::type &blockOffset, typename MultiArrayShape<N>::type &blockShape, MultiArrayView<N, T, UnstridedArrayTag> &array, const hid_t datatype, const int numBandsOfType)
1342  {
1343  //Prepare to read without using HDF5ImportInfo
1344  //ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName) ;
1345  hssize_t dimensions = getDatasetDimensions(datasetName);
1346 
1347  std::string errorMessage ("HDF5File::readBlock(): Unable to open dataset '" + datasetName + "'.");
1348  HDF5Handle datasetHandle (getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
1349 
1350 
1351  int offset = (numBandsOfType > 1);
1352 
1353  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
1354  "readHDF5_block(): Array dimension disagrees with data dimension.");
1355 
1356  vigra_precondition(blockShape == array.shape(),
1357  "readHDF5_block(): Array shape disagrees with block size.");
1358 
1359  // hyperslab parameters for position, size, ...
1360  hsize_t boffset [N];
1361  hsize_t bshape [N];
1362  hsize_t bones [N];
1363 
1364  for(int i = 0; i < N; i++){
1365  // virgra and hdf5 use different indexing
1366  boffset[i] = blockOffset[N-1-i];
1367  //bshape[i] = blockShape[i];
1368  bshape[i] = blockShape[N-1-i];
1369  //boffset[i] = blockOffset[N-1-i];
1370  bones[i] = 1;
1371  }
1372 
1373  // create a target dataspace in memory with the shape of the desired block
1374  HDF5Handle memspace_handle (H5Screate_simple(N,bshape,NULL),&H5Sclose,"Unable to create target dataspace");
1375 
1376  // get file dataspace and select the desired block
1377  HDF5Handle dataspaceHandle (H5Dget_space(datasetHandle),&H5Sclose,"Unable to get dataspace");
1378  H5Sselect_hyperslab(dataspaceHandle, H5S_SELECT_SET, boffset, bones, bones, bshape);
1379 
1380  // now read the data
1381  H5Dread( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, array.data() ); // .data() possible since void pointer!
1382  }
1383 
1384 }; /* class HDF5File */
1385 
1386 
1387 
1388 
1389 
1390 
1391 
1392 namespace detail {
1393 
1394 template <class Shape>
1395 inline void
1396 selectHyperslabs(HDF5Handle & mid1, HDF5Handle & mid2, Shape const & shape, int & counter, const int elements, const int numBandsOfType)
1397 {
1398  // select hyperslab in HDF5 file
1399  hsize_t shapeHDF5[2];
1400  shapeHDF5[0] = 1;
1401  shapeHDF5[1] = elements;
1402  hsize_t startHDF5[2];
1403  startHDF5[0] = 0;
1404  startHDF5[1] = counter * numBandsOfType * shape[0]; // we have to reserve space for the pixel type channel(s)
1405  hsize_t strideHDF5[2];
1406  strideHDF5[0] = 1;
1407  strideHDF5[1] = 1;
1408  hsize_t countHDF5[2];
1409  countHDF5[0] = 1;
1410  countHDF5[1] = numBandsOfType * shape[0];
1411  hsize_t blockHDF5[2];
1412  blockHDF5[0] = 1;
1413  blockHDF5[1] = 1;
1414  mid1 = HDF5Handle(H5Screate_simple(2, shapeHDF5, NULL),
1415  &H5Sclose, "unable to create hyperslabs.");
1416  H5Sselect_hyperslab(mid1, H5S_SELECT_SET, startHDF5, strideHDF5, countHDF5, blockHDF5);
1417  // select hyperslab in input data object
1418  hsize_t shapeData[2];
1419  shapeData[0] = 1;
1420  shapeData[1] = numBandsOfType * shape[0];
1421  hsize_t startData[2];
1422  startData[0] = 0;
1423  startData[1] = 0;
1424  hsize_t strideData[2];
1425  strideData[0] = 1;
1426  strideData[1] = 1;
1427  hsize_t countData[2];
1428  countData[0] = 1;
1429  countData[1] = numBandsOfType * shape[0];
1430  hsize_t blockData[2];
1431  blockData[0] = 1;
1432  blockData[1] = 1;
1433  mid2 = HDF5Handle(H5Screate_simple(2, shapeData, NULL),
1434  &H5Sclose, "unable to create hyperslabs.");
1435  H5Sselect_hyperslab(mid2, H5S_SELECT_SET, startData, strideData, countData, blockData);
1436 }
1437 
1438 template <class DestIterator, class Shape, class T>
1439 inline void
1440 readHDF5Impl(DestIterator d, Shape const & shape, const hid_t dataset_id, const hid_t datatype, ArrayVector<T> & buffer, int & counter, const int elements, const int numBandsOfType, MetaInt<0>)
1441 {
1442  HDF5Handle mid1, mid2;
1443 
1444  // select hyperslabs
1445  selectHyperslabs(mid1, mid2, shape, counter, elements, numBandsOfType);
1446 
1447  // read from hdf5
1448  H5Dread(dataset_id, datatype, mid2, mid1, H5P_DEFAULT, buffer.data());
1449 
1450  // increase counter
1451  counter++;
1452 
1453 
1454  //std::cout << "numBandsOfType: " << numBandsOfType << std::endl;
1455  DestIterator dend = d + shape[0];
1456  int k = 0;
1457  for(; d < dend; ++d, k++)
1458  {
1459  *d = buffer[k];
1460  //std::cout << buffer[k] << "| ";
1461  }
1462 
1463 }
1464 
1465 template <class DestIterator, class Shape, class T, int N>
1466 void
1467 readHDF5Impl(DestIterator d, Shape const & shape, const hid_t dataset_id, const hid_t datatype, ArrayVector<T> & buffer, int & counter, const int elements, const int numBandsOfType, MetaInt<N>)
1468 {
1469  DestIterator dend = d + shape[N];
1470  for(; d < dend; ++d)
1471  {
1472  readHDF5Impl(d.begin(), shape, dataset_id, datatype, buffer, counter, elements, numBandsOfType, MetaInt<N-1>());
1473  }
1474 }
1475 
1476 } // namespace detail
1477 
1478  /** \brief Read the data specified by the given \ref vigra::HDF5ImportInfo object
1479  and write the into the given 'array'.
1480 
1481  The array must have the correct number of dimensions and shape for the dataset
1482  represented by 'info'. When the element type of 'array' differs from the stored element
1483  type, HDF5 will convert the type on the fly (except when the HDF5 version is 1.6 or below,
1484  in which case an error will result). Multi-channel element types (i.e. \ref vigra::RGBValue
1485  and \ref vigra::TinyVector) are recognized and handled correctly.
1486 
1487  <b> Declaration:</b>
1488 
1489  \code
1490  namespace vigra {
1491  template<unsigned int N, class T, class StrideTag>
1492  void
1493  readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, StrideTag> array);
1494  }
1495  \endcode
1496 
1497  <b> Usage:</b>
1498 
1499  <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br>
1500  Namespace: vigra
1501 
1502  \code
1503 
1504  HDF5ImportInfo info(filename, dataset_name);
1505  vigra_precondition(info.numDimensions() == 3, "Dataset must be 3-dimensional.");
1506 
1507  MultiArrayShape<3>::type shape(info.shape().begin());
1508  MultiArray<3, int> array(shape);
1509 
1510  readHDF5(info, array);
1511  \endcode
1512 */
1513 doxygen_overloaded_function(template <...> void readHDF5)
1514 
1515 // scalar and unstrided target multi array
1516 template<unsigned int N, class T>
1517 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, UnstridedArrayTag> array) // scalar
1518 {
1519  readHDF5(info, array, detail::getH5DataType<T>(), 1);
1520 }
1521 
1522 // non-scalar (TinyVector) and unstrided target multi array
1523 template<unsigned int N, class T, int SIZE>
1524 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> array)
1525 {
1526  readHDF5(info, array, detail::getH5DataType<T>(), SIZE);
1527 }
1528 
1529 // non-scalar (RGBValue) and unstrided target multi array
1530 template<unsigned int N, class T>
1531 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> array)
1532 {
1533  readHDF5(info, array, detail::getH5DataType<T>(), 3);
1534 }
1535 
1536 // unstrided target multi array
1537 template<unsigned int N, class T>
1538 void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, UnstridedArrayTag> array, const hid_t datatype, const int numBandsOfType)
1539 {
1540  int offset = (numBandsOfType > 1);
1541 
1542  //std::cout << "offset: " << offset << ", N: " << N << ", dims: " << info.numDimensions() << std::endl;
1543  vigra_precondition(( (N + offset ) == info.numDimensions()), // the object in the HDF5 file may have one additional dimension which we then interpret as the pixel type bands
1544  "readHDF5(): Array dimension disagrees with HDF5ImportInfo.numDimensions().");
1545 
1546  typename MultiArrayShape<N>::type shape;
1547  for(int k=offset; k<info.numDimensions(); ++k) {
1548  shape[k-offset] = info.shapeOfDimension(k);
1549  }
1550 
1551  vigra_precondition(shape == array.shape(),
1552  "readHDF5(): Array shape disagrees with HDF5ImportInfo.");
1553 
1554  // simply read in the data as is
1555  H5Dread( info.getDatasetHandle(), datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data() ); // .data() possible since void pointer!
1556 }
1557 
1558 // scalar and strided target multi array
1559 template<unsigned int N, class T>
1560 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, StridedArrayTag> array) // scalar
1561 {
1562  readHDF5(info, array, detail::getH5DataType<T>(), 1);
1563 }
1564 
1565 // non-scalar (TinyVector) and strided target multi array
1566 template<unsigned int N, class T, int SIZE>
1567 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, TinyVector<T, SIZE>, StridedArrayTag> array)
1568 {
1569  readHDF5(info, array, detail::getH5DataType<T>(), SIZE);
1570 }
1571 
1572 // non-scalar (RGBValue) and strided target multi array
1573 template<unsigned int N, class T>
1574 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, RGBValue<T>, StridedArrayTag> array)
1575 {
1576  readHDF5(info, array, detail::getH5DataType<T>(), 3);
1577 }
1578 
1579 // strided target multi array
1580 template<unsigned int N, class T>
1581 void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, StridedArrayTag> array, const hid_t datatype, const int numBandsOfType)
1582 {
1583  int offset = (numBandsOfType > 1);
1584 
1585  //std::cout << "offset: " << offset << ", N: " << N << ", dims: " << info.numDimensions() << std::endl;
1586  vigra_precondition(( (N + offset ) == info.numDimensions()), // the object in the HDF5 file may have one additional dimension which we then interpret as the pixel type bands
1587  "readHDF5(): Array dimension disagrees with HDF5ImportInfo.numDimensions().");
1588 
1589  typename MultiArrayShape<N>::type shape;
1590  for(int k=offset; k<info.numDimensions(); ++k) {
1591  shape[k-offset] = info.shapeOfDimension(k);
1592  }
1593 
1594  vigra_precondition(shape == array.shape(),
1595  "readHDF5(): Array shape disagrees with HDF5ImportInfo.");
1596 
1597  //Get the data
1598  int counter = 0;
1599  int elements = numBandsOfType;
1600  for(unsigned int i=0;i<N;++i)
1601  elements *= shape[i];
1602  ArrayVector<T> buffer(shape[0]);
1603  detail::readHDF5Impl(array.traverser_begin(), shape, info.getDatasetHandle(), datatype, buffer, counter, elements, numBandsOfType, vigra::MetaInt<N-1>());
1604 }
1605 
1606 inline hid_t openGroup(hid_t parent, std::string group_name)
1607 {
1608  //std::cout << group_name << std::endl;
1609  size_t last_slash = group_name.find_last_of('/');
1610  if (last_slash == std::string::npos || last_slash != group_name.size() - 1)
1611  group_name = group_name + '/';
1612  std::string::size_type begin = 0, end = group_name.find('/');
1613  int ii = 0;
1614  while (end != std::string::npos)
1615  {
1616  std::string group(group_name.begin()+begin, group_name.begin()+end);
1617  hid_t prev_parent = parent;
1618  parent = H5Gopen(prev_parent, group.c_str(), H5P_DEFAULT);
1619 
1620  if(ii != 0) H5Gclose(prev_parent);
1621  if(parent < 0) return parent;
1622  ++ii;
1623  begin = end + 1;
1624  end = group_name.find('/', begin);
1625  }
1626  return parent;
1627 }
1628 
1629 inline hid_t createGroup(hid_t parent, std::string group_name)
1630 {
1631  if(group_name.size() == 0 ||*group_name.rbegin() != '/')
1632  group_name = group_name + '/';
1633  if(group_name == "/")
1634  return H5Gopen(parent, group_name.c_str(), H5P_DEFAULT);
1635 
1636  std::string::size_type begin = 0, end = group_name.find('/');
1637  int ii = 0;
1638  while (end != std::string::npos)
1639  {
1640  std::string group(group_name.begin()+begin, group_name.begin()+end);
1641  hid_t prev_parent = parent;
1642 
1643  if(H5LTfind_dataset(parent, group.c_str()) == 0)
1644  {
1645  parent = H5Gcreate(prev_parent, group.c_str(), H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
1646  } else {
1647  parent = H5Gopen(prev_parent, group.c_str(), H5P_DEFAULT);
1648  }
1649 
1650  if(ii != 0) H5Gclose(prev_parent);
1651  if(parent < 0) return parent;
1652  ++ii;
1653  begin = end + 1;
1654  end = group_name.find('/', begin);
1655  }
1656  return parent;
1657 }
1658 
1659 inline void deleteDataset(hid_t parent, std::string dataset_name)
1660 {
1661  // delete existing data and create new dataset
1662  if(H5LTfind_dataset(parent, dataset_name.c_str()))
1663  {
1664  //std::cout << "dataset already exists" << std::endl;
1665 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6)
1666  if(H5Gunlink(parent, dataset_name.c_str()) < 0)
1667  {
1668  vigra_postcondition(false, "writeToHDF5File(): Unable to delete existing data.");
1669  }
1670 #else
1671  if(H5Ldelete(parent, dataset_name.c_str(), H5P_DEFAULT ) < 0)
1672  {
1673  vigra_postcondition(false, "createDataset(): Unable to delete existing data.");
1674  }
1675 #endif
1676  }
1677 }
1678 
1679 inline hid_t createFile(std::string filePath, bool append_ = true)
1680 {
1681  FILE * pFile;
1682  pFile = fopen ( filePath.c_str(), "r" );
1683  hid_t file_id;
1684  if ( pFile == NULL )
1685  {
1686  file_id = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
1687  }
1688  else if(append_)
1689  {
1690  fclose( pFile );
1691  file_id = H5Fopen(filePath.c_str(), H5F_ACC_RDWR, H5P_DEFAULT);
1692  }
1693  else
1694  {
1695  fclose(pFile);
1696  std::remove(filePath.c_str());
1697  file_id = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
1698  }
1699  return file_id;
1700 }
1701 
1702 template<unsigned int N, class T, class Tag>
1703 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)
1704 {
1705  std::string path_name(pathInFile), group_name, data_set_name, message;
1706  std::string::size_type delimiter = path_name.rfind('/');
1707 
1708  //create or open file
1709  file_handle = HDF5Handle(createFile(filePath), &H5Fclose,
1710  "createDataset(): unable to open output file.");
1711 
1712  // get the groupname and the filename
1713  if(delimiter == std::string::npos)
1714  {
1715  group_name = "/";
1716  data_set_name = path_name;
1717  }
1718  else
1719  {
1720  group_name = std::string(path_name.begin(), path_name.begin()+delimiter);
1721  data_set_name = std::string(path_name.begin()+delimiter+1, path_name.end());
1722  }
1723 
1724  // create all groups
1725  HDF5Handle group(createGroup(file_handle, group_name), &H5Gclose,
1726  "createDataset(): Unable to create and open group. generic v");
1727 
1728  // delete the dataset if it already exists
1729  deleteDataset(group, data_set_name);
1730 
1731  // create dataspace
1732  // add an extra dimension in case that the data is non-scalar
1733  HDF5Handle dataspace_handle;
1734  if(numBandsOfType > 1) {
1735  // invert dimensions to guarantee c-order
1736  hsize_t shape_inv[N+1]; // one additional dimension for pixel type channel(s)
1737  for(unsigned int k=0; k<N; ++k) {
1738  shape_inv[N-1-k] = array.shape(k); // the channels (eg of an RGB image) are represented by the first dimension (before inversion)
1739  //std::cout << shape_inv[N-k] << " (" << N << ")";
1740  }
1741  shape_inv[N] = numBandsOfType;
1742 
1743  // create dataspace
1744  dataspace_handle = HDF5Handle(H5Screate_simple(N+1, shape_inv, NULL),
1745  &H5Sclose, "createDataset(): unable to create dataspace for non-scalar data.");
1746  } else {
1747  // invert dimensions to guarantee c-order
1748  hsize_t shape_inv[N];
1749  for(unsigned int k=0; k<N; ++k)
1750  shape_inv[N-1-k] = array.shape(k);
1751 
1752  // create dataspace
1753  dataspace_handle = HDF5Handle(H5Screate_simple(N, shape_inv, NULL),
1754  &H5Sclose, "createDataset(): unable to create dataspace for scalar data.");
1755  }
1756 
1757  //alloc memory for dataset.
1758  dataset_handle = HDF5Handle(H5Dcreate(group,
1759  data_set_name.c_str(),
1760  datatype,
1761  dataspace_handle,
1762  H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT),
1763  &H5Dclose, "createDataset(): unable to create dataset.");
1764 }
1765 
1766 
1767 
1768 namespace detail {
1769 
1770 template <class DestIterator, class Shape, class T>
1771 inline void
1772 writeHDF5Impl(DestIterator d, Shape const & shape, const hid_t dataset_id, const hid_t datatype, ArrayVector<T> & buffer, int & counter, const int elements, const int numBandsOfType, MetaInt<0>)
1773 {
1774  DestIterator dend = d + (typename DestIterator::difference_type)shape[0];
1775  int k = 0;
1776  //std::cout << "new:" << std::endl;
1777  for(; d < dend; ++d, k++)
1778  {
1779  buffer[k] = *d;
1780  //std::cout << buffer[k] << " ";
1781  }
1782  //std::cout << std::endl;
1783  HDF5Handle mid1, mid2;
1784 
1785  // select hyperslabs
1786  selectHyperslabs(mid1, mid2, shape, counter, elements, numBandsOfType);
1787 
1788  // write to hdf5
1789  H5Dwrite(dataset_id, datatype, mid2, mid1, H5P_DEFAULT, buffer.data());
1790  // increase counter
1791  counter++;
1792 }
1793 
1794 template <class DestIterator, class Shape, class T, int N>
1795 void
1796 writeHDF5Impl(DestIterator d, Shape const & shape, const hid_t dataset_id, const hid_t datatype, ArrayVector<T> & buffer, int & counter, const int elements, const int numBandsOfType, MetaInt<N>)
1797 {
1798  DestIterator dend = d + (typename DestIterator::difference_type)shape[N];
1799  for(; d < dend; ++d)
1800  {
1801  writeHDF5Impl(d.begin(), shape, dataset_id, datatype, buffer, counter, elements, numBandsOfType, MetaInt<N-1>());
1802  }
1803 }
1804 
1805 } // namespace detail
1806 
1807  /** \brief Store array data in an HDF5 file.
1808 
1809  The number of dimensions, shape and element type of the stored dataset is automatically
1810  determined from the properties of the given \a array. Strided arrays are stored in an
1811  unstrided way, i.e. in contiguous scan-order. Multi-channel element types
1812  (i.e. \ref vigra::RGBValue and \ref vigra::TinyVector) are recognized and handled correctly
1813  (in particular, the will form the innermost dimension of the stored dataset).
1814  \a pathInFile may contain '/'-separated group names, but must end with the name
1815  of the dataset to be created.
1816 
1817  <b> Declaration:</b>
1818 
1819  \code
1820  namespace vigra {
1821  template<unsigned int N, class T, class StrideTag>
1822  void
1823  writeHDF5(const char* filePath, const char* pathInFile,
1824  MultiArrayView<N, T, StrideTag>const & array);
1825  }
1826  \endcode
1827 
1828  <b> Usage:</b>
1829 
1830  <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br>
1831  Namespace: vigra
1832 
1833  \code
1834  MultiArrayShape<3>::type shape(100, 200, 20);
1835  MultiArray<3, int> array(shape);
1836  ... // fill array with data
1837 
1838  writeHDF5("mydata.h5", "/group1/my_dataset", array);
1839  \endcode
1840 */
1841 doxygen_overloaded_function(template <...> void writeHDF5)
1842 
1843 // scalar and unstrided multi arrays
1844 template<unsigned int N, class T>
1845 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, UnstridedArrayTag> & array) // scalar
1846 {
1847  writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 1);
1848 }
1849 
1850 // non-scalar (TinyVector) and unstrided multi arrays
1851 template<unsigned int N, class T, int SIZE>
1852 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array)
1853 {
1854  writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), SIZE);
1855 }
1856 
1857 // non-scalar (RGBValue) and unstrided multi arrays
1858 template<unsigned int N, class T>
1859 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array)
1860 {
1861  writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 3);
1862 }
1863 
1864 // unstrided multi arrays
1865 template<unsigned int N, class T>
1866 void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, UnstridedArrayTag> & array, const hid_t datatype, const int numBandsOfType)
1867 {
1868  HDF5Handle file_handle;
1869  HDF5Handle dataset_handle;
1870  createDataset(filePath, pathInFile, array, datatype, numBandsOfType, file_handle, dataset_handle);
1871 
1872  // Write the data to the HDF5 dataset as is
1873  H5Dwrite( dataset_handle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data()); // .data() possible since void pointer!
1874 
1875  H5Fflush(file_handle, H5F_SCOPE_GLOBAL);
1876 }
1877 
1878 
1879 // scalar and strided multi arrays
1880 template<unsigned int N, class T>
1881 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, StridedArrayTag> & array) // scalar
1882 {
1883  writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 1);
1884 }
1885 
1886 // non-scalar (TinyVector) and strided multi arrays
1887 template<unsigned int N, class T, int SIZE>
1888 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, TinyVector<T, SIZE>, StridedArrayTag> & array)
1889 {
1890  writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), SIZE);
1891 }
1892 
1893 // non-scalar (RGBValue) and strided multi arrays
1894 template<unsigned int N, class T>
1895 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, RGBValue<T>, StridedArrayTag> & array)
1896 {
1897  writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 3);
1898 }
1899 
1900 // strided multi arrays
1901 template<unsigned int N, class T>
1902 void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, StridedArrayTag> & array, const hid_t datatype, const int numBandsOfType)
1903 {
1904  HDF5Handle file_handle;
1905  HDF5Handle dataset_handle;
1906  createDataset(filePath, pathInFile, array, datatype, numBandsOfType, file_handle, dataset_handle);
1907 
1909  vigra::TinyVector<int,N> stride;
1910  int elements = numBandsOfType;
1911  for(unsigned int k=0; k<N; ++k)
1912  {
1913  shape[k] = array.shape(k);
1914  stride[k] = array.stride(k);
1915  elements *= (int)shape[k];
1916  }
1917  int counter = 0;
1918 
1919  ArrayVector<T> buffer((int)array.shape(0));
1920  detail::writeHDF5Impl(array.traverser_begin(), shape, dataset_handle, datatype, buffer, counter, elements, numBandsOfType, vigra::MetaInt<N-1>());
1921 
1922  H5Fflush(file_handle, H5F_SCOPE_GLOBAL);
1923 
1924 }
1925 
1926 namespace detail
1927 {
1928 struct MaxSizeFnc
1929 {
1930  size_t size;
1931 
1932  MaxSizeFnc()
1933  : size(0)
1934  {}
1935 
1936  void operator()(std::string const & in)
1937  {
1938  size = in.size() > size ?
1939  in.size() :
1940  size;
1941  }
1942 };
1943 }
1944 
1945 
1946 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR == 8) || DOXYGEN
1947 /** Write a numeric MultiArray as an attribute with name \a name
1948  of the dataset specified by the handle \a loc.
1949 
1950  <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br>
1951  Namespace: vigra
1952 */
1953 template<size_t N, class T, class C>
1954 void writeHDF5Attr(hid_t loc,
1955  const char* name,
1956  MultiArrayView<N, T, C> const & array)
1957 {
1958  if(H5Aexists(loc, name) > 0)
1959  H5Adelete(loc, name);
1960 
1961  ArrayVector<hsize_t> shape(array.shape().begin(),
1962  array.shape().end());
1963  HDF5Handle
1964  dataspace_handle(H5Screate_simple(N, shape.data(), NULL),
1965  &H5Sclose,
1966  "writeToHDF5File(): unable to create dataspace.");
1967 
1968  HDF5Handle attr(H5Acreate(loc,
1969  name,
1970  detail::getH5DataType<T>(),
1971  dataspace_handle,
1972  H5P_DEFAULT ,H5P_DEFAULT ),
1973  &H5Aclose,
1974  "writeHDF5Attr: unable to create Attribute");
1975 
1976  //copy data - since attributes are small - who cares!
1977  ArrayVector<T> buffer;
1978  for(int ii = 0; ii < array.size(); ++ii)
1979  buffer.push_back(array[ii]);
1980  H5Awrite(attr, detail::getH5DataType<T>(), buffer.data());
1981 }
1982 
1983 /** Write a string MultiArray as an attribute with name \a name
1984  of the dataset specified by the handle \a loc.
1985 
1986  <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br>
1987  Namespace: vigra
1988 */
1989 template<size_t N, class C>
1990 void writeHDF5Attr(hid_t loc,
1991  const char* name,
1992  MultiArrayView<N, std::string, C> const & array)
1993 {
1994  if(H5Aexists(loc, name) > 0)
1995  H5Adelete(loc, name);
1996 
1997  ArrayVector<hsize_t> shape(array.shape().begin(),
1998  array.shape().end());
1999  HDF5Handle
2000  dataspace_handle(H5Screate_simple(N, shape.data(), NULL),
2001  &H5Sclose,
2002  "writeToHDF5File(): unable to create dataspace.");
2003 
2004  HDF5Handle atype(H5Tcopy (H5T_C_S1),
2005  &H5Tclose,
2006  "writeToHDF5File(): unable to create type.");
2007 
2008  detail::MaxSizeFnc max_size;
2009  max_size = std::for_each(array.data(),array.data()+ array.size(), max_size);
2010  H5Tset_size (atype, max_size.size);
2011 
2012  HDF5Handle attr(H5Acreate(loc,
2013  name,
2014  atype,
2015  dataspace_handle,
2016  H5P_DEFAULT ,H5P_DEFAULT ),
2017  &H5Aclose,
2018  "writeHDF5Attr: unable to create Attribute");
2019 
2020  std::string buf ="";
2021  for(int ii = 0; ii < array.size(); ++ii)
2022  {
2023  buf = buf + array[ii]
2024  + std::string(max_size.size - array[ii].size(), ' ');
2025  }
2026  H5Awrite(attr, atype, buf.c_str());
2027 }
2028 
2029 /** Write a numeric ArrayVectorView as an attribute with name \a name
2030  of the dataset specified by the handle \a loc.
2031 
2032  <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br>
2033  Namespace: vigra
2034 */
2035 template<class T>
2036 inline void writeHDF5Attr( hid_t loc,
2037  const char* name,
2038  ArrayVectorView<T> & array)
2039 {
2040  writeHDF5Attr(loc, name,
2042  array.data()));
2043 }
2044 
2045 /** write an Attribute given a file and a path in the file.
2046  * the path in the file should have the format
2047  * [attribute] or /[subgroups/]dataset.attribute or
2048  * /[subgroups/]group.attribute.
2049  * The attribute is written to the root group, a dataset or a subgroup
2050  * respectively
2051  */
2052 template<class Arr>
2053 inline void writeHDF5Attr( std::string filePath,
2054  std::string pathInFile,
2055  Arr & ar)
2056 {
2057  std::string path_name(pathInFile), group_name, data_set_name, message, attr_name;
2058  std::string::size_type delimiter = path_name.rfind('/');
2059 
2060  //create or open file
2061  HDF5Handle file_id(createFile(filePath), &H5Fclose,
2062  "writeToHDF5File(): unable to open output file.");
2063 
2064  // get the groupname and the filename
2065  if(delimiter == std::string::npos)
2066  {
2067  group_name = "/";
2068  data_set_name = path_name;
2069  }
2070 
2071  else
2072  {
2073  group_name = std::string(path_name.begin(), path_name.begin()+delimiter);
2074  data_set_name = std::string(path_name.begin()+delimiter+1, path_name.end());
2075  }
2076  delimiter = data_set_name.rfind('.');
2077  if(delimiter == std::string::npos)
2078  {
2079  attr_name = path_name;
2080  data_set_name = "/";
2081  }
2082  else
2083  {
2084  attr_name = std::string(data_set_name.begin()+delimiter+1, data_set_name.end());
2085  data_set_name = std::string(data_set_name.begin(), data_set_name.begin()+delimiter);
2086  }
2087 
2088  HDF5Handle group(openGroup(file_id, group_name), &H5Gclose,
2089  "writeToHDF5File(): Unable to create and open group. attr ver");
2090 
2091  if(data_set_name != "/")
2092  {
2093  HDF5Handle dset(H5Dopen(group, data_set_name.c_str(), H5P_DEFAULT), &H5Dclose,
2094  "writeHDF5Attr():unable to open dataset");
2095  writeHDF5Attr(hid_t(dset), attr_name.c_str(), ar);
2096  }
2097  else
2098  {
2099  writeHDF5Attr(hid_t(group), attr_name.c_str(), ar);
2100  }
2101 
2102 }
2103 #endif
2104 
2105 //@}
2106 
2107 } // namespace vigra
2108 
2109 #endif // VIGRA_HDF5IMPEX_HXX

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

html generated using doxygen and Python
vigra 1.7.1 (Thu Jun 14 2012)