36 #ifndef VIGRA_HDF5IMPEX_HXX 
   37 #define VIGRA_HDF5IMPEX_HXX 
   43 #define H5Gcreate_vers 2 
   44 #define H5Gopen_vers 2 
   45 #define H5Dopen_vers 2 
   46 #define H5Dcreate_vers 2 
   47 #define H5Acreate_vers 2 
   48 #define H5Eset_auto_vers 2 
   49 #define H5Eget_auto_vers 2 
   53 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6) 
   55 #   define H5Gopen(a, b, c) H5Gopen(a, b) 
   58 #  define H5Gcreate(a, b, c, d, e) H5Gcreate(a, b, 1) 
   61 #  define H5Dopen(a, b, c) H5Dopen(a, b) 
   64 #  define H5Dcreate(a, b, c, d, e, f, g) H5Dcreate(a, b, c, d, f) 
   67 #  define H5Acreate(a, b, c, d, e, f) H5Acreate(a, b, c, d, e) 
   69 # ifndef H5Pset_obj_track_times 
   70 #  define H5Pset_obj_track_times(a, b) do {} while (0) 
   78 #include "multi_array.hxx" 
   79 #include "multi_iterator_coupled.hxx" 
   80 #include "multi_impex.hxx" 
  102 inline bool isHDF5(
char const * filename)
 
  105     return _access(filename, 0) != -1 && H5Fis_hdf5(filename);
 
  107     return access(filename, F_OK) == 0 && H5Fis_hdf5(filename);
 
  132     H5E_auto1_t old_func1_;
 
  133     H5E_auto2_t old_func2_;
 
  134     void *old_client_data_;
 
  135     int error_handler_version_;
 
  144     , old_client_data_(0)
 
  145     , error_handler_version_(-1)
 
  147         if(H5Eget_auto2(H5E_DEFAULT, &old_func2_, &old_client_data_) >= 0)
 
  150             H5Eset_auto2(H5E_DEFAULT, NULL, NULL);
 
  151             error_handler_version_ = 2;
 
  153         else if(H5Eget_auto1(&old_func1_, &old_client_data_) >= 0)
 
  158             H5Eset_auto1(NULL, NULL);
 
  159             error_handler_version_ = 1;
 
  165         if(error_handler_version_ == 1)
 
  166             H5Eset_auto1(old_func1_, old_client_data_);
 
  167         else if(error_handler_version_ == 2)
 
  168             H5Eset_auto2(H5E_DEFAULT, old_func2_, old_client_data_);
 
  213     typedef herr_t (*Destructor)(hid_t);
 
  217     Destructor destructor_;
 
  248     HDF5Handle(hid_t h, Destructor destructor, 
const char * error_message)
 
  250       destructor_(destructor)
 
  253             vigra_fail(error_message);
 
  261     : handle_( h.handle_ ),
 
  262       destructor_(h.destructor_)
 
  273         if(h.handle_ != handle_)
 
  277             destructor_ = h.destructor_;
 
  299         if(handle_ && destructor_)
 
  300             res = (*destructor_)(handle_);
 
  321     void reset(hid_t h, Destructor destructor, 
const char * error_message)
 
  324             vigra_fail(error_message);
 
  329             destructor_ = destructor;
 
  339         std::swap(handle_, h.handle_);
 
  340         std::swap(destructor_, h.destructor_);
 
  358     operator hid_t()
 const 
  367         return handle_ == h.handle_;
 
  381         return handle_ != h.handle_;
 
  434     typedef herr_t (*Destructor)(hid_t);
 
  438     Destructor destructor_;
 
  473       destructor_(destructor),
 
  477             vigra_fail(error_message);
 
  479             refcount_ = 
new size_t(1);
 
  486     : handle_( h.handle_ ),
 
  487       destructor_(h.destructor_),
 
  488       refcount_(h.refcount_)
 
  500         if(h.handle_ != handle_)
 
  504             destructor_ = h.destructor_;
 
  505             refcount_ = h.refcount_;
 
  536                     res = (*destructor_)(handle_);
 
  550     void reset(hid_t h, Destructor destructor, 
const char * error_message)
 
  553             vigra_fail(error_message);
 
  558             destructor_ = destructor;
 
  560                 refcount_ = 
new size_t(1);
 
  588         std::swap(handle_, h.handle_);
 
  589         std::swap(destructor_, h.destructor_);
 
  590         std::swap(refcount_, h.refcount_);
 
  608     operator hid_t()
 const 
  617         return handle_ == h.handle_;
 
  631         return handle_ != h.handle_;
 
  685     enum PixelType { UINT8, UINT16, UINT32, UINT64,
 
  686                      INT8, INT16, INT32, INT64,
 
  699     VIGRA_EXPORT 
HDF5ImportInfo( 
const char* filePath, 
const char* pathInFile );
 
  701     VIGRA_EXPORT ~HDF5ImportInfo();
 
  705     VIGRA_EXPORT 
const std::string& 
getFilePath() 
const;
 
  778     VIGRA_EXPORT PixelType 
pixelType() 
const;
 
  782     std::string m_filename, m_path, m_pixeltype;
 
  783     hssize_t m_dimensions;
 
  791 struct HDF5TypeTraits
 
  793     static hid_t getH5DataType()
 
  795         std::runtime_error(
"getH5DataType(): invalid type");
 
  799     static int numberOfBands()
 
  801         std::runtime_error(
"numberOfBands(): invalid type");
 
  807 inline hid_t getH5DataType()
 
  809     return HDF5TypeTraits<T>::getH5DataType();
 
  812 #define VIGRA_H5_DATATYPE(type, h5type) \ 
  814 struct HDF5TypeTraits<type> \ 
  816     static hid_t getH5DataType() \ 
  820     static int numberOfBands() \ 
  824     typedef type value_type; \ 
  827 struct HDF5TypeTraits<TinyVector<type, M> > \ 
  829     static hid_t getH5DataType() \ 
  833     static int numberOfBands() \ 
  837     typedef type value_type; \ 
  840 struct HDF5TypeTraits<RGBValue<type> > \ 
  842     static hid_t getH5DataType() \ 
  846     static int numberOfBands() \ 
  850     typedef type value_type; \ 
  853 VIGRA_H5_DATATYPE(
char, H5T_NATIVE_CHAR)
 
  854 VIGRA_H5_DATATYPE(
signed char, H5T_NATIVE_SCHAR)
 
  855 VIGRA_H5_DATATYPE(
unsigned char, H5T_NATIVE_UCHAR)
 
  856 VIGRA_H5_DATATYPE(
signed short, H5T_NATIVE_SHORT)
 
  857 VIGRA_H5_DATATYPE(
unsigned short, H5T_NATIVE_USHORT)
 
  858 VIGRA_H5_DATATYPE(
signed int, H5T_NATIVE_INT)
 
  859 VIGRA_H5_DATATYPE(
unsigned int, H5T_NATIVE_UINT)
 
  860 VIGRA_H5_DATATYPE(
signed long, H5T_NATIVE_LONG)
 
  861 VIGRA_H5_DATATYPE(
unsigned long, H5T_NATIVE_ULONG)
 
  862 VIGRA_H5_DATATYPE(
signed long long, H5T_NATIVE_LLONG)
 
  863 VIGRA_H5_DATATYPE(
unsigned long long, H5T_NATIVE_ULLONG)
 
  864 VIGRA_H5_DATATYPE(
float, H5T_NATIVE_FLOAT)
 
  865 VIGRA_H5_DATATYPE(
double, H5T_NATIVE_DOUBLE)
 
  866 VIGRA_H5_DATATYPE(
long double, H5T_NATIVE_LDOUBLE)
 
  870 struct HDF5TypeTraits<
char*>
 
  872     static hid_t getH5DataType()
 
  874         hid_t stringtype = H5Tcopy (H5T_C_S1);
 
  875         H5Tset_size(stringtype, H5T_VARIABLE);
 
  879     static int numberOfBands()
 
  886 struct HDF5TypeTraits<const char*>
 
  888     static hid_t getH5DataType()
 
  890         hid_t stringtype = H5Tcopy (H5T_C_S1);
 
  891         H5Tset_size(stringtype, H5T_VARIABLE);
 
  895     static int numberOfBands()
 
  901 #undef VIGRA_H5_DATATYPE 
  905 inline hid_t getH5DataType<FFTWComplex<float> >()
 
  907     hid_t complex_id = H5Tcreate (H5T_COMPOUND, 
sizeof (FFTWComplex<float>));
 
  908     H5Tinsert (complex_id, 
"real", 0, H5T_NATIVE_FLOAT);
 
  909     H5Tinsert (complex_id, 
"imaginary", 
sizeof(
float), H5T_NATIVE_FLOAT);
 
  914 inline hid_t getH5DataType<FFTWComplex<double> >()
 
  916     hid_t complex_id = H5Tcreate (H5T_COMPOUND, 
sizeof (FFTWComplex<double>));
 
  917     H5Tinsert (complex_id, 
"real", 0, H5T_NATIVE_DOUBLE);
 
  918     H5Tinsert (complex_id, 
"imaginary", 
sizeof(
double), H5T_NATIVE_DOUBLE);
 
  927 void HDF5_ls_insert(
void*, 
const std::string &);
 
  932 VIGRA_EXPORT H5O_type_t HDF5_get_type(hid_t, 
const char*);
 
  933 extern "C" VIGRA_EXPORT herr_t HDF5_ls_inserter_callback(hid_t, 
const char*, 
const H5L_info_t*, 
void*);
 
  934 extern "C" VIGRA_EXPORT herr_t HDF5_listAttributes_inserter_callback(hid_t, 
const char*, 
const H5A_info_t*, 
void*);
 
  991         virtual void insert(
const std::string &) = 0;
 
  992         virtual ~ls_closure() {}
 
  996     struct lsOpData : 
public ls_closure
 
  998         std::vector<std::string> & objects;
 
  999         lsOpData(std::vector<std::string> & o) : objects(o) {}
 
 1000         void insert(
const std::string & x)
 
 1002             objects.push_back(x);
 
 1007     template<
class Container>
 
 1008     struct ls_container_data : 
public ls_closure
 
 1010         Container & objects;
 
 1011         ls_container_data(Container & o) : objects(o) {}
 
 1012         void insert(
const std::string & x)
 
 1014             objects.insert(std::string(x));
 
 1021     friend void HDF5_ls_insert(
void*, 
const std::string &);
 
 1036         ReadOnly = OpenReadOnly, 
 
 1056     : track_time(track_creation_times ? 1 : 0),
 
 1065     explicit HDF5File(std::string filePath, 
OpenMode mode = ReadOnly, 
bool track_creation_times = 
false)
 
 1066         : track_time(track_creation_times ? 1 : 0)
 
 1068         open(filePath, mode);
 
 1076     explicit HDF5File(
char const * filePath, 
OpenMode mode = ReadOnly, 
bool track_creation_times = 
false)
 
 1077         : track_time(track_creation_times ? 1 : 0)
 
 1079         open(std::string(filePath), mode);
 
 1099                       const std::string & pathname = 
"",
 
 1100                       bool read_only = 
false)
 
 1101     : fileHandle_(fileHandle),
 
 1102       read_only_(read_only)
 
 1109         cGroupHandle_ = 
HDF5Handle(openCreateGroup_(pathname), &H5Gclose,
 
 1110                 "HDF5File(fileHandle, pathname): Failed to open group");
 
 1113         hbool_t track_times_tmp;
 
 1114         HDF5Handle plist_id(H5Fget_create_plist(fileHandle_), &H5Pclose,
 
 1115                 "HDF5File(fileHandle, pathname): Failed to open file creation property list");
 
 1116         herr_t status = H5Pget_obj_track_times(plist_id, &track_times_tmp );
 
 1117         vigra_postcondition(status >= 0,
 
 1118             "HDF5File(fileHandle, pathname): cannot access track time attribute");
 
 1119         track_time = track_times_tmp;
 
 1127     : fileHandle_(other.fileHandle_),
 
 1128       track_time(other.track_time),
 
 1129       read_only_(other.read_only_)
 
 1131         cGroupHandle_ = 
