Fatras/Digitization/ChannelizerTests.cpp

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 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 https://mozilla.org/MPL/2.0/.
0008
0009 #include <boost/test/unit_test.hpp>
0010
0011 #include "Acts/Definitions/Units.hpp"
0012 #include "Acts/Surfaces/CurvilinearSurface.hpp"
0013 #include "Acts/Surfaces/Surface.hpp"
0014 #include "Acts/Utilities/IAxis.hpp"
0015 #include "Acts/Utilities/IMultiAxis.hpp"
0016 #include "ActsFatras/Digitization/Channelizer.hpp"
0017 #include "ActsFatras/Digitization/SurfaceMask.hpp"
0018
0019 #include <memory>
0020 #include <numeric>
0021
0022 using namespace Acts;
0023 using namespace Acts::UnitLiterals;
0024 using namespace ActsFatras;
0025
0026 struct Helper {
0027   std::shared_ptr<Surface> surface;
0028   std::shared_ptr<const IMultiAxis> segmentation;
0029
0030   GeometryContext gctx = GeometryContext::dangerouslyDefaultConstruct();
0031   double thickness = 125_um;
0032   Vector3 driftDir = Vector3::Zero();
0033
0034   ActsFatras::Channelizer channelizer;
0035
0036   Helper() {
0037     surface =
0038         CurvilinearSurface(Vector3::Zero(), Vector3{0.0, 0.0, 1.0}).surface();
0039
0040     float pitchSize = 50_um;
0041     float min = -200_um;
0042     float max = 200_um;
0043     int bins = static_cast<int>((max - min) / pitchSize);
0044     const auto axisX = IAxis::createEquidistant(
0045         AxisBoundaryType::Bound, min, max, bins, AxisDirection::AxisX);
0046     const auto axisY = IAxis::createEquidistant(
0047         AxisBoundaryType::Bound, min, max, bins, AxisDirection::AxisY);
0048     segmentation = IMultiAxis::create(*axisX, *axisY);
0049   }
0050
0051   auto channelize(const Vector3 &pos3, const Vector3 &dir3) const {
0052     Vector4 pos4 = Vector4::Zero();
0053     pos4.segment<3>(ePos0) = pos3;
0054     Vector4 mom4 = Vector4::Zero();
0055     mom4.segment<3>(eMom0) = dir3;
0056     ActsFatras::Hit hit({}, {}, pos4, mom4, mom4);
0057     auto res = channelizer.channelize(hit, *surface, gctx, driftDir,
0058                                       *segmentation, thickness);
0059     BOOST_REQUIRE(res.ok());
0060     return *res;
0061   }
0062 };
0063
0064 namespace ActsTests {
0065
0066 BOOST_AUTO_TEST_SUITE(DigitizationSuite)
0067
0068 BOOST_AUTO_TEST_CASE(test_upright_particle) {
0069   Helper helper;
0070
0071   Vector3 pos3 = Vector3{10_um, 10_um, 0.0};
0072   Vector3 dir3 = Vector3{0.0, 0.0, helper.thickness}.normalized();
0073
0074   auto segments = helper.channelize(pos3, dir3);
0075
0076   BOOST_CHECK_EQUAL(segments.size(), 1);
0077   BOOST_CHECK_CLOSE(segments[0].activation, helper.thickness, 1.e-8);
0078 }
0079
0080 BOOST_AUTO_TEST_CASE(test_tilted_particle) {
0081   Helper helper;
0082
0083   const double disp = 10_um;
0084
0085   Vector3 hitPosition = Vector3{10_um, 10_um, 0.0};
0086   Vector3 hitDirection = Vector3({disp, 0.0, helper.thickness}).normalized();
0087
0088   auto segments = helper.channelize(hitPosition, hitDirection);
0089
0090   std::cout << "Segments:\n";
0091   for (const auto &seg : segments) {
0092     std::cout << " - (" << seg.bin[0] << ", " << seg.bin[1]
0093               << "), activation: " << seg.activation << "\n";
0094   }
0095
0096   BOOST_CHECK_EQUAL(segments.size(), 1);
0097   BOOST_CHECK_CLOSE(segments[0].activation, std::hypot(disp, helper.thickness),
0098                     1.e-8);
0099 }
0100
0101 BOOST_AUTO_TEST_CASE(test_more_tilted_particle) {
0102   Helper helper;
0103
0104   const double disp = 50_um;
0105
0106   Vector3 hitPosition = Vector3{10_um, 10_um, 0.0};
0107   Vector3 hitDirection = Vector3{disp, 0.0, helper.thickness}.normalized();
0108
0109   auto segments = helper.channelize(hitPosition, hitDirection);
0110
0111   BOOST_CHECK_EQUAL(segments.size(), 2);
0112   auto sum =
0113       std::accumulate(segments.begin(), segments.end(), 0.0,
0114                       [](double s, auto seg) { return s + seg.activation; });
0115   BOOST_CHECK_CLOSE(sum, std::hypot(disp, helper.thickness), 1.e-8);
0116 }
0117
0118 // This should go directly up on the segment border
0119 BOOST_AUTO_TEST_CASE(test_pathological_upright_particle) {
0120   Helper helper;
0121
0122   Vector3 hitPosition = Vector3{0.0, 10_um, 0.0};
0123   Vector3 hitDirection = Vector3{0.0, 0.0, helper.thickness}.normalized();
0124
0125   auto segments = helper.channelize(hitPosition, hitDirection);
0126
0127   BOOST_CHECK_EQUAL(segments.size(), 1);
0128   BOOST_CHECK_CLOSE(segments[0].activation, helper.thickness, 1.e-8);
0129 }
0130
0131 // This should go directly up on the segment border
0132 // TODO why does this does not activate both cells with half of the path???
0133 BOOST_AUTO_TEST_CASE(test_pathological_tilted_particle) {
0134   Helper helper;
0135
0136   double disp = 2.0_um;
0137
0138   Vector3 hitPosition = Vector3{-0.5 * disp, 10_um, 0.0};
0139   Vector3 hitDirection = Vector3{disp, 0.0, helper.thickness}.normalized();
0140
0141   auto segments = helper.channelize(hitPosition, hitDirection);
0142
0143   std::cout << "Segments:\n";
0144   for (const auto &seg : segments) {
0145     std::cout << " - (" << seg.bin[0] << ", " << seg.bin[1]
0146               << "), activation: " << seg.activation << "\n";
0147   }
0148
0149   BOOST_CHECK_EQUAL(segments.size(), 2);
0150   auto sum =
0151       std::accumulate(segments.begin(), segments.end(), 0.0,
0152                       [](double s, auto seg) { return s + seg.activation; });
0153   BOOST_CHECK_CLOSE(sum, std::hypot(disp, helper.thickness), 1.e-8);
0154 }
0155
0156 BOOST_AUTO_TEST_SUITE_END()
0157
0158 }  // namespace ActsTests