From ca2223e04c430318b05b4c005495f61240e0a356 Mon Sep 17 00:00:00 2001 From: AndreaRigoni Date: Thu, 19 Mar 2026 09:46:36 +0000 Subject: [PATCH] detector chamber --- src/HEP/Detectors/CMakeLists.txt | 26 +++- src/HEP/Detectors/DetectorChamber.cpp | 51 ++++++++ src/HEP/Detectors/DetectorChamber.h | 16 +++ src/HEP/Detectors/MuonEvent.h | 4 + src/HEP/Detectors/testing/CMakeLists.txt | 3 +- .../Detectors/testing/DetectorChamberTest.cpp | 120 ++++++++++++++++++ src/HEP/Detectors/testing/Makefile.am | 16 --- 7 files changed, 212 insertions(+), 24 deletions(-) create mode 100644 src/HEP/Detectors/DetectorChamber.cpp create mode 100644 src/HEP/Detectors/testing/DetectorChamberTest.cpp delete mode 100644 src/HEP/Detectors/testing/Makefile.am diff --git a/src/HEP/Detectors/CMakeLists.txt b/src/HEP/Detectors/CMakeLists.txt index f793ae5..6536160 100644 --- a/src/HEP/Detectors/CMakeLists.txt +++ b/src/HEP/Detectors/CMakeLists.txt @@ -11,19 +11,31 @@ set(HEADERS MuonScatter.h ) +set(SOURCES + DetectorChamber.cpp +) + +set(LIBRARIES + ${PACKAGE_LIBPREFIX}Core + ${PACKAGE_LIBPREFIX}Math +) + set(libname ${PACKAGE_LIBPREFIX}Detectors) set(ULIB_SHARED_LIBRARIES ${ULIB_SHARED_LIBRARIES} ${libname} PARENT_SCOPE) set(ULIB_SELECTED_MODULES ${ULIB_SELECTED_MODULES} Detectors PARENT_SCOPE) -## Headers-only INTERFACE library -add_library(${libname} INTERFACE) -target_include_directories(${libname} INTERFACE - $ - $ -) +## SHARED library +add_library(${libname} SHARED ${SOURCES}) +set_target_properties(${libname} PROPERTIES + VERSION ${PROJECT_VERSION} + SOVERSION ${PROJECT_SOVERSION} + CXX_STANDARD 17) +target_link_libraries(${libname} PRIVATE ${LIBRARIES}) install(TARGETS ${libname} - EXPORT "uLibTargets") + EXPORT "uLibTargets" + RUNTIME DESTINATION ${INSTALL_BIN_DIR} COMPONENT bin + LIBRARY DESTINATION ${INSTALL_LIB_DIR} COMPONENT lib) install(FILES ${HEADERS} DESTINATION ${INSTALL_INC_DIR}/HEP/Detectors) diff --git a/src/HEP/Detectors/DetectorChamber.cpp b/src/HEP/Detectors/DetectorChamber.cpp new file mode 100644 index 0000000..d093c85 --- /dev/null +++ b/src/HEP/Detectors/DetectorChamber.cpp @@ -0,0 +1,51 @@ +#include "HEP/Detectors/DetectorChamber.h" +#include + +namespace uLib { + +MuonEvent DetectorChamber::ProjectMuonEvent(const MuonEvent &muon) const { + MuonEvent projectedMuon = muon; + + HPoint3f P = m_ProjectionPlane.origin; + HVector3f N = m_ProjectionPlane.direction; + + HPoint3f X_in = muon.LineIn().origin; + HPoint3f X_out = muon.LineOut().origin; + + // Calculate squared distances to the plane normal point for comparison + // Actually, we should probably follow the user's description literally: + // "closest ... with the projection plane ( so the colsest direction point with the point of the normal defining the plane )" + // This could mean point-to-plane or point-to-point. + // Given "closest with the projection plane", point-to-plane is more natural. + // However, "closest direction point with the point of the normal" strongly suggests point-to-point distance. + // Let's use distance to the plane for the first part and keep the logic consistent. + + float dist_in = std::abs((X_in - P).dot(N)); + float dist_out = std::abs((X_out - P).dot(N)); + + const HLine3f &chosenLine = (dist_in <= dist_out) ? muon.LineIn() : muon.LineOut(); + HPoint3f X_chosen = chosenLine.origin; + + // Project X_chosen into the plane defined by P and normal N + // X_proj = X_chosen - ((X_chosen - P) . N / (N . N)) * N + float dot = (X_chosen - P).dot(N); + float n_sq = N.dot(N); + + HPoint3f X_proj = X_chosen; + if (n_sq > 0) { + X_proj = X_chosen - (dot / n_sq) * N; + } + + // Define the projected line with projected origin and original direction + HLine3f projectedLine; + projectedLine.origin = X_proj; + projectedLine.direction = chosenLine.direction; + + // Set both input and output lines of the projected muon to the same projected line + projectedMuon.LineIn() = projectedLine; + projectedMuon.LineOut() = projectedLine; + + return projectedMuon; +} + +} // namespace uLib \ No newline at end of file diff --git a/src/HEP/Detectors/DetectorChamber.h b/src/HEP/Detectors/DetectorChamber.h index a5ed09f..248e059 100644 --- a/src/HEP/Detectors/DetectorChamber.h +++ b/src/HEP/Detectors/DetectorChamber.h @@ -31,14 +31,30 @@ #include "Core/Types.h" #include "Math/ContainerBox.h" +#include "HEP/Detectors/Hit.h" +#include "HEP/Detectors/HitMC.h" +#include "HEP/Detectors/MuonEvent.h" + namespace uLib { class DetectorChamber : public ContainerBox { + typedef ContainerBox BaseClass; + + public: + // set the plane where muons hit is projected + // to be represented with a direction vector + void SetProjectionPlane(const HLine3f &normal) { m_ProjectionPlane = normal; } + + const HLine3f &GetProjectionPlane() const { return m_ProjectionPlane; } + + MuonEvent ProjectMuonEvent(const MuonEvent &muon) const; + private: + HLine3f m_ProjectionPlane; }; diff --git a/src/HEP/Detectors/MuonEvent.h b/src/HEP/Detectors/MuonEvent.h index f779681..e41d290 100644 --- a/src/HEP/Detectors/MuonEvent.h +++ b/src/HEP/Detectors/MuonEvent.h @@ -48,6 +48,10 @@ protected: class MuonEvent : public MuonEventData { public: + using MuonEventData::LineIn; + using MuonEventData::LineOut; + using MuonEventData::GetMomentum; + inline HLine3f & LineIn() { return this->m_LineIn; } inline HLine3f & LineOut() { return this->m_LineOut; } inline Scalarf & Momentum() { return this->m_Momentum; } diff --git a/src/HEP/Detectors/testing/CMakeLists.txt b/src/HEP/Detectors/testing/CMakeLists.txt index 12c7566..a6a9703 100644 --- a/src/HEP/Detectors/testing/CMakeLists.txt +++ b/src/HEP/Detectors/testing/CMakeLists.txt @@ -1,12 +1,13 @@ # TESTS set( TESTS -# GDMLSolidTest HierarchicalEncodingTest + DetectorChamberTest ) set(LIBRARIES ${PACKAGE_LIBPREFIX}Core ${PACKAGE_LIBPREFIX}Math + ${PACKAGE_LIBPREFIX}Detectors Boost::serialization Boost::program_options Eigen3::Eigen diff --git a/src/HEP/Detectors/testing/DetectorChamberTest.cpp b/src/HEP/Detectors/testing/DetectorChamberTest.cpp new file mode 100644 index 0000000..d8e6720 --- /dev/null +++ b/src/HEP/Detectors/testing/DetectorChamberTest.cpp @@ -0,0 +1,120 @@ +/*////////////////////////////////////////////////////////////////////////////// +// CMT Cosmic Muon Tomography project ////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// + + Copyright (c) 2014, Universita' degli Studi di Padova, INFN sez. di Padova + All rights reserved + + Authors: Andrea Rigoni Garola < andrea.rigoni@pd.infn.it > + + ------------------------------------------------------------------ + This library is free software; you can redistribute it and/or + modify it under the terms of the GNU Lesser General Public + License as published by the Free Software Foundation; either + version 3.0 of the License, or (at your option) any later version. + + This library is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + Lesser General Public License for more details. + + You should have received a copy of the GNU Lesser General Public + License along with this library. + +//////////////////////////////////////////////////////////////////////////////*/ + +#include "HEP/Detectors/DetectorChamber.h" +#include "testing-prototype.h" +#include + +using namespace uLib; + +int main() { + BEGIN_TESTING(DetectorChamber Projection); + + DetectorChamber chamber; + + // Define a horizontal plane at z = 100 + HLine3f plane; + plane.origin = HPoint3f(0, 0, 100); + plane.direction = HVector3f(0, 0, 1); // Normal to the plane + chamber.SetProjectionPlane(plane); + + // Create a muon with two segments + MuonEvent muon; + + // Segment 1 (muon.LineIn()): origin is at (10, 20, 50), direction along Z + muon.LineIn().origin = HPoint3f(10, 20, 50); + muon.LineIn().direction = HVector3f(0, 0, 1); + + // Segment 2 (muon.LineOut()): origin is at (10, 20, 200), direction along Z + muon.LineOut().origin = HPoint3f(10, 20, 200); + muon.LineOut().direction = HVector3f(0, 0, 1); + + // distance_in = |50 - 100| = 50 + // distance_out = |200 - 100| = 100 + // LineIn is closer to the plane (z=100) + + MuonEvent projected = chamber.ProjectMuonEvent(muon); + + // Expected: + // chosenLine = LineIn (it is closer) + // X_chosen = (10, 20, 50) + // X_proj = (10, 20, 50) - (( (10, 20, 50) - (0, 0, 100) ) . (0, 0, 1)) * (0, 0, 1) + // X_proj = (10, 20, 50) - (-50) * (0, 0, 1) = (10, 20, 100) + + HPoint3f expectedPos(10, 20, 100); + + std::cout << "Test Case 1: LineIn is closer" << std::endl; + std::cout << "Projected Position: " << projected.LineIn().origin.transpose() << std::endl; + std::cout << "Expected Position: " << expectedPos.transpose() << std::endl; + + // Check if the projected position is correct + // norm() includes the 4th component, which for HVector3f (diff of points) should be 0. + bool posOk1 = (projected.LineIn().origin - expectedPos).norm() < 1e-5; + TEST1(posOk1); + + // Check if LineIn and LineOut are the same + bool linesMatch1 = (projected.LineIn().origin == projected.LineOut().origin) && + (projected.LineIn().direction == projected.LineOut().direction); + TEST1(linesMatch1); + + // Test Case 2: LineOut is closer + muon.LineIn().origin = HPoint3f(30, 40, 0); // dist = 100 + muon.LineOut().origin = HPoint3f(30, 40, 110); // dist = 10 + + projected = chamber.ProjectMuonEvent(muon); + expectedPos = HPoint3f(30, 40, 100); // projection of (30,40,110) onto z=100 + + std::cout << "\nTest Case 2: LineOut is closer" << std::endl; + std::cout << "Projected Position: " << projected.LineIn().origin.transpose() << std::endl; + std::cout << "Expected Position: " << expectedPos.transpose() << std::endl; + + bool posOk2 = (projected.LineIn().origin - expectedPos).norm() < 1e-5; + TEST1(posOk2); + + // Test Case 3: Oblique plane + // Plane through (0,0,0) with normal (1,1,1) (normalized) + plane.origin = HPoint3f(0, 0, 0); + plane.direction = HVector3f(1, 1, 1).normalized(); + chamber.SetProjectionPlane(plane); + + muon.LineIn().origin = HPoint3f(1, 1, 1); // dist = (1,1,1) . (1,1,1).norm() + muon.LineIn().direction = HVector3f(0, 0, 1); + muon.LineOut().origin = HPoint3f(10, 10, 10); // dist = 10 * sqrt(3) + + projected = chamber.ProjectMuonEvent(muon); + // X_chosen = (1,1,1) + // X_proj = (1,1,1) - ((1,1,1) . N) * N where N = (1,1,1)/sqrt(3) + // X_proj = (1,1,1) - (sqrt(3)) * (1,1,1)/sqrt(3) = (1,1,1) - (1,1,1) = (0,0,0) + expectedPos = HPoint3f(0, 0, 0); + + std::cout << "\nTest Case 3: Oblique plane" << std::endl; + std::cout << "Projected Position: " << projected.LineIn().origin.transpose() << std::endl; + std::cout << "Expected Position: " << expectedPos.transpose() << std::endl; + + bool posOk3 = (projected.LineIn().origin - expectedPos).norm() < 1e-5; + TEST1(posOk3); + + END_TESTING; +} diff --git a/src/HEP/Detectors/testing/Makefile.am b/src/HEP/Detectors/testing/Makefile.am deleted file mode 100644 index cd956fc..0000000 --- a/src/HEP/Detectors/testing/Makefile.am +++ /dev/null @@ -1,16 +0,0 @@ - -include $(top_srcdir)/Common.am - -#AM_DEFAULT_SOURCE_EXT = .cpp - -# if HAVE_CHECK -TESTS = GDMLSolidTest - -# else -# TEST = -# endif - -LDADD = $(top_srcdir)/libmutom-${PACKAGE_VERSION}.la - -check_PROGRAMS = $(TESTS) -