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

vigra/hdf5impex.hxx VIGRA

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

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

html generated using doxygen and Python
vigra 1.9.0 (Tue Nov 6 2012)