Files and structure for python interface

This commit is contained in:
Erik Fröjdh
2024-10-28 11:22:12 +01:00
parent 5d643dc133
commit a4fb217e3f
26 changed files with 2238 additions and 43 deletions

View File

@ -4,15 +4,20 @@ set(SourceFiles
${CMAKE_CURRENT_SOURCE_DIR}/defs.cpp
${CMAKE_CURRENT_SOURCE_DIR}/Dtype.cpp
${CMAKE_CURRENT_SOURCE_DIR}/Frame.cpp
${CMAKE_CURRENT_SOURCE_DIR}/File.cpp
${CMAKE_CURRENT_SOURCE_DIR}/NumpyFile.cpp
${CMAKE_CURRENT_SOURCE_DIR}/RawFile.cpp
${CMAKE_CURRENT_SOURCE_DIR}/SubFile.cpp
${CMAKE_CURRENT_SOURCE_DIR}/NumpyHelpers.cpp
)
add_library(aare_core STATIC ${SourceFiles})
target_include_directories(aare_core PUBLIC ${CMAKE_SOURCE_DIR}/include)
target_link_libraries(aare_core PUBLIC fmt::fmt PRIVATE aare_compiler_flags )
target_link_libraries(aare_core PUBLIC fmt::fmt PRIVATE aare_compiler_flags nlohmann_json::nlohmann_json)
if (AARE_PYTHON_BINDINGS)
set_property(TARGET aare_core PROPERTY POSITION_INDEPENDENT_CODE ON)
set_property(TARGET aare_core PROPERTY POSITION_INDEPENDENT_CODE ON)
endif()
if(AARE_TESTS)
@ -25,11 +30,13 @@ if(AARE_TESTS)
${CMAKE_CURRENT_SOURCE_DIR}/NDView.test.cpp
${CMAKE_CURRENT_SOURCE_DIR}/ClusterFinder.test.cpp
${CMAKE_CURRENT_SOURCE_DIR}/Pedestal.test.cpp
${CMAKE_CURRENT_SOURCE_DIR}/NumpyFile.test.cpp
${CMAKE_CURRENT_SOURCE_DIR}/NumpyHelpers.test.cpp
${CMAKE_CURRENT_SOURCE_DIR}/RawFile.test.cpp
# ${CMAKE_CURRENT_SOURCE_DIR}/test/CircularFifo.test.cpp
# ${CMAKE_CURRENT_SOURCE_DIR}/test/wrappers.test.cpp
# ${CMAKE_CURRENT_SOURCE_DIR}/test/Transforms.test.cpp
)
target_sources(tests PRIVATE ${TestSources} )
target_link_libraries(tests PRIVATE aare_core)
endif()

70
src/File.cpp Normal file
View File

@ -0,0 +1,70 @@
#include "aare/File.hpp"
#include "aare/NumpyFile.hpp"
#include "aare/RawFile.hpp"
#include <fmt/format.h>
namespace aare {
File::File(const std::filesystem::path &fname, const std::string &mode, const FileConfig &cfg)
: file_impl(nullptr), is_npy(true) {
if (mode != "r" && mode != "w" && mode != "a") {
throw std::invalid_argument("Unsupported file mode");
}
if ((mode == "r" || mode == "a") && !std::filesystem::exists(fname)) {
throw std::runtime_error(fmt::format("File does not exist: {}", fname.string()));
}
if (fname.extension() == ".raw" || fname.extension() == ".json") {
// aare::logger::debug("Loading raw file");
file_impl = new RawFile(fname, mode, cfg);
is_npy = false;
}
// check if extension is numpy
else if (fname.extension() == ".npy") {
// aare::logger::debug("Loading numpy file");
file_impl = new NumpyFile(fname, mode, cfg);
} else {
throw std::runtime_error("Unsupported file type");
}
}
void File::write(Frame &frame, sls_detector_header header) {
if (is_npy) {
// aare::logger::info("ignoring header for npy file");
dynamic_cast<NumpyFile *>(file_impl)->write(frame);
} else {
dynamic_cast<RawFile *>(file_impl)->write(frame, header);
}
}
Frame File::read() { return file_impl->read(); }
size_t File::total_frames() const { return file_impl->total_frames(); }
std::vector<Frame> File::read(size_t n_frames) { return file_impl->read(n_frames); }
void File::read_into(std::byte *image_buf) { file_impl->read_into(image_buf); }
void File::read_into(std::byte *image_buf, size_t n_frames) { file_impl->read_into(image_buf, n_frames); }
size_t File::frame_number(size_t frame_index) { return file_impl->frame_number(frame_index); }
size_t File::bytes_per_frame() { return file_impl->bytes_per_frame(); }
size_t File::pixels_per_frame() { return file_impl->pixels_per_frame(); }
void File::seek(size_t frame_number) { file_impl->seek(frame_number); }
size_t File::tell() const { return file_impl->tell(); }
size_t File::rows() const { return file_impl->rows(); }
size_t File::cols() const { return file_impl->cols(); }
size_t File::bitdepth() const { return file_impl->bitdepth(); }
void File::set_total_frames(size_t total_frames) { return file_impl->set_total_frames(total_frames); }
File::~File() { delete file_impl; }
DetectorType File::detector_type() const { return file_impl->detector_type(); }
xy File::geometry() const {
if (is_npy) {
return {1, 1};
}
return reinterpret_cast<RawFile *>(file_impl)->geometry();
}
Frame File::iread(size_t frame_number) { return file_impl->iread(frame_number); }
File::File(File &&other) noexcept : file_impl(other.file_impl), is_npy(other.is_npy) { other.file_impl = nullptr; }
// write move assignment operator
} // namespace aare

200
src/NumpyFile.cpp Normal file
View File

