merged with developer

.gitea/workflows/rh8-native.yml

@@ -1,18 +1,24 @@
 name: Build on RHEL8
 
 on:
+  push:
   workflow_dispatch:
 
 permissions:
   contents: read
 
 jobs:
-  buildh:
+  build:
     runs-on: "ubuntu-latest"
     container:
       image: gitea.psi.ch/images/rhel8-developer-gitea-actions
     steps:
-      - uses: actions/checkout@v4
+      # workaround until actions/checkout@v4 is available for RH8
+      # - uses: actions/checkout@v4
+      - name: Clone repository
+        run: |
+          echo Cloning ${{ github.ref_name }}
+          git clone https://${{secrets.GITHUB_TOKEN}}@gitea.psi.ch/${{ github.repository }}.git --branch=${{ github.ref_name }} .
 
 
       - name: Install dependencies
@@ -22,7 +28,7 @@ jobs:
       - name: Build library 
         run: |
             mkdir build && cd build
-            cmake .. -DAARE_PYTHON_BINDINGS=ON -DAARE_TESTS=ON
+            cmake .. -DAARE_PYTHON_BINDINGS=ON -DAARE_TESTS=ON -DPython_FIND_VIRTUALENV=FIRST
             make -j 2 
 
       - name: C++ unit tests

.gitea/workflows/rh9-native.yml

@@ -8,7 +8,7 @@ permissions:
   contents: read
 
 jobs:
-  buildh:
+  build:
     runs-on: "ubuntu-latest"
     container:
       image: gitea.psi.ch/images/rhel9-developer-gitea-actions

.github/workflows/build_wheel.yml

@@ -0,0 +1,64 @@
+name: Build wheel
+
+on:
+    workflow_dispatch:
+    pull_request:
+    push:
+      branches:
+        - main
+    release:
+      types:
+        - published
+  
+
+jobs:
+  build_wheels:
+    name: Build wheels on ${{ matrix.os }}
+    runs-on: ${{ matrix.os }}
+    strategy:
+      matrix:
+        os: [ubuntu-latest,]
+
+    steps:
+      - uses: actions/checkout@v4
+
+      - name: Build wheels
+        run: pipx run cibuildwheel==2.23.0
+
+      - uses: actions/upload-artifact@v4
+        with:
+          name: cibw-wheels-${{ matrix.os }}-${{ strategy.job-index }}
+          path: ./wheelhouse/*.whl
+    
+  build_sdist:
+    name: Build source distribution
+    runs-on: ubuntu-latest
+    steps:
+    - uses: actions/checkout@v4
+
+    - name: Build sdist
+      run: pipx run build --sdist
+
+    - uses: actions/upload-artifact@v4
+      with:
+        name: cibw-sdist
+        path: dist/*.tar.gz
+
+  upload_pypi:
+    needs: [build_wheels, build_sdist]
+    runs-on: ubuntu-latest
+    environment: pypi
+    permissions:
+        id-token: write
+    if: github.event_name == 'release' && github.event.action == 'published'
+    # or, alternatively, upload to PyPI on every tag starting with 'v' (remove on: release above to use this)
+    # if: github.event_name == 'push' && startsWith(github.ref, 'refs/tags/v')
+    steps:
+        - uses: actions/download-artifact@v4
+          with:
+            # unpacks all CIBW artifacts into dist/
+            pattern: cibw-*
+            path: dist
+            merge-multiple: true
+
+        - uses: pypa/gh-action-pypi-publish@release/v1

.gitignore

@@ -17,7 +17,8 @@ Testing/
 ctbDict.cpp
 ctbDict.h
 
-
+wheelhouse/
+dist/
 
 *.pyc
 */__pycache__/*

CMakeLists.txt

@@ -1,4 +1,4 @@
-cmake_minimum_required(VERSION 3.14)
+cmake_minimum_required(VERSION 3.15)
 
 project(aare 
     VERSION 1.0.0
@@ -11,6 +11,14 @@ set(CMAKE_CXX_STANDARD 17)
 set(CMAKE_CXX_STANDARD_REQUIRED ON)
 set(CMAKE_CXX_EXTENSIONS OFF)
 
+execute_process(
+        COMMAND git log -1 --format=%h
+        WORKING_DIRECTORY ${CMAKE_CURRENT_LIST_DIR}
+        OUTPUT_VARIABLE GIT_HASH
+        OUTPUT_STRIP_TRAILING_WHITESPACE
+        )
+message(STATUS "Building from git hash: ${GIT_HASH}")
+
 if (${CMAKE_VERSION} VERSION_GREATER "3.24")
     cmake_policy(SET CMP0135 NEW) #Fetch content download timestamp
 endif()
@@ -381,7 +389,9 @@ set(SourceFiles
     ${CMAKE_CURRENT_SOURCE_DIR}/src/RawFile.cpp
     ${CMAKE_CURRENT_SOURCE_DIR}/src/RawSubFile.cpp
     ${CMAKE_CURRENT_SOURCE_DIR}/src/RawMasterFile.cpp
+
     ${CMAKE_CURRENT_SOURCE_DIR}/src/utils/task.cpp
+    ${CMAKE_CURRENT_SOURCE_DIR}/src/utils/ifstream_helpers.cpp
 )
 
 add_library(aare_core STATIC ${SourceFiles})
@@ -415,6 +425,7 @@ if(AARE_TESTS)
     set(TestSources
             ${CMAKE_CURRENT_SOURCE_DIR}/src/algorithm.test.cpp
             ${CMAKE_CURRENT_SOURCE_DIR}/src/defs.test.cpp
+            ${CMAKE_CURRENT_SOURCE_DIR}/src/decode.test.cpp
             ${CMAKE_CURRENT_SOURCE_DIR}/src/Dtype.test.cpp
             ${CMAKE_CURRENT_SOURCE_DIR}/src/Frame.test.cpp
             ${CMAKE_CURRENT_SOURCE_DIR}/src/geo_helpers.test.cpp

conda-recipe/meta.yaml

@@ -1,6 +1,7 @@
 package:
   name: aare
-  version: 2025.4.1 #TODO! how to not duplicate this?
+  version: 2025.4.22 #TODO! how to not duplicate this?
+
 
 
 
@@ -38,6 +39,7 @@ requirements:
   run:
     - python {{python}}
     - numpy {{ numpy }}
+    - matplotlib
 
 
 test:

include/aare/ClusterFile.hpp

@@ -160,9 +160,8 @@ class ClusterFile {
     }
 
     /**
-     * @brief Set the gain map to use when reading clusters. If
-     * set the gain map will be applied to the clusters that
-     * pass ROI and noise_map selection.
+     * @brief Set the gain map to use when reading clusters. If set the gain map will be applied
+     * to the clusters that pass ROI and noise_map selection. The gain map is expected to be in ADU/energy.
      */
     void set_gain_map(const NDView<double, 2> gain_map) {
         m_gain_map = GainMap(gain_map);

include/aare/NDView.hpp

@@ -184,4 +184,9 @@ std::ostream& operator <<(std::ostream& os, const NDView<T, Ndim>& arr){
 }
 
 
+template <typename T>
+NDView<T,1> make_view(std::vector<T>& vec){
+    return NDView<T,1>(vec.data(), {static_cast<int64_t>(vec.size())});
+}
+
 } // namespace aare

