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