HDF5Handle(openCreateGroup_(other.currentGroupName_()), &H5Gclose,
 
 1132                                    "HDF5File(HDF5File const &): Failed to open group.");
 
 1155             fileHandle_ = other.fileHandle_;
 
 1156             cGroupHandle_ = 
HDF5Handle(openCreateGroup_(other.currentGroupName_()), &H5Gclose,
 
 1157                                        "HDF5File::operator=(): Failed to open group.");
 
 1158             track_time = other.track_time;
 
 1159             read_only_ = other.read_only_;
 
 1164     int file_use_count()
 const 
 1171         return fileHandle_ != (hid_t)0;
 
 1174     bool isReadOnly()
 const 
 1179     void setReadOnly(
bool stat=
true)
 
 1191         std::string errorMessage = 
"HDF5File.open(): Could not open or create file '" + filePath + 
"'.";
 
 1192         fileHandle_ = 
HDF5HandleShared(createFile_(filePath, mode), &H5Fclose, errorMessage.c_str());
 
 1193         cGroupHandle_ = 
HDF5Handle(openCreateGroup_(
"/"), &H5Gclose, 
"HDF5File.open(): Failed to open root group.");
 
 1194         setReadOnly(mode == OpenReadOnly);
 
 1201         bool success = cGroupHandle_.
close() >= 0 && fileHandle_.
close() >= 0;
 
 1202         vigra_postcondition(success, 
"HDF5File.close() failed.");
 
 1209         std::string message = 
"HDF5File::root(): Could not open group '/'.";
 
 1210         cGroupHandle_ = 
HDF5Handle(H5Gopen(fileHandle_, 
"/", H5P_DEFAULT),&H5Gclose,message.c_str());
 
 1216     inline void cd(std::string groupName)
 
 1227         std::string groupName = currentGroupName_();
 
 1230         if(groupName == 
"/"){
 
 1234         size_t lastSlash = groupName.find_last_of(
'/');
 
 1236         std::string parentGroup (groupName.begin(), groupName.begin()+lastSlash+1);
 
 1249         std::string groupName = currentGroupName_();
 
 1251         for(
int i = 0; i<levels; i++)
 
 1256                 if(groupName != currentGroupName_())
 
 1268     inline void mkdir(std::string groupName)
 
 1270         vigra_precondition(!isReadOnly(),
 
 1271             "HDF5File::mkdir(): file is read-only.");
 
 1273         std::string message = 
"HDF5File::mkdir(): Could not create group '" + groupName + 
"'.\n";
 
 1278         HDF5Handle(openCreateGroup_(groupName.c_str()),&H5Gclose,message.c_str());
 
 1285     inline void cd_mk(std::string groupName)
 
 1287         vigra_precondition(!isReadOnly(),
 
 1288             "HDF5File::cd_mk(): file is read-only.");
 
 1290         std::string  message = 
"HDF5File::cd_mk(): Could not create group '" + groupName + 
"'.";
 
 1295         cGroupHandle_ = 
HDF5Handle(openCreateGroup_(groupName.c_str()),&H5Gclose,message.c_str());
 
 1299     void ls_H5Literate(ls_closure & data)
 const 
 1301         H5Literate(cGroupHandle_, H5_INDEX_NAME, H5_ITER_NATIVE, NULL,
 
 1302                    HDF5_ls_inserter_callback, static_cast<void*>(&data));
 
 1310     inline std::vector<std::string> 
ls()
 const 
 1312         std::vector<std::string> list;
 
 1313         lsOpData data(list);
 
 1314         ls_H5Literate(data);
 
 1331     template<
class Container>
 
 1332     void ls(Container & cont)
 const 
 1334         ls_container_data<Container> data(cont);
 
 1335         ls_H5Literate(data);
 
 1340     inline std::string 
pwd()
 const 
 1342         return currentGroupName_();
 
 1358         return (H5Lexists(fileHandle_, datasetName.c_str(), H5P_DEFAULT) > 0);
 
 1369         return getDatasetDimensions_(datasetHandle);
 
 1372     hssize_t getDatasetDimensions_(hid_t dataset)
 const 
 1374         std::string errorMessage = 
"HDF5File::getDatasetDimensions(): Unable to access dataspace.";
 
 1375         HDF5Handle dataspaceHandle(H5Dget_space(dataset), &H5Sclose, errorMessage.c_str());
 
 1378         return H5Sget_simple_extent_ndims(dataspaceHandle);
 
 1400         std::string errorMessage = 
"HDF5File::getDatasetShape(): Unable to open dataset '" + datasetName + 
"'.";
 
 1401         HDF5Handle datasetHandle = 
HDF5Handle(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
 
 1403         errorMessage = 
"HDF5File::getDatasetShape(): Unable to access dataspace.";
 
 1404         HDF5Handle dataspaceHandle(H5Dget_space(datasetHandle), &H5Sclose, errorMessage.c_str());
 
 1410         H5Sget_simple_extent_dims(dataspaceHandle, shape.
data(), maxdims.
data());
 
 1413         std::reverse(shape.
begin(), shape.
end());
 
 1436         std::string errorMessage = 
"HDF5File::getChunkShape(): Unable to open dataset '" + datasetName + 
"'.";
 
 1437         HDF5Handle datasetHandle = 
HDF5Handle(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
 
 1439         errorMessage = 
"HDF5File::getChunkShape(): Unable to access dataspace.";
 
 1440         HDF5Handle dataspaceHandle(H5Dget_space(datasetHandle), &H5Sclose, errorMessage.c_str());
 
 1441         HDF5Handle properties(H5Dget_create_plist(datasetHandle),
 
 1442                              &H5Pclose, 
"HDF5File::read(): failed to get property list");
 
 1449         H5Pget_chunk(properties, dimensions, shape.
data());
 
 1452         std::reverse(shape.
begin(), shape.
end());
 
 1477         hid_t datatype = H5Dget_type(datasetHandle);
 
 1478         H5T_class_t dataclass = H5Tget_class(datatype);
 
 1479         size_t datasize  = H5Tget_size(datatype);
 
 1480         H5T_sign_t datasign  = H5Tget_sign(datatype);
 
 1482         if(dataclass == H5T_FLOAT)
 
 1486             else if(datasize == 8)
 
 1489         else if(dataclass == H5T_INTEGER)
 
 1491             if(datasign == H5T_SGN_NONE)
 
 1495                 else if(datasize == 2)
 
 1497                 else if(datasize == 4)
 
 1499                 else if(datasize == 8)
 
 1506                 else if(datasize == 2)
 
 1508                 else if(datasize == 4)
 
 1510                 else if(datasize == 8)
 
 1521         std::string errorMessage = 
"HDF5File::getDatasetHandle(): Unable to open dataset '" + datasetName + 
"'.";
 
 1529         std::string errorMessage = 
"HDF5File::getDatasetHandle(): Unable to open dataset '" + datasetName + 
"'.";
 
 1536                               std::string function_name = 
"HDF5File::getGroupHandle()")
 
 1538         std::string errorMessage = function_name + 
": Group '" + group_name + 
"' not found.";
 
 1544         vigra_precondition(group_name == 
"/" || H5Lexists(fileHandle_, group_name.c_str(), H5P_DEFAULT) != 0,
 
 1545                            errorMessage.c_str());
 
 1548         return HDF5Handle(openCreateGroup_(group_name), &H5Gclose, 
"Internal error");
 
 1552     void ls_H5Aiterate(std::string 
const & group_or_dataset, ls_closure & data)
 const 
 1554         H5O_type_t h5_type = get_object_type_(group_or_dataset);
 
 1555         vigra_precondition(h5_type == H5O_TYPE_GROUP || h5_type == H5O_TYPE_DATASET,
 
 1556             "HDF5File::listAttributes(): object \"" + group_or_dataset + 
"\" is neither a group nor a dataset.");
 
 1558         HDF5Handle object_handle(h5_type == H5O_TYPE_GROUP
 
 1559                                      ? const_cast<HDF5File*>(
this)->openCreateGroup_(group_or_dataset)
 
 1560                                      : getDatasetHandle_(group_or_dataset),
 
 1561                                  h5_type == H5O_TYPE_GROUP
 
 1564                                  "HDF5File::listAttributes(): unable to open object.");
 
 1566         H5Aiterate2(object_handle, H5_INDEX_NAME, H5_ITER_NATIVE, &n,
 
 1567                     HDF5_listAttributes_inserter_callback, static_cast<void*>(&data));
 
 1575     inline std::vector<std::string> 
listAttributes(std::string 
const & group_or_dataset)
 const 
 1577         std::vector<std::string> list;
 
 1578         lsOpData data(list);
 
 1579         ls_H5Aiterate(group_or_dataset, data);
 
 1589     template<
class Container>
 
 1590     void listAttributes(std::string 
const & group_or_dataset, Container & container)
 const 
 1592         ls_container_data<Container> data(container);
 
 1593         ls_H5Aiterate(group_or_dataset, data);
 
 1600         std::string message = 
"HDF5File::getAttributeHandle(): Attribute '" + attribute_name + 
"' not found.";
 
 1602                           &H5Aclose, message.c_str());
 
 1610     template<
unsigned int N, 
class T, 
class Str
ide>
 
 1612                                std::string attribute_name,
 
 1618         write_attribute_(object_name, attribute_name, array, detail::getH5DataType<T>(), 1);
 
 1621     template<
unsigned int N, 
class T, 
int SIZE, 
class Str
ide>
 
 1623                                std::string attributeName,
 
 1629         write_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), SIZE);
 
 1632     template<
unsigned int N, 
class T, 
class Str
ide>
 
 1634                                std::string attributeName,
 
 1635                                const MultiArrayView<N, RGBValue<T>, Stride> & array)
 
 1640         write_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), 3);
 
 1646     inline void writeAttribute(std::string object_name, std::string attribute_name, 
char data)
 
 1647         { writeAtomicAttribute(object_name,attribute_name,data); }
 
 1648     inline void writeAttribute(std::string datasetName, std::string attributeName, 
signed char data)
 
 1649         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1650     inline void writeAttribute(std::string datasetName, std::string attributeName, 
signed short data)
 
 1651         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1652     inline void writeAttribute(std::string datasetName, std::string attributeName, 
signed int data)
 
 1653         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1654     inline void writeAttribute(std::string datasetName, std::string attributeName, 
signed long data)
 
 1655         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1656     inline void writeAttribute(std::string datasetName, std::string attributeName, 
signed long long data)
 
 1657         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1658     inline void writeAttribute(std::string datasetName, std::string attributeName, 
unsigned char data)
 
 1659         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1660     inline void writeAttribute(std::string datasetName, std::string attributeName, 
unsigned short data)
 
 1661         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1662     inline void writeAttribute(std::string datasetName, std::string attributeName, 
unsigned int data)
 
 1663         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1664     inline void writeAttribute(std::string datasetName, std::string attributeName, 
unsigned long data)
 
 1665         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1666     inline void writeAttribute(std::string datasetName, std::string attributeName, 
unsigned long long data)
 
 1667         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1668     inline void writeAttribute(std::string datasetName, std::string attributeName, 
float data)
 
 1669         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1670     inline void writeAttribute(std::string datasetName, std::string attributeName, 
double data)
 
 1671         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1672     inline void writeAttribute(std::string datasetName, std::string attributeName, 
long double data)
 
 1673         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1674     inline void writeAttribute(std::string datasetName, std::string attributeName, 
const char* data)
 
 1675         { writeAtomicAttribute(datasetName,attributeName,data); }
 
 1676     inline void writeAttribute(std::string datasetName, std::string attributeName, std::string 
const & data)
 
 1677         { writeAtomicAttribute(datasetName,attributeName,data.c_str()); }
 
 1684         htri_t exists = H5Aexists_by_name(fileHandle_, obj_path.c_str(),
 
 1685                                           attribute_name.c_str(), H5P_DEFAULT);
 
 1686         vigra_precondition(exists >= 0, 
"HDF5File::existsAttribute(): " 
 1687                                         "object '" + object_name + 
"' " 
 1697     template<
unsigned int N, 
class T, 
class Str
ide>
 
 1699                               std::string attribute_name,
 
 1705         read_attribute_(object_name, attribute_name, array, detail::getH5DataType<T>(), 1);
 
 1708     template<
unsigned int N, 
class T, 
int SIZE, 
class Str
ide>
 
 1710                               std::string attributeName,
 
 1716         read_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), SIZE);
 
 1719     template<
