[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]
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 # include <H5LT.h> 00066 #else 00067 # include <hdf5_hl.h> 00068 #endif 00069 00070 #include "impex.hxx" 00071 #include "multi_array.hxx" 00072 #include "multi_impex.hxx" 00073 00074 namespace vigra { 00075 00076 /** \addtogroup VigraHDF5Impex Import/Export of Images and Arrays in HDF5 Format 00077 00078 Supports arrays with arbitrary element types and arbitrary many dimensions. 00079 See the <a href="http://www.hdfgroup.org/HDF5/">HDF5 Website</a> for more 00080 information on the HDF5 file format. 00081 **/ 00082 //@{ 00083 00084 /** \brief Wrapper for hid_t objects. 00085 00086 Newly created or opened HDF5 handles are usually stored as objects of type 'hid_t'. When the handle 00087 is no longer needed, the appropriate close function must be called. However, if a function is 00088 aborted by an exception, this is difficult to ensure. Class HDF5Handle is a smart pointer that 00089 solves this problem by calling the close function in the destructor (This is analogous to how 00090 std::auto_ptr calls 'delete' on the contained pointer). A pointer to the close function must be 00091 passed to the constructor, along with an error message that is raised when creation/opening fails. 00092 00093 Since HDF5Handle objects are convertible to hid_t, they can be used in the code in place 00094 of the latter. 00095 00096 <b>Usage:</b> 00097 00098 \code 00099 HDF5Handle file_id(H5Fopen(filename, H5F_ACC_RDWR, H5P_DEFAULT), 00100 &H5Fclose, 00101 "Error message."); 00102 00103 ... // use file_id in the same way as a plain hid_t object 00104 \endcode 00105 00106 <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br> 00107 Namespace: vigra 00108 **/ 00109 class HDF5Handle 00110 { 00111 public: 00112 typedef herr_t (*Destructor)(hid_t); 00113 00114 private: 00115 hid_t handle_; 00116 Destructor destructor_; 00117 00118 public: 00119 00120 /** \brief Default constuctor. 00121 Creates a NULL handle. 00122 **/ 00123 HDF5Handle() 00124 : handle_( 0 ), 00125 destructor_(0) 00126 {} 00127 00128 /** \brief Create a wrapper for a hid_t object. 00129 00130 The hid_t object \a h is assumed to be the return value of an open or create function. 00131 It will be closed with the given close function \a destructor as soon as this 00132 HDF5Handle is destructed, except when \a destructor is a NULL pointer (in which 00133 case nothing happens at destruction time). If \a h has a value that indicates 00134 failed opening or creation (by HDF5 convention, this means if it is a negative number), 00135 an exception is raised by calling <tt>vigra_fail(error_message)</tt>. 00136 00137 <b>Usage:</b> 00138 00139 \code 00140 HDF5Handle file_id(H5Fopen(filename, H5F_ACC_RDWR, H5P_DEFAULT), 00141 &H5Fclose, 00142 "Error message."); 00143 00144 ... // use file_id in the same way 00145 \endcode 00146 **/ 00147 HDF5Handle(hid_t h, Destructor destructor, const char * error_message) 00148 : handle_( h ), 00149 destructor_(destructor) 00150 { 00151 if(handle_ < 0) 00152 vigra_fail(error_message); 00153 } 00154 00155 /** \brief Copy constructor. 00156 Hands over ownership of the RHS handle (analogous to std::auto_ptr). 00157 **/ 00158 HDF5Handle(HDF5Handle const & h) 00159 : handle_( h.handle_ ), 00160 destructor_(h.destructor_) 00161 { 00162 const_cast<HDF5Handle &>(h).handle_ = 0; 00163 } 00164 00165 /** \brief Assignment. 00166 Calls close() for the LHS handle and hands over ownership of the 00167 RHS handle (analogous to std::auto_ptr). 00168 **/ 00169 HDF5Handle & operator=(HDF5Handle const & h) 00170 { 00171 if(h.handle_ != handle_) 00172 { 00173 close(); 00174 handle_ = h.handle_; 00175 destructor_ = h.destructor_; 00176 const_cast<HDF5Handle &>(h).handle_ = 0; 00177 } 00178 return *this; 00179 } 00180 00181 /** \brief Destreuctor. 00182 Calls close() for the contained handle. 00183 **/ 00184 ~HDF5Handle() 00185 { 00186 close(); 00187 } 00188 00189 /** \brief Explicitly call the stored function (if one has been stored within 00190 this object) for the contained handle and set the handle to NULL. 00191 **/ 00192 herr_t close() 00193 { 00194 herr_t res = 1; 00195 if(handle_ && destructor_) 00196 res = (*destructor_)(handle_); 00197 handle_ = 0; 00198 return res; 00199 } 00200 00201 /** \brief Get a temporary hid_t object for the contained handle. 00202 Do not call a close function on the return value - a crash will be likely 00203 otherwise. 00204 **/ 00205 hid_t get() const 00206 { 00207 return handle_; 00208 } 00209 00210 /** \brief Convert to a plain hid_t object. 00211 00212 This function ensures that hid_t objects can be transparently replaced with 00213 HDF5Handle objects in user code. Do not call a close function on the return 00214 value - a crash will be likely otherwise. 00215 **/ 00216 operator hid_t() const 00217 { 00218 return handle_; 00219 } 00220 00221 /** \brief Equality comparison of the contained handle. 00222 **/ 00223 bool operator==(HDF5Handle const & h) const 00224 { 00225 return handle_ == h.handle_; 00226 } 00227 00228 /** \brief Equality comparison of the contained handle. 00229 **/ 00230 bool operator==(hid_t h) const 00231 { 00232 return handle_ == h; 00233 } 00234 00235 /** \brief Unequality comparison of the contained handle. 00236 **/ 00237 bool operator!=(HDF5Handle const & h) const 00238 { 00239 return handle_ != h.handle_; 00240 } 00241 00242 /** \brief Unequality comparison of the contained handle. 00243 **/ 00244 bool operator!=(hid_t h) const 00245 { 00246 return handle_ != h; 00247 } 00248 }; 00249 00250 00251 /********************************************************/ 00252 /* */ 00253 /* HDF5ImportInfo */ 00254 /* */ 00255 /********************************************************/ 00256 00257 /** \brief Argument object for the function readHDF5(). 00258 00259 See \ref readHDF5() for a usage example. This object must be 00260 used to read an image or array from a HDF5 file 00261 and enquire about its properties. 00262 00263 <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br> 00264 Namespace: vigra 00265 **/ 00266 class HDF5ImportInfo 00267 { 00268 public: 00269 enum PixelType { UINT8, UINT16, UINT32, UINT64, 00270 INT8, INT16, INT32, INT64, 00271 FLOAT, DOUBLE }; 00272 00273 /** Construct HDF5ImportInfo object. 00274 00275 The dataset \a pathInFile in the HDF5 file \a filename is accessed to 00276 read its properties. \a pathInFile may contain '/'-separated group 00277 names, but must end with the name of the desired dataset: 00278 00279 \code 00280 HDF5ImportInfo info(filename, "/group1/group2/my_dataset"); 00281 \endcode 00282 **/ 00283 VIGRA_EXPORT HDF5ImportInfo( const char* filePath, const char* pathInFile ); 00284 00285 VIGRA_EXPORT ~HDF5ImportInfo(); 00286 00287 /** Get the filename of this HDF5 object. 00288 **/ 00289 VIGRA_EXPORT const std::string& getFilePath() const; 00290 00291 /** Get the dataset's full name in the HDF5 file. 00292 **/ 00293 VIGRA_EXPORT const std::string& getPathInFile() const; 00294 00295 /** Get a handle to the file represented by this info object. 00296 **/ 00297 VIGRA_EXPORT hid_t getH5FileHandle() const; 00298 00299 /** Get a handle to the dataset represented by this info object. 00300 **/ 00301 VIGRA_EXPORT hid_t getDatasetHandle() const; 00302 00303 /** Get the number of dimensions of the dataset represented by this info object. 00304 **/ 00305 VIGRA_EXPORT MultiArrayIndex numDimensions() const; 00306 00307 /** Get the shape of the dataset represented by this info object. 00308 **/ 00309 VIGRA_EXPORT ArrayVector<hsize_t> const & shape() const 00310 { 00311 return m_dims; 00312 } 00313 00314 /** Get the shape (length) of the dataset along dimension \a dim. 00315 **/ 00316 VIGRA_EXPORT MultiArrayIndex shapeOfDimension(const int dim) const; 00317 00318 /** Query the pixel type of the dataset. 00319 00320 Possible values are: 00321 <DL> 00322 <DT>"UINT8"<DD> 8-bit unsigned integer (unsigned char) 00323 <DT>"INT16"<DD> 16-bit signed integer (short) 00324 <DT>"UINT16"<DD> 16-bit unsigned integer (unsigned short) 00325 <DT>"INT32"<DD> 32-bit signed integer (long) 00326 <DT>"UINT32"<DD> 32-bit unsigned integer (unsigned long) 00327 <DT>"FLOAT"<DD> 32-bit floating point (float) 00328 <DT>"DOUBLE"<DD> 64-bit floating point (double) 00329 </DL> 00330 **/ 00331 VIGRA_EXPORT const char * getPixelType() const; 00332 00333 /** Query the pixel type of the dataset. 00334 00335 Same as getPixelType(), but the result is returned as a 00336 ImageImportInfo::PixelType enum. This is useful to implement 00337 a switch() on the pixel type. 00338 00339 Possible values are: 00340 <DL> 00341 <DT>UINT8<DD> 8-bit unsigned integer (unsigned char) 00342 <DT>INT16<DD> 16-bit signed integer (short) 00343 <DT>UINT16<DD> 16-bit unsigned integer (unsigned short) 00344 <DT>INT32<DD> 32-bit signed integer (long) 00345 <DT>UINT32<DD> 32-bit unsigned integer (unsigned long) 00346 <DT>FLOAT<DD> 32-bit floating point (float) 00347 <DT>DOUBLE<DD> 64-bit floating point (double) 00348 </DL> 00349 **/ 00350 VIGRA_EXPORT PixelType pixelType() const; 00351 00352 private: 00353 HDF5Handle m_file_handle, m_dataset_handle; 00354 std::string m_filename, m_path, m_pixeltype; 00355 hssize_t m_dimensions; 00356 ArrayVector<hsize_t> m_dims; 00357 }; 00358 00359 00360 namespace detail { 00361 00362 template<class type> 00363 inline hid_t getH5DataType() 00364 { 00365 std::runtime_error("getH5DataType(): invalid type"); 00366 return 0; 00367 } 00368 00369 #define VIGRA_H5_DATATYPE(type, h5type) \ 00370 template<> \ 00371 inline hid_t getH5DataType<type>() \ 00372 { return h5type;} 00373 VIGRA_H5_DATATYPE(char, H5T_NATIVE_CHAR) 00374 VIGRA_H5_DATATYPE(Int8, H5T_NATIVE_INT8) 00375 VIGRA_H5_DATATYPE(Int16, H5T_NATIVE_INT16) 00376 VIGRA_H5_DATATYPE(Int32, H5T_NATIVE_INT32) 00377 VIGRA_H5_DATATYPE(Int64, H5T_NATIVE_INT64) 00378 VIGRA_H5_DATATYPE(UInt8, H5T_NATIVE_UINT8) 00379 VIGRA_H5_DATATYPE(UInt16, H5T_NATIVE_UINT16) 00380 VIGRA_H5_DATATYPE(UInt32, H5T_NATIVE_UINT32) 00381 VIGRA_H5_DATATYPE(UInt64, H5T_NATIVE_UINT64) 00382 VIGRA_H5_DATATYPE(float, H5T_NATIVE_FLOAT) 00383 VIGRA_H5_DATATYPE(double, H5T_NATIVE_DOUBLE) 00384 VIGRA_H5_DATATYPE(long double, H5T_NATIVE_LDOUBLE) 00385 00386 #undef VIGRA_H5_DATATYPE 00387 00388 } // namespace detail 00389 00390 00391 00392 00393 /********************************************************/ 00394 /* */ 00395 /* HDF5File */ 00396 /* */ 00397 /********************************************************/ 00398 00399 00400 /** \brief Access to HDF5 files 00401 00402 HDF5File proviedes a convenient way of accessing data in HDF5 files. vigra::MultiArray 00403 structures of any dimension can be stored to / loaded from HDF5 files. Typical 00404 HDF5 features like subvolume access, chunks and data compression are available, 00405 string attributes can be attached to any dataset or group. Group- or dataset-handles 00406 are encapsulated in the class and managed automatically. The internal file-system like 00407 structure can be accessed by functions like "cd()" or "mkdir()". 00408 00409 00410 <b>Example:</b> 00411 Write the MultiArray out_multi_array to file. Change the current directory to 00412 "/group" and read in the same MultiArray as in_multi_array. 00413 \code 00414 HDF5File file("/path/to/file",HDF5File::New); 00415 file.mkdir("group"); 00416 file.write("/group/dataset", out_multi_array); 00417 00418 file.cd("/group"); 00419 file.read("dataset", in_multi_array); 00420 00421 \endcode 00422 00423 <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br> 00424 Namespace: vigra 00425 **/ 00426 class HDF5File 00427 { 00428 private: 00429 HDF5Handle fileHandle_; 00430 00431 // current group handle 00432 HDF5Handle cGroupHandle_; 00433 00434 00435 // datastructure to hold a list of dataset and group names 00436 struct lsOpData 00437 { 00438 std::vector<std::string> objects; 00439 }; 00440 00441 // operator function, used by H5Literate 00442 static herr_t opFunc (hid_t loc_id, const char *name, const H5L_info_t *info, void *operator_data) 00443 { 00444 // get information about object 00445 H5O_info_t infobuf; 00446 H5Oget_info_by_name (loc_id, name, &infobuf, H5P_DEFAULT); 00447 00448 // add name to list, if object is a dataset or a group 00449 if(infobuf.type == H5O_TYPE_GROUP) 00450 { 00451 (*(struct lsOpData *) operator_data).objects.push_back(std::string(name)+"/"); 00452 } 00453 if(infobuf.type == H5O_TYPE_DATASET) 00454 { 00455 (*(struct lsOpData *) operator_data).objects.push_back(std::string(name)); 00456 } 00457 00458 return 0; 00459 } 00460 00461 public: 00462 /** \brief Set how a file is opened. 00463 OpenMode::New creates a new file. If the file already exists, overwrite it. 00464 00465 OpenMode::Open opens a file for reading/writing. The file will be created, 00466 if nescessary. 00467 **/ 00468 enum OpenMode { 00469 New, // Create new empty file (existing file will be deleted). 00470 Open // Open file. Create if not existing. 00471 }; 00472 00473 00474 /** \brief Create a HDF5File object. 00475 00476 Creates or opens HDF5 file at position filename. The current group is set 00477 to "/". 00478 **/ 00479 HDF5File(std::string filename, OpenMode mode) 00480 { 00481 std::string errorMessage = "HDF5File: Could not create file '" + filename + "'."; 00482 fileHandle_ = HDF5Handle(createFile_(filename, mode), &H5Fclose, errorMessage.c_str()); 00483 cGroupHandle_ = HDF5Handle(openCreateGroup_("/"), &H5Gclose, "HDF5File(): Failed to open root group."); 00484 } 00485 00486 00487 /** \brief Destructor to make sure that all data is flushed before closing the file. 00488 */ 00489 ~HDF5File() 00490 { 00491 //Write everything to disk before closing 00492 H5Fflush(fileHandle_, H5F_SCOPE_GLOBAL); 00493 } 00494 00495 00496 /** \brief Change current group to "/". 00497 */ 00498 void root() 00499 { 00500 std::string message = "HDF5File::root(): Could not open group '/'."; 00501 cGroupHandle_ = HDF5Handle(H5Gopen(fileHandle_, "/", H5P_DEFAULT),&H5Gclose,message.c_str()); 00502 } 00503 00504 00505 /** \brief Change the current group. 00506 If the first character is a "/", the path will be interpreted as absolute path, 00507 otherwise it will be interpreted as path relative to the current group. 00508 */ 00509 void cd(std::string groupName) 00510 { 00511 std::string message = "HDF5File::cd(): Could not open group '" + groupName + "'.\n"; 00512 if(groupName == "/") 00513 { 00514 cGroupHandle_ = HDF5Handle(openCreateGroup_("/"),&H5Gclose,message.c_str()); 00515 return; 00516 } 00517 else if(groupName =="..") 00518 { 00519 cd_up(); 00520 } 00521 else{ 00522 if(relativePath_(groupName)) 00523 { 00524 if (H5Lexists(cGroupHandle_, groupName.c_str(), H5P_DEFAULT) == 0) 00525 { 00526 std::cerr << message; 00527 return; 00528 } 00529 cGroupHandle_ = HDF5Handle(openCreateGroup_(groupName),&H5Gclose,message.c_str()); 00530 } 00531 else 00532 { 00533 if (H5Lexists(fileHandle_, groupName.c_str(), H5P_DEFAULT) == 0) 00534 { 00535 std::cerr << message; 00536 return; 00537 } 00538 cGroupHandle_ = HDF5Handle(openCreateGroup_(groupName),&H5Gclose,message.c_str()); 00539 } 00540 } 00541 } 00542 00543 /** \brief Change the current group to its parent group. 00544 returns true if successful, false otherwise. 00545 */ 00546 bool cd_up() 00547 { 00548 std::string groupName = currentGroupName_(); 00549 00550 //do not try to move up if we already in "/" 00551 if(groupName == "/"){ 00552 return false; 00553 } 00554 00555 size_t lastSlash = groupName.find_last_of('/'); 00556 00557 std::string parentGroup (groupName.begin(), groupName.begin()+lastSlash+1); 00558 00559 cd(parentGroup); 00560 00561 return true; 00562 } 00563 void cd_up(int levels) 00564 { 00565 for(int i = 0; i<levels; i++) 00566 cd_up(); 00567 } 00568 00569 00570 /** \brief Create a new group. 00571 If the first character is a "/", the path will be interpreted as absolute path, 00572 otherwise it will be interpreted as path relative to the current group. 00573 */ 00574 void mkdir(std::string groupName) 00575 { 00576 hid_t handle = openCreateGroup_(groupName.c_str()); 00577 if (handle != cGroupHandle_){ 00578 H5Gclose(handle); 00579 } 00580 } 00581 00582 /** \brief Change the current group; create it if nescessary. 00583 If the first character is a "/", the path will be interpreted as absolute path, 00584 otherwise it will be interpreted as path relative to the current group. 00585 */ 00586 void cd_mk(std::string groupName) 00587 { 00588 std::string message = "HDF5File::cd_mk(): Could not create group '" + groupName + "'."; 00589 cGroupHandle_ = HDF5Handle(openCreateGroup_(groupName.c_str()),&H5Gclose,message.c_str()); 00590 } 00591 00592 00593 00594 /** \brief List the content of the current group. 00595 The function returns a vector of strings holding the entries of the current 00596 group. Only datasets and groups are listed, other objects (e.g. datatypes) 00597 are ignored. Group names always have an ending "/". 00598 */ 00599 std::vector<std::string> ls() 00600 { 00601 lsOpData data; 00602 H5Literate(cGroupHandle_,H5_INDEX_NAME,H5_ITER_NATIVE,NULL, &opFunc, (void *) &data); 00603 00604 return data.objects; 00605 } 00606 00607 00608 /** \brief Get the path of the current group. 00609 */ 00610 std::string pwd() 00611 { 00612 return currentGroupName_(); 00613 } 00614 00615 /** \brief Get the name of the associated file. 00616 */ 00617 std::string filename() 00618 { 00619 return fileName_(); 00620 } 00621 00622 /** \brief Get the number of dimensions of a certain dataset 00623 If the first character is a "/", the path will be interpreted as absolute path, 00624 otherwise it will be interpreted as path relative to the current group. 00625 */ 00626 hssize_t getDatasetDimensions(std::string datasetName) 00627 { 00628 //Open dataset and dataspace 00629 std::string errorMessage = "HDF5File::getDatasetDimensions(): Unable to open dataset '" + datasetName + "'."; 00630 HDF5Handle datasetHandle = HDF5Handle(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str()); 00631 00632 errorMessage = "HDF5File::getDatasetDimensions(): Unable to access dataspace."; 00633 HDF5Handle dataspaceHandle(H5Dget_space(datasetHandle), &H5Sclose, errorMessage.c_str()); 00634 00635 //return dimension information 00636 return H5Sget_simple_extent_ndims(dataspaceHandle); 00637 } 00638 00639 /** \brief Get the shape of each dimension of a certain dataset. 00640 Normally, this function is called after determining the dimension of the 00641 dataset using \ref getDatasetDimensions(). 00642 If the first character is a "/", the path will be interpreted as absolute path, 00643 otherwise it will be interpreted as path relative to the current group. 00644 */ 00645 ArrayVector<hsize_t> getDatasetShape(std::string datasetName) 00646 { 00647 //Open dataset and dataspace 00648 std::string errorMessage = "HDF5File::getDatasetShape(): Unable to open dataset '" + datasetName + "'."; 00649 HDF5Handle datasetHandle = HDF5Handle(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str()); 00650 00651 errorMessage = "HDF5File::getDatasetShape(): Unable to access dataspace."; 00652 HDF5Handle dataspaceHandle(H5Dget_space(datasetHandle), &H5Sclose, errorMessage.c_str()); 00653 00654 //get dimension information 00655 ArrayVector<hsize_t>::size_type dimensions = H5Sget_simple_extent_ndims(dataspaceHandle); 00656 00657 ArrayVector<hsize_t> shape(dimensions); 00658 ArrayVector<hsize_t> maxdims(dimensions); 00659 H5Sget_simple_extent_dims(dataspaceHandle, shape.data(), maxdims.data()); 00660 00661 00662 // invert the dimensions to guarantee c-order 00663 ArrayVector<hsize_t> shape_inv(dimensions); 00664 for(ArrayVector<hsize_t>::size_type i=0; i<shape.size(); i++) { 00665 shape_inv[i] = shape[dimensions-1-i]; 00666 } 00667 00668 return shape_inv; 00669 } 00670 00671 /** \brief Attach a string attribute to an existing object. 00672 The attribute can be attached to datasets and groups. The string may have arbitrary length. 00673 */ 00674 void setAttribute(std::string datasetName, std::string attributeName, std::string text) 00675 { 00676 std::string groupname = SplitString(datasetName).first(); 00677 std::string setname = SplitString(datasetName).last(); 00678 00679 std::string errorMessage ("HDF5File::setAttribute(): Unable to open group '" + groupname + "'."); 00680 HDF5Handle groupHandle (openCreateGroup_(groupname), &H5Gclose, errorMessage.c_str()); 00681 00682 H5LTset_attribute_string(groupHandle,setname.c_str(), attributeName.c_str(),text.c_str()); 00683 } 00684 00685 00686 /** \brief Get an attribute string of an object. 00687 */ 00688 std::string getAttribute(std::string datasetName, std::string attributeName) 00689 { 00690 typedef ArrayVector<char>::size_type SizeType; 00691 if (H5Lexists(fileHandle_, datasetName.c_str(), H5P_DEFAULT) == 0) 00692 { 00693 std::cerr << "HDF5File::getAttribute(): Dataset '" << datasetName << "' does not exist.\n"; 00694 return "error - dataset not found"; 00695 } 00696 00697 std::string groupname = SplitString(datasetName).first(); 00698 std::string setname = SplitString(datasetName).last(); 00699 00700 std::string errorMessage ("HDF5File::setAttribute(): Unable to open group '" + groupname + "'."); 00701 HDF5Handle groupHandle (openCreateGroup_(groupname), &H5Gclose, errorMessage.c_str()); 00702 00703 // get the size of the attribute 00704 HDF5Handle AttrHandle (H5Aopen_by_name(groupHandle,setname.c_str(),attributeName.c_str(),H5P_DEFAULT, H5P_DEFAULT),&H5Aclose, "HDF5File::getAttribute(): Unable to open attribute."); 00705 SizeType len = (SizeType)H5Aget_storage_size(AttrHandle); 00706 00707 //read the attribute 00708 ArrayVector<char> text (len+1,0); 00709 H5LTget_attribute_string(groupHandle, setname.c_str(), attributeName.c_str(), text.begin()); 00710 00711 return std::string(text.begin()); 00712 } 00713 00714 00715 // Writing data 00716 00717 /** \brief Write multi arrays. 00718 Chunks can be activated by setting \code iChunkSize = size; //size > 0 \endcode . 00719 The chunks will be hypercubes with edge length size. 00720 00721 Compression can be activated by setting \code compression = parameter; // 0 < parameter <= 9 \endcode 00722 where 0 stands for no compression and 9 for maximum compression. 00723 00724 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 00725 otherwise it will be interpreted as path relative to the current group. 00726 */ 00727 template<unsigned int N, class T> 00728 inline void write(std::string datasetName, const MultiArrayView<N, T, UnstridedArrayTag> & array, int iChunkSize = 0, int compression = 0) 00729 { 00730 typename MultiArrayShape<N>::type chunkSize; 00731 for(int i = 0; i < N; i++){ 00732 chunkSize[i] = iChunkSize; 00733 } 00734 write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize, compression); 00735 } 00736 00737 /** \brief Write multi arrays. 00738 Chunks can be activated by providing a MultiArrayShape as chunkSize. 00739 chunkSize must have equal dimension as array. 00740 00741 Compression can be activated by setting \code compression = parameter; // 0 < parameter <= 9 \endcode 00742 where 0 stands for no compression and 9 for maximum compression. 00743 00744 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 00745 otherwise it will be interpreted as path relative to the current group. 00746 */ 00747 template<unsigned int N, class T> 00748 inline void write(std::string datasetName, const MultiArrayView<N, T, UnstridedArrayTag> & array, typename MultiArrayShape<N>::type chunkSize, int compression = 0) 00749 { 00750 write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize, compression); 00751 } 00752 00753 /** \brief Write a multi array into a larger volume. 00754 blockOffset determines the position, where array is written. 00755 00756 Chunks can be activated by providing a MultiArrayShape as chunkSize. 00757 chunkSize must have equal dimension as array. 00758 00759 Compression can be activated by setting \code compression = parameter; // 0 < parameter <= 9 \endcode 00760 where 0 stands for no compression and 9 for maximum compression. 00761 00762 If the first character of datasetName 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 template<unsigned int N, class T> 00766 inline void writeBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, const MultiArrayView<N, T, UnstridedArrayTag> & array) 00767 { 00768 writeBlock_(datasetName, blockOffset, array, detail::getH5DataType<T>(), 1); 00769 } 00770 00771 // non-scalar (TinyVector) and unstrided multi arrays 00772 template<unsigned int N, class T, int SIZE> 00773 inline void write(std::string datasetName, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array, int iChunkSize = 0, int compression = 0) 00774 { 00775 typename MultiArrayShape<N>::type chunkSize; 00776 for(int i = 0; i < N; i++){ 00777 chunkSize[i] = iChunkSize; 00778 } 00779 write_(datasetName, array, detail::getH5DataType<T>(), SIZE, chunkSize, compression); 00780 } 00781 00782 template<unsigned int N, class T, int SIZE> 00783 inline void write(std::string datasetName, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array, typename MultiArrayShape<N>::type chunkSize, int compression = 0) 00784 { 00785 write_(datasetName, array, detail::getH5DataType<T>(), SIZE, chunkSize, compression); 00786 } 00787 00788 template<unsigned int N, class T, int SIZE> 00789 inline void writeBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array) 00790 { 00791 writeBlock_(datasetName, blockOffset, array, detail::getH5DataType<T>(), SIZE); 00792 } 00793 00794 // non-scalar (RGBValue) and unstrided multi arrays 00795 template<unsigned int N, class T> 00796 inline void write(std::string datasetName, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array, int iChunkSize = 0, int compression = 0) 00797 { 00798 typename MultiArrayShape<N>::type chunkSize; 00799 for(int i = 0; i < N; i++){ 00800 chunkSize[i] = iChunkSize; 00801 } 00802 write_(datasetName, array, detail::getH5DataType<T>(), 3, chunkSize, compression); 00803 } 00804 00805 template<unsigned int N, class T> 00806 inline void write(std::string datasetName, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array, typename MultiArrayShape<N>::type chunkSize, int compression = 0) 00807 { 00808 write_(datasetName, array, detail::getH5DataType<T>(), 3, chunkSize, compression); 00809 } 00810 00811 template<unsigned int N, class T> 00812 inline void writeBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array) 00813 { 00814 writeBlock_(datasetName, blockOffset, array, detail::getH5DataType<T>(), 3); 00815 } 00816 00817 /** \brief Write single value as dataset. 00818 This functions allows to write data of atomic datatypes (int, long, double) 00819 as a dataset in the HDF5 file. So it is not nescessary to create a MultiArray 00820 of size 1 to write a single number. 00821 00822 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 00823 otherwise it will be interpreted as path relative to the current group. 00824 */ 00825 template<class T> 00826 inline void writeAtomic(std::string datasetName, const T data) 00827 { 00828 typename MultiArrayShape<1>::type chunkSize; 00829 chunkSize[0] = 0; 00830 MultiArray<1,T> array(MultiArrayShape<1>::type(1)); 00831 array[0] = data; 00832 write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize,0); 00833 } 00834 00835 00836 /* Reading data. */ 00837 00838 /** \brief Read data into a multi array. 00839 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 00840 otherwise it will be interpreted as path relative to the current group. 00841 */ 00842 template<unsigned int N, class T> 00843 inline void read(std::string datasetName, MultiArrayView<N, T, UnstridedArrayTag> array) 00844 { 00845 read_(datasetName, array, detail::getH5DataType<T>(), 1); 00846 } 00847 00848 /** \brief Read a block of data into s multi array. 00849 This function allows to read a small block out of a larger volume stored 00850 in a HDF5 dataset. 00851 00852 blockOffset determines the position of the block. 00853 blockSize determines the size in each dimension of the block. 00854 00855 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 00856 otherwise it will be interpreted as path relative to the current group. 00857 */ 00858 template<unsigned int N, class T> 00859 inline void readBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, typename MultiArrayShape<N>::type blockShape, MultiArrayView<N, T, UnstridedArrayTag> array) 00860 { 00861 readBlock_(datasetName, blockOffset, blockShape, array, detail::getH5DataType<T>(), 1); 00862 } 00863 00864 // non-scalar (TinyVector) and unstrided target multi array 00865 template<unsigned int N, class T, int SIZE> 00866 inline void read(std::string datasetName, MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> array) 00867 { 00868 read_(datasetName, array, detail::getH5DataType<T>(), SIZE); 00869 } 00870 00871 template<unsigned int N, class T, int SIZE> 00872 inline void readBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, typename MultiArrayShape<N>::type blockShape, MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> array) 00873 { 00874 readBlock_(datasetName, blockOffset, blockShape, array, detail::getH5DataType<T>(), SIZE); 00875 } 00876 00877 // non-scalar (RGBValue) and unstrided target multi array 00878 template<unsigned int N, class T> 00879 inline void read(std::string datasetName, MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> array) 00880 { 00881 read_(datasetName, array, detail::getH5DataType<T>(), 3); 00882 } 00883 00884 template<unsigned int N, class T> 00885 inline void readBlock(std::string datasetName, typename MultiArrayShape<N>::type blockOffset, typename MultiArrayShape<N>::type blockShape, MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> array) 00886 { 00887 readBlock_(datasetName, blockOffset, blockShape, array, detail::getH5DataType<T>(), 3); 00888 } 00889 00890 /** \brief Read a single value. 00891 This functions allows to read a single datum of atomic datatype (int, long, double) 00892 from the HDF5 file. So it is not nescessary to create a MultiArray 00893 of size 1 to read a single number. 00894 00895 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 00896 otherwise it will be interpreted as path relative to the current group. 00897 */ 00898 template<class T> 00899 inline void readAtomic(std::string datasetName, T & data) 00900 { 00901 MultiArray<1,T> array(MultiArrayShape<1>::type(1)); 00902 read_(datasetName, array, detail::getH5DataType<T>(), 1); 00903 data = array[0]; 00904 } 00905 00906 00907 /** \brief Create a new dataset. 00908 This function can be used to create a dataset filled with a default value, 00909 for example before writing data into it using \ref writeBlock(). 00910 Attention: only atomic datatypes are provided. For spectral data, add an 00911 dimension (case RGB: add one dimension of size 3). 00912 00913 shape determines the dimension and the size of the dataset. 00914 00915 Chunks can be activated by providing a MultiArrayShape as chunkSize. 00916 chunkSize must have equal dimension as array. 00917 00918 Compression can be activated by setting \code compression = parameter; // 0 < parameter <= 9 \endcode 00919 where 0 stands for no compression and 9 for maximum compression. 00920 00921 If the first character of datasetName is a "/", the path will be interpreted as absolute path, 00922 otherwise it will be interpreted as path relative to the current group. 00923 */ 00924 template<unsigned int N, class T> 00925 inline void createDataset(std::string datasetName, typename MultiArrayShape<N>::type shape, T init, int iChunkSize = 0, int compressionParameter = 0) 00926 { 00927 typename MultiArrayShape<N>::type chunkSize; 00928 for(int i = 0; i < N; i++){ 00929 chunkSize[i] = iChunkSize; 00930 } 00931 createDataset<N,T>(datasetName, shape, init, chunkSize, compressionParameter); 00932 } 00933 00934 template<unsigned int N, class T> 00935 inline void createDataset(std::string datasetName, typename MultiArrayShape<N>::type shape, T init, typename MultiArrayShape<N>::type chunkSize, int compressionParameter = 0) 00936 { 00937 std::string groupname = SplitString(datasetName).first(); 00938 std::string setname = SplitString(datasetName).last(); 00939 00940 hid_t parent = openCreateGroup_(groupname); 00941 00942 // delete the dataset if it already exists 00943 deleteDataset_(parent, setname); 00944 00945 // create dataspace 00946 // add an extra dimension in case that the data is non-scalar 00947 HDF5Handle dataspaceHandle; 00948 00949 // invert dimensions to guarantee c-order 00950 hsize_t shape_inv[N]; 00951 for(unsigned int k=0; k<N; ++k) 00952 shape_inv[N-1-k] = shape[k]; 00953 00954 // create dataspace 00955 dataspaceHandle = HDF5Handle(H5Screate_simple(N, shape_inv, NULL), 00956 &H5Sclose, "HDF5File::createDataset(): unable to create dataspace for scalar data."); 00957 00958 // set fill value 00959 HDF5Handle plist ( H5Pcreate(H5P_DATASET_CREATE), &H5Pclose, "HDF5File::createDataset(): unable to create property list." ); 00960 H5Pset_fill_value(plist,detail::getH5DataType<T>(), &init); 00961 00962 // enable chunks 00963 if(chunkSize[0] > 0) 00964 { 00965 hsize_t cSize [N]; 00966 for(int i = 0; i<N; i++) 00967 { 00968 cSize[i] = chunkSize[N-1-i]; 00969 } 00970 H5Pset_chunk (plist, N, cSize); 00971 } 00972 00973 // enable compression 00974 if(compressionParameter > 0) 00975 { 00976 H5Pset_deflate(plist, compressionParameter); 00977 } 00978 00979 //create the dataset. 00980 HDF5Handle datasetHandle ( H5Dcreate(parent, setname.c_str(), detail::getH5DataType<T>(), dataspaceHandle, H5P_DEFAULT, plist, H5P_DEFAULT), 00981 &H5Dclose, "HDF5File::createDataset(): unable to create dataset."); 00982 if(parent != cGroupHandle_) 00983 H5Gclose(parent); 00984 } 00985 00986 /** \brief Immediately write all data to disk 00987 */ 00988 void flushToDisk() 00989 { 00990 H5Fflush(fileHandle_, H5F_SCOPE_GLOBAL); 00991 } 00992 00993 00994 private: 00995 00996 /* Simple extension of std::string for splitting into two parts 00997 * 00998 * Strings (in particular: file/dataset paths) will be split into two 00999 * parts. The split is made at the last occurance of the delimiter. 01000 * 01001 * For example, "/path/to/some/file" will be split (delimiter = "/") into 01002 * first() = "/path/to/some" and last() = "file". 01003 */ 01004 class SplitString: public std::string { 01005 public: 01006 SplitString(std::string &sstring): std::string(sstring) {}; 01007 01008 // return the part of the string before the delimiter 01009 std::string first(char delimiter = '/') 01010 { 01011 size_t last = find_last_of(delimiter); 01012 if(last == std::string::npos) // delimiter not found --> no first 01013 return ""; 01014 01015 return std::string(begin(), begin()+last+1); 01016 } 01017 01018 // return the part of the string after the delimiter 01019 std::string last(char delimiter = '/') 01020 { 01021 size_t last = find_last_of(delimiter); 01022 if(last == std::string::npos) // delimiter not found --> only last 01023 return std::string(*this); 01024 return std::string(begin()+last+1, end()); 01025 } 01026 }; 01027 01028 inline bool relativePath_(std::string & path) 01029 { 01030 std::string::size_type pos = path.find('/') ; 01031 if(pos == 0) 01032 return false; 01033 01034 return true; 01035 } 01036 01037 01038 std::string currentGroupName_() 01039 { 01040 int len = H5Iget_name(cGroupHandle_,NULL,1000); 01041 ArrayVector<char> name (len+1,0); 01042 H5Iget_name(cGroupHandle_,name.begin(),len+1); 01043 01044 return std::string(name.begin()); 01045 } 01046 01047 std::string fileName_() 01048 { 01049 int len = H5Fget_name(fileHandle_,NULL,1000); 01050 ArrayVector<char> name (len+1,0); 01051 H5Fget_name(fileHandle_,name.begin(),len+1); 01052 01053 return std::string(name.begin()); 01054 } 01055 01056 01057 inline hid_t createFile_(std::string filePath, OpenMode mode = Open) 01058 { 01059 // try to open file 01060 FILE * pFile; 01061 pFile = fopen ( filePath.c_str(), "r" ); 01062 hid_t fileId; 01063 01064 // check if opening was successful (= file exists) 01065 if ( pFile == NULL ) 01066 { 01067 fileId = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT); 01068 } 01069 else if(mode == Open) 01070 { 01071 fclose( pFile ); 01072 fileId = H5Fopen(filePath.c_str(), H5F_ACC_RDWR, H5P_DEFAULT); 01073 } 01074 else 01075 { 01076 fclose(pFile); 01077 std::remove(filePath.c_str()); 01078 fileId = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT); 01079 } 01080 return fileId; 01081 } 01082 01083 01084 // open a group and subgroups. Create if nescessary. 01085 inline hid_t openCreateGroup_(std::string groupName) 01086 { 01087 // check if groupName is absolute or relative 01088 hid_t parent; 01089 if(relativePath_(groupName)) 01090 { 01091 parent = cGroupHandle_; 01092 } 01093 else 01094 { 01095 // open root group 01096 parent = H5Gopen(fileHandle_, "/", H5P_DEFAULT); 01097 if(groupName == "/") 01098 { 01099 return parent; 01100 } 01101 01102 // remove leading / 01103 groupName = std::string(groupName.begin()+1, groupName.end()); 01104 } 01105 01106 01107 // check if the groupName has the right format 01108 if( groupName.size() != 0 && *groupName.rbegin() != '/') 01109 { 01110 groupName = groupName + '/'; 01111 } 01112 01113 01114 //create subgroups one by one 01115 std::string::size_type begin = 0, end = groupName.find('/'); 01116 int ii = 0; 01117 while (end != std::string::npos) 01118 { 01119 std::string group(groupName.begin()+begin, groupName.begin()+end); 01120 hid_t prevParent = parent; 01121 01122 if(H5LTfind_dataset(parent, group.c_str()) == 0) 01123 { 01124 parent = H5Gcreate(prevParent, group.c_str(), H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); 01125 } else { 01126 parent = H5Gopen(prevParent, group.c_str(), H5P_DEFAULT); 01127 } 01128 01129 if(ii != 0) 01130 { 01131 H5Gclose(prevParent); 01132 } 01133 if(parent < 0) 01134 { 01135 return parent; 01136 } 01137 ++ii; 01138 begin = end + 1; 01139 end = groupName.find('/', begin); 01140 } 01141 01142 01143 return parent; 01144 } 01145 01146 01147 inline void deleteDataset_(hid_t parent, std::string datasetName) 01148 { 01149 // delete existing data and create new dataset 01150 if(H5LTfind_dataset(parent, datasetName.c_str())) 01151 { 01152 //std::cout << "dataset already exists" << std::endl; 01153 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6) 01154 if(H5Gunlink(parent, datasetName.c_str()) < 0) 01155 { 01156 vigra_postcondition(false, "HDF5File::deleteDataset_(): Unable to delete existing data."); 01157 } 01158 #else 01159 if(H5Ldelete(parent, datasetName.c_str(), H5P_DEFAULT ) < 0) 01160 { 01161 vigra_postcondition(false, "HDF5File::deleteDataset_(): Unable to delete existing data."); 01162 } 01163 #endif 01164 } 01165 } 01166 01167 hid_t getDatasetHandle_(std::string datasetName) 01168 { 01169 std::string groupname = SplitString(datasetName).first(); 01170 std::string setname = SplitString(datasetName).last(); 01171 01172 if(relativePath_(datasetName)) 01173 { 01174 if (H5Lexists(cGroupHandle_, datasetName.c_str(), H5P_DEFAULT) != 1) 01175 { 01176 std::cerr << "HDF5File::getDatasetHandle_(): Dataset '" << datasetName << "' does not exist.\n"; 01177 return -1; 01178 } 01179 //Open parent group 01180 hid_t groupHandle = openCreateGroup_(groupname); 01181 01182 hid_t datasetHandle = H5Dopen(groupHandle, setname.c_str(), H5P_DEFAULT); 01183 01184 if(groupHandle != cGroupHandle_){ 01185 H5Gclose(groupHandle); 01186 } 01187 01188 //return dataset handle 01189 return datasetHandle; 01190 } 01191 else 01192 { 01193 if (H5Lexists(fileHandle_, datasetName.c_str(), H5P_DEFAULT) != 1) 01194 { 01195 std::cerr << "HDF5File::getDatasetHandle_(): Dataset '" << datasetName << "' does not exist.\n"; 01196 return -1; 01197 } 01198 01199 //Open parent group 01200 hid_t groupHandle = openCreateGroup_(groupname); 01201 01202 hid_t datasetHandle = H5Dopen(groupHandle, setname.c_str(), H5P_DEFAULT); 01203 01204 if(groupHandle != cGroupHandle_){ 01205 H5Gclose(groupHandle); 01206 } 01207 01208 //return dataset handle 01209 return datasetHandle; 01210 } 01211 } 01212 01213 01214 // unstrided multi arrays 01215 template<unsigned int N, class T> 01216 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) 01217 { 01218 std::string groupname = SplitString(datasetName).first(); 01219 std::string setname = SplitString(datasetName).last(); 01220 01221 // shape of the array. Add one dimension, if array contains non-scalars. 01222 ArrayVector<hsize_t> shape(N + (numBandsOfType > 1),0); 01223 for(int i = 0; i < N; i++){ 01224 shape[N-1-i] = array.shape(i); // reverse order 01225 } 01226 01227 if(numBandsOfType > 1) 01228 shape[N] = numBandsOfType; 01229 01230 HDF5Handle dataspace ( H5Screate_simple(N + (numBandsOfType > 1), shape.begin(), NULL), &H5Sclose, "HDF5File::write(): Can not create dataspace."); 01231 01232 // create and open group: 01233 std::string errorMessage ("HDF5File::write(): can not create group '" + groupname + "'."); 01234 hid_t groupHandle = openCreateGroup_(groupname); 01235 if(groupHandle <= 0) 01236 { 01237 std::cerr << errorMessage << "\n"; 01238 } 01239 01240 // delete dataset, if it already exists 01241 deleteDataset_(groupHandle, setname.c_str()); 01242 01243 // set up properties list 01244 HDF5Handle plist ( H5Pcreate(H5P_DATASET_CREATE), &H5Pclose, "HDF5File::write(): unable to create property list." ); 01245 01246 // enable chunks 01247 if(chunkSize[0] > 0) 01248 { 01249 ArrayVector<hsize_t> cSize(N + (numBandsOfType > 1),0); 01250 for(int i = 0; i<N; i++) 01251 { 01252 cSize[i] = chunkSize[N-1-i]; 01253 } 01254 if(numBandsOfType > 1) 01255 cSize[N] = numBandsOfType; 01256 01257 H5Pset_chunk (plist, N + (numBandsOfType > 1), cSize.begin()); 01258 } 01259 01260 // enable compression 01261 if(compressionParameter > 0) 01262 { 01263 H5Pset_deflate(plist, compressionParameter); 01264 } 01265 01266 // create dataset 01267 HDF5Handle datasetHandle (H5Dcreate(groupHandle, setname.c_str(), datatype, dataspace,H5P_DEFAULT, plist, H5P_DEFAULT), &H5Dclose, "HDF5File::write(): Can not create dataset."); 01268 01269 // Write the data to the HDF5 dataset as is 01270 H5Dwrite( datasetHandle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data()); 01271 01272 if(groupHandle != cGroupHandle_) 01273 { 01274 H5Gclose(groupHandle); 01275 } 01276 } 01277 01278 // unstrided target multi array 01279 template<unsigned int N, class T> 01280 void read_(std::string datasetName, MultiArrayView<N, T, UnstridedArrayTag> array, const hid_t datatype, const int numBandsOfType) 01281 { 01282 //Prepare to read without using HDF5ImportInfo 01283 ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName); 01284 hssize_t dimensions = getDatasetDimensions(datasetName); 01285 01286 std::string errorMessage ("HDF5File::read(): Unable to open dataset '" + datasetName + "'."); 01287 HDF5Handle datasetHandle (getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str()); 01288 01289 int offset = (numBandsOfType > 1); 01290 01291 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 01292 "HDF5File::read(): Array dimension disagrees with dataset dimension."); 01293 01294 typename MultiArrayShape<N>::type shape; 01295 for(int k=offset; k< MultiArrayIndex(dimensions); ++k) { 01296 shape[k-offset] = MultiArrayIndex(dimshape[k]); 01297 } 01298 01299 vigra_precondition(shape == array.shape(), 01300 "HDF5File::read(): Array shape disagrees with dataset shape."); 01301 01302 // simply read in the data as is 01303 H5Dread( datasetHandle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data() ); // .data() possible since void pointer! 01304 } 01305 01306 template<unsigned int N, class T> 01307 void writeBlock_(std::string datasetName, typename MultiArrayShape<N>::type &blockOffset, const MultiArrayView<N, T, UnstridedArrayTag> & array, const hid_t datatype, const int numBandsOfType) 01308 { 01309 // open dataset if it exists 01310 std::string errorMessage = "HDF5File::writeBlock(): Error opening dataset '" + datasetName + "'."; 01311 HDF5Handle datasetHandle (getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str()); 01312 01313 01314 // hyperslab parameters for position, size, ... 01315 hsize_t boffset [N]; 01316 hsize_t bshape [N]; 01317 hsize_t bones [N]; 01318 01319 for(int i = 0; i < N; i++){ 01320 boffset[i] = blockOffset[N-1-i]; 01321 bshape[i] = array.size(N-1-i); 01322 bones[i] = 1; 01323 } 01324 01325 // create a target dataspace in memory with the shape of the desired block 01326 HDF5Handle memspace_handle (H5Screate_simple(N,bshape,NULL),&H5Sclose,"Unable to get origin dataspace"); 01327 01328 // get file dataspace and select the desired block 01329 HDF5Handle dataspaceHandle (H5Dget_space(datasetHandle),&H5Sclose,"Unable to create target dataspace"); 01330 H5Sselect_hyperslab(dataspaceHandle, H5S_SELECT_SET, boffset, bones, bones, bshape); 01331 01332 // Write the data to the HDF5 dataset as is 01333 H5Dwrite( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, array.data()); // .data() possible since void pointer! 01334 } 01335 01336 01337 01338 // unstrided target multi array 01339 // read a block of a HDF5 dataset into a MultiArray 01340 template<unsigned int N, class T> 01341 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) 01342 { 01343 //Prepare to read without using HDF5ImportInfo 01344 //ArrayVector<hsize_t> dimshape = getDatasetShape(datasetName) ; 01345 hssize_t dimensions = getDatasetDimensions(datasetName); 01346 01347 std::string errorMessage ("HDF5File::readBlock(): Unable to open dataset '" + datasetName + "'."); 01348 HDF5Handle datasetHandle (getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str()); 01349 01350 01351 int offset = (numBandsOfType > 1); 01352 01353 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 01354 "readHDF5_block(): Array dimension disagrees with data dimension."); 01355 01356 vigra_precondition(blockShape == array.shape(), 01357 "readHDF5_block(): Array shape disagrees with block size."); 01358 01359 // hyperslab parameters for position, size, ... 01360 hsize_t boffset [N]; 01361 hsize_t bshape [N]; 01362 hsize_t bones [N]; 01363 01364 for(int i = 0; i < N; i++){ 01365 // virgra and hdf5 use different indexing 01366 boffset[i] = blockOffset[N-1-i]; 01367 //bshape[i] = blockShape[i]; 01368 bshape[i] = blockShape[N-1-i]; 01369 //boffset[i] = blockOffset[N-1-i]; 01370 bones[i] = 1; 01371 } 01372 01373 // create a target dataspace in memory with the shape of the desired block 01374 HDF5Handle memspace_handle (H5Screate_simple(N,bshape,NULL),&H5Sclose,"Unable to create target dataspace"); 01375 01376 // get file dataspace and select the desired block 01377 HDF5Handle dataspaceHandle (H5Dget_space(datasetHandle),&H5Sclose,"Unable to get dataspace"); 01378 H5Sselect_hyperslab(dataspaceHandle, H5S_SELECT_SET, boffset, bones, bones, bshape); 01379 01380 // now read the data 01381 H5Dread( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, array.data() ); // .data() possible since void pointer! 01382 } 01383 01384 }; /* class HDF5File */ 01385 01386 01387 01388 01389 01390 01391 01392 namespace detail { 01393 01394 template <class Shape> 01395 inline void 01396 selectHyperslabs(HDF5Handle & mid1, HDF5Handle & mid2, Shape const & shape, int & counter, const int elements, const int numBandsOfType) 01397 { 01398 // select hyperslab in HDF5 file 01399 hsize_t shapeHDF5[2]; 01400 shapeHDF5[0] = 1; 01401 shapeHDF5[1] = elements; 01402 hsize_t startHDF5[2]; 01403 startHDF5[0] = 0; 01404 startHDF5[1] = counter * numBandsOfType * shape[0]; // we have to reserve space for the pixel type channel(s) 01405 hsize_t strideHDF5[2]; 01406 strideHDF5[0] = 1; 01407 strideHDF5[1] = 1; 01408 hsize_t countHDF5[2]; 01409 countHDF5[0] = 1; 01410 countHDF5[1] = numBandsOfType * shape[0]; 01411 hsize_t blockHDF5[2]; 01412 blockHDF5[0] = 1; 01413 blockHDF5[1] = 1; 01414 mid1 = HDF5Handle(H5Screate_simple(2, shapeHDF5, NULL), 01415 &H5Sclose, "unable to create hyperslabs."); 01416 H5Sselect_hyperslab(mid1, H5S_SELECT_SET, startHDF5, strideHDF5, countHDF5, blockHDF5); 01417 // select hyperslab in input data object 01418 hsize_t shapeData[2]; 01419 shapeData[0] = 1; 01420 shapeData[1] = numBandsOfType * shape[0]; 01421 hsize_t startData[2]; 01422 startData[0] = 0; 01423 startData[1] = 0; 01424 hsize_t strideData[2]; 01425 strideData[0] = 1; 01426 strideData[1] = 1; 01427 hsize_t countData[2]; 01428 countData[0] = 1; 01429 countData[1] = numBandsOfType * shape[0]; 01430 hsize_t blockData[2]; 01431 blockData[0] = 1; 01432 blockData[1] = 1; 01433 mid2 = HDF5Handle(H5Screate_simple(2, shapeData, NULL), 01434 &H5Sclose, "unable to create hyperslabs."); 01435 H5Sselect_hyperslab(mid2, H5S_SELECT_SET, startData, strideData, countData, blockData); 01436 } 01437 01438 template <class DestIterator, class Shape, class T> 01439 inline void 01440 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>) 01441 { 01442 HDF5Handle mid1, mid2; 01443 01444 // select hyperslabs 01445 selectHyperslabs(mid1, mid2, shape, counter, elements, numBandsOfType); 01446 01447 // read from hdf5 01448 H5Dread(dataset_id, datatype, mid2, mid1, H5P_DEFAULT, buffer.data()); 01449 01450 // increase counter 01451 counter++; 01452 01453 01454 //std::cout << "numBandsOfType: " << numBandsOfType << std::endl; 01455 DestIterator dend = d + shape[0]; 01456 int k = 0; 01457 for(; d < dend; ++d, k++) 01458 { 01459 *d = buffer[k]; 01460 //std::cout << buffer[k] << "| "; 01461 } 01462 01463 } 01464 01465 template <class DestIterator, class Shape, class T, int N> 01466 void 01467 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>) 01468 { 01469 DestIterator dend = d + shape[N]; 01470 for(; d < dend; ++d) 01471 { 01472 readHDF5Impl(d.begin(), shape, dataset_id, datatype, buffer, counter, elements, numBandsOfType, MetaInt<N-1>()); 01473 } 01474 } 01475 01476 } // namespace detail 01477 01478 /** \brief Read the data specified by the given \ref vigra::HDF5ImportInfo object 01479 and write the into the given 'array'. 01480 01481 The array must have the correct number of dimensions and shape for the dataset 01482 represented by 'info'. When the element type of 'array' differs from the stored element 01483 type, HDF5 will convert the type on the fly (except when the HDF5 version is 1.6 or below, 01484 in which case an error will result). Multi-channel element types (i.e. \ref vigra::RGBValue 01485 and \ref vigra::TinyVector) are recognized and handled correctly. 01486 01487 <b> Declaration:</b> 01488 01489 \code 01490 namespace vigra { 01491 template<unsigned int N, class T, class StrideTag> 01492 void 01493 readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, StrideTag> array); 01494 } 01495 \endcode 01496 01497 <b> Usage:</b> 01498 01499 <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br> 01500 Namespace: vigra 01501 01502 \code 01503 01504 HDF5ImportInfo info(filename, dataset_name); 01505 vigra_precondition(info.numDimensions() == 3, "Dataset must be 3-dimensional."); 01506 01507 MultiArrayShape<3>::type shape(info.shape().begin()); 01508 MultiArray<3, int> array(shape); 01509 01510 readHDF5(info, array); 01511 \endcode 01512 */ 01513 doxygen_overloaded_function(template <...> void readHDF5) 01514 01515 // scalar and unstrided target multi array 01516 template<unsigned int N, class T> 01517 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, UnstridedArrayTag> array) // scalar 01518 { 01519 readHDF5(info, array, detail::getH5DataType<T>(), 1); 01520 } 01521 01522 // non-scalar (TinyVector) and unstrided target multi array 01523 template<unsigned int N, class T, int SIZE> 01524 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> array) 01525 { 01526 readHDF5(info, array, detail::getH5DataType<T>(), SIZE); 01527 } 01528 01529 // non-scalar (RGBValue) and unstrided target multi array 01530 template<unsigned int N, class T> 01531 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> array) 01532 { 01533 readHDF5(info, array, detail::getH5DataType<T>(), 3); 01534 } 01535 01536 // unstrided target multi array 01537 template<unsigned int N, class T> 01538 void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, UnstridedArrayTag> array, const hid_t datatype, const int numBandsOfType) 01539 { 01540 int offset = (numBandsOfType > 1); 01541 01542 //std::cout << "offset: " << offset << ", N: " << N << ", dims: " << info.numDimensions() << std::endl; 01543 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 01544 "readHDF5(): Array dimension disagrees with HDF5ImportInfo.numDimensions()."); 01545 01546 typename MultiArrayShape<N>::type shape; 01547 for(int k=offset; k<info.numDimensions(); ++k) { 01548 shape[k-offset] = info.shapeOfDimension(k); 01549 } 01550 01551 vigra_precondition(shape == array.shape(), 01552 "readHDF5(): Array shape disagrees with HDF5ImportInfo."); 01553 01554 // simply read in the data as is 01555 H5Dread( info.getDatasetHandle(), datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data() ); // .data() possible since void pointer! 01556 } 01557 01558 // scalar and strided target multi array 01559 template<unsigned int N, class T> 01560 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, StridedArrayTag> array) // scalar 01561 { 01562 readHDF5(info, array, detail::getH5DataType<T>(), 1); 01563 } 01564 01565 // non-scalar (TinyVector) and strided target multi array 01566 template<unsigned int N, class T, int SIZE> 01567 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, TinyVector<T, SIZE>, StridedArrayTag> array) 01568 { 01569 readHDF5(info, array, detail::getH5DataType<T>(), SIZE); 01570 } 01571 01572 // non-scalar (RGBValue) and strided target multi array 01573 template<unsigned int N, class T> 01574 inline void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, RGBValue<T>, StridedArrayTag> array) 01575 { 01576 readHDF5(info, array, detail::getH5DataType<T>(), 3); 01577 } 01578 01579 // strided target multi array 01580 template<unsigned int N, class T> 01581 void readHDF5(const HDF5ImportInfo &info, MultiArrayView<N, T, StridedArrayTag> array, const hid_t datatype, const int numBandsOfType) 01582 { 01583 int offset = (numBandsOfType > 1); 01584 01585 //std::cout << "offset: " << offset << ", N: " << N << ", dims: " << info.numDimensions() << std::endl; 01586 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 01587 "readHDF5(): Array dimension disagrees with HDF5ImportInfo.numDimensions()."); 01588 01589 typename MultiArrayShape<N>::type shape; 01590 for(int k=offset; k<info.numDimensions(); ++k) { 01591 shape[k-offset] = info.shapeOfDimension(k); 01592 } 01593 01594 vigra_precondition(shape == array.shape(), 01595 "readHDF5(): Array shape disagrees with HDF5ImportInfo."); 01596 01597 //Get the data 01598 int counter = 0; 01599 int elements = numBandsOfType; 01600 for(unsigned int i=0;i<N;++i) 01601 elements *= shape[i]; 01602 ArrayVector<T> buffer(shape[0]); 01603 detail::readHDF5Impl(array.traverser_begin(), shape, info.getDatasetHandle(), datatype, buffer, counter, elements, numBandsOfType, vigra::MetaInt<N-1>()); 01604 } 01605 01606 inline hid_t openGroup(hid_t parent, std::string group_name) 01607 { 01608 //std::cout << group_name << std::endl; 01609 size_t last_slash = group_name.find_last_of('/'); 01610 if (last_slash == std::string::npos || last_slash != group_name.size() - 1) 01611 group_name = group_name + '/'; 01612 std::string::size_type begin = 0, end = group_name.find('/'); 01613 int ii = 0; 01614 while (end != std::string::npos) 01615 { 01616 std::string group(group_name.begin()+begin, group_name.begin()+end); 01617 hid_t prev_parent = parent; 01618 parent = H5Gopen(prev_parent, group.c_str(), H5P_DEFAULT); 01619 01620 if(ii != 0) H5Gclose(prev_parent); 01621 if(parent < 0) return parent; 01622 ++ii; 01623 begin = end + 1; 01624 end = group_name.find('/', begin); 01625 } 01626 return parent; 01627 } 01628 01629 inline hid_t createGroup(hid_t parent, std::string group_name) 01630 { 01631 if(group_name.size() == 0 ||*group_name.rbegin() != '/') 01632 group_name = group_name + '/'; 01633 if(group_name == "/") 01634 return H5Gopen(parent, group_name.c_str(), H5P_DEFAULT); 01635 01636 std::string::size_type begin = 0, end = group_name.find('/'); 01637 int ii = 0; 01638 while (end != std::string::npos) 01639 { 01640 std::string group(group_name.begin()+begin, group_name.begin()+end); 01641 hid_t prev_parent = parent; 01642 01643 if(H5LTfind_dataset(parent, group.c_str()) == 0) 01644 { 01645 parent = H5Gcreate(prev_parent, group.c_str(), H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); 01646 } else { 01647 parent = H5Gopen(prev_parent, group.c_str(), H5P_DEFAULT); 01648 } 01649 01650 if(ii != 0) H5Gclose(prev_parent); 01651 if(parent < 0) return parent; 01652 ++ii; 01653 begin = end + 1; 01654 end = group_name.find('/', begin); 01655 } 01656 return parent; 01657 } 01658 01659 inline void deleteDataset(hid_t parent, std::string dataset_name) 01660 { 01661 // delete existing data and create new dataset 01662 if(H5LTfind_dataset(parent, dataset_name.c_str())) 01663 { 01664 //std::cout << "dataset already exists" << std::endl; 01665 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6) 01666 if(H5Gunlink(parent, dataset_name.c_str()) < 0) 01667 { 01668 vigra_postcondition(false, "writeToHDF5File(): Unable to delete existing data."); 01669 } 01670 #else 01671 if(H5Ldelete(parent, dataset_name.c_str(), H5P_DEFAULT ) < 0) 01672 { 01673 vigra_postcondition(false, "createDataset(): Unable to delete existing data."); 01674 } 01675 #endif 01676 } 01677 } 01678 01679 inline hid_t createFile(std::string filePath, bool append_ = true) 01680 { 01681 FILE * pFile; 01682 pFile = fopen ( filePath.c_str(), "r" ); 01683 hid_t file_id; 01684 if ( pFile == NULL ) 01685 { 01686 file_id = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT); 01687 } 01688 else if(append_) 01689 { 01690 fclose( pFile ); 01691 file_id = H5Fopen(filePath.c_str(), H5F_ACC_RDWR, H5P_DEFAULT); 01692 } 01693 else 01694 { 01695 fclose(pFile); 01696 std::remove(filePath.c_str()); 01697 file_id = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT); 01698 } 01699 return file_id; 01700 } 01701 01702 template<unsigned int N, class T, class Tag> 01703 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) 01704 { 01705 std::string path_name(pathInFile), group_name, data_set_name, message; 01706 std::string::size_type delimiter = path_name.rfind('/'); 01707 01708 //create or open file 01709 file_handle = HDF5Handle(createFile(filePath), &H5Fclose, 01710 "createDataset(): unable to open output file."); 01711 01712 // get the groupname and the filename 01713 if(delimiter == std::string::npos) 01714 { 01715 group_name = "/"; 01716 data_set_name = path_name; 01717 } 01718 else 01719 { 01720 group_name = std::string(path_name.begin(), path_name.begin()+delimiter); 01721 data_set_name = std::string(path_name.begin()+delimiter+1, path_name.end()); 01722 } 01723 01724 // create all groups 01725 HDF5Handle group(createGroup(file_handle, group_name), &H5Gclose, 01726 "createDataset(): Unable to create and open group. generic v"); 01727 01728 // delete the dataset if it already exists 01729 deleteDataset(group, data_set_name); 01730 01731 // create dataspace 01732 // add an extra dimension in case that the data is non-scalar 01733 HDF5Handle dataspace_handle; 01734 if(numBandsOfType > 1) { 01735 // invert dimensions to guarantee c-order 01736 hsize_t shape_inv[N+1]; // one additional dimension for pixel type channel(s) 01737 for(unsigned int k=0; k<N; ++k) { 01738 shape_inv[N-1-k] = array.shape(k); // the channels (eg of an RGB image) are represented by the first dimension (before inversion) 01739 //std::cout << shape_inv[N-k] << " (" << N << ")"; 01740 } 01741 shape_inv[N] = numBandsOfType; 01742 01743 // create dataspace 01744 dataspace_handle = HDF5Handle(H5Screate_simple(N+1, shape_inv, NULL), 01745 &H5Sclose, "createDataset(): unable to create dataspace for non-scalar data."); 01746 } else { 01747 // invert dimensions to guarantee c-order 01748 hsize_t shape_inv[N]; 01749 for(unsigned int k=0; k<N; ++k) 01750 shape_inv[N-1-k] = array.shape(k); 01751 01752 // create dataspace 01753 dataspace_handle = HDF5Handle(H5Screate_simple(N, shape_inv, NULL), 01754 &H5Sclose, "createDataset(): unable to create dataspace for scalar data."); 01755 } 01756 01757 //alloc memory for dataset. 01758 dataset_handle = HDF5Handle(H5Dcreate(group, 01759 data_set_name.c_str(), 01760 datatype, 01761 dataspace_handle, 01762 H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT), 01763 &H5Dclose, "createDataset(): unable to create dataset."); 01764 } 01765 01766 01767 01768 namespace detail { 01769 01770 template <class DestIterator, class Shape, class T> 01771 inline void 01772 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>) 01773 { 01774 DestIterator dend = d + (typename DestIterator::difference_type)shape[0]; 01775 int k = 0; 01776 //std::cout << "new:" << std::endl; 01777 for(; d < dend; ++d, k++) 01778 { 01779 buffer[k] = *d; 01780 //std::cout << buffer[k] << " "; 01781 } 01782 //std::cout << std::endl; 01783 HDF5Handle mid1, mid2; 01784 01785 // select hyperslabs 01786 selectHyperslabs(mid1, mid2, shape, counter, elements, numBandsOfType); 01787 01788 // write to hdf5 01789 H5Dwrite(dataset_id, datatype, mid2, mid1, H5P_DEFAULT, buffer.data()); 01790 // increase counter 01791 counter++; 01792 } 01793 01794 template <class DestIterator, class Shape, class T, int N> 01795 void 01796 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>) 01797 { 01798 DestIterator dend = d + (typename DestIterator::difference_type)shape[N]; 01799 for(; d < dend; ++d) 01800 { 01801 writeHDF5Impl(d.begin(), shape, dataset_id, datatype, buffer, counter, elements, numBandsOfType, MetaInt<N-1>()); 01802 } 01803 } 01804 01805 } // namespace detail 01806 01807 /** \brief Store array data in an HDF5 file. 01808 01809 The number of dimensions, shape and element type of the stored dataset is automatically 01810 determined from the properties of the given \a array. Strided arrays are stored in an 01811 unstrided way, i.e. in contiguous scan-order. Multi-channel element types 01812 (i.e. \ref vigra::RGBValue and \ref vigra::TinyVector) are recognized and handled correctly 01813 (in particular, the will form the innermost dimension of the stored dataset). 01814 \a pathInFile may contain '/'-separated group names, but must end with the name 01815 of the dataset to be created. 01816 01817 <b> Declaration:</b> 01818 01819 \code 01820 namespace vigra { 01821 template<unsigned int N, class T, class StrideTag> 01822 void 01823 writeHDF5(const char* filePath, const char* pathInFile, 01824 MultiArrayView<N, T, StrideTag>const & array); 01825 } 01826 \endcode 01827 01828 <b> Usage:</b> 01829 01830 <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br> 01831 Namespace: vigra 01832 01833 \code 01834 MultiArrayShape<3>::type shape(100, 200, 20); 01835 MultiArray<3, int> array(shape); 01836 ... // fill array with data 01837 01838 writeHDF5("mydata.h5", "/group1/my_dataset", array); 01839 \endcode 01840 */ 01841 doxygen_overloaded_function(template <...> void writeHDF5) 01842 01843 // scalar and unstrided multi arrays 01844 template<unsigned int N, class T> 01845 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, UnstridedArrayTag> & array) // scalar 01846 { 01847 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 1); 01848 } 01849 01850 // non-scalar (TinyVector) and unstrided multi arrays 01851 template<unsigned int N, class T, int SIZE> 01852 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, TinyVector<T, SIZE>, UnstridedArrayTag> & array) 01853 { 01854 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), SIZE); 01855 } 01856 01857 // non-scalar (RGBValue) and unstrided multi arrays 01858 template<unsigned int N, class T> 01859 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, RGBValue<T>, UnstridedArrayTag> & array) 01860 { 01861 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 3); 01862 } 01863 01864 // unstrided multi arrays 01865 template<unsigned int N, class T> 01866 void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, UnstridedArrayTag> & array, const hid_t datatype, const int numBandsOfType) 01867 { 01868 HDF5Handle file_handle; 01869 HDF5Handle dataset_handle; 01870 createDataset(filePath, pathInFile, array, datatype, numBandsOfType, file_handle, dataset_handle); 01871 01872 // Write the data to the HDF5 dataset as is 01873 H5Dwrite( dataset_handle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data()); // .data() possible since void pointer! 01874 01875 H5Fflush(file_handle, H5F_SCOPE_GLOBAL); 01876 } 01877 01878 01879 // scalar and strided multi arrays 01880 template<unsigned int N, class T> 01881 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, StridedArrayTag> & array) // scalar 01882 { 01883 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 1); 01884 } 01885 01886 // non-scalar (TinyVector) and strided multi arrays 01887 template<unsigned int N, class T, int SIZE> 01888 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, TinyVector<T, SIZE>, StridedArrayTag> & array) 01889 { 01890 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), SIZE); 01891 } 01892 01893 // non-scalar (RGBValue) and strided multi arrays 01894 template<unsigned int N, class T> 01895 inline void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, RGBValue<T>, StridedArrayTag> & array) 01896 { 01897 writeHDF5(filePath, pathInFile, array, detail::getH5DataType<T>(), 3); 01898 } 01899 01900 // strided multi arrays 01901 template<unsigned int N, class T> 01902 void writeHDF5(const char* filePath, const char* pathInFile, const MultiArrayView<N, T, StridedArrayTag> & array, const hid_t datatype, const int numBandsOfType) 01903 { 01904 HDF5Handle file_handle; 01905 HDF5Handle dataset_handle; 01906 createDataset(filePath, pathInFile, array, datatype, numBandsOfType, file_handle, dataset_handle); 01907 01908 vigra::TinyVector<int,N> shape; 01909 vigra::TinyVector<int,N> stride; 01910 int elements = numBandsOfType; 01911 for(unsigned int k=0; k<N; ++k) 01912 { 01913 shape[k] = array.shape(k); 01914 stride[k] = array.stride(k); 01915 elements *= (int)shape[k]; 01916 } 01917 int counter = 0; 01918 01919 ArrayVector<T> buffer((int)array.shape(0)); 01920 detail::writeHDF5Impl(array.traverser_begin(), shape, dataset_handle, datatype, buffer, counter, elements, numBandsOfType, vigra::MetaInt<N-1>()); 01921 01922 H5Fflush(file_handle, H5F_SCOPE_GLOBAL); 01923 01924 } 01925 01926 namespace detail 01927 { 01928 struct MaxSizeFnc 01929 { 01930 size_t size; 01931 01932 MaxSizeFnc() 01933 : size(0) 01934 {} 01935 01936 void operator()(std::string const & in) 01937 { 01938 size = in.size() > size ? 01939 in.size() : 01940 size; 01941 } 01942 }; 01943 } 01944 01945 01946 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR == 8) || DOXYGEN 01947 /** Write a numeric MultiArray as an attribute with name \a name 01948 of the dataset specified by the handle \a loc. 01949 01950 <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br> 01951 Namespace: vigra 01952 */ 01953 template<size_t N, class T, class C> 01954 void writeHDF5Attr(hid_t loc, 01955 const char* name, 01956 MultiArrayView<N, T, C> const & array) 01957 { 01958 if(H5Aexists(loc, name) > 0) 01959 H5Adelete(loc, name); 01960 01961 ArrayVector<hsize_t> shape(array.shape().begin(), 01962 array.shape().end()); 01963 HDF5Handle 01964 dataspace_handle(H5Screate_simple(N, shape.data(), NULL), 01965 &H5Sclose, 01966 "writeToHDF5File(): unable to create dataspace."); 01967 01968 HDF5Handle attr(H5Acreate(loc, 01969 name, 01970 detail::getH5DataType<T>(), 01971 dataspace_handle, 01972 H5P_DEFAULT ,H5P_DEFAULT ), 01973 &H5Aclose, 01974 "writeHDF5Attr: unable to create Attribute"); 01975 01976 //copy data - since attributes are small - who cares! 01977 ArrayVector<T> buffer; 01978 for(int ii = 0; ii < array.size(); ++ii) 01979 buffer.push_back(array[ii]); 01980 H5Awrite(attr, detail::getH5DataType<T>(), buffer.data()); 01981 } 01982 01983 /** Write a string MultiArray as an attribute with name \a name 01984 of the dataset specified by the handle \a loc. 01985 01986 <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br> 01987 Namespace: vigra 01988 */ 01989 template<size_t N, class C> 01990 void writeHDF5Attr(hid_t loc, 01991 const char* name, 01992 MultiArrayView<N, std::string, C> const & array) 01993 { 01994 if(H5Aexists(loc, name) > 0) 01995 H5Adelete(loc, name); 01996 01997 ArrayVector<hsize_t> shape(array.shape().begin(), 01998 array.shape().end()); 01999 HDF5Handle 02000 dataspace_handle(H5Screate_simple(N, shape.data(), NULL), 02001 &H5Sclose, 02002 "writeToHDF5File(): unable to create dataspace."); 02003 02004 HDF5Handle atype(H5Tcopy (H5T_C_S1), 02005 &H5Tclose, 02006 "writeToHDF5File(): unable to create type."); 02007 02008 detail::MaxSizeFnc max_size; 02009 max_size = std::for_each(array.data(),array.data()+ array.size(), max_size); 02010 H5Tset_size (atype, max_size.size); 02011 02012 HDF5Handle attr(H5Acreate(loc, 02013 name, 02014 atype, 02015 dataspace_handle, 02016 H5P_DEFAULT ,H5P_DEFAULT ), 02017 &H5Aclose, 02018 "writeHDF5Attr: unable to create Attribute"); 02019 02020 std::string buf =""; 02021 for(int ii = 0; ii < array.size(); ++ii) 02022 { 02023 buf = buf + array[ii] 02024 + std::string(max_size.size - array[ii].size(), ' '); 02025 } 02026 H5Awrite(attr, atype, buf.c_str()); 02027 } 02028 02029 /** Write a numeric ArrayVectorView as an attribute with name \a name 02030 of the dataset specified by the handle \a loc. 02031 02032 <b>\#include</b> <<a href="hdf5impex_8hxx_source.html">vigra/hdf5impex.hxx</a>><br> 02033 Namespace: vigra 02034 */ 02035 template<class T> 02036 inline void writeHDF5Attr( hid_t loc, 02037 const char* name, 02038 ArrayVectorView<T> & array) 02039 { 02040 writeHDF5Attr(loc, name, 02041 MultiArrayView<1, T>(MultiArrayShape<1>::type(array.size()), 02042 array.data())); 02043 } 02044 02045 /** write an Attribute given a file and a path in the file. 02046 * the path in the file should have the format 02047 * [attribute] or /[subgroups/]dataset.attribute or 02048 * /[subgroups/]group.attribute. 02049 * The attribute is written to the root group, a dataset or a subgroup 02050 * respectively 02051 */ 02052 template<class Arr> 02053 inline void writeHDF5Attr( std::string filePath, 02054 std::string pathInFile, 02055 Arr & ar) 02056 { 02057 std::string path_name(pathInFile), group_name, data_set_name, message, attr_name; 02058 std::string::size_type delimiter = path_name.rfind('/'); 02059 02060 //create or open file 02061 HDF5Handle file_id(createFile(filePath), &H5Fclose, 02062 "writeToHDF5File(): unable to open output file."); 02063 02064 // get the groupname and the filename 02065 if(delimiter == std::string::npos) 02066 { 02067 group_name = "/"; 02068 data_set_name = path_name; 02069 } 02070 02071 else 02072 { 02073 group_name = std::string(path_name.begin(), path_name.begin()+delimiter); 02074 data_set_name = std::string(path_name.begin()+delimiter+1, path_name.end()); 02075 } 02076 delimiter = data_set_name.rfind('.'); 02077 if(delimiter == std::string::npos) 02078 { 02079 attr_name = path_name; 02080 data_set_name = "/"; 02081 } 02082 else 02083 { 02084 attr_name = std::string(data_set_name.begin()+delimiter+1, data_set_name.end()); 02085 data_set_name = std::string(data_set_name.begin(), data_set_name.begin()+delimiter); 02086 } 02087 02088 HDF5Handle group(openGroup(file_id, group_name), &H5Gclose, 02089 "writeToHDF5File(): Unable to create and open group. attr ver"); 02090 02091 if(data_set_name != "/") 02092 { 02093 HDF5Handle dset(H5Dopen(group, data_set_name.c_str(), H5P_DEFAULT), &H5Dclose, 02094 "writeHDF5Attr():unable to open dataset"); 02095 writeHDF5Attr(hid_t(dset), attr_name.c_str(), ar); 02096 } 02097 else 02098 { 02099 writeHDF5Attr(hid_t(group), attr_name.c_str(), ar); 02100 } 02101 02102 } 02103 #endif 02104 02105 //@} 02106 02107 } // namespace vigra 02108 02109 #endif // VIGRA_HDF5IMPEX_HXX
© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de) |
html generated using doxygen and Python
|