detectors/aare
25
0
Fork 0
mirror of https://github.com/slsdetectorgroup/aare.git synced 2025-06-05 04:20:40 +02:00
Code Issues Actions Packages Releases Activity

tests and fix

Browse Source
This commit is contained in:
froejdh_e 2025-04-14 16:38:25 +02:00
parent 3f753ec900
commit 7c93632605
13 changed files with 359 additions and 117 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
include/aare/ClusterFile.hpp
View File

@ -282,7 +282,7 @@ template <typename ClusterType, typename Enable>
ClusterVector<ClusterType>
ClusterFile<ClusterType, Enable>::read_clusters_with_cut(size_t n_clusters) {
ClusterVector<ClusterType> clusters;
clusters.resize(n_clusters);
clusters.reserve(n_clusters);
// if there are photons left from previous frame read them first
if (m_num_left) {

8
include/aare/ClusterVector.hpp
View File

@ -384,6 +384,14 @@ class ClusterVector<Cluster<T, ClusterSizeX, ClusterSizeY, CoordType>> {
void set_frame_number(uint64_t frame_number) {
m_frame_number = frame_number;
}
std::vector<T> sum() {
std::vector<T> sums(m_data.size());
for (size_t i = 0; i < m_data.size(); i++) {
sums[i] = m_data[i].sum();
}
return sums;
}
};
} // namespace aare

4
python/CMakeLists.txt
View File

@ -28,6 +28,9 @@ target_link_libraries(_aare PRIVATE aare_core aare_compiler_flags)
set( PYTHON_FILES
aare/__init__.py
aare/CtbRawFile.py
aare/ClusterFinder.py
aare/ClusterVector.py
aare/func.py
aare/RawFile.py
aare/transform.py
@ -35,6 +38,7 @@ set( PYTHON_FILES
aare/utils.py
)
# Copy the python files to the build directory
foreach(FILE ${PYTHON_FILES})
configure_file(${FILE} ${CMAKE_BINARY_DIR}/${FILE} )

14
python/aare/ClusterFinder.py Normal file
View File

@ -0,0 +1,14 @@
from ._aare import ClusterFinder_Cluster3x3i
import numpy as np
def ClusterFinder(image_size, cluster_size, n_sigma=5, dtype = np.int32, capacity = 1024):
"""
Factory function to create a ClusterFinder object. Provides a cleaner syntax for
the templated ClusterFinder in C++.
"""
if dtype == np.int32 and cluster_size == (3,3):
return ClusterFinder_Cluster3x3i(image_size, n_sigma = n_sigma, capacity=capacity)
else:
#TODO! add the other formats
raise ValueError(f"Unsupported dtype: {dtype}. Only np.int32 is supported.")

11
python/aare/ClusterVector.py Normal file
View File

@ -0,0 +1,11 @@
from ._aare import ClusterVector_Cluster3x3i
import numpy as np
def ClusterVector(cluster_size, dtype = np.int32):
if dtype == np.int32 and cluster_size == (3,3):
return ClusterVector_Cluster3x3i()
else:
raise ValueError(f"Unsupported dtype: {dtype}. Only np.int32 is supported.")

4
python/aare/__init__.py
View File

@ -11,8 +11,12 @@ from ._aare import ROI
# from ._aare import ClusterFinderMT, ClusterCollector, ClusterFileSink, ClusterVector_i
from .ClusterFinder import ClusterFinder
from .ClusterVector import ClusterVector
from ._aare import fit_gaus, fit_pol1
from ._aare import Interpolator
from ._aare import calculate_eta2
from .CtbRawFile import CtbRawFile
from .RawFile import RawFile
from .ScanParameters import ScanParameters

103
python/src/bind_ClusterVector.hpp Normal file
View File

@ -0,0 +1,103 @@
#include "aare/ClusterCollector.hpp"
#include "aare/ClusterFileSink.hpp"
#include "aare/ClusterFinder.hpp"
#include "aare/ClusterFinderMT.hpp"
#include "aare/ClusterVector.hpp"
#include "aare/NDView.hpp"
#include "aare/Pedestal.hpp"
#include "np_helper.hpp"
#include <cstdint>
#include <filesystem>
#include <pybind11/pybind11.h>
#include <pybind11/stl.h>
#include <pybind11/stl_bind.h>
namespace py = pybind11;
using pd_type = double;
using namespace aare;
#pragma GCC diagnostic push
#pragma GCC diagnostic ignored "-Wunused-parameter"
template <typename Type, uint8_t ClusterSizeX, uint8_t ClusterSizeY,
typename CoordType = uint16_t>
void define_ClusterVector(py::module &m, const std::string &typestr) {
using ClusterType =
Cluster<Type, ClusterSizeX, ClusterSizeY, CoordType, void>;
auto class_name = fmt::format("ClusterVector_{}", typestr);
py::class_<ClusterVector<
Cluster<Type, ClusterSizeX, ClusterSizeY, CoordType, void>, void>>(
m, class_name.c_str(),
py::buffer_protocol())
.def(py::init()) // TODO change!!!
.def("push_back",
[](ClusterVector<ClusterType> &self, const ClusterType &cluster) {
self.push_back(cluster);
})
.def("sum", [](ClusterVector<ClusterType> &self) {
auto *vec = new std::vector<Type>(self.sum());
return return_vector(vec);
})
.def_property_readonly("size", &ClusterVector<ClusterType>::size)
.def("item_size", &ClusterVector<ClusterType>::item_size)
.def_property_readonly("fmt",
[typestr](ClusterVector<ClusterType> &self) {
return fmt_format<ClusterType>;
})
.def_property_readonly("cluster_size_x",
&ClusterVector<ClusterType>::cluster_size_x)
.def_property_readonly("cluster_size_y",
&ClusterVector<ClusterType>::cluster_size_y)
.def_property_readonly("capacity",
&ClusterVector<ClusterType>::capacity)
.def_property("frame_number", &ClusterVector<ClusterType>::frame_number,
&ClusterVector<ClusterType>::set_frame_number)
.def_buffer(
[typestr](ClusterVector<ClusterType> &self) -> py::buffer_info {
return py::buffer_info(
self.data(), /* Pointer to buffer */
self.item_size(), /* Size of one scalar */
fmt_format<ClusterType>, /* Format descriptor */
1, /* Number of dimensions */
{self.size()}, /* Buffer dimensions */
{self.item_size()} /* Strides (in bytes) for each index */
);
});
// Free functions using ClusterVector
m.def("hitmap",
[](std::array<size_t, 2> image_size, ClusterVector<ClusterType> &cv) {
// Create a numpy array to hold the hitmap
// The shape of the array is (image_size[0], image_size[1])
// note that the python array is passed as [row, col] which
// is the opposite of the clusters [x,y]
py::array_t<int32_t> hitmap(image_size);
auto r = hitmap.mutable_unchecked<2>();
// Initialize hitmap to 0
for (py::ssize_t i = 0; i < r.shape(0); i++)
for (py::ssize_t j = 0; j < r.shape(1); j++)
r(i, j) = 0;
// Loop over the clusters and increment the hitmap
// Skip out of bound clusters
for (const auto& cluster : cv) {
auto x = cluster.x;
auto y = cluster.y;
if(x<image_size[1] && y<image_size[0])
r(cluster.y, cluster.x) += 1;
}
return hitmap;
});
}

66
python/src/cluster.hpp
View File

@ -64,53 +64,7 @@ void define_cluster(py::module &m, const std::string &typestr) {
*/
}
template <typename Type, uint8_t ClusterSizeX, uint8_t ClusterSizeY,
typename CoordType = uint16_t>
void define_cluster_vector(py::module &m, const std::string &typestr) {
using ClusterType =
Cluster<Type, ClusterSizeX, ClusterSizeY, CoordType, void>;
auto class_name = fmt::format("ClusterVector_{}", typestr);
py::class_<ClusterVector<
Cluster<Type, ClusterSizeX, ClusterSizeY, CoordType, void>, void>>(
m, class_name.c_str(),
py::buffer_protocol())
.def(py::init()) // TODO change!!!
.def("push_back",
[](ClusterVector<ClusterType> &self, const ClusterType &cluster) {
self.push_back(cluster);
})
// implement push_back
.def_property_readonly("size", &ClusterVector<ClusterType>::size)
.def("item_size", &ClusterVector<ClusterType>::item_size)
.def_property_readonly("fmt",
[typestr](ClusterVector<ClusterType> &self) {
return fmt_format<ClusterType>;
})
.def_property_readonly("cluster_size_x",
&ClusterVector<ClusterType>::cluster_size_x)
.def_property_readonly("cluster_size_y",
&ClusterVector<ClusterType>::cluster_size_y)
.def_property_readonly("capacity",
&ClusterVector<ClusterType>::capacity)
.def_property("frame_number", &ClusterVector<ClusterType>::frame_number,
&ClusterVector<ClusterType>::set_frame_number)
.def_buffer(
[typestr](ClusterVector<ClusterType> &self) -> py::buffer_info {
return py::buffer_info(
self.data(), /* Pointer to buffer */
self.item_size(), /* Size of one scalar */
fmt_format<ClusterType>, /* Format descriptor */
1, /* Number of dimensions */
{self.size()}, /* Buffer dimensions */
{self.item_size()} /* Strides (in bytes) for each index */
);
});
}
template <typename T, uint8_t ClusterSizeX, uint8_t ClusterSizeY,
typename CoordType = uint16_t>
@ -252,25 +206,5 @@ void define_cluster_finder_bindings(py::module &m, const std::string &typestr) {
},
py::arg(), py::arg("frame_number") = 0);
m.def("hitmap",
[](std::array<size_t, 2> image_size, ClusterVector<ClusterType> &cv) {
py::array_t<int32_t> hitmap(image_size);
auto r = hitmap.mutable_unchecked<2>();
// Initialize hitmap to 0
for (py::ssize_t i = 0; i < r.shape(0); i++)
for (py::ssize_t j = 0; j < r.shape(1); j++)
r(i, j) = 0;
size_t stride = cv.item_size();
auto ptr = cv.data();
for (size_t i = 0; i < cv.size(); i++) {
auto x = *reinterpret_cast<int16_t *>(ptr);
auto y = *reinterpret_cast<int16_t *>(ptr + sizeof(int16_t));
r(y, x) += 1;
ptr += stride;
}
return hitmap;
});
}
#pragma GCC diagnostic pop