@ -0,0 +1,200 @@
#include "aare/NumpyFile.hpp"
#include "aare/NumpyHelpers.hpp"
namespace aare {
NumpyFile::NumpyFile(const std::filesystem::path &fname, const std::string &mode, FileConfig cfg) {
// TODO! add opts to constructor
m_fname = fname;
m_mode = mode;
if (mode == "r") {
fp = fopen(m_fname.string().c_str(), "rb");
if (!fp) {
throw std::runtime_error(fmt::format("Could not open: {} for reading", m_fname.string()));
}
load_metadata();
} else if (mode == "w") {
m_bitdepth = cfg.dtype.bitdepth();
m_rows = cfg.rows;
m_cols = cfg.cols;
m_header = {cfg.dtype, false, {cfg.rows, cfg.cols}};
m_header.shape = {0, cfg.rows, cfg.cols};
fp = fopen(m_fname.string().c_str(), "wb");
if (!fp) {
throw std::runtime_error(fmt::format("Could not open: {} for reading", m_fname.string()));
}
initial_header_len = aare::NumpyHelpers::write_header(std::filesystem::path(m_fname.c_str()), m_header);
}
m_pixels_per_frame = std::accumulate(m_header.shape.begin() + 1, m_header.shape.end(), 1, std::multiplies<>());
m_bytes_per_frame = m_header.dtype.bitdepth() / 8 * m_pixels_per_frame;
}
void NumpyFile::write(Frame &frame) { write_impl(frame.data(), frame.bytes()); }
void NumpyFile::write_impl(void *data, uint64_t size) {
if (fp == nullptr) {
throw std::runtime_error("File not open");
}
if (!(m_mode == "w" || m_mode == "a")) {
throw std::invalid_argument("File not open for writing");
}
if (fseek(fp, 0, SEEK_END))
throw std::runtime_error("Could not seek to end of file");
size_t const rc = fwrite(data, size, 1, fp);
if (rc != 1) {
throw std::runtime_error("Error writing frame to file");
}
m_header.shape[0]++;
}
Frame NumpyFile::get_frame(size_t frame_number) {
Frame frame(m_header.shape[1], m_header.shape[2], m_header.dtype);
get_frame_into(frame_number, frame.data());
return frame;
}
void NumpyFile::get_frame_into(size_t frame_number, std::byte *image_buf) {
if (fp == nullptr) {
throw std::runtime_error("File not open");
}
if (frame_number > m_header.shape[0]) {
throw std::invalid_argument("Frame number out of range");
}
if (fseek(fp, header_size + frame_number * m_bytes_per_frame, SEEK_SET)) // NOLINT
throw std::runtime_error("Could not seek to frame");
size_t const rc = fread(image_buf, m_bytes_per_frame, 1, fp);
if (rc != 1) {
throw std::runtime_error("Error reading frame from file");
}
}
size_t NumpyFile::pixels_per_frame() { return m_pixels_per_frame; };
size_t NumpyFile::bytes_per_frame() { return m_bytes_per_frame; };
std::vector<Frame> NumpyFile::read(size_t n_frames) {
// TODO: implement this in a more efficient way
std::vector<Frame> frames;
for (size_t i = 0; i < n_frames; i++) {
frames.push_back(get_frame(current_frame));
current_frame++;
}
return frames;
}
void NumpyFile::read_into(std::byte *image_buf, size_t n_frames) {
// TODO: implement this in a more efficient way
for (size_t i = 0; i < n_frames; i++) {
get_frame_into(current_frame++, image_buf);
image_buf += m_bytes_per_frame;
}
}
NumpyFile::~NumpyFile() noexcept {
if (m_mode == "w" || m_mode == "a") {
// determine number of frames
if (fseek(fp, 0, SEEK_END)) {
std::cout << "Could not seek to end of file" << std::endl;
}
size_t const file_size = ftell(fp);
size_t const data_size = file_size - initial_header_len;
size_t const n_frames = data_size / m_bytes_per_frame;
// update number of frames in header (first element of shape)
m_header.shape[0] = n_frames;
if (fseek(fp, 0, SEEK_SET)) {
std::cout << "Could not seek to beginning of file" << std::endl;
}
// create string stream to contain header
std::stringstream ss;
aare::NumpyHelpers::write_header(ss, m_header);
std::string const header_str = ss.str();
// write header
size_t const rc = fwrite(header_str.c_str(), header_str.size(), 1, fp);
if (rc != 1) {
std::cout << "Error writing header to numpy file in destructor" << std::endl;
}
}
if (fp != nullptr) {
if (fclose(fp)) {
std::cout << "Error closing file" << std::endl;
}
}
}
void NumpyFile::load_metadata() {
// read magic number
std::array<char, 6> tmp{};
size_t rc = fread(tmp.data(), tmp.size(), 1, fp);
if (rc != 1) {
throw std::runtime_error("Error reading magic number");
}
if (tmp != aare::NumpyHelpers::magic_str) {
for (auto item : tmp)
fmt::print("{}, ", static_cast<int>(item));
fmt::print("\n");
throw std::runtime_error("Not a numpy file");
}
// read version
rc = fread(reinterpret_cast<char *>(&major_ver_), sizeof(major_ver_), 1, fp);
rc += fread(reinterpret_cast<char *>(&minor_ver_), sizeof(minor_ver_), 1, fp);
if (rc != 2) {
throw std::runtime_error("Error reading numpy version");
}
if (major_ver_ == 1) {
header_len_size = 2;
} else if (major_ver_ == 2) {
header_len_size = 4;
} else {
throw std::runtime_error("Unsupported numpy version");
}
// read header length
rc = fread(reinterpret_cast<char *>(&header_len), header_len_size, 1, fp);
if (rc != 1) {
throw std::runtime_error("Error reading header length");
}
header_size = aare::NumpyHelpers::magic_string_length + 2 + header_len_size + header_len;
if (header_size % 16 != 0) {
fmt::print("Warning: header length is not a multiple of 16\n");
}
// read header
std::string header(header_len, '\0');
rc = fread(header.data(), header_len, 1, fp);
if (rc != 1) {
throw std::runtime_error("Error reading header");
}
// parse header
std::vector<std::string> const keys{"descr", "fortran_order", "shape"};
auto dict_map = aare::NumpyHelpers::parse_dict(header, keys);
if (dict_map.empty())
throw std::runtime_error("invalid dictionary in header");
std::string const descr_s = dict_map["descr"];
std::string const fortran_s = dict_map["fortran_order"];
std::string const shape_s = dict_map["shape"];
std::string const descr = aare::NumpyHelpers::parse_str(descr_s);
aare::Dtype const dtype = aare::NumpyHelpers::parse_descr(descr);
// convert literal Python bool to C++ bool
bool const fortran_order = aare::NumpyHelpers::parse_bool(fortran_s);
// parse the shape tuple
auto shape_v = aare::NumpyHelpers::parse_tuple(shape_s);
std::vector<size_t> shape;
for (const auto &item : shape_v) {
auto dim = static_cast<size_t>(std::stoul(item));
shape.push_back(dim);
}
m_header = {dtype, fortran_order, shape};
}
} // namespace aare

