HighFive 2.3.1
HighFive - Header-only C++ HDF5 interface
Loading...
Searching...
No Matches
H5Slice_traits_misc.hpp
Go to the documentation of this file.
1/*
2 * Copyright (c), 2017, Adrien Devresse <adrien.devresse@epfl.ch>
3 *
4 * Distributed under the Boost Software License, Version 1.0.
5 * (See accompanying file LICENSE_1_0.txt or copy at
6 * http://www.boost.org/LICENSE_1_0.txt)
7 *
8 */
9#ifndef H5SLICE_TRAITS_MISC_HPP
10#define H5SLICE_TRAITS_MISC_HPP
11
12#include <algorithm>
13#include <cassert>
14#include <functional>
15#include <numeric>
16#include <sstream>
17#include <string>
18
19#ifdef H5_USE_BOOST
20// starting Boost 1.64, serialization header must come before ublas
21#include <boost/multi_array.hpp>
22#include <boost/numeric/ublas/matrix.hpp>
23#include <boost/serialization/vector.hpp>
24#endif
25
26#include <H5Dpublic.h>
27#include <H5Ppublic.h>
28
29#include "H5ReadWrite_misc.hpp"
30#include "H5Converter_misc.hpp"
31
32namespace HighFive {
33
34namespace details {
35
36// map the correct reference to the dataset depending of the layout
37// dataset -> itself
38// subselection -> parent dataset
39inline const DataSet& get_dataset(const Selection& sel) {
40 return sel.getDataset();
41}
42
43inline const DataSet& get_dataset(const DataSet& ds) {
44 return ds;
45}
46
47// map the correct memspace identifier depending of the layout
48// dataset -> entire memspace
49// selection -> resolve space id
50inline hid_t get_memspace_id(const Selection& ptr) {
51 return ptr.getMemSpace().getId();
52}
53
54inline hid_t get_memspace_id(const DataSet&) {
55 return H5S_ALL;
56}
57} // namespace details
58
59inline ElementSet::ElementSet(std::initializer_list<std::size_t> list)
60 : _ids(list) {}
61
62inline ElementSet::ElementSet(std::initializer_list<std::vector<std::size_t>> list)
63 : ElementSet(std::vector<std::vector<std::size_t>>(list)) {}
64
65inline ElementSet::ElementSet(const std::vector<std::size_t>& element_ids)
66 : _ids(element_ids) {}
67
68inline ElementSet::ElementSet(const std::vector<std::vector<std::size_t>>& element_ids) {
69 for (const auto& vec : element_ids) {
70 std::copy(vec.begin(), vec.end(), std::back_inserter(_ids));
71 }
73
74template <typename Derivate>
75inline Selection SliceTraits<Derivate>::select(const std::vector<size_t>& offset,
76 const std::vector<size_t>& count,
77 const std::vector<size_t>& stride) const {
78 // hsize_t type conversion
79 // TODO : normalize hsize_t type in HighFive namespace
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());
87
88 DataSpace space = slice.getSpace().clone();
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");
93 }
94
95 return Selection(DataSpace(count), space, details::get_dataset(slice));
96}
97
98template <typename Derivate>
99inline Selection SliceTraits<Derivate>::select(const std::vector<size_t>& columns) const {
100 const auto& slice = static_cast<const Derivate&>(*this);
101 const DataSpace& space = slice.getSpace();
102 const DataSet& dataset = details::get_dataset(slice);
103 std::vector<size_t> dims = space.getDimensions();
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);
108
109 H5Sselect_none(space.getId());
111 for (const auto& column : columns) {
112 offsets[offsets.size() - 1] = column;
113
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");
117 }
118 }
119
120 dims[dims.size() - 1] = columns.size();
121 return Selection(DataSpace(dims), space, dataset);
123
124// no data conversion on 64bits platforms
125template <typename T>
126typename std::enable_if<std::is_same<std::size_t, T>::value>::type
128 typename std::vector<T>&,
129 const std::size_t,
130 const std::vector<std::size_t>& element_ids) {
131 data = reinterpret_cast<const T*>(&(element_ids[0]));
132}
133
134// data conversion on 32bits platforms
135template <typename T>
136typename std::enable_if<!std::is_same<std::size_t, T>::value>::type
138 typename std::vector<T>& raw_elements,
139 const std::size_t length,
140 const std::vector<std::size_t>& element_ids) {
141 raw_elements.resize(length);
142 std::copy(element_ids.begin(), element_ids.end(), raw_elements.begin());
143 data = raw_elements.data();
144}
145
146template <typename Derivate>
148 const auto& slice = static_cast<const Derivate&>(*this);
149 const hsize_t* data = nullptr;
150 const DataSpace space = slice.getSpace().clone();
151 const std::size_t length = elements._ids.size();
152 if (length % space.getNumberDimensions() != 0) {
153 throw DataSpaceException("Number of coordinates in elements picking "
154 "should be a multiple of the dimensions.");
155 }
156 const std::size_t num_elements = length / space.getNumberDimensions();
157 std::vector<hsize_t> raw_elements;
158
159 // optimised at compile time
160 access_with_conversion<>(data, raw_elements, length, elements._ids);
161
162 if (H5Sselect_elements(space.getId(), H5S_SELECT_SET, num_elements, data) < 0) {
163 HDF5ErrMapper::ToException<DataSpaceException>("Unable to select elements");
164 }
165
166 return Selection(DataSpace(num_elements), space, details::get_dataset(slice));
167}
168
169template <typename Derivate>
170template <typename T>
171inline void SliceTraits<Derivate>::read(T& array) const {
172 const auto& slice = static_cast<const Derivate&>(*this);
173 const DataSpace& mem_space = slice.getMemSpace();
174 const details::BufferInfo<T> buffer_info(slice.getDataType(),
175 [slice]() -> std::string { return details::get_dataset(slice).getPath(); });
176
177 if (!details::checkDimensions(mem_space, buffer_info.n_dimensions)) {
178 std::ostringstream ss;
179 ss << "Impossible to read DataSet of dimensions "
180 << mem_space.getNumberDimensions() << " into arrays of dimensions "
181 << buffer_info.n_dimensions;
182 throw DataSpaceException(ss.str());
183 }
184 details::data_converter<T> converter(mem_space);
185 read(converter.transform_read(array), buffer_info.data_type);
186 // re-arrange results
187 converter.process_result(array);
188}
189
190
191template <typename Derivate>
192template <typename T>
193inline void SliceTraits<Derivate>::read(T* array, const DataType& dtype) const {
194 static_assert(!std::is_const<T>::value,
195 "read() requires a non-const structure to read data into");
196 const auto& slice = static_cast<const Derivate&>(*this);
197 using element_type = typename details::inspector<T>::base_type;
198
199 // Auto-detect mem datatype if not provided
200 const DataType& mem_datatype =
201 dtype.empty() ? create_and_check_datatype<element_type>() : dtype;
202
203 if (H5Dread(details::get_dataset(slice).getId(),
204 mem_datatype.getId(),
205 details::get_memspace_id(slice),
206 slice.getSpace().getId(), H5P_DEFAULT, static_cast<void*>(array)) < 0) {
207 HDF5ErrMapper::ToException<DataSetException>("Error during HDF5 Read: ");
208 }
209}
210
211
212template <typename Derivate>
213template <typename T>
214inline void SliceTraits<Derivate>::write(const T& buffer) {
215 const auto& slice = static_cast<const Derivate&>(*this);
216 const DataSpace& mem_space = slice.getMemSpace();
217 const details::BufferInfo<T> buffer_info(slice.getDataType(),
218 [slice]() -> std::string { return details::get_dataset(slice).getPath(); });
219
220 if (!details::checkDimensions(mem_space, buffer_info.n_dimensions)) {
221 std::ostringstream ss;
222 ss << "Impossible to write buffer of dimensions " << buffer_info.n_dimensions
223 << " into dataset of dimensions " << mem_space.getNumberDimensions();
224 throw DataSpaceException(ss.str());
225 }
226 details::data_converter<T> converter(mem_space);
227 write_raw(converter.transform_write(buffer), buffer_info.data_type);
228}
229
230
231template <typename Derivate>
232template <typename T>
233inline void SliceTraits<Derivate>::write_raw(const T* buffer, const DataType& dtype) {
234 using element_type = typename details::inspector<T>::base_type;
235 const auto& slice = static_cast<const Derivate&>(*this);
236 const auto& mem_datatype =
237 dtype.empty() ? create_and_check_datatype<element_type>() : dtype;
238
239 if (H5Dwrite(details::get_dataset(slice).getId(),
240 mem_datatype.getId(),
241 details::get_memspace_id(slice),
242 slice.getSpace().getId(), H5P_DEFAULT,
243 static_cast<const void*>(buffer)) < 0) {
244 HDF5ErrMapper::ToException<DataSetException>("Error during HDF5 Write: ");
245 }
246}
247
248} // namespace HighFive
249
250#endif // H5SLICE_TRAITS_MISC_HPP
Class representing a dataset.
Definition: H5DataSet.hpp:31
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:94
std::vector< size_t > getDimensions() const
getDimensions
Definition: H5Dataspace_misc.hpp:103
DataSpace clone() const
Definition: H5Dataspace_misc.hpp:86
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:171
void write(const T &buffer)
Definition: H5Slice_traits_misc.hpp:214
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:233
Definition: H5_definitions.hpp:15
std::enable_if< std::is_same< std::size_t, T >::value >::type access_with_conversion(const T *&data, typename std::vector< T > &, const std::size_t, const std::vector< std::size_t > &element_ids)
Definition: H5Slice_traits_misc.hpp:127