From bd2aff778e5b8307401c72aa826dc38b3a8fff82 Mon Sep 17 00:00:00 2001 From: Luis Antonio Obis Aparicio Date: Thu, 1 Jun 2023 23:18:44 +0200 Subject: [PATCH 01/10] experimental new readout alternative --- inc/TRestDetectorExperimentalReadout.h | 28 ++++++++ inc/TRestDetectorExperimentalReadoutModule.h | 42 +++++++++++ inc/TRestDetectorExperimentalReadoutPixel.h | 71 +++++++++++++++++++ src/TRestDetectorExperimentalReadout.cxx | 5 ++ ...TRestDetectorExperimentalReadoutModule.cxx | 5 ++ src/TRestDetectorExperimentalReadoutPixel.cxx | 5 ++ 6 files changed, 156 insertions(+) create mode 100644 inc/TRestDetectorExperimentalReadout.h create mode 100644 inc/TRestDetectorExperimentalReadoutModule.h create mode 100644 inc/TRestDetectorExperimentalReadoutPixel.h create mode 100644 src/TRestDetectorExperimentalReadout.cxx create mode 100644 src/TRestDetectorExperimentalReadoutModule.cxx create mode 100644 src/TRestDetectorExperimentalReadoutPixel.cxx diff --git a/inc/TRestDetectorExperimentalReadout.h b/inc/TRestDetectorExperimentalReadout.h new file mode 100644 index 00000000..751e8e73 --- /dev/null +++ b/inc/TRestDetectorExperimentalReadout.h @@ -0,0 +1,28 @@ +// +// Created by lobis on 01-Jun-23. +// + +#ifndef REST_TRESTDETECTOREXPERIMENTALREADOUT_H +#define REST_TRESTDETECTOREXPERIMENTALREADOUT_H + +#include + +#include + +#include "TRestDetectorExperimentalReadoutModule.h" + +class TRestDetectorExperimentalReadout : public TNamed { + private: + std::map fModules; + // std::vector fModules; + + public: + TRestDetectorExperimentalReadout() = default; + ~TRestDetectorExperimentalReadout() override = default; + + // size_t GetNumberOfModules() { return fModules.size(); } + + ClassDef(TRestDetectorExperimentalReadout, 1); +}; + +#endif // REST_TRESTDETECTOREXPERIMENTALREADOUT_H diff --git a/inc/TRestDetectorExperimentalReadoutModule.h b/inc/TRestDetectorExperimentalReadoutModule.h new file mode 100644 index 00000000..59f9fb6d --- /dev/null +++ b/inc/TRestDetectorExperimentalReadoutModule.h @@ -0,0 +1,42 @@ +// +// Created by lobis on 01-Jun-23. +// + +#ifndef REST_TRESTDETECTOREXPERIMENTALREADOUTMODULE_H +#define REST_TRESTDETECTOREXPERIMENTALREADOUTMODULE_H + +#include + +#include "TRestDetectorExperimentalReadoutPixel.h" + +class TRestDetectorExperimentalReadoutModule { + private: + // std::map fPixels; + + TVector3 fPosition; // position of the module in 3D space + TVector3 fNormal; // normal vector to the module + double fHeight; // height of the module (drift distance) + std::pair fCoordinateAxes; + + public: + + TRestDetectorExperimentalReadoutModule(const TVector3& position, double height, const TVector3& normal, + double rotationInDegrees = 0) { + fPosition = position; + fHeight = height; + fNormal = normal; + + fCoordinateAxes = {{1, 0, 0}, {0, 1, 0}}; + + // rotate the coordinate axes by rotationInDegrees around the normal vector + fCoordinateAxes.first.Rotate(rotationInDegrees * TMath::RadToDeg(), fNormal); + fCoordinateAxes.second.Rotate(rotationInDegrees * TMath::RadToDeg(), fNormal); + } + + ~TRestDetectorExperimentalReadoutModule() = default; + TRestDetectorExperimentalReadoutModule() = default; + + ClassDef(TRestDetectorExperimentalReadoutModule, 1); +}; + +#endif // REST_TRESTDETECTOREXPERIMENTALREADOUTMODULE_H diff --git a/inc/TRestDetectorExperimentalReadoutPixel.h b/inc/TRestDetectorExperimentalReadoutPixel.h new file mode 100644 index 00000000..e8de6baa --- /dev/null +++ b/inc/TRestDetectorExperimentalReadoutPixel.h @@ -0,0 +1,71 @@ +// +// Created by lobis on 01-Jun-23. +// + +#ifndef REST_TRESTDETECTOREXPERIMENTALREADOUTPIXEL_H +#define REST_TRESTDETECTOREXPERIMENTALREADOUTPIXEL_H + +#include +#include + +#include +#include +#include + +typedef unsigned short id; + +class TRestDetectorExperimentalReadoutPixel { + private: + std::vector fVertices; + + void InitializeVertices(const std::vector& vertices) { + // at least three non-aligned vertices + if (vertices.size() < 3) { + throw std::invalid_argument("Pixel must have at least three vertices"); + } + // check the vertices form a convex polygon (TODO) + + fVertices = vertices; + } + + static std::vector GetRectangularVertices(const TVector2& center, + const std::pair& size) { + // rectangular pixel + std::vector vertices = {{center.X() - size.first / 2, center.Y() - size.second / 2}, + {center.X() + size.first / 2, center.Y() - size.second / 2}, + {center.X() + size.first / 2, center.Y() + size.second / 2}, + {center.X() - size.first / 2, center.Y() + size.second / 2}}; + + return vertices; + } + + public: + explicit TRestDetectorExperimentalReadoutPixel(const std::vector& vertices) { InitializeVertices(vertices); } + + TRestDetectorExperimentalReadoutPixel(const TVector2& center, const std::pair& size) { + // rectangular pixel + InitializeVertices(GetRectangularVertices(center, size)); + } + + TRestDetectorExperimentalReadoutPixel(const TVector2& center, double size) { + // square pixel + InitializeVertices(GetRectangularVertices(center, {size, size})); + } + + ~TRestDetectorExperimentalReadoutPixel() = default; + + TVector2 GetCenter() const { + TVector2 center(0, 0); + for (auto& vertex : fVertices) { + center += vertex; + } + center *= 1. / fVertices.size(); + return center; + } + + std::vector GetVertices() const { return fVertices; } + + ClassDef(TRestDetectorExperimentalReadoutPixel, 1); +}; + +#endif // REST_TRESTDETECTOREXPERIMENTALREADOUTPIXEL_H diff --git a/src/TRestDetectorExperimentalReadout.cxx b/src/TRestDetectorExperimentalReadout.cxx new file mode 100644 index 00000000..d5da5c4d --- /dev/null +++ b/src/TRestDetectorExperimentalReadout.cxx @@ -0,0 +1,5 @@ +#include "TRestDetectorExperimentalReadout.h" + +using namespace std; + +ClassImp(TRestDetectorExperimentalReadout); diff --git a/src/TRestDetectorExperimentalReadoutModule.cxx b/src/TRestDetectorExperimentalReadoutModule.cxx new file mode 100644 index 00000000..192be6b2 --- /dev/null +++ b/src/TRestDetectorExperimentalReadoutModule.cxx @@ -0,0 +1,5 @@ +#include "TRestDetectorExperimentalReadoutModule.h" + +using namespace std; + +ClassImp(TRestDetectorExperimentalReadoutModule); diff --git a/src/TRestDetectorExperimentalReadoutPixel.cxx b/src/TRestDetectorExperimentalReadoutPixel.cxx new file mode 100644 index 00000000..91c3ab80 --- /dev/null +++ b/src/TRestDetectorExperimentalReadoutPixel.cxx @@ -0,0 +1,5 @@ +#include "TRestDetectorExperimentalReadoutPixel.h" + +using namespace std; + +ClassImp(TRestDetectorExperimentalReadoutPixel); From 837b5684383bda1b4ee7bceab714e88c4de8e37c Mon Sep 17 00:00:00 2001 From: Luis Antonio Obis Aparicio Date: Fri, 2 Jun 2023 13:21:36 +0200 Subject: [PATCH 02/10] complex hull calculations --- inc/TRestDetectorExperimentalReadout.h | 6 +- inc/TRestDetectorExperimentalReadoutModule.h | 27 ++++- inc/TRestDetectorExperimentalReadoutPixel.h | 70 ++++++++++++- test/src/TRestDetectorReadout.cxx | 102 +++++++++++++++++++ 4 files changed, 195 insertions(+), 10 deletions(-) create mode 100644 test/src/TRestDetectorReadout.cxx diff --git a/inc/TRestDetectorExperimentalReadout.h b/inc/TRestDetectorExperimentalReadout.h index 751e8e73..e5d1e3cd 100644 --- a/inc/TRestDetectorExperimentalReadout.h +++ b/inc/TRestDetectorExperimentalReadout.h @@ -20,7 +20,11 @@ class TRestDetectorExperimentalReadout : public TNamed { TRestDetectorExperimentalReadout() = default; ~TRestDetectorExperimentalReadout() override = default; - // size_t GetNumberOfModules() { return fModules.size(); } + size_t GetNumberOfModules() { return fModules.size(); } + + void AddModule(const TRestDetectorExperimentalReadoutModule& module) { + fModules.insert({fModules.size(), module}); + } ClassDef(TRestDetectorExperimentalReadout, 1); }; diff --git a/inc/TRestDetectorExperimentalReadoutModule.h b/inc/TRestDetectorExperimentalReadoutModule.h index 59f9fb6d..fd6843fb 100644 --- a/inc/TRestDetectorExperimentalReadoutModule.h +++ b/inc/TRestDetectorExperimentalReadoutModule.h @@ -11,7 +11,7 @@ class TRestDetectorExperimentalReadoutModule { private: - // std::map fPixels; + std::vector fPixels; TVector3 fPosition; // position of the module in 3D space TVector3 fNormal; // normal vector to the module @@ -19,7 +19,6 @@ class TRestDetectorExperimentalReadoutModule { std::pair fCoordinateAxes; public: - TRestDetectorExperimentalReadoutModule(const TVector3& position, double height, const TVector3& normal, double rotationInDegrees = 0) { fPosition = position; @@ -33,8 +32,28 @@ class TRestDetectorExperimentalReadoutModule { fCoordinateAxes.second.Rotate(rotationInDegrees * TMath::RadToDeg(), fNormal); } - ~TRestDetectorExperimentalReadoutModule() = default; - TRestDetectorExperimentalReadoutModule() = default; + virtual ~TRestDetectorExperimentalReadoutModule() = default; + TRestDetectorExperimentalReadoutModule() = default; + + size_t GetNumberOfPixels() const { return fPixels.size(); } + + std::vector GetPixels() const { return fPixels; } + + void AddPixel(const TRestDetectorExperimentalReadoutPixel& pixel, id pixelId = 0) { + fPixels.push_back(pixel); + } + + std::vector GetConvexHull() const { + std::vector vertices; + + for (auto& pixel : fPixels) { + for (auto& vertex : pixel.GetVertices()) { + vertices.push_back(vertex); + } + } + + return TRestDetectorExperimentalReadoutPixel::findConvexHull(vertices); + } ClassDef(TRestDetectorExperimentalReadoutModule, 1); }; diff --git a/inc/TRestDetectorExperimentalReadoutPixel.h b/inc/TRestDetectorExperimentalReadoutPixel.h index e8de6baa..90fb261a 100644 --- a/inc/TRestDetectorExperimentalReadoutPixel.h +++ b/inc/TRestDetectorExperimentalReadoutPixel.h @@ -8,6 +8,7 @@ #include #include +#include #include #include #include @@ -23,9 +24,9 @@ class TRestDetectorExperimentalReadoutPixel { if (vertices.size() < 3) { throw std::invalid_argument("Pixel must have at least three vertices"); } - // check the vertices form a convex polygon (TODO) - - fVertices = vertices; + // compute convex hull + std::vector convexHull = findConvexHull(vertices); + fVertices = convexHull; } static std::vector GetRectangularVertices(const TVector2& center, @@ -39,8 +40,32 @@ class TRestDetectorExperimentalReadoutPixel { return vertices; } + // Helper function to calculate the cross product of two vectors + static double crossProduct(const TVector2& A, const TVector2& B, const TVector2& C) { + return (B.X() - A.X()) * (C.Y() - A.Y()) - (B.Y() - A.Y()) * (C.X() - A.X()); + } + + // Helper function to calculate the squared distance between two points + static double squaredDistance(const TVector2& A, const TVector2& B) { + double dx = B.X() - A.X(); + double dy = B.Y() - A.Y(); + return dx * dx + dy * dy; + } + + // Comparator function to compare two points based on their polar angle with respect to the anchor point + static bool comparePoints(const TVector2& A, const TVector2& B, const TVector2& anchor) { + double cross = crossProduct(anchor, A, B); + if (cross == 0) { + // If two points have the same polar angle, choose the one that is closer to the anchor point + return squaredDistance(anchor, A) < squaredDistance(anchor, B); + } + return cross > 0; + } + public: - explicit TRestDetectorExperimentalReadoutPixel(const std::vector& vertices) { InitializeVertices(vertices); } + explicit TRestDetectorExperimentalReadoutPixel(const std::vector& vertices) { + InitializeVertices(vertices); + } TRestDetectorExperimentalReadoutPixel(const TVector2& center, const std::pair& size) { // rectangular pixel @@ -52,7 +77,8 @@ class TRestDetectorExperimentalReadoutPixel { InitializeVertices(GetRectangularVertices(center, {size, size})); } - ~TRestDetectorExperimentalReadoutPixel() = default; + TRestDetectorExperimentalReadoutPixel() = default; + virtual ~TRestDetectorExperimentalReadoutPixel() = default; TVector2 GetCenter() const { TVector2 center(0, 0); @@ -65,6 +91,40 @@ class TRestDetectorExperimentalReadoutPixel { std::vector GetVertices() const { return fVertices; } + // Function to find the convex hull of a set of points using the Graham's scan algorithm + static std::vector findConvexHull(const std::vector& _points) { + std::vector points = _points; + // Find the bottommost point as the anchor point + TVector2 anchor = points[0]; + for (const TVector2& point : points) { + if (point.Y() < anchor.Y() || (point.Y() == anchor.Y() && point.X() < anchor.X())) { + anchor = point; + } + } + + // Sort the points based on their polar angle with respect to the anchor point + std::sort(points.begin(), points.end(), + [&](const TVector2& A, const TVector2& B) { return comparePoints(A, B, anchor); }); + + // Build the convex hull + std::vector convexHull; + convexHull.push_back(points[0]); + convexHull.push_back(points[1]); + int hullSize = 2; + + for (int i = 2; i < points.size(); i++) { + while (hullSize >= 2 && + crossProduct(convexHull[hullSize - 2], convexHull[hullSize - 1], points[i]) <= 0) { + convexHull.pop_back(); + hullSize--; + } + convexHull.push_back(points[i]); + hullSize++; + } + + return convexHull; + } + ClassDef(TRestDetectorExperimentalReadoutPixel, 1); }; diff --git a/test/src/TRestDetectorReadout.cxx b/test/src/TRestDetectorReadout.cxx new file mode 100644 index 00000000..a0bee2e6 --- /dev/null +++ b/test/src/TRestDetectorReadout.cxx @@ -0,0 +1,102 @@ + +#include +#include +#include +#include +#include + +#include + +namespace fs = std::filesystem; + +using namespace std; + +const auto filesPath = fs::path(__FILE__).parent_path().parent_path() / "files"; + +std::tuple GetGraphData(std::vector vertices) { + const auto n = vertices.size() + 1; + + double x[n]; + double y[n]; + + for (int i = 0; i < vertices.size(); i++) { + x[i] = vertices[i].X(); + y[i] = vertices[i].Y(); + } + + x[n - 1] = vertices[0].X(); + y[n - 1] = vertices[0].Y(); + + return {n, x, y}; +} + +TEST(TRestDetectorReadout, Initialize) { + auto readout = TRestDetectorExperimentalReadout(); + + auto module = TRestDetectorExperimentalReadoutModule(); + readout.AddModule(module); + + double pixelSize = 1.0; + int moduleSizeInPixels = 10; + + for (int i = 0; i < moduleSizeInPixels; i++) { + for (int j = 0; j < moduleSizeInPixels; j++) { + TVector2 center = {i * pixelSize, j * pixelSize}; + auto pixel = TRestDetectorExperimentalReadoutPixel(center, pixelSize); + module.AddPixel(pixel); + } + } + + auto canvas = new TCanvas("canvas", "Shape Canvas", 800, 600); + + // first draw the convex hull for the module + const auto hullVertices = module.GetConvexHull(); + + TGraph* graph; + { + const auto& vertices = hullVertices; + const auto n = vertices.size() + 1; + double x[n]; + double y[n]; + + for (int i = 0; i < n - 1; i++) { + x[i] = vertices[i].X(); + y[i] = vertices[i].Y(); + } + x[n - 1] = vertices[0].X(); + y[n - 1] = vertices[0].Y(); + graph = new TGraph(n, x, y); + } + graph->Draw(); + + // print pixel info + for (int i = 0; i < module.GetNumberOfPixels(); i++) { + const auto pixel = module.GetPixels()[i]; + cout << "Pixel X: " << pixel.GetCenter().X() << " Y: " << pixel.GetCenter().Y() << endl; + const auto vertices = pixel.GetVertices(); + + const auto n = vertices.size() + 1; + double x[n]; + double y[n]; + + for (int i = 0; i < n - 1; i++) { + x[i] = vertices[i].X(); + y[i] = vertices[i].Y(); + } + x[n - 1] = vertices[0].X(); + y[n - 1] = vertices[0].Y(); + + for (int i = 0; i < n; i++) { + cout << "i: " << i << " X: " << x[i] << " Y: " << y[i] << endl; + } + + auto graph = new TGraph(vertices.size() + 1, x, y); + + graph->SetLineColor(kRed); + graph->SetLineWidth(2); + + graph->Draw("SAME"); + } + + canvas->Update(); +} From ff49290d3300fc15fdd1c07f0827cfef56086ae7 Mon Sep 17 00:00:00 2001 From: Luis Antonio Obis Aparicio Date: Mon, 5 Jun 2023 10:44:07 +0200 Subject: [PATCH 03/10] cleanup --- inc/TRestDetectorExperimentalReadoutModule.h | 31 +----- inc/TRestDetectorExperimentalReadoutPixel.h | 102 ++++-------------- ...TRestDetectorExperimentalReadoutModule.cxx | 32 ++++++ src/TRestDetectorExperimentalReadoutPixel.cxx | 99 +++++++++++++++++ test/src/TRestDetectorReadout.cxx | 2 +- 5 files changed, 155 insertions(+), 111 deletions(-) diff --git a/inc/TRestDetectorExperimentalReadoutModule.h b/inc/TRestDetectorExperimentalReadoutModule.h index fd6843fb..7738ac38 100644 --- a/inc/TRestDetectorExperimentalReadoutModule.h +++ b/inc/TRestDetectorExperimentalReadoutModule.h @@ -20,40 +20,19 @@ class TRestDetectorExperimentalReadoutModule { public: TRestDetectorExperimentalReadoutModule(const TVector3& position, double height, const TVector3& normal, - double rotationInDegrees = 0) { - fPosition = position; - fHeight = height; - fNormal = normal; + double rotationInDegrees = 0); - fCoordinateAxes = {{1, 0, 0}, {0, 1, 0}}; - - // rotate the coordinate axes by rotationInDegrees around the normal vector - fCoordinateAxes.first.Rotate(rotationInDegrees * TMath::RadToDeg(), fNormal); - fCoordinateAxes.second.Rotate(rotationInDegrees * TMath::RadToDeg(), fNormal); - } - - virtual ~TRestDetectorExperimentalReadoutModule() = default; + // needed for root TRestDetectorExperimentalReadoutModule() = default; + virtual ~TRestDetectorExperimentalReadoutModule() = default; size_t GetNumberOfPixels() const { return fPixels.size(); } std::vector GetPixels() const { return fPixels; } - void AddPixel(const TRestDetectorExperimentalReadoutPixel& pixel, id pixelId = 0) { - fPixels.push_back(pixel); - } - - std::vector GetConvexHull() const { - std::vector vertices; - - for (auto& pixel : fPixels) { - for (auto& vertex : pixel.GetVertices()) { - vertices.push_back(vertex); - } - } + void InsertPixel(const TRestDetectorExperimentalReadoutPixel& pixel); - return TRestDetectorExperimentalReadoutPixel::findConvexHull(vertices); - } + std::vector GetConvexHull() const; ClassDef(TRestDetectorExperimentalReadoutModule, 1); }; diff --git a/inc/TRestDetectorExperimentalReadoutPixel.h b/inc/TRestDetectorExperimentalReadoutPixel.h index 90fb261a..1a820179 100644 --- a/inc/TRestDetectorExperimentalReadoutPixel.h +++ b/inc/TRestDetectorExperimentalReadoutPixel.h @@ -19,64 +19,19 @@ class TRestDetectorExperimentalReadoutPixel { private: std::vector fVertices; - void InitializeVertices(const std::vector& vertices) { - // at least three non-aligned vertices - if (vertices.size() < 3) { - throw std::invalid_argument("Pixel must have at least three vertices"); - } - // compute convex hull - std::vector convexHull = findConvexHull(vertices); - fVertices = convexHull; - } + void InitializeVertices(const std::vector& vertices); static std::vector GetRectangularVertices(const TVector2& center, - const std::pair& size) { - // rectangular pixel - std::vector vertices = {{center.X() - size.first / 2, center.Y() - size.second / 2}, - {center.X() + size.first / 2, center.Y() - size.second / 2}, - {center.X() + size.first / 2, center.Y() + size.second / 2}, - {center.X() - size.first / 2, center.Y() + size.second / 2}}; - - return vertices; - } - - // Helper function to calculate the cross product of two vectors - static double crossProduct(const TVector2& A, const TVector2& B, const TVector2& C) { - return (B.X() - A.X()) * (C.Y() - A.Y()) - (B.Y() - A.Y()) * (C.X() - A.X()); - } - - // Helper function to calculate the squared distance between two points - static double squaredDistance(const TVector2& A, const TVector2& B) { - double dx = B.X() - A.X(); - double dy = B.Y() - A.Y(); - return dx * dx + dy * dy; - } - - // Comparator function to compare two points based on their polar angle with respect to the anchor point - static bool comparePoints(const TVector2& A, const TVector2& B, const TVector2& anchor) { - double cross = crossProduct(anchor, A, B); - if (cross == 0) { - // If two points have the same polar angle, choose the one that is closer to the anchor point - return squaredDistance(anchor, A) < squaredDistance(anchor, B); - } - return cross > 0; - } + const std::pair& size); public: - explicit TRestDetectorExperimentalReadoutPixel(const std::vector& vertices) { - InitializeVertices(vertices); - } + explicit TRestDetectorExperimentalReadoutPixel(const std::vector& vertices); - TRestDetectorExperimentalReadoutPixel(const TVector2& center, const std::pair& size) { - // rectangular pixel - InitializeVertices(GetRectangularVertices(center, size)); - } + TRestDetectorExperimentalReadoutPixel(const TVector2& center, const std::pair& size); - TRestDetectorExperimentalReadoutPixel(const TVector2& center, double size) { - // square pixel - InitializeVertices(GetRectangularVertices(center, {size, size})); - } + TRestDetectorExperimentalReadoutPixel(const TVector2& center, double size); + // needed for root TRestDetectorExperimentalReadoutPixel() = default; virtual ~TRestDetectorExperimentalReadoutPixel() = default; @@ -91,41 +46,20 @@ class TRestDetectorExperimentalReadoutPixel { std::vector GetVertices() const { return fVertices; } - // Function to find the convex hull of a set of points using the Graham's scan algorithm - static std::vector findConvexHull(const std::vector& _points) { - std::vector points = _points; - // Find the bottommost point as the anchor point - TVector2 anchor = points[0]; - for (const TVector2& point : points) { - if (point.Y() < anchor.Y() || (point.Y() == anchor.Y() && point.X() < anchor.X())) { - anchor = point; - } - } + ClassDef(TRestDetectorExperimentalReadoutPixel, 1); +}; - // Sort the points based on their polar angle with respect to the anchor point - std::sort(points.begin(), points.end(), - [&](const TVector2& A, const TVector2& B) { return comparePoints(A, B, anchor); }); - - // Build the convex hull - std::vector convexHull; - convexHull.push_back(points[0]); - convexHull.push_back(points[1]); - int hullSize = 2; - - for (int i = 2; i < points.size(); i++) { - while (hullSize >= 2 && - crossProduct(convexHull[hullSize - 2], convexHull[hullSize - 1], points[i]) <= 0) { - convexHull.pop_back(); - hullSize--; - } - convexHull.push_back(points[i]); - hullSize++; - } +namespace readout { - return convexHull; - } +// Helper function to calculate the cross product of two vectors +double crossProduct(const TVector2& A, const TVector2& B, const TVector2& C); - ClassDef(TRestDetectorExperimentalReadoutPixel, 1); -}; +// Helper function to calculate the squared distance between two points +double squaredDistance(const TVector2& A, const TVector2& B); +// Comparator function to compare two points based on their polar angle with respect to the anchor point +bool comparePoints(const TVector2& A, const TVector2& B, const TVector2& anchor); +// Function to find the convex hull of a set of points using the Graham's scan algorithm +std::vector findConvexHull(const std::vector& _points); +} // namespace readout #endif // REST_TRESTDETECTOREXPERIMENTALREADOUTPIXEL_H diff --git a/src/TRestDetectorExperimentalReadoutModule.cxx b/src/TRestDetectorExperimentalReadoutModule.cxx index 192be6b2..96c33053 100644 --- a/src/TRestDetectorExperimentalReadoutModule.cxx +++ b/src/TRestDetectorExperimentalReadoutModule.cxx @@ -1,5 +1,37 @@ #include "TRestDetectorExperimentalReadoutModule.h" using namespace std; +using namespace readout; ClassImp(TRestDetectorExperimentalReadoutModule); + +TRestDetectorExperimentalReadoutModule::TRestDetectorExperimentalReadoutModule(const TVector3& position, + double height, + const TVector3& normal, + double rotationInDegrees) { + fPosition = position; + fHeight = height; + fNormal = normal; + + fCoordinateAxes = {{1, 0, 0}, {0, 1, 0}}; + + // rotate the coordinate axes by rotationInDegrees around the normal vector + fCoordinateAxes.first.Rotate(rotationInDegrees * TMath::RadToDeg(), fNormal); + fCoordinateAxes.second.Rotate(rotationInDegrees * TMath::RadToDeg(), fNormal); +} + +void TRestDetectorExperimentalReadoutModule::InsertPixel(const TRestDetectorExperimentalReadoutPixel& pixel) { + fPixels.push_back(pixel); +} + +std::vector TRestDetectorExperimentalReadoutModule::GetConvexHull() const { + std::vector vertices; + + for (auto& pixel : fPixels) { + for (auto& vertex : pixel.GetVertices()) { + vertices.push_back(vertex); + } + } + + return findConvexHull(vertices); +} \ No newline at end of file diff --git a/src/TRestDetectorExperimentalReadoutPixel.cxx b/src/TRestDetectorExperimentalReadoutPixel.cxx index 91c3ab80..2f748b9e 100644 --- a/src/TRestDetectorExperimentalReadoutPixel.cxx +++ b/src/TRestDetectorExperimentalReadoutPixel.cxx @@ -1,5 +1,104 @@ #include "TRestDetectorExperimentalReadoutPixel.h" using namespace std; +using namespace readout; ClassImp(TRestDetectorExperimentalReadoutPixel); + +TRestDetectorExperimentalReadoutPixel::TRestDetectorExperimentalReadoutPixel( + const std::vector& vertices) { + InitializeVertices(vertices); +} + +TRestDetectorExperimentalReadoutPixel::TRestDetectorExperimentalReadoutPixel( + const TVector2& center, const std::pair& size) { + // rectangular pixel + InitializeVertices(GetRectangularVertices(center, size)); +} + +TRestDetectorExperimentalReadoutPixel::TRestDetectorExperimentalReadoutPixel(const TVector2& center, + double size) { + // square pixel + InitializeVertices(GetRectangularVertices(center, {size, size})); +} + +void TRestDetectorExperimentalReadoutPixel::InitializeVertices(const std::vector& vertices) { + // at least three non-aligned vertices + if (vertices.size() < 3) { + throw std::invalid_argument("Pixel must have at least three vertices"); + } + // compute convex hull + std::vector convexHull = findConvexHull(vertices); + fVertices = convexHull; +} + +std::vector TRestDetectorExperimentalReadoutPixel::GetRectangularVertices( + const TVector2& center, const std::pair& size) { + // rectangular pixel + std::vector vertices = {{center.X() - size.first / 2, center.Y() - size.second / 2}, + {center.X() + size.first / 2, center.Y() - size.second / 2}, + {center.X() + size.first / 2, center.Y() + size.second / 2}, + {center.X() - size.first / 2, center.Y() + size.second / 2}}; + + return vertices; +} + +// tools + +namespace readout { +// Helper function to calculate the cross product of two vectors +double crossProduct(const TVector2& A, const TVector2& B, const TVector2& C) { + return (B.X() - A.X()) * (C.Y() - A.Y()) - (B.Y() - A.Y()) * (C.X() - A.X()); +} + +// Helper function to calculate the squared distance between two points +double squaredDistance(const TVector2& A, const TVector2& B) { + double dx = B.X() - A.X(); + double dy = B.Y() - A.Y(); + return dx * dx + dy * dy; +} + +// Comparator function to compare two points based on their polar angle with respect to the anchor point +bool comparePoints(const TVector2& A, const TVector2& B, const TVector2& anchor) { + double cross = crossProduct(anchor, A, B); + if (cross == 0) { + // If two points have the same polar angle, choose the one that is closer to the anchor point + return squaredDistance(anchor, A) < squaredDistance(anchor, B); + } + return cross > 0; +} + +// Function to find the convex hull of a set of points using the Graham's scan algorithm +std::vector findConvexHull(const std::vector& _points) { + std::vector points = _points; + // Find the bottommost point as the anchor point + TVector2 anchor = points[0]; + for (const TVector2& point : points) { + if (point.Y() < anchor.Y() || (point.Y() == anchor.Y() && point.X() < anchor.X())) { + anchor = point; + } + } + + // Sort the points based on their polar angle with respect to the anchor point + std::sort(points.begin(), points.end(), + [&](const TVector2& A, const TVector2& B) { return comparePoints(A, B, anchor); }); + + // Build the convex hull + std::vector convexHull; + convexHull.push_back(points[0]); + convexHull.push_back(points[1]); + int hullSize = 2; + + for (int i = 2; i < points.size(); i++) { + while (hullSize >= 2 && + crossProduct(convexHull[hullSize - 2], convexHull[hullSize - 1], points[i]) <= 0) { + convexHull.pop_back(); + hullSize--; + } + convexHull.push_back(points[i]); + hullSize++; + } + + return convexHull; +} +} // namespace readout \ No newline at end of file diff --git a/test/src/TRestDetectorReadout.cxx b/test/src/TRestDetectorReadout.cxx index a0bee2e6..4f2578f7 100644 --- a/test/src/TRestDetectorReadout.cxx +++ b/test/src/TRestDetectorReadout.cxx @@ -43,7 +43,7 @@ TEST(TRestDetectorReadout, Initialize) { for (int j = 0; j < moduleSizeInPixels; j++) { TVector2 center = {i * pixelSize, j * pixelSize}; auto pixel = TRestDetectorExperimentalReadoutPixel(center, pixelSize); - module.AddPixel(pixel); + module.InsertPixel(pixel); } } From af4873c456bda3879955cf94d827420ee4ac65e3 Mon Sep 17 00:00:00 2001 From: Luis Antonio Obis Aparicio Date: Mon, 5 Jun 2023 13:41:19 +0200 Subject: [PATCH 04/10] added auxiliary methods --- inc/TRestDetectorExperimentalReadoutModule.h | 30 +++- inc/TRestDetectorExperimentalReadoutPixel.h | 17 +- ...TRestDetectorExperimentalReadoutModule.cxx | 160 +++++++++++++++++- src/TRestDetectorExperimentalReadoutPixel.cxx | 41 +++++ test/src/TRestDetectorReadout.cxx | 22 ++- 5 files changed, 247 insertions(+), 23 deletions(-) diff --git a/inc/TRestDetectorExperimentalReadoutModule.h b/inc/TRestDetectorExperimentalReadoutModule.h index 7738ac38..bc399178 100644 --- a/inc/TRestDetectorExperimentalReadoutModule.h +++ b/inc/TRestDetectorExperimentalReadoutModule.h @@ -6,9 +6,12 @@ #define REST_TRESTDETECTOREXPERIMENTALREADOUTMODULE_H #include +#include #include "TRestDetectorExperimentalReadoutPixel.h" +class KDTree; + class TRestDetectorExperimentalReadoutModule { private: std::vector fPixels; @@ -18,21 +21,42 @@ class TRestDetectorExperimentalReadoutModule { double fHeight; // height of the module (drift distance) std::pair fCoordinateAxes; + // auxiliary data structures + // std::unique_ptr fKDTree; + KDTree* fKDTree = nullptr; + double fSearchRadius = 0; + std::vector fConvexHull; + + void BuildKDTree(); + std::vector ComputeConvexHull(); + public: TRestDetectorExperimentalReadoutModule(const TVector3& position, double height, const TVector3& normal, double rotationInDegrees = 0); // needed for root TRestDetectorExperimentalReadoutModule() = default; - virtual ~TRestDetectorExperimentalReadoutModule() = default; + virtual ~TRestDetectorExperimentalReadoutModule() { + // delete fKDTree; + } size_t GetNumberOfPixels() const { return fPixels.size(); } std::vector GetPixels() const { return fPixels; } - void InsertPixel(const TRestDetectorExperimentalReadoutPixel& pixel); + void SetPixels(const std::vector&); + + std::vector GetConvexHull() const { return fConvexHull; } + + std::vector GetPixelsInPoint(const TVector2& point) const; + + double GetZ(const TVector3& point) const; + + bool IsInside(const TVector2& point) const; + bool IsInside(const TVector3& point) const; - std::vector GetConvexHull() const; + TVector2 TransformToModuleCoordinates(const TVector3& point) const; + TVector3 TransformToAbsoluteCoordinates(const TVector2& point) const; ClassDef(TRestDetectorExperimentalReadoutModule, 1); }; diff --git a/inc/TRestDetectorExperimentalReadoutPixel.h b/inc/TRestDetectorExperimentalReadoutPixel.h index 1a820179..0a9b645d 100644 --- a/inc/TRestDetectorExperimentalReadoutPixel.h +++ b/inc/TRestDetectorExperimentalReadoutPixel.h @@ -18,6 +18,8 @@ typedef unsigned short id; class TRestDetectorExperimentalReadoutPixel { private: std::vector fVertices; + TVector2 fCenter; + double fRadius; void InitializeVertices(const std::vector& vertices); @@ -35,17 +37,13 @@ class TRestDetectorExperimentalReadoutPixel { TRestDetectorExperimentalReadoutPixel() = default; virtual ~TRestDetectorExperimentalReadoutPixel() = default; - TVector2 GetCenter() const { - TVector2 center(0, 0); - for (auto& vertex : fVertices) { - center += vertex; - } - center *= 1. / fVertices.size(); - return center; - } + TVector2 GetCenter() const { return fCenter; } + double GetRadius() const { return fRadius; } std::vector GetVertices() const { return fVertices; } + bool IsInside(const TVector2& point) const; + ClassDef(TRestDetectorExperimentalReadoutPixel, 1); }; @@ -61,5 +59,8 @@ double squaredDistance(const TVector2& A, const TVector2& B); bool comparePoints(const TVector2& A, const TVector2& B, const TVector2& anchor); // Function to find the convex hull of a set of points using the Graham's scan algorithm std::vector findConvexHull(const std::vector& _points); + +bool IsPointInsideConvexHull(const TVector2& point, const std::vector& convexHull); + } // namespace readout #endif // REST_TRESTDETECTOREXPERIMENTALREADOUTPIXEL_H diff --git a/src/TRestDetectorExperimentalReadoutModule.cxx b/src/TRestDetectorExperimentalReadoutModule.cxx index 96c33053..493afd01 100644 --- a/src/TRestDetectorExperimentalReadoutModule.cxx +++ b/src/TRestDetectorExperimentalReadoutModule.cxx @@ -20,11 +20,24 @@ TRestDetectorExperimentalReadoutModule::TRestDetectorExperimentalReadoutModule(c fCoordinateAxes.second.Rotate(rotationInDegrees * TMath::RadToDeg(), fNormal); } -void TRestDetectorExperimentalReadoutModule::InsertPixel(const TRestDetectorExperimentalReadoutPixel& pixel) { - fPixels.push_back(pixel); +void TRestDetectorExperimentalReadoutModule::SetPixels( + const std::vector& pixels) { + fPixels = pixels; + fConvexHull = ComputeConvexHull(); + + // search radius is the maximum radius of all pixels + double searchRadius = 0; + for (auto& pixel : fPixels) { + if (pixel.GetRadius() > searchRadius) { + searchRadius = pixel.GetRadius(); + } + } + fSearchRadius = searchRadius; + + BuildKDTree(); } -std::vector TRestDetectorExperimentalReadoutModule::GetConvexHull() const { +vector TRestDetectorExperimentalReadoutModule::ComputeConvexHull() { std::vector vertices; for (auto& pixel : fPixels) { @@ -34,4 +47,143 @@ std::vector TRestDetectorExperimentalReadoutModule::GetConvexHull() co } return findConvexHull(vertices); -} \ No newline at end of file +} + +// KD Tree + +struct KDNode { + int index; // Index of the point in the original collection + KDNode* left; + KDNode* right; + + KDNode(int _index) : index(_index), left(nullptr), right(nullptr) {} +}; + +class KDTree { + private: + KDNode* root; + std::vector& points; + + KDNode* buildTree(std::vector& indices, int start, int end, int depth) { + if (start > end) { + return nullptr; + } + + int axis = depth % 2; // Split alternately on x and y coordinates + int mid = (start + end) / 2; + + // Sort the indices based on the current splitting axis + if (axis == 0) { + std::sort(indices.begin() + start, indices.begin() + end + 1, + [this](int a, int b) { return points[a].X() < points[b].X(); }); + } else { + std::sort(indices.begin() + start, indices.begin() + end + 1, + [this](int a, int b) { return points[a].Y() < points[b].Y(); }); + } + + KDNode* node = new KDNode(indices[mid]); + + node->left = buildTree(indices, start, mid - 1, depth + 1); + node->right = buildTree(indices, mid + 1, end, depth + 1); + + return node; + } + + void query(KDNode* node, const TVector2& target, double distance, std::vector& result) { + if (node == nullptr) { + return; + } + + double dist = (points[node->index] - target).Mod(); + if (dist <= distance) { + result.push_back(node->index); + } + + int axis = points[node->index].X() < target.X() + ? 0 + : 1; // Determine the axis to search based on target position + + if (axis == 0) { + query(node->left, target, distance, result); + + if (target.X() + distance > points[node->index].X()) { + query(node->right, target, distance, result); + } + } else { + query(node->right, target, distance, result); + + if (target.X() - distance < points[node->index].X()) { + query(node->left, target, distance, result); + } + } + } + + public: + KDTree(std::vector& _points) : points(_points) { + std::vector indices(points.size()); + for (int i = 0; i < points.size(); ++i) { + indices[i] = i; + } + root = buildTree(indices, 0, points.size() - 1, 0); + } + + std::vector queryIndices(const TVector2& target, double distance) { + std::vector result; + query(root, target, distance, result); + return result; + } +}; + +void TRestDetectorExperimentalReadoutModule::BuildKDTree() { + std::vector centers; + centers.reserve(fPixels.size()); + for (const auto& pixel : fPixels) { + centers.push_back(pixel.GetCenter()); + } + + // auto tree = KDTree(centers); + + // delete fKDTree; + fKDTree = new KDTree(centers); +} + +std::vector TRestDetectorExperimentalReadoutModule::GetPixelsInPoint( + const TVector2& point) const { + auto indices = fKDTree->queryIndices(point, fSearchRadius); + cout << "Found " << indices.size() << " pixels in point" << endl; + + std::vector pixels; + pixels.reserve(indices.size()); + for (auto& index : indices) { + pixels.push_back(fPixels[index]); + } + + return pixels; +} + +bool TRestDetectorExperimentalReadoutModule::IsInside(const TVector2& point) const { + return IsPointInsideConvexHull(point, fConvexHull); +} + +bool TRestDetectorExperimentalReadoutModule::IsInside(const TVector3& point) const { + const double z = GetZ(point); + if (z < 0 || z > fHeight) { + return false; + } + return IsInside(TransformToModuleCoordinates(point)); +} + +TVector2 TRestDetectorExperimentalReadoutModule::TransformToModuleCoordinates(const TVector3& point) const { + return { + (point - fPosition).Dot(fCoordinateAxes.first), + (point - fPosition).Dot(fCoordinateAxes.second), + }; +} + +TVector3 TRestDetectorExperimentalReadoutModule::TransformToAbsoluteCoordinates(const TVector2& point) const { + return fPosition + point.X() * fCoordinateAxes.first + point.Y() * fCoordinateAxes.second; +} + +double TRestDetectorExperimentalReadoutModule::GetZ(const TVector3& point) const { + (point - fPosition).Dot(fNormal); +} diff --git a/src/TRestDetectorExperimentalReadoutPixel.cxx b/src/TRestDetectorExperimentalReadoutPixel.cxx index 2f748b9e..2c4212e5 100644 --- a/src/TRestDetectorExperimentalReadoutPixel.cxx +++ b/src/TRestDetectorExperimentalReadoutPixel.cxx @@ -30,6 +30,24 @@ void TRestDetectorExperimentalReadoutPixel::InitializeVertices(const std::vector // compute convex hull std::vector convexHull = findConvexHull(vertices); fVertices = convexHull; + + // compute center and radius + TVector2 center(0, 0); + for (auto& vertex : fVertices) { + center += vertex; + } + center *= 1. / double(fVertices.size()); + fCenter = center; + + // the radius is the maximum distance from the center to any vertex + double radius = 0; + for (auto& vertex : fVertices) { + double distance = (vertex - center).Mod(); + if (distance > radius) { + radius = distance; + } + } + fRadius = radius; } std::vector TRestDetectorExperimentalReadoutPixel::GetRectangularVertices( @@ -43,6 +61,10 @@ std::vector TRestDetectorExperimentalReadoutPixel::GetRectangularVerti return vertices; } +bool TRestDetectorExperimentalReadoutPixel::IsInside(const TVector2& point) const { + return readout::IsPointInsideConvexHull(point, fVertices); +} + // tools namespace readout { @@ -101,4 +123,23 @@ std::vector findConvexHull(const std::vector& _points) { return convexHull; } + +bool IsPointInsideConvexHull(const TVector2& point, const std::vector& convexHull) { + // Iterate through each edge of the convex hull + for (size_t i = 0; i < convexHull.size(); ++i) { + const TVector2& vertex1 = convexHull[i]; + const TVector2& vertex2 = convexHull[(i + 1) % convexHull.size()]; + + // Calculate the cross product of the point with the current edge + double crossProductValue = readout::crossProduct(vertex1, vertex2, point); + + // If the cross product is negative for any edge, the point is outside the convex hull + if (crossProductValue < 0.0) { + return false; + } + } + + // If the point is not outside any edge, it is inside the convex hull + return true; +} } // namespace readout \ No newline at end of file diff --git a/test/src/TRestDetectorReadout.cxx b/test/src/TRestDetectorReadout.cxx index 4f2578f7..b29af73e 100644 --- a/test/src/TRestDetectorReadout.cxx +++ b/test/src/TRestDetectorReadout.cxx @@ -37,16 +37,21 @@ TEST(TRestDetectorReadout, Initialize) { readout.AddModule(module); double pixelSize = 1.0; - int moduleSizeInPixels = 10; + int pixelsPerSide = 25; - for (int i = 0; i < moduleSizeInPixels; i++) { - for (int j = 0; j < moduleSizeInPixels; j++) { + std::vector pixels; + for (int i = 0; i < pixelsPerSide; i++) { + for (int j = 0; j < pixelsPerSide; j++) { TVector2 center = {i * pixelSize, j * pixelSize}; - auto pixel = TRestDetectorExperimentalReadoutPixel(center, pixelSize); - module.InsertPixel(pixel); + pixels.emplace_back(center, pixelSize); } } + module.SetPixels(pixels); + + module.GetPixelsInPoint({0, 0}); + return; + auto canvas = new TCanvas("canvas", "Shape Canvas", 800, 600); // first draw the convex hull for the module @@ -72,7 +77,7 @@ TEST(TRestDetectorReadout, Initialize) { // print pixel info for (int i = 0; i < module.GetNumberOfPixels(); i++) { const auto pixel = module.GetPixels()[i]; - cout << "Pixel X: " << pixel.GetCenter().X() << " Y: " << pixel.GetCenter().Y() << endl; + // cout << "Pixel X: " << pixel.GetCenter().X() << " Y: " << pixel.GetCenter().Y() << endl; const auto vertices = pixel.GetVertices(); const auto n = vertices.size() + 1; @@ -87,7 +92,7 @@ TEST(TRestDetectorReadout, Initialize) { y[n - 1] = vertices[0].Y(); for (int i = 0; i < n; i++) { - cout << "i: " << i << " X: " << x[i] << " Y: " << y[i] << endl; + // cout << "i: " << i << " X: " << x[i] << " Y: " << y[i] << endl; } auto graph = new TGraph(vertices.size() + 1, x, y); @@ -96,7 +101,8 @@ TEST(TRestDetectorReadout, Initialize) { graph->SetLineWidth(2); graph->Draw("SAME"); + canvas->Update(); } - + gPad->SetFixedAspectRatio(); canvas->Update(); } From 6e429c6a715c766538b18edc75fe9125365e83b6 Mon Sep 17 00:00:00 2001 From: Luis Antonio Obis Aparicio Date: Tue, 6 Jun 2023 10:24:01 +0200 Subject: [PATCH 05/10] add some documentation --- inc/TRestDetectorExperimentalReadout.h | 60 +++++++--- inc/TRestDetectorExperimentalReadoutModule.h | 94 +++++++++++++-- inc/TRestDetectorExperimentalReadoutPixel.h | 111 +++++++++++++++--- src/TRestDetectorExperimentalReadout.cxx | 21 ++++ ...TRestDetectorExperimentalReadoutModule.cxx | 4 +- 5 files changed, 245 insertions(+), 45 deletions(-) diff --git a/inc/TRestDetectorExperimentalReadout.h b/inc/TRestDetectorExperimentalReadout.h index e5d1e3cd..039f556a 100644 --- a/inc/TRestDetectorExperimentalReadout.h +++ b/inc/TRestDetectorExperimentalReadout.h @@ -1,6 +1,7 @@ -// -// Created by lobis on 01-Jun-23. -// +/** +* @file TRestDetectorExperimentalReadout.h +* @brief Defines the TRestDetectorExperimentalReadout class. +*/ #ifndef REST_TRESTDETECTOREXPERIMENTALREADOUT_H #define REST_TRESTDETECTOREXPERIMENTALREADOUT_H @@ -11,22 +12,53 @@ #include "TRestDetectorExperimentalReadoutModule.h" +/** +* @class TRestDetectorExperimentalReadout +* @brief Represents an experimental readout detector. +* @inherits TNamed +* @author Luis Obis +*/ class TRestDetectorExperimentalReadout : public TNamed { - private: - std::map fModules; - // std::vector fModules; + private: + std::vector fModules; /**< The modules of the detector. */ - public: - TRestDetectorExperimentalReadout() = default; - ~TRestDetectorExperimentalReadout() override = default; + public: + /** + * @brief Constructs a TRestDetectorExperimentalReadout object. + */ + TRestDetectorExperimentalReadout() = default; - size_t GetNumberOfModules() { return fModules.size(); } + /** + * @brief Destructor for TRestDetectorExperimentalReadout. + */ + ~TRestDetectorExperimentalReadout() override = default; - void AddModule(const TRestDetectorExperimentalReadoutModule& module) { - fModules.insert({fModules.size(), module}); - } + /** + * @brief Returns the number of modules in the detector. + * @return The number of modules. + */ + size_t GetNumberOfModules() { return fModules.size(); } - ClassDef(TRestDetectorExperimentalReadout, 1); + /** + * @brief Adds a module to the detector. + * @param module The module to add. + */ + void AddModule(const TRestDetectorExperimentalReadoutModule& module) { fModules.push_back(module); } + + /** + * @brief Returns the modules of the detector. + * @return The modules. + */ + std::vector GetModules() const; + + /** + * @brief Returns the modules that contain a given point. + * @param point The point to check. + * @return The modules containing the point. + */ + std::vector GetModulesInPoint(const TVector3& point) const; + + ClassDef(TRestDetectorExperimentalReadout, 1); }; #endif // REST_TRESTDETECTOREXPERIMENTALREADOUT_H diff --git a/inc/TRestDetectorExperimentalReadoutModule.h b/inc/TRestDetectorExperimentalReadoutModule.h index bc399178..ce5c48ab 100644 --- a/inc/TRestDetectorExperimentalReadoutModule.h +++ b/inc/TRestDetectorExperimentalReadoutModule.h @@ -1,36 +1,56 @@ -// -// Created by lobis on 01-Jun-23. -// +/** + * @file TRestDetectorExperimentalReadoutModule.h + * @brief Defines the TRestDetectorExperimentalReadoutModule class. + */ #ifndef REST_TRESTDETECTOREXPERIMENTALREADOUTMODULE_H #define REST_TRESTDETECTOREXPERIMENTALREADOUTMODULE_H #include -#include #include "TRestDetectorExperimentalReadoutPixel.h" class KDTree; +/** + * @class TRestDetectorExperimentalReadoutModule + * @brief Represents a module in an experimental readout detector. + * @author Luis Obis + */ class TRestDetectorExperimentalReadoutModule { private: - std::vector fPixels; + std::vector fPixels; /**< The pixels of the module. */ - TVector3 fPosition; // position of the module in 3D space - TVector3 fNormal; // normal vector to the module - double fHeight; // height of the module (drift distance) - std::pair fCoordinateAxes; + TVector3 fPosition; /**< The position of the module in 3D space. */ + TVector3 fNormal; /**< The normal vector to the module. */ + double fHeight; /**< The height of the module (drift distance). */ + std::pair fCoordinateAxes; /**< The coordinate axes of the module. */ // auxiliary data structures // std::unique_ptr fKDTree; - KDTree* fKDTree = nullptr; - double fSearchRadius = 0; - std::vector fConvexHull; + KDTree* fKDTree = nullptr; /**< The KDTree used for spatial queries. */ + double fSearchRadius = 0; /**< The search radius for spatial queries. */ + std::vector fConvexHull; /**< The convex hull of the module. */ + /** + * @brief Builds the KDTree for the module. + */ void BuildKDTree(); + + /** + * @brief Computes the convex hull of the module. + * @return The convex hull of the module. + */ std::vector ComputeConvexHull(); public: + /** + * @brief Constructs a TRestDetectorExperimentalReadoutModule object with the specified parameters. + * @param position The position of the module in 3D space. + * @param height The height of the module (drift distance). + * @param normal The normal vector to the module. + * @param rotationInDegrees The rotation angle of the module in degrees (default: 0). + */ TRestDetectorExperimentalReadoutModule(const TVector3& position, double height, const TVector3& normal, double rotationInDegrees = 0); @@ -40,22 +60,70 @@ class TRestDetectorExperimentalReadoutModule { // delete fKDTree; } + /** + * @brief Returns the number of pixels in the module. + * @return The number of pixels. + */ size_t GetNumberOfPixels() const { return fPixels.size(); } + /** + * @brief Returns the pixels of the module. + * @return The pixels. + */ std::vector GetPixels() const { return fPixels; } - void SetPixels(const std::vector&); + /** + * @brief Sets the pixels of the module. + * @param pixels The pixels to set. + */ + void SetPixels(const std::vector& pixels); + /** + * @brief Returns the convex hull of the module. + * @return The convex hull. + */ std::vector GetConvexHull() const { return fConvexHull; } + /** + * @brief Returns the pixels that contain a given point. + * @param point The point to check. + * @return The pixels containing the point. + */ std::vector GetPixelsInPoint(const TVector2& point) const; + /** + * @brief Returns the Z coordinate of a given point inside the module. + * @param point The point to check. + * @return The Z coordinate of the point. + */ double GetZ(const TVector3& point) const; + /** + * @brief Checks if a given point is inside the module in the module's local XY plane. + * @param point The point to check. + * @return True if the point is inside the module, false otherwise. + */ bool IsInside(const TVector2& point) const; + + /** + * @brief Checks if a given point is inside the module in 3D space. + * @param point The point to check. + * @return True if the point is inside the module, false otherwise. + */ bool IsInside(const TVector3& point) const; + /** + * @brief Transforms a point from 3D space to the module's local XY plane coordinates. + * @param point The point to transform. + * @return The transformed point in the local XY plane. + */ TVector2 TransformToModuleCoordinates(const TVector3& point) const; + + /** + * @brief Transforms a point from the module's local XY plane coordinates to 3D space. + * @param point The point to transform. + * @return The transformed point in 3D space. + */ TVector3 TransformToAbsoluteCoordinates(const TVector2& point) const; ClassDef(TRestDetectorExperimentalReadoutModule, 1); diff --git a/inc/TRestDetectorExperimentalReadoutPixel.h b/inc/TRestDetectorExperimentalReadoutPixel.h index 0a9b645d..7e83f1ec 100644 --- a/inc/TRestDetectorExperimentalReadoutPixel.h +++ b/inc/TRestDetectorExperimentalReadoutPixel.h @@ -1,6 +1,7 @@ -// -// Created by lobis on 01-Jun-23. -// +/** + * @file TRestDetectorExperimentalReadoutPixel.h + * @brief Defines the TRestDetectorExperimentalReadoutPixel class. + */ #ifndef REST_TRESTDETECTOREXPERIMENTALREADOUTPIXEL_H #define REST_TRESTDETECTOREXPERIMENTALREADOUTPIXEL_H @@ -13,54 +14,132 @@ #include #include -typedef unsigned short id; +/** + * @class TRestDetectorExperimentalReadoutPixel + * @brief Represents a pixel in an experimental readout detector. + * @author Luis Obis + +*/ class TRestDetectorExperimentalReadoutPixel { private: - std::vector fVertices; - TVector2 fCenter; - double fRadius; + std::vector fVertices; /**< The vertices of the pixel. */ + TVector2 fCenter; /**< The center of the pixel in mm. A pixel is guaranteed to be inside a circle with + center fCenter and radius fRadius. */ + double fRadius; /**< The radius of the pixel in mm. A pixel is guaranteed to be inside a circle with + center fCenter and radius fRadius. */ + + /** + @brief Initializes the vertices of the pixel. + @param vertices The vertices of the pixel. + */ void InitializeVertices(const std::vector& vertices); + /** + @brief Calculates the vertices of a rectangular pixel. + @param center The center of the pixel. + @param size The size of the pixel (width, height). + @return The vertices of the rectangular pixel. + */ static std::vector GetRectangularVertices(const TVector2& center, const std::pair& size); public: + /** + + @brief Constructs a TRestDetectorExperimentalReadoutPixel object with given vertices. + @param vertices The vertices of the pixel. + */ explicit TRestDetectorExperimentalReadoutPixel(const std::vector& vertices); + /** + @brief Constructs a rectangular TRestDetectorExperimentalReadoutPixel object. + @param center The center of the pixel. + @param size The size of the pixel (width, height). + */ TRestDetectorExperimentalReadoutPixel(const TVector2& center, const std::pair& size); + /** + @brief Constructs a square TRestDetectorExperimentalReadoutPixel object. + @param center The center of the pixel. + @param size The size of the pixel (width and height). + */ TRestDetectorExperimentalReadoutPixel(const TVector2& center, double size); - // needed for root TRestDetectorExperimentalReadoutPixel() = default; virtual ~TRestDetectorExperimentalReadoutPixel() = default; + /** + + @brief Returns the center of the pixel. + @return The center of the pixel. + */ TVector2 GetCenter() const { return fCenter; } + /** + + @brief Returns the radius of the pixel. + @return The radius of the pixel. + */ double GetRadius() const { return fRadius; } + /** + @brief Returns the vertices of the pixel. + @return The vertices of the pixel. + */ std::vector GetVertices() const { return fVertices; } + /** + @brief Checks if a given point is inside the pixel. + @param point The point to check. + @return True if the point is inside the pixel, false otherwise. + */ bool IsInside(const TVector2& point) const; - ClassDef(TRestDetectorExperimentalReadoutPixel, 1); }; namespace readout { -// Helper function to calculate the cross product of two vectors +/** + +@brief Calculates the triple product of three vectors. +@param A The first vector. +@param B The second vector. +@param C The third vector. +@return The cross product of the vectors. +*/ double crossProduct(const TVector2& A, const TVector2& B, const TVector2& C); +/** -// Helper function to calculate the squared distance between two points +@brief Calculates the squared distance between two points. +@param A The first point. +@param B The second point. +@return The squared distance between the points. +*/ double squaredDistance(const TVector2& A, const TVector2& B); - -// Comparator function to compare two points based on their polar angle with respect to the anchor point +/** + +@brief Comparator function to compare two points based on their polar angle with respect to the anchor point. +@param A The first point. +@param B The second point. +@param anchor The anchor point. +@return True if point A is smaller than point B, false otherwise. +*/ bool comparePoints(const TVector2& A, const TVector2& B, const TVector2& anchor); -// Function to find the convex hull of a set of points using the Graham's scan algorithm +/** + +@brief Finds the convex hull of a set of points using the Graham's scan algorithm. +@param _points The input points. +@return The convex hull of the points. +*/ std::vector findConvexHull(const std::vector& _points); +/** +@brief Checks if a given point is inside a convex hull. +@param point The point to check. +@param convexHull The convex hull. +@return True if the point is inside the convex hull, false otherwise. +*/ bool IsPointInsideConvexHull(const TVector2& point, const std::vector& convexHull); - } // namespace readout -#endif // REST_TRESTDETECTOREXPERIMENTALREADOUTPIXEL_H +#endif // REST_TRESTDETECTOREXPERIMENTALREADOUTPIXEL_H \ No newline at end of file diff --git a/src/TRestDetectorExperimentalReadout.cxx b/src/TRestDetectorExperimentalReadout.cxx index d5da5c4d..3d146931 100644 --- a/src/TRestDetectorExperimentalReadout.cxx +++ b/src/TRestDetectorExperimentalReadout.cxx @@ -3,3 +3,24 @@ using namespace std; ClassImp(TRestDetectorExperimentalReadout); + +std::vector +TRestDetectorExperimentalReadout::GetModulesInPoint(const TVector3& point) const { + std::vector modules; + for (const auto& module : fModules) { + if (module.IsInside(point)) { + modules.push_back(&module); + } + } + return modules; +} + +std::vector TRestDetectorExperimentalReadout::GetModules() + const { + std::vector modules; + modules.reserve(fModules.size()); + for (const auto& module : fModules) { + modules.push_back(&module); + } + return modules; +} diff --git a/src/TRestDetectorExperimentalReadoutModule.cxx b/src/TRestDetectorExperimentalReadoutModule.cxx index 493afd01..b69a4c22 100644 --- a/src/TRestDetectorExperimentalReadoutModule.cxx +++ b/src/TRestDetectorExperimentalReadoutModule.cxx @@ -121,10 +121,10 @@ class KDTree { public: KDTree(std::vector& _points) : points(_points) { std::vector indices(points.size()); - for (int i = 0; i < points.size(); ++i) { + for (size_t i = 0; i < points.size(); ++i) { indices[i] = i; } - root = buildTree(indices, 0, points.size() - 1, 0); + root = buildTree(indices, 0, int(points.size()) - 1, 0); } std::vector queryIndices(const TVector2& target, double distance) { From 9470f2134a967a323da481b07b1f881439ba7113 Mon Sep 17 00:00:00 2001 From: Luis Antonio Obis Aparicio Date: Tue, 6 Jun 2023 11:37:06 +0200 Subject: [PATCH 06/10] update geometrical transformations --- inc/TRestDetectorExperimentalReadoutModule.h | 35 ++++++----- ...TRestDetectorExperimentalReadoutModule.cxx | 58 ++++++++++++++----- 2 files changed, 62 insertions(+), 31 deletions(-) diff --git a/inc/TRestDetectorExperimentalReadoutModule.h b/inc/TRestDetectorExperimentalReadoutModule.h index ce5c48ab..c86375fe 100644 --- a/inc/TRestDetectorExperimentalReadoutModule.h +++ b/inc/TRestDetectorExperimentalReadoutModule.h @@ -21,13 +21,14 @@ class TRestDetectorExperimentalReadoutModule { private: std::vector fPixels; /**< The pixels of the module. */ - TVector3 fPosition; /**< The position of the module in 3D space. */ - TVector3 fNormal; /**< The normal vector to the module. */ - double fHeight; /**< The height of the module (drift distance). */ - std::pair fCoordinateAxes; /**< The coordinate axes of the module. */ - + TVector3 fPosition = {0, 0, 0}; /**< The position of the module in 3D space. */ + TVector3 fNormal = {0, 0, 1}; /**< The normal vector to the module. */ + double fHeight = 1.0; /**< The height of the module (drift distance). */ + std::pair fCoordinateAxes = {{1, 0, 0}, + {0, 1, 0}}; /**< The coordinate axes of the module. */ + double fRotation = + 0; /**< The rotation angle of the module with respect to the normal vector in degrees. */ // auxiliary data structures - // std::unique_ptr fKDTree; KDTree* fKDTree = nullptr; /**< The KDTree used for spatial queries. */ double fSearchRadius = 0; /**< The search radius for spatial queries. */ std::vector fConvexHull; /**< The convex hull of the module. */ @@ -43,23 +44,25 @@ class TRestDetectorExperimentalReadoutModule { */ std::vector ComputeConvexHull(); - public: - /** - * @brief Constructs a TRestDetectorExperimentalReadoutModule object with the specified parameters. - * @param position The position of the module in 3D space. - * @param height The height of the module (drift distance). - * @param normal The normal vector to the module. - * @param rotationInDegrees The rotation angle of the module in degrees (default: 0). - */ - TRestDetectorExperimentalReadoutModule(const TVector3& position, double height, const TVector3& normal, - double rotationInDegrees = 0); + void UpdateAxes(); + public: // needed for root TRestDetectorExperimentalReadoutModule() = default; virtual ~TRestDetectorExperimentalReadoutModule() { // delete fKDTree; } + void SetPosition(const TVector3& position) { fPosition = position; } + void SetNormal(const TVector3& normal); + void SetHeight(double height); + void SetRotation(double rotationInDegrees); + + TVector3 GetPosition() const { return fPosition; } + TVector3 GetNormal() const { return fNormal; } + double GetHeight() const { return fHeight; } + std::pair GetCoordinateAxes() const { return fCoordinateAxes; } + double GetRotation() const { return fRotation; } /** * @brief Returns the number of pixels in the module. * @return The number of pixels. diff --git a/src/TRestDetectorExperimentalReadoutModule.cxx b/src/TRestDetectorExperimentalReadoutModule.cxx index b69a4c22..09aad684 100644 --- a/src/TRestDetectorExperimentalReadoutModule.cxx +++ b/src/TRestDetectorExperimentalReadoutModule.cxx @@ -5,21 +5,6 @@ using namespace readout; ClassImp(TRestDetectorExperimentalReadoutModule); -TRestDetectorExperimentalReadoutModule::TRestDetectorExperimentalReadoutModule(const TVector3& position, - double height, - const TVector3& normal, - double rotationInDegrees) { - fPosition = position; - fHeight = height; - fNormal = normal; - - fCoordinateAxes = {{1, 0, 0}, {0, 1, 0}}; - - // rotate the coordinate axes by rotationInDegrees around the normal vector - fCoordinateAxes.first.Rotate(rotationInDegrees * TMath::RadToDeg(), fNormal); - fCoordinateAxes.second.Rotate(rotationInDegrees * TMath::RadToDeg(), fNormal); -} - void TRestDetectorExperimentalReadoutModule::SetPixels( const std::vector& pixels) { fPixels = pixels; @@ -187,3 +172,46 @@ TVector3 TRestDetectorExperimentalReadoutModule::TransformToAbsoluteCoordinates( double TRestDetectorExperimentalReadoutModule::GetZ(const TVector3& point) const { (point - fPosition).Dot(fNormal); } + +void TRestDetectorExperimentalReadoutModule::SetNormal(const TVector3& normal) { + fNormal = normal.Unit(); + UpdateAxes(); +} + +void TRestDetectorExperimentalReadoutModule::SetRotation(double rotationInDegrees) { + fRotation = rotationInDegrees; + UpdateAxes(); +} + +void TRestDetectorExperimentalReadoutModule::SetHeight(double height) { + if (height <= 0) { + cerr << "TRestDetectorExperimentalReadoutModule::SetHeight: Height must be positive" << endl; + return; + } + fHeight = height; +} + +void TRestDetectorExperimentalReadoutModule::UpdateAxes() { + // idempotent + + fCoordinateAxes.first = {1, 0, 0}; + fCoordinateAxes.second = {0, 1, 0}; + + // Check if fNormal is different from the original normal + const TVector3 originalNormal = {0, 0, 1}; + if (fNormal != originalNormal) { + // Calculate the rotation axis by taking the cross product between the original normal and fNormal + TVector3 rotationAxis = originalNormal.Cross(fNormal); + + // Calculate the rotation angle using the dot product between the original normal and fNormal + double rotationAngle = acos(originalNormal.Dot(fNormal) / (originalNormal.Mag() * fNormal.Mag())); + + // Rotate the axes around the rotation axis by the rotation angle + fCoordinateAxes.first.Rotate(rotationAngle, rotationAxis); + fCoordinateAxes.second.Rotate(rotationAngle, rotationAxis); + } + + // rotate around normal by rotation angle + fCoordinateAxes.first.Rotate(fRotation * TMath::DegToRad(), fNormal); + fCoordinateAxes.second.Rotate(fRotation * TMath::DegToRad(), fNormal); +} From d013c75ae0e269d5a3190bb0b1690e6b69e6bfca Mon Sep 17 00:00:00 2001 From: Luis Antonio Obis Aparicio Date: Tue, 6 Jun 2023 11:54:33 +0200 Subject: [PATCH 07/10] update setters and getters --- inc/TRestDetectorExperimentalReadout.h | 87 ++++++++++--------- inc/TRestDetectorExperimentalReadoutModule.h | 3 +- src/TRestDetectorExperimentalReadout.cxx | 25 ++++-- ...TRestDetectorExperimentalReadoutModule.cxx | 14 ++- test/src/TRestDetectorReadout.cxx | 2 +- 5 files changed, 78 insertions(+), 53 deletions(-) diff --git a/inc/TRestDetectorExperimentalReadout.h b/inc/TRestDetectorExperimentalReadout.h index 039f556a..37d81d44 100644 --- a/inc/TRestDetectorExperimentalReadout.h +++ b/inc/TRestDetectorExperimentalReadout.h @@ -1,7 +1,7 @@ /** -* @file TRestDetectorExperimentalReadout.h -* @brief Defines the TRestDetectorExperimentalReadout class. -*/ + * @file TRestDetectorExperimentalReadout.h + * @brief Defines the TRestDetectorExperimentalReadout class. + */ #ifndef REST_TRESTDETECTOREXPERIMENTALREADOUT_H #define REST_TRESTDETECTOREXPERIMENTALREADOUT_H @@ -13,52 +13,57 @@ #include "TRestDetectorExperimentalReadoutModule.h" /** -* @class TRestDetectorExperimentalReadout -* @brief Represents an experimental readout detector. -* @inherits TNamed -* @author Luis Obis -*/ + * @class TRestDetectorExperimentalReadout + * @brief Represents an experimental readout detector. + * @inherits TNamed + * @author Luis Obis + */ class TRestDetectorExperimentalReadout : public TNamed { - private: - std::vector fModules; /**< The modules of the detector. */ + private: + std::vector fModules; /**< The modules of the detector. */ - public: - /** - * @brief Constructs a TRestDetectorExperimentalReadout object. - */ - TRestDetectorExperimentalReadout() = default; + public: + /** + * @brief Constructs a TRestDetectorExperimentalReadout object. + */ + TRestDetectorExperimentalReadout() = default; - /** - * @brief Destructor for TRestDetectorExperimentalReadout. - */ - ~TRestDetectorExperimentalReadout() override = default; + /** + * @brief Destructor for TRestDetectorExperimentalReadout. + */ + ~TRestDetectorExperimentalReadout() override = default; - /** - * @brief Returns the number of modules in the detector. - * @return The number of modules. - */ - size_t GetNumberOfModules() { return fModules.size(); } + /** + * @brief Returns the number of modules in the detector. + * @return The number of modules. + */ + size_t GetNumberOfModules() { return fModules.size(); } - /** - * @brief Adds a module to the detector. - * @param module The module to add. - */ - void AddModule(const TRestDetectorExperimentalReadoutModule& module) { fModules.push_back(module); } + /** + * @brief Adds a module to the detector. + * @param module The module to add. + */ + void AddModule(const TRestDetectorExperimentalReadoutModule& module) { fModules.push_back(module); } - /** - * @brief Returns the modules of the detector. - * @return The modules. - */ - std::vector GetModules() const; + /** + * @brief Returns the modules of the detector. + * @return The modules. + */ + std::vector GetModules() const; - /** - * @brief Returns the modules that contain a given point. - * @param point The point to check. - * @return The modules containing the point. - */ - std::vector GetModulesInPoint(const TVector3& point) const; + /** + * @brief Returns the modules that contain a given point. + * @param point The point to check. + * @return The modules containing the point. + */ + std::vector GetModulesForPoint( + const TVector3& point) const; - ClassDef(TRestDetectorExperimentalReadout, 1); + std::vector GetPixelsForPoint(const TVector3& point) const; + + std::vector GetChannelsForPoint(const TVector3& point) const; + + ClassDef(TRestDetectorExperimentalReadout, 1); }; #endif // REST_TRESTDETECTOREXPERIMENTALREADOUT_H diff --git a/inc/TRestDetectorExperimentalReadoutModule.h b/inc/TRestDetectorExperimentalReadoutModule.h index c86375fe..66b143ef 100644 --- a/inc/TRestDetectorExperimentalReadoutModule.h +++ b/inc/TRestDetectorExperimentalReadoutModule.h @@ -92,7 +92,8 @@ class TRestDetectorExperimentalReadoutModule { * @param point The point to check. * @return The pixels containing the point. */ - std::vector GetPixelsInPoint(const TVector2& point) const; + std::vector GetPixelsForPoint(const TVector2& point) const; + std::vector GetPixelsForPoint(const TVector3& point) const; /** * @brief Returns the Z coordinate of a given point inside the module. diff --git a/src/TRestDetectorExperimentalReadout.cxx b/src/TRestDetectorExperimentalReadout.cxx index 3d146931..052371eb 100644 --- a/src/TRestDetectorExperimentalReadout.cxx +++ b/src/TRestDetectorExperimentalReadout.cxx @@ -4,9 +4,9 @@ using namespace std; ClassImp(TRestDetectorExperimentalReadout); -std::vector -TRestDetectorExperimentalReadout::GetModulesInPoint(const TVector3& point) const { - std::vector modules; +vector TRestDetectorExperimentalReadout::GetModulesForPoint( + const TVector3& point) const { + vector modules; for (const auto& module : fModules) { if (module.IsInside(point)) { modules.push_back(&module); @@ -15,12 +15,25 @@ TRestDetectorExperimentalReadout::GetModulesInPoint(const TVector3& point) const return modules; } -std::vector TRestDetectorExperimentalReadout::GetModules() - const { - std::vector modules; +vector TRestDetectorExperimentalReadout::GetModules() const { + vector modules; modules.reserve(fModules.size()); for (const auto& module : fModules) { modules.push_back(&module); } return modules; } + +vector TRestDetectorExperimentalReadout::GetPixelsForPoint( + const TVector3& point) const { + vector pixels; + for (const auto& module : GetModulesForPoint(point)) { + module->GetPixelsForPoint(point); + } +} + +vector TRestDetectorExperimentalReadout::GetChannelsForPoint(const TVector3& point) const { + vector channels; + + return channels; +} diff --git a/src/TRestDetectorExperimentalReadoutModule.cxx b/src/TRestDetectorExperimentalReadoutModule.cxx index 09aad684..0be99293 100644 --- a/src/TRestDetectorExperimentalReadoutModule.cxx +++ b/src/TRestDetectorExperimentalReadoutModule.cxx @@ -132,15 +132,15 @@ void TRestDetectorExperimentalReadoutModule::BuildKDTree() { fKDTree = new KDTree(centers); } -std::vector TRestDetectorExperimentalReadoutModule::GetPixelsInPoint( - const TVector2& point) const { +std::vector +TRestDetectorExperimentalReadoutModule::GetPixelsForPoint(const TVector2& point) const { auto indices = fKDTree->queryIndices(point, fSearchRadius); cout << "Found " << indices.size() << " pixels in point" << endl; - std::vector pixels; + std::vector pixels; pixels.reserve(indices.size()); for (auto& index : indices) { - pixels.push_back(fPixels[index]); + pixels.push_back(&fPixels[index]); } return pixels; @@ -215,3 +215,9 @@ void TRestDetectorExperimentalReadoutModule::UpdateAxes() { fCoordinateAxes.first.Rotate(fRotation * TMath::DegToRad(), fNormal); fCoordinateAxes.second.Rotate(fRotation * TMath::DegToRad(), fNormal); } + +std::vector +TRestDetectorExperimentalReadoutModule::GetPixelsForPoint(const TVector3& point) const { + const TVector2 pointInModuleCoordinates = TransformToModuleCoordinates(point); + return GetPixelsForPoint(pointInModuleCoordinates); +} diff --git a/test/src/TRestDetectorReadout.cxx b/test/src/TRestDetectorReadout.cxx index b29af73e..e521b3a5 100644 --- a/test/src/TRestDetectorReadout.cxx +++ b/test/src/TRestDetectorReadout.cxx @@ -49,7 +49,7 @@ TEST(TRestDetectorReadout, Initialize) { module.SetPixels(pixels); - module.GetPixelsInPoint({0, 0}); + module.GetPixelsForPoint({0, 0}); return; auto canvas = new TCanvas("canvas", "Shape Canvas", 800, 600); From 479056e2eb76686da2b6c449fa90238b22c03811 Mon Sep 17 00:00:00 2001 From: Luis Antonio Obis Aparicio Date: Tue, 6 Jun 2023 12:23:36 +0200 Subject: [PATCH 08/10] added channels --- inc/TRestDetectorExperimentalReadoutModule.h | 5 ++++- inc/TRestDetectorExperimentalReadoutPixel.h | 12 ++++++++--- src/TRestDetectorExperimentalReadout.cxx | 4 ++++ ...TRestDetectorExperimentalReadoutModule.cxx | 21 +++++++++++++++++++ src/TRestDetectorExperimentalReadoutPixel.cxx | 12 ++++++----- 5 files changed, 45 insertions(+), 9 deletions(-) diff --git a/inc/TRestDetectorExperimentalReadoutModule.h b/inc/TRestDetectorExperimentalReadoutModule.h index 66b143ef..0b25393e 100644 --- a/inc/TRestDetectorExperimentalReadoutModule.h +++ b/inc/TRestDetectorExperimentalReadoutModule.h @@ -32,7 +32,8 @@ class TRestDetectorExperimentalReadoutModule { KDTree* fKDTree = nullptr; /**< The KDTree used for spatial queries. */ double fSearchRadius = 0; /**< The search radius for spatial queries. */ std::vector fConvexHull; /**< The convex hull of the module. */ - + std::map> + fChannelToPixels; /**< A map from channel to pixels. */ /** * @brief Builds the KDTree for the module. */ @@ -94,6 +95,8 @@ class TRestDetectorExperimentalReadoutModule { */ std::vector GetPixelsForPoint(const TVector2& point) const; std::vector GetPixelsForPoint(const TVector3& point) const; + std::vector GetPixelsForChannel( + unsigned short channel) const; /** * @brief Returns the Z coordinate of a given point inside the module. diff --git a/inc/TRestDetectorExperimentalReadoutPixel.h b/inc/TRestDetectorExperimentalReadoutPixel.h index 7e83f1ec..c158b952 100644 --- a/inc/TRestDetectorExperimentalReadoutPixel.h +++ b/inc/TRestDetectorExperimentalReadoutPixel.h @@ -29,6 +29,7 @@ class TRestDetectorExperimentalReadoutPixel { double fRadius; /**< The radius of the pixel in mm. A pixel is guaranteed to be inside a circle with center fCenter and radius fRadius. */ + unsigned short fChannel; /**< The channel of the pixel. */ /** @brief Initializes the vertices of the pixel. @@ -51,21 +52,23 @@ class TRestDetectorExperimentalReadoutPixel { @brief Constructs a TRestDetectorExperimentalReadoutPixel object with given vertices. @param vertices The vertices of the pixel. */ - explicit TRestDetectorExperimentalReadoutPixel(const std::vector& vertices); + explicit TRestDetectorExperimentalReadoutPixel(const std::vector& vertices, + unsigned short channel = 0); /** @brief Constructs a rectangular TRestDetectorExperimentalReadoutPixel object. @param center The center of the pixel. @param size The size of the pixel (width, height). */ - TRestDetectorExperimentalReadoutPixel(const TVector2& center, const std::pair& size); + TRestDetectorExperimentalReadoutPixel(const TVector2& center, const std::pair& size, + unsigned short channel = 0); /** @brief Constructs a square TRestDetectorExperimentalReadoutPixel object. @param center The center of the pixel. @param size The size of the pixel (width and height). */ - TRestDetectorExperimentalReadoutPixel(const TVector2& center, double size); + TRestDetectorExperimentalReadoutPixel(const TVector2& center, double size, unsigned short channel = 0); // needed for root TRestDetectorExperimentalReadoutPixel() = default; virtual ~TRestDetectorExperimentalReadoutPixel() = default; @@ -95,6 +98,9 @@ class TRestDetectorExperimentalReadoutPixel { @return True if the point is inside the pixel, false otherwise. */ bool IsInside(const TVector2& point) const; + + unsigned short GetChannel() const { return fChannel; } + ClassDef(TRestDetectorExperimentalReadoutPixel, 1); }; diff --git a/src/TRestDetectorExperimentalReadout.cxx b/src/TRestDetectorExperimentalReadout.cxx index 052371eb..02d817bc 100644 --- a/src/TRestDetectorExperimentalReadout.cxx +++ b/src/TRestDetectorExperimentalReadout.cxx @@ -35,5 +35,9 @@ vector TRestDetectorExperimentalRe vector TRestDetectorExperimentalReadout::GetChannelsForPoint(const TVector3& point) const { vector channels; + for (const auto& pixel : GetPixelsForPoint(point)) { + channels.push_back(pixel->GetChannel()); + } + return channels; } diff --git a/src/TRestDetectorExperimentalReadoutModule.cxx b/src/TRestDetectorExperimentalReadoutModule.cxx index 0be99293..a3356e7c 100644 --- a/src/TRestDetectorExperimentalReadoutModule.cxx +++ b/src/TRestDetectorExperimentalReadoutModule.cxx @@ -10,6 +10,12 @@ void TRestDetectorExperimentalReadoutModule::SetPixels( fPixels = pixels; fConvexHull = ComputeConvexHull(); + // build map of pixels by channel (a pixel channel cannot be updated so this will remain current) + fChannelToPixels.clear(); + for (auto& pixel : fPixels) { + fChannelToPixels[pixel.GetChannel()].push_back(&pixel); + } + // search radius is the maximum radius of all pixels double searchRadius = 0; for (auto& pixel : fPixels) { @@ -221,3 +227,18 @@ TRestDetectorExperimentalReadoutModule::GetPixelsForPoint(const TVector3& point) const TVector2 pointInModuleCoordinates = TransformToModuleCoordinates(point); return GetPixelsForPoint(pointInModuleCoordinates); } + +std::vector +TRestDetectorExperimentalReadoutModule::GetPixelsForChannel(unsigned short channel) const { + std::vector pixels; + // if channel is not in map, return + if (fChannelToPixels.find(channel) == fChannelToPixels.end()) { + return pixels; + } + + pixels.reserve(fChannelToPixels.at(channel).size()); + for (const auto& pixel : fChannelToPixels.at(channel)) { + pixels.push_back(pixel); + } + return pixels; +} diff --git a/src/TRestDetectorExperimentalReadoutPixel.cxx b/src/TRestDetectorExperimentalReadoutPixel.cxx index 2c4212e5..239e4a91 100644 --- a/src/TRestDetectorExperimentalReadoutPixel.cxx +++ b/src/TRestDetectorExperimentalReadoutPixel.cxx @@ -6,20 +6,22 @@ using namespace readout; ClassImp(TRestDetectorExperimentalReadoutPixel); TRestDetectorExperimentalReadoutPixel::TRestDetectorExperimentalReadoutPixel( - const std::vector& vertices) { + const std::vector& vertices, unsigned short channel) + : fChannel(channel) { InitializeVertices(vertices); } TRestDetectorExperimentalReadoutPixel::TRestDetectorExperimentalReadoutPixel( - const TVector2& center, const std::pair& size) { + const TVector2& center, const std::pair& size, unsigned short channel) + : TRestDetectorExperimentalReadoutPixel(GetRectangularVertices(center, size), channel) { // rectangular pixel - InitializeVertices(GetRectangularVertices(center, size)); } TRestDetectorExperimentalReadoutPixel::TRestDetectorExperimentalReadoutPixel(const TVector2& center, - double size) { + double size, + unsigned short channel) + : TRestDetectorExperimentalReadoutPixel(GetRectangularVertices(center, {size, size}), channel) { // square pixel - InitializeVertices(GetRectangularVertices(center, {size, size})); } void TRestDetectorExperimentalReadoutPixel::InitializeVertices(const std::vector& vertices) { From 5842ad09799ee8547a6ad8f22fc018f91baf4e32 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Tue, 6 Jun 2023 10:29:45 +0000 Subject: [PATCH 09/10] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- inc/TRestDetectorExperimentalReadoutPixel.h | 2 +- src/TRestDetectorExperimentalReadoutPixel.cxx | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/inc/TRestDetectorExperimentalReadoutPixel.h b/inc/TRestDetectorExperimentalReadoutPixel.h index c158b952..97dea0c5 100644 --- a/inc/TRestDetectorExperimentalReadoutPixel.h +++ b/inc/TRestDetectorExperimentalReadoutPixel.h @@ -148,4 +148,4 @@ std::vector findConvexHull(const std::vector& _points); */ bool IsPointInsideConvexHull(const TVector2& point, const std::vector& convexHull); } // namespace readout -#endif // REST_TRESTDETECTOREXPERIMENTALREADOUTPIXEL_H \ No newline at end of file +#endif // REST_TRESTDETECTOREXPERIMENTALREADOUTPIXEL_H diff --git a/src/TRestDetectorExperimentalReadoutPixel.cxx b/src/TRestDetectorExperimentalReadoutPixel.cxx index 239e4a91..c27df273 100644 --- a/src/TRestDetectorExperimentalReadoutPixel.cxx +++ b/src/TRestDetectorExperimentalReadoutPixel.cxx @@ -144,4 +144,4 @@ bool IsPointInsideConvexHull(const TVector2& point, const std::vector& // If the point is not outside any edge, it is inside the convex hull return true; } -} // namespace readout \ No newline at end of file +} // namespace readout From 2ad0ca80718e88b1e7fc17dcfd9883dda8e5d9d7 Mon Sep 17 00:00:00 2001 From: Luis Antonio Obis Aparicio Date: Tue, 6 Jun 2023 14:20:59 +0200 Subject: [PATCH 10/10] added channel class --- inc/TRestDetectorExperimentalReadoutChannel.h | 30 +++++++++++++++++++ ...RestDetectorExperimentalReadoutChannel.cxx | 6 ++++ 2 files changed, 36 insertions(+) create mode 100644 inc/TRestDetectorExperimentalReadoutChannel.h create mode 100644 src/TRestDetectorExperimentalReadoutChannel.cxx diff --git a/inc/TRestDetectorExperimentalReadoutChannel.h b/inc/TRestDetectorExperimentalReadoutChannel.h new file mode 100644 index 00000000..c0388764 --- /dev/null +++ b/inc/TRestDetectorExperimentalReadoutChannel.h @@ -0,0 +1,30 @@ +/** + * @file TRestDetectorExperimentalReadoutModule.h + * @brief Defines the TRestDetectorExperimentalReadoutModule class. + */ + +#ifndef REST_TRESTDETECTOREXPERIMENTALREADOUTCHANNEL_H +#define REST_TRESTDETECTOREXPERIMENTALREADOUTCHANNEL_H + +#include + +#include "TRestDetectorExperimentalReadoutPixel.h" + +/** + * @class TRestDetectorExperimentalReadoutChannel + * @brief Represents a channel in an experimental readout detector. + * @author Luis Obis + */ +class TRestDetectorExperimentalReadoutChannel { + private: + std::vector fPixels; /**< The pixels of the channel. */ + unsigned short fChannel; /**< The channel of the pixel. */ + + public: + TRestDetectorExperimentalReadoutChannel() = default; + ~TRestDetectorExperimentalReadoutChannel() = default; + + ClassDef(TRestDetectorExperimentalReadoutChannel, 1); +}; + +#endif // REST_TRESTDETECTOREXPERIMENTALREADOUTCHANNEL_H diff --git a/src/TRestDetectorExperimentalReadoutChannel.cxx b/src/TRestDetectorExperimentalReadoutChannel.cxx new file mode 100644 index 00000000..f113732d --- /dev/null +++ b/src/TRestDetectorExperimentalReadoutChannel.cxx @@ -0,0 +1,6 @@ +#include "TRestDetectorExperimentalReadoutChannel.h" + +using namespace std; +using namespace readout; + +ClassImp(TRestDetectorExperimentalReadoutChannel);