50
src/NumpyFile.test.cpp Normal file
View File

@ -0,0 +1,50 @@
#include "aare/NumpyFile.hpp"
#include "aare/NDArray.hpp"
#include <catch2/catch_test_macros.hpp>
#include "test_config.hpp"
using aare::Dtype;
using aare::NumpyFile;
TEST_CASE("Read a 1D numpy file with int32 data type") {
auto fpath = test_data_path() / "numpy" / "test_1d_int32.npy";
REQUIRE(std::filesystem::exists(fpath));
NumpyFile f(fpath);
// we know the file contains 10 elements of np.int32 containing values 0-9
REQUIRE(f.dtype() == Dtype::INT32);
REQUIRE(f.shape() == std::vector<size_t>{10});
// use the load function to read the full file into a NDArray
auto data = f.load<int32_t, 1>();
for (int32_t i = 0; i < 10; i++) {
REQUIRE(data(i) == i);
}
}
TEST_CASE("Read a 3D numpy file with np.double data type") {
auto fpath = test_data_path() / "numpy" / "test_3d_double.npy";
REQUIRE(std::filesystem::exists(fpath));
NumpyFile f(fpath);
// we know the file contains 10 elements of np.int32 containing values 0-9
REQUIRE(f.dtype() == Dtype::DOUBLE);
REQUIRE(f.shape() == std::vector<size_t>{3, 2, 5});
// use the load function to read the full file into a NDArray
// numpy code to generate the array
// arr2[0,0,0] = 1.0
// arr2[0,0,1] = 2.0
// arr2[0,1,0] = 72.0
// arr2[2,0,4] = 63.0
auto data = f.load<double, 3>();
REQUIRE(data(0, 0, 0) == 1.0);
REQUIRE(data(0, 0, 1) == 2.0);
REQUIRE(data(0, 1, 0) == 72.0);
REQUIRE(data(2, 0, 4) == 63.0);
}

269
src/NumpyHelpers.cpp Normal file
View File

