Fatras/Digitization/ChannelizerTests.cpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2023 CERN for the benefit of the Acts project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
0008
0009 #include <boost/test/unit_test.hpp>
0010
0011 #include "Acts/Definitions/Units.hpp"
0012 #include "Acts/Surfaces/PlaneSurface.hpp"
0013 #include "Acts/Utilities/BinUtility.hpp"
0014 #include "ActsFatras/Digitization/Channelizer.hpp"
0015 #include "ActsFatras/Digitization/PlanarSurfaceDrift.hpp"
0016 #include "ActsFatras/Digitization/PlanarSurfaceMask.hpp"
0017
0018 #include <numeric>
0019
0020 using namespace Acts::UnitLiterals;
0021
0022 struct Helper {
0023   std::shared_ptr<Acts::Surface> surface;
0024   Acts::BinUtility segmentation;
0025
0026   Acts::GeometryContext gctx{};
0027   double thickness = 125_um;
0028   Acts::Vector3 driftDir = Acts::Vector3::Zero();
0029
0030   ActsFatras::Channelizer channelizer;
0031
0032   Helper() {
0033     surface = Acts::Surface::makeShared<Acts::PlaneSurface>(
0034         Acts::Vector3::Zero(), Acts::Vector3{0.0, 0.0, 1.0});
0035
0036     float pitchSize = 50_um;
0037     float min = -200_um;
0038     float max = 200_um;
0039     int bins = static_cast<int>((max - min) / pitchSize);
0040     segmentation = Acts::BinUtility(bins, min, max, Acts::BinningOption::open,
0041                                     Acts::BinningValue::binX);
0042     segmentation += Acts::BinUtility(bins, min, max, Acts::BinningOption::open,
0043                                      Acts::BinningValue::binY);
0044   }
0045
0046   auto channelize(const Acts::Vector3 &pos3, const Acts::Vector3 &dir3) const {
0047     Acts::Vector4 pos4 = Acts::Vector4::Zero();
0048     pos4.segment<3>(Acts::ePos0) = pos3;
0049     Acts::Vector4 mom4 = Acts::Vector4::Zero();
0050     mom4.segment<3>(Acts::eMom0) = dir3;
0051     ActsFatras::Hit hit({}, {}, pos4, mom4, mom4);
0052     auto res = channelizer.channelize(hit, *surface, gctx, driftDir,
0053                                       segmentation, thickness);
0054     BOOST_REQUIRE(res.ok());
0055     return *res;
0056   }
0057 };
0058
0059 BOOST_AUTO_TEST_CASE(test_upright_particle) {
0060   Helper helper;
0061
0062   Acts::Vector3 pos3 = Acts::Vector3{10_um, 10_um, 0.0};
0063   Acts::Vector3 dir3 = Acts::Vector3{0.0, 0.0, helper.thickness}.normalized();
0064
0065   auto segments = helper.channelize(pos3, dir3);
0066
0067   BOOST_CHECK_EQUAL(segments.size(), 1);
0068   BOOST_CHECK_CLOSE(segments[0].activation, helper.thickness, 1.e-8);
0069 }
0070
0071 BOOST_AUTO_TEST_CASE(test_tilted_particle) {
0072   Helper helper;
0073
0074   const double disp = 10_um;
0075
0076   Acts::Vector3 hitPosition = Acts::Vector3{10_um, 10_um, 0.0};
0077   Acts::Vector3 hitDirection =
0078       Acts::Vector3({disp, 0.0, helper.thickness}).normalized();
0079
0080   auto segments = helper.channelize(hitPosition, hitDirection);
0081
0082   std::cout << "Segments:\n";
0083   for (const auto &seg : segments) {
0084     std::cout << " - (" << seg.bin[0] << ", " << seg.bin[1]
0085               << "), activation: " << seg.activation << "\n";
0086   }
0087
0088   BOOST_CHECK_EQUAL(segments.size(), 1);
0089   BOOST_CHECK_CLOSE(segments[0].activation, std::hypot(disp, helper.thickness),
0090                     1.e-8);
0091 }
0092
0093 BOOST_AUTO_TEST_CASE(test_more_tilted_particle) {
0094   Helper helper;
0095
0096   const double disp = 50_um;
0097
0098   Acts::Vector3 hitPosition = Acts::Vector3{10_um, 10_um, 0.0};
0099   Acts::Vector3 hitDirection =
0100       Acts::Vector3{disp, 0.0, helper.thickness}.normalized();
0101
0102   auto segments = helper.channelize(hitPosition, hitDirection);
0103
0104   BOOST_CHECK_EQUAL(segments.size(), 2);
0105   auto sum =
0106       std::accumulate(segments.begin(), segments.end(), 0.0,
0107                       [](double s, auto seg) { return s + seg.activation; });
0108   BOOST_CHECK_CLOSE(sum, std::hypot(disp, helper.thickness), 1.e-8);
0109 }
0110
0111 // This should go directly up on the segment border
0112 BOOST_AUTO_TEST_CASE(test_pathological_upright_particle) {
0113   Helper helper;
0114
0115   Acts::Vector3 hitPosition = Acts::Vector3{0.0, 10_um, 0.0};
0116   Acts::Vector3 hitDirection =
0117       Acts::Vector3{0.0, 0.0, helper.thickness}.normalized();
0118
0119   auto segments = helper.channelize(hitPosition, hitDirection);
0120
0121   BOOST_CHECK_EQUAL(segments.size(), 1);
0122   BOOST_CHECK_CLOSE(segments[0].activation, helper.thickness, 1.e-8);
0123 }
0124
0125 // This should go directly up on the segment border
0126 // TODO why does this does not activate both cells with half of the path???
0127 BOOST_AUTO_TEST_CASE(test_pathological_tilted_particle) {
0128   Helper helper;
0129
0130   double disp = 2.0_um;
0131
0132   Acts::Vector3 hitPosition = Acts::Vector3{-0.5 * disp, 10_um, 0.0};
0133   Acts::Vector3 hitDirection =
0134       Acts::Vector3{disp, 0.0, helper.thickness}.normalized();
0135
0136   auto segments = helper.channelize(hitPosition, hitDirection);
0137
0138   std::cout << "Segments:\n";
0139   for (const auto &seg : segments) {
0140     std::cout << " - (" << seg.bin[0] << ", " << seg.bin[1]
0141               << "), activation: " << seg.activation << "\n";
0142   }
0143
0144   BOOST_CHECK_EQUAL(segments.size(), 2);
0145   auto sum =
0146       std::accumulate(segments.begin(), segments.end(), 0.0,
0147                       [](double s, auto seg) { return s + seg.activation; });
0148   BOOST_CHECK_CLOSE(sum, std::hypot(disp, helper.thickness), 1.e-8);
0149 }