[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]
vigra/hdf5impex.hxx | ![]() |
00001 /************************************************************************/ 00002 /* */ 00003 /* Copyright 2009 by Michael Hanselmann and Ullrich Koethe */ 00004 /* */ 00005 /* This file is part of the VIGRA computer vision library. */ 00006 /* The VIGRA Website is */ 00007 /* http://hci.iwr.uni-heidelberg.de/vigra/ */ 00008 /* Please direct questions, bug reports, and contributions to */ 00009 /* ullrich.koethe@iwr.uni-heidelberg.de or */ 00010 /* vigra@informatik.uni-hamburg.de */ 00011 /* */ 00012 /* Permission is hereby granted, free of charge, to any person */ 00013 /* obtaining a copy of this software and associated documentation */ 00014 /* files (the "Software"), to deal in the Software without */ 00015 /* restriction, including without limitation the rights to use, */ 00016 /* copy, modify, merge, publish, distribute, sublicense, and/or */ 00017 /* sell copies of the Software, and to permit persons to whom the */ 00018 /* Software is furnished to do so, subject to the following */ 00019 /* conditions: */ 00020 /* */ 00021 /* The above copyright notice and this permission notice shall be */ 00022 /* included in all copies or substantial portions of the */ 00023 /* Software. */ 00024 /* */ 00025 /* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND */ 00026 /* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES */ 00027 /* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND */ 00028 /* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT */ 00029 /* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, */ 00030 /* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING */ 00031 /* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR */ 00032 /* OTHER DEALINGS IN THE SOFTWARE. */ 00033 /* */ 00034 /************************************************************************/ 00035 00036 #ifndef VIGRA_HDF5IMPEX_HXX 00037 #define VIGRA_HDF5IMPEX_HXX 00038 00039 #include <string> 00040 00041 #define H5Gcreate_vers 2 00042 #define H5Gopen_vers 2 00043 #define H5Dopen_vers 2 00044 #define H5Dcreate_vers 2 00045 #define H5Acreate_vers 2 00046 00047 #include <hdf5.h> 00048 00049 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6) 00050 # ifndef H5Gopen 00051 # define H5Gopen(a, b, c) H5Gopen(a, b) 00052 # endif 00053 # ifndef H5Gcreate 00054 # define H5Gcreate(a, b, c, d, e) H5Gcreate(a, b, 1) 00055 # endif 00056 # ifndef H5Dopen 00057 # define H5Dopen(a, b, c) H5Dopen(a, b) 00058 # endif 00059 # ifndef H5Dcreate 00060 # define H5Dcreate(a, b, c, d, e, f, g) H5Dcreate(a, b, c, d, f) 00061 # endif 00062 # ifndef H5Acreate 00063 # define H5Acreate(a, b, c, d, e, f) H5Acreate(a, b, c, d, e) 00064 # endif 00065 # ifndef H5Pset_obj_track_times 00066 # define H5Pset_obj_track_times(a, b) do {} while (0) 00067 # endif 00068 # include <H5LT.h> 00069 #else 00070 # include <hdf5_hl.h> 00071 #endif 00072 00073 #include "impex.hxx" 00074 #include "multi_array.hxx" 00075 #include "multi_impex.hxx" 00076 #include "utilities.hxx" 00077 #include "error.hxx" 00078 00079 namespace vigra { 00080 00081 /** \addtogroup VigraHDF5Impex Import/Export of Images and Arrays in HDF5 Format 00082 00083 Supports arrays with arbitrary element types and arbitrary many dimensions. 00084 See the <a href="http://www.hdfgroup.org/HDF5/">HDF5 Website</a> for more 00085 information on the HDF5 file format. 00086 */ 00087 //@{ 00088 00089 /** \brief Wrapper for hid_t objects. 00090 00091 Newly created or opened HDF5 handles are usually stored as objects of type 'hid_t'. When the handle 00092 is no longer needed, the appropriate close function must be called. However, if a function is 00093 aborted by an exception, this is difficult to ensure. Class HDF5Handle is a smart pointer that 00094 solves this problem by calling the close function in the destructor (This is analogous to how 00095 std::auto_ptr calls 'delete' on the contained pointer). A pointer to the close function must be 00096 passed to the constructor, along with an error message that is raised when creation/opening fails. 00097 00098 Since HDF5Handle objects are convertible to hid_t, they can be used in the code in place 00099 of the latter. 00100 00101 <b>Usage:</b> 00102 00103 \code 00104 HDF5Handle file_id(H5Fopen(filename, H5F_ACC_RDWR, H5P_DEFAULT), 00105 &H5Fclose, 00106 "Error message."); 00107 00108 ... // use file_id in the same way as a plain hid_t object 00109 \endcode 00110 00111 <b>\#include</b> <vigra/hdf5impex.hxx><br> 00112 Namespace: vigra 00113 */ 00114 class HDF5Handle 00115 { 00116 public: 00117 typedef herr_t (*Destructor)(hid_t); 00118 00119 private: 00120 hid_t handle_; 00121 Destructor destructor_; 00122 00123 public: 00124 00125 /** \brief Default constructor. 00126 Creates a NULL handle. 00127 **/ 00128 HDF5Handle() 00129 : handle_( 0 ), 00130 destructor_(0) 00131 {} 00132 00133 /** \brief Create a wrapper for a hid_t object. 00134 00135 The hid_t object \a h is assumed to be the return value of an open or create function. 00136 It will be closed with the given close function \a destructor as soon as this 00137 HDF5Handle is destructed, except when \a destructor is a NULL pointer (in which 00138 case nothing happens at destruction time). If \a h has a value that indicates 00139 failed opening or creation (by HDF5 convention, this means if it is a negative number), 00140 an exception is raised by calling <tt>vigra_fail(error_message)</tt>. 00141 00142 <b>Usage:</b> 00143 00144 \code 00145 HDF5Handle file_id(H5Fopen(filename, H5F_ACC_RDWR, H5P_DEFAULT), 00146 &H5Fclose, 00147 "Error message."); 00148 00149 ... // use file_id in the same way 00150 \endcode 00151 */ 00152 HDF5Handle(hid_t h, Destructor destructor, const char * error_message) 00153 : handle_( h ), 00154 destructor_(destructor) 00155 { 00156 if(handle_ < 0) 00157 vigra_fail(error_message); 00158 } 00159 00160 /** \brief Copy constructor. 00161 Hands over ownership of the RHS handle (analogous to std::auto_ptr). 00162 */ 00163 HDF5Handle(HDF5Handle const & h) 00164 : handle_( h.handle_ ), 00165 destructor_(h.destructor_) 00166 { 00167 const_cast<HDF5Handle &>(h).handle_ = 0; 00168 } 00169 00170 /** \brief Assignment. 00171 Calls close() for the LHS handle and hands over ownership of the 00172 RHS handle (analogous to std::auto_ptr). 00173 */ 00174 HDF5Handle & operator=(HDF5Handle const & h) 00175 { 00176 if(h.handle_ != handle_) 00177 { 00178 close(); 00179 handle_ = h.handle_; 00180 destructor_ = h.destructor_; 00181 const_cast<HDF5Handle &>(h).handle_ = 0; 00182 } 00183 return *this; 00184 } 00185 00186 /** \brief Destructor. 00187 Calls close() for the contained handle. 00188 */ 00189 ~HDF5Handle() 00190 { 00191 close(); 00192 } 00193 00194 /** \brief Explicitly call the stored function (if one has been stored within 00195 this object) for the contained handle and set the handle to NULL. 00196 */ 00197 herr_t close() 00198 { 00199 herr_t res = 1; 00200 if(handle_ && destructor_) 00201 res = (*destructor_)(handle_); 00202 handle_ = 0; 00203 return res; 00204 } 00205 00206 /** \brief Get a temporary hid_t object for the contained handle. 00207 Do not call a close function on the return value - a crash will be likely 00208 otherwise. 00209 */ 00210 hid_t get() const 00211 { 00212 return handle_; 00213 } 00214 00215 /** \brief Convert to a plain hid_t object. 00216 00217 This function ensures that hid_t objects can be transparently replaced with 00218 HDF5Handle objects in user code. Do not call a close function on the return 00219 value - a crash will be likely otherwise. 00220 */ 00221 operator hid_t() const 00222 { 00223 return handle_; 00224 } 00225 00226 /** \brief Equality comparison of the contained handle. 00227 */ 00228 bool operator==(HDF5Handle const & h) const 00229 { 00230 return handle_ == h.handle_; 00231 } 00232 00233 /** \brief Equality comparison of the contained handle. 00234 */ 00235 bool operator==(hid_t h) const 00236 { 00237 return handle_ == h; 00238 } 00239 00240 /** \brief Inequality comparison of the contained handle. 00241 */ 00242 bool operator!=(HDF5Handle const & h) const 00243 { 00244 return handle_ != h.handle_; 00245 } 00246 00247 /** \brief Inequality comparison of the contained handle. 00248 */ 00249 bool operator!=(hid_t h) const 00250 { 00251 return handle_ != h; 00252 } 00253 }; 00254 00255 00256 /********************************************************/ 00257 /* */ 00258 /* HDF5ImportInfo */ 00259 /* */ 00260 /********************************************************/ 00261 00262 /** \brief Argument object for the function readHDF5(). 00263 00264 See \ref readHDF5() for a usage example. This object must be 00265 used to read an image or array from an HDF5 file 00266 and enquire about its properties. 00267 00268 <b>\#include</b> <vigra/hdf5impex.hxx><br> 00269 Namespace: vigra 00270 */ 00271 class HDF5ImportInfo 00272 { 00273 public: 00274 enum PixelType { UINT8, UINT16, UINT32, UINT64, 00275 INT8, INT16, INT32, INT64, 00276 FLOAT, DOUBLE }; 00277 00278 /** Construct HDF5ImportInfo object. 00279 00280 The dataset \a pathInFile in the HDF5 file \a filename is accessed to 00281 read its properties. \a pathInFile may contain '/'-separated group 00282 names, but must end with the name of the desired dataset: 00283 00284 \code 00285 HDF5ImportInfo info(filename, "/group1/group2/my_dataset"); 00286 \endcode 00287 */ 00288 VIGRA_EXPORT HDF5ImportInfo( const char* filePath, const char* pathInFile ); 00289 00290 VIGRA_EXPORT ~HDF5ImportInfo(); 00291 00292 /** Get the filename of this HDF5 object. 00293 */ 00294 VIGRA_EXPORT const std::string& getFilePath() const; 00295 00296 /** Get the dataset's full name in the HDF5 file. 00297 */ 00298 VIGRA_EXPORT const std::string& getPathInFile() const; 00299 00300 /** Get a handle to the file represented by this info object. 00301 */ 00302 VIGRA_EXPORT hid_t getH5FileHandle() const; 00303 00304 /** Get a handle to the dataset represented by this info object. 00305 */ 00306 VIGRA_EXPORT hid_t getDatasetHandle() const; 00307 00308 /** Get the number of dimensions of the dataset represented by this info object. 00309 */ 00310 VIGRA_EXPORT MultiArrayIndex numDimensions() const; 00311 00312 /** Get the shape of the dataset represented by this info object. 00313 */ 00314 VIGRA_EXPORT ArrayVector<hsize_t> const & shape() const 00315 { 00316 return m_dims; 00317 } 00318 00319 /** Get the shape (length) of the dataset along dimension \a dim. 00320 */ 00321 VIGRA_EXPORT MultiArrayIndex shapeOfDimension(const int dim) const; 00322 00323 /** Query the pixel type of the dataset. 00324 00325 Possible values are: 00326 <DL> 00327 <DT>"UINT8"<DD> 8-bit unsigned integer (unsigned char) 00328 <DT>"INT16"<DD> 16-bit signed integer (short) 00329 <DT>"UINT16"<DD> 16-bit unsigned integer (unsigned short) 00330 <DT>"INT32"<DD> 32-bit signed integer (long) 00331 <DT>"UINT32"<DD> 32-bit unsigned integer (unsigned long) 00332 <DT>"FLOAT"<DD> 32-bit floating point (float) 00333 <DT>"DOUBLE"<DD> 64-bit floating point (double) 00334 </DL> 00335 */ 00336 VIGRA_EXPORT const char * getPixelType() const; 00337 00338 /** Query the pixel type of the dataset. 00339 00340 Same as getPixelType(), but the result is returned as a 00341 ImageImportInfo::PixelType enum. This is useful to implement 00342 a switch() on the pixel type. 00343 00344 Possible values are: 00345 <DL> 00346 <DT>UINT8<DD> 8-bit unsigned integer (unsigned char) 00347 <DT>INT16<DD> 16-bit signed integer (short) 00348 <DT>UINT16<DD> 16-bit unsigned integer (unsigned short) 00349 <DT>INT32<DD> 32-bit signed integer (long) 00350 <DT>UINT32<DD> 32-bit unsigned integer (unsigned long) 00351 <DT>FLOAT<DD> 32-bit floating point (float) 00352 <DT>DOUBLE<DD> 64-bit floating point (double) 00353 </DL> 00354 */ 00355 VIGRA_EXPORT PixelType pixelType() const; 00356 00357 private: 00358 HDF5Handle m_file_handle, m_dataset_handle; 00359 std::string m_filename, m_path, m_pixeltype; 00360 hssize_t m_dimensions; 00361 ArrayVector<hsize_t> m_dims; 00362 }; 00363 00364 00365 namespace detail { 00366 00367 template<class type> 00368 inline hid_t getH5DataType() 00369 { 00370 std::runtime_error("getH5DataType(): invalid type"); 00371 return 0; 00372 } 00373 00374 #define VIGRA_H5_DATATYPE(type, h5type) \ 00375 template<> \ 00376 inline hid_t getH5DataType<type>() \ 00377 { return h5type;} 00378 00379 VIGRA_H5_DATATYPE(char, H5T_NATIVE_CHAR) 00380 VIGRA_H5_DATATYPE(float, H5T_NATIVE_FLOAT) 00381 VIGRA_H5_DATATYPE(double, H5T_NATIVE_DOUBLE) 00382 VIGRA_H5_DATATYPE(long double, H5T_NATIVE_LDOUBLE) 00383 00384 // char arrays with flexible length require 'handcrafted' H5 datatype 00385 template<> 00386 inline hid_t getH5DataType<char*>() 00387 { 00388 hid_t stringtype = H5Tcopy (H5T_C_S1); 00389 H5Tset_size(stringtype, H5T_VARIABLE); 00390 return stringtype; 00391 } 00392 template<> 00393 inline hid_t getH5DataType<const char*>() 00394 { 00395 hid_t stringtype = H5Tcopy (H5T_C_S1); 00396 H5Tset_size(stringtype, H5T_VARIABLE); 00397 return stringtype; 00398 } 00399 #undef VIGRA_H5_DATATYPE 00400 00401 #define VIGRA_H5_SIGNED_DATATYPE(type) \ 00402 template<> \ 00403 inline hid_t getH5DataType<type>() \ 00404 { static hid_t types[] = {0, H5T_NATIVE_INT8, H5T_NATIVE_INT16, 0, H5T_NATIVE_INT32, 0,0,0,H5T_NATIVE_INT64}; \ 00405 return types[sizeof(type)];} 00406 00407 VIGRA_H5_SIGNED_DATATYPE(signed char) 00408 VIGRA_H5_SIGNED_DATATYPE(signed short) 00409 VIGRA_H5_SIGNED_DATATYPE(signed int) 00410 VIGRA_H5_SIGNED_DATATYPE(signed long) 00411 VIGRA_H5_SIGNED_DATATYPE(signed long long) 00412 00413 #undef VIGRA_H5_SIGNED_DATATYPE 00414 00415 #define VIGRA_H5_UNSIGNED_DATATYPE(type) \ 00416 template<> \ 00417 inline hid_t getH5DataType<type>() \ 00418 { static hid_t types[] = {0, H5T_NATIVE_UINT8, H5T_NATIVE_UINT16, 0, H5T_NATIVE_UINT32, 0,0,0,H5T_NATIVE_UINT64}; \ 00419 return types[sizeof(type)];} 00420 00421 VIGRA_H5_UNSIGNED_DATATYPE(unsigned char) 00422 VIGRA_H5_UNSIGNED_DATATYPE(unsigned short) 00423 VIGRA_H5_UNSIGNED_DATATYPE(unsigned int) 00424 VIGRA_H5_UNSIGNED_DATATYPE(unsigned long) 00425 VIGRA_H5_UNSIGNED_DATATYPE(unsigned long long) 00426 00427 #undef VIGRA_H5_UNSIGNED_DATATYPE 00428 00429 #if 0 00430 template<> 00431 inline hid_t getH5DataType<FFTWComplex<float> >() 00432 { 00433 hid_t complex_id = H5Tcreate (H5T_COMPOUND, sizeof (FFTWComplex<float>)); 00434 H5Tinsert (complex_id, "real", 0, H5T_NATIVE_FLOAT); 00435 H5Tinsert (complex_id, "imaginary", sizeof(float), H5T_NATIVE_FLOAT); 00436 return complex_id; 00437 } 00438 00439 template<> 00440 inline hid_t getH5DataType<FFTWComplex<double> >() 00441 { 00442 hid_t complex_id = H5Tcreate (H5T_COMPOUND, sizeof (FFTWComplex<double>)); 00443 H5Tinsert (complex_id, "real", 0, H5T_NATIVE_DOUBLE); 00444 H5Tinsert (complex_id, "imaginary", sizeof(double), H5T_NATIVE_DOUBLE); 00445 return complex_id; 00446 } 00447 #endif 00448 00449 00450 } // namespace detail 00451 00452 // helper friend function for callback HDF5_ls_inserter_callback() 00453 void HDF5_ls_insert(void*, const std::string &); 00454 // callback function for ls(), called via HDF5File::H5Literate() 00455 // see http://www.parashift.com/c++-faq-lite/pointers-to-members.html#faq-33.2 00456 // for as to why. 00457 00458 VIGRA_EXPORT H5O_type_t HDF5_get_type(hid_t, const char*); 00459 extern "C" VIGRA_EXPORT herr_t HDF5_ls_inserter_callback(hid_t, const char*, const H5L_info_t*, void*); 00460 00461 /********************************************************/ 00462 /* */ 00463 /* HDF5File */ 00464 /* */ 00465 /********************************************************/ 00466 00467 00468 /** \brief Access to HDF5 files 00469 00470 HDF5File provides a convenient way of accessing data in HDF5 files. vigra::MultiArray 00471 structures of any dimension can be stored to / loaded from HDF5 files. Typical 00472 HDF5 features like subvolume access, chunks and data compression are available, 00473 string attributes can be attached to any dataset or group. Group- or dataset-handles 00474 are encapsulated in the class and managed automatically. The internal file-system like 00475 structure can be accessed by functions like "cd()" or "mkdir()". 00476 00477 00478 <b>Example:</b> 00479 Write the MultiArray out_multi_array to file. Change the current directory to 00480 "/group" and read in the same MultiArray as in_multi_array. 00481 \code 00482 HDF5File file("/path/to/file",HDF5File::New); 00483 file.mkdir("group"); 00484 file.write("/group/dataset", out_multi_array); 00485 00486 file.cd("/group"); 00487 file.read("dataset", in_multi_array); 00488 00489 \endcode 00490 00491 <b>\#include</b> <vigra/hdf5impex.hxx><br> 00492 Namespace: vigra 00493 */ 00494 class HDF5File 00495 { 00496 private: 00497 HDF5Handle fileHandle_; 00498 00499 // current group handle 00500 HDF5Handle cGroupHandle_; 00501 00502 // time tagging of datasets, turned off (= 0) by default. 00503 int track_time; 00504 00505 // helper class for ls() 00506 struct ls_closure 00507 { 00508 virtual void insert(const std::string &) = 0; 00509 virtual ~ls_closure() {} 00510 }; 00511 // datastructure to hold a list of dataset and group names 00512 struct lsOpData : public ls_closure 00513 { 00514 std::vector<std::string> & objects; 00515 lsOpData(std::vector<std::string> & o) : objects(o) {} 00516 void insert(const std::string & x) 00517 { 00518 objects.push_back(x); 00519 } 00520 }; 00521 // (associative-)container closure 00522 template<class Container> 00523 struct ls_container_data : public ls_closure 00524 { 00525 Container & objects; 00526 ls_container_data(Container & o) : objects(o) {} 00527 void insert(const std::string & x) 00528 { 00529 objects.insert(std::string(x)); 00530 } 00531 }; 00532 00533 public: 00534 00535 // helper for callback HDF5_ls_inserter_callback(), used by ls() 00536 friend void HDF5_ls_insert(void*, const std::string &); 00537 00538 /** \brief Set how a file is opened. 00539 OpenMode::New creates a new file. If the file already exists, overwrite it. 00540 00541 OpenMode::Open opens a file for reading/writing. The file will be created, 00542 if necessary. 00543 */ 00544 enum OpenMode { 00545 New, // Create new empty file (existing file will be deleted). 00546 Open // Open file. Create if not existing. 00547 }; 00548 00549 00550 00551 00552 /** \brief Create an HDF5File object. 00553 00554 Creates or opens HDF5 file at position filename. The current group is set 00555 to "/". 00556 */ 00557 HDF5File(std::string filename, OpenMode mode, int track_creation_times = 0) 00558 : track_time(track_creation_times) 00559 { 00560 std::string errorMessage = "HDF5File: Could not create file '" + filename + "'."; 00561 fileHandle_ = HDF5Handle(createFile_(filename, mode), &H5Fclose, errorMessage.c_str()); 00562 cGroupHandle_ = HDF5Handle(openCreateGroup_("/"), &H5Gclose, "HDF5File(): Failed to open root group."); 00563 } 00564 00565 00566 00567 00568 /** \brief Destructor to make sure that all data is flushed before closing the file. 00569 */ 00570 ~HDF5File() 00571 { 00572 //Write everything to disk before closing 00573 H5Fflush(fileHandle_, H5F_SCOPE_GLOBAL); 00574 } 00575 00576 00577 00578 00579 /** \brief Change current group to "/". 00580 */ 00581 inline void root() 00582 { 00583 std::string message = "HDF5File::root(): Could not open group '/'."; 00584 cGroupHandle_ = HDF5Handle(H5Gopen(fileHandle_, "/", H5P_DEFAULT),&H5Gclose,message.c_str()); 00585 } 00586 00587 00588 00589 00590 /** \brief Change the current group. 00591 Both absolute and relative group names are allowed. 00592 */ 00593 inline void cd(std::string groupName) 00594 { 00595 std::string message = "HDF5File::cd(): Could not open group '" + groupName + "'.\n"; 00596 00597 // make groupName clean 00598 groupName = get_absolute_path(groupName); 00599 00600 if(groupName == "/"){ 00601 cGroupHandle_ = HDF5Handle(openCreateGroup_("/"),&H5Gclose,message.c_str()); 00602 return; 00603 } 00604 else{ 00605 if (H5Lexists(fileHandle_, groupName.c_str(), H5P_DEFAULT) == 0) 00606 { 00607 std::cerr << message; 00608 return; 00609 } 00610 cGroupHandle_ = HDF5Handle(openCreateGroup_(groupName),&H5Gclose,message.c_str()); 00611 } 00612 } 00613 00614 00615 00616 00617 /** \brief Change the current group to its parent group. 00618 returns true if successful, false otherwise. 00619 */ 00620 inline bool cd_up() 00621 { 00622 std::string groupName = currentGroupName_(); 00623 00624 //do not try to move up if we already in "/" 00625 if(groupName == "/"){ 00626 return false; 00627 } 00628 00629 size_t lastSlash = groupName.find_last_of('/'); 00630 00631 std::string parentGroup (groupName.begin(), groupName.begin()+lastSlash+1); 00632 00633 cd(parentGroup); 00634 00635 return true; 00636 } 00637 inline void cd_up(int levels) 00638 { 00639 for(int i = 0; i<levels; i++) 00640 cd_up(); 00641 } 00642 00643 00644 00645 00646 /** \brief Create a new group. 00647 If the first character is a "/", the path will be interpreted as absolute path, 00648 otherwise it will be interpreted as path relative to the current group. 00649 */ 00650 inline void mkdir(std::string groupName) 00651 { 00652 // make groupName clean 00653 groupName = get_absolute_path(groupName); 00654 00655 hid_t handle = openCreateGroup_(groupName.c_str()); 00656 if (handle != cGroupHandle_){ 00657 H5Gclose(handle); 00658 } 00659 } 00660 00661 00662 00663 00664 /** \brief Change the current group; create it if necessary. 00665 If the first character is a "/", the path will be interpreted as absolute path, 00666 otherwise it will be interpreted as path relative to the current group. 00667 */ 00668 inline void cd_mk(std::string groupName) 00669 { 00670 // make groupName clean 00671 groupName = get_absolute_path(groupName); 00672 00673 std::string message = "HDF5File::cd_mk(): Could not create group '" + groupName + "'."; 00674 cGroupHandle_ = HDF5Handle(openCreateGroup_(groupName.c_str()),&H5Gclose,message.c_str()); 00675 } 00676 00677 // helper function for the various ls() variants. 00678 void ls_H5Literate(ls_closure & data) 00679 { 00680 H5Literate(cGroupHandle_, H5_INDEX_NAME, H5_ITER_NATIVE, NULL, 00681 HDF5_ls_inserter_callback, static_cast<void*>(&data)); 00682 } 00683 00684 /** \brief List the contents of the current group. 00685 The function returns a vector of strings holding the entries of the 00686 current group. Only datasets and groups are listed, other objects 00687 (e.g. datatypes) are ignored. Group names always have an ending "/". 00688 */ 00689 /** 00690 * 00691 */ 00692 inline std::vector<std::string> ls() 00693 { 00694 std::vector<std::string> list; 00695 lsOpData data(list); 00696 ls_H5Literate(data); 00697 return list; 00698 } 00699 00700 /** \brief List the contents of the current group into a container-like 00701 object via insert(). Only datasets and groups are inserted, other 00702 objects (e.g., datatypes) are ignored. Group names always have an 00703 ending "/". 00704 00705 The argument cont is presumably an associative container, however, 00706 only its member function <tt>cont.insert(std::string)</tt> will be 00707 called. 00708 \param cont reference to a container supplying a member function 00709 <tt>insert(const i_type &)</tt>, where <tt>i_type</tt> 00710 is convertible to <tt>std::string</tt>. 00711 */ 00712 template<class Container> 00713 void ls(Container & cont) 00714 { 00715 ls_container_data<Container> data(cont); 00716 ls_H5Literate(data); 00717 } 00718 00719 00720 /** \brief Get the path of the current group. 00721 */ 00722 inline std::string pwd() 00723 { 00724 return currentGroupName_(); 00725 } 00726 00727 00728 /** \brief Get the name of the associated file. 00729 */ 00730 inline std::string filename() 00731 { 00732 return fileName_(); 00733 } 00734 00735 00736 /** \brief Get the number of dimensions of a certain dataset 00737 If the first character is a "/", the path will be interpreted as absolute path, 00738 otherwise it will be interpreted as path relative to the current group. 00739 */ 00740 inline hssize_t getDatasetDimensions(std::string datasetName) 00741 { 00742 // make datasetName clean 00743 datasetName = get_absolute_path(datasetName); 00744 00745 //Open dataset and dataspace 00746 std::string errorMessage = "HDF5File::getDatasetDimensions(): Unable to open dataset '" + datasetName + "'."; 00747 HDF5Handle datasetHandle = HDF5Handle(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str()); 00748 00749 errorMessage = "HDF5File::getDatasetDimensions(): Unable to access dataspace."; 00750 HDF5Handle dataspaceHandle(H5Dget_space(datasetHandle), &H5Sclose, errorMessage.c_str()); 00751 00752 //return dimension information 00753 return H5Sget_simple_extent_ndims(dataspaceHandle); 00754 } 00755 00756 00757 00758 00759 /** \brief Get the shape of each dimension of a certain dataset. 00760 Normally, this function is called after determining the dimension of the 00761 dataset using \ref getDatasetDimensions(). 00762 If the first character is a "/", the path will be interpreted as absolute path, 00763 otherwise it will be interpreted as path relative to the current group. 00764 */ 00765 inline ArrayVector<hsize_t> getDatasetShape(std::string datasetName) 00766 { 00767 // make datasetName clean 00768 datasetName = get_absolute_path(datasetName); 00769 00770 //Open dataset and dataspace 00771 std::string errorMessage = "HDF5File::getDatasetShape(): Unable to open dataset '" + datasetName + "'."; 00772 HDF5Handle datasetHandle = HDF5Handle(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str()); 00773 00774 errorMessage = "HDF5File::getDatasetShape(): Unable to access dataspace."; 00775 HDF5Handle dataspaceHandle(H5Dget_space(datasetHandle), &H5Sclose, errorMessage.c_str()); 00776 00777 //get dimension information 00778 ArrayVector<hsize_t>::size_type dimensions = H5Sget_simple_extent_ndims(dataspaceHandle); 00779 00780 ArrayVector<hsize_t> shape(dimensions); 00781 ArrayVector<hsize_t> maxdims(dimensions); 00782 H5Sget_simple_extent_dims(dataspaceHandle, shape.data(), maxdims.data()); 00783 00784 // invert the dimensions to guarantee c-order 00785 ArrayVector<hsize_t> shape_inv(dimensions); 00786 for(ArrayVector<hsize_t>::size_type i=0; i<shape.size(); i++) { 00787 shape_inv[i] = shape[dimensions-1-i]; 00788 } 00789 00790 return shape_inv; 00791 } 00792 00793 00794 00795 00796 /** \brief Obtain the HDF5 handle of a dataset. 00797 */ 00798 inline hid_t getDatasetHandle(std::string dataset_name) 00799 { 00800 return getDatasetHandle_(dataset_name); 00801 } 00802 00803 00804 00805 00806 /** \brief Obtain the HDF5 handle of a group. 00807 */ 00808 inline hid_t getGroupHandle(std::string group_name) 00809 { 00810 // make group_name clean 00811 group_name = get_absolute_path(group_name); 00812 00813 // group must exist 00814 vigra_precondition((H5Lexists(fileHandle_, group_name.c_str(), H5P_DEFAULT) == 1), "Error: Group '" + group_name + "' does not exist."); 00815 00816 //open group and return group handle 00817 return openCreateGroup_( group_name); 00818 } 00819 00820 00821 00822 00823 /** \brief Obtain the HDF5 handle of a attribute. 00824 */ 00825 inline hid_t getAttributeHandle(std::string dataset_name, std::string attribute_name) 00826 { 00827 HDF5Handle dset (getDatasetHandle_(dataset_name),&H5Dclose,std::string("Error: Dataset '"+dataset_name+"' not found.").c_str()); 00828 00829 return H5Aopen(dset, attribute_name.c_str(),H5P_DEFAULT); 00830 } 00831 00832 00833 00834 00835 /* Writing Attributes */ 00836 00837 /** \brief Write MultiArray Attributes. 00838 * In contrast to datasets, subarray access, chunks and compression are not available. 00839 */ 00840 template<unsigned int N, class T> 00841 inline void writeAttribute(std::string object_name, std::string attribute_name, const MultiArrayView<N, T, UnstridedArrayTag> & array) 00842 { 00843 // make object_name clean 00844 object_name = get_absolute_path(object_name); 00845 00846 write_attribute_(object_name, attribute_name, array, detail::getH5DataType<T>(), 1); 00847 } 00848 00849 00850 00851 template<unsigned int N, class T, int SIZE> 00852 inline void writeAttribute(std::string datasetName, std::string attributeName, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array) 00853 { 00854 // make datasetName clean 00855 datasetName = get_absolute_path(datasetName); 00856 00857 write_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), SIZE); 00858 } 00859 00860 00861 00862 template<unsigned int N, class T> 00863 inline void writeAttribute(std::string datasetName, std::string attributeName, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array) 00864 { 00865 // make datasetName clean 00866 datasetName = get_absolute_path(datasetName); 00867 00868 write_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), 3); 00869 } 00870 00871 00872 00873 00874 /** \brief Write a single value. 00875 Specialization of the write function for simple datatypes 00876 */ 00877 inline void writeAttribute(std::string object_name, std::string attribute_name, char data) 00878 { writeAtomicAttribute(object_name,attribute_name,data); } 00879 inline void writeAttribute(std::string datasetName, std::string attributeName, signed char data) 00880 { writeAtomicAttribute(datasetName,attributeName,data); } 00881 inline void writeAttribute(std::string datasetName, std::string attributeName, signed short data) 00882 { writeAtomicAttribute(datasetName,attributeName,data); } 00883 inline void writeAttribute(std::string datasetName, std::string attributeName, signed int data) 00884 { writeAtomicAttribute(datasetName,attributeName,data); } 00885 inline void writeAttribute(std::string datasetName, std::string attributeName, signed long data) 00886 { writeAtomicAttribute(datasetName,attributeName,data); } 00887 inline void writeAttribute(std::string datasetName, std::string attributeName, signed long long data) 00888 { writeAtomicAttribute(datasetName,attributeName,data); } 00889 inline void writeAttribute(std::string datasetName, std::string attributeName, unsigned char data) 00890 { writeAtomicAttribute(datasetName,attributeName,data); } 00891 inline void writeAttribute(std::string datasetName, std::string attributeName, unsigned short data) 00892 { writeAtomicAttribute(datasetName,attributeName,data); } 00893 inline void writeAttribute(std::string datasetName, std::string attributeName, unsigned int data) 00894 { writeAtomicAttribute(datasetName,attributeName,data); } 00895 inline void writeAttribute(std::string datasetName, std::string attributeName, unsigned long data) 00896 { writeAtomicAttribute(datasetName,attributeName,data); } 00897 inline void writeAttribute(std::string datasetName, std::string attributeName, unsigned long long data) 00898 { writeAtomicAttribute(datasetName,attributeName,data); } 00899 inline void writeAttribute(std::string datasetName, std::string attributeName, float data) 00900 { writeAtomicAttribute(datasetName,attributeName,data); } 00901 inline void writeAttribute(std::string datasetName, std::string attributeName, double data) 00902 { writeAtomicAttribute(datasetName,attributeName,data); } 00903 inline void writeAttribute(std::string datasetName, std::string attributeName, long double data) 00904 { writeAtomicAttribute(datasetName,attributeName,data); } 00905 inline void writeAttribute(std::string datasetName, std::string attributeName, const char* data) 00906 { writeAtomicAttribute(datasetName,attributeName,data); } 00907 inline void writeAttribute(std::string datasetName, std::string attributeName, std::string const & data) 00908 { writeAtomicAttribute(datasetName,attributeName,data.c_str()); } 00909 00910 /** \brief Test if attribute exists. 00911 00912 */ 00913 bool existsAttribute(std::string object_name, std::string attribute_name) 00914 { 00915 std::string obj_path = get_absolute_path(object_name); 00916 htri_t exists = H5Aexists_by_name(fileHandle_, obj_path.c_str(), 00917 attribute_name.c_str(), H5P_DEFAULT); 00918 vigra_precondition(exists >= 0, "HDF5File::existsAttribute(): " 00919 "object \"" + object_name + "\" " 00920 "not found."); 00921 return exists != 0; 00922 } 00923 00924 // Reading Attributes 00925 00926 /** \brief Read MultiArray Attributes. 00927 * In contrast to datasets, subarray access is not available. 00928 */ 00929 template<unsigned int N, class T> 00930 inline void readAttribute(std::string object_name, std::string attribute_name, const MultiArrayView<N, T, UnstridedArrayTag> & array) 00931 { 00932 // make object_name clean 00933 object_name = get_absolute_path(object_name); 00934 00935 read_attribute_(object_name, attribute_name, array, detail::getH5DataType<T>(), 1); 00936 } 00937 00938 00939 00940 template<unsigned int N, class T, int SIZE> 00941 inline void readAttribute(std::string datasetName, std::string attributeName, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array) 00942 { 00943 // make datasetName clean 00944 datasetName = get_absolute_path(datasetName); 00945 00946 read_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), SIZE); 00947 } 00948 00949 00950 00951 template<unsigned int N, class T> 00952 inline void readAttribute(std::string datasetName, std::string attributeName, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array) 00953 { 00954 // make datasetName clean 00955 datasetName = get_absolute_path(datasetName); 00956 00957 read_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), 3); 00958 } 00959 00960 00961 00962 00963 /** \brief Read a single value. 00964 Specialization of the read function for simple datatypes 00965 */ 00966 inline void readAttribute(std::string object_name, std::string attribute_name, char &data) 00967 { readAtomicAttribute(object_name,attribute_name,data); } 00968 inline void readAttribute(std::string datasetName, std::string attributeName, signed char &data) 00969 { readAtomicAttribute(datasetName,attributeName,data); } 00970 inline void readAttribute(std::string datasetName, std::string attributeName, signed short &data) 00971 { readAtomicAttribute(datasetName,attributeName,data); } 00972 inline void readAttribute(std::string datasetName, std::string attributeName, signed int &data) 00973 { readAtomicAttribute(datasetName,attributeName,data); } 00974 inline void readAttribute(std::string datasetName, std::string attributeName, signed long &data) 00975 { readAtomicAttribute(datasetName,attributeName,data); } 00976 inline void readAttribute(std::string datasetName, std::string attributeName, signed long long &data) 00977 { readAtomicAttribute(datasetName,attributeName,data); } 00978 inline void readAttribute(std::string datasetName, std::string attributeName, unsigned char &data) 00979 { readAtomicAttribute(datasetName,attributeName,data); } 00980 inline void readAttribute(std::string datasetName, std::string attributeName, unsigned short &data) 00981 { readAtomicAttribute(datasetName,attributeName,data); } 00982 inline void readAttribute(std::string datasetName, std::string attributeName, unsigned int &data) 00983 { readAtomicAttribute(datasetName,attributeName,data); } 00984 inline void readAttribute(std::string datasetName, std::string attributeName, unsigned long &data) 00985 { readAtomicAttribute(datasetName,attributeName,data); } 00986 inline void readAttribute(std::string datasetName, std::string attributeName, unsigned long long &data) 00987 { readAtomicAttribute(datasetName,attributeName,data); } 00988 inline void readAttribute(std::string datasetName, std::string attributeName, float &data) 00989 { readAtomicAttribute(datasetName,attributeName,data); } 00990 inline void readAttribute(std::string datasetName, std::string attributeName, double &data) 00991 { readAtomicAttribute(datasetName,attributeName,data); } 00992 inline void readAttribute(std::string datasetName, std::string attributeName, long double &data) 00993 { readAtomicAttribute(datasetName,attributeName,data); } 00994 inline void readAttribute(std::string datasetName, std::string attributeName, std::string &data) 00995 { readAtomicAttribute(datasetName,attributeName,data); } 00996 00997 00998 00999 01000 // Writing data 01001 01002 /** \brief Write multi arrays. 01003 01004 Chunks can be activated by setting 01005 \code iChunkSize = size; //size > 0 01006 \endcode . 01007 The chunks will be hypercubes with edge length size. 01008 01009 Compression can be activated by setting 01010 \code compression = parameter; // 0 < parameter <= 9 01011 \endcode 01012 where 0 stands for no compression and 9 for maximum compression. 01013 01014 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 01015 otherwise it will be interpreted as path relative to the current group. 01016 */ 01017 template<unsigned int N, class T> 01018 inline void write(std::string datasetName, const MultiArrayView<N, T, UnstridedArrayTag> & array, int iChunkSize = 0, int compression = 0) 01019 { 01020 // make datasetName clean 01021 datasetName = get_absolute_path(datasetName); 01022 01023 typename MultiArrayShape<N>::type chunkSize; 01024 for(int i = 0; i < N; i++){ 01025 chunkSize[i] = iChunkSize; 01026 } 01027 write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize, compression); 01028 } 01029 01030 01031 01032 /** \brief Write multi arrays. 01033 Chunks can be activated by providing a MultiArrayShape as chunkSize. 01034 chunkSize must have equal dimension as array. 01035 01036 Compression can be activated by setting 01037 \code compression = parameter; // 0 < parameter <= 9 01038 \endcode 01039 where 0 stands for no compression and 9 for maximum compression. 01040 01041 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 01042 otherwise it will be interpreted as path relative to the current group. 01043 */ 01044 template<unsigned int N, class T> 01045 inline void write(std::string datasetName, const MultiArrayView<N, T, UnstridedArrayTag> & array, typename MultiArrayShape<N>::type chunkSize, int compression = 0) 01046 { 01047 // make datasetName clean 01048 datasetName = get_absolute_path(datasetName); 01049 01050 write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize, compression); 01051 } 01052 01053 01054 01055 /** \brief Write a multi array into a larger volume. 01056 blockOffset determines the position, where array is written. 01057 01058 Chunks can be activated by providing a MultiArrayShape as chunkSize. 01059 chunkSize must have equal dimension as array. 01060 01061 Compression can be activated by setting 01062 \code compression = parameter; // 0 < parameter <= 9 01063 \endcode 01064 where 0 stands for no compression and 9 for maximum compression. 01065 01066 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 01067 otherwise it will be interpreted as path relative to the current group. 01068 */ 01069 template<unsigned int N, class T> 01070 inline void writeBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, const MultiArrayView<N, T, UnstridedArrayTag> & array) 01071 { 01072 // make datasetName clean 01073 datasetName = get_absolute_path(datasetName); 01074 01075 writeBlock_(datasetName, blockOffset, array, detail::getH5DataType<T>(), 1); 01076 } 01077 01078 01079 01080 01081 // non-scalar (TinyVector) and unstrided multi arrays 01082 template<unsigned int N, class T, int SIZE> 01083 inline void write(std::string datasetName, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array, int iChunkSize = 0, int compression = 0) 01084 { 01085 // make datasetName clean 01086 datasetName = get_absolute_path(datasetName); 01087 01088 typename MultiArrayShape<N>::type chunkSize; 01089 for(int i = 0; i < N; i++){ 01090 chunkSize[i] = iChunkSize; 01091 } 01092 write_(datasetName, array, detail::getH5DataType<T>(), SIZE, chunkSize, compression); 01093 } 01094 01095 01096 01097 template<unsigned int N, class T, int SIZE> 01098 inline void write(std::string datasetName, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array, typename MultiArrayShape<N>::type chunkSize, int compression = 0) 01099 { 01100 // make datasetName clean 01101 datasetName = get_absolute_path(datasetName); 01102 01103 write_(datasetName, array, detail::getH5DataType<T>(), SIZE, chunkSize, compression); 01104 } 01105 01106 01107 /** \brief Write array vectors. 01108 01109 Compression can be activated by setting 01110 \code compression = parameter; // 0 < parameter <= 9 01111 \endcode 01112 where 0 stands for no compression and 9 for maximum compression. 01113 01114 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 01115 otherwise it will be interpreted as path relative to the current group. 01116 */ 01117 template<class T> 01118 void write(const std::string & datasetName, 01119 const ArrayVectorView<T> & array, 01120 int compression = 0) 01121 { 01122 // convert to a (trivial) MultiArrayView and forward. 01123 MultiArrayShape<1>::type shape(array.size()); 01124 const MultiArrayView<1, T> m_array(shape, const_cast<T*>(array.data())); 01125 write(datasetName, m_array, compression); 01126 } 01127 01128 01129 template<unsigned int N, class T, int SIZE> 01130 inline void writeBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array) 01131 { 01132 // make datasetName clean 01133 datasetName = get_absolute_path(datasetName); 01134 01135 writeBlock_(datasetName, blockOffset, array, detail::getH5DataType<T>(), SIZE); 01136 } 01137 01138 01139 01140 01141 // non-scalar (RGBValue) and unstrided multi arrays 01142 template<unsigned int N, class T> 01143 inline void write(std::string datasetName, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array, int iChunkSize = 0, int compression = 0) 01144 { 01145 // make datasetName clean 01146 datasetName = get_absolute_path(datasetName); 01147 01148 typename MultiArrayShape<N>::type chunkSize; 01149 for(int i = 0; i < N; i++){ 01150 chunkSize[i] = iChunkSize; 01151 } 01152 write_(datasetName, array, detail::getH5DataType<T>(), 3, chunkSize, compression); 01153 } 01154 01155 01156 01157 template<unsigned int N, class T> 01158 inline void write(std::string datasetName, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array, typename MultiArrayShape<N>::type chunkSize, int compression = 0) 01159 { 01160 // make datasetName clean 01161 datasetName = get_absolute_path(datasetName); 01162 01163 write_(datasetName, array, detail::getH5DataType<T>(), 3, chunkSize, compression); 01164 } 01165 01166 01167 01168 template<unsigned int N, class T> 01169 inline void writeBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array) 01170 { 01171 // make datasetName clean 01172 datasetName = get_absolute_path(datasetName); 01173 01174 writeBlock_(datasetName, blockOffset, array, detail::getH5DataType<T>(), 3); 01175 } 01176 01177 01178 01179 01180 /** \brief Write a single value. 01181 Specialization of the write function for simple datatypes 01182 */ 01183 inline void write(std::string datasetName, char data) { writeAtomic(datasetName,data); } 01184 inline void write(std::string datasetName, signed char data) { writeAtomic(datasetName,data); } 01185 inline void write(std::string datasetName, signed short data) { writeAtomic(datasetName,data); } 01186 inline void write(std::string datasetName, signed int data) { writeAtomic(datasetName,data); } 01187 inline void write(std::string datasetName, signed long data) { writeAtomic(datasetName,data); } 01188 inline void write(std::string datasetName, signed long long data) { writeAtomic(datasetName,data); } 01189 inline void write(std::string datasetName, unsigned char data) { writeAtomic(datasetName,data); } 01190 inline void write(std::string datasetName, unsigned short data) { writeAtomic(datasetName,data); } 01191 inline void write(std::string datasetName, unsigned int data) { writeAtomic(datasetName,data); } 01192 inline void write(std::string datasetName, unsigned long data) { writeAtomic(datasetName,data); } 01193 inline void write(std::string datasetName, unsigned long long data) { writeAtomic(datasetName,data); } 01194 inline void write(std::string datasetName, float data) { writeAtomic(datasetName,data); } 01195 inline void write(std::string datasetName, double data) { writeAtomic(datasetName,data); } 01196 inline void write(std::string datasetName, long double data) { writeAtomic(datasetName,data); } 01197 inline void write(std::string datasetName, const char* data) { writeAtomic(datasetName,data); } 01198 inline void write(std::string datasetName, std::string const & data) { writeAtomic(datasetName,data.c_str()); } 01199 01200 01201 01202 01203 01204 // Reading data 01205 01206 /** \brief Read data into a multi array. 01207 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 01208 otherwise it will be interpreted as path relative to the current group. 01209 */ 01210 template<unsigned int N, class T> 01211 inline void read(std::string datasetName, MultiArrayView<N, T, UnstridedArrayTag> & array) 01212 { 01213 // make datasetName clean 01214 datasetName = get_absolute_path(datasetName); 01215 01216 read_(datasetName, array, detail::getH5DataType<T>(), 1); 01217 } 01218 01219 01220 01221 /** \brief Read data into a MultiArray. Resize MultiArray to the correct size. 01222 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 01223 otherwise it will be interpreted as path relative to the current group. 01224 */ 01225 template<unsigned int N, class T> 01226 inline void readAndResize(std::string datasetName, MultiArray<N, T> & array) 01227 { 01228 // make datasetName clean 01229 datasetName = get_absolute_path(datasetName); 01230 01231 // get dataset dimension 01232 ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName); 01233 hssize_t dimensions = getDatasetDimensions(datasetName); 01234 01235 // check if dimensions are correct 01236 vigra_precondition((N == MultiArrayIndex(dimensions)), // the object in the HDF5 file may have one additional dimension which we then interpret as the pixel type bands 01237 "HDF5File::readAndResize(): Array dimension disagrees with dataset dimension."); 01238 typename MultiArrayShape<N>::type shape; 01239 01240 for(int k=0; k< MultiArrayIndex(dimensions); ++k) { 01241 shape[k] = MultiArrayIndex(dimshape[k]); 01242 } 01243 01244 // reshape target MultiArray 01245 array.reshape(shape); 01246 01247 read_(datasetName, array, detail::getH5DataType<T>(), 1); 01248 } 01249 01250 /** \brief Read data into an array vector. 01251 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 01252 otherwise it will be interpreted as path relative to the current group. 01253 */ 01254 template<class T> 01255 inline void read(const std::string & datasetName, ArrayVectorView<T> & array) 01256 { 01257 // convert to a (trivial) MultiArrayView and forward. 01258 MultiArrayShape<1>::type shape(array.size()); 01259 MultiArrayView<1, T> m_array(shape, (array.data())); 01260 read(datasetName, m_array); 01261 } 01262 01263 /** \brief Read data into an array vector. Resize the array vector to the correct size. 01264 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 01265 otherwise it will be interpreted as path relative to the current group. 01266 */ 01267 template<class T> 01268 inline void readAndResize(std::string datasetName, 01269 ArrayVector<T> & array) 01270 { 01271 // make dataset name clean 01272 datasetName = get_absolute_path(datasetName); 01273 01274 // get dataset dimension 01275 ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName); 01276 hssize_t dimensions = getDatasetDimensions(datasetName); 01277 01278 // check if dimensions are correct 01279 vigra_precondition((1 == MultiArrayIndex(dimensions)), 01280 "HDF5File::readAndResize(): Array dimension disagrees with Dataset dimension must equal one for vigra::ArrayVector."); 01281 // resize target array vector 01282 array.resize((typename ArrayVector<T>::size_type)dimshape[0]); 01283 // convert to a (trivial) MultiArrayView and forward. 01284 MultiArrayShape<1>::type shape(array.size()); 01285 MultiArrayView<1, T> m_array(shape, (array.data())); 01286 01287 read_(datasetName, m_array, detail::getH5DataType<T>(), 1); 01288 } 01289 01290 /** \brief Read a block of data into a multi array. 01291 This function allows to read a small block out of a larger volume stored 01292 in an HDF5 dataset. 01293 01294 blockOffset determines the position of the block. 01295 blockSize determines the size in each dimension of the block. 01296 01297 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 01298 otherwise it will be interpreted as path relative to the current group. 01299 */ 01300 template<unsigned int N, class T> 01301 inline void readBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, typename MultiArrayShape<N>::type blockShape, MultiArrayView<N, T, UnstridedArrayTag> & array) 01302 { 01303 // make datasetName clean 01304 datasetName = get_absolute_path(datasetName); 01305 01306 readBlock_(datasetName, blockOffset, blockShape, array, detail::getH5DataType<T>(), 1); 01307 } 01308 01309 01310 01311 01312 // non-scalar (TinyVector) and unstrided target MultiArrayView 01313 template<unsigned int N, class T, int SIZE> 01314 inline void read(std::string datasetName, MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array) 01315 { 01316 // make datasetName clean 01317 datasetName = get_absolute_path(datasetName); 01318 01319 read_(datasetName, array, detail::getH5DataType<T>(), SIZE); 01320 } 01321 01322 01323 01324 // non-scalar (TinyVector) MultiArray 01325 template<unsigned int N, class T, int SIZE> 01326 inline void readAndResize(std::string datasetName, MultiArray<N, TinyVector<T, SIZE> > & array) 01327 { 01328 // make datasetName clean 01329 datasetName = get_absolute_path(datasetName); 01330 01331 // get dataset dimension 01332 ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName); 01333 hssize_t dimensions = getDatasetDimensions(datasetName); 01334 01335 // check if dimensions are correct 01336 vigra_precondition(((N+1) == MultiArrayIndex(dimensions)), // the object in the HDF5 file may have one additional dimension which we then interpret as the pixel type bands 01337 "HDF5File::readAndResize(): Array dimension disagrees with dataset dimension."); 01338 typename MultiArrayShape<N>::type shape; 01339 01340 for(int k=1; k< MultiArrayIndex(dimensions); ++k) { 01341 shape[k-1] = MultiArrayIndex(dimshape[k]); 01342 } 01343 01344 // reshape target MultiArray 01345 array.reshape(shape); 01346 01347 read_(datasetName, array, detail::getH5DataType<T>(), SIZE); 01348 } 01349 01350 01351 01352 template<unsigned int N, class T, int SIZE> 01353 inline void readBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, typename MultiArrayShape<N>::type blockShape, MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array) 01354 { 01355 // make datasetName clean 01356 datasetName = get_absolute_path(datasetName); 01357 01358 readBlock_(datasetName, blockOffset, blockShape, array, detail::getH5DataType<T>(), SIZE); 01359 } 01360 01361 01362 01363 01364 // non-scalar (RGBValue) and unstrided target MultiArrayView 01365 template<unsigned int N, class T> 01366 inline void read(std::string datasetName, MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array) 01367 { 01368 // make datasetName clean 01369 datasetName = get_absolute_path(datasetName); 01370 01371 read_(datasetName, array, detail::getH5DataType<T>(), 3); 01372 } 01373 01374 01375 01376 // non-scalar (RGBValue) MultiArray 01377 template<unsigned int N, class T> 01378 inline void readAndResize(std::string datasetName, MultiArray<N, RGBValue<T> > & array) 01379 { 01380 // make datasetName clean 01381 datasetName = get_absolute_path(datasetName); 01382 01383 // get dataset dimension 01384 ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName); 01385 hssize_t dimensions = getDatasetDimensions(datasetName); 01386 01387 // check if dimensions are correct 01388 vigra_precondition(((N+1) == MultiArrayIndex(dimensions)), // the object in the HDF5 file may have one additional dimension which we then interpret as the pixel type bands 01389 "HDF5File::readAndResize(): Array dimension disagrees with dataset dimension."); 01390 typename MultiArrayShape<N>::type shape; 01391 01392 for(int k=1; k< MultiArrayIndex(dimensions); ++k) { 01393 shape[k-1] = MultiArrayIndex(dimshape[k]); 01394 } 01395 01396 // reshape target MultiArray 01397 array.reshape(shape); 01398 01399 read_(datasetName, array, detail::getH5DataType<T>(), 3); 01400 } 01401 01402 01403 01404 template<unsigned int N, class T> 01405 inline void readBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, typename MultiArrayShape<N>::type blockShape, MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array) 01406 { 01407 // make datasetName clean 01408 datasetName = get_absolute_path(datasetName); 01409 01410 readBlock_(datasetName, blockOffset, blockShape, array, detail::getH5DataType<T>(), 3); 01411 } 01412 01413 01414 01415 01416 /** \brief Read a single value. 01417 Specialization of the read function for simple datatypes 01418 */ 01419 inline void read(std::string datasetName, char &data) { readAtomic(datasetName,data); } 01420 inline void read(std::string datasetName, signed char &data) { readAtomic(datasetName,data); } 01421 inline void read(std::string datasetName, signed short &data) { readAtomic(datasetName,data); } 01422 inline void read(std::string datasetName, signed int &data) { readAtomic(datasetName,data); } 01423 inline void read(std::string datasetName, signed long &data) { readAtomic(datasetName,data); } 01424 inline void read(std::string datasetName, signed long long &data) { readAtomic(datasetName,data); } 01425 inline void read(std::string datasetName, unsigned char &data) { readAtomic(datasetName,data); } 01426 inline void read(std::string datasetName, unsigned short &data) { readAtomic(datasetName,data); } 01427 inline void read(std::string datasetName, unsigned int &data) { readAtomic(datasetName,data); } 01428 inline void read(std::string datasetName, unsigned long &data) { readAtomic(datasetName,data); } 01429 inline void read(std::string datasetName, unsigned long long &data) { readAtomic(datasetName,data); } 01430 inline void read(std::string datasetName, float &data) { readAtomic(datasetName,data); } 01431 inline void read(std::string datasetName, double &data) { readAtomic(datasetName,data); } 01432 inline void read(std::string datasetName, long double &data) { readAtomic(datasetName,data); } 01433 inline void read(std::string datasetName, std::string &data) { readAtomic(datasetName,data); } 01434 01435 01436 01437 01438 01439 /** \brief Create a new dataset. 01440 This function can be used to create a dataset filled with a default value, 01441 for example before writing data into it using \ref writeBlock(). 01442 Attention: only atomic datatypes are provided. For spectral data, add an 01443 dimension (case RGB: add one dimension of size 3). 01444 01445 shape determines the dimension and the size of the dataset. 01446 01447 Chunks can be activated by providing a MultiArrayShape as chunkSize. 01448 chunkSize must have equal dimension as array. 01449 01450 Compression can be activated by setting 01451 \code compression = parameter; // 0 < parameter <= 9 01452 \endcode 01453 where 0 stands for no compression and 9 for maximum compression. 01454 01455 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 01456 otherwise it will be interpreted as path relative to the current group. 01457 */ 01458 template<unsigned int N, class T> 01459 inline void createDataset(std::string datasetName, typename MultiArrayShape<N>::type shape, T init = T(), int iChunkSize = 0, int compressionParameter = 0) 01460 { 01461 // make datasetName clean 01462 datasetName = get_absolute_path(datasetName); 01463 01464 typename MultiArrayShape<N>::type chunkSize; 01465 for(int i = 0; i < N; i++){ 01466 chunkSize[i] = iChunkSize; 01467 } 01468 createDataset<N,T>(datasetName, shape, init, chunkSize, compressionParameter); 01469 } 01470 01471 01472 01473 template<unsigned int N, class T> 01474 inline void createDataset(std::string datasetName, typename MultiArrayShape<N>::type shape, T init, typename MultiArrayShape<N>::type chunkSize, int compressionParameter = 0) 01475 { 01476 // make datasetName clean 01477 datasetName = get_absolute_path(datasetName); 01478 01479 std::string groupname = SplitString(datasetName).first(); 01480 std::string setname = SplitString(datasetName).last(); 01481 01482 hid_t parent = openCreateGroup_(groupname); 01483 01484 // delete the dataset if it already exists 01485 deleteDataset_(parent, setname); 01486 01487 // create dataspace 01488 // add an extra dimension in case that the data is non-scalar 01489 HDF5Handle dataspaceHandle; 01490 01491 // invert dimensions to guarantee c-order 01492 hsize_t shape_inv[N]; 01493 for(unsigned int k=0; k<N; ++k) 01494 shape_inv[N-1-k] = shape[k]; 01495 01496 // create dataspace 01497 dataspaceHandle = HDF5Handle(H5Screate_simple(N, shape_inv, NULL), 01498 &H5Sclose, "HDF5File::createDataset(): unable to create dataspace for scalar data."); 01499 01500 // set fill value 01501 HDF5Handle plist ( H5Pcreate(H5P_DATASET_CREATE), &H5Pclose, "HDF5File::createDataset(): unable to create property list." ); 01502 H5Pset_fill_value(plist,detail::getH5DataType<T>(), &init); 01503 01504 // turn off time tagging of datasets by default. 01505 H5Pset_obj_track_times(plist, track_time); 01506 01507 // enable chunks 01508 if(chunkSize[0] > 0) 01509 { 01510 hsize_t cSize [N]; 01511 for(int i = 0; i<N; i++) 01512 { 01513 cSize[i] = chunkSize[N-1-i]; 01514 } 01515 H5Pset_chunk (plist, N, cSize); 01516 } 01517 01518 // enable compression 01519 if(compressionParameter > 0) 01520 { 01521 H5Pset_deflate(plist, compressionParameter); 01522 } 01523 01524 //create the dataset. 01525 HDF5Handle datasetHandle ( H5Dcreate(parent, setname.c_str(), detail::getH5DataType<T>(), dataspaceHandle, H5P_DEFAULT, plist, H5P_DEFAULT), 01526 &H5Dclose, "HDF5File::createDataset(): unable to create dataset."); 01527 if(parent != cGroupHandle_) 01528 H5Gclose(parent); 01529 } 01530 01531 01532 01533 01534 /** \brief Immediately write all data to disk 01535 */ 01536 inline void flushToDisk() 01537 { 01538 H5Fflush(fileHandle_, H5F_SCOPE_GLOBAL); 01539 } 01540 01541 01542 private: 01543 01544 /* Simple extension of std::string for splitting into two parts 01545 * 01546 * Strings (in particular: file/dataset paths) will be split into two 01547 * parts. The split is made at the last occurrence of the delimiter. 01548 * 01549 * For example, "/path/to/some/file" will be split (delimiter = "/") into 01550 * first() = "/path/to/some" and last() = "file". 01551 */ 01552 class SplitString: public std::string { 01553 public: 01554 SplitString(std::string &sstring): std::string(sstring) {}; 01555 01556 // return the part of the string before the delimiter 01557 std::string first(char delimiter = '/') 01558 { 01559 size_t last = find_last_of(delimiter); 01560 if(last == std::string::npos) // delimiter not found --> no first 01561 return ""; 01562 01563 return std::string(begin(), begin()+last+1); 01564 } 01565 01566 // return the part of the string after the delimiter 01567 std::string last(char delimiter = '/') 01568 { 01569 size_t last = find_last_of(delimiter); 01570 if(last == std::string::npos) // delimiter not found --> only last 01571 return std::string(*this); 01572 return std::string(begin()+last+1, end()); 01573 } 01574 }; 01575 01576 public: 01577 01578 /** \brief takes any path and converts it into an absolute path 01579 in the current file. 01580 01581 Elements like "." and ".." are treated as expected. 01582 Links are not supported or resolved. 01583 */ 01584 inline std::string get_absolute_path(std::string path) { 01585 // check for empty input or "." and return the current folder 01586 if(path.length() == 0 || path == "."){ 01587 return currentGroupName_(); 01588 } 01589 01590 std::string str; 01591 // convert to absolute path 01592 if(relativePath_(path)){ 01593 std::string cname = currentGroupName_(); 01594 if (cname == "/") 01595 str = currentGroupName_()+path; 01596 else 01597 str = currentGroupName_()+"/"+path; 01598 }else{ 01599 str = path; 01600 } 01601 01602 // cut out "./" 01603 std::string::size_type startpos = 0; 01604 while(str.find(std::string("./"), startpos) != std::string::npos){ 01605 std::string::size_type pos = str.find(std::string("./"), startpos); 01606 startpos = pos+1; 01607 // only cut if "./" is not part of "../" (see below) 01608 if(str.substr(pos-1,3) != "../"){ 01609 // cut out part of the string 01610 str = str.substr(0,pos) + str.substr(pos+2,str.length()-pos-2); 01611 startpos = pos; 01612 } 01613 } 01614 01615 // cut out pairs of "bla/../" 01616 while(str.find(std::string("..")) != std::string::npos){ 01617 std::string::size_type pos = str.find(std::string("..")); 01618 01619 // find first slash after ".." 01620 std::string::size_type end = str.find("/",pos); 01621 if(end != std::string::npos){ 01622 // also include slash 01623 end++; 01624 }else{ 01625 // no "/" after ".." --> this is a group, add a "/" 01626 str = str + "/"; 01627 end = str.length(); 01628 } 01629 01630 // find first slash before ".." 01631 std::string::size_type prev_slash = str.rfind("/",pos); 01632 // if the root slash is the first before ".." --> Error 01633 vigra_invariant(prev_slash != 0 && prev_slash != std::string::npos, 01634 "Error parsing path: "+str); 01635 // find second slash before ".." 01636 std::string::size_type begin = str.rfind("/",prev_slash-1); 01637 01638 // cut out part of the string 01639 str = str.substr(0,begin+1) + str.substr(end,str.length()-end); 01640 } 01641 01642 return str; 01643 } 01644 01645 private: 01646 01647 /* checks if the given path is a relative path. 01648 */ 01649 inline bool relativePath_(std::string & path) 01650 { 01651 std::string::size_type pos = path.find('/') ; 01652 if(pos == 0) 01653 return false; 01654 01655 return true; 01656 } 01657 01658 01659 01660 01661 /* return the name of the current group 01662 */ 01663 inline std::string currentGroupName_() 01664 { 01665 int len = H5Iget_name(cGroupHandle_,NULL,1000); 01666 ArrayVector<char> name (len+1,0); 01667 H5Iget_name(cGroupHandle_,name.begin(),len+1); 01668 01669 return std::string(name.begin()); 01670 } 01671 01672 01673 01674 01675 /* return the name of the current file 01676 */ 01677 inline std::string fileName_() 01678 { 01679 int len = H5Fget_name(fileHandle_,NULL,1000); 01680 ArrayVector<char> name (len+1,0); 01681 H5Fget_name(fileHandle_,name.begin(),len+1); 01682 01683 return std::string(name.begin()); 01684 } 01685 01686 01687 01688 01689 /* create an empty file and open is 01690 */ 01691 inline hid_t createFile_(std::string filePath, OpenMode mode = Open) 01692 { 01693 // try to open file 01694 FILE * pFile; 01695 pFile = fopen ( filePath.c_str(), "r" ); 01696 hid_t fileId; 01697 01698 // check if opening was successful (= file exists) 01699 if ( pFile == NULL ) 01700 { 01701 fileId = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT); 01702 } 01703 else if(mode == Open) 01704 { 01705 fclose( pFile ); 01706 fileId = H5Fopen(filePath.c_str(), H5F_ACC_RDWR, H5P_DEFAULT); 01707 } 01708 else 01709 { 01710 fclose(pFile); 01711 std::remove(filePath.c_str()); 01712 fileId = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT); 01713 } 01714 return fileId; 01715 } 01716 01717 01718 01719 01720 /* open a group and subgroups. Create if necessary. 01721 */ 01722 inline hid_t openCreateGroup_(std::string groupName) 01723 { 01724 // make groupName clean 01725 groupName = get_absolute_path(groupName); 01726 01727 // open root group 01728 hid_t parent = H5Gopen(fileHandle_, "/", H5P_DEFAULT); 01729 if(groupName == "/") 01730 { 01731 return parent; 01732 } 01733 01734 // remove leading / 01735 groupName = std::string(groupName.begin()+1, groupName.end()); 01736 01737 // check if the groupName has finishing slash 01738 if( groupName.size() != 0 && *groupName.rbegin() != '/') 01739 { 01740 groupName = groupName + '/'; 01741 } 01742 01743 //open or create subgroups one by one 01744 std::string::size_type begin = 0, end = groupName.find('/'); 01745 int ii = 0; 01746 while (end != std::string::npos) 01747 { 01748 std::string group(groupName.begin()+begin, groupName.begin()+end); 01749 hid_t prevParent = parent; 01750 01751 if(H5LTfind_dataset(parent, group.c_str()) == 0) 01752 { 01753 parent = H5Gcreate(prevParent, group.c_str(), H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); 01754 } else { 01755 parent = H5Gopen(prevParent, group.c_str(), H5P_DEFAULT); 01756 } 01757 01758 if(ii != 0) 01759 { 01760 H5Gclose(prevParent); 01761 } 01762 if(parent < 0) 01763 { 01764 return parent; 01765 } 01766 ++ii; 01767 begin = end + 1; 01768 end = groupName.find('/', begin); 01769 } 01770 01771 return parent; 01772 } 01773 01774 01775 01776 01777 /* delete a dataset by unlinking it from the file structure. This does not 01778 delete the data! 01779 */ 01780 inline void deleteDataset_(hid_t parent, std::string datasetName) 01781 { 01782 // delete existing data and create new dataset 01783 if(H5LTfind_dataset(parent, datasetName.c_str())) 01784 { 01785 01786 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6) 01787 if(H5Gunlink(parent, datasetName.c_str()) < 0) 01788 { 01789 vigra_postcondition(false, "HDF5File::deleteDataset_(): Unable to delete existing data."); 01790 } 01791 #else 01792 if(H5Ldelete(parent, datasetName.c_str(), H5P_DEFAULT ) < 0) 01793 { 01794 vigra_postcondition(false, "HDF5File::deleteDataset_(): Unable to delete existing data."); 01795 } 01796 #endif 01797 } 01798 } 01799 01800 01801 01802 01803 /* get the handle of a dataset specified by a string 01804 */ 01805 inline hid_t getDatasetHandle_(std::string datasetName) 01806 { 01807 // make datasetName clean 01808 datasetName = get_absolute_path(datasetName); 01809 01810 std::string groupname = SplitString(datasetName).first(); 01811 std::string setname = SplitString(datasetName).last(); 01812 01813 if (H5Lexists(fileHandle_, datasetName.c_str(), H5P_DEFAULT) <= 0) 01814 { 01815 std::cerr << "HDF5File::getDatasetHandle_(): Dataset '" << datasetName << "' does not exist.\n"; 01816 return -1; 01817 } 01818 01819 //Open parent group 01820 hid_t groupHandle = openCreateGroup_(groupname); 01821 01822 hid_t datasetHandle = H5Dopen(groupHandle, setname.c_str(), H5P_DEFAULT); 01823 01824 if(groupHandle != cGroupHandle_) 01825 { 01826 H5Gclose(groupHandle); 01827 } 01828 01829 //return dataset handle 01830 return datasetHandle; 01831 01832 } 01833 01834 01835 /* get the type of an object specified by a string 01836 */ 01837 H5O_type_t get_object_type_(std::string name) 01838 { 01839 name = get_absolute_path(name); 01840 std::string group_name = SplitString(name).first(); 01841 std::string object_name = SplitString(name).last(); 01842 if (!object_name.size()) 01843 return H5O_TYPE_GROUP; 01844 01845 htri_t exists = H5Lexists(fileHandle_, name.c_str(), H5P_DEFAULT); 01846 vigra_precondition(exists > 0, "HDF5File::get_object_type_(): " 01847 "object \"" + name + "\" " 01848 "not found."); 01849 // open parent group 01850 hid_t group_handle = openCreateGroup_(group_name); 01851 H5O_type_t h5_type = HDF5_get_type(group_handle, name.c_str()); 01852 if (group_handle != cGroupHandle_) 01853 { 01854 H5Gclose(group_handle); 01855 } 01856 return h5_type; 01857 } 01858 01859 /* low-level write function to write vigra MultiArray data as an attribute 01860 */ 01861 template<unsigned int N, class T> 01862 void write_attribute_(std::string name, const std::string & attribute_name, 01863 const MultiArrayView<N, T, UnstridedArrayTag> & array, 01864 const hid_t datatype, const int numBandsOfType) 01865 { 01866 01867 // shape of the array. Add one dimension, if array contains non-scalars. 01868 ArrayVector<hsize_t> shape(N + (numBandsOfType > 1),0); 01869 for(unsigned int i = 0; i < N; i++){ 01870 shape[N-1-i] = array.shape(i); // reverse order 01871 } 01872 01873 if(numBandsOfType > 1) 01874 shape[N] = numBandsOfType; 01875 01876 HDF5Handle dataspace(H5Screate_simple(N + (numBandsOfType > 1), 01877 shape.begin(), NULL), 01878 &H5Sclose, "HDF5File::writeAttribute(): Can not" 01879 " create dataspace."); 01880 01881 std::string errorMessage ("HDF5File::writeAttribute(): can not find " 01882 "object '" + name + "'."); 01883 01884 H5O_type_t h5_type = get_object_type_(name); 01885 bool is_group = h5_type == H5O_TYPE_GROUP; 01886 if (!is_group && h5_type != H5O_TYPE_DATASET) 01887 vigra_precondition(0, "HDF5File::writeAttribute(): object \"" 01888 + name + "\" is neither a group nor a " 01889 "dataset."); 01890 // get parent object handle 01891 HDF5Handle object_handle(is_group 01892 ? openCreateGroup_(name) 01893 : getDatasetHandle_(name), 01894 is_group 01895 ? &H5Gclose 01896 : &H5Dclose, 01897 errorMessage.c_str()); 01898 // create / open attribute 01899 bool exists = existsAttribute(name, attribute_name); 01900 HDF5Handle attributeHandle(exists 01901 ? H5Aopen(object_handle, 01902 attribute_name.c_str(), 01903 H5P_DEFAULT) 01904 : H5Acreate(object_handle, 01905 attribute_name.c_str(), datatype, 01906 dataspace, H5P_DEFAULT, 01907 H5P_DEFAULT), 01908 &H5Aclose, 01909 "HDF5File::writeAttribute(): Can not create" 01910 " attribute."); 01911 01912 // Write the data to the HDF5 object 01913 H5Awrite(attributeHandle, datatype, array.data()); 01914 } 01915 01916 01917 01918 01919 /* Write single value attribute 01920 This function allows to write data of atomic datatypes (int, long, double) 01921 as an attribute in the HDF5 file. So it is not necessary to create a MultiArray 01922 of size 1 to write a single number. 01923 */ 01924 template<class T> 01925 inline void writeAtomicAttribute(std::string datasetName, std::string attributeName, const T data) 01926 { 01927 // make datasetName clean 01928 datasetName = get_absolute_path(datasetName); 01929 01930 typename MultiArrayShape<1>::type chunkSize; 01931 chunkSize[0] = 0; 01932 MultiArray<1,T> array(MultiArrayShape<1>::type(1)); 01933 array[0] = data; 01934 write_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), 1); 01935 } 01936 01937 01938 01939 01940 /* low-level read function to write vigra MultiArray data from attributes 01941 */ 01942 template<unsigned int N, class T> 01943 inline void read_attribute_(std::string datasetName, std::string attributeName, MultiArrayView<N, T, UnstridedArrayTag> array, const hid_t datatype, const int numBandsOfType) 01944 { 01945 std::string dataset_path = get_absolute_path(datasetName); 01946 // open Attribute handle 01947 std::string message = "Error: could not get handle for attribute '"+attributeName+"'' of object '"+dataset_path+"'."; 01948 HDF5Handle attr_handle (H5Aopen_by_name(fileHandle_,dataset_path.c_str(),attributeName.c_str(),H5P_DEFAULT,H5P_DEFAULT),&H5Aclose, message.c_str()); 01949 01950 // get Attribute dataspace 01951 message = "Error: could not get dataspace for attribute '"+attributeName+"'' of object '"+dataset_path+"'."; 01952 HDF5Handle attr_dataspace_handle (H5Aget_space(attr_handle),&H5Sclose,message.c_str()); 01953 01954 // obtain Attribute shape 01955 int dims = H5Sget_simple_extent_ndims(attr_dataspace_handle); 01956 ArrayVector<hsize_t> shape_inv(dims); 01957 H5Sget_simple_extent_dims(attr_dataspace_handle, shape_inv.data(), NULL); 01958 01959 // invert the dimensions to guarantee c-order 01960 ArrayVector<hsize_t> dimshape(dims); 01961 for(ArrayVector<hsize_t>::size_type i=0; i<shape_inv.size(); i++) { 01962 dimshape[i] = shape_inv[dims-1-i]; 01963 } 01964 01965 int offset = (numBandsOfType > 1); 01966 message = "Error: Array dimension disagrees with dataset dimension."; 01967 // the object in the HDF5 file may have one additional dimension which we then interpret as the pixel type bands 01968 vigra_precondition( ( (N + offset ) == MultiArrayIndex(dims)), message); 01969 01970 typename MultiArrayShape<N>::type shape; 01971 for(int k=offset; k< MultiArrayIndex(dims); ++k) { 01972 shape[k-offset] = MultiArrayIndex(dimshape[k]); 01973 } 01974 01975 message = "Error: Array shape disagrees with dataset shape"; 01976 vigra_precondition(shape == array.shape(),message); 01977 01978 // simply read in the data as is 01979 H5Aread( attr_handle, datatype, array.data()); 01980 } 01981 01982 01983 01984 01985 /* Read a single value attribute. 01986 This functions allows to read a single value attribute of atomic datatype (int, long, double) 01987 from the HDF5 file. So it is not necessary to create a MultiArray 01988 of size 1 to read a single number. 01989 */ 01990 template<class T> 01991 inline void readAtomicAttribute(std::string datasetName, std::string attributeName, T & data) 01992 { 01993 // make datasetName clean 01994 datasetName = get_absolute_path(datasetName); 01995 01996 MultiArray<1,T> array(MultiArrayShape<1>::type(1)); 01997 read_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), 1); 01998 data = array[0]; 01999 } 02000 02001 02002 02003 inline void readAtomicAttribute(std::string datasetName, std::string attributeName, std::string & data) 02004 { 02005 // make datasetName clean 02006 datasetName = get_absolute_path(datasetName); 02007 02008 MultiArray<1,const char *> array(MultiArrayShape<1>::type(1)); 02009 read_attribute_(datasetName, attributeName, array, detail::getH5DataType<const char *>(), 1); 02010 data = std::string(array[0]); 02011 } 02012 02013 02014 02015 02016 /* low-level write function to write vigra unstrided MultiArray data 02017 */ 02018 template<unsigned int N, class T> 02019 inline 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) 02020 { 02021 std::string groupname = SplitString(datasetName).first(); 02022 std::string setname = SplitString(datasetName).last(); 02023 02024 // shape of the array. Add one dimension, if array contains non-scalars. 02025 ArrayVector<hsize_t> shape(N + (numBandsOfType > 1),0); 02026 for(unsigned int i = 0; i < N; i++){ 02027 shape[N-1-i] = array.shape(i); // reverse order 02028 } 02029 02030 if(numBandsOfType > 1) 02031 shape[N] = numBandsOfType; 02032 02033 HDF5Handle dataspace ( H5Screate_simple(N + (numBandsOfType > 1), shape.begin(), NULL), &H5Sclose, "HDF5File::write(): Can not create dataspace."); 02034 02035 // create and open group: 02036 std::string errorMessage ("HDF5File::write(): can not create group '" + groupname + "'."); 02037 hid_t groupHandle = openCreateGroup_(groupname); 02038 if(groupHandle <= 0) 02039 { 02040 std::cerr << errorMessage << "\n"; 02041 } 02042 02043 // delete dataset, if it already exists 02044 deleteDataset_(groupHandle, setname.c_str()); 02045 02046 // set up properties list 02047 HDF5Handle plist ( H5Pcreate(H5P_DATASET_CREATE), &H5Pclose, "HDF5File::write(): unable to create property list." ); 02048 02049 // turn off time tagging of datasets by default. 02050 H5Pset_obj_track_times(plist, track_time); 02051 02052 // enable chunks 02053 if(chunkSize[0] > 0) 02054 { 02055 ArrayVector<hsize_t> cSize(N + (numBandsOfType > 1),0); 02056 for(unsigned int i = 0; i<N; i++) 02057 { 02058 cSize[i] = chunkSize[N-1-i]; 02059 } 02060 if(numBandsOfType > 1) 02061 cSize[N] = numBandsOfType; 02062 02063 H5Pset_chunk (plist, N + (numBandsOfType > 1), cSize.begin()); 02064 } 02065 02066 // enable compression 02067 if(compressionParameter > 0) 02068 { 02069 H5Pset_deflate(plist, compressionParameter); 02070 } 02071 02072 // create dataset 02073 HDF5Handle datasetHandle (H5Dcreate(groupHandle, setname.c_str(), datatype, dataspace,H5P_DEFAULT, plist, H5P_DEFAULT), &H5Dclose, "HDF5File::write(): Can not create dataset."); 02074 02075 // Write the data to the HDF5 dataset as is 02076 herr_t write_status = H5Dwrite(datasetHandle, datatype, H5S_ALL, 02077 H5S_ALL, H5P_DEFAULT, array.data()); 02078 vigra_precondition(write_status >= 0, "HDF5File::write_(): write to " 02079 "dataset \"" + datasetName + "\" " 02080 "failed."); 02081 02082 if(groupHandle != cGroupHandle_) 02083 { 02084 H5Gclose(groupHandle); 02085 } 02086 } 02087 02088 02089 02090 02091 /* Write single value as dataset. 02092 This functions allows to write data of atomic datatypes (int, long, double) 02093 as a dataset in the HDF5 file. So it is not necessary to create a MultiArray 02094 of size 1 to write a single number. 02095 02096 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 02097 otherwise it will be interpreted as path relative to the current group. 02098 */ 02099 template<class T> 02100 inline void writeAtomic(std::string datasetName, const T data) 02101 { 02102 // make datasetName clean 02103 datasetName = get_absolute_path(datasetName); 02104 02105 typename MultiArrayShape<1>::type chunkSize; 02106 chunkSize[0] = 0; 02107 MultiArray<1,T> array(MultiArrayShape<1>::type(1)); 02108 array[0] = data; 02109 write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize,0); 02110 } 02111 02112 02113 02114 02115 /* low-level read function to read vigra unstrided MultiArray data 02116 */ 02117 template<unsigned int N, class T> 02118 inline void read_(std::string datasetName, MultiArrayView<N, T, UnstridedArrayTag> array, const hid_t datatype, const int numBandsOfType) 02119 { 02120 //Prepare to read without using HDF5ImportInfo 02121 ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName); 02122 hssize_t dimensions = getDatasetDimensions(datasetName); 02123 02124 std::string errorMessage ("HDF5File::read(): Unable to open dataset '" + datasetName + "'."); 02125 HDF5Handle datasetHandle (getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str()); 02126 02127 int offset = (numBandsOfType > 1); 02128 02129 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 02130 "HDF5File::read(): Array dimension disagrees with dataset dimension."); 02131 02132 typename MultiArrayShape<N>::type shape; 02133 for(int k=offset; k< MultiArrayIndex(dimensions); ++k) { 02134 shape[k-offset] = MultiArrayIndex(dimshape[k]); 02135 } 02136 02137 vigra_precondition(shape == array.shape(), 02138 "HDF5File::read(): Array shape disagrees with dataset shape."); 02139 02140 // simply read in the data as is 02141 H5Dread( datasetHandle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data() ); // .data() possible since void pointer! 02142 } 02143 02144 02145 02146 02147 /* Read a single value. 02148 This functions allows to read a single datum of atomic datatype (int, long, double) 02149 from the HDF5 file. So it is not necessary to create a MultiArray 02150 of size 1 to read a single number. 02151 02152 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 02153 otherwise it will be interpreted as path relative to the current group. 02154 */ 02155 template<class T> 02156 inline void readAtomic(std::string datasetName, T & data) 02157 { 02158 // make datasetName clean 02159 datasetName = get_absolute_path(datasetName); 02160 02161 MultiArray<1,T> array(MultiArrayShape<1>::type(1)); 02162 read_(datasetName, array, detail::getH5DataType<T>(), 1); 02163 data = array[0]; 02164 } 02165 inline void readAtomic(std::string datasetName, std::string & data) 02166 { 02167 // make datasetName clean 02168 datasetName = get_absolute_path(datasetName); 02169 02170 MultiArray<1,const char *> array(MultiArrayShape<1>::type(1)); 02171 read_(datasetName, array, detail::getH5DataType<const char *>(), 1); 02172 data = std::string(array[0]); 02173 } 02174 02175 02176 02177 02178 /* low-level write function to write vigra unstrided MultiArray data into a sub-block of a dataset 02179 */ 02180 template<unsigned int N, class T> 02181 inline void writeBlock_(std::string datasetName, typename MultiArrayShape<N>::type &blockOffset, const MultiArrayView<N, T, UnstridedArrayTag> & array, const hid_t datatype, const int numBandsOfType) 02182 { 02183 // open dataset if it exists 02184 std::string errorMessage = "HDF5File::writeBlock(): Error opening dataset '" + datasetName + "'."; 02185 HDF5Handle datasetHandle (getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str()); 02186 02187 02188 // hyperslab parameters for position, size, ... 02189 hsize_t boffset [N]; 02190 hsize_t bshape [N]; 02191 hsize_t bones [N]; 02192 02193 for(int i = 0; i < N; i++){ 02194 boffset[i] = blockOffset[N-1-i]; 02195 bshape[i] = array.size(N-1-i); 02196 bones[i] = 1; 02197 } 02198 02199 // create a target dataspace in memory with the shape of the desired block 02200 HDF5Handle memspace_handle (H5Screate_simple(N,bshape,NULL),&H5Sclose,"Unable to get origin dataspace"); 02201 02202 // get file dataspace and select the desired block 02203 HDF5Handle dataspaceHandle (H5Dget_space(datasetHandle),&H5Sclose,"Unable to create target dataspace"); 02204 H5Sselect_hyperslab(dataspaceHandle, H5S_SELECT_SET, boffset, bones, bones, bshape); 02205 02206 // Write the data to the HDF5 dataset as is 02207 H5Dwrite( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, array.data()); // .data() possible since void pointer! 02208 } 02209 02210 02211 02212 02213 /* low-level read function to read vigra unstrided MultiArray data from a sub-block of a dataset 02214 */ 02215 template<unsigned int N, class T> 02216 inline 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) 02217 { 02218 //Prepare to read without using HDF5ImportInfo 02219 //ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName) ; 02220 hssize_t dimensions = getDatasetDimensions(datasetName); 02221 02222 std::string errorMessage ("HDF5File::readBlock(): Unable to open dataset '" + datasetName + "'."); 02223 HDF5Handle datasetHandle (getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str()); 02224 02225 02226 int offset = (numBandsOfType > 1); 02227 02228 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 02229 "readHDF5_block(): Array dimension disagrees with data dimension."); 02230 02231 vigra_precondition(blockShape == array.shape(), 02232 "readHDF5_block(): Array shape disagrees with block size."); 02233 02234 // hyperslab parameters for position, size, ... 02235 hsize_t boffset [N]; 02236 hsize_t bshape [N]; 02237 hsize_t bones [N]; 02238 02239 for(int i = 0; i < N; i++){ 02240 // vigra and hdf5 use different indexing 02241 boffset[i] = blockOffset[N-1-i]; 02242 //bshape[i] = blockShape[i]; 02243 bshape[i] = blockShape[N-1-i]; 02244 //boffset[i] = blockOffset[N-1-i]; 02245 bones[i] = 1; 02246 } 02247 02248 // create a target dataspace in memory with the shape of the desired block 02249 HDF5Handle memspace_handle (H5Screate_simple(N,bshape,NULL),&H5Sclose,"Unable to create target dataspace"); 02250 02251 // get file dataspace and select the desired block 02252 HDF5Handle dataspaceHandle (H5Dget_space(datasetHandle),&H5Sclose,"Unable to get dataspace"); 02253 H5Sselect_hyperslab(dataspaceHandle, H5S_SELECT_SET, boffset, bones, bones, bshape); 02254 02255 // now read the data 02256 H5Dread( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, array.data() ); // .data() possible since void pointer! 02257 } 02258 02259 }; /* class HDF5File */ 02260 02261 namespace detail { 02262 02263 template <class Shape> 02264 inline void 02265 selectHyperslabs(HDF5Handle & mid1, HDF5Handle & mid2, Shape const & shape, int & counter, const int elements, const int numBandsOfType) 02266 { 02267 // select hyperslab in HDF5 file 02268 hsize_t shapeHDF5[2]; 02269 shapeHDF5[0] = 1; 02270 shapeHDF5[1] = elements; 02271 hsize_t startHDF5[2]; 02272 startHDF5[0] = 0; 02273 startHDF5[1] = counter * numBandsOfType * shape[0]; // we have to reserve space for the pixel type channel(s) 02274 hsize_t strideHDF5[2]; 02275 strideHDF5[0] = 1; 02276 strideHDF5[1] = 1; 02277 hsize_t countHDF5[2]; 02278 countHDF5[0] = 1; 02279 countHDF5[1] = numBandsOfType * shape[0]; 02280 hsize_t blockHDF5[2]; 02281 blockHDF5[0] = 1; 02282 blockHDF5[1] = 1; 02283 mid1 = HDF5Handle(H5Screate_simple(2, shapeHDF5, NULL), 02284 &H5Sclose, "unable to create hyperslabs."); 02285 H5Sselect_hyperslab(mid1, H5S_SELECT_SET, startHDF5, strideHDF5, countHDF5, blockHDF5); 02286 // select hyperslab in input data object 02287 hsize_t shapeData[2]; 02288 shapeData[0] = 1; 02289 shapeData[1] = numBandsOfType * shape[0]; 02290 hsize_t startData[2]; 02291 startData[0] = 0; 02292 startData[1] = 0; 02293 hsize_t strideData[2]; 02294 strideData[0] = 1; 02295 strideData[1] = 1; 02296 hsize_t countData[2]; 02297 countData[0] = 1; 02298 countData[1] = numBandsOfType * shape[0]; 02299 hsize_t blockData[2]; 02300 blockData[0] = 1; 02301 blockData[1] = 1; 02302 mid2 = HDF5Handle(H5Screate_simple(2, shapeData, NULL), 02303 &H5Sclose, "unable to create hyperslabs."); 02304 H5Sselect_hyperslab(mid2, H5S_SELECT_SET, startData, strideData, countData, blockData); 02305 } 02306 02307 template <class DestIterator, class Shape, class T> 02308 inline void 02309 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>) 02310 { 02311 HDF5Handle mid1, mid2; 02312 02313 // select hyperslabs 02314 selectHyperslabs(mid1, mid2, shape, counter, elements, numBandsOfType); 02315 02316 // read from hdf5 02317 herr_t read_status = H5Dread(dataset_id, datatype, mid2, mid1, H5P_DEFAULT, buffer.data()); 02318 vigra_precondition(read_status >= 0, "readHDF5Impl(): read from dataset failed."); 02319 02320 // increase counter 02321 counter++; 02322 02323 //std::cout << "numBandsOfType: " << numBandsOfType << std::endl; 02324 DestIterator dend = d + shape[0]; 02325 int k = 0; 02326 for(; d < dend; ++d, k++) 02327 { 02328 *d = buffer[k]; 02329 //std::cout << buffer[k] << "| "; 02330 } 02331 02332 } 02333 02334 template <class DestIterator, class Shape, class T, int N> 02335 void 02336 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>) 02337 { 02338 DestIterator dend = d + shape[N]; 02339 for(; d < dend; ++d) 02340 { 02341 readHDF5Impl(d.begin(), shape, dataset_id, datatype, buffer, counter, elements, numBandsOfType, MetaInt<N-1>()); 02342 } 02343 } 02344 02345 } // namespace detail 02346 02347 /** \brief Read the data specified by the given \ref vigra::HDF5ImportInfo object 02348 and write the into the given 'array'. 02349 02350 The array must have the correct number of dimensions and shape for the dataset 02351 represented by 'info'. When the element type of 'array' differs from the stored element 02352 type, HDF5 will convert the type on the fly (except when the HDF5 version is 1.6 or below, 02353 in which case an error will result). Multi-channel element types (i.e. \ref vigra::RGBValue 02354 and \ref vigra::TinyVector) are recognized and handled correctly. 02355 02356 <b> Declaration:</b> 02357 02358 \code 02359 namespace vigra { 02360 template<unsigned int N, class T, class StrideTag> 02361 void 02362 readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, StrideTag> array); 02363 } 02364 \endcode 02365 02366 <b> Usage:</b> 02367 02368 <b>\#include</b> <vigra/hdf5impex.hxx><br> 02369 Namespace: vigra 02370 02371 \code 02372 02373 HDF5ImportInfo info(filename, dataset_name); 02374 vigra_precondition(info.numDimensions() == 3, "Dataset must be 3-dimensional."); 02375 02376 MultiArrayShape<3>::type shape(info.shape().begin()); 02377 MultiArray<3, int> array(shape); 02378 02379 readHDF5(info, array); 02380 \endcode 02381 */ 02382 doxygen_overloaded_function(template <...> void readHDF5) 02383 02384 // scalar and unstrided target multi array 02385 template<unsigned int N, class T> 02386 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, UnstridedArrayTag> array) // scalar 02387 { 02388 readHDF5(info, array, detail::getH5DataType<T>(), 1); 02389 } 02390 02391 // non-scalar (TinyVector) and unstrided target multi array 02392 template<unsigned int N, class T, int SIZE> 02393 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> array) 02394 { 02395 readHDF5(info, array, detail::getH5DataType<T>(), SIZE); 02396 } 02397 02398 // non-scalar (RGBValue) and unstrided target multi array 02399 template<unsigned int N, class T> 02400 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> array) 02401 { 02402 readHDF5(info, array, detail::getH5DataType<T>(), 3); 02403 } 02404 02405 // unstrided target multi array 02406 template<unsigned int N, class T> 02407 void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, UnstridedArrayTag> array, const hid_t datatype, const int numBandsOfType) 02408 { 02409 int offset = (numBandsOfType > 1); 02410 02411 //std::cout << "offset: " << offset << ", N: " << N << ", dims: " << info.numDimensions() << std::endl; 02412 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 02413 "readHDF5(): Array dimension disagrees with HDF5ImportInfo.numDimensions()."); 02414 02415 typename MultiArrayShape<N>::type shape; 02416 for(int k=offset; k<info.numDimensions(); ++k) { 02417 shape[k-offset] = info.shapeOfDimension(k); 02418 } 02419 02420 vigra_precondition(shape == array.shape(), 02421 "readHDF5(): Array shape disagrees with HDF5ImportInfo."); 02422 02423 // simply read in the data as is 02424 H5Dread( info.getDatasetHandle(), datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data() ); // .data() possible since void pointer! 02425 } 02426 02427 // scalar and strided target multi array 02428 template<unsigned int N, class T> 02429 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, StridedArrayTag> array) // scalar 02430 { 02431 readHDF5(info, array, detail::getH5DataType<T>(), 1); 02432 } 02433 02434 // non-scalar (TinyVector) and strided target multi array 02435 template<unsigned int N, class T, int SIZE> 02436 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, TinyVector<T, SIZE>, StridedArrayTag> array) 02437 { 02438 readHDF5(info, array, detail::getH5DataType<T>(), SIZE); 02439 } 02440 02441 // non-scalar (RGBValue) and strided target multi array 02442 template<unsigned int N, class T> 02443 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, RGBValue<T>, StridedArrayTag> array) 02444 { 02445 readHDF5(info, array, detail::getH5DataType<T>(), 3); 02446 } 02447 02448 // strided target multi array 02449 template<unsigned int N, class T> 02450 void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, StridedArrayTag> array, const hid_t datatype, const int numBandsOfType) 02451 { 02452 int offset = (numBandsOfType > 1); 02453 02454 //std::cout << "offset: " << offset << ", N: " << N << ", dims: " << info.numDimensions() << std::endl; 02455 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 02456 "readHDF5(): Array dimension disagrees with HDF5ImportInfo.numDimensions()."); 02457 02458 typename MultiArrayShape<N>::type shape; 02459 for(int k=offset; k<info.numDimensions(); ++k) { 02460 shape[k-offset] = info.shapeOfDimension(k); 02461 } 02462 02463 vigra_precondition(shape == array.shape(), 02464 "readHDF5(): Array shape disagrees with HDF5ImportInfo."); 02465 02466 //Get the data 02467 int counter = 0; 02468 int elements = numBandsOfType; 02469 for(unsigned int i=0;i<N;++i) 02470 elements *= shape[i]; 02471 ArrayVector<T> buffer(shape[0]); 02472 detail::readHDF5Impl(array.traverser_begin(), shape, info.getDatasetHandle(), datatype, buffer, counter, elements, numBandsOfType, vigra::MetaInt<N-1>()); 02473 } 02474 02475 inline hid_t openGroup(hid_t parent, std::string group_name) 02476 { 02477 //std::cout << group_name << std::endl; 02478 size_t last_slash = group_name.find_last_of('/'); 02479 if (last_slash == std::string::npos || last_slash != group_name.size() - 1) 02480 group_name = group_name + '/'; 02481 std::string::size_type begin = 0, end = group_name.find('/'); 02482 int ii = 0; 02483 while (end != std::string::npos) 02484 { 02485 std::string group(group_name.begin()+begin, group_name.begin()+end); 02486 hid_t prev_parent = parent; 02487 parent = H5Gopen(prev_parent, group.c_str(), H5P_DEFAULT); 02488 02489 if(ii != 0) H5Gclose(prev_parent); 02490 if(parent < 0) return parent; 02491 ++ii; 02492 begin = end + 1; 02493 end = group_name.find('/', begin); 02494 } 02495 return parent; 02496 } 02497 02498 inline hid_t createGroup(hid_t parent, std::string group_name) 02499 { 02500 if(group_name.size() == 0 ||*group_name.rbegin() != '/') 02501 group_name = group_name + '/'; 02502 if(group_name == "/") 02503 return H5Gopen(parent, group_name.c_str(), H5P_DEFAULT); 02504 02505 std::string::size_type begin = 0, end = group_name.find('/'); 02506 int ii = 0; 02507 while (end != std::string::npos) 02508 { 02509 std::string group(group_name.begin()+begin, group_name.begin()+end); 02510 hid_t prev_parent = parent; 02511 02512 if(H5LTfind_dataset(parent, group.c_str()) == 0) 02513 { 02514 parent = H5Gcreate(prev_parent, group.c_str(), H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); 02515 } else { 02516 parent = H5Gopen(prev_parent, group.c_str(), H5P_DEFAULT); 02517 } 02518 02519 if(ii != 0) H5Gclose(prev_parent); 02520 if(parent < 0) return parent; 02521 ++ii; 02522 begin = end + 1; 02523 end = group_name.find('/', begin); 02524 } 02525 return parent; 02526 } 02527 02528 inline void deleteDataset(hid_t parent, std::string dataset_name) 02529 { 02530 // delete existing data and create new dataset 02531 if(H5LTfind_dataset(parent, dataset_name.c_str())) 02532 { 02533 //std::cout << "dataset already exists" << std::endl; 02534 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6) 02535 if(H5Gunlink(parent, dataset_name.c_str()) < 0) 02536 { 02537 vigra_postcondition(false, "writeToHDF5File(): Unable to delete existing data."); 02538 } 02539 #else 02540 if(H5Ldelete(parent, dataset_name.c_str(), H5P_DEFAULT ) < 0) 02541 { 02542 vigra_postcondition(false, "createDataset(): Unable to delete existing data."); 02543 } 02544 #endif 02545 } 02546 } 02547 02548 inline hid_t createFile(std::string filePath, bool append_ = true) 02549 { 02550 FILE * pFile; 02551 pFile = fopen ( filePath.c_str(), "r" ); 02552 hid_t file_id; 02553 if ( pFile == NULL ) 02554 { 02555 file_id = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT); 02556 } 02557 else if(append_) 02558 { 02559 fclose( pFile ); 02560 file_id = H5Fopen(filePath.c_str(), H5F_ACC_RDWR, H5P_DEFAULT); 02561 } 02562 else 02563 { 02564 fclose(pFile); 02565 std::remove(filePath.c_str()); 02566 file_id = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT); 02567 } 02568 return file_id; 02569 } 02570 02571 template<unsigned int N, class T, class Tag> 02572 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) 02573 { 02574 std::string path_name(pathInFile), group_name, data_set_name, message; 02575 std::string::size_type delimiter = path_name.rfind('/'); 02576 02577 //create or open file 02578 file_handle = HDF5Handle(createFile(filePath), &H5Fclose, 02579 "createDataset(): unable to open output file."); 02580 02581 // get the groupname and the filename 02582 if(delimiter == std::string::npos) 02583 { 02584 group_name = "/"; 02585 data_set_name = path_name; 02586 } 02587 else 02588 { 02589 group_name = std::string(path_name.begin(), path_name.begin()+delimiter); 02590 data_set_name = std::string(path_name.begin()+delimiter+1, path_name.end()); 02591 } 02592 02593 // create all groups 02594 HDF5Handle group(createGroup(file_handle, group_name), &H5Gclose, 02595 "createDataset(): Unable to create and open group. generic v"); 02596 02597 // delete the dataset if it already exists 02598 deleteDataset(group, data_set_name); 02599 02600 // create dataspace 02601 // add an extra dimension in case that the data is non-scalar 02602 HDF5Handle dataspace_handle; 02603 if(numBandsOfType > 1) { 02604 // invert dimensions to guarantee c-order 02605 hsize_t shape_inv[N+1]; // one additional dimension for pixel type channel(s) 02606 for(unsigned int k=0; k<N; ++k) { 02607 shape_inv[N-1-k] = array.shape(k); // the channels (eg of an RGB image) are represented by the first dimension (before inversion) 02608 //std::cout << shape_inv[N-k] << " (" << N << ")"; 02609 } 02610 shape_inv[N] = numBandsOfType; 02611 02612 // create dataspace 02613 dataspace_handle = HDF5Handle(H5Screate_simple(N+1, shape_inv, NULL), 02614 &H5Sclose, "createDataset(): unable to create dataspace for non-scalar data."); 02615 } else { 02616 // invert dimensions to guarantee c-order 02617 hsize_t shape_inv[N]; 02618 for(unsigned int k=0; k<N; ++k) 02619 shape_inv[N-1-k] = array.shape(k); 02620 02621 // create dataspace 02622 dataspace_handle = HDF5Handle(H5Screate_simple(N, shape_inv, NULL), 02623 &H5Sclose, "createDataset(): unable to create dataspace for scalar data."); 02624 } 02625 02626 //alloc memory for dataset. 02627 dataset_handle = HDF5Handle(H5Dcreate(group, 02628 data_set_name.c_str(), 02629 datatype, 02630 dataspace_handle, 02631 H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT), 02632 &H5Dclose, "createDataset(): unable to create dataset."); 02633 } 02634 02635 02636 02637 namespace detail { 02638 02639 template <class DestIterator, class Shape, class T> 02640 inline void 02641 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>) 02642 { 02643 DestIterator dend = d + (typename DestIterator::difference_type)shape[0]; 02644 int k = 0; 02645 //std::cout << "new:" << std::endl; 02646 for(; d < dend; ++d, k++) 02647 { 02648 buffer[k] = *d; 02649 //std::cout << buffer[k] << " "; 02650 } 02651 //std::cout << std::endl; 02652 HDF5Handle mid1, mid2; 02653 02654 // select hyperslabs 02655 selectHyperslabs(mid1, mid2, shape, counter, elements, numBandsOfType); 02656 02657 // write to hdf5 02658 H5Dwrite(dataset_id, datatype, mid2, mid1, H5P_DEFAULT, buffer.data()); 02659 // increase counter 02660 counter++; 02661 } 02662 02663 template <class DestIterator, class Shape, class T, int N> 02664 void 02665 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>) 02666 { 02667 DestIterator dend = d + (typename DestIterator::difference_type)shape[N]; 02668 for(; d < dend; ++d) 02669 { 02670 writeHDF5Impl(d.begin(), shape, dataset_id, datatype, buffer, counter, elements, numBandsOfType, MetaInt<N-1>()); 02671 } 02672 } 02673 02674 } // namespace detail 02675 02676 /** \brief Store array data in an HDF5 file. 02677 02678 The number of dimensions, shape and element type of the stored dataset is automatically 02679 determined from the properties of the given \a array. Strided arrays are stored in an 02680 unstrided way, i.e. in contiguous scan-order. Multi-channel element types 02681 (i.e. \ref vigra::RGBValue and \ref vigra::TinyVector) are recognized and handled correctly 02682 (in particular, the will form the innermost dimension of the stored dataset). 02683 \a pathInFile may contain '/'-separated group names, but must end with the name 02684 of the dataset to be created. 02685 02686 <b> Declaration:</b> 02687 02688 \code 02689 namespace vigra { 02690 template<unsigned int N, class T, class StrideTag> 02691 void 02692 writeHDF5(const char* filePath, const char* pathInFile, 02693 MultiArrayView<N, T, StrideTag>const & array); 02694 } 02695 \endcode 02696 02697 <b> Usage:</b> 02698 02699 <b>\#include</b> <vigra/hdf5impex.hxx><br> 02700 Namespace: vigra 02701 02702 \code 02703 MultiArrayShape<3>::type shape(100, 200, 20); 02704 MultiArray<3, int> array(shape); 02705 ... // fill array with data 02706 02707 writeHDF5("mydata.h5", "/group1/my_dataset", array); 02708 \endcode 02709 */ 02710 doxygen_overloaded_function(template <...> void writeHDF5) 02711 02712 // scalar and unstrided multi arrays 02713 template<unsigned int N, class T> 02714 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, UnstridedArrayTag> & array) // scalar 02715 { 02716 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 1); 02717 } 02718 02719 // non-scalar (TinyVector) and unstrided multi arrays 02720 template<unsigned int N, class T, int SIZE> 02721 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array) 02722 { 02723 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), SIZE); 02724 } 02725 02726 // non-scalar (RGBValue) and unstrided multi arrays 02727 template<unsigned int N, class T> 02728 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array) 02729 { 02730 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 3); 02731 } 02732 02733 // non-scalar (FFTWComplex) and unstrided multi arrays 02734 template<unsigned int N, class T> 02735 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, FFTWComplex<T>, UnstridedArrayTag> & array) 02736 { 02737 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 2); 02738 } 02739 02740 // unstrided multi arrays 02741 template<unsigned int N, class T> 02742 void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, UnstridedArrayTag> & array, const hid_t datatype, const int numBandsOfType) 02743 { 02744 HDF5Handle file_handle; 02745 HDF5Handle dataset_handle; 02746 createDataset(filePath, pathInFile, array, datatype, numBandsOfType, file_handle, dataset_handle); 02747 02748 // Write the data to the HDF5 dataset as is 02749 H5Dwrite( dataset_handle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data()); // .data() possible since void pointer! 02750 02751 H5Fflush(file_handle, H5F_SCOPE_GLOBAL); 02752 } 02753 02754 02755 // scalar and strided multi arrays 02756 template<unsigned int N, class T> 02757 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, StridedArrayTag> & array) // scalar 02758 { 02759 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 1); 02760 } 02761 02762 // non-scalar (TinyVector) and strided multi arrays 02763 template<unsigned int N, class T, int SIZE> 02764 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, TinyVector<T, SIZE>, StridedArrayTag> & array) 02765 { 02766 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), SIZE); 02767 } 02768 02769 // non-scalar (RGBValue) and strided multi arrays 02770 template<unsigned int N, class T> 02771 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, RGBValue<T>, StridedArrayTag> & array) 02772 { 02773 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 3); 02774 } 02775 02776 // non-scalar (RGBValue) and strided multi arrays 02777 template<unsigned int N, class T> 02778 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, FFTWComplex<T>, StridedArrayTag> & array) 02779 { 02780 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 2); 02781 } 02782 02783 // strided multi arrays 02784 template<unsigned int N, class T> 02785 void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, StridedArrayTag> & array, const hid_t datatype, const int numBandsOfType) 02786 { 02787 HDF5Handle file_handle; 02788 HDF5Handle dataset_handle; 02789 createDataset(filePath, pathInFile, array, datatype, numBandsOfType, file_handle, dataset_handle); 02790 02791 vigra::TinyVector<int,N> shape; 02792 vigra::TinyVector<int,N> stride; 02793 int elements = numBandsOfType; 02794 for(unsigned int k=0; k<N; ++k) 02795 { 02796 shape[k] = array.shape(k); 02797 stride[k] = array.stride(k); 02798 elements *= (int)shape[k]; 02799 } 02800 int counter = 0; 02801 02802 ArrayVector<T> buffer((int)array.shape(0)); 02803 detail::writeHDF5Impl(array.traverser_begin(), shape, dataset_handle, datatype, buffer, counter, elements, numBandsOfType, vigra::MetaInt<N-1>()); 02804 02805 H5Fflush(file_handle, H5F_SCOPE_GLOBAL); 02806 02807 } 02808 02809 namespace detail 02810 { 02811 struct MaxSizeFnc 02812 { 02813 size_t size; 02814 02815 MaxSizeFnc() 02816 : size(0) 02817 {} 02818 02819 void operator()(std::string const & in) 02820 { 02821 size = in.size() > size ? 02822 in.size() : 02823 size; 02824 } 02825 }; 02826 } 02827 02828 02829 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR == 8) || DOXYGEN 02830 /** Write a numeric MultiArray as an attribute with name \a name 02831 of the dataset specified by the handle \a loc. 02832 02833 <b>\#include</b> <vigra/hdf5impex.hxx><br> 02834 Namespace: vigra 02835 */ 02836 template<size_t N, class T, class C> 02837 void writeHDF5Attr(hid_t loc, 02838 const char* name, 02839 MultiArrayView<N, T, C> const & array) 02840 { 02841 if(H5Aexists(loc, name) > 0) 02842 H5Adelete(loc, name); 02843 02844 ArrayVector<hsize_t> shape(array.shape().begin(), 02845 array.shape().end()); 02846 HDF5Handle 02847 dataspace_handle(H5Screate_simple(N, shape.data(), NULL), 02848 &H5Sclose, 02849 "writeToHDF5File(): unable to create dataspace."); 02850 02851 HDF5Handle attr(H5Acreate(loc, 02852 name, 02853 detail::getH5DataType<T>(), 02854 dataspace_handle, 02855 H5P_DEFAULT ,H5P_DEFAULT ), 02856 &H5Aclose, 02857 "writeHDF5Attr: unable to create Attribute"); 02858 02859 //copy data - since attributes are small - who cares! 02860 ArrayVector<T> buffer; 02861 for(int ii = 0; ii < array.size(); ++ii) 02862 buffer.push_back(array[ii]); 02863 H5Awrite(attr, detail::getH5DataType<T>(), buffer.data()); 02864 } 02865 02866 /** Write a string MultiArray as an attribute with name \a name 02867 of the dataset specified by the handle \a loc. 02868 02869 <b>\#include</b> <vigra/hdf5impex.hxx><br> 02870 Namespace: vigra 02871 */ 02872 template<size_t N, class C> 02873 void writeHDF5Attr(hid_t loc, 02874 const char* name, 02875 MultiArrayView<N, std::string, C> const & array) 02876 { 02877 if(H5Aexists(loc, name) > 0) 02878 H5Adelete(loc, name); 02879 02880 ArrayVector<hsize_t> shape(array.shape().begin(), 02881 array.shape().end()); 02882 HDF5Handle 02883 dataspace_handle(H5Screate_simple(N, shape.data(), NULL), 02884 &H5Sclose, 02885 "writeToHDF5File(): unable to create dataspace."); 02886 02887 HDF5Handle atype(H5Tcopy (H5T_C_S1), 02888 &H5Tclose, 02889 "writeToHDF5File(): unable to create type."); 02890 02891 detail::MaxSizeFnc max_size; 02892 max_size = std::for_each(array.data(),array.data()+ array.size(), max_size); 02893 H5Tset_size (atype, max_size.size); 02894 02895 HDF5Handle attr(H5Acreate(loc, 02896 name, 02897 atype, 02898 dataspace_handle, 02899 H5P_DEFAULT ,H5P_DEFAULT ), 02900 &H5Aclose, 02901 "writeHDF5Attr: unable to create Attribute"); 02902 02903 std::string buf =""; 02904 for(int ii = 0; ii < array.size(); ++ii) 02905 { 02906 buf = buf + array[ii] 02907 + std::string(max_size.size - array[ii].size(), ' '); 02908 } 02909 H5Awrite(attr, atype, buf.c_str()); 02910 } 02911 02912 /** Write a numeric ArrayVectorView as an attribute with name \a name 02913 of the dataset specified by the handle \a loc. 02914 02915 <b>\#include</b> <vigra/hdf5impex.hxx><br> 02916 Namespace: vigra 02917 */ 02918 template<class T> 02919 inline void writeHDF5Attr( hid_t loc, 02920 const char* name, 02921 ArrayVectorView<T> & array) 02922 { 02923 writeHDF5Attr(loc, name, 02924 MultiArrayView<1, T>(MultiArrayShape<1>::type(array.size()), 02925 array.data())); 02926 } 02927 02928 /** write an Attribute given a file and a path in the file. 02929 the path in the file should have the format 02930 [attribute] or /[subgroups/]dataset.attribute or 02931 /[subgroups/]group.attribute. 02932 The attribute is written to the root group, a dataset or a subgroup 02933 respectively 02934 */ 02935 template<class Arr> 02936 inline void writeHDF5Attr( std::string filePath, 02937 std::string pathInFile, 02938 Arr & ar) 02939 { 02940 std::string path_name(pathInFile), group_name, data_set_name, message, attr_name; 02941 std::string::size_type delimiter = path_name.rfind('/'); 02942 02943 //create or open file 02944 HDF5Handle file_id(createFile(filePath), &H5Fclose, 02945 "writeToHDF5File(): unable to open output file."); 02946 02947 // get the groupname and the filename 02948 if(delimiter == std::string::npos) 02949 { 02950 group_name = "/"; 02951 data_set_name = path_name; 02952 } 02953 02954 else 02955 { 02956 group_name = std::string(path_name.begin(), path_name.begin()+delimiter); 02957 data_set_name = std::string(path_name.begin()+delimiter+1, path_name.end()); 02958 } 02959 delimiter = data_set_name.rfind('.'); 02960 if(delimiter == std::string::npos) 02961 { 02962 attr_name = path_name; 02963 data_set_name = "/"; 02964 } 02965 else 02966 { 02967 attr_name = std::string(data_set_name.begin()+delimiter+1, data_set_name.end()); 02968 data_set_name = std::string(data_set_name.begin(), data_set_name.begin()+delimiter); 02969 } 02970 02971 HDF5Handle group(openGroup(file_id, group_name), &H5Gclose, 02972 "writeToHDF5File(): Unable to create and open group. attr ver"); 02973 02974 if(data_set_name != "/") 02975 { 02976 HDF5Handle dset(H5Dopen(group, data_set_name.c_str(), H5P_DEFAULT), &H5Dclose, 02977 "writeHDF5Attr():unable to open dataset"); 02978 writeHDF5Attr(hid_t(dset), attr_name.c_str(), ar); 02979 } 02980 else 02981 { 02982 writeHDF5Attr(hid_t(group), attr_name.c_str(), ar); 02983 } 02984 02985 } 02986 #endif 02987 02988 //@} 02989 02990 } // namespace vigra 02991 02992 #endif // VIGRA_HDF5IMPEX_HXX
© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de) |
html generated using doxygen and Python
|