Core/Material/MaterialInteractionAssignmentTests.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/Algebra.hpp"
0012 #include "Acts/Geometry/GeometryContext.hpp"
0013 #include "Acts/Geometry/GeometryIdentifier.hpp"
0014 #include "Acts/Material/Material.hpp"
0015 #include "Acts/Material/MaterialInteraction.hpp"
0016 #include "Acts/Material/MaterialInteractionAssignment.hpp"
0017 #include "Acts/Material/MaterialSlab.hpp"
0018 #include "Acts/Material/interface/IAssignmentFinder.hpp"
0019 #include "Acts/Propagator/SurfaceCollector.hpp"
0020 #include "Acts/Surfaces/CylinderSurface.hpp"
0021 #include "Acts/Utilities/Enumerate.hpp"
0022 #include "Acts/Utilities/VectorHelpers.hpp"
0023
0024 #include <limits>
0025
0026 namespace Acts::Test {
0027
0028 auto tContext = GeometryContext();
0029
0030 BOOST_AUTO_TEST_SUITE(MaterialInteractionAssignmentSuite)
0031
0032 BOOST_AUTO_TEST_CASE(AssignToClosest) {
0033   // Create a vector of surfaces
0034   std::vector<std::shared_ptr<Surface>> surfaces = {
0035       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 20.0, 100.0),
0036       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 30.0, 100.0),
0037       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 50.0,
0038                                            100.0)};
0039
0040   for (auto [is, surface] : enumerate(surfaces)) {
0041     surface->assignGeometryId(GeometryIdentifier().setSensitive(is + 1));
0042   }
0043
0044   std::vector<IAssignmentFinder::SurfaceAssignment> intersectedSurfaces = {
0045       {surfaces[0].get(), {20., 0., 0.}, {1., 0., 0.}},
0046       {surfaces[1].get(), {30., 0., 0.}, {1., 0., 0.}},
0047       {surfaces[2].get(), {50., 0., 0.}, {1., 0., 0.}}};
0048
0049   // Create a material
0050   Material material = Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0);
0051
0052   std::vector<MaterialInteraction> materialInteractions;
0053   materialInteractions.reserve(50);
0054   // Generate some material interactions
0055   for (unsigned int i = 1; i < 50; ++i) {
0056     MaterialInteraction materialInteraction;
0057     materialInteraction.materialSlab = MaterialSlab(material, 0.1);
0058     materialInteraction.position = Vector3{i * 1.0, 0, 0};
0059     materialInteraction.direction = Vector3{1.0, 0., 0.};
0060     materialInteractions.push_back(materialInteraction);
0061   }
0062
0063   MaterialInteractionAssignment::Options options;
0064
0065   // Assign the material interactions to the surface hits
0066   auto [assigned, unassigned, surfacesLeft] =
0067       MaterialInteractionAssignment::assign(tContext, materialInteractions,
0068                                             intersectedSurfaces, options);
0069
0070   // Check that the material interaction was assigned
0071   BOOST_CHECK_EQUAL(assigned.size(), materialInteractions.size());
0072   BOOST_CHECK_EQUAL(unassigned.size(), 0u);
0073   BOOST_CHECK_EQUAL(surfacesLeft.size(), 0u);
0074
0075   // Check that it is assigned to the closest surface always
0076   for (const auto& mi : assigned) {
0077     double minDistance = std::numeric_limits<double>::max();
0078     const Surface* closestSurface = nullptr;
0079     for (const auto& [surface, position, direction] : intersectedSurfaces) {
0080       double distance = (mi.position - position).norm();
0081       if (distance < minDistance) {
0082         minDistance = distance;
0083         closestSurface = surface;
0084       }
0085     }
0086     BOOST_CHECK_EQUAL(mi.surface, closestSurface);
0087   }
0088 }
0089
0090 BOOST_AUTO_TEST_CASE(AssignToClosest_withGlobalVeto) {
0091   // Create a vector of surfaces
0092   std::vector<std::shared_ptr<Surface>> surfaces = {
0093       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 20.0, 100.0),
0094       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 30.0, 100.0),
0095       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 50.0,
0096                                            100.0)};
0097
0098   for (auto [is, surface] : enumerate(surfaces)) {
0099     surface->assignGeometryId(GeometryIdentifier().setSensitive(is + 1));
0100   }
0101
0102   std::vector<IAssignmentFinder::SurfaceAssignment> intersectedSurfaces = {
0103       {surfaces[0].get(), {20., 0., 0.}, {1., 0., 0.}},
0104       {surfaces[1].get(), {30., 0., 0.}, {1., 0., 0.}},
0105       {surfaces[2].get(), {50., 0., 0.}, {1., 0., 0.}}};
0106
0107   // Create a material
0108   Material material = Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0);
0109
0110   std::vector<MaterialInteraction> materialInteractions;
0111   materialInteractions.reserve(50);
0112   // Generate some material interactions
0113   for (unsigned int i = 1; i < 50; ++i) {
0114     MaterialInteraction materialInteraction;
0115     materialInteraction.materialSlab = MaterialSlab(material, 0.1);
0116     materialInteraction.position = Vector3{i * 1.0, 0, 0};
0117     materialInteraction.direction = Vector3{1.0, 0., 0.};
0118     materialInteractions.push_back(materialInteraction);
0119   }
0120
0121   // Veto everything above 40 mm
0122   struct RadialVeto {
0123     double rMax = 40.0;
0124     bool operator()(const MaterialInteraction& mi) const {
0125       return VectorHelpers::perp(mi.position) > rMax;
0126     }
0127   };
0128   MaterialInteractionAssignment::Options options;
0129   options.globalVetos.push_back(RadialVeto{40});
0130
0131   // Assign the material interactions to the surface hits
0132   auto [assigned, unassigned, surfacesLeft] =
0133       MaterialInteractionAssignment::assign(tContext, materialInteractions,
0134                                             intersectedSurfaces, options);
0135
0136   // Check that the material interaction was assigned
0137   BOOST_CHECK_EQUAL(assigned.size(), 40u);
0138   BOOST_CHECK_EQUAL(unassigned.size(), 9u);
0139   BOOST_CHECK_EQUAL(surfacesLeft.size(), 1u);
0140 }
0141
0142 BOOST_AUTO_TEST_CASE(AssignToClosest_withLocalVeto) {
0143   // Create a vector of surfaces
0144   std::vector<std::shared_ptr<Surface>> surfaces = {
0145       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 20.0, 100.0),
0146       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 30.0, 100.0),
0147       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 50.0,
0148                                            100.0)};
0149
0150   for (auto [is, surface] : enumerate(surfaces)) {
0151     surface->assignGeometryId(GeometryIdentifier().setSensitive(is + 1));
0152   }
0153
0154   std::vector<IAssignmentFinder::SurfaceAssignment> intersectedSurfaces = {
0155       {surfaces[0].get(), {20., 0., 0.}, {1., 0., 0.}},
0156       {surfaces[1].get(), {30., 0., 0.}, {1., 0., 0.}},
0157       {surfaces[2].get(), {50., 0., 0.}, {1., 0., 0.}}};
0158
0159   // Create a material
0160   Material material = Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0);
0161
0162   std::vector<MaterialInteraction> materialInteractions;
0163   materialInteractions.reserve(50);
0164   // Generate some material interactions
0165   for (unsigned int i = 1; i < 50; ++i) {
0166     MaterialInteraction materialInteraction;
0167     materialInteraction.materialSlab = MaterialSlab(material, 0.1);
0168     materialInteraction.position = Vector3{i * 1.0, 0, 0};
0169     materialInteraction.direction = Vector3{1.0, 0., 0.};
0170     materialInteractions.push_back(materialInteraction);
0171   }
0172
0173   // Veto in a specific one
0174   struct VetoThisOne {
0175     bool operator()(
0176         const MaterialInteraction& /*m*/,
0177         const IAssignmentFinder::SurfaceAssignment& /*suggestedAssignment*/)
0178         const {
0179       return true;
0180     }
0181   };
0182
0183   // We assign this to
0184   std::vector<
0185       std::pair<GeometryIdentifier, MaterialInteractionAssignment::LocalVeto>>
0186       localVetoVector = {{GeometryIdentifier().setSensitive(2), VetoThisOne{}}};
0187   GeometryHierarchyMap<MaterialInteractionAssignment::LocalVeto> localVetos(
0188       localVetoVector);
0189   MaterialInteractionAssignment::Options options;
0190   options.localVetos = localVetos;
0191
0192   // Assign the material interactions to the surface hits
0193   auto [assigned, unassigned, surfacesLeft] =
0194       MaterialInteractionAssignment::assign(tContext, materialInteractions,
0195                                             intersectedSurfaces, options);
0196
0197   // Check that the material interaction was assigned
0198   BOOST_CHECK_EQUAL(assigned.size(), 34u);
0199   BOOST_CHECK_EQUAL(unassigned.size(), 15u);
0200   BOOST_CHECK_EQUAL(surfacesLeft.size(), 1u);
0201 }
0202
0203 BOOST_AUTO_TEST_CASE(AssignToClosest_withReassignment) {
0204   // Create a vector of surfaces
0205   std::vector<std::shared_ptr<Surface>> surfaces = {
0206       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 20.0, 100.0),
0207       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 30.0, 100.0),
0208       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 50.0,
0209                                            100.0)};
0210
0211   for (auto [is, surface] : enumerate(surfaces)) {
0212     surface->assignGeometryId(GeometryIdentifier().setSensitive(is + 1));
0213   }
0214
0215   std::vector<IAssignmentFinder::SurfaceAssignment> intersectedSurfaces = {
0216       {surfaces[0].get(), {20., 0., 0.}, {1., 0., 0.}},
0217       {surfaces[1].get(), {30., 0., 0.}, {1., 0., 0.}},
0218       {surfaces[2].get(), {50., 0., 0.}, {1., 0., 0.}}};
0219
0220   // Create a material
0221   Material material = Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0);
0222
0223   std::vector<MaterialInteraction> materialInteractions;
0224   materialInteractions.reserve(50);
0225   // Generate some material interactions
0226   for (unsigned int i = 1; i < 50; ++i) {
0227     MaterialInteraction materialInteraction;
0228     materialInteraction.materialSlab = MaterialSlab(material, 0.1);
0229     materialInteraction.position = Vector3{i * 1.0, 0, 0};
0230     materialInteraction.direction = Vector3{1.0, 0., 0.};
0231     materialInteractions.push_back(materialInteraction);
0232   }
0233
0234   // Veto in a specific one
0235   struct ReAssignToNeighbor {
0236     void operator()(
0237         MaterialInteraction& m,
0238         const IAssignmentFinder::SurfaceAssignment& /*suggestedAssignment*/,
0239         const IAssignmentFinder::SurfaceAssignment& n) const {
0240       auto [surface, position, direction] = n;
0241       m.surface = surface;
0242       m.position = position;
0243       m.direction = direction;
0244       m.intersectionID = surface->geometryId();
0245       return;
0246     }
0247   };
0248
0249   // We assign this to
0250   std::vector<std::pair<GeometryIdentifier,
0251                         MaterialInteractionAssignment::ReAssignment>>
0252       reassignmentVector = {
0253           {GeometryIdentifier().setSensitive(2), ReAssignToNeighbor{}}};
0254   GeometryHierarchyMap<MaterialInteractionAssignment::ReAssignment>
0255       reassignments(reassignmentVector);
0256   MaterialInteractionAssignment::Options options;
0257   options.reAssignments = reassignments;
0258
0259   // Assign the material interactions to the surface hits
0260   auto [assigned, unassigned, surfaceLeft] =
0261       MaterialInteractionAssignment::assign(tContext, materialInteractions,
0262                                             intersectedSurfaces, options);
0263
0264   // Check that the material interaction was assigned
0265   BOOST_CHECK_EQUAL(assigned.size(), 49u);
0266   BOOST_CHECK_EQUAL(unassigned.size(), 0u);
0267   BOOST_CHECK_EQUAL(surfaceLeft.size(), 1u);
0268
0269   // Check that the geoid with number 2 never shows up
0270   for (const auto& mi : assigned) {
0271     BOOST_CHECK_NE(mi.intersectionID, GeometryIdentifier().setSensitive(2));
0272   }
0273 }
0274
0275 BOOST_AUTO_TEST_CASE(AssignWithPathLength) {
0276   auto surface =
0277       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 20.0, 100.0);
0278   surface->assignGeometryId(GeometryIdentifier().setSensitive(1));
0279
0280   Vector3 position = {20., 10., 0.};
0281   Vector3 direction = position.normalized();
0282
0283   IAssignmentFinder::SurfaceAssignment surfaceHit{surface.get(), position,
0284                                                   direction};
0285
0286   Material material = Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0);
0287
0288   MaterialInteraction materialInteraction;
0289   materialInteraction.materialSlab = MaterialSlab(material, 0.1);
0290   materialInteraction.position = position + 0.5 * direction;
0291   materialInteraction.direction = direction;
0292
0293   MaterialInteractionAssignment::Options options;
0294
0295   auto [assigned, unassigned, surfaceLeft] =
0296       MaterialInteractionAssignment::assign(tContext, {materialInteraction},
0297                                             {surfaceHit}, options);
0298
0299   // Check that the material interaction was assigned
0300   BOOST_CHECK_EQUAL(assigned.size(), 1u);
0301   BOOST_CHECK_EQUAL(unassigned.size(), 0u);
0302   BOOST_CHECK_EQUAL(surfaceLeft.size(), 0u);
0303
0304   // Check that the path correction is set
0305   BOOST_CHECK_NE(assigned[0].pathCorrection, 0.);
0306 }
0307
0308 BOOST_AUTO_TEST_SUITE_END()
0309
0310 }  // namespace Acts::Test