HighFive 2.7.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#pragma once
10
11#include <algorithm>
12#include <cassert>
13#include <functional>
14#include <numeric>
15#include <sstream>
16#include <string>
17
18#include <H5Dpublic.h>
19#include <H5Ppublic.h>
20
21#include "H5ReadWrite_misc.hpp"
22#include "H5Converter_misc.hpp"
23
24namespace HighFive {
25
26namespace details {
27
28// map the correct reference to the dataset depending of the layout
29// dataset -> itself
30// subselection -> parent dataset
31inline const DataSet& get_dataset(const Selection& sel) {
32 return sel.getDataset();
33}
34
35inline const DataSet& get_dataset(const DataSet& ds) {
36 return ds;
37}
38
39// map the correct memspace identifier depending of the layout
40// dataset -> entire memspace
41// selection -> resolve space id
42inline hid_t get_memspace_id(const Selection& ptr) {
43 return ptr.getMemSpace().getId();
44}
45
46inline hid_t get_memspace_id(const DataSet&) {
47 return H5S_ALL;
48}
49} // namespace details
50
51inline ElementSet::ElementSet(std::initializer_list<std::size_t> list)
52 : _ids(list) {}
53
54inline ElementSet::ElementSet(std::initializer_list<std::vector<std::size_t>> list)
55 : ElementSet(std::vector<std::vector<std::size_t>>(list)) {}
56
57inline ElementSet::ElementSet(const std::vector<std::size_t>& element_ids)
58 : _ids(element_ids) {}
59
60inline ElementSet::ElementSet(const std::vector<std::vector<std::size_t>>& element_ids) {
61 for (const auto& vec: element_ids) {
62 std::copy(vec.begin(), vec.end(), std::back_inserter(_ids));
63 }
64}
65
66template <typename Derivate>
68 const DataSpace& memspace) const {
69 // Note: The current limitation are that memspace must describe a
70 // packed memspace.
71 //
72 // The reason for this is that we're unable to unpack general
73 // hyperslabs when the memory is not contiguous, e.g.
74 // `std::vector<std::vector<double>>`.
75 const auto& slice = static_cast<const Derivate&>(*this);
76 auto filespace = hyperslab.apply(slice.getSpace());
77
78 return detail::make_selection(memspace, filespace, details::get_dataset(slice));
79}
80
81template <typename Derivate>
82inline Selection SliceTraits<Derivate>::select(const HyperSlab& hyper_slab) const {
83 const auto& slice = static_cast<const Derivate&>(*this);
84 auto filespace = slice.getSpace();
85 filespace = hyper_slab.apply(filespace);
86
87 auto n_elements = H5Sget_select_npoints(filespace.getId());
88 auto memspace = DataSpace(std::array<size_t, 1>{size_t(n_elements)});
89
90 return detail::make_selection(memspace, filespace, details::get_dataset(slice));
91}
92
93
94template <typename Derivate>
95inline Selection SliceTraits<Derivate>::select(const std::vector<size_t>& offset,
96 const std::vector<size_t>& count,
97 const std::vector<size_t>& stride,
98 const std::vector<size_t>& block) const {
99 auto slab = HyperSlab(RegularHyperSlab(offset, count, stride, block));
100 auto memspace = DataSpace(count);
101 return select_impl(slab, memspace);
102}
103
104template <typename Derivate>
105inline Selection SliceTraits<Derivate>::select(const std::vector<size_t>& columns) const {
106 const auto& slice = static_cast<const Derivate&>(*this);
107 const DataSpace& space = slice.getSpace();
108 std::vector<size_t> dims = space.getDimensions();
109
110 std::vector<size_t> counts = dims;
111 counts.back() = 1;
112
113 std::vector<size_t> offsets(dims.size(), 0);
114
115 HyperSlab slab;
116 for (const auto& column: columns) {
117 offsets.back() = column;
118 slab |= RegularHyperSlab(offsets, counts);
119 }
120
121 std::vector<size_t> memdims = dims;
122 memdims.back() = columns.size();
123
124 return select_impl(slab, DataSpace(memdims));
125}
126
127template <typename Derivate>
129 const auto& slice = static_cast<const Derivate&>(*this);
130 const hsize_t* data = nullptr;
131 const DataSpace space = slice.getSpace().clone();
132 const std::size_t length = elements._ids.size();
133 if (length % space.getNumberDimensions() != 0) {
134 throw DataSpaceException(
135 "Number of coordinates in elements picking "
136 "should be a multiple of the dimensions.");
137 }
138 const std::size_t num_elements = length / space.getNumberDimensions();
139 std::vector<hsize_t> raw_elements;
140
141 // optimised at compile time
142 // switch for data conversion on 32bits platforms
143 if (std::is_same<std::size_t, hsize_t>::value) {
144 // `if constexpr` can't be used, thus a reinterpret_cast is needed.
145 data = reinterpret_cast<const hsize_t*>(&(elements._ids[0]));
146 } else {
147 raw_elements.resize(length);
148 std::copy(elements._ids.begin(), elements._ids.end(), raw_elements.begin());
149 data = raw_elements.data();
150 }
151
152 if (H5Sselect_elements(space.getId(), H5S_SELECT_SET, num_elements, data) < 0) {
153 HDF5ErrMapper::ToException<DataSpaceException>("Unable to select elements");
154 }
155
156 return detail::make_selection(DataSpace(num_elements), space, details::get_dataset(slice));
157}
158
159
160template <typename Derivate>
161template <typename T>
162inline T SliceTraits<Derivate>::read(const DataTransferProps& xfer_props) const {
163 T array;
164 read(array, xfer_props);
165 return array;
166}
167
168
169template <typename Derivate>
170template <typename T>
171inline void SliceTraits<Derivate>::read(T& array, const DataTransferProps& xfer_props) const {
172 const auto& slice = static_cast<const Derivate&>(*this);
173 const DataSpace& mem_space = slice.getMemSpace();
174
175 const details::BufferInfo<T> buffer_info(
176 slice.getDataType(),
177 [&slice]() -> std::string { return details::get_dataset(slice).getPath(); },
178 details::BufferInfo<T>::Operation::read);
179
180 if (!details::checkDimensions(mem_space, buffer_info.n_dimensions)) {
181 std::ostringstream ss;
182 ss << "Impossible to read DataSet of dimensions " << mem_space.getNumberDimensions()
183 << " into arrays of dimensions " << buffer_info.n_dimensions;
184 throw DataSpaceException(ss.str());
185 }
186 auto dims = mem_space.getDimensions();
187
188 if (mem_space.getElementCount() == 0) {
189 auto effective_dims = details::squeezeDimensions(dims,
190 details::inspector<T>::recursive_ndim);
191
192 details::inspector<T>::prepare(array, effective_dims);
193 return;
194 }
195
196 auto r = details::data_converter::get_reader<T>(dims, array);
197 read(r.get_pointer(), buffer_info.data_type, xfer_props);
198 // re-arrange results
199 r.unserialize();
200 auto t = create_datatype<typename details::inspector<T>::base_type>();
201 auto c = t.getClass();
202 if (c == DataTypeClass::VarLen || t.isVariableStr()) {
203#if H5_VERSION_GE(1, 12, 0)
204 // This one have been created in 1.12.0
205 (void) H5Treclaim(t.getId(), mem_space.getId(), xfer_props.getId(), r.get_pointer());
206#else
207 // This one is deprecated since 1.12.0
208 (void) H5Dvlen_reclaim(t.getId(), mem_space.getId(), xfer_props.getId(), r.get_pointer());
209#endif
210 }
211}
212
213
214template <typename Derivate>
215template <typename T>
216inline void SliceTraits<Derivate>::read(T* array,
217 const DataType& dtype,
218 const DataTransferProps& xfer_props) const {
219 static_assert(!std::is_const<T>::value,
220 "read() requires a non-const structure to read data into");
221 const auto& slice = static_cast<const Derivate&>(*this);
222 using element_type = typename details::inspector<T>::base_type;
223
224 // Auto-detect mem datatype if not provided
225 const DataType& mem_datatype = dtype.empty() ? create_and_check_datatype<element_type>()
226 : dtype;
227
228 if (H5Dread(details::get_dataset(slice).getId(),
229 mem_datatype.getId(),
230 details::get_memspace_id(slice),
231 slice.getSpace().getId(),
232 xfer_props.getId(),
233 static_cast<void*>(array)) < 0) {
234 HDF5ErrMapper::ToException<DataSetException>("Error during HDF5 Read.");
235 }
236}
237
238
239template <typename Derivate>
240template <typename T>
241inline void SliceTraits<Derivate>::write(const T& buffer, const DataTransferProps& xfer_props) {
242 const auto& slice = static_cast<const Derivate&>(*this);
243 const DataSpace& mem_space = slice.getMemSpace();
244
245 if (mem_space.getElementCount() == 0) {
246 return;
247 }
248
249 const details::BufferInfo<T> buffer_info(
250 slice.getDataType(),
251 [&slice]() -> std::string { return details::get_dataset(slice).getPath(); },
252 details::BufferInfo<T>::Operation::write);
253
254 if (!details::checkDimensions(mem_space, buffer_info.n_dimensions)) {
255 std::ostringstream ss;
256 ss << "Impossible to write buffer of dimensions "
257 << details::format_vector(mem_space.getDimensions())
258 << " into dataset with n = " << buffer_info.n_dimensions << " dimensions.";
259 throw DataSpaceException(ss.str());
260 }
261 auto w = details::data_converter::serialize<T>(buffer);
262 write_raw(w.get_pointer(), buffer_info.data_type, xfer_props);
263}
264
265
266template <typename Derivate>
267template <typename T>
268inline void SliceTraits<Derivate>::write_raw(const T* buffer,
269 const DataType& dtype,
270 const DataTransferProps& xfer_props) {
271 using element_type = typename details::inspector<T>::base_type;
272 const auto& slice = static_cast<const Derivate&>(*this);
273 const auto& mem_datatype = dtype.empty() ? create_and_check_datatype<element_type>() : dtype;
274
275 if (H5Dwrite(details::get_dataset(slice).getId(),
276 mem_datatype.getId(),
277 details::get_memspace_id(slice),
278 slice.getSpace().getId(),
279 xfer_props.getId(),
280 static_cast<const void*>(buffer)) < 0) {
281 HDF5ErrMapper::ToException<DataSetException>("Error during HDF5 Write: ");
282 }
284
285
286} // namespace HighFive
Exception specific to HighFive DataSpace interface.
Definition H5Exception.hpp:112
Class representing the space (dimensions) of a dataset.
Definition H5DataSpace.hpp:25
size_t getNumberDimensions() const
getNumberDimensions
Definition H5Dataspace_misc.hpp:93
size_t getElementCount() const
getElementCount
Definition H5Dataspace_misc.hpp:112
std::vector< size_t > getDimensions() const
getDimensions
Definition H5Dataspace_misc.hpp:102
DataSpace clone() const
Definition H5Dataspace_misc.hpp:85
HDF5 Data Type.
Definition H5DataType.hpp:54
bool empty() const noexcept
Check the DataType was default constructed. Such value might represent auto-detection of the datatype...
Definition H5DataType_misc.hpp:35
Definition H5Slice_traits.hpp:21
ElementSet(std::initializer_list< std::size_t > list)
Create a list of points of N-dimension for selection.
Definition H5Slice_traits_misc.hpp:51
Definition H5Slice_traits.hpp:120
DataSpace apply(const DataSpace &space_) const
Definition H5Slice_traits.hpp:173
hid_t getId() const noexcept
getId
Definition H5Object_misc.hpp:65
HDF5 property Lists.
Definition H5PropertyList.hpp:79
Selection: represent a view on a slice/part of a dataset.
Definition H5Selection.hpp:27
DataSpace getSpace() const noexcept
getSpace
Definition H5Selection_misc.hpp:20
T read(const DataTransferProps &xfer_props=DataTransferProps()) const
Definition H5Slice_traits_misc.hpp:162
void write(const T &buffer, const DataTransferProps &xfer_props=DataTransferProps())
Definition H5Slice_traits_misc.hpp:241
void write_raw(const T *buffer, const DataType &dtype=DataType(), const DataTransferProps &xfer_props=DataTransferProps())
Definition H5Slice_traits_misc.hpp:268
Selection select(const HyperSlab &hyperslab) const
Select an hyperslab in the current Slice/Dataset.
Definition H5Slice_traits_misc.hpp:82
Selection select_impl(const HyperSlab &hyperslab, const DataSpace &memspace) const
Definition H5Slice_traits_misc.hpp:67
Definition H5_definitions.hpp:15
Definition H5Slice_traits.hpp:72