unsigned int N, 
class T, 
class Str
ide>
 
 1721                               std::string attributeName,
 
 1722                               MultiArrayView<N, RGBValue<T>, Stride> array)
 
 1727         read_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), 3);
 
 1733     inline void readAttribute(std::string object_name, std::string attribute_name, 
char &data)
 
 1734         { readAtomicAttribute(object_name,attribute_name,data); }
 
 1735     inline void readAttribute(std::string datasetName, std::string attributeName, 
signed char &data)
 
 1736         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1737     inline void readAttribute(std::string datasetName, std::string attributeName, 
signed short &data)
 
 1738         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1739     inline void readAttribute(std::string datasetName, std::string attributeName, 
signed int &data)
 
 1740         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1741     inline void readAttribute(std::string datasetName, std::string attributeName, 
signed long &data)
 
 1742         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1743     inline void readAttribute(std::string datasetName, std::string attributeName, 
signed long long &data)
 
 1744         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1745     inline void readAttribute(std::string datasetName, std::string attributeName, 
unsigned char &data)
 
 1746         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1747     inline void readAttribute(std::string datasetName, std::string attributeName, 
unsigned short &data)
 
 1748         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1749     inline void readAttribute(std::string datasetName, std::string attributeName, 
unsigned int &data)
 
 1750         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1751     inline void readAttribute(std::string datasetName, std::string attributeName, 
unsigned long &data)
 
 1752         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1753     inline void readAttribute(std::string datasetName, std::string attributeName, 
unsigned long long &data)
 
 1754         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1755     inline void readAttribute(std::string datasetName, std::string attributeName, 
float &data)
 
 1756         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1757     inline void readAttribute(std::string datasetName, std::string attributeName, 
double &data)
 
 1758         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1759     inline void readAttribute(std::string datasetName, std::string attributeName, 
long double &data)
 
 1760         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1761     inline void readAttribute(std::string datasetName, std::string attributeName, std::string &data)
 
 1762         { readAtomicAttribute(datasetName,attributeName,data); }
 
 1790     template<
unsigned int N, 
class T, 
class Str
ide>
 
 1791     inline void write(std::string datasetName,
 
 1793                       int iChunkSize = 0, 
int compression = 0)
 
 1799         for(
unsigned int i = 0; i < N; i++){
 
 1800             chunkSize[i] = iChunkSize;
 
 1802         write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize, compression);
 
 1822     template<
unsigned int N, 
class T, 
class Str
ide>
 
 1823     inline void write(std::string datasetName,
 
 1830         write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize, compression);
 
 1844     template<
unsigned int N, 
class T, 
class Str
ide>
 
 1851         typedef detail::HDF5TypeTraits<T> TypeTraits;
 
 1852         writeBlock_(datasetName, blockOffset, array,
 
 1853                     TypeTraits::getH5DataType(), TypeTraits::numberOfBands());
 
 1856     template<
unsigned int N, 
class T, 
class Str
ide>
 
 1861         typedef detail::HDF5TypeTraits<T> TypeTraits;
 
 1862         return writeBlock_(dataset, blockOffset, array,
 
 1863                            TypeTraits::getH5DataType(), TypeTraits::numberOfBands());
 
 1867     template<
unsigned int N, 
class T, 
int SIZE, 
class Str
ide>
 
 1868     inline void write(std::string datasetName,
 
 1869                       const MultiArrayView<N, TinyVector<T, SIZE>, Stride> & array,
 
 1870                       int iChunkSize = 0, 
int compression = 0)
 
 1875         typename MultiArrayShape<N>::type chunkSize;
 
 1876         for(
unsigned i = 0; i < N; i++){
 
 1877             chunkSize[i] = iChunkSize;
 
 1879         write_(datasetName, array, detail::getH5DataType<T>(), SIZE, chunkSize, compression);
 
 1882     template<
unsigned int N, 
class T, 
int SIZE, 
class Str
ide>
 
 1883     inline void write(std::string datasetName,
 
 1884                       const MultiArrayView<N, TinyVector<T, SIZE>, Stride> & array,
 
 1885                       typename MultiArrayShape<N>::type chunkSize, 
int compression = 0)
 
 1890         write_(datasetName, array, detail::getH5DataType<T>(), SIZE, chunkSize, compression);
 
 1904     void write(
const std::string & datasetName,
 
 1906                       int compression = 0)
 
 1911         write(datasetName, m_array, compression);
 
 1915     template<
unsigned int N, 
class T, 
class Str
ide>
 
 1916     inline void write(std::string datasetName,
 
 1918                       int iChunkSize = 0, 
int compression = 0)
 
 1924         for(
unsigned i = 0; i < N; i++){
 
 1925             chunkSize[i] = iChunkSize;
 
 1927         write_(datasetName, array, detail::getH5DataType<T>(), 3, chunkSize, compression);
 
 1930     template<
unsigned int N, 
class T, 
class Str
ide>
 
 1931     inline void write(std::string datasetName,
 
 1932                       const MultiArrayView<N, RGBValue<T>, Stride> & array,
 
 1933                       typename MultiArrayShape<N>::type chunkSize, 
int compression = 0)
 
 1938         write_(datasetName, array, detail::getH5DataType<T>(), 3, chunkSize, compression);
 
 1944     inline void write(std::string datasetName, 
char data) { writeAtomic(datasetName,data); }
 
 1945     inline void write(std::string datasetName, 
signed char data) { writeAtomic(datasetName,data); }
 
 1946     inline void write(std::string datasetName, 
signed short data) { writeAtomic(datasetName,data); }
 
 1947     inline void write(std::string datasetName, 
signed int data) { writeAtomic(datasetName,data); }
 
 1948     inline void write(std::string datasetName, 
signed long data) { writeAtomic(datasetName,data); }
 
 1949     inline void write(std::string datasetName, 
signed long long data) { writeAtomic(datasetName,data); }
 
 1950     inline void write(std::string datasetName, 
unsigned char data) { writeAtomic(datasetName,data); }
 
 1951     inline void write(std::string datasetName, 
unsigned short data) { writeAtomic(datasetName,data); }
 
 1952     inline void write(std::string datasetName, 
unsigned int data) { writeAtomic(datasetName,data); }
 
 1953     inline void write(std::string datasetName, 
unsigned long data) { writeAtomic(datasetName,data); }
 
 1954     inline void write(std::string datasetName, 
unsigned long long data) { writeAtomic(datasetName,data); }
 
 1955     inline void write(std::string datasetName, 
float data) { writeAtomic(datasetName,data); }
 
 1956     inline void write(std::string datasetName, 
double data) { writeAtomic(datasetName,data); }
 
 1957     inline void write(std::string datasetName, 
long double data) { writeAtomic(datasetName,data); }
 
 1958     inline void write(std::string datasetName, 
const char* data) { writeAtomic(datasetName,data); }
 
 1959     inline void write(std::string datasetName, std::string 
const & data) { writeAtomic(datasetName,data.c_str()); }
 
 1972     template<
unsigned int N, 
class T, 
class Str
ide>
 
 1978         read_(datasetName, array, detail::getH5DataType<T>(), 1);
 
 1990     template<
unsigned int N, 
class T, 
class Alloc>
 
 2001             "HDF5File::readAndResize(): Array dimension disagrees with dataset dimension.");
 
 2005         for(
int k=0; k < static_cast<int>(dimshape.
size()); ++k)
 
 2006             shape[k] = static_cast<MultiArrayIndex>(dimshape[k]);
 
 2009         read_(datasetName, array, detail::getH5DataType<T>(), 1);
 
 2022         read(datasetName, m_array);
 
 2041             "HDF5File::readAndResize(): Array dimension disagrees with Dataset dimension must equal one for vigra::ArrayVector.");
 
 2049         read_(datasetName, m_array, detail::getH5DataType<T>(), 1);
 
 2067     template<
unsigned int N, 
class T, 
class Str
ide>
 
 2075         typedef detail::HDF5TypeTraits<T> TypeTraits;
 
 2076         readBlock_(datasetName, blockOffset, blockShape, array,
 
 2077                    TypeTraits::getH5DataType(), TypeTraits::numberOfBands());
 
 2080     template<
unsigned int N, 
class T, 
class Str
ide>
 
 2086         typedef detail::HDF5TypeTraits<T> TypeTraits;
 
 2087         return readBlock_(dataset, blockOffset, blockShape, array,
 
 2088                           TypeTraits::getH5DataType(), TypeTraits::numberOfBands());
 
 2092     template<
unsigned int N, 
class T, 
int SIZE, 
class Str
ide>
 
 2093     inline void read(std::string datasetName, MultiArrayView<N, TinyVector<T, SIZE>, Stride> array)
 
 2098         read_(datasetName, array, detail::getH5DataType<T>(), SIZE);
 
 2102     template<
unsigned int N, 
class T, 
int SIZE, 
class Alloc>
 
 2103     inline void readAndResize(std::string datasetName, MultiArray<N, TinyVector<T, SIZE>, Alloc> & array)
 
 2113                            SIZE == dimshape[0], 
 
 2114             "HDF5File::readAndResize(): Array dimension disagrees with dataset dimension.");
 
 2117         typename MultiArrayShape<N>::type shape;
 
 2118         for(
int k=1; k < static_cast<int>(dimshape.size()); ++k)
 
 2119             shape[k-1] = static_cast<MultiArrayIndex>(dimshape[k]);
 
 2120         array.reshape(shape);
 
 2122         read_(datasetName, array, detail::getH5DataType<T>(), SIZE);
 
 2126     template<
unsigned int N, 
class T, 
class Str
ide>
 
 2127     inline void read(std::string datasetName, MultiArrayView<N, RGBValue<T>, Stride> array)
 
 2132         read_(datasetName, array, detail::getH5DataType<T>(), 3);
 
 2136     template<
unsigned int N, 
class T, 
class Alloc>
 
 2137     inline void readAndResize(std::string datasetName, MultiArray<N, RGBValue<T>, Alloc> & array)
 
 2148             "HDF5File::readAndResize(): Array dimension disagrees with dataset dimension.");
 
 2151         typename MultiArrayShape<N>::type shape;
 
 2152         for(
int k=1; k < static_cast<int>(dimshape.size()); ++k)
 
 2153             shape[k-1] = static_cast<MultiArrayIndex>(dimshape[k]);
 
 2154         array.reshape(shape);
 
 2156         read_(datasetName, array, detail::getH5DataType<T>(), 3);
 
 2162     inline void read(std::string datasetName, 
char &data) { readAtomic(datasetName,data); }
 
 2163     inline void read(std::string datasetName, 
signed char &data) { readAtomic(datasetName,data); }
 
 2164     inline void read(std::string datasetName, 
signed short &data) { readAtomic(datasetName,data); }
 
 2165     inline void read(std::string datasetName, 
signed int &data) { readAtomic(datasetName,data); }
 
 2166     inline void read(std::string datasetName, 
signed long &data) { readAtomic(datasetName,data); }
 
 2167     inline void read(std::string datasetName, 
signed long long &data) { readAtomic(datasetName,data); }
 
 2168     inline void read(std::string datasetName, 
unsigned char &data) { readAtomic(datasetName,data); }
 
 2169     inline void read(std::string datasetName, 
unsigned short &data) { readAtomic(datasetName,data); }
 
 2170     inline void read(std::string datasetName, 
unsigned int &data) { readAtomic(datasetName,data); }
 
 2171     inline void read(std::string datasetName, 
unsigned long &data) { readAtomic(datasetName,data); }
 
 2172     inline void read(std::string datasetName, 
unsigned long long &data) { readAtomic(datasetName,data); }
 
 2173     inline void read(std::string datasetName, 
float &data) { readAtomic(datasetName,data); }
 
 2174     inline void read(std::string datasetName, 
double &data) { readAtomic(datasetName,data); }
 
 2175     inline void read(std::string datasetName, 
long double &data) { readAtomic(datasetName,data); }
 
 2176     inline void read(std::string datasetName, std::string &data) { readAtomic(datasetName,data); }
 
 2202     template<
int N, 
class T>
 
 2205                   TinyVector<MultiArrayIndex, N> 
const & shape,
 
 2206                   typename detail::HDF5TypeTraits<T>::value_type init =
 
 2207                                                      typename detail::HDF5TypeTraits<T>::value_type(),
 
 2209                   TinyVector<MultiArrayIndex, N> 
const & chunkSize = TinyVector<MultiArrayIndex, N>(),
 
 2211                   TinyVector<MultiArrayIndex, N> 
const & chunkSize = (TinyVector<MultiArrayIndex, N>()),
 
 2213                   int compressionParameter = 0);
 
 2216     template<
int N, 
class T>
 
 2219                   TinyVector<MultiArrayIndex, N> 
const & shape,
 
 2222                   int compressionParameter = 0)
 
 2224         typename MultiArrayShape<N>::type chunkSize;
 
 2225         for(
int i = 0; i < N; i++){
 
 2226             chunkSize[i] = iChunkSize;
 
 2228         return this->
template createDataset<N, T>(datasetName, shape, init,
 
 2229                                                   chunkSize, compressionParameter);
 
 2237             H5Fflush(fileHandle_, H5F_SCOPE_GLOBAL);
 
 2250     class SplitString: 
public std::string {
 
 2252         SplitString(std::string &sstring): std::string(sstring) {};
 
 2255         std::string first(
char delimiter = 
'/')
 
 2257             size_t lastPos = find_last_of(delimiter);
 
 2258             if(lastPos == std::string::npos) 
 
 2261             return std::string(begin(), begin()+lastPos+1);
 
 2265         std::string last(
char delimiter = 
'/')
 
 2267             size_t lastPos = find_last_of(delimiter);
 
 2268             if(lastPos == std::string::npos) 
 
 2269                 return std::string(*
this);
 
 2270             return std::string(begin()+lastPos+1, end());
 
 2274     template <
class Shape>
 
 2275     ArrayVector<hsize_t>
 
 2276     defineChunks(Shape chunks, Shape 
const & shape, 
int numBands, 
int compression = 0)
 
 2278         if(
prod(chunks) > 0)
 
 2280             ArrayVector<hsize_t> res(chunks.begin(), chunks.end());
 
 2282                 res.insert(res.begin(), 
static_cast<hsize_t
>(numBands));
 
 2285         else if(compression > 0)
 
 2288             chunks = min(detail::ChunkShape<Shape::static_size>::defaultShape(), shape);
 
 2289             ArrayVector<hsize_t> res(chunks.begin(), chunks.end());
 
 2291                 res.insert(res.begin(), 
static_cast<hsize_t
>(numBands));
 
 2296             return ArrayVector<hsize_t>();
 
 2310         if(path.length() == 0 || path == 
"."){
 
 2311             return currentGroupName_();
 
 2316         if(relativePath_(path)){
 
 2317             std::string cname = currentGroupName_();
 
 2319                 str = currentGroupName_()+path;
 
 2321                 str = currentGroupName_()+
"/"+path;
 
 2327         std::string::size_type startpos = 0;
 
 2328         while(str.find(std::string(
"./"), startpos) != std::string::npos){
 
 2329             std::string::size_type pos = str.find(std::string(
"./"), startpos);
 
 2332             if(str.substr(pos-1,3) != 
"../"){
 
 2334                 str = str.substr(0,pos) + str.substr(pos+2,str.length()-pos-2);
 
 2340         while(str.find(std::string(
"..")) != std::string::npos){
 
 2341             std::string::size_type pos = str.find(std::string(
".."));
 
 2344             std::string::size_type end = str.find(
"/",pos);
 
 2345             if(end != std::string::npos){
 
 2355             std::string::size_type prev_slash = str.rfind(
"/",pos);
 
 2357             vigra_invariant(prev_slash != 0 && prev_slash != std::string::npos,
 
 2358                             "Error parsing path: "+str);
 
 2360             std::string::size_type begin = str.rfind(
"/",prev_slash-1);
 
 2363             str = str.substr(0,begin+1) + str.substr(end,str.length()-end);
 
 2373     inline bool relativePath_(std::string & path)
 const 
 2375         std::string::size_type pos = path.find(
'/') ;
 
 2384     inline std::string currentGroupName_()
 const 
 2386         int len = H5Iget_name(cGroupHandle_,NULL,1000);
 
 2387         ArrayVector<char> name (len+1,0);
 
 2388         H5Iget_name(cGroupHandle_,name.begin(),len+1);
 
 2390         return std::string(name.begin());
 
 2395     inline std::string fileName_()
 const 
 2397         int len = H5Fget_name(fileHandle_,NULL,1000);
 
 2398         ArrayVector<char> name (len+1,0);
 
 2399         H5Fget_name(fileHandle_,name.begin(),len+1);
 
 2401         return std::string(name.begin());
 
 2406     inline hid_t createFile_(std::string filePath, 
OpenMode mode = Open)
 
 2410         pFile = fopen ( filePath.c_str(), 
"r" );
 
 2414         if ( pFile == NULL )
 
 2416             vigra_precondition(mode != OpenReadOnly,
 
 2417                 "HDF5File::open(): cannot open non-existing file in read-only mode.");
 
 2418             fileId = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
 
 2423             if(mode == OpenReadOnly)
 
 2425                 fileId = H5Fopen(filePath.c_str(), H5F_ACC_RDONLY, H5P_DEFAULT);
 
 2427             else if(mode == New)
 
 2429                 std::remove(filePath.c_str());
 
 2430                 fileId = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
 
 2434                 fileId = H5Fopen(filePath.c_str(), H5F_ACC_RDWR, H5P_DEFAULT);
 
 2445     hid_t openGroup_(std::string groupName)
 const 
 2447         return const_cast<HDF5File *
>(
this)->openCreateGroup_(groupName, 
false);
 
 2458     hid_t openCreateGroup_(std::string groupName, 
bool create = 
true)
 
 2464         hid_t parent = H5Gopen(fileHandle_, 
"/", H5P_DEFAULT);
 
 2465         if(groupName == 
"/")
 
 2471         groupName = std::string(groupName.begin()+1, groupName.end());
 
 2474         if( groupName.size() != 0 && *groupName.rbegin() != 
'/')
 
 2476             groupName = groupName + 
'/';
 
 2482         HDF5DisableErrorOutput disable_error;
 
 2485         std::string::size_type begin = 0, end = groupName.find(
'/');
 
 2486         while (end != std::string::npos)
 
 2488             std::string group(groupName.begin()+begin, groupName.begin()+end);
 
 2490             hid_t prevParent = parent;
 
 2491             parent = H5Gopen(prevParent, group.c_str(), H5P_DEFAULT);
 
 2492             if(parent < 0 && create) 
 
 2493                 parent = H5Gcreate(prevParent, group.c_str(), H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
 
 2494             H5Gclose(prevParent);
 
 2500             end = groupName.find(
'/', begin);
 
 2509     inline void deleteDataset_(hid_t parent, std::string datasetName)
 
 2512         if(H5LTfind_dataset(parent, datasetName.c_str()))
 
 2515     #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6) 
 2516             if(H5Gunlink(parent, datasetName.c_str()) < 0)
 
 2518                 vigra_postcondition(
false, 
"HDF5File::deleteDataset_(): Unable to delete existing data.");
 
 2521             if(H5Ldelete(parent, datasetName.c_str(), H5P_DEFAULT ) < 0)
 
 2523                 vigra_postcondition(
false, 
"HDF5File::deleteDataset_(): Unable to delete existing data.");
 
 2531     hid_t getDatasetHandle_(std::string datasetName)
 const 
 2536         std::string groupname = SplitString(datasetName).first();
 
 2537         std::string setname = SplitString(datasetName).last();
 
 2539         if(H5Lexists(fileHandle_, datasetName.c_str(), H5P_DEFAULT) <= 0)
 
 2541             std::cerr << 
"HDF5File::getDatasetHandle_(): Dataset '" << datasetName << 
"' does not exist.\n";
 
 2546         HDF5Handle groupHandle(openGroup_(groupname), &H5Gclose, 
"HDF5File::getDatasetHandle_(): Internal error");
 
 2548         return H5Dopen(groupHandle, setname.c_str(), H5P_DEFAULT);
 
 2553     H5O_type_t get_object_type_(std::string name)
 const 
 2556         std::string group_name = SplitString(name).first();
 
 2557         std::string object_name = SplitString(name).last();
 
 2558         if (!object_name.size())
 
 2559             return H5O_TYPE_GROUP;
 
 2561         htri_t exists = H5Lexists(fileHandle_, name.c_str(), H5P_DEFAULT);
 
 2562         vigra_precondition(exists > 0,  
"HDF5File::get_object_type_(): " 
 2563                                         "object \"" + name + 
"\" " 
 2566         HDF5Handle group_handle(openGroup_(group_name), &H5Gclose, 
"Internal error");
 
 2567         return HDF5_get_type(group_handle, name.c_str());
 
 2572     template<
unsigned int N, 
class T, 
class Str
ide>
 
 2573     void write_attribute_(std::string name,
 
 2574                           const std::string & attribute_name,
 
 2575                           const MultiArrayView<N, T, Stride> & array,
 
 2576                           const hid_t datatype,
 
 2577                           const int numBandsOfType);
 
 2585     inline void writeAtomicAttribute(std::string datasetName, std::string attributeName, 
const T data)
 
 2590         typename MultiArrayShape<1>::type chunkSize;
 
 2592         MultiArray<1,T> array(MultiArrayShape<1>::type(1));
 
 2594         write_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), 1);
 
 2599     template<
unsigned int N, 
class T, 
class Str
ide>
 
 2600     void read_attribute_(std::string datasetName,
 
 2601                          std::string attributeName,
 
 2602                          MultiArrayView<N, T, Stride> array,
 
 2603                          const hid_t datatype, 
const int numBandsOfType);
 
 2611     inline void readAtomicAttribute(std::string datasetName, std::string attributeName, T & data)
 
 2616         MultiArray<1,T> array(MultiArrayShape<1>::type(1));
 
 2617         read_attribute_(datasetName, attributeName, array, detail::getH5DataType<T>(), 1);
 
 2621     inline void readAtomicAttribute(std::string datasetName, std::string attributeName, std::string & data)
 
 2626         MultiArray<1,const char *> array(MultiArrayShape<1>::type(1));
 
 2627         read_attribute_(datasetName, attributeName, array, detail::getH5DataType<const char *>(), 1);
 
 2628         data = std::string(array[0]);
 
 2633     template<
unsigned int N, 
class T, 
class Str
ide>
 
 2634     void write_(std::string &datasetName,
 
 2635                        const MultiArrayView<N, T, Stride> & array,
 
 2636                        const hid_t datatype,
 
 2637                        const int numBandsOfType,
 
 2638                        typename MultiArrayShape<N>::type &chunkSize,
 
 2639                        int compressionParameter = 0);
 
 2650     inline void writeAtomic(std::string datasetName, 
const T data)
 
 2655         typename MultiArrayShape<1>::type chunkSize;
 
 2657         MultiArray<1,T> array(MultiArrayShape<1>::type(1));
 
 2659         write_(datasetName, array, detail::getH5DataType<T>(), 1, chunkSize,0);
 
 2664     template<
unsigned int N, 
class T, 
class Str
ide>
 
 2665     void read_(std::string datasetName,
 
 2666                       MultiArrayView<N, T, Stride> array,
 
 2667                       const hid_t datatype, 
const int numBandsOfType);
 
 2678     inline void readAtomic(std::string datasetName, T & data)
 
 2683         MultiArray<1,T> array(MultiArrayShape<1>::type(1));
 
 2684         read_(datasetName, array, detail::getH5DataType<T>(), 1);
 
 2688     inline void readAtomic(std::string datasetName, std::string & data)
 
 2693         MultiArray<1,const char *> array(MultiArrayShape<1>::type(1));
 
 2694         read_(datasetName, array, detail::getH5DataType<const char *>(), 1);
 
 2695         data = std::string(array[0]);
 
 2700     template<
unsigned int N, 
class T, 
class Str
ide>
 
 2701     void writeBlock_(std::string datasetName,
 
 2702                      typename MultiArrayShape<N>::type &blockOffset,
 
 2703                      const MultiArrayView<N, T, Stride> & array,
 
 2704                      const hid_t datatype,
 
 2705                      const int numBandsOfType)
 
 2708         std::string errorMessage = 
"HDF5File::writeBlock(): Error opening dataset '" + datasetName + 
"'.";
 
 2709         HDF5HandleShared dataset(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
 
 2710         herr_t status = writeBlock_(dataset, blockOffset, array, datatype, numBandsOfType);
 
 2711         vigra_postcondition(status >= 0,
 
 2712             "HDF5File::writeBlock(): write to dataset '" + datasetName + 
"' via H5Dwrite() failed.");
 
 2719     template<
unsigned int N, 
class T, 
class Str
ide>
 
 2720     herr_t writeBlock_(HDF5HandleShared dataset,
 
 2721                        typename MultiArrayShape<N>::type &blockOffset,
 
 2722                        const MultiArrayView<N, T, Stride> & array,
 
 2723                        const hid_t datatype,
 
 2724                        const int numBandsOfType);
 
 2730     template<
unsigned int N, 
class T, 
class Str
ide>
 
 2731     void readBlock_(std::string datasetName,
 
 2732                     typename MultiArrayShape<N>::type &blockOffset,
 
 2733                     typename MultiArrayShape<N>::type &blockShape,
 
 2734                     MultiArrayView<N, T, Stride> array,
 
 2735                     const hid_t datatype, 
const int numBandsOfType)
 
 2737         std::string errorMessage (
"HDF5File::readBlock(): Unable to open dataset '" + datasetName + 
"'.");
 
 2738         HDF5HandleShared dataset(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
 
 2739         herr_t status = readBlock_(dataset, blockOffset, blockShape, array, datatype, numBandsOfType);
 
 2740         vigra_postcondition(status >= 0,
 
 2741             "HDF5File::readBlock(): read from dataset '" + datasetName + 
"' via H5Dread() failed.");
 
 2749     template<
unsigned int N, 
class T, 
class Str
ide>
 
 2750     herr_t readBlock_(HDF5HandleShared dataset,
 
 2751                       typename MultiArrayShape<N>::type &blockOffset,
 
 2752                       typename MultiArrayShape<N>::type &blockShape,
 
 2753                       MultiArrayView<N, T, Stride> array,
 
 2754                       const hid_t datatype, 
const int numBandsOfType);
 
 2759 template<
int N, 
class T>
 
 2763                         typename detail::HDF5TypeTraits<T>::value_type init,
 
 2765                          int compressionParameter)
 
 2767     vigra_precondition(!isReadOnly(),
 
 2768         "HDF5File::createDataset(): file is read-only.");
 
 2773     std::string groupname = SplitString(datasetName).first();
 
 2774     std::string setname = SplitString(datasetName).last();
 
 2776     hid_t parent = openCreateGroup_(groupname);
 
 2779     deleteDataset_(parent, setname);
 
 2783     typedef detail::HDF5TypeTraits<T> TypeTraits;
 
 2785     if(TypeTraits::numberOfBands() > 1)
 
 2787         shape_inv.resize(N+1);
 
 2788         shape_inv[N] = TypeTraits::numberOfBands();
 
 2792         shape_inv.resize(N);
 
 2794     for(
int k=0; k<N; ++k)
 
 2795         shape_inv[N-1-k] = shape[k];
 
 2799     dataspaceHandle = 
HDF5Handle(H5Screate_simple(shape_inv.
size(), shape_inv.
data(), NULL),
 
 2800                                 &H5Sclose, 
"HDF5File::createDataset(): unable to create dataspace for scalar data.");
 
 2803     HDF5Handle plist ( H5Pcreate(H5P_DATASET_CREATE), &H5Pclose, 
"HDF5File::createDataset(): unable to create property list." );
 
 2804     H5Pset_fill_value(plist, TypeTraits::getH5DataType(), &init);
 
 2807     H5Pset_obj_track_times(plist, track_time);
 
 2810     ArrayVector<hsize_t> chunks(defineChunks(chunkSize, shape, TypeTraits::numberOfBands(), compressionParameter));
 
 2811     if(chunks.
size() > 0)
 
 2813         std::reverse(chunks.
begin(), chunks.
end());
 
 2814         H5Pset_chunk (plist, chunks.
size(), chunks.
begin());
 
 2818     if(compressionParameter > 0)
 
 2820         H5Pset_deflate(plist, compressionParameter);
 
 2825                                              TypeTraits::getH5DataType(),
 
 2826                                              dataspaceHandle, H5P_DEFAULT, plist, H5P_DEFAULT),
 
 2828                                    "HDF5File::createDataset(): unable to create dataset.");
 
 2829     if(parent != cGroupHandle_)
 
 2832     return datasetHandle;
 
 2837 template<
unsigned int N, 
class T, 
class Str
ide>
 
 2838 void HDF5File::write_(std::string &datasetName,
 
 2840                       const hid_t datatype,
 
 2841                       const int numBandsOfType,
 
 2843                       int compressionParameter)
 
 2845     vigra_precondition(!isReadOnly(),
 
 2846         "HDF5File::write(): file is read-only.");
 
 2848     std::string groupname = SplitString(datasetName).first();
 
 2849     std::string setname = SplitString(datasetName).last();
 
 2853     std::reverse(shape.begin(), shape.end());
 
 2855     if(numBandsOfType > 1)
 
 2856         shape.push_back(numBandsOfType);
 
 2858     HDF5Handle dataspace(H5Screate_simple(shape.size(), shape.begin(), NULL), &H5Sclose,
 
 2859                          "HDF5File::write(): Can not create dataspace.");
 
 2862     std::string errorMessage (
"HDF5File::write(): can not create group '" + groupname + 
"'.");
 
 2863     HDF5Handle groupHandle(openCreateGroup_(groupname), &H5Gclose, errorMessage.c_str());
 
 2866     deleteDataset_(groupHandle, setname.c_str());
 
 2869     HDF5Handle plist(H5Pcreate(H5P_DATASET_CREATE), &H5Pclose,
 
 2870                      "HDF5File::write(): unable to create property list." );
 
 2873     H5Pset_obj_track_times(plist, track_time);
 
 2877     if(chunks.size() > 0)
 
 2879         std::reverse(chunks.begin(), chunks.end());
 
 2880         H5Pset_chunk (plist, chunks.size(), chunks.begin());
 
 2884     if(compressionParameter > 0)
 
 2886         H5Pset_deflate(plist, compressionParameter);
 
 2890     HDF5Handle datasetHandle(H5Dcreate(groupHandle, setname.c_str(), datatype, dataspace,H5P_DEFAULT, plist, H5P_DEFAULT),
 
 2891                              &H5Dclose, 
"HDF5File::write(): Can not create dataset.");
 
 2897         status = H5Dwrite(datasetHandle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.
data());
 
 2907         int offset = numBandsOfType > 1 ? 1 : 0;
 
 2908         std::reverse(shape.begin(), shape.end());
 
 2909         if(chunks.size() > 0)
 
 2912             std::reverse(chunks.begin(), chunks.end());
 
 2917             ArrayVector<hsize_t>(shape.size(), 1).swap(chunks);
 
 2918             chunks[0] = numBandsOfType;
 
 2920             for(
unsigned int k=0; k<N;  ++k)
 
 2922                 chunks[k+offset] = array.
shape(k);
 
 2923                 prod *= array.
shape(k);
 
 2929         ArrayVector<hsize_t> null(shape.size(), 0),
 
 2930                              start(shape.size(), 0),
 
 2931                              count(shape.size(), 1);
 
 2933         count[N-1-offset] = numBandsOfType;
 
 2935         typedef typename MultiArrayShape<N>::type Shape;
 
 2936         Shape chunkCount, chunkMaxShape;
 
 2937         for(
unsigned int k=offset; k<chunks.size(); ++k)
 
 2939             chunkMaxShape[k-offset] = chunks[k];
 
 2943         typename CoupledIteratorType<N>::type chunkIter = createCoupledIterator(chunkCount),
 
 2944                                               chunkEnd  = chunkIter.getEndIterator();
 
 2945         for(; chunkIter != chunkEnd; ++chunkIter)
 
 2947             Shape chunkStart(chunkIter.point() * chunkMaxShape),
 
 2948                   chunkStop(min(chunkStart + chunkMaxShape, array.
shape()));
 
 2949             MultiArray<N, T> buffer(array.
subarray(chunkStart, chunkStop));
 
 2951             for(
unsigned int k=0; k<N; ++k)
 
 2953                 start[N-1-k] = chunkStart[k];
 
 2954                 count[N-1-k] = buffer.shape(k);
 
 2959                 count[N] = numBandsOfType;
 
 2961             HDF5Handle filespace(H5Dget_space(datasetHandle),
 
 2962                                  &H5Sclose, 
"HDF5File::write(): unable to create hyperslabs.");
 
 2963             status = H5Sselect_hyperslab(filespace, H5S_SELECT_SET, start.data(), NULL, count.data(), NULL);
 
 2967             HDF5Handle dataspace2(H5Screate_simple(count.size(), count.data(), NULL),
 
 2968                                  &H5Sclose, 
"HDF5File::write(): unable to create hyperslabs.");
 
 2969             status = H5Sselect_hyperslab(dataspace2, H5S_SELECT_SET, null.data(), NULL, count.data(), NULL);
 
 2973             status = H5Dwrite(datasetHandle, datatype, dataspace2, filespace, H5P_DEFAULT, buffer.data());
 
 2978     vigra_postcondition(status >= 0,
 
 2979         "HDF5File::write(): write to dataset '" + datasetName + 
"' via H5Dwrite() failed.");
 
 2984 template<
unsigned int N, 
class T, 
class Str
ide>
 
 2985 herr_t HDF5File::writeBlock_(HDF5HandleShared datasetHandle,
 
 2986                              typename MultiArrayShape<N>::type &blockOffset,
 
 2987                              const MultiArrayView<N, T, Stride> & array,
 
 2988                              const hid_t datatype,
 
 2989                              const int numBandsOfType)
 
 2991     vigra_precondition(!isReadOnly(),
 
 2992         "HDF5File::writeBlock(): file is read-only.");
 
 2994     ArrayVector<hsize_t> boffset, bshape, bones(N+1, 1);
 
 2995     hssize_t dimensions = getDatasetDimensions_(datasetHandle);
 
 2996     if(numBandsOfType > 1)
 
 2998         vigra_precondition(N+1 == dimensions,
 
 2999             "HDF5File::readBlock(): Array dimension disagrees with data dimension.");
 
 3001         boffset.resize(N+1);
 
 3002         bshape[N] = numBandsOfType;
 
 3007         vigra_precondition(N == dimensions,
 
 3008             "HDF5File::readBlock(): Array dimension disagrees with data dimension.");
 
 3013     for(
unsigned i = 0; i < N; ++i)
 
 3016         bshape[N-1-i] = array.shape(i);
 
 3017         boffset[N-1-i] = blockOffset[i];
 
 3021     HDF5Handle memspace_handle (H5Screate_simple(bshape.size(), bshape.data(), NULL),
 
 3023                                 "Unable to get origin dataspace");
 
 3026     HDF5Handle dataspaceHandle (H5Dget_space(datasetHandle),&H5Sclose,
"Unable to create target dataspace");
 
 3027     H5Sselect_hyperslab(dataspaceHandle, H5S_SELECT_SET,
 
 3028                         boffset.data(), bones.data(), bones.data(), bshape.data());
 
 3031     if(array.isUnstrided())
 
 3034         status = H5Dwrite( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, array.data());
 
 3039         MultiArray<N, T> buffer(array);
 
 3040         status = H5Dwrite( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, buffer.data());
 
 3047 template<
unsigned int N, 
class T, 
class Str
ide>
 
 3048 void HDF5File::write_attribute_(std::string name,
 
 3049                                 const std::string & attribute_name,
 
 3050                                 const MultiArrayView<N, T, Stride> & array,
 
 3051                                 const hid_t datatype,
 
 3052                                 const int numBandsOfType)
 
 3054     vigra_precondition(!isReadOnly(),
 
 3055         "HDF5File::writeAttribute(): file is read-only.");
 
 3058     ArrayVector<hsize_t> shape(array.shape().begin(), array.shape().end());
 
 3059     std::reverse(shape.begin(), shape.end());
 
 3060     if(numBandsOfType > 1)
 
 3061         shape.push_back(numBandsOfType);
 
 3063     HDF5Handle dataspace(H5Screate_simple(shape.size(),
 
 3064                                           shape.begin(), NULL),
 
 3065                          &H5Sclose, 
"HDF5File::writeAttribute(): Can not" 
 3066                                     " create dataspace.");
 
 3068     std::string errorMessage (
"HDF5File::writeAttribute(): can not find " 
 3069                               "object '" + name + 
"'.");
 
 3071     H5O_type_t h5_type = get_object_type_(name);
 
 3072     bool is_group = h5_type == H5O_TYPE_GROUP;
 
 3073     if (!is_group && h5_type != H5O_TYPE_DATASET)
 
 3074         vigra_precondition(0, 
"HDF5File::writeAttribute(): object \"" 
 3075                                + name + 
"\" is neither a group nor a " 
 3078     HDF5Handle object_handle(is_group
 
 3079                                  ? openCreateGroup_(name)
 
 3080                                  : getDatasetHandle_(name),
 
 3084                              errorMessage.c_str());
 
 3087     HDF5Handle attributeHandle(exists
 
 3088                                ? H5Aopen(object_handle,
 
 3089                                          attribute_name.c_str(),
 
 3091                                : H5Acreate(object_handle,
 
 3092                                            attribute_name.c_str(), datatype,
 
 3093                                            dataspace, H5P_DEFAULT,
 
 3096                                "HDF5File::writeAttribute(): Can not create" 
 3099     if(array.isUnstrided())
 
 3102         status = H5Awrite(attributeHandle, datatype, array.data());
 
 3108         MultiArray<N, T> buffer(array);
 
 3109         status = H5Awrite(attributeHandle, datatype, buffer.data());
 
 3111     vigra_postcondition(status >= 0,
 
 3112         "HDF5File::writeAttribute(): write to attribute '" + attribute_name + 
"' via H5Awrite() failed.");
 
 3117 template<
unsigned int N, 
class T, 
class Str
ide>
 
 3118 void HDF5File::read_(std::string datasetName,
 
 3119                      MultiArrayView<N, T, Stride> array,
 
 3120                      const hid_t datatype, 
const int numBandsOfType)
 
 3125     std::string errorMessage (
"HDF5File::read(): Unable to open dataset '" + datasetName + 
"'.");
 
 3126     HDF5Handle datasetHandle(getDatasetHandle_(datasetName), &H5Dclose, errorMessage.c_str());
 
 3130     int offset = (numBandsOfType > 1)
 
 3135         "HDF5File::read(): Array dimension disagrees with dataset dimension.");
 
 3137     typename MultiArrayShape<N>::type shape;
 
 3138     for(
int k=offset; k < (int)dimshape.size(); ++k)
 
 3141     vigra_precondition(shape == array.shape(),
 
 3142                        "HDF5File::read(): Array shape disagrees with dataset shape.");
 
 3144         vigra_precondition(dimshape[0] == static_cast<hsize_t>(numBandsOfType),
 
 3145                            "HDF5File::read(): Band count doesn't match destination array compound type.");
 
 3148     if(array.isUnstrided())
 
 3151         status = H5Dread(datasetHandle, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, array.data() );
 
 3157         ArrayVector<hsize_t> null(dimshape.size(), 0),
 
 3158                              chunks(dimshape.size(), 1),
 
 3159                              start(dimshape.size(), 0),
 
 3160                              count(dimshape.size(), 1);
 
 3162         HDF5Handle properties(H5Dget_create_plist(datasetHandle),
 
 3163                              &H5Pclose, 
"HDF5File::read(): failed to get property list");
 
 3164         if(H5D_CHUNKED == H5Pget_layout(properties))
 
 3167             H5Pget_chunk(properties, static_cast<int>(chunks.size()), chunks.data());
 
 3168             std::reverse(chunks.begin(), chunks.end());
 
 3173             chunks[0] = numBandsOfType;
 
 3175             for(
unsigned int k=0; k<N; ++k)
 
 3177                 chunks[k+offset] = array.shape(k);
 
 3178                 prod *= array.shape(k);
 
 3184         count[N-1-offset] = 
static_cast<hsize_t
>(numBandsOfType);
 
 3186         typedef typename MultiArrayShape<N>::type Shape;
 
 3187         Shape chunkCount, chunkMaxShape;
 
 3188         for(
unsigned int k=offset; k<chunks.size(); ++k)
 
 3190             chunkMaxShape[k-offset] = chunks[k];
 
 3194         typename CoupledIteratorType<N>::type chunkIter = createCoupledIterator(chunkCount),
 
 3195                                               chunkEnd  = chunkIter.getEndIterator();
 
 3196         for(; chunkIter != chunkEnd; ++chunkIter)
 
 3198             Shape chunkStart(chunkIter.point() * chunkMaxShape),
 
 3199                   chunkStop(min(chunkStart + chunkMaxShape, array.shape()));
 
 3200             MultiArray<N, T> buffer(chunkStop - chunkStart);
 
 3202             for(
unsigned int k=0; k<N; ++k)
 
 3204                 start[N-1-k] = chunkStart[k];
 
 3205                 count[N-1-k] = buffer.shape(k);
 
 3210                 count[N] = numBandsOfType;
 
 3212             HDF5Handle filespace(H5Dget_space(datasetHandle),
 
 3213                                  &H5Sclose, 
"HDF5File::read(): unable to create hyperslabs.");
 
 3214             status = H5Sselect_hyperslab(filespace, H5S_SELECT_SET, start.data(), NULL, count.data(), NULL);
 
 3218             HDF5Handle dataspace(H5Screate_simple(count.size(), count.data(), NULL),
 
 3219                                  &H5Sclose, 
"HDF5File::read(): unable to create hyperslabs.");
 
 3220             status = H5Sselect_hyperslab(dataspace, H5S_SELECT_SET, null.data(), NULL, count.data(), NULL);
 
 3224             status = H5Dread(datasetHandle, datatype, dataspace, filespace, H5P_DEFAULT, buffer.data());
 
 3228             array.subarray(chunkStart, chunkStop) = buffer;
 
 3231     vigra_postcondition(status >= 0,
 
 3232         "HDF5File::read(): read from dataset '" + datasetName + 
"' via H5Dread() failed.");
 
 3237 template<
unsigned int N, 
class T, 
class Str
ide>
 
 3238 herr_t HDF5File::readBlock_(HDF5HandleShared datasetHandle,
 
 3239                             typename MultiArrayShape<N>::type &blockOffset,
 
 3240                             typename MultiArrayShape<N>::type &blockShape,
 
 3241                             MultiArrayView<N, T, Stride> array,
 
 3242                             const hid_t datatype, 
const int numBandsOfType)
 
 3244     vigra_precondition(blockShape == array.shape(),
 
 3245          "HDF5File::readBlock(): Array shape disagrees with block size.");
 
 3247     ArrayVector<hsize_t> boffset, bshape, bones(N+1, 1);
 
 3248     hssize_t dimensions = getDatasetDimensions_(datasetHandle);
 
 3249     if(numBandsOfType > 1)
 
 3251         vigra_precondition(N+1 == dimensions,
 
 3252             "HDF5File::readBlock(): Array dimension disagrees with data dimension.");
 
 3254         boffset.resize(N+1);
 
 3255         bshape[N] = numBandsOfType;
 
 3260         vigra_precondition(N == dimensions,
 
 3261             "HDF5File::readBlock(): Array dimension disagrees with data dimension.");
 
 3266     for(
unsigned i = 0; i < N; ++i)
 
 3269         bshape[N-1-i] = blockShape[i];
 
 3270         boffset[N-1-i] = blockOffset[i];
 
 3274     HDF5Handle memspace_handle(H5Screate_simple(bshape.size(), bshape.data(), NULL),
 
 3276                                "Unable to create target dataspace");
 
 3279     HDF5Handle dataspaceHandle(H5Dget_space(datasetHandle), &H5Sclose,
 
 3280                                "Unable to get dataspace");
 
 3281     H5Sselect_hyperslab(dataspaceHandle, H5S_SELECT_SET,
 
 3282                         boffset.data(), bones.data(), bones.data(), bshape.data());
 
 3285     if(array.isUnstrided())
 
 3288         status = H5Dread( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, array.data());
 
 3293         MultiArray<N, T> buffer(array.shape());
 
 3294         status = H5Dread( datasetHandle, datatype, memspace_handle, dataspaceHandle, H5P_DEFAULT, buffer.data());
 
 3304 template<
unsigned int N, 
class T, 
class Str
ide>
 
 3305 void HDF5File::read_attribute_(std::string datasetName,
 
 3306                                std::string attributeName,
 
 3307                                MultiArrayView<N, T, Stride> array,
 
 3308                                const hid_t datatype, 
const int numBandsOfType)
 
 3312     std::string message = 
"HDF5File::readAttribute(): could not get handle for attribute '"+attributeName+
"'' of object '"+dataset_path+
"'.";
 
 3313     HDF5Handle attr_handle (H5Aopen_by_name(fileHandle_,dataset_path.c_str(),attributeName.c_str(),H5P_DEFAULT,H5P_DEFAULT),&H5Aclose, message.c_str());
 
 3316     message = 
"HDF5File::readAttribute(): could not get dataspace for attribute '"+attributeName+
"'' of object '"+dataset_path+
"'.";
 
 3317     HDF5Handle attr_dataspace_handle (H5Aget_space(attr_handle),&H5Sclose,message.c_str());
 
 3320     int raw_dims = H5Sget_simple_extent_ndims(attr_dataspace_handle);
 
 3321     int dims = std::max(raw_dims, 1); 
 
 3322     ArrayVector<hsize_t> dimshape(dims);
 
 3324         H5Sget_simple_extent_dims(attr_dataspace_handle, dimshape.data(), NULL);
 
 3329     std::reverse(dimshape.begin(), dimshape.end());
 
 3331     int offset = (numBandsOfType > 1)
 
 3334     message = 
"HDF5File::readAttribute(): Array dimension disagrees with dataset dimension.";
 
 3338     for(
int k=offset; k < (int)dimshape.size(); ++k)
 
 3339         vigra_precondition(array.shape()[k-offset] == (
MultiArrayIndex)dimshape[k],
 
 3340                            "HDF5File::readAttribute(): Array shape disagrees with dataset shape");
 
 3343     if(array.isUnstrided())
 
 3346         status = H5Aread( attr_handle, datatype, array.data());
 
 3352         MultiArray<N, T> buffer(array.shape());
 
 3353         status = H5Aread( attr_handle, datatype, buffer.data() );
 
 3358     vigra_postcondition(status >= 0,
 
 3359         "HDF5File::readAttribute(): read from attribute '" + attributeName + 
"' via H5Aread() failed.");
 
 3401 template<
unsigned int N, 
class T, 
class Str
ideTag>
 
 3402 inline void readHDF5(
const HDF5ImportInfo &info, MultiArrayView<N, T, StrideTag> array)
 
 3404     HDF5File file(info.getFilePath(), HDF5File::OpenReadOnly);
 
 3405     file.read(info.getPathInFile(), array);
 
 3408 inline hid_t openGroup(hid_t parent, std::string group_name)
 
 3411     size_t last_slash = group_name.find_last_of(
'/');
 
 3412     if (last_slash == std::string::npos || last_slash != group_name.size() - 1)
 
 3413         group_name = group_name + 
'/';
 
 3414     std::string::size_type begin = 0, end = group_name.find(
'/');
 
 3416     while (end != std::string::npos)
 
 3418         std::string group(group_name.begin()+begin, group_name.begin()+end);
 
 3419         hid_t prev_parent = parent;
 
 3420         parent = H5Gopen(prev_parent, group.c_str(), H5P_DEFAULT);
 
 3422         if(ii != 0)     H5Gclose(prev_parent);
 
 3423         if(parent < 0)  
return parent;
 
 3426         end = group_name.find(
'/', begin);
 
 3431 inline hid_t createGroup(hid_t parent, std::string group_name)
 
 3433     if(group_name.size() == 0 ||*group_name.rbegin() != 
'/')
 
 3434         group_name = group_name + 
'/';
 
 3435     if(group_name == 
"/")
 
 3436         return H5Gopen(parent, group_name.c_str(), H5P_DEFAULT);
 
 3438     std::string::size_type begin = 0, end = group_name.find(
'/');
 
 3440     while (end != std::string::npos)
 
 3442         std::string group(group_name.begin()+begin, group_name.begin()+end);
 
 3443         hid_t prev_parent = parent;
 
 3445         if(H5LTfind_dataset(parent, group.c_str()) == 0)
 
 3447             parent = H5Gcreate(prev_parent, group.c_str(), H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
 
 3449             parent = H5Gopen(prev_parent, group.c_str(), H5P_DEFAULT);
 
 3452         if(ii != 0)     H5Gclose(prev_parent);
 
 3453         if(parent < 0)  
return parent;
 
 3456         end = group_name.find(
'/', begin);
 
 3461 inline void deleteDataset(hid_t parent, std::string dataset_name)
 
 3464     if(H5LTfind_dataset(parent, dataset_name.c_str()))
 
 3467 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR <= 6) 
 3468         if(H5Gunlink(parent, dataset_name.c_str()) < 0)
 
 3470             vigra_postcondition(
false, 
"writeToHDF5File(): Unable to delete existing data.");
 
 3473         if(H5Ldelete(parent, dataset_name.c_str(), H5P_DEFAULT ) < 0)
 
 3475             vigra_postcondition(
false, 
"createDataset(): Unable to delete existing data.");
 
 3481 inline hid_t createFile(std::string filePath, 
bool append_ = 
true)
 
 3484     pFile = fopen ( filePath.c_str(), 
"r" );
 
 3486     if ( pFile == NULL )
 
 3488         file_id = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
 
 3493         file_id = H5Fopen(filePath.c_str(), H5F_ACC_RDWR, H5P_DEFAULT);
 
 3498         std::remove(filePath.c_str());
 
 3499         file_id = H5Fcreate(filePath.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
 
 3504 template<
unsigned int N, 
class T, 
class Tag>
 
 3505 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)
 
 3507     std::string path_name(pathInFile), group_name, data_set_name, message;
 
 3508     std::string::size_type delimiter = path_name.rfind(
'/');
 
 3511     file_handle = HDF5Handle(createFile(filePath), &H5Fclose,
 
 3512                        "createDataset(): unable to open output file.");
 
 3515     if(delimiter == std::string::npos)
 
 3518         data_set_name = path_name;
 
 3522         group_name = std::string(path_name.begin(), path_name.begin()+delimiter);
 
 3523         data_set_name = std::string(path_name.begin()+delimiter+1, path_name.end());
 
 3527     HDF5Handle group(createGroup(file_handle, group_name), &H5Gclose,
 
 3528                      "createDataset(): Unable to create and open group. generic v");
 
 3531     deleteDataset(group, data_set_name);
 
 3535     HDF5Handle dataspace_handle;
 
 3536     if(numBandsOfType > 1) {
 
 3538         hsize_t shape_inv[N+1]; 
 
 3539         for(
unsigned int k=0; k<N; ++k) {
 
 3540             shape_inv[N-1-k] = array.shape(k);  
 
 3543         shape_inv[N] = numBandsOfType;
 
 3546         dataspace_handle = HDF5Handle(H5Screate_simple(N+1, shape_inv, NULL),
 
 3547                                     &H5Sclose, 
"createDataset(): unable to create dataspace for non-scalar data.");
 
 3550         hsize_t shape_inv[N];
 
 3551         for(
unsigned int k=0; k<N; ++k)
 
 3552             shape_inv[N-1-k] = array.shape(k);
 
 3555         dataspace_handle = HDF5Handle(H5Screate_simple(N, shape_inv, NULL),
 
 3556                                     &H5Sclose, 
"createDataset(): unable to create dataspace for scalar data.");
 
 3560     dataset_handle = HDF5Handle(H5Dcreate(group,
 
 3561                                         data_set_name.c_str(),
 
 3564                                         H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT),
 
 3565                               &H5Dclose, 
"createDataset(): unable to create dataset.");
 
 3608 template<
unsigned int N, 
class T, 
class Str
ideTag>
 
 3609 inline void writeHDF5(
const char* filePath, 
const char* pathInFile, 
const MultiArrayView<N, T, StrideTag> & array)
 
 3611     HDF5File file(filePath, HDF5File::Open);
 
 3612     file.write(pathInFile, array);
 
 3625     void operator()(std::string 
const & in)
 
 3627         size = in.size() > size ?
 
 3635 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR == 8) || DOXYGEN 
 3642 template<
size_t N, 
class T, 
class C>
 
 3647     if(H5Aexists(loc, name) > 0)
 
 3648         H5Adelete(loc, name);
 
 3651                                array.
shape().end());
 
 3653         dataspace_handle(H5Screate_simple(N, shape.data(), NULL),
 
 3655                          "writeToHDF5File(): unable to create dataspace.");
 
 3659                               detail::getH5DataType<T>(),
 
 3661                               H5P_DEFAULT ,H5P_DEFAULT ),
 
 3663                     "writeHDF5Attr: unable to create Attribute");
 
 3667     for(
int ii = 0; ii < array.
size(); ++ii)
 
 3668         buffer.push_back(array[ii]);
 
 3669     H5Awrite(attr, detail::getH5DataType<T>(), buffer.
data());
 
 3678 template<
size_t N, 
class C>
 
 3683     if(H5Aexists(loc, name) > 0)
 
 3684         H5Adelete(loc, name);
 
 3687                                array.
shape().end());
 
 3689         dataspace_handle(H5Screate_simple(N, shape.data(), NULL),
 
 3691                          "writeToHDF5File(): unable to create dataspace.");
 
 3695                      "writeToHDF5File(): unable to create type.");
 
 3697     detail::MaxSizeFnc max_size;
 
 3698     max_size = std::for_each(array.
data(),array.
data()+ array.
size(), max_size);
 
 3699     H5Tset_size (atype, max_size.size);
 
 3705                               H5P_DEFAULT ,H5P_DEFAULT ),
 
 3707                     "writeHDF5Attr: unable to create Attribute");
 
 3709     std::string buf =
"";
 
 3710     for(
int ii = 0; ii < array.
size(); ++ii)
 
 3712         buf = buf + array[ii]
 
 3713                   + std::string(max_size.size - array[ii].
size(), 
' ');
 
 3715     H5Awrite(attr, atype, buf.c_str());
 
 3743                             std::string pathInFile,
 
 3746     std::string path_name(pathInFile), group_name, data_set_name, message, attr_name;
 
 3747     std::string::size_type delimiter = path_name.rfind(
'/');
 
 3750     HDF5Handle file_id(createFile(filePath), &H5Fclose,
 
 3751                        "writeToHDF5File(): unable to open output file.");
 
 3754     if(delimiter == std::string::npos)
 
 3757         data_set_name = path_name;
 
 3762         group_name = std::string(path_name.begin(), path_name.begin()+delimiter);
 
 3763         data_set_name = std::string(path_name.begin()+delimiter+1, path_name.end());
 
 3765     delimiter = data_set_name.rfind(
'.');
 
 3766     if(delimiter == std::string::npos)
 
 3768         attr_name = path_name;
 
 3769         data_set_name = 
"/";
 
 3773         attr_name = std::string(data_set_name.begin()+delimiter+1, data_set_name.end());
 
 3774         data_set_name = std::string(data_set_name.begin(), data_set_name.begin()+delimiter);
 
 3777     HDF5Handle group(openGroup(file_id, group_name), &H5Gclose,
 
 3778                      "writeToHDF5File(): Unable to create and open group. attr ver");
 
 3780     if(data_set_name != 
"/")
 
 3782         HDF5Handle dset(H5Dopen(group, data_set_name.c_str(), H5P_DEFAULT), &H5Dclose,
 
 3783                         "writeHDF5Attr():unable to open dataset");
 
 3798 #endif // VIGRA_HDF5IMPEX_HXX 
HDF5File & operator=(HDF5File const &other)
Assign a HDF5File object. 
Definition: hdf5impex.hxx:1150
void read(const std::string &datasetName, ArrayVectorView< T > array)
Read data into an array vector. If the first character of datasetName is a "/", the path will be inte...
Definition: hdf5impex.hxx:2017
HDF5Handle()
Default constructor. Creates a NULL handle. 
Definition: hdf5impex.hxx:224
void write(const std::string &datasetName, const ArrayVectorView< T > &array, int compression=0)
Write array vectors. 
Definition: hdf5impex.hxx:1904
bool cd_up(int levels)
Change the current group to its parent group. Returns true if successful, false otherwise. If unsuccessful, the group will not change. 
Definition: hdf5impex.hxx:1247
std::vector< std::string > listAttributes(std::string const &group_or_dataset) const 
List the attribute names of the given group or dataset. 
Definition: hdf5impex.hxx:1575
bool operator==(HDF5Handle const &h) const 
Equality comparison of the contained handle. 
Definition: hdf5impex.hxx:365
Temporarily disable HDF5's native error output. 
Definition: hdf5impex.hxx:130
ArrayVector< hsize_t > getDatasetShape(std::string datasetName) const 
Get the shape of each dimension of a certain dataset. 
Definition: hdf5impex.hxx:1394
Wrapper for unique hid_t objects. 
Definition: hdf5impex.hxx:210
MultiArrayIndex numDimensions() const 
bool operator!=(hid_t h) const 
Inequality comparison of the contained handle. 
Definition: hdf5impex.hxx:636
HDF5HandleShared()
Default constructor. Creates a NULL handle. 
Definition: hdf5impex.hxx:446
const difference_type & shape() const 
Definition: multi_array.hxx:1648
void readAttribute(std::string object_name, std::string attribute_name, MultiArrayView< N, T, Stride > array)
Read MultiArray Attributes. In contrast to datasets, subarray access is not available. 
Definition: hdf5impex.hxx:1698
void cd_mk(std::string groupName)
Change the current group; create it if necessary. If the first character is a "/", the path will be interpreted as absolute path, otherwise it will be interpreted as path relative to the current group. 
Definition: hdf5impex.hxx:1285
hssize_t getDatasetDimensions(std::string datasetName) const 
Get the number of dimensions of a certain dataset If the first character is a "/", the path will be interpreted as absolute path, otherwise it will be interpreted as path relative to the current group. 
Definition: hdf5impex.hxx:1365
Definition: array_vector.hxx:76
const_iterator begin() const 
Definition: array_vector.hxx:223
pointer data() const 
Definition: multi_array.hxx:1898
HDF5ImportInfo(const char *filePath, const char *pathInFile)
HDF5Handle(hid_t h, Destructor destructor, const char *error_message)
Create a wrapper for a hid_t object. 
Definition: hdf5impex.hxx:248
std::string get_absolute_path(std::string path) const 
takes any path and converts it into an absolute path in the current file. 
Definition: hdf5impex.hxx:2308
void swap(HDF5HandleShared &h)
Swap the contents of two handle wrappers. 
Definition: hdf5impex.hxx:586
void reshape(const difference_type &shape)
Definition: multi_array.hxx:2861
HDF5HandleShared getDatasetHandleShared(std::string const &datasetName) const 
Obtain a shared HDF5 handle of a dataset. 
Definition: hdf5impex.hxx:1527
ArrayVector< hsize_t > getChunkShape(std::string datasetName) const 
Get the shape of chunks along each dimension of a certain dataset. 
Definition: hdf5impex.hxx:1430
void mkdir(std::string groupName)
Create a new group. If the first character is a "/", the path will be interpreted as absolute path...
Definition: hdf5impex.hxx:1268
bool existsAttribute(std::string object_name, std::string attribute_name)
Test if attribute exists. 
Definition: hdf5impex.hxx:1681
const std::string & getFilePath() const 
std::ptrdiff_t MultiArrayIndex
Definition: multi_fwd.hxx:60
void writeHDF5Attr(hid_t loc, const char *name, MultiArrayView< N, T, C > const &array)
Definition: hdf5impex.hxx:3643
HDF5HandleShared createDataset(std::string datasetName, TinyVector< MultiArrayIndex, N > const &shape, typename detail::HDF5TypeTraits< T >::value_type init=typename detail::HDF5TypeTraits< T >::value_type(), TinyVector< MultiArrayIndex, N > const &chunkSize=(TinyVector< MultiArrayIndex, N >()), int compressionParameter=0)
Create a new dataset. This function can be used to create a dataset filled with a default value init...
Definition: hdf5impex.hxx:2761
bool operator==(hid_t h) const 
Equality comparison of the contained handle. 
Definition: hdf5impex.hxx:372
hid_t release()
Return the contained handle and set the wrapper to NULL without calling close(). 
Definition: hdf5impex.hxx:309
const std::string & getPathInFile() const 
herr_t close()
Close the handle if this is the unique (i.e. last) owner. 
Definition: hdf5impex.hxx:527
void ls(Container &cont) const 
List the contents of the current group into a container-like object via insert(). ...
Definition: hdf5impex.hxx:1332
void write(std::string datasetName, const MultiArrayView< N, T, Stride > &array, typename MultiArrayShape< N >::type chunkSize, int compression=0)
Write multi arrays. Chunks can be activated by providing a MultiArrayShape as chunkSize. chunkSize must have equal dimension as array. 
Definition: hdf5impex.hxx:1823
Argument object for the function readHDF5(). 
Definition: hdf5impex.hxx:682
bool isHDF5(char const *filename)
Check if given filename refers to a HDF5 file. 
Definition: hdf5impex.hxx:102
void reset(hid_t h, Destructor destructor, const char *error_message)
Reset the handle to a new value. 
Definition: hdf5impex.hxx:550
HDF5File(HDF5File const &other)
Copy a HDF5File object. 
Definition: hdf5impex.hxx:1126
HDF5HandleShared & operator=(HDF5HandleShared const &h)
Assignment. Call close() for the present LHS handle and share ownership with the RHS handle (analogou...
Definition: hdf5impex.hxx:498
HDF5HandleShared(HDF5HandleShared const &h)
Copy constructor. Shares ownership with the RHS handle (analogous to std::shared_ptr). 
Definition: hdf5impex.hxx:485
bool operator!=(hid_t h) const 
Inequality comparison of the contained handle. 
Definition: hdf5impex.hxx:386
bool operator!=(HDF5Handle const &h) const 
Inequality comparison of the contained handle. 
Definition: hdf5impex.hxx:379
difference_type_1 size() const 
Definition: multi_array.hxx:1641
Wrapper for shared hid_t objects. 
Definition: hdf5impex.hxx:431
std::string filename() const 
Get the name of the associated file. 
Definition: hdf5impex.hxx:1347
HDF5File()
Default constructor. 
Definition: hdf5impex.hxx:1047
HDF5Handle & operator=(HDF5Handle const &h)
Assignment. Calls close() for the LHS handle and hands over ownership of the RHS handle (analogous to...
Definition: hdf5impex.hxx:271
~HDF5File()
The destructor flushes and closes the file. 
Definition: hdf5impex.hxx:1137
NumericTraits< V >::Promote prod(TinyVectorBase< V, SIZE, D1, D2 > const &l)
product of the vector's elements 
Definition: tinyvector.hxx:2097
bool operator==(HDF5HandleShared const &h) const 
Equality comparison of the contained handle. 
Definition: hdf5impex.hxx:615
HDF5File(HDF5HandleShared const &fileHandle, const std::string &pathname="", bool read_only=false)
Initialize an HDF5File object from HDF5 file handle. 
Definition: hdf5impex.hxx:1098
HDF5Handle getAttributeHandle(std::string dataset_name, std::string attribute_name) const 
Obtain the HDF5 handle of a attribute. 
Definition: hdf5impex.hxx:1598
void readAndResize(std::string datasetName, MultiArray< N, T, Alloc > &array)
Read data into a MultiArray. Resize MultiArray to the correct size. If the first character of dataset...
Definition: hdf5impex.hxx:1991
bool operator==(hid_t h) const 
Equality comparison of the contained handle. 
Definition: hdf5impex.hxx:622
HDF5File(char const *filePath, OpenMode mode=ReadOnly, bool track_creation_times=false)
Open or create an HDF5File object. 
Definition: hdf5impex.hxx:1076
void readHDF5(...)
Read the data specified by the given vigra::HDF5ImportInfo object and write the into the given 'array...
void readAndResize(std::string datasetName, ArrayVector< T > &array)
Read data into an array vector. Resize the array vector to the correct size. If the first character o...
Definition: hdf5impex.hxx:2030
~HDF5Handle()
Destructor. Calls close() for the contained handle. 
Definition: hdf5impex.hxx:286
bool operator!=(HDF5HandleShared const &h) const 
Inequality comparison of the contained handle. 
Definition: hdf5impex.hxx:629
void read(std::string datasetName, char &data)
Read a single value. Specialization of the read function for simple datatypes. 
Definition: hdf5impex.hxx:2162
doxygen_overloaded_function(template<...> void separableConvolveBlockwise) template< unsigned int N
Separated convolution on ChunkedArrays. 
void open(std::string filePath, OpenMode mode)
Open or create the given file in the given mode and set the group to "/". If another file is currentl...
Definition: hdf5impex.hxx:1187
PixelType pixelType() const 
OpenMode
Set how a file is opened. 
Definition: hdf5impex.hxx:1031
void readAttribute(std::string object_name, std::string attribute_name, char &data)
Read a single value. Specialization of the read function for simple datatypes. 
Definition: hdf5impex.hxx:1733
void writeBlock(std::string datasetName, typename MultiArrayShape< N >::type blockOffset, const MultiArrayView< N, T, Stride > &array)
Write a multi array into a larger volume. blockOffset determines the position, where array is written...
Definition: hdf5impex.hxx:1845
void flushToDisk()
Immediately write all data to disk. 
Definition: hdf5impex.hxx:2234
HDF5Handle(HDF5Handle const &h)
Copy constructor. 
Definition: hdf5impex.hxx:260
void reset(hid_t h, Destructor destructor, const char *error_message)
Reset the wrapper to a new handle. 
Definition: hdf5impex.hxx:321
Class for fixed size vectors.This class contains an array of size SIZE of the specified VALUETYPE...
Definition: accessor.hxx:940
void writeHDF5(...)
Store array data in an HDF5 file. 
MultiArrayIndex shapeOfDimension(const int dim) const 
HDF5Handle getGroupHandle(std::string group_name, std::string function_name="HDF5File::getGroupHandle()")
Obtain the HDF5 handle of a group (create the group if it doesn't exist). 
Definition: hdf5impex.hxx:1535
std::string pwd() const 
Get the path of the current group. 
Definition: hdf5impex.hxx:1340
void read(std::string datasetName, MultiArrayView< N, T, Stride > array)
Read data into a multi array. If the first character of datasetName is a "/", the path will be interp...
Definition: hdf5impex.hxx:1973
void listAttributes(std::string const &group_or_dataset, Container &container) const 
Insert the attribute names of the given group or dataset into the given container by calling containe...
Definition: hdf5impex.hxx:1590
bool unique() const 
Check if this is the unique owner of the conatined handle. 
Definition: hdf5impex.hxx:577
void write(std::string datasetName, char data)
Write a single value. Specialization of the write function for simple datatypes. 
Definition: hdf5impex.hxx:1944
const char * getPixelType() const 
void writeAttribute(std::string object_name, std::string attribute_name, char data)
Write a single value. Specialization of the write function for simple datatypes. 
Definition: hdf5impex.hxx:1646
void writeAttribute(std::string object_name, std::string attribute_name, const MultiArrayView< N, T, Stride > &array)
Write MultiArray Attributes. In contrast to datasets, subarray access, chunks and compression are not...
Definition: hdf5impex.hxx:1611
herr_t close()
Explicitly call the stored destructor (if one has been registered in the constructor) for the contain...
Definition: hdf5impex.hxx:296
void write(std::string datasetName, const MultiArrayView< N, T, Stride > &array, int iChunkSize=0, int compression=0)
Write multi arrays. 
Definition: hdf5impex.hxx:1791
image import and export functions 
bool isUnstrided(unsigned int dimension=N-1) const 
Definition: multi_array.hxx:1286
void close()
Close the current file. 
Definition: hdf5impex.hxx:1199
bool cd_up()
Change the current group to its parent group. Returns true if successful, false otherwise. If unsuccessful, the group will not change. 
Definition: hdf5impex.hxx:1225
Base class for, and view to, vigra::MultiArray. 
Definition: multi_array.hxx:704
std::string getDatasetType(std::string const &datasetName) const 
Definition: hdf5impex.hxx:1473
bool existsDataset(std::string datasetName) const 
Check if given datasetName exists. 
Definition: hdf5impex.hxx:1354
const_iterator end() const 
Definition: array_vector.hxx:237
const_pointer data() const 
Definition: array_vector.hxx:209
ArrayVector< hsize_t > const & shape() const 
size_type size() const 
Definition: array_vector.hxx:358
MultiArrayView subarray(difference_type p, difference_type q) const 
Definition: multi_array.hxx:1528
void readBlock(std::string datasetName, typename MultiArrayShape< N >::type blockOffset, typename MultiArrayShape< N >::type blockShape, MultiArrayView< N, T, Stride > array)
Read a block of data into a multi array. This function allows to read a small block out of a larger v...
Definition: hdf5impex.hxx:2068
Class for a single RGB value. 
Definition: accessor.hxx:938
int ceil(FixedPoint< IntBits, FracBits > v)
rounding up. 
Definition: fixedpoint.hxx:675
hid_t getDatasetHandle() const 
HDF5HandleShared(hid_t h, Destructor destructor, const char *error_message)
Create a shared wrapper for a plain hid_t object. 
Definition: hdf5impex.hxx:471
HDF5File(std::string filePath, OpenMode mode=ReadOnly, bool track_creation_times=false)
Open or create an HDF5File object. 
Definition: hdf5impex.hxx:1065
hid_t getH5FileHandle() const 
~HDF5HandleShared()
Destructor (calls close()) 
Definition: hdf5impex.hxx:514
HDF5Handle getDatasetHandle(std::string const &datasetName) const 
Obtain the HDF5 handle of a dataset. 
Definition: hdf5impex.hxx:1519
void root()
Change current group to "/". 
Definition: hdf5impex.hxx:1207
HDF5File(bool track_creation_times)
Construct with time tagging of datasets enabled. 
Definition: hdf5impex.hxx:1055
void cd(std::string groupName)
Change the current group. Both absolute and relative group names are allowed. 
Definition: hdf5impex.hxx:1216
size_t use_count() const 
Get the number of owners of the contained handle. 
Definition: hdf5impex.hxx:566
void swap(HDF5Handle &h)
Swap the contents of two handle wrappers. 
Definition: hdf5impex.hxx:337
Access to HDF5 files. 
Definition: hdf5impex.hxx:974
std::vector< std::string > ls() const 
List the contents of the current group. The function returns a vector of strings holding the entries ...
Definition: hdf5impex.hxx:1310