@ -0,0 +1,269 @@
/*
28-03-2024 modified by: Bechir Braham <bechir.braham@psi.ch>
Copyright 2017-2023 Leon Merten Lohse
Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files (the "Software"), to deal
in the Software without restriction, including without limitation the rights
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in
all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
SOFTWARE.
*/
#include "aare/NumpyHelpers.hpp"
#include <iterator>
namespace aare {
std::string NumpyHeader::to_string() const {
std::stringstream sstm;
sstm << "dtype: " << dtype.to_string() << ", fortran_order: " << fortran_order << ' ';
sstm << "shape: (";
for (auto item : shape)
sstm << item << ',';
sstm << ')';
return sstm.str();
}
namespace NumpyHelpers {
std::unordered_map<std::string, std::string> parse_dict(std::string in, const std::vector<std::string> &keys) {
std::unordered_map<std::string, std::string> map;
if (keys.empty())
return map;
in = trim(in);
// unwrap dictionary
if ((in.front() == '{') && (in.back() == '}'))
in = in.substr(1, in.length() - 2);
else
throw std::runtime_error("Not a Python dictionary.");
std::vector<std::pair<size_t, std::string>> positions;
for (auto const &key : keys) {
size_t const pos = in.find("'" + key + "'");
if (pos == std::string::npos)
throw std::runtime_error("Missing '" + key + "' key.");
std::pair<size_t, std::string> const position_pair{pos, key};
positions.push_back(position_pair);
}
// sort by position in dict
std::sort(positions.begin(), positions.end());
for (size_t i = 0; i < positions.size(); ++i) {
std::string raw_value;
size_t const begin{positions[i].first};
size_t end{std::string::npos};
std::string const key = positions[i].second;
if (i + 1 < positions.size())
end = positions[i + 1].first;
raw_value = in.substr(begin, end - begin);
raw_value = trim(raw_value);
if (raw_value.back() == ',')
raw_value.pop_back();
map[key] = get_value_from_map(raw_value);
}
return map;
}
aare::Dtype parse_descr(std::string typestring) {
if (typestring.length() < 3) {
throw std::runtime_error("invalid typestring (length)");
}
constexpr char little_endian_char = '<';
constexpr char big_endian_char = '>';
constexpr char no_endian_char = '|';
constexpr std::array<char, 3> endian_chars = {little_endian_char, big_endian_char, no_endian_char};
constexpr std::array<char, 4> numtype_chars = {'f', 'i', 'u', 'c'};
const char byteorder_c = typestring[0];
const char kind_c = typestring[1];
std::string const itemsize_s = typestring.substr(2);
if (!in_array(byteorder_c, endian_chars)) {
throw std::runtime_error("invalid typestring (byteorder)");
}
if (!in_array(kind_c, numtype_chars)) {
throw std::runtime_error("invalid typestring (kind)");
}
if (!is_digits(itemsize_s)) {
throw std::runtime_error("invalid typestring (itemsize)");
}
return aare::Dtype(typestring);
}
bool parse_bool(const std::string &in) {
if (in == "True")
return true;
if (in == "False")
return false;
throw std::runtime_error("Invalid python boolean.");
}
std::string get_value_from_map(const std::string &mapstr) {
size_t const sep_pos = mapstr.find_first_of(':');
if (sep_pos == std::string::npos)
return "";
std::string const tmp = mapstr.substr(sep_pos + 1);
return trim(tmp);
}
bool is_digits(const std::string &str) { return std::all_of(str.begin(), str.end(), ::isdigit); }
std::vector<std::string> parse_tuple(std::string in) {
std::vector<std::string> v;
const char separator = ',';
in = trim(in);
if ((in.front() == '(') && (in.back() == ')'))
in = in.substr(1, in.length() - 2);
else
throw std::runtime_error("Invalid Python tuple.");
std::istringstream iss(in);
for (std::string token; std::getline(iss, token, separator);) {
v.push_back(token);
}
return v;
}
std::string trim(const std::string &str) {
const std::string whitespace = " \t\n";
auto begin = str.find_first_not_of(whitespace);
if (begin == std::string::npos)
return "";
auto end = str.find_last_not_of(whitespace);
return str.substr(begin, end - begin + 1);
}
std::string parse_str(const std::string &in) {
if ((in.front() == '\'') && (in.back() == '\''))
return in.substr(1, in.length() - 2);
throw std::runtime_error("Invalid python string.");
}
void write_magic(std::ostream &ostream, int version_major, int version_minor) {
ostream.write(magic_str.data(), magic_string_length);
ostream.put(static_cast<char>(version_major));
ostream.put(static_cast<char>(version_minor));
}
template <typename T> inline std::string write_tuple(const std::vector<T> &v) {
if (v.empty())
return "()";
std::ostringstream ss;
ss.imbue(std::locale("C"));
if (v.size() == 1) {
ss << "(" << v.front() << ",)";
} else {
const std::string delimiter = ", ";
// v.size() > 1
ss << "(";
// for (size_t i = 0; i < v.size() - 1; ++i) {
// ss << v[i] << delimiter;
// }
// ss << v.back();
std::copy(v.begin(), v.end() - 1, std::ostream_iterator<T>(ss, ", "));
ss << v.back();
ss << ")";
}
return ss.str();
}
inline std::string write_boolean(bool b) {
if (b)
return "True";
return "False";
}
inline std::string write_header_dict(const std::string &descr, bool fortran_order, const std::vector<size_t> &shape) {
std::string const s_fortran_order = write_boolean(fortran_order);
std::string const shape_s = write_tuple(shape);
return "{'descr': '" + descr + "', 'fortran_order': " + s_fortran_order + ", 'shape': " + shape_s + ", }";
}
size_t write_header(const std::filesystem::path &fname, const NumpyHeader &header) {
std::ofstream out(fname, std::ios::binary | std::ios::out);
return write_header(out, header);
}
size_t write_header(std::ostream &out, const NumpyHeader &header) {
std::string const header_dict = write_header_dict(header.dtype.to_string(), header.fortran_order, header.shape);
size_t length = magic_string_length + 2 + 2 + header_dict.length() + 1;
int version_major = 1;
int version_minor = 0;
if (length >= static_cast<size_t>(255) * 255) {
length = magic_string_length + 2 + 4 + header_dict.length() + 1;
version_major = 2;
version_minor = 0;
}
size_t const padding_len = 16 - length % 16;
std::string const padding(padding_len, ' ');
// write magic
write_magic(out, version_major, version_minor);
// write header length
if (version_major == 1 && version_minor == 0) {
auto header_len = static_cast<uint16_t>(header_dict.length() + padding.length() + 1);
std::array<uint8_t, 2> header_len_le16{static_cast<uint8_t>((header_len >> 0) & 0xff),
static_cast<uint8_t>((header_len >> 8) & 0xff)};
out.write(reinterpret_cast<char *>(header_len_le16.data()), 2);
} else {
auto header_len = static_cast<uint32_t>(header_dict.length() + padding.length() + 1);
std::array<uint8_t, 4> header_len_le32{
static_cast<uint8_t>((header_len >> 0) & 0xff), static_cast<uint8_t>((header_len >> 8) & 0xff),
static_cast<uint8_t>((header_len >> 16) & 0xff), static_cast<uint8_t>((header_len >> 24) & 0xff)};
out.write(reinterpret_cast<char *>(header_len_le32.data()), 4);
}
out << header_dict << padding << '\n';
return length;
}
} // namespace NumpyHelpers
} // namespace aare

62
src/NumpyHelpers.test.cpp Normal file
View File

@ -0,0 +1,62 @@
#include "aare/NumpyHelpers.hpp" //Is this really a public header?
#include <catch2/catch_test_macros.hpp>
using namespace aare::NumpyHelpers;
TEST_CASE("is_digits with a few standard cases") {
REQUIRE(is_digits(""));
REQUIRE(is_digits("123"));
REQUIRE(is_digits("0"));
REQUIRE_FALSE(is_digits("hej123"));
REQUIRE_FALSE(is_digits("a"));
REQUIRE_FALSE(is_digits(" "));
REQUIRE_FALSE(is_digits("abcdef"));
}
TEST_CASE("Check for quotes and return stripped string") {
REQUIRE(parse_str("'hej'") == "hej");
REQUIRE(parse_str("'hej hej'") == "hej hej");
REQUIRE(parse_str("''") == "");
}
TEST_CASE("parsing a string without quotes throws") { REQUIRE_THROWS(parse_str("hej")); }
TEST_CASE("trim whitespace") {
REQUIRE(trim(" hej ") == "hej");
REQUIRE(trim("hej") == "hej");
REQUIRE(trim(" hej") == "hej");
REQUIRE(trim("hej ") == "hej");
REQUIRE(trim(" ") == "");
REQUIRE(trim(" \thej hej ") == "hej hej");
}
TEST_CASE("parse data type descriptions") {
REQUIRE(parse_descr("<i1") == aare::Dtype::INT8);
REQUIRE(parse_descr("<i2") == aare::Dtype::INT16);
REQUIRE(parse_descr("<i4") == aare::Dtype::INT32);
REQUIRE(parse_descr("<i8") == aare::Dtype::INT64);
REQUIRE(parse_descr("<u1") == aare::Dtype::UINT8);
REQUIRE(parse_descr("<u2") == aare::Dtype::UINT16);
REQUIRE(parse_descr("<u4") == aare::Dtype::UINT32);
REQUIRE(parse_descr("<u8") == aare::Dtype::UINT64);
REQUIRE(parse_descr("<f4") == aare::Dtype::FLOAT);
REQUIRE(parse_descr("<f8") == aare::Dtype::DOUBLE);
}
TEST_CASE("is element in array") {
REQUIRE(in_array(1, std::array<int, 3>{1, 2, 3}));
REQUIRE_FALSE(in_array(4, std::array<int, 3>{1, 2, 3}));
REQUIRE(in_array(1, std::array<int, 1>{1}));
REQUIRE_FALSE(in_array(1, std::array<int, 0>{}));
}
TEST_CASE("Parse numpy dict") {
std::string in = "{'descr': '<f4', 'fortran_order': False, 'shape': (3, 4)}";
std::vector<std::string> keys{"descr", "fortran_order", "shape"};
auto map = parse_dict(in, keys);
REQUIRE(map["descr"] == "'<f4'");
REQUIRE(map["fortran_order"] == "False");
REQUIRE(map["shape"] == "(3, 4)");
}

