9 #ifndef H5SLICE_TRAITS_MISC_HPP
10 #define H5SLICE_TRAITS_MISC_HPP
21 #include <boost/multi_array.hpp>
22 #include <boost/numeric/ublas/matrix.hpp>
23 #include <boost/serialization/vector.hpp>
26 #include <H5Dpublic.h>
27 #include <H5Ppublic.h>
39 inline const DataSet& get_dataset(
const Selection& sel) {
40 return sel.getDataset();
43 inline const DataSet& get_dataset(
const DataSet& ds) {
50 inline hid_t get_memspace_id(
const Selection& ptr) {
51 return ptr.getMemSpace().getId();
54 inline hid_t get_memspace_id(
const DataSet&) {
63 :
ElementSet(std::vector<std::vector<std::size_t>>(list)) {}
66 : _ids(element_ids) {}
69 for (
const auto& vec : element_ids) {
70 std::copy(vec.begin(), vec.end(), std::back_inserter(_ids));
74 template <
typename Derivate>
76 const std::vector<size_t>& count,
77 const std::vector<size_t>& stride)
const {
80 const auto& slice =
static_cast<const Derivate&
>(*this);
81 std::vector<hsize_t> offset_local(offset.size());
82 std::vector<hsize_t> count_local(count.size());
83 std::vector<hsize_t> stride_local(stride.size());
84 std::copy(offset.begin(), offset.end(), offset_local.begin());
85 std::copy(count.begin(), count.end(), count_local.begin());
86 std::copy(stride.begin(), stride.end(), stride_local.begin());
89 if (H5Sselect_hyperslab(space.
getId(), H5S_SELECT_SET, offset_local.data(),
90 stride.empty() ? NULL : stride_local.data(),
91 count_local.data(), NULL) < 0) {
92 HDF5ErrMapper::ToException<DataSpaceException>(
"Unable to select hyperslap");
98 template <
typename Derivate>
100 const auto& slice =
static_cast<const Derivate&
>(*this);
102 const DataSet& dataset = details::get_dataset(slice);
104 std::vector<hsize_t> counts(dims.size());
105 std::copy(dims.begin(), dims.end(), counts.begin());
106 counts[dims.size() - 1] = 1;
107 std::vector<hsize_t> offsets(dims.size(), 0);
109 H5Sselect_none(space.
getId());
111 for (
const auto& column : columns) {
112 offsets[offsets.size() - 1] = column;
114 if (H5Sselect_hyperslab(space.
getId(), H5S_SELECT_OR, offsets.data(), 0,
115 counts.data(), 0) < 0) {
116 HDF5ErrMapper::ToException<DataSpaceException>(
"Unable to select hyperslap");
120 dims[dims.size() - 1] = columns.size();
124 template <
typename Derivate>
126 const auto& slice =
static_cast<const Derivate&
>(*this);
127 const hsize_t* data =
nullptr;
129 const std::size_t length = elements._ids.size();
132 "should be a multiple of the dimensions.");
135 std::vector<hsize_t> raw_elements;
139 if (std::is_same<std::size_t, hsize_t>::value) {
141 data =
reinterpret_cast<const hsize_t*
>(&(elements._ids[0]));
143 raw_elements.resize(length);
144 std::copy(elements._ids.begin(), elements._ids.end(), raw_elements.begin());
145 data = raw_elements.data();
148 if (H5Sselect_elements(space.
getId(), H5S_SELECT_SET, num_elements, data) < 0) {
149 HDF5ErrMapper::ToException<DataSpaceException>(
"Unable to select elements");
156 template <
typename Derivate>
157 template <
typename T>
159 const auto& slice =
static_cast<const Derivate&
>(*this);
160 const DataSpace& mem_space = slice.getMemSpace();
161 const details::BufferInfo<T> buffer_info(slice.getDataType());
163 if (!details::checkDimensions(mem_space, buffer_info.n_dimensions)) {
164 std::ostringstream ss;
165 ss <<
"Impossible to read DataSet of dimensions "
167 << buffer_info.n_dimensions;
170 details::data_converter<T> converter(mem_space);
171 read(converter.transform_read(array), buffer_info.data_type);
173 converter.process_result(array);
177 template <
typename Derivate>
178 template <
typename T>
180 static_assert(!std::is_const<T>::value,
181 "read() requires a non-const structure to read data into");
182 const auto& slice =
static_cast<const Derivate&
>(*this);
183 using element_type =
typename details::type_of_array<T>::type;
187 dtype.
empty() ? create_and_check_datatype<element_type>() : dtype;
189 if (H5Dread(details::get_dataset(slice).getId(),
190 mem_datatype.
getId(),
191 details::get_memspace_id(slice),
192 slice.getSpace().getId(), H5P_DEFAULT,
static_cast<void*
>(array)) < 0) {
193 HDF5ErrMapper::ToException<DataSetException>(
"Error during HDF5 Read: ");
198 template <
typename Derivate>
199 template <
typename T>
201 const auto& slice =
static_cast<const Derivate&
>(*this);
202 const DataSpace& mem_space = slice.getMemSpace();
203 const details::BufferInfo<T> buffer_info(slice.getDataType());
205 if (!details::checkDimensions(mem_space, buffer_info.n_dimensions)) {
206 std::ostringstream ss;
207 ss <<
"Impossible to write buffer of dimensions " << buffer_info.n_dimensions
211 details::data_converter<T> converter(mem_space);
212 write_raw(converter.transform_write(buffer), buffer_info.data_type);
216 template <
typename Derivate>
217 template <
typename T>
219 using element_type =
typename details::type_of_array<T>::type;
220 const auto& slice =
static_cast<const Derivate&
>(*this);
221 const auto& mem_datatype =
222 dtype.
empty() ? create_and_check_datatype<element_type>() : dtype;
224 if (H5Dwrite(details::get_dataset(slice).getId(),
225 mem_datatype.getId(),
226 details::get_memspace_id(slice),
227 slice.getSpace().getId(), H5P_DEFAULT,
228 static_cast<const void*
>(buffer)) < 0) {
229 HDF5ErrMapper::ToException<DataSetException>(
"Error during HDF5 Write: ");
Class representing a dataset.
Definition: H5DataSet.hpp:29
Exception specific to HighFive DataSpace interface.
Definition: H5Exception.hpp:99
Class representing the space (dimensions) of a dataset.
Definition: H5DataSpace.hpp:37
size_t getNumberDimensions() const
getNumberDimensions
Definition: H5Dataspace_misc.hpp:90
std::vector< size_t > getDimensions() const
getDimensions
Definition: H5Dataspace_misc.hpp:99
DataSpace clone() const
Definition: H5Dataspace_misc.hpp:82
HDF5 Data Type.
Definition: H5DataType.hpp:42
bool empty() const noexcept
Check the DataType was default constructed. Such value might represent auto-detection of the datatype...
Definition: H5DataType_misc.hpp:28
Definition: H5Slice_traits.hpp:20
ElementSet(std::initializer_list< std::size_t > list)
Create a list of points of N-dimension for selection.
Definition: H5Slice_traits_misc.hpp:59
hid_t getId() const noexcept
getId
Definition: H5Object_misc.hpp:55
Selection: represent a view on a slice/part of a dataset.
Definition: H5Selection.hpp:23
void read(T &array) const
Definition: H5Slice_traits_misc.hpp:158
void write(const T &buffer)
Definition: H5Slice_traits_misc.hpp:200
Selection select(const std::vector< size_t > &offset, const std::vector< size_t > &count, const std::vector< size_t > &stride=std::vector< size_t >()) const
Select a region in the current Slice/Dataset of count points at offset separated by stride....
Definition: H5Slice_traits_misc.hpp:75
void write_raw(const T *buffer, const DataType &dtype=DataType())
Definition: H5Slice_traits_misc.hpp:218
Definition: H5_definitions.hpp:15