17
python/src/module.cpp
View File

@ -1,4 +1,9 @@
// Files with bindings to the different classes
//New style file naming
#include "bind_ClusterVector.hpp"
//TODO! migrate the other names
#include "cluster.hpp"
#include "cluster_file.hpp"
#include "ctb_raw_file.hpp"
@ -39,12 +44,12 @@ PYBIND11_MODULE(_aare, m) {
define_cluster_file_io_bindings<float, 2, 2, uint16_t>(m, "Cluster2x2f");
define_cluster_file_io_bindings<double, 2, 2, uint16_t>(m, "Cluster2x2d");
define_cluster_vector<int, 3, 3, uint16_t>(m, "Cluster3x3i");
define_cluster_vector<double, 3, 3, uint16_t>(m, "Cluster3x3d");
define_cluster_vector<float, 3, 3, uint16_t>(m, "Cluster3x3f");
define_cluster_vector<int, 2, 2, uint16_t>(m, "Cluster2x2i");
define_cluster_vector<double, 2, 2, uint16_t>(m, "Cluster2x2d");
define_cluster_vector<float, 2, 2, uint16_t>(m, "Cluster2x2f");
define_ClusterVector<int, 3, 3, uint16_t>(m, "Cluster3x3i");
define_ClusterVector<double, 3, 3, uint16_t>(m, "Cluster3x3d");
define_ClusterVector<float, 3, 3, uint16_t>(m, "Cluster3x3f");
define_ClusterVector<int, 2, 2, uint16_t>(m, "Cluster2x2i");
define_ClusterVector<double, 2, 2, uint16_t>(m, "Cluster2x2d");
define_ClusterVector<float, 2, 2, uint16_t>(m, "Cluster2x2f");
define_cluster_finder_bindings<int, 3, 3, uint16_t>(m, "Cluster3x3i");
define_cluster_finder_bindings<double, 3, 3, uint16_t>(m, "Cluster3x3d");

59
python/tests/test_Cluster.py
View File

@ -1,12 +1,12 @@
import pytest
import numpy as np
import aare._aare as aare
from aare import _aare #import the C++ module
from conftest import test_data_path
def test_cluster_vector_can_be_converted_to_numpy():
cv = aare.ClusterVector_Cluster3x3i()
cv = _aare.ClusterVector_Cluster3x3i()
arr = np.array(cv, copy=False)
assert arr.shape == (0,) # 4 for x, y, size, energy and 9 for the cluster data
@ -14,24 +14,23 @@ def test_cluster_vector_can_be_converted_to_numpy():
def test_ClusterVector():
"""Test ClusterVector"""
clustervector = aare.ClusterVector_Cluster3x3i()
clustervector = _aare.ClusterVector_Cluster3x3i()
assert clustervector.cluster_size_x == 3
assert clustervector.cluster_size_y == 3
assert clustervector.item_size() == 4+9*4
assert clustervector.frame_number == 0
assert clustervector.capacity == 1024
assert clustervector.size == 0
cluster = aare.Cluster3x3i(0,0,np.ones(9, dtype=np.int32))
cluster = _aare.Cluster3x3i(0,0,np.ones(9, dtype=np.int32))
clustervector.push_back(cluster)
assert clustervector.size == 1
with pytest.raises(TypeError): # Or use the appropriate exception type
clustervector.push_back(aare.Cluster2x2i(0,0,np.ones(4, dtype=np.int32)))
clustervector.push_back(_aare.Cluster2x2i(0,0,np.ones(4, dtype=np.int32)))
with pytest.raises(TypeError):
clustervector.push_back(aare.Cluster3x3f(0,0,np.ones(9, dtype=np.float32)))
clustervector.push_back(_aare.Cluster3x3f(0,0,np.ones(9, dtype=np.float32)))
def test_Interpolator():
"""Test Interpolator"""
@ -41,13 +40,13 @@ def test_Interpolator():
ybins = np.linspace(0, 5, 30, dtype=np.float64)
etacube = np.zeros(shape=[30, 30, 20], dtype=np.float64)
interpolator = aare.Interpolator(etacube, xbins, ybins, ebins)
interpolator = _aare.Interpolator(etacube, xbins, ybins, ebins)
assert interpolator.get_ietax().shape == (30,30,20)
assert interpolator.get_ietay().shape == (30,30,20)
clustervector = aare.ClusterVector_Cluster3x3i()
clustervector = _aare.ClusterVector_Cluster3x3i()
cluster = aare.Cluster3x3i(0,0, np.ones(9, dtype=np.int32))
cluster = _aare.Cluster3x3i(0,0, np.ones(9, dtype=np.int32))
clustervector.push_back(cluster)
interpolated_photons = interpolator.interpolate(clustervector)
@ -58,9 +57,9 @@ def test_Interpolator():
assert interpolated_photons[0]["y"] == -1
assert interpolated_photons[0]["energy"] == 4 #eta_sum = 4, dx, dy = -1,-1 m_ietax = 0, m_ietay = 0
clustervector = aare.ClusterVector_Cluster2x2i()
clustervector = _aare.ClusterVector_Cluster2x2i()
cluster = aare.Cluster2x2i(0,0, np.ones(4, dtype=np.int32))
cluster = _aare.Cluster2x2i(0,0, np.ones(4, dtype=np.int32))
clustervector.push_back(cluster)
interpolated_photons = interpolator.interpolate(clustervector)
@ -71,28 +70,15 @@ def test_Interpolator():
assert interpolated_photons[0]["y"] == 0
assert interpolated_photons[0]["energy"] == 4
@pytest.mark.files
def test_cluster_file(test_data_path):
"""Test ClusterFile"""
cluster_file = aare.ClusterFile_Cluster3x3i(test_data_path / "clust/single_frame_97_clustrers.clust")
clustervector = cluster_file.read_clusters(10) #conversion does not work
cluster_file.close()
assert clustervector.size == 10
###reading with wrong file
with pytest.raises(TypeError):
cluster_file = aare.ClusterFile_Cluster2x2i(test_data_path / "clust/single_frame_97_clustrers.clust")
cluster_file.close()
def test_calculate_eta():
"""Calculate Eta"""
clusters = aare.ClusterVector_Cluster3x3i()
clusters.push_back(aare.Cluster3x3i(0,0, np.ones(9, dtype=np.int32)))
clusters.push_back(aare.Cluster3x3i(0,0, np.array([1,1,1,2,2,2,3,3,3])))
clusters = _aare.ClusterVector_Cluster3x3i()
clusters.push_back(_aare.Cluster3x3i(0,0, np.ones(9, dtype=np.int32)))
clusters.push_back(_aare.Cluster3x3i(0,0, np.array([1,1,1,2,2,2,3,3,3])))
eta2 = aare.calculate_eta2(clusters)
eta2 = _aare.calculate_eta2(clusters)
assert eta2.shape == (2,2)
assert eta2[0,0] == 0.5
@ -103,7 +89,7 @@ def test_calculate_eta():
def test_cluster_finder():
"""Test ClusterFinder"""
clusterfinder = aare.ClusterFinder_Cluster3x3i([100,100])
clusterfinder = _aare.ClusterFinder_Cluster3x3i([100,100])
#frame = np.random.rand(100,100)
frame = np.zeros(shape=[100,100])
@ -115,18 +101,7 @@ def test_cluster_finder():
assert clusters.size == 0
#TODO dont understand behavior
def test_cluster_collector():
"""Test ClusterCollector"""
clusterfinder = aare.ClusterFinderMT_Cluster3x3i([100,100]) #TODO: no idea what the data is in InputQueue not zero
clustercollector = aare.ClusterCollector_Cluster3x3i(clusterfinder)
cluster_vectors = clustercollector.steal_clusters()
assert len(cluster_vectors) == 1 #single thread execution
assert cluster_vectors[0].size == 0 #

64
python/tests/test_ClusterFile.py Normal file
View File

@ -0,0 +1,64 @@
import pytest
import numpy as np
import boost_histogram as bh
import time
from pathlib import Path
import pickle
from aare import ClusterFile
from conftest import test_data_path
@pytest.mark.files
def test_cluster_file(test_data_path):
"""Test ClusterFile"""
f = ClusterFile(test_data_path / "clust/single_frame_97_clustrers.clust")
cv = f.read_clusters(10) #conversion does not work
assert cv.frame_number == 135
assert cv.size == 10
#Known data
#frame_number, num_clusters [135] 97
#[ 1 200] [0 1 2 3 4 5 6 7 8]
#[ 2 201] [ 9 10 11 12 13 14 15 16 17]
#[ 3 202] [18 19 20 21 22 23 24 25 26]
#[ 4 203] [27 28 29 30 31 32 33 34 35]
#[ 5 204] [36 37 38 39 40 41 42 43 44]
#[ 6 205] [45 46 47 48 49 50 51 52 53]
#[ 7 206] [54 55 56 57 58 59 60 61 62]
#[ 8 207] [63 64 65 66 67 68 69 70 71]
#[ 9 208] [72 73 74 75 76 77 78 79 80]
#[ 10 209] [81 82 83 84 85 86 87 88 89]
#conversion to numpy array
arr = np.array(cv, copy = False)
assert arr.size == 10
for i in range(10):
assert arr[i]['x'] == i+1
@pytest.mark.files
def test_read_clusters_and_fill_histogram(test_data_path):
# Create the histogram
n_bins = 100
xmin = -100
xmax = 1e4
hist_aare = bh.Histogram(bh.axis.Regular(n_bins, xmin, xmax))
fname = test_data_path / "clust/beam_En700eV_-40deg_300V_10us_d0_f0_100.clust"
#Read clusters and fill the histogram with pixel values
with ClusterFile(fname, chunk_size = 10000) as f:
for clusters in f:
arr = np.array(clusters, copy = False)
hist_aare.fill(arr['data'].flat)
#Load the histogram from the pickle file
with open(fname.with_suffix('.pkl'), 'rb') as f:
hist_py = pickle.load(f)
#Compare the two histograms
assert hist_aare == hist_py

54
python/tests/test_ClusterVector.py Normal file
View File

@ -0,0 +1,54 @@
import pytest
import numpy as np
import boost_histogram as bh
import time
from pathlib import Path
import pickle
from aare import ClusterFile
from aare import _aare
from conftest import test_data_path
def test_create_cluster_vector():
cv = _aare.ClusterVector_Cluster3x3i()
assert cv.cluster_size_x == 3
assert cv.cluster_size_y == 3
assert cv.size == 0
def test_push_back_on_cluster_vector():
cv = _aare.ClusterVector_Cluster2x2i()
assert cv.cluster_size_x == 2
assert cv.cluster_size_y == 2
assert cv.size == 0
cluster = _aare.Cluster2x2i(19, 22, np.ones(4, dtype=np.int32))
cv.push_back(cluster)
assert cv.size == 1
arr = np.array(cv, copy=False)
assert arr[0]['x'] == 19
assert arr[0]['y'] == 22
def test_make_a_hitmap_from_cluster_vector():
cv = _aare.ClusterVector_Cluster3x3i()
# Push back 4 clusters with different positions
cv.push_back(_aare.Cluster3x3i(0, 0, np.ones(9, dtype=np.int32)))
cv.push_back(_aare.Cluster3x3i(1, 1, np.ones(9, dtype=np.int32)))
cv.push_back(_aare.Cluster3x3i(1, 1, np.ones(9, dtype=np.int32)))
cv.push_back(_aare.Cluster3x3i(2, 2, np.ones(9, dtype=np.int32)))
ref = np.zeros((5, 5), dtype=np.int32)
ref[0,0] = 1
ref[1,1] = 2
ref[2,2] = 1
img = _aare.hitmap((5,5), cv)
# print(img)
# print(ref)
assert (img == ref).all()

70
src/CalculateEta.test.cpp
View File

@ -37,7 +37,7 @@ auto get_test_parameters() {
Eta2<int>{3. / 5, 4. / 6, 1, 11}));
}
TEST_CASE("compute_largest_2x2_subcluster", "[.eta_calculation]") {
TEST_CASE("compute_largest_2x2_subcluster", "[eta_calculation]") {
auto [cluster, expected_eta] = get_test_parameters();
auto [sum, index] = std::visit(
@ -47,7 +47,7 @@ TEST_CASE("compute_largest_2x2_subcluster", "[.eta_calculation]") {
CHECK(expected_eta.sum == sum);
}
TEST_CASE("calculate_eta2", "[.eta_calculation]") {
TEST_CASE("calculate_eta2", "[eta_calculation]") {
auto [cluster, expected_eta] = get_test_parameters();
@ -60,3 +60,69 @@ TEST_CASE("calculate_eta2", "[.eta_calculation]") {
CHECK(eta.c == expected_eta.c);
CHECK(eta.sum == expected_eta.sum);
}
//3x3 cluster layout (in case of cBottomLeft etc corner):
// 6, 7, 8
// 3, 4, 5
// 0, 1, 2
TEST_CASE("Calculate eta2 for a 3x3 int32 cluster with the largest 2x2 sum in the bottom left",
"[eta_calculation]") {
// Create a 3x3 cluster
Cluster<int32_t, 3, 3> cl;
cl.x = 0;
cl.y = 0;
cl.data[0] = 30;
cl.data[1] = 23;
cl.data[2] = 5;
cl.data[3] = 20;
cl.data[4] = 50;
cl.data[5] = 3;
cl.data[6] = 8;
cl.data[7] = 2;
cl.data[8] = 3;
// 8, 2, 3
// 20, 50, 3
// 30, 23, 5
auto eta = calculate_eta2(cl);
CHECK(eta.c == corner::cBottomLeft);
CHECK(eta.x == 50.0 / (20 + 50)); // 4/(3+4)
CHECK(eta.y == 50.0 / (23 + 50)); // 4/(1+4)
CHECK(eta.sum == 30+23+20+50);
}
TEST_CASE("Calculate eta2 for a 3x3 int32 cluster with the largest 2x2 sum in the top left",
"[eta_calculation]") {
// Create a 3x3 cluster
Cluster<int32_t, 3, 3> cl;
cl.x = 0;
cl.y = 0;
cl.data[0] = 8;
cl.data[1] = 12;
cl.data[2] = 5;
cl.data[3] = 77;
cl.data[4] = 80;
cl.data[5] = 3;
cl.data[6] = 82;
cl.data[7] = 91;
cl.data[8] = 3;
// 82, 91, 3
// 77, 80, 3
// 8, 12, 5
auto eta = calculate_eta2(cl);
CHECK(eta.c == corner::cTopLeft);
CHECK(eta.x == 80. / (77 + 80)); // 4/(3+4)
CHECK(eta.y == 91.0 / (91 + 80)); // 7/(7+4)
CHECK(eta.sum == 77+80+82+91);
}