include/aare/RawSubFile.hpp

@@ -22,7 +22,7 @@ class RawSubFile {
     size_t m_rows{};
     size_t m_cols{};
     size_t m_bytes_per_frame{};
-    size_t n_frames{};
+    size_t m_num_frames{};
     uint32_t m_pos_row{};
     uint32_t m_pos_col{};
  
@@ -53,6 +53,7 @@ class RawSubFile {
     size_t tell();
 
     void read_into(std::byte *image_buf, DetectorHeader *header = nullptr);
+    void read_into(std::byte *image_buf, size_t n_frames, DetectorHeader *header= nullptr);
     void get_part(std::byte *buffer, size_t frame_index);
     
     void read_header(DetectorHeader *header);
@@ -66,6 +67,8 @@ class RawSubFile {
     size_t pixels_per_frame() const { return m_rows * m_cols; }
     size_t bytes_per_pixel() const { return m_bitdepth / bits_per_byte; }
 
+    size_t frames_in_file() const { return m_num_frames; }
+
 private:
   template <typename T>
   void read_with_map(std::byte *image_buf);

include/aare/decode.hpp

@@ -1,6 +1,7 @@
 #pragma once
 
 #include <cstdint>
+#include <vector>
 #include <aare/NDView.hpp>
 namespace aare {
 
@@ -10,4 +11,16 @@ uint16_t adc_sar_04_decode64to16(uint64_t input);
 void adc_sar_05_decode64to16(NDView<uint64_t, 2> input, NDView<uint16_t,2> output);
 void adc_sar_04_decode64to16(NDView<uint64_t, 2> input, NDView<uint16_t,2> output);
 
-} // namespace aare
+
+/**
+ * @brief Apply custom weights to a 16-bit input value. Will sum up weights[i]**i
+ * for each bit i that is set in the input value.
+ * @throws std::out_of_range if weights.size() < 16
+ * @param input 16-bit input value
+ * @param weights vector of weights, size must be less than or equal to 16
+ */
+double apply_custom_weights(uint16_t input, const NDView<double, 1> weights);
+
+void apply_custom_weights(NDView<uint16_t, 1> input, NDView<double, 1> output, const NDView<double, 1> weights);
+
+} // namespace aare

include/aare/utils/ifstream_helpers.hpp

@@ -0,0 +1,12 @@
+#pragma once
+
+#include <fstream>
+#include <string>
+namespace aare {
+
+/** 
+ * @brief Get the error message from an ifstream object
+*/
+std::string ifstream_error_msg(std::ifstream &ifs);
+
+} // namespace aare

pyproject.toml

@@ -4,19 +4,31 @@ build-backend = "scikit_build_core.build"
 
 [project]
 name = "aare"
-version = "2025.4.1"
+version = "2025.4.22"
+requires-python = ">=3.11"
+dependencies = [
+    "numpy",
+    "matplotlib",
+]
+
+
+[tool.cibuildwheel]
+
+build = "cp{311,312,313}-manylinux_x86_64"
 
 
 
 
 [tool.scikit-build]
-cmake.verbose = true
+build.verbose = true
+cmake.build-type = "Release"
+install.components = ["python"]
 
 [tool.scikit-build.cmake.define]
 AARE_PYTHON_BINDINGS = "ON"
-AARE_SYSTEM_LIBRARIES = "ON"
 AARE_INSTALL_PYTHONEXT = "ON"
 
+
 [tool.pytest.ini_options]
 markers = [
     "files: marks tests that need additional data (deselect with '-m \"not files\"')",

python/CMakeLists.txt

@@ -1,12 +1,13 @@
 
-find_package (Python 3.10 COMPONENTS Interpreter Development REQUIRED)
+find_package (Python 3.10 COMPONENTS Interpreter Development.Module REQUIRED)
+set(PYBIND11_FINDPYTHON ON) # Needed for RH8
 
 # Download or find pybind11 depending on configuration
 if(AARE_FETCH_PYBIND11)
     FetchContent_Declare(
         pybind11
         GIT_REPOSITORY https://github.com/pybind/pybind11
-        GIT_TAG        v2.13.0
+        GIT_TAG        v2.13.6
     )
     FetchContent_MakeAvailable(pybind11)
 else()
@@ -62,10 +63,16 @@ endforeach(FILE ${PYTHON_EXAMPLES})
 
 
 if(AARE_INSTALL_PYTHONEXT)
-    install(TARGETS _aare
+    install(
+        TARGETS _aare
         EXPORT "${TARGETS_EXPORT_NAME}"
         LIBRARY DESTINATION aare
+        COMPONENT python
     )
 
-    install(FILES ${PYTHON_FILES} DESTINATION aare)
+    install(
+        FILES ${PYTHON_FILES} 
+        DESTINATION aare
+        COMPONENT python
+        )
 endif()

python/aare/__init__.py

@@ -18,6 +18,10 @@ from .ClusterVector import ClusterVector
 from ._aare import fit_gaus, fit_pol1
 from ._aare import Interpolator
 from ._aare import calculate_eta2
+
+
+from ._aare import apply_custom_weights
+
 from .CtbRawFile import CtbRawFile
 from .RawFile import RawFile
 from .ScanParameters import ScanParameters

python/src/ctb_raw_file.hpp

@@ -10,6 +10,8 @@
 #include "aare/decode.hpp"
 // #include "aare/fClusterFileV2.hpp"
 
+#include "np_helper.hpp"
+
 #include <cstdint>
 #include <filesystem>
 #include <pybind11/iostream.h>
@@ -65,35 +67,54 @@ m.def("adc_sar_04_decode64to16", [](py::array_t<uint8_t> input) {
     return output;
 });
 
- py::class_<CtbRawFile>(m, "CtbRawFile")
-        .def(py::init<const std::filesystem::path &>())
-        .def("read_frame",
-             [](CtbRawFile &self) {
-                 size_t image_size = self.image_size_in_bytes();
-                 py::array image;
-                 std::vector<ssize_t> shape;
-                 shape.reserve(2);
-                 shape.push_back(1);
-                 shape.push_back(image_size);
+m.def(
+    "apply_custom_weights",
+    [](py::array_t<uint16_t, py::array::c_style | py::array::forcecast> &input,
+       py::array_t<double, py::array::c_style | py::array::forcecast>
+           &weights) {
+        
 
-                 py::array_t<DetectorHeader> header(1);
+        // Create new array with same shape as the input array (uninitialized values)
+        py::buffer_info buf = input.request();
+        py::array_t<double> output(buf.shape);
 
-                 // always read bytes
-                 image = py::array_t<uint8_t>(shape);
+        // Use NDViews to call into the C++ library
+        auto weights_view = make_view_1d(weights);
+        NDView<uint16_t, 1> input_view(input.mutable_data(), {input.size()});
+        NDView<double, 1> output_view(output.mutable_data(), {output.size()});
 
-                 self.read_into(
-                     reinterpret_cast<std::byte *>(image.mutable_data()),
-                     header.mutable_data());
+        apply_custom_weights(input_view, output_view, weights_view);
+        return output;
+    });
 
-                 return py::make_tuple(header, image);
-             })
-        .def("seek", &CtbRawFile::seek)
-        .def("tell", &CtbRawFile::tell)
-        .def("master", &CtbRawFile::master)
+py::class_<CtbRawFile>(m, "CtbRawFile")
+    .def(py::init<const std::filesystem::path &>())
+    .def("read_frame",
+         [](CtbRawFile &self) {
+             size_t image_size = self.image_size_in_bytes();
+             py::array image;
+             std::vector<ssize_t> shape;
+             shape.reserve(2);
+             shape.push_back(1);
+             shape.push_back(image_size);
 
-        .def_property_readonly("image_size_in_bytes",
-                               &CtbRawFile::image_size_in_bytes)
+             py::array_t<DetectorHeader> header(1);
 
-        .def_property_readonly("frames_in_file", &CtbRawFile::frames_in_file);
+             // always read bytes
+             image = py::array_t<uint8_t>(shape);
 
-}
+             self.read_into(reinterpret_cast<std::byte *>(image.mutable_data()),
+                            header.mutable_data());
+
+             return py::make_tuple(header, image);
+         })
+    .def("seek", &CtbRawFile::seek)
+    .def("tell", &CtbRawFile::tell)
+    .def("master", &CtbRawFile::master)
+
+    .def_property_readonly("image_size_in_bytes",
+                           &CtbRawFile::image_size_in_bytes)
+
+    .def_property_readonly("frames_in_file", &CtbRawFile::frames_in_file);
+
+}

python/src/file.hpp

@@ -20,6 +20,9 @@
 namespace py = pybind11;
 using namespace ::aare;
 
+
+
+
 //Disable warnings for unused parameters, as we ignore some
 //in the __exit__ method
 #pragma GCC diagnostic push
@@ -214,36 +217,9 @@ void define_file_io_bindings(py::module &m) {
     
 
 
-    py::class_<RawSubFile>(m, "RawSubFile")
-        .def(py::init<const std::filesystem::path &, DetectorType, size_t,
-                      size_t, size_t>())
-        .def_property_readonly("bytes_per_frame", &RawSubFile::bytes_per_frame)
-        .def_property_readonly("pixels_per_frame",
-                               &RawSubFile::pixels_per_frame)
-        .def("seek", &RawSubFile::seek)
-        .def("tell", &RawSubFile::tell)
-        .def_property_readonly("rows", &RawSubFile::rows)
-        .def_property_readonly("cols", &RawSubFile::cols)
-        .def("read_frame",
-             [](RawSubFile &self) {
-                 const uint8_t item_size = self.bytes_per_pixel();
-                 py::array image;
-                 std::vector<ssize_t> shape;
-                 shape.reserve(2);
-                 shape.push_back(self.rows());
-                 shape.push_back(self.cols());
-                 if (item_size == 1) {
-                     image = py::array_t<uint8_t>(shape);
-                 } else if (item_size == 2) {
-                     image = py::array_t<uint16_t>(shape);
-                 } else if (item_size == 4) {
-                     image = py::array_t<uint32_t>(shape);
-                 }
-                 fmt::print("item_size: {} rows: {} cols: {}\n", item_size, self.rows(), self.cols());
-                 self.read_into(
-                     reinterpret_cast<std::byte *>(image.mutable_data()));
-                 return image;
-             });
+
+
+
 
 #pragma GCC diagnostic pop
     // py::class_<ClusterHeader>(m, "ClusterHeader")

python/src/module.cpp

@@ -15,6 +15,7 @@
 #include "raw_file.hpp"
 #include "raw_master_file.hpp"
 #include "var_cluster.hpp"
+#include "raw_sub_file.hpp"
 
 #include "jungfrau_data_file.hpp"
 
@@ -27,6 +28,7 @@ namespace py = pybind11;
 PYBIND11_MODULE(_aare, m) {
     define_file_io_bindings(m);
     define_raw_file_io_bindings(m);
+    define_raw_sub_file_io_bindings(m);
     define_ctb_raw_file_io_bindings(m);
     define_raw_master_file_bindings(m);
     define_var_cluster_finder_bindings(m);

python/src/raw_sub_file.hpp

@@ -0,0 +1,110 @@
+#include "aare/CtbRawFile.hpp"
+#include "aare/File.hpp"
+#include "aare/Frame.hpp"
+#include "aare/RawFile.hpp"
+#include "aare/RawMasterFile.hpp"
+#include "aare/RawSubFile.hpp"
+
+#include "aare/defs.hpp"
+// #include "aare/fClusterFileV2.hpp"
+
+#include <cstdint>
+#include <filesystem>
+#include <pybind11/iostream.h>
+#include <pybind11/numpy.h>
+#include <pybind11/pybind11.h>
+#include <pybind11/stl.h>
+#include <pybind11/stl/filesystem.h>
+#include <string>
+
+namespace py = pybind11;
+using namespace ::aare;
+
+auto read_frame_from_RawSubFile(RawSubFile &self) {
+    py::array_t<DetectorHeader> header(1);
+    const uint8_t item_size = self.bytes_per_pixel();
+    std::vector<ssize_t> shape{static_cast<ssize_t>(self.rows()),
+        static_cast<ssize_t>(self.cols())};
+    
+    py::array image;
+    if (item_size == 1) {
+        image = py::array_t<uint8_t>(shape);
+    } else if (item_size == 2) {
+        image = py::array_t<uint16_t>(shape);
+    } else if (item_size == 4) {
+        image = py::array_t<uint32_t>(shape);
+    }
+    self.read_into(reinterpret_cast<std::byte *>(image.mutable_data()),
+                   header.mutable_data());
+
+    return py::make_tuple(header, image);
+}
+
+auto read_n_frames_from_RawSubFile(RawSubFile &self, size_t n_frames) {
+    py::array_t<DetectorHeader> header(n_frames);
+    const uint8_t item_size = self.bytes_per_pixel();
+    std::vector<ssize_t> shape{
+        static_cast<ssize_t>(n_frames), 
+        static_cast<ssize_t>(self.rows()),
+        static_cast<ssize_t>(self.cols())
+    };
+    
+    py::array image;
+    if (item_size == 1) {
+        image = py::array_t<uint8_t>(shape);
+    } else if (item_size == 2) {
+        image = py::array_t<uint16_t>(shape);
+    } else if (item_size == 4) {
+        image = py::array_t<uint32_t>(shape);
+    }
+    self.read_into(reinterpret_cast<std::byte *>(image.mutable_data()), n_frames,
+                   header.mutable_data());
+
+    return py::make_tuple(header, image);
+}
+
+
+//Disable warnings for unused parameters, as we ignore some
+//in the __exit__ method
+#pragma GCC diagnostic push
+#pragma GCC diagnostic ignored "-Wunused-parameter"
+
+void define_raw_sub_file_io_bindings(py::module &m) {
+    py::class_<RawSubFile>(m, "RawSubFile")
+        .def(py::init<const std::filesystem::path &, DetectorType, size_t,
+                      size_t, size_t>())
+        .def_property_readonly("bytes_per_frame", &RawSubFile::bytes_per_frame)
+        .def_property_readonly("pixels_per_frame",
+                               &RawSubFile::pixels_per_frame)
+        .def_property_readonly("bytes_per_pixel", &RawSubFile::bytes_per_pixel)                    
+        .def("seek", &RawSubFile::seek)
+        .def("tell", &RawSubFile::tell)
+        .def_property_readonly("rows", &RawSubFile::rows)
+        .def_property_readonly("cols", &RawSubFile::cols)
+        .def_property_readonly("frames_in_file", &RawSubFile::frames_in_file)
+        .def("read_frame", &read_frame_from_RawSubFile)
+        .def("read_n", &read_n_frames_from_RawSubFile)
+        .def("read", [](RawSubFile &self){
+            self.seek(0);
+            auto n_frames = self.frames_in_file();
+            return read_n_frames_from_RawSubFile(self, n_frames);
+        })
+        .def("__enter__", [](RawSubFile &self) { return &self; })
+        .def("__exit__",
+             [](RawSubFile &self,
+                const std::optional<pybind11::type> &exc_type,
+                const std::optional<pybind11::object> &exc_value,
+                const std::optional<pybind11::object> &traceback) {
+             })
+        .def("__iter__", [](RawSubFile &self) { return &self; })
+        .def("__next__", [](RawSubFile &self) {
+            try {
+                return read_frame_from_RawSubFile(self);
+            } catch (std::runtime_error &e) {
+                throw py::stop_iteration();
+            }
+        });
+
+}
+
+#pragma GCC diagnostic pop

python/tests/test_RawSubFile.py

@@ -0,0 +1,36 @@
+import pytest
+import numpy as np
+from aare import RawSubFile, DetectorType
+
+
+@pytest.mark.files
+def test_read_a_jungfrau_RawSubFile(test_data_path):
+    with RawSubFile(test_data_path / "raw/jungfrau/jungfrau_single_d0_f1_0.raw", DetectorType.Jungfrau, 512, 1024, 16) as f:
+        assert f.frames_in_file == 3
+
+        headers, frames = f.read()
+
+    assert headers.size == 3
+    assert frames.shape == (3, 512, 1024)
+    
+    # Frame numbers in this file should be 4, 5, 6
+    for i,h in zip(range(4,7,1), headers):
+        assert h["frameNumber"] == i
+
+    # Compare to canned data using numpy
+    data = np.load(test_data_path / "raw/jungfrau/jungfrau_single_0.npy")
+    assert np.all(data[3:6] == frames)
+
+@pytest.mark.files
+def test_iterate_over_a_jungfrau_RawSubFile(test_data_path):
+
+    data = np.load(test_data_path / "raw/jungfrau/jungfrau_single_0.npy")
+
+    with RawSubFile(test_data_path / "raw/jungfrau/jungfrau_single_d0_f0_0.raw", DetectorType.Jungfrau, 512, 1024, 16) as f:
+        i = 0
+        for header, frame in f:
+            assert header["frameNumber"] == i+1
+            assert np.all(frame == data[i])
+            i += 1
+        assert i == 3
+        assert header["frameNumber"] == 3

src/ClusterFile.test.cpp

@@ -27,6 +27,7 @@ TEST_CASE("Read one frame from a cluster file", "[.files]") {
                      std::begin(expected_cluster_data)));
 }
 
+
 TEST_CASE("Read one frame using ROI", "[.files]") {
     // We know that the frame has 97 clusters
     auto fpath = test_data_path() / "clust" / "single_frame_97_clustrers.clust";
@@ -60,6 +61,8 @@ TEST_CASE("Read one frame using ROI", "[.files]") {
                      std::begin(expected_cluster_data)));
 }
 
+
+
 TEST_CASE("Read clusters from single frame file", "[.files]") {
 
     //  frame_number, num_clusters   [135] 97
@@ -162,6 +165,7 @@ TEST_CASE("Read clusters from single frame file", "[.files]") {
     // [ 97 296] [864 865 866 867 868 869 870 871 872]
 
     auto fpath = test_data_path() / "clust" / "single_frame_97_clustrers.clust";
+
     REQUIRE(std::filesystem::exists(fpath));
 
     SECTION("Read fewer clusters than available") {

src/NDView.test.cpp

@@ -3,6 +3,7 @@
 
 #include <iostream>
 #include <vector>
+#include <numeric>
 
 using aare::NDView;
 using aare::Shape;
@@ -21,10 +22,8 @@ TEST_CASE("Element reference 1D") {
 }
 
 TEST_CASE("Element reference 2D") {
-    std::vector<int> vec;
-    for (int i = 0; i != 12; ++i) {
-        vec.push_back(i);
-    }
+    std::vector<int> vec(12);
+    std::iota(vec.begin(), vec.end(), 0);
 
     NDView<int, 2> data(vec.data(), Shape<2>{3, 4});
     REQUIRE(vec.size() == static_cast<size_t>(data.size()));
@@ -58,10 +57,8 @@ TEST_CASE("Element reference 3D") {
 }
 
 TEST_CASE("Plus and miuns with single value") {
-    std::vector<int> vec;
-    for (int i = 0; i != 12; ++i) {
-        vec.push_back(i);
-    }
+    std::vector<int> vec(12);
+    std::iota(vec.begin(), vec.end(), 0);
     NDView<int, 2> data(vec.data(), Shape<2>{3, 4});
     data += 5;
     int i = 0;
@@ -116,10 +113,8 @@ TEST_CASE("elementwise assign") {
 }
 
 TEST_CASE("iterators") {
-    std::vector<int> vec;
-    for (int i = 0; i != 12; ++i) {
-        vec.push_back(i);
-    }
+    std::vector<int> vec(12);
+    std::iota(vec.begin(), vec.end(), 0);
     NDView<int, 1> data(vec.data(), Shape<1>{12});
     int i = 0;
     for (const auto item : data) {
@@ -167,27 +162,31 @@ TEST_CASE("divide with another span") {
 }
 
 TEST_CASE("Retrieve shape") {
-    std::vector<int> vec;
-    for (int i = 0; i != 12; ++i) {
-        vec.push_back(i);
-    }
+    std::vector<int> vec(12);
+    std::iota(vec.begin(), vec.end(), 0);
     NDView<int, 2> data(vec.data(), Shape<2>{3, 4});
     REQUIRE(data.shape()[0] == 3);
     REQUIRE(data.shape()[1] == 4);
 }
 
 TEST_CASE("compare two views") {
-    std::vector<int> vec1;
-    for (int i = 0; i != 12; ++i) {
-        vec1.push_back(i);
-    }
+    std::vector<int> vec1(12);
+    std::iota(vec1.begin(), vec1.end(), 0);
     NDView<int, 2> view1(vec1.data(), Shape<2>{3, 4});
 
-    std::vector<int> vec2;
-    for (int i = 0; i != 12; ++i) {
-        vec2.push_back(i);
-    }
+    std::vector<int> vec2(12);
+    std::iota(vec2.begin(), vec2.end(), 0);
     NDView<int, 2> view2(vec2.data(), Shape<2>{3, 4});
 
     REQUIRE((view1 == view2));
+}
+
+
+TEST_CASE("Create a view over a vector"){
+    std::vector<int> vec(12);
+    std::iota(vec.begin(), vec.end(), 0);
+    auto v = aare::make_view(vec);
+    REQUIRE(v.shape()[0] == 12);
+    REQUIRE(v[0] == 0);
+    REQUIRE(v[11] == 11);
 }

src/RawSubFile.cpp

@@ -1,9 +1,12 @@
 #include "aare/RawSubFile.hpp"
 #include "aare/PixelMap.hpp"
+#include "aare/utils/ifstream_helpers.hpp"
 #include <cstring> // memcpy
 #include <fmt/core.h>
 #include <iostream>
 
+
+
 namespace aare {
 
 RawSubFile::RawSubFile(const std::filesystem::path &fname,
@@ -20,7 +23,7 @@ RawSubFile::RawSubFile(const std::filesystem::path &fname,
     }
 
     if (std::filesystem::exists(fname)) {
-        n_frames = std::filesystem::file_size(fname) /
+        m_num_frames = std::filesystem::file_size(fname) /
                    (sizeof(DetectorHeader) + rows * cols * bitdepth / 8);
     } else {
         throw std::runtime_error(
@@ -35,7 +38,7 @@ RawSubFile::RawSubFile(const std::filesystem::path &fname,
     }
 
 #ifdef AARE_VERBOSE
-    fmt::print("Opened file: {} with {} frames\n", m_fname.string(), n_frames);
+    fmt::print("Opened file: {} with {} frames\n", m_fname.string(), m_num_frames);
     fmt::print("m_rows: {}, m_cols: {}, m_bitdepth: {}\n", m_rows, m_cols,
                m_bitdepth);
     fmt::print("file size: {}\n", std::filesystem::file_size(fname));
@@ -43,8 +46,8 @@ RawSubFile::RawSubFile(const std::filesystem::path &fname,
 }
 
 void RawSubFile::seek(size_t frame_index) {
-    if (frame_index >= n_frames) {
-        throw std::runtime_error(LOCATION + fmt::format("Frame index {} out of range in a file with {} frames", frame_index, n_frames));
+    if (frame_index >= m_num_frames) {
+        throw std::runtime_error(LOCATION + fmt::format("Frame index {} out of range in a file with {} frames", frame_index, m_num_frames));
     }
     m_file.seekg((sizeof(DetectorHeader) + bytes_per_frame()) * frame_index);
 }
@@ -60,6 +63,10 @@ void RawSubFile::read_into(std::byte *image_buf, DetectorHeader *header) {
         m_file.seekg(sizeof(DetectorHeader), std::ios::cur);
     }
 
+    if (m_file.fail()){
+        throw std::runtime_error(LOCATION + ifstream_error_msg(m_file));
+    }
+
     // TODO! expand support for different bitdepths
     if (m_pixel_map) {
         // read into a temporary buffer and then copy the data to the buffer
@@ -79,8 +86,24 @@ void RawSubFile::read_into(std::byte *image_buf, DetectorHeader *header) {
         // read directly into the buffer
         m_file.read(reinterpret_cast<char *>(image_buf), bytes_per_frame());
     }
+
+    if (m_file.fail()){
+        throw std::runtime_error(LOCATION + ifstream_error_msg(m_file));
+    }
 }
 
+void RawSubFile::read_into(std::byte *image_buf, size_t n_frames, DetectorHeader *header) {
+    for (size_t i = 0; i < n_frames; i++) {
+        read_into(image_buf, header);
+        image_buf += bytes_per_frame();
+        if (header) {
+            ++header;
+        }
+    }
+}
+
+
+
 template <typename T>
 void RawSubFile::read_with_map(std::byte *image_buf) {
     auto part_buffer = new std::byte[bytes_per_frame()];

src/decode.cpp

@@ -1,5 +1,5 @@
 #include "aare/decode.hpp"
-
+#include <cmath>
 namespace aare {
 
 uint16_t adc_sar_05_decode64to16(uint64_t input){
@@ -22,6 +22,10 @@ uint16_t adc_sar_05_decode64to16(uint64_t input){
 }
 
 void adc_sar_05_decode64to16(NDView<uint64_t, 2> input, NDView<uint16_t,2> output){
+    if(input.shape() != output.shape()){
+        throw std::invalid_argument(LOCATION + " input and output shapes must match");
+    }
+
     for(int64_t i = 0; i < input.shape(0); i++){
         for(int64_t j = 0; j < input.shape(1); j++){
             output(i,j) = adc_sar_05_decode64to16(input(i,j));
@@ -49,6 +53,9 @@ uint16_t adc_sar_04_decode64to16(uint64_t input){
 }
 
 void adc_sar_04_decode64to16(NDView<uint64_t, 2> input, NDView<uint16_t,2> output){
+    if(input.shape() != output.shape()){
+        throw std::invalid_argument(LOCATION + " input and output shapes must match");
+    }
     for(int64_t i = 0; i < input.shape(0); i++){
         for(int64_t j = 0; j < input.shape(1); j++){
             output(i,j) = adc_sar_04_decode64to16(input(i,j));
@@ -56,6 +63,40 @@ void adc_sar_04_decode64to16(NDView<uint64_t, 2> input, NDView<uint16_t,2> outpu
     }
 }
 
+double apply_custom_weights(uint16_t input, const NDView<double, 1> weights) {
+    if(weights.size() > 16){
+        throw std::invalid_argument("weights size must be less than or equal to 16");
+    }
+
+    double result = 0.0;
+    for (ssize_t i = 0; i < weights.size(); ++i) {
+        result += ((input >> i) & 1) * std::pow(weights[i], i);
+    }
+    return result;
+
+}
+
+void apply_custom_weights(NDView<uint16_t, 1> input, NDView<double, 1> output, const NDView<double,1> weights) {
+    if(input.shape() != output.shape()){
+        throw std::invalid_argument(LOCATION + " input and output shapes must match");
+    }
+
+    //Calculate weights to avoid repeatedly calling std::pow
+    std::vector<double> weights_powers(weights.size());
+    for (ssize_t i = 0; i < weights.size(); ++i) {
+        weights_powers[i] = std::pow(weights[i], i);
+    }
+
+    // Apply custom weights to each element in the input array
+    for (ssize_t i = 0; i < input.shape(0); i++) {            
+        double result = 0.0;
+        for (size_t bit_index = 0; bit_index < weights_powers.size(); ++bit_index) {
+            result += ((input(i) >> bit_index) & 1) * weights_powers[bit_index];
+        }
+        output(i) = result;
+    }
+}
+
 
 
 } // namespace aare

src/decode.test.cpp

@@ -0,0 +1,80 @@
+#include "aare/decode.hpp"
+
+#include <catch2/matchers/catch_matchers_floating_point.hpp>
+#include <catch2/catch_test_macros.hpp>
+#include "aare/NDArray.hpp"
+using Catch::Matchers::WithinAbs;
+#include <vector>
+
+TEST_CASE("test_adc_sar_05_decode64to16"){
+    uint64_t input = 0;
+    uint16_t output = aare::adc_sar_05_decode64to16(input);
+    CHECK(output == 0);
+
+
+    // bit 29 on th input is bit 0 on the output
+    input = 1UL << 29;
+    output = aare::adc_sar_05_decode64to16(input);
+    CHECK(output == 1);
+
+    // test all bits by iteratting through the bitlist
+    std::vector<int> bitlist = {29, 19, 28, 18, 31, 21, 27, 20, 24, 23, 25, 22};
+    for (size_t i = 0; i < bitlist.size(); i++) {
+        input = 1UL << bitlist[i];
+        output = aare::adc_sar_05_decode64to16(input);
+        CHECK(output == (1 << i));
+    }
+
+    
+    // test a few "random" values
+    input = 0;
+    input |= (1UL << 29);
+    input |= (1UL << 19);
+    input |= (1UL << 28);
+    output = aare::adc_sar_05_decode64to16(input);
+    CHECK(output == 7UL);
+
+
+    input = 0;
+    input |= (1UL << 18);
+    input |= (1UL << 27);
+    input |= (1UL << 25);
+    output = aare::adc_sar_05_decode64to16(input);
+    CHECK(output == 1096UL);
+
+    input = 0;
+    input |= (1UL << 25);
+    input |= (1UL << 22);
+    output = aare::adc_sar_05_decode64to16(input);
+    CHECK(output == 3072UL);
+ }
+
+
+ TEST_CASE("test_apply_custom_weights") {
+
+    uint16_t input = 1;
+    aare::NDArray<double, 1> weights_data({3}, 0.0);
+    weights_data(0) = 1.7;
+    weights_data(1) = 2.1;
+    weights_data(2) = 1.8;
+
+    auto weights = weights_data.view();
+
+
+    double output = aare::apply_custom_weights(input, weights);
+    CHECK_THAT(output, WithinAbs(1.0, 0.001));
+
+    input = 1 << 1;
+    output = aare::apply_custom_weights(input, weights);
+    CHECK_THAT(output, WithinAbs(2.1, 0.001));
+
+
+    input = 1 << 2;
+    output = aare::apply_custom_weights(input, weights);
+    CHECK_THAT(output, WithinAbs(3.24, 0.001));
+
+    input = 0b111;
+    output = aare::apply_custom_weights(input, weights);
+    CHECK_THAT(output, WithinAbs(6.34, 0.001));
+
+ }

src/utils/ifstream_helpers.cpp

@@ -0,0 +1,18 @@
+#include "aare/utils/ifstream_helpers.hpp"
+
+namespace aare {
+
+std::string ifstream_error_msg(std::ifstream &ifs) {
+    std::ios_base::iostate state = ifs.rdstate();
+    if (state & std::ios_base::eofbit) {
+        return " End of file reached";
+    } else if (state & std::ios_base::badbit) {
+        return " Bad file stream";
+    } else if (state & std::ios_base::failbit) {
+        return " File read failed";
+    }else{
+        return " Unknown/no error";
+    }
+}
+
+} // namespace aare