9#include <Poco/Logger.h>
15#include <boost/numeric/conversion/cast.hpp>
26const auto g_log = &Poco::Logger::get(
"H5Util");
28const std::string NX_ATTR_CLASS(
"NX_class");
29const std::string CAN_SAS_ATTR_CLASS(
"canSAS_class");
36template <
typename NumT> DataType
getType() {
throw DataTypeIException(); }
38template <> MANTID_NEXUS_DLL DataType
getType<float>() {
return PredType::NATIVE_FLOAT; }
40template <> MANTID_NEXUS_DLL DataType
getType<double>() {
return PredType::NATIVE_DOUBLE; }
42template <> MANTID_NEXUS_DLL DataType
getType<int8_t>() {
return PredType::NATIVE_INT8; }
44template <> MANTID_NEXUS_DLL DataType
getType<uint8_t>() {
return PredType::NATIVE_UINT8; }
46template <> MANTID_NEXUS_DLL DataType
getType<int16_t>() {
return PredType::NATIVE_INT16; }
48template <> MANTID_NEXUS_DLL DataType
getType<uint16_t>() {
return PredType::NATIVE_UINT16; }
50template <> MANTID_NEXUS_DLL DataType
getType<int32_t>() {
return PredType::NATIVE_INT32; }
52template <> MANTID_NEXUS_DLL DataType
getType<uint32_t>() {
return PredType::NATIVE_UINT32; }
54template <> MANTID_NEXUS_DLL DataType
getType<int64_t>() {
return PredType::NATIVE_INT64; }
56template <> MANTID_NEXUS_DLL DataType
getType<uint64_t>() {
return PredType::NATIVE_UINT64; }
58template <> MANTID_NEXUS_DLL DataType
getType<char>() {
return PredType::C_S1; }
62template <> MANTID_NEXUS_DLL DataType
getType<bool>() {
return PredType::NATIVE_HBOOL; }
67template <
typename InT,
typename OutT>
68inline void narrowingVectorCast(std::vector<InT>
const &input, std::vector<OutT> &output) {
69 output.resize(input.size());
70 std::transform(input.cbegin(), input.cend(), output.begin(),
71 [](
const auto &a) { return boost::numeric_cast<OutT>(a); });
74template <
typename InT,
typename OutT, Narrowing narrow>
inline void throwIfUnallowedNarrowing() {
76 std::string msg =
"H5Util: narrowing from " + getType<InT>().fromClass();
77 msg +=
" to " + getType<OutT>().fromClass() +
" is prevented.";
78 throw DataTypeIException(msg);
90#define RUN_H5UTIL_FUNCTION(dataType, func_name, ...) \
91 if (dataType == PredType::NATIVE_FLOAT) { \
92 func_name<float, OutT, narrow>(__VA_ARGS__); \
93 } else if (dataType == PredType::NATIVE_DOUBLE) { \
94 func_name<double, OutT, narrow>(__VA_ARGS__); \
95 } else if (dataType == PredType::NATIVE_INT8) { \
96 func_name<int8_t, OutT, narrow>(__VA_ARGS__); \
97 } else if (dataType == PredType::NATIVE_UINT8) { \
98 func_name<uint8_t, OutT, narrow>(__VA_ARGS__); \
99 } else if (dataType == PredType::NATIVE_INT16) { \
100 func_name<int16_t, OutT, narrow>(__VA_ARGS__); \
101 } else if (dataType == PredType::NATIVE_UINT16) { \
102 func_name<uint16_t, OutT, narrow>(__VA_ARGS__); \
103 } else if (dataType == PredType::NATIVE_INT32) { \
104 func_name<int32_t, OutT, narrow>(__VA_ARGS__); \
105 } else if (dataType == PredType::NATIVE_UINT32) { \
106 func_name<uint32_t, OutT, narrow>(__VA_ARGS__); \
107 } else if (dataType == PredType::NATIVE_INT64) { \
108 func_name<int64_t, OutT, narrow>(__VA_ARGS__); \
109 } else if (dataType == PredType::NATIVE_UINT64) { \
110 func_name<uint64_t, OutT, narrow>(__VA_ARGS__); \
111 } else if (dataType == PredType::C_S1) { \
112 func_name<char, OutT, narrow>(__VA_ARGS__); \
114 std::string msg = "H5Util: error in narrowing, unknown H5Cpp PredType "; \
115 msg += dataType.fromClass() + "."; \
116 throw DataTypeIException(msg); \
120 H5::FileAccPropList access_plist;
121 access_plist.setFcloseDegree(H5F_CLOSE_STRONG);
125bool isHdf5(std::string
const &filename) {
130 htri_t ret_value = H5Fis_accessible(filename.c_str(), plist);
135 else if (ret_value == 0)
139 throw H5::FileIException(
"H5File::isHdf5",
"H5Fis_accessible returned negative value");
145 return DataSpace(1, dims);
150template <
typename NumT> H5::DataSet writeScalarDataSet(
Group &
group,
const std::string &
name,
const NumT &
value) {
151 static_assert(std::is_integral<NumT>::value || std::is_floating_point<NumT>::value,
152 "The writeNumAttribute function only accepts integral of "
153 "floating point values.");
154 auto dataType = getType<NumT>();
156 H5::DataSet data =
group.createDataSet(
name, dataType, dataSpace);
157 data.write(&
value, dataType);
162H5::DataSet writeScalarDataSet<std::string>(
Group &
group,
const std::string &
name,
const std::string &
value) {
163 StrType dataType(0,
value.length() + 1);
165 H5::DataSet data =
group.createDataSet(
name, dataType, dataSpace);
166 data.write(
value, dataType);
201 DSetCreatPropList propList;
202 hsize_t chunk_dims[1] = {length};
203 propList.setChunk(1, chunk_dims);
204 propList.setDeflate(deflateLevel);
209 StrType attrType(0, H5T_VARIABLE);
210 DataSpace attrSpace(H5S_SCALAR);
211 auto attribute =
object.createAttribute(
name, attrType, attrSpace);
212 attribute.write(attrType,
value);
215template <
typename NumT>
217 static_assert(std::is_integral<NumT>::value || std::is_floating_point<NumT>::value,
218 "The writeNumAttribute function only accepts integral or "
219 "floating point values.");
220 auto attrType = getType<NumT>();
221 DataSpace attrSpace(H5S_SCALAR);
223 auto attribute =
object.createAttribute(
name, attrType, attrSpace);
225 std::array<NumT, 1> valueArray = {{
value}};
226 attribute.write(attrType, valueArray.data());
229template <
typename NumT>
231 static_assert(std::is_integral<NumT>::value || std::is_floating_point<NumT>::value,
232 "The writeNumAttribute function only accepts integral of "
233 "floating point values.");
234 auto attrType = getType<NumT>();
237 auto attribute =
object.createAttribute(
name, attrType, attrSpace);
238 attribute.write(attrType,
value.data());
245 const std::map<std::string, std::string> &attributes) {
247 for (
const auto &attribute : attributes) {
253 DataType dataType(getType<NumT>());
258 auto data =
group.createDataSet(
name, dataType, dataSpace, propList);
259 data.write(values.data(), dataType);
266std::string
readString(H5::H5File &file,
const std::string &address) {
268 auto data = file.openDataSet(address);
270 }
catch (
const H5::FileIException &e) {
273 }
catch (
const H5::GroupIException &e) {
281 const auto data =
group.openDataSet(
name);
283 }
catch (
const H5::GroupIException &e) {
291 dataset.read(
value, dataset.getDataType(), dataset.getSpace());
293#if H5_VERSION_GE(1, 14, 0)
296 if (
const auto pos =
value.rfind(
'\0'); pos != std::string::npos)
310 DataSet dataset =
group.openDataSet(
name);
317 std::vector<std::string> result;
318 DataSpace dataspace = dataset.getSpace();
319 DataType datatype = dataset.getDataType();
321 dataspace.getSimpleExtentDims(dims,
nullptr);
322 result.reserve(dims[0]);
324 rdata =
new char *[dims[0]];
325 dataset.read(rdata, datatype);
327 for (
size_t i = 0; i < dims[0]; ++i)
328 result.emplace_back(std::string(rdata[i]));
330 dataset.vlenReclaim(rdata, datatype, dataspace);
337bool hasAttribute(
const H5::H5Object &
object,
const char *attributeName) {
338 const htri_t
exists = H5Aexists(
object.getId(), attributeName);
342void readStringAttribute(
const H5::H5Object &
object,
const std::string &attributeName, std::string &output) {
343 const auto attribute =
object.openAttribute(attributeName);
344 attribute.read(attribute.getDataType(), output);
348template <
typename OutT, Narrowing narrow>
351 DataSet dataset =
group.openDataSet(
name);
352 readArray1DCoerce<OutT, narrow>(dataset, output);
353 }
catch (
const H5::GroupIException &) {
355 }
catch (
const H5::DataTypeIException &) {
360template <
typename OutT, Narrowing narrow>
362 std::vector<OutT> result;
364 DataSet dataset =
group.openDataSet(
name);
365 readArray1DCoerce<OutT, narrow>(dataset, result);
366 }
catch (
const H5::GroupIException &) {
368 }
catch (
const H5::DataTypeIException &) {
376template <
typename InT,
typename OutT, Narrowing narrow>
377void convertingVectorRead(
const DataSet &dataset,
const DataType &dataType, std::vector<OutT> &output,
378 const DataSpace &memspace,
const DataSpace &filespace) {
379 throwIfUnallowedNarrowing<InT, OutT, narrow>();
381 std::size_t dataSize = filespace.getSelectNpoints();
382 if constexpr (std::is_same_v<InT, OutT>) {
383 output.resize(dataSize);
384 dataset.read(output.data(), dataType, memspace, filespace);
386 std::vector<InT> temp(dataSize);
387 dataset.read(temp.data(), dataType, memspace, filespace);
388 narrowingVectorCast<InT, OutT>(temp, output);
392template <
typename InT,
typename OutT, Narrowing narrow>
393void convertingNumArrayAttributeRead(
Attribute &attribute, DataType
const &dataType, std::vector<OutT> &
value) {
394 throwIfUnallowedNarrowing<InT, OutT, narrow>();
396 std::size_t dataSize = (attribute.getSpace()).getSelectNpoints();
397 if constexpr (std::is_same_v<InT, OutT>) {
398 value.resize(dataSize);
399 attribute.read(dataType,
value.data());
401 std::vector<InT> temp(dataSize);
402 attribute.read(dataType, temp.data());
403 narrowingVectorCast<InT, OutT>(temp,
value);
407template <
typename InT,
typename OutT, Narrowing narrow>
408void convertingScalarRead(
Attribute &attribute,
const DataType &dataType, OutT &
value) {
409 throwIfUnallowedNarrowing<InT, OutT, narrow>();
411 if constexpr (std::is_same_v<InT, OutT>) {
412 attribute.read(dataType, &
value);
415 attribute.read(dataType, &temp);
416 value = boost::numeric_cast<OutT>(temp);
426template <
typename OutT, Narrowing narrow>
428 auto attribute =
object.openAttribute(attributeName);
429 auto dataType = attribute.getDataType();
440template <
typename OutT, Narrowing narrow>
442 auto attribute =
object.openAttribute(attributeName);
443 auto dataType = attribute.getDataType();
445 std::vector<OutT>
value;
455template <
typename OutT, Narrowing narrow>
456void readArray1DCoerce(
const H5::DataSet &dataset, std::vector<OutT> &output,
const size_t length,
457 const size_t offset) {
458 DataSpace filespace = dataset.getSpace();
459 const auto length_actual =
static_cast<size_t>(filespace.getSelectNpoints());
461 if (offset >= length_actual && offset != 0) {
462 std::stringstream msg;
463 msg <<
"Tried to read offset=" << offset <<
" into array that is only lenght=" << length_actual <<
" long";
464 throw std::runtime_error(msg.str());
469 const hsize_t rankedextent[1] = {
470 static_cast<hsize_t>(std::min(length, length_actual - offset))};
472 if (rankedextent[0] < length_actual)
473 filespace.selectHyperslab(H5S_SELECT_SET, rankedextent, rankedoffset);
476 DataSpace memspace(1, rankedextent);
479 const DataType dataType = dataset.getDataType();
480 RUN_H5UTIL_FUNCTION(dataType, convertingVectorRead, dataset, dataType, output, memspace, filespace);
484bool groupExists(H5::H5Object
const &h5,
const std::string &groupAddress) {
488 h5.openGroup(groupAddress);
489 }
catch (
const H5::Exception &) {
496bool keyHasValue(H5::H5Object
const &h5,
const std::string &key,
const std::string &
value) {
501 attr.read(attr.getDataType(), value_);
504 }
catch (
const H5::Exception &) {
510void copyGroup(H5::H5Object &dest,
const std::string &destGroupAddress, H5::H5Object &src,
511 const std::string &srcGroupAddress) {
514 throw std::invalid_argument(std::string(
"H5Util::copyGroup: source group '") + srcGroupAddress +
515 "' must exist and " +
"destination group '" + destGroupAddress +
"' must not exist.");
524 if (H5Pset_create_intermediate_group(lcpl, 1) < 0)
525 throw std::runtime_error(
"H5Util::copyGroup: 'H5Pset_create_intermediate_group' error return.");
527 if (H5Ocopy(src.getId(), srcGroupAddress.c_str(), dest.getId(), destGroupAddress.c_str(), H5P_DEFAULT, lcpl) < 0)
528 throw std::runtime_error(
"H5Util::copyGroup: 'H5Ocopy' error return.");
537 if (H5Ldelete(h5.getId(), target.c_str(), H5P_DEFAULT) < 0)
538 throw std::runtime_error(
"H5Util::deleteObjectLink: 'H5Ldelete' error return.");
548 const double &
value);
550 const int8_t &
value);
552 const uint8_t &
value);
554 const int16_t &
value);
556 const uint16_t &
value);
558 const int32_t &
value);
560 const uint32_t &
value);
562 const int64_t &
value);
564 const uint64_t &
value);
569 const std::vector<float> &
value);
571 const std::vector<double> &
value);
573 const std::vector<int8_t> &
value);
575 const std::vector<uint8_t> &
value);
577 const std::vector<int16_t> &
value);
579 const std::vector<uint16_t> &
value);
581 const std::vector<int32_t> &
value);
583 const std::vector<uint32_t> &
value);
585 const std::vector<int64_t> &
value);
587 const std::vector<uint64_t> &
value);
589 const std::vector<char> &
value);
600template MANTID_NEXUS_DLL
double readNumAttributeCoerce(
const H5::H5Object &
object,
const std::string &attributeName);
601template MANTID_NEXUS_DLL int8_t
readNumAttributeCoerce(
const H5::H5Object &
object,
const std::string &attributeName);
602template MANTID_NEXUS_DLL uint8_t
readNumAttributeCoerce(
const H5::H5Object &
object,
const std::string &attributeName);
603template MANTID_NEXUS_DLL int16_t
readNumAttributeCoerce(
const H5::H5Object &
object,
const std::string &attributeName);
604template MANTID_NEXUS_DLL uint16_t
readNumAttributeCoerce(
const H5::H5Object &
object,
const std::string &attributeName);
605template MANTID_NEXUS_DLL int32_t
readNumAttributeCoerce(
const H5::H5Object &
object,
const std::string &attributeName);
606template MANTID_NEXUS_DLL uint32_t
readNumAttributeCoerce(
const H5::H5Object &
object,
const std::string &attributeName);
607template MANTID_NEXUS_DLL int64_t
readNumAttributeCoerce(
const H5::H5Object &
object,
const std::string &attributeName);
608template MANTID_NEXUS_DLL uint64_t
readNumAttributeCoerce(
const H5::H5Object &
object,
const std::string &attributeName);
609template MANTID_NEXUS_DLL
char readNumAttributeCoerce(
const H5::H5Object &
object,
const std::string &attributeName);
620 const std::string &attributeName);
622 const std::string &attributeName);
624 const std::string &attributeName);
626 const std::string &attributeName);
628 const std::string &attributeName);
630 const std::string &attributeName);
641 const std::vector<float> &values);
643 const std::vector<double> &values);
645 const std::vector<int32_t> &values);
647 const std::vector<uint32_t> &values);
649 const std::vector<int64_t> &values);
651 const std::vector<uint64_t> &values);
656template MANTID_NEXUS_DLL
void
658 const std::map<std::string, std::string> &attributes);
659template MANTID_NEXUS_DLL
void
661 const std::map<std::string, std::string> &attributes);
662template MANTID_NEXUS_DLL
void
664 const std::map<std::string, std::string> &attributes);
665template MANTID_NEXUS_DLL
void
667 const std::map<std::string, std::string> &attributes);
668template MANTID_NEXUS_DLL
void
670 const std::map<std::string, std::string> &attributes);
671template MANTID_NEXUS_DLL
void
673 const std::map<std::string, std::string> &attributes);
674template MANTID_NEXUS_DLL
void
676 const std::map<std::string, std::string> &attributes);
682 std::vector<float> &output);
684 std::vector<double> &output);
686 std::vector<int32_t> &output);
688 std::vector<uint32_t> &output);
690 std::vector<int64_t> &output);
692 std::vector<uint64_t> &output);
701template MANTID_NEXUS_DLL
void readArray1DCoerce(
const DataSet &dataset, std::vector<float> &output,
702 const size_t length,
const size_t offset);
703template MANTID_NEXUS_DLL
void readArray1DCoerce(
const DataSet &dataset, std::vector<double> &output,
704 const size_t length,
const size_t offset);
705template MANTID_NEXUS_DLL
void readArray1DCoerce(
const DataSet &dataset, std::vector<int8_t> &output,
706 const size_t length,
const size_t offset);
707template MANTID_NEXUS_DLL
void readArray1DCoerce(
const DataSet &dataset, std::vector<uint8_t> &output,
708 const size_t length,
const size_t offset);
709template MANTID_NEXUS_DLL
void readArray1DCoerce(
const DataSet &dataset, std::vector<int16_t> &output,
710 const size_t length,
const size_t offset);
711template MANTID_NEXUS_DLL
void readArray1DCoerce(
const DataSet &dataset, std::vector<uint16_t> &output,
712 const size_t length,
const size_t offset);
713template MANTID_NEXUS_DLL
void readArray1DCoerce(
const DataSet &dataset, std::vector<int32_t> &output,
714 const size_t length,
const size_t offset);
715template MANTID_NEXUS_DLL
void readArray1DCoerce(
const DataSet &dataset, std::vector<uint32_t> &output,
716 const size_t length,
const size_t offset);
717template MANTID_NEXUS_DLL
void readArray1DCoerce(
const DataSet &dataset, std::vector<int64_t> &output,
718 const size_t length,
const size_t offset);
719template MANTID_NEXUS_DLL
void readArray1DCoerce(
const DataSet &dataset, std::vector<uint64_t> &output,
720 const size_t length,
const size_t offset);
721template MANTID_NEXUS_DLL
void readArray1DCoerce(
const DataSet &dataset, std::vector<char> &output,
const size_t length,
722 const size_t offset);
double value
The value of the point.
#define RUN_H5UTIL_FUNCTION(dataType, func_name,...)
This macro matches the dataType in terms of PredType, with the actual primitive datatype.
#define UNUSED_ARG(x)
Function arguments are sometimes unused in certain implmentations but are required for documentation ...
Attribute is a non-fitting parameter.
The class Group represents a set of symmetry operations (or symmetry group).
void information(const std::string &msg)
Logs at information level.
A wrapper class for managing HDF5 object handles (hid_t).
void reset(hid_t const id=INVALID_ID)
Close the existing ID and replace with the new ID; or, set to invalid.
Kernel::Logger g_log("ExperimentInfo")
static logger object
bool exists(Nexus::File &file, const std::string &name)
Operations for help with narrowing casts.
MANTID_NEXUS_DLL DataType getType< uint64_t >()
MANTID_NEXUS_DLL void readStringAttribute(const H5::H5Object &object, const std::string &attributeName, std::string &output)
MANTID_NEXUS_DLL bool groupExists(H5::H5Object const &h5, const std::string &groupAddress)
Test if a group already exists within an HDF5 file or parent group.
MANTID_NEXUS_DLL std::string readString(H5::H5File &file, const std::string &address)
MANTID_NEXUS_DLL DataType getType< bool >()
void writeArray1D(H5::Group &group, const std::string &name, const std::vector< NumT > &values)
MANTID_NEXUS_DLL void writeStrAttribute(const H5::H5Object &object, const std::string &name, const std::string &value)
MANTID_NEXUS_DLL bool hasAttribute(const H5::H5Object &object, const char *attributeName)
MANTID_NEXUS_DLL H5::FileAccPropList defaultFileAcc()
Default file access is H5F_CLOSE_STRONG.
void writeNumAttribute(const H5::H5Object &object, const std::string &name, const NumT &value)
MANTID_NEXUS_DLL bool isHdf5(std::string const &filename)
Determine if a given file can be opened with HDF5 using the CORRECT file access level (H5F_CLOSE_STRO...
MANTID_NEXUS_DLL H5::DataSpace getDataSpace(const size_t length)
Create a 1D data-space to hold data of length.
MANTID_NEXUS_DLL void deleteObjectLink(H5::H5Object &h5, const std::string &target)
Delete a target link for a group or dataset from a parent group.
MANTID_NEXUS_DLL DataType getType< uint8_t >()
MANTID_NEXUS_DLL DataType getType< int64_t >()
MANTID_NEXUS_DLL DataType getType< char >()
std::vector< NumT > readNumArrayAttributeCoerce(const H5::H5Object &object, const std::string &attributeName)
Read a numerical array from an attribute, coerced to type of OutT.
MANTID_NEXUS_DLL DataType getType< std::string >()
MANTID_NEXUS_DLL DataType getType< uint32_t >()
MANTID_NEXUS_DLL H5::Group createGroupCanSAS(H5::Group &group, const std::string &name, const std::string &nxtype, const std::string &cstype)
NumT readNumAttributeCoerce(const H5::H5Object &object, const std::string &attributeName)
Read a single quantity from an attribute, coerced to type of OutT.
MANTID_NEXUS_DLL bool keyHasValue(H5::H5Object const &h5, const std::string &key, const std::string &value)
Test if an attribute is present and has a specific string value for an HDF5 group or dataset.
MANTID_NEXUS_DLL H5::Group createGroupNXS(H5::H5File &file, const std::string &name, const std::string &nxtype)
MANTID_NEXUS_DLL DataType getType< int8_t >()
void readArray1DCoerce(const H5::Group &group, const std::string &name, std::vector< NumT > &output)
MANTID_NEXUS_DLL DataType getType< float >()
MANTID_NEXUS_DLL DataType getType< double >()
MANTID_NEXUS_DLL void copyGroup(H5::H5Object &dest, const std::string &destGroupAddress, H5::H5Object &src, const std::string &srcGroupAddress)
Copy a group and all of its contents, between the same or different HDF5 files or groups.
H5::DataType getType()
Convert a primitive type to the appropriate H5::DataType.
MANTID_NEXUS_DLL DataType getType< uint16_t >()
MANTID_NEXUS_DLL std::vector< std::string > readStringVector(H5::Group &, const std::string &)
MANTID_NEXUS_DLL void write(H5::Group &group, const std::string &name, const std::string &value)
void writeScalarDataSetWithStrAttributes(H5::Group &group, const std::string &name, const T &value, const std::map< std::string, std::string > &attributes)
MANTID_NEXUS_DLL DataType getType< int32_t >()
MANTID_NEXUS_DLL H5::DSetCreatPropList setCompressionAttributes(const std::size_t length, const int deflateLevel=6)
Sets up the chunking and compression rate.
MANTID_NEXUS_DLL DataType getType< int16_t >()