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().withSensitive(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().withSensitive(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().withSensitive(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 = {
0187           {GeometryIdentifier().withSensitive(2), VetoThisOne{}}};
0188   GeometryHierarchyMap<MaterialInteractionAssignment::LocalVeto> localVetos(
0189       localVetoVector);
0190   MaterialInteractionAssignment::Options options;
0191   options.localVetos = localVetos;
0192
0193   // Assign the material interactions to the surface hits
0194   auto [assigned, unassigned, surfacesLeft] =
0195       MaterialInteractionAssignment::assign(tContext, materialInteractions,
0196                                             intersectedSurfaces, options);
0197
0198   // Check that the material interaction was assigned
0199   BOOST_CHECK_EQUAL(assigned.size(), 34u);
0200   BOOST_CHECK_EQUAL(unassigned.size(), 15u);
0201   BOOST_CHECK_EQUAL(surfacesLeft.size(), 1u);
0202 }
0203
0204 BOOST_AUTO_TEST_CASE(AssignToClosest_withReassignment) {
0205   // Create a vector of surfaces
0206   std::vector<std::shared_ptr<Surface>> surfaces = {
0207       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 20.0, 100.0),
0208       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 30.0, 100.0),
0209       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 50.0,
0210                                            100.0)};
0211
0212   for (auto [is, surface] : enumerate(surfaces)) {
0213     surface->assignGeometryId(GeometryIdentifier().withSensitive(is + 1));
0214   }
0215
0216   std::vector<IAssignmentFinder::SurfaceAssignment> intersectedSurfaces = {
0217       {surfaces[0].get(), {20., 0., 0.}, {1., 0., 0.}},
0218       {surfaces[1].get(), {30., 0., 0.}, {1., 0., 0.}},
0219       {surfaces[2].get(), {50., 0., 0.}, {1., 0., 0.}}};
0220
0221   // Create a material
0222   Material material = Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0);
0223
0224   std::vector<MaterialInteraction> materialInteractions;
0225   materialInteractions.reserve(50);
0226   // Generate some material interactions
0227   for (unsigned int i = 1; i < 50; ++i) {
0228     MaterialInteraction materialInteraction;
0229     materialInteraction.materialSlab = MaterialSlab(material, 0.1);
0230     materialInteraction.position = Vector3{i * 1.0, 0, 0};
0231     materialInteraction.direction = Vector3{1.0, 0., 0.};
0232     materialInteractions.push_back(materialInteraction);
0233   }
0234
0235   // Veto in a specific one
0236   struct ReAssignToNeighbor {
0237     void operator()(
0238         MaterialInteraction& m,
0239         const IAssignmentFinder::SurfaceAssignment& /*suggestedAssignment*/,
0240         const IAssignmentFinder::SurfaceAssignment& n) const {
0241       auto [surface, position, direction] = n;
0242       m.surface = surface;
0243       m.position = position;
0244       m.direction = direction;
0245       m.intersectionID = surface->geometryId();
0246       return;
0247     }
0248   };
0249
0250   // We assign this to
0251   std::vector<std::pair<GeometryIdentifier,
0252                         MaterialInteractionAssignment::ReAssignment>>
0253       reassignmentVector = {
0254           {GeometryIdentifier().withSensitive(2), ReAssignToNeighbor{}}};
0255   GeometryHierarchyMap<MaterialInteractionAssignment::ReAssignment>
0256       reassignments(reassignmentVector);
0257   MaterialInteractionAssignment::Options options;
0258   options.reAssignments = reassignments;
0259
0260   // Assign the material interactions to the surface hits
0261   auto [assigned, unassigned, surfaceLeft] =
0262       MaterialInteractionAssignment::assign(tContext, materialInteractions,
0263                                             intersectedSurfaces, options);
0264
0265   // Check that the material interaction was assigned
0266   BOOST_CHECK_EQUAL(assigned.size(), 49u);
0267   BOOST_CHECK_EQUAL(unassigned.size(), 0u);
0268   BOOST_CHECK_EQUAL(surfaceLeft.size(), 1u);
0269
0270   // Check that the geoid with number 2 never shows up
0271   for (const auto& mi : assigned) {
0272     BOOST_CHECK_NE(mi.intersectionID, GeometryIdentifier().withSensitive(2));
0273   }
0274 }
0275
0276 BOOST_AUTO_TEST_CASE(AssignWithPathLength) {
0277   auto surface =
0278       Surface::makeShared<CylinderSurface>(Transform3::Identity(), 20.0, 100.0);
0279   surface->assignGeometryId(GeometryIdentifier().withSensitive(1));
0280
0281   Vector3 position = {20., 10., 0.};
0282   Vector3 direction = position.normalized();
0283
0284   IAssignmentFinder::SurfaceAssignment surfaceHit{surface.get(), position,
0285                                                   direction};
0286
0287   Material material = Material::fromMolarDensity(1.0, 2.0, 3.0, 4.0, 5.0);
0288
0289   MaterialInteraction materialInteraction;
0290   materialInteraction.materialSlab = MaterialSlab(material, 0.1);
0291   materialInteraction.position = position + 0.5 * direction;
0292   materialInteraction.direction = direction;
0293
0294   MaterialInteractionAssignment::Options options;
0295
0296   auto [assigned, unassigned, surfaceLeft] =
0297       MaterialInteractionAssignment::assign(tContext, {materialInteraction},
0298                                             {surfaceHit}, options);
0299
0300   // Check that the material interaction was assigned
0301   BOOST_CHECK_EQUAL(assigned.size(), 1u);
0302   BOOST_CHECK_EQUAL(unassigned.size(), 0u);
0303   BOOST_CHECK_EQUAL(surfaceLeft.size(), 0u);
0304
0305   // Check that the path correction is set
0306   BOOST_CHECK_NE(assigned[0].pathCorrection, 0.);
0307 }
0308
0309 BOOST_AUTO_TEST_SUITE_END()
0310
0311 }  // namespace Acts::Test