428
src/RawFile.cpp Normal file
View File

@ -0,0 +1,428 @@
#include "aare/RawFile.hpp"
#include "aare/defs.hpp"
#include "aare/json.hpp"
#include <fmt/format.h>
#include <nlohmann/json.hpp>
using json = nlohmann::json;
namespace aare {
RawFile::RawFile(const std::filesystem::path &fname, const std::string &mode, const FileConfig &config) {
m_mode = mode;
m_fname = fname;
if (mode == "r" || mode == "r+") {
if (config != FileConfig()) {
// aare::logger::warn(
// "In read mode it is not necessary to provide a config, the provided config will be ignored");
}
parse_fname();
parse_metadata();
find_number_of_subfiles();
find_geometry();
open_subfiles();
} else if (mode == "w" || mode == "w+") {
if (std::filesystem::exists(fname)) {
// handle mode w as w+ (no overrwriting)
throw std::runtime_error(LOCATION + "File already exists");
}
parse_config(config);
parse_fname();
write_master_file();
n_subfiles = 1;
n_subfile_parts = 1;
subfile_cols = m_cols;
subfile_rows = m_rows;
open_subfiles();
} else {
throw std::runtime_error(LOCATION + "Unsupported mode");
}
}
void RawFile::parse_config(const FileConfig &config) {
m_bitdepth = config.dtype.bitdepth();
m_total_frames = config.total_frames;
m_rows = config.rows;
m_cols = config.cols;
m_type = config.detector_type;
max_frames_per_file = config.max_frames_per_file;
m_geometry = config.geometry;
version = config.version;
subfile_rows = config.geometry.row;
subfile_cols = config.geometry.col;
if (m_geometry != aare::xy{1, 1}) {
throw std::runtime_error(LOCATION + "Only geometry {1,1} files are supported for writing");
}
}
void RawFile::write_master_file() {
if (m_ext != ".json") {
throw std::runtime_error(LOCATION + "only json master files are supported for writing");
}
std::ofstream ofs(master_fname(), std::ios::binary);
std::string ss;
ss.reserve(1024);
ss += "{\n\t";
aare::write_str(ss, "Version", version);
ss += "\n\t";
aare::write_digit(ss, "Total Frames", m_total_frames);
ss += "\n\t";
aare::write_str(ss, "Detector Type", toString(m_type));
ss += "\n\t";
aare::write_str(ss, "Geometry", m_geometry.to_string());
ss += "\n\t";
uint64_t img_size = (m_cols * m_rows) / (static_cast<size_t>(m_geometry.col * m_geometry.row));
img_size *= m_bitdepth;
aare::write_digit(ss, "Image Size in bytes", img_size);
ss += "\n\t";
aare::write_digit(ss, "Max Frames Per File", max_frames_per_file);
ss += "\n\t";
aare::write_digit(ss, "Dynamic Range", m_bitdepth);
ss += "\n\t";
const aare::xy pixels = {static_cast<uint32_t>(m_rows / m_geometry.row),
static_cast<uint32_t>(m_cols / m_geometry.col)};
aare::write_str(ss, "Pixels", pixels.to_string());
ss += "\n\t";
aare::write_digit(ss, "Number of rows", m_rows);
ss += "\n\t";
const std::string tmp = "{\n"
" \"Frame Number\": \"8 bytes\",\n"
" \"Exposure Length\": \"4 bytes\",\n"
" \"Packet Number\": \"4 bytes\",\n"
" \"Bunch Id\": \"8 bytes\",\n"
" \"Timestamp\": \"8 bytes\",\n"
" \"Module Id\": \"2 bytes\",\n"
" \"Row\": \"2 bytes\",\n"
" \"Column\": \"2 bytes\",\n"
" \"Reserved\": \"2 bytes\",\n"
" \"Debug\": \"4 bytes\",\n"
" \"RoundRNumber\": \"2 bytes\",\n"
" \"DetType\": \"1 byte\",\n"
" \"Version\": \"1 byte\",\n"
" \"Packet Mask\": \"64 bytes\"\n"
" }";
ss += "\"Frame Header Format\":" + tmp + "\n";
ss += "}";
ofs << ss;
ofs.close();
}
void RawFile::open_subfiles() {
if (m_mode == "r")
for (size_t i = 0; i != n_subfiles; ++i) {
auto v = std::vector<SubFile *>(n_subfile_parts);
for (size_t j = 0; j != n_subfile_parts; ++j) {
v[j] = new SubFile(data_fname(i, j), m_type, subfile_rows, subfile_cols, m_bitdepth);
}
subfiles.push_back(v);
}
else {
auto v = std::vector<SubFile *>(n_subfile_parts); // only one subfile is implemented
v[0] = new SubFile(data_fname(0, 0), m_type, m_rows, m_cols, m_bitdepth, "w");
subfiles.push_back(v);
}
}
sls_detector_header RawFile::read_header(const std::filesystem::path &fname) {
sls_detector_header h{};
FILE *fp = fopen(fname.string().c_str(), "r");
if (!fp)
throw std::runtime_error(fmt::format("Could not open: {} for reading", fname.string()));
size_t const rc = fread(reinterpret_cast<char *>(&h), sizeof(h), 1, fp);
if (rc != 1)
throw std::runtime_error(LOCATION + "Could not read header from file");
if (fclose(fp)) {
throw std::runtime_error(LOCATION + "Could not close file");
}
return h;
}
bool RawFile::is_master_file(const std::filesystem::path &fpath) {
std::string const stem = fpath.stem().string();
return stem.find("_master_") != std::string::npos;
}
void RawFile::find_number_of_subfiles() {
int n_mod = 0;
while (std::filesystem::exists(data_fname(++n_mod, 0)))
;
n_subfiles = n_mod;
}
inline std::filesystem::path RawFile::data_fname(size_t mod_id, size_t file_id) {
return this->m_base_path / fmt::format("{}_d{}_f{}_{}.raw", this->m_base_name, file_id, mod_id, this->m_findex);
}
inline std::filesystem::path RawFile::master_fname() {
return this->m_base_path / fmt::format("{}_master_{}{}", this->m_base_name, this->m_findex, this->m_ext);
}
void RawFile::find_geometry() {
uint16_t r{};
uint16_t c{};
for (size_t i = 0; i < n_subfile_parts; i++) {
for (size_t j = 0; j != n_subfiles; ++j) {
auto h = this->read_header(data_fname(j, i));
r = std::max(r, h.row);
c = std::max(c, h.column);
positions.push_back({h.row, h.column});
}
}
r++;
c++;
m_rows = (r * subfile_rows);
m_cols = (c * subfile_cols);
m_rows += static_cast<size_t>((r - 1) * cfg.module_gap_row);
}
void RawFile::parse_metadata() {
if (m_ext == ".raw") {
parse_raw_metadata();
if (m_bitdepth == 0) {
switch (m_type) {
case DetectorType::Eiger:
m_bitdepth = 32;
break;
default:
m_bitdepth = 16;
}
}
} else if (m_ext == ".json") {
parse_json_metadata();
} else {
throw std::runtime_error(LOCATION + "Unsupported file type");
}
n_subfile_parts = static_cast<size_t>(m_geometry.row) * m_geometry.col;
}
void RawFile::parse_json_metadata() {
std::ifstream ifs(master_fname());
json j;
ifs >> j;
double v = j["Version"];
version = fmt::format("{:.1f}", v);
m_type = StringTo<DetectorType>(j["Detector Type"].get<std::string>());
timing_mode = StringTo<TimingMode>(j["Timing Mode"].get<std::string>());
m_total_frames = j["Frames in File"];
subfile_rows = j["Pixels"]["y"];
subfile_cols = j["Pixels"]["x"];
max_frames_per_file = j["Max Frames Per File"];
try {
m_bitdepth = j.at("Dynamic Range");
} catch (const json::out_of_range &e) {
m_bitdepth = 16;
}
// only Eiger had quad
if (m_type == DetectorType::Eiger) {
quad = (j["Quad"] == 1);
}
m_geometry = {j["Geometry"]["y"], j["Geometry"]["x"]};
}
void RawFile::parse_raw_metadata() {
std::ifstream ifs(master_fname());
for (std::string line; std::getline(ifs, line);) {
if (line == "#Frame Header")
break;
auto pos = line.find(':');
auto key_pos = pos;
while (key_pos != std::string::npos && std::isspace(line[--key_pos]))
;
if (key_pos != std::string::npos) {
auto key = line.substr(0, key_pos + 1);
auto value = line.substr(pos + 2);
// do the actual parsing
if (key == "Version") {
version = value;
} else if (key == "TimeStamp") {
} else if (key == "Detector Type") {
m_type = StringTo<DetectorType>(value);
} else if (key == "Timing Mode") {
timing_mode = StringTo<TimingMode>(value);
} else if (key == "Pixels") {
// Total number of pixels cannot be found yet looking at
// submodule
pos = value.find(',');
subfile_cols = std::stoi(value.substr(1, pos));
subfile_rows = std::stoi(value.substr(pos + 1));
} else if (key == "Total Frames") {
m_total_frames = std::stoi(value);
} else if (key == "Dynamic Range") {
m_bitdepth = std::stoi(value);
} else if (key == "Quad") {
quad = (value == "1");
} else if (key == "Max Frames Per File") {
max_frames_per_file = std::stoi(value);
} else if (key == "Geometry") {
pos = value.find(',');
m_geometry = {static_cast<uint32_t>(std::stoi(value.substr(1, pos))),
static_cast<uint32_t>(std::stoi(value.substr(pos + 1)))};
}
}
}
}
void RawFile::parse_fname() {
bool wrong_format = false;
m_base_path = m_fname.parent_path().string();
m_base_name = m_fname.stem().string();
m_ext = m_fname.extension().string();
try {
auto pos = m_base_name.rfind('_');
m_findex = std::stoi(m_base_name.substr(pos + 1));
} catch (const std::invalid_argument &e) {
m_findex = 0;
wrong_format = true;
}
auto pos = m_base_name.find("_master_");
if (pos != std::string::npos) {
m_base_name.erase(pos);
wrong_format = true;
}
if (wrong_format && (m_mode == "w+" || m_mode == "w")) {
// aare::logger::warn("Master Filename", m_fname, "is not in the correct format");
// aare::logger::warn("using", master_fname(), "as the master file");
}
}
Frame RawFile::get_frame(size_t frame_index) {
auto f = Frame(this->m_rows, this->m_cols, Dtype::from_bitdepth(this->m_bitdepth));
std::byte *frame_buffer = f.data();
get_frame_into(frame_index, frame_buffer);
return f;
}
void RawFile::get_frame_into(size_t frame_index, std::byte *frame_buffer) {
if (frame_index > this->m_total_frames) {
throw std::runtime_error(LOCATION + "Frame number out of range");
}
std::vector<size_t> frame_numbers(this->n_subfile_parts);
std::vector<size_t> frame_indices(this->n_subfile_parts, frame_index);
if (n_subfile_parts != 1) {
for (size_t part_idx = 0; part_idx != this->n_subfile_parts; ++part_idx) {
auto subfile_id = frame_index / this->max_frames_per_file;
frame_numbers[part_idx] =
this->subfiles[subfile_id][part_idx]->frame_number(frame_index % this->max_frames_per_file);
}
// 1. if frame number vector is the same break
while (std::adjacent_find(frame_numbers.begin(), frame_numbers.end(), std::not_equal_to<>()) !=
frame_numbers.end()) {
// 2. find the index of the minimum frame number,
auto min_frame_idx =
std::distance(frame_numbers.begin(), std::min_element(frame_numbers.begin(), frame_numbers.end()));
// 3. increase its index and update its respective frame number
frame_indices[min_frame_idx]++;
// 4. if we can't increase its index => throw error
if (frame_indices[min_frame_idx] >= this->m_total_frames) {
throw std::runtime_error(LOCATION + "Frame number out of range");
}
auto subfile_id = frame_indices[min_frame_idx] / this->max_frames_per_file;
frame_numbers[min_frame_idx] = this->subfiles[subfile_id][min_frame_idx]->frame_number(
frame_indices[min_frame_idx] % this->max_frames_per_file);
}
}
if (this->m_geometry.col == 1) {
// get the part from each subfile and copy it to the frame
for (size_t part_idx = 0; part_idx != this->n_subfile_parts; ++part_idx) {
auto corrected_idx = frame_indices[part_idx];
auto subfile_id = corrected_idx / this->max_frames_per_file;
auto part_offset = this->subfiles[subfile_id][part_idx]->bytes_per_part();
this->subfiles[subfile_id][part_idx]->get_part(frame_buffer + part_idx * part_offset,
corrected_idx % this->max_frames_per_file);
}
} else {
// create a buffer that will hold a the frame part
auto bytes_per_part = this->subfile_rows * this->subfile_cols * this->m_bitdepth / 8;
auto *part_buffer = new std::byte[bytes_per_part];
for (size_t part_idx = 0; part_idx != this->n_subfile_parts; ++part_idx) {
auto corrected_idx = frame_indices[part_idx];
auto subfile_id = corrected_idx / this->max_frames_per_file;
this->subfiles[subfile_id][part_idx]->get_part(part_buffer, corrected_idx % this->max_frames_per_file);
for (size_t cur_row = 0; cur_row < (this->subfile_rows); cur_row++) {
auto irow = cur_row + (part_idx / this->m_geometry.col) * this->subfile_rows;
auto icol = (part_idx % this->m_geometry.col) * this->subfile_cols;
auto dest = (irow * this->m_cols + icol);
dest = dest * this->m_bitdepth / 8;
memcpy(frame_buffer + dest, part_buffer + cur_row * this->subfile_cols * this->m_bitdepth / 8,
this->subfile_cols * this->m_bitdepth / 8);
}
}
delete[] part_buffer;
}
}
void RawFile::write(Frame &frame, sls_detector_header header) {
if (m_mode == "r") {
throw std::runtime_error(LOCATION + "File is open in read mode");
}
size_t const subfile_id = this->current_frame / this->max_frames_per_file;
for (size_t part_idx = 0; part_idx != this->n_subfile_parts; ++part_idx) {
this->subfiles[subfile_id][part_idx]->write_part(frame.data(), header,
this->current_frame % this->max_frames_per_file);
}
this->current_frame++;
}
std::vector<Frame> RawFile::read(size_t n_frames) {
// TODO: implement this in a more efficient way
std::vector<Frame> frames;
for (size_t i = 0; i < n_frames; i++) {
frames.push_back(this->get_frame(this->current_frame));
this->current_frame++;
}
return frames;
}
void RawFile::read_into(std::byte *image_buf, size_t n_frames) {
// TODO: implement this in a more efficient way
for (size_t i = 0; i < n_frames; i++) {
this->get_frame_into(this->current_frame++, image_buf);
image_buf += this->bytes_per_frame();
}
}
size_t RawFile::frame_number(size_t frame_index) {
if (frame_index > this->m_total_frames) {
throw std::runtime_error(LOCATION + "Frame number out of range");
}
size_t const subfile_id = frame_index / this->max_frames_per_file;
return this->subfiles[subfile_id][0]->frame_number(frame_index % this->max_frames_per_file);
}
RawFile::~RawFile() noexcept {
// update master file
if (m_mode == "w" || m_mode == "w+" || m_mode == "r+") {
try {
write_master_file();
} catch (...) {
// aare::logger::warn(LOCATION + "Could not update master file");
}
}
for (auto &vec : subfiles) {
for (auto *subfile : vec) {
delete subfile;
}
}
}
} // namespace aare

114
src/RawFile.test.cpp Normal file
View File

@ -0,0 +1,114 @@
#include "aare/File.hpp"
#include <catch2/catch_test_macros.hpp>
#include <filesystem>
#include "test_config.hpp"
using aare::File;
TEST_CASE("Read number of frames from a jungfrau raw file") {
auto fpath = test_data_path() / "jungfrau" / "jungfrau_single_master_0.json";
REQUIRE(std::filesystem::exists(fpath));
File f(fpath, "r");
REQUIRE(f.total_frames() == 10);
}
TEST_CASE("Read frame numbers from a jungfrau raw file") {
auto fpath = test_data_path() / "jungfrau" / "jungfrau_single_master_0.json";
REQUIRE(std::filesystem::exists(fpath));
File f(fpath, "r");
// we know this file has 10 frames with frame numbers 1 to 10
// f0 1,2,3
// f1 4,5,6
// f2 7,8,9
// f3 10
for (size_t i = 0; i < 10; i++) {
CHECK(f.frame_number(i) == i + 1);
}
}
TEST_CASE("Read data from a jungfrau 500k single port raw file") {
auto fpath = test_data_path() / "jungfrau" / "jungfrau_single_master_0.json";
REQUIRE(std::filesystem::exists(fpath));
File f(fpath, "r");
// we know this file has 10 frames with pixel 0,0 being: 2123, 2051, 2109, 2117, 2089, 2095, 2072, 2126, 2097, 2102
std::vector<uint16_t> pixel_0_0 = {2123, 2051, 2109, 2117, 2089, 2095, 2072, 2126, 2097, 2102};
for (size_t i = 0; i < 10; i++) {
auto frame = f.read();
CHECK(frame.rows() == 512);
CHECK(frame.cols() == 1024);
CHECK(frame.view<uint16_t>()(0, 0) == pixel_0_0[i]);
}
}
TEST_CASE("Read frame numbers from a raw file") {
auto fpath = test_data_path() / "eiger" / "eiger_500k_16bit_master_0.json";
REQUIRE(std::filesystem::exists(fpath));
// we know this file has 3 frames with frame numbers 14, 15, 16
std::vector<size_t> frame_numbers = {14, 15, 16};
File f(fpath, "r");
for (size_t i = 0; i < 3; i++) {
CHECK(f.frame_number(i) == frame_numbers[i]);
}
}
TEST_CASE("Compare reading from a numpy file with a raw file") {
auto fpath_raw = test_data_path() / "jungfrau" / "jungfrau_single_master_0.json";
REQUIRE(std::filesystem::exists(fpath_raw));
auto fpath_npy = test_data_path() / "jungfrau" / "jungfrau_single_0.npy";
REQUIRE(std::filesystem::exists(fpath_npy));
File raw(fpath_raw, "r");
File npy(fpath_npy, "r");
CHECK(raw.total_frames() == 10);
CHECK(npy.total_frames() == 10);
for (size_t i = 0; i < 10; ++i) {
auto raw_frame = raw.read();
auto npy_frame = npy.read();
CHECK((raw_frame.view<uint16_t>() == npy_frame.view<uint16_t>()));
}
}
TEST_CASE("Read multipart files") {
auto fpath = test_data_path() / "jungfrau" / "jungfrau_double_master_0.json";
REQUIRE(std::filesystem::exists(fpath));
File f(fpath, "r");
// we know this file has 10 frames check read_multiport.py for the values
std::vector<uint16_t> pixel_0_0 = {2099, 2121, 2108, 2084, 2084, 2118, 2066, 2108, 2112, 2116};
std::vector<uint16_t> pixel_0_1 = {2842, 2796, 2865, 2798, 2805, 2817, 2852, 2789, 2792, 2833};
std::vector<uint16_t> pixel_255_1023 = {2149, 2037, 2115, 2102, 2118, 2090, 2036, 2071, 2073, 2142};
std::vector<uint16_t> pixel_511_1023 = {3231, 3169, 3167, 3162, 3168, 3160, 3171, 3171, 3169, 3171};
std::vector<uint16_t> pixel_1_0 = {2748, 2614, 2665, 2629, 2618, 2630, 2631, 2634, 2577, 2598};
for (size_t i = 0; i < 10; i++) {
auto frame = f.read();
CHECK(frame.rows() == 512);
CHECK(frame.cols() == 1024);
CHECK(frame.view<uint16_t>()(0, 0) == pixel_0_0[i]);
CHECK(frame.view<uint16_t>()(0, 1) == pixel_0_1[i]);
CHECK(frame.view<uint16_t>()(1, 0) == pixel_1_0[i]);
CHECK(frame.view<uint16_t>()(255, 1023) == pixel_255_1023[i]);
CHECK(frame.view<uint16_t>()(511, 1023) == pixel_511_1023[i]);
}
}
TEST_CASE("Read file with unordered frames") {
auto fpath = test_data_path() / "mythen" / "scan242_master_3.raw";
REQUIRE(std::filesystem::exists(fpath));
File f(fpath, "r");
REQUIRE_THROWS((f.read()));
}

71
src/SubFile.cpp Normal file
View File

@ -0,0 +1,71 @@
#include "aare/SubFile.hpp"
#include <cstring> // memcpy
#include <fmt/core.h>
#include <iostream>
namespace aare {
SubFile::SubFile(const std::filesystem::path &fname, DetectorType detector, size_t rows, size_t cols, size_t bitdepth,
const std::string &mode)
: m_bitdepth(bitdepth), m_fname(fname), m_rows(rows), m_cols(cols), m_mode(mode) {
if (std::filesystem::exists(fname)) {
n_frames = std::filesystem::file_size(fname) / (sizeof(sls_detector_header) + rows * cols * bitdepth / 8);
} else {
n_frames = 0;
}
if (mode == "r") {
fp = fopen(m_fname.string().c_str(), "rb");
} else {
// if file exists, open in read/write mode (without truncating the file)
// if file does not exist, open in write mode
if (std::filesystem::exists(fname)) {
fp = fopen(m_fname.string().c_str(), "r+b");
} else {
fp = fopen(m_fname.string().c_str(), "wb");
}
}
if (fp == nullptr) {
throw std::runtime_error(LOCATION + "Could not open file for writing");
}
}
size_t SubFile::get_part(std::byte *buffer, size_t frame_index) {
if (frame_index >= n_frames) {
throw std::runtime_error("Frame number out of range");
}
fseek(fp, sizeof(sls_detector_header) + (sizeof(sls_detector_header) + bytes_per_part()) * frame_index, // NOLINT
SEEK_SET);
return fread(buffer, this->bytes_per_part(), 1, this->fp);
}
size_t SubFile::write_part(std::byte *buffer, sls_detector_header header, size_t frame_index) {
if (frame_index > n_frames) {
throw std::runtime_error("Frame number out of range");
}
fseek(fp, static_cast<int64_t>((sizeof(sls_detector_header) + bytes_per_part()) * frame_index), SEEK_SET);
auto wc = fwrite(reinterpret_cast<char *>(&header), sizeof(header), 1, fp);
wc += fwrite(buffer, bytes_per_part(), 1, fp);
return wc;
}
size_t SubFile::frame_number(size_t frame_index) {
sls_detector_header h{};
fseek(fp, (sizeof(sls_detector_header) + bytes_per_part()) * frame_index, SEEK_SET); // NOLINT
size_t const rc = fread(reinterpret_cast<char *>(&h), sizeof(h), 1, fp);
if (rc != 1)
throw std::runtime_error(LOCATION + "Could not read header from file");
return h.frameNumber;
}
SubFile::~SubFile() {
if (fp) {
fclose(fp);
}
}
} // namespace aare