Acts/Material/PropagatorMaterialAssigner.hpp

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 #pragma once
0010
0011 #include "Acts/Geometry/GeometryContext.hpp"
0012 #include "Acts/MagneticField/MagneticFieldContext.hpp"
0013 #include "Acts/Material/interface/IAssignmentFinder.hpp"
0014 #include "Acts/Propagator/ActorList.hpp"
0015 #include "Acts/Propagator/SurfaceCollector.hpp"
0016 #include "Acts/Surfaces/Surface.hpp"
0017 #include "Acts/Utilities/Logger.hpp"
0018 #include "Acts/Utilities/VectorHelpers.hpp"
0019
0020 #include <map>
0021 #include <utility>
0022 #include <vector>
0023
0024 namespace Acts {
0025
0026 /// @brief selector for finding surface
0027 ///
0028 /// This is to be used with a SurfaceCollector<>
0029 struct MaterialSurfaceIdentifier {
0030   /// check if the surface has material
0031   /// @param sf Surface to check for material
0032   /// @return True if the surface has material assigned to it
0033   bool operator()(const Surface& sf) const {
0034     return (sf.surfaceMaterial() != nullptr);
0035   }
0036 };
0037
0038 /// An Interaction volume collector with unique counting
0039 struct InteractionVolumeCollector {
0040   /// Simple result struct that holds the collected volumes
0041   ///
0042   /// @note the map is to avoid double counting as this is
0043   /// called in an action list
0044   struct this_result {
0045     /// Map of collected volume assignments by geometry identifier
0046     std::map<GeometryIdentifier, IAssignmentFinder::VolumeAssignment> collected;
0047   };
0048
0049   /// Type alias for volume collection result
0050   using result_type = this_result;
0051
0052   /// Collector action for the ActionList of the Propagator
0053   /// It checks if the propagator state has a current volume,
0054   /// in which case the action is performed:
0055   /// - it records the volume given the configuration
0056   ///
0057   /// @tparam propagator_state_t is the type of Propagator state
0058   /// @tparam stepper_t Type of the stepper used for the propagation
0059   /// @tparam navigator_t Type of the navigator used for the propagation
0060   ///
0061   /// @param [in,out] state is the mutable stepper state object
0062   /// @param [in] stepper The stepper in use
0063   /// @param [in] navigator The navigator in use
0064   /// @param [in,out] result is the mutable result object
0065   /// @return Result object indicating success or failure
0066   template <typename propagator_state_t, typename stepper_t,
0067             typename navigator_t>
0068   Result<void> act(propagator_state_t& state, const stepper_t& stepper,
0069                    const navigator_t& navigator, result_type& result,
0070                    const Logger& /*logger*/) const {
0071     // Retrieve the current volume
0072     auto currentVolume = navigator.currentVolume(state.navigation);
0073
0074     // The current volume has been assigned by the navigator
0075     if (currentVolume != nullptr) {
0076       auto collIt = result.collected.find(currentVolume->geometryId());
0077       // Check if the volume has been collected and if it has material
0078       if (collIt == result.collected.end() &&
0079           currentVolume->volumeMaterial() != nullptr) {
0080         Vector3 entryPosition = stepper.position(state.stepping);
0081         Vector3 exitPosition = entryPosition;
0082         IAssignmentFinder::VolumeAssignment vAssignment{
0083             InteractionVolume(currentVolume), entryPosition, exitPosition};
0084         result.collected.emplace(currentVolume->geometryId(), vAssignment);
0085       } else if (collIt != result.collected.end()) {
0086         // Update the exit position
0087         (collIt->second).exit = stepper.position(state.stepping);
0088       }
0089     }
0090     return Result<void>::success();
0091   }
0092 };
0093
0094 /// @class PropagatorMaterialAssigner
0095 ///
0096 /// A Propagator based material assigner that uses the navigation and
0097 /// transport of the propagator to assign the material to the surface
0098 /// or the volume.
0099 ///
0100 /// @note eventual navigation problems would affect he material mapping
0101 template <typename propagator_t>
0102 class PropagatorMaterialAssigner final : public IAssignmentFinder {
0103  public:
0104   /// @brief  Construct with propagator
0105   /// @param propagator
0106   explicit PropagatorMaterialAssigner(propagator_t propagator)
0107       : m_propagator(std::move(propagator)) {}
0108
0109   /// @brief Method for generating assignment candidates for the
0110   /// material interaction assignment to surfaces or volumes
0111   ///
0112   /// @param gctx is the geometry context
0113   /// @param mctx is the magnetic field context
0114   /// @param position is the position of the initial ray
0115   /// @param direction is the direction of initial ray
0116   ///
0117   /// @return a vector of Surface Assignments and Volume Assignments
0118   std::pair<std::vector<IAssignmentFinder::SurfaceAssignment>,
0119             std::vector<IAssignmentFinder::VolumeAssignment>>
0120   assignmentCandidates(const GeometryContext& gctx,
0121                        const MagneticFieldContext& mctx,
0122                        const Vector3& position,
0123                        const Vector3& direction) const final {
0124     // Return container
0125     std::pair<std::vector<IAssignmentFinder::SurfaceAssignment>,
0126               std::vector<IAssignmentFinder::VolumeAssignment>>
0127         candidates;
0128
0129     using VectorHelpers::makeVector4;
0130     // Neutral curvilinear parameters
0131     BoundTrackParameters start = BoundTrackParameters::createCurvilinear(
0132         makeVector4(position, 0), direction, 1, std::nullopt,
0133         ParticleHypothesis::geantino());
0134
0135     // Prepare Action list and abort list
0136     using MaterialSurfaceCollector =
0137         SurfaceCollector<MaterialSurfaceIdentifier>;
0138     using ActorList = ActorList<MaterialSurfaceCollector,
0139                                 InteractionVolumeCollector, EndOfWorldReached>;
0140     using PropagatorOptions =
0141         typename propagator_t::template Options<ActorList>;
0142
0143     PropagatorOptions options(gctx, mctx);
0144
0145     const auto& result = m_propagator.propagate(start, options).value();
0146
0147     // The surface collection results
0148     auto scResult =
0149         result.template get<MaterialSurfaceCollector::result_type>();
0150     auto mappingSurfaces = scResult.collected;
0151
0152     for (auto& mSurface : mappingSurfaces) {
0153       candidates.first.push_back(IAssignmentFinder::SurfaceAssignment{
0154           mSurface.surface, mSurface.position, direction});
0155     }
0156
0157     // The volume collection results
0158     auto vcResult =
0159         result.template get<InteractionVolumeCollector::result_type>();
0160     for (const auto& [geoId, vIntersection] : vcResult.collected) {
0161       candidates.second.push_back(vIntersection);
0162     }
0163
0164     // Return the candidates
0165     return candidates;
0166   }
0167
0168  private:
0169   propagator_t m_propagator;
0170 };
0171
0172 }  // namespace Acts