Page MenuHomec4science

hdf5_eigen.inl
No OneTemporary

File Metadata

Created
Sun, Jul 28, 16:23

hdf5_eigen.inl

/* =============================================================================
Copyright (c) 2014 - 2016
F. Georget <fabieng@princeton.edu> Princeton University
All rights reserved.
Redistribution and use in source and binary forms, with or without modification,
are permitted provided that the following conditions are met:
1. Redistributions of source code must retain the above copyright notice,
this list of conditions and the following disclaimer.
2. Redistributions in binary form must reproduce the above copyright notice,
this list of conditions and the following disclaimer in the documentation
and/or other materials provided with the distribution.
3. Neither the name of the copyright holder nor the names of its
contributors may be used to endorse or promote products derived from this
software without specific prior written permission.
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR
ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. *
============================================================================= */
#include "hdf5_eigen.hpp" // for syntax coloring only
//! \file hdf5_eigen.inl
//! \brief implementation of saving a eigen array into hdf5 file
#include "specmicp_hdf5.hpp"
#include "specmicp_common/compat.hpp"
#include <H5Cpp.h>
#include <type_traits>
namespace specmicp {
namespace io {
namespace internal {
//! \brief Implementation of how to save an Eigen matrix in HDF5 format
//!
//! \internal
template <typename Derived>
class SPECMICP_DLL_LOCAL EigenMatrixHDF5Saver
{
public:
EigenMatrixHDF5Saver() {}
static std::unique_ptr<H5::DataSet> save(
HDF5File& file,
const std::string& name,
const std::string& section,
const Eigen::DenseBase<Derived>& matrix
);
static std::unique_ptr<H5::DataSet> save(
H5::Group& file,
const std::string& name,
const Eigen::DenseBase<Derived>& matrix
);
private:
//! \brief Algorithms for a column-major matrix
struct ColumnDispatcher {
static std::unique_ptr<H5::DataSet> save(
H5::Group& group,
const std::string& name,
const Eigen::DenseBase<Derived>& matrix
);
};
//! \brief Algorithms for a row-major matrix
struct RowDispatcher {
static std::unique_ptr<H5::DataSet> save(
H5::Group& group,
const std::string& name,
const Eigen::DenseBase<Derived>& matrix
);
};
//! \brief Algorithms for a column vector
struct VectorDispatcher {
static std::unique_ptr<H5::DataSet> save(
H5::Group& group,
const std::string& name,
const Eigen::DenseBase<Derived>& matrix
);
};
};
//! \brief Implementation of how to read an Eigen matrix in HDF5 format
//!
//! \internal
template <typename Derived>
class SPECMICP_DLL_LOCAL EigenMatrixHDF5Reader
{
public:
static int read(
const HDF5File& file,
const std::string& name,
const std::string& section,
Eigen::DenseBase<Derived>& matrix
);
static int read(
const H5::Group& group,
const std::string& name,
Eigen::DenseBase<Derived>& matrix
);
static int read(
const H5::DataSet& dataset,
Eigen::DenseBase<Derived>& matrix
);
private:
//! \brief Algorithms for a column-major matrix
struct ColumnDispatcher {
static int read(
const H5::DataSet& dataset,
Eigen::DenseBase<Derived>& m_matrix
);
};
//! \brief Algorithms for a row-major matrix
struct RowDispatcher {
static int read(
const H5::DataSet& dataset,
Eigen::DenseBase<Derived>& m_matrix
);
};
//! \brief Algorithms for a column vector
struct VectorDispatcher {
static int read(
const H5::DataSet& dataset,
Eigen::DenseBase<Derived>& m_matrix
);
};
};
} //end namespace internal
// implementation
template <typename Derived>
std::unique_ptr<H5::DataSet> save_eigen_matrix(
HDF5File& file,
const std::string& name,
const std::string& section,
const Eigen::DenseBase<Derived>& matrix
)
{
return internal::EigenMatrixHDF5Saver<Derived>::save(
file, name, section, matrix);
}
// implementation
template <typename Derived>
std::unique_ptr<H5::DataSet> save_eigen_matrix(
H5::Group& group,
const std::string& name,
const Eigen::DenseBase<Derived>& matrix
)
{
return internal::EigenMatrixHDF5Saver<Derived>::save(group, name, matrix);
}
template <typename Derived>
int read_eigen_matrix(
const HDF5File& file,
const std::string& name,
const std::string& section,
Eigen::DenseBase<Derived>& matrix
)
{
return internal::EigenMatrixHDF5Reader<Derived>::read(
file, name, section, matrix);
}
template <typename Derived>
int read_eigen_matrix(
const H5::Group& group,
const std::string& name,
Eigen::DenseBase<Derived>& matrix
)
{
return internal::EigenMatrixHDF5Reader<Derived>::read(
group, name, matrix);
}
template <typename Derived>
int read_eigen_matrix(
const H5::DataSet& dataset,
Eigen::DenseBase<Derived>& matrix
)
{
return internal::EigenMatrixHDF5Reader<Derived>::read(
dataset, matrix);
}
// ========================
//
// Implementation
//
// ========================
namespace internal {
// HDF5 type
// ========
// Return a HDF5 type from the type of the matrix
//! \brief Return the HDF5 DataType from a type
template <typename T>
struct DataType;
template <>
struct DataType<double>
{
static inline H5::PredType get() {
return H5::PredType::NATIVE_DOUBLE;
}
};
template <>
struct DataType<float>
{
static inline H5::PredType get() {
return H5::PredType::NATIVE_FLOAT;
}
};
template <>
struct DataType<int>
{
static inline H5::PredType get() {
return H5::PredType::NATIVE_INT;
}
};
template <>
struct DataType<long int>
{
static inline H5::PredType get() {
return H5::PredType::NATIVE_LONG;
}
};
//! \brief Return a HDF5 type from the type of the matrix
template <typename Derived>
H5::PredType get_datatype()
{
return DataType<typename Eigen::DenseBase<Derived>::Scalar>::get();
}
// Save functions
// ===============
template <typename Derived>
std::unique_ptr<H5::DataSet> EigenMatrixHDF5Saver<Derived>::save(
HDF5File& file,
const std::string& name,
const std::string& section,
const Eigen::DenseBase<Derived>& matrix
)
{
auto group = file.open_group(section);
return EigenMatrixHDF5Saver<Derived>::save(*group.get(), name, matrix);
}
template <typename Derived>
std::unique_ptr<H5::DataSet> EigenMatrixHDF5Saver<Derived>::save(
H5::Group& group,
const std::string& name,
const Eigen::DenseBase<Derived>& matrix
)
{
// this awful looking line just check the compile-time information of the
// matrix and call the correct algorithm (*Dispatcher::save())
// if vector
return std::conditional<Eigen::DenseBase<Derived>::ColsAtCompileTime == 1,
// then
VectorDispatcher,
// else if row major matrix
typename std::conditional<Eigen::DenseBase<Derived>::IsRowMajor,
// then
RowDispatcher,
// else (colum major matrix)
ColumnDispatcher
>::type
>::type::save(group, name, matrix);
}
template <typename Derived>
std::unique_ptr<H5::DataSet>
EigenMatrixHDF5Saver<Derived>::ColumnDispatcher::save(
H5::Group& group,
const std::string& name,
const Eigen::DenseBase<Derived>& matrix
)
{
H5::DSetCreatPropList plist;
H5::PredType scalar_type = get_datatype<Derived>();
// assume column major for now
auto rows = static_cast<hsize_t>(matrix.rows());
auto cols = static_cast<hsize_t>(matrix.cols());
auto stride = static_cast<hsize_t>(matrix.outerStride()); // should be rows for a column major
assert(stride == rows);
// file space
// read a col
hsize_t fdim[] = {rows, cols};
H5::DataSpace fspace(2, fdim);
hsize_t fstride[] = {1, cols};
hsize_t fcount[] = {1, 1};
hsize_t fblock[] = {1, cols};
// matrix space
// write a col
hsize_t mdim[] = {cols, stride};
H5::DataSpace mspace(2, mdim);
hsize_t mstride[] = {stride, 1};
hsize_t mcount[] = {1, 1};
hsize_t mblock[] = {cols, 1};
// where to write
std::unique_ptr<H5::DataSet> dataset = make_unique<H5::DataSet>(
group.createDataSet(name, scalar_type, fspace, plist));
// column-major to row-major transformation
for (hsize_t row=0; row<fdim[1]; ++row)
{
hsize_t fpos[] = {row, 0};
hsize_t mpos[] = {0, row};
fspace.selectHyperslab(H5S_SELECT_SET, fcount, fpos, fstride, fblock);
mspace.selectHyperslab(H5S_SELECT_SET, mcount, mpos, mstride, mblock);
dataset->write(matrix.derived().data(), scalar_type, mspace, fspace);
}
return dataset;
}
template <typename Derived>
std::unique_ptr<H5::DataSet>
EigenMatrixHDF5Saver<Derived>::RowDispatcher::save(
H5::Group& group,
const std::string& name,
const Eigen::DenseBase<Derived>& matrix
)
{
H5::DSetCreatPropList plist;
H5::PredType scalar_type = get_datatype<Derived>();
// assume column major for now
auto rows = static_cast<hsize_t>(matrix.rows());
auto cols = static_cast<hsize_t>(matrix.cols());
// file space
// read a col
hsize_t fdim[] = {rows, cols};
H5::DataSpace fspace(2, fdim);
// where to write
std::unique_ptr<H5::DataSet> dataset = make_unique<H5::DataSet>(
group.createDataSet(name, scalar_type, fspace, plist));
dataset->write(matrix.derived().data(), scalar_type, fspace);
return dataset;
}
template <typename Derived>
std::unique_ptr<H5::DataSet>
EigenMatrixHDF5Saver<Derived>::VectorDispatcher::save(
H5::Group& group,
const std::string& name,
const Eigen::DenseBase<Derived>& matrix
)
{
H5::DSetCreatPropList plist;
H5::PredType scalar_type = get_datatype<Derived>();
// assume column major for now
auto rows = static_cast<hsize_t>(matrix.rows());
// file space
// read a col
hsize_t fdim[] = {rows};
H5::DataSpace fspace(1, fdim);
// where to write
std::unique_ptr<H5::DataSet> dataset = make_unique<H5::DataSet>(
group.createDataSet(name, scalar_type, fspace, plist));
dataset->write(matrix.derived().data(), scalar_type, fspace);
return dataset;
}
// Read function
// =============
template <typename Derived>
int EigenMatrixHDF5Reader<Derived>::read(
const HDF5File& file,
const std::string& name,
const std::string& section,
Eigen::DenseBase<Derived>& matrix
)
{
auto dataset = file.open_dataset(file.complete_name(name, section));
return read(*dataset.get(), matrix);
}
template <typename Derived>
int EigenMatrixHDF5Reader<Derived>::read(
const H5::Group& group,
const std::string& name,
Eigen::DenseBase<Derived>& matrix
)
{
auto dataset = group.openDataSet(name);
return read(dataset, matrix);
}
template <typename Derived>
int EigenMatrixHDF5Reader<Derived>::read(
const H5::DataSet& dataset,
Eigen::DenseBase<Derived>& matrix
)
{
// this awful looking line just check the compile-time information of the
// matrix and call the correct algorithm (*Dispatcher::save())
// if vector
return std::conditional<Eigen::DenseBase<Derived>::ColsAtCompileTime == 1,
// then
VectorDispatcher,
// else if row major matrix
typename std::conditional<Eigen::DenseBase<Derived>::IsRowMajor,
// then
RowDispatcher,
// else (colum major matrix)
ColumnDispatcher
>::type
>::type::read(dataset, matrix);
}
template <typename Derived>
int
EigenMatrixHDF5Reader<Derived>::ColumnDispatcher::read(
const H5::DataSet& dataset,
Eigen::DenseBase<Derived>& matrix
)
{
// read in row major
Eigen::Matrix<
typename Eigen::DenseBase<Derived>::Scalar,
Eigen::Dynamic, Eigen::Dynamic,
Eigen::RowMajor> row_maj_matrix;
H5::DataSpace dataspace = dataset.getSpace();
int rank = dataspace.getSimpleExtentNdims();
if (rank != 2)
{
return -1;
}
hsize_t dim[2];
dataspace.getSimpleExtentDims(dim);
row_maj_matrix.resize(dim[0], dim[1]);
H5::PredType scalar_type = get_datatype<Derived>();
dataset.read(row_maj_matrix.data(), scalar_type);
// then transform
matrix = row_maj_matrix;
return 0;
}
template <typename Derived>
int
EigenMatrixHDF5Reader<Derived>::RowDispatcher::read(
const H5::DataSet& dataset,
Eigen::DenseBase<Derived>& matrix
)
{
H5::DataSpace dataspace = dataset.getSpace();
int rank = dataspace.getSimpleExtentNdims();
if (rank != 2)
{
return -1;
}
hsize_t dim[2];
dataspace.getSimpleExtentDims(dim);
matrix.derived().resize(dim[0], dim[1]);
H5::PredType scalar_type = get_datatype<Derived>();
dataset.read(matrix.derived().data(), scalar_type);
return 1;
}
template <typename Derived>
int
EigenMatrixHDF5Reader<Derived>::VectorDispatcher::read(
const H5::DataSet& dataset,
Eigen::DenseBase<Derived>& matrix
)
{
H5::DataSpace dataspace = dataset.getSpace();
int rank = dataspace.getSimpleExtentNdims();
if (rank != 1)
{
return -1;
}
hsize_t dim[1];
dataspace.getSimpleExtentDims(dim);
matrix.derived().resize(dim[0], 1);
H5::PredType scalar_type = get_datatype<Derived>();
dataset.read(matrix.derived().data(), scalar_type);
return 1;
}
} //end namespace internal
} //end namespace io
} //end namespace specmicp

Event Timeline