Acts/Propagator/MaterialInteractor.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/Definitions/Units.hpp"
0012 #include "Acts/Geometry/TrackingVolume.hpp"
0013 #include "Acts/Material/MaterialInteraction.hpp"
0014 #include "Acts/Material/MaterialSlab.hpp"
0015 #include "Acts/Propagator/Propagator.hpp"
0016 #include "Acts/Propagator/detail/PointwiseMaterialInteraction.hpp"
0017 #include "Acts/Propagator/detail/VolumeMaterialInteraction.hpp"
0018 #include "Acts/Surfaces/Surface.hpp"
0019
0020 #include <sstream>
0021
0022 namespace Acts {
0023
0024 /// Material interactor propagator action.
0025 ///
0026 /// Apply material interactions at a surface and update the track state.
0027 struct MaterialInteractor {
0028   /// Whether to consider multiple scattering.
0029   bool multipleScattering = true;
0030   /// Whether to consider energy loss.
0031   bool energyLoss = true;
0032   /// Whether to record all material interactions.
0033   bool recordInteractions = false;
0034   /// Whether to add or remove noise.
0035   NoiseUpdateMode noiseUpdateMode = NoiseUpdateMode::addNoise;
0036
0037   /// Type alias for material interaction result
0038   using result_type = RecordedMaterial;
0039
0040   /// @brief Interaction with detector material for the ActionList
0041   /// of the Propagator
0042   ///
0043   /// It checks if the state has a current surface, in which case
0044   /// the action is performed: the covariance is transported to the position,
0045   /// multiple scattering and energy loss is applied  according to the
0046   /// configuration.
0047   ///
0048   /// @tparam propagator_state_t is the type of Propagator state
0049   /// @tparam stepper_t Type of the stepper of the propagation
0050   /// @tparam navigator_t Type of the navigator of the propagation
0051   ///
0052   /// @param state is the mutable propagator state object
0053   /// @param stepper The stepper in use
0054   /// @param navigator The navigator in use
0055   /// @param result is the mutable result state object
0056   /// @param logger a logger instance
0057   template <typename propagator_state_t, typename stepper_t,
0058             typename navigator_t>
0059   void act(propagator_state_t& state, const stepper_t& stepper,
0060            const navigator_t& navigator, result_type& result,
0061            const Logger& logger) const {
0062     if (state.stage == PropagatorStage::postPropagation) {
0063       return;
0064     }
0065
0066     // Do nothing if nothing is what is requested.
0067     if (!(multipleScattering || energyLoss || recordInteractions)) {
0068       return;
0069     }
0070
0071     // Handle surface material
0072
0073     // Note that start and target surface conditions are handled in the
0074     // interaction code
0075     const Surface* surface = navigator.currentSurface(state.navigation);
0076
0077     // We only have material interactions if there is potential material
0078     if (surface && surface->surfaceMaterial()) {
0079       ACTS_VERBOSE("MaterialInteractor | " << "Found material on surface "
0080                                            << surface->geometryId());
0081
0082       // Prepare relevant input particle properties
0083       detail::PointwiseMaterialInteraction interaction(surface, state, stepper);
0084
0085       // Determine the effective traversed material and its properties
0086       // Material exists but it's not real, i.e. vacuum; there is nothing to do
0087       if (interaction.evaluateMaterialSlab(state, navigator)) {
0088         // Evaluate the material effects
0089         interaction.evaluatePointwiseMaterialInteraction(multipleScattering,
0090                                                          energyLoss);
0091
0092         if (energyLoss) {
0093           using namespace UnitLiterals;
0094           ACTS_VERBOSE("MaterialInteractor | "
0095                        << interaction.slab << " absPdg=" << interaction.absPdg
0096                        << " mass=" << interaction.mass / 1_MeV << "MeV"
0097                        << " momentum=" << interaction.momentum / 1_GeV << "GeV"
0098                        << " energyloss=" << interaction.Eloss / 1_MeV << "MeV");
0099         }
0100
0101         // To integrate process noise, we need to transport
0102         // the covariance to the current position in space
0103         if (interaction.performCovarianceTransport) {
0104           stepper.transportCovarianceToCurvilinear(state.stepping);
0105         }
0106         // Apply the material interactions
0107         interaction.updateState(state, stepper, noiseUpdateMode);
0108
0109         // Record the result
0110         recordResult(interaction, result);
0111       }
0112     }
0113
0114     // Handle volume material
0115
0116     // In case of Volume material update the result of the previous step
0117     if (!result.materialInteractions.empty() &&
0118         !result.materialInteractions.back().volume.empty() &&
0119         result.materialInteractions.back().updatedVolumeStep == false) {
0120       updateResult(state, stepper, result);
0121     }
0122
0123     auto volume = navigator.currentVolume(state.navigation);
0124
0125     // We only have material interactions if there is potential material
0126     if (volume && volume->volumeMaterial()) {
0127       ACTS_VERBOSE("MaterialInteractor | " << "Found material in volume "
0128                                            << volume->geometryId());
0129
0130       // Prepare relevant input particle properties
0131       detail::VolumeMaterialInteraction interaction(volume, state, stepper);
0132       // Determine the effective traversed material and its properties
0133       // Material exists but it's not real, i.e. vacuum; there is nothing to do
0134       if (interaction.evaluateMaterialSlab(state, navigator)) {
0135         // Record the result
0136         recordResult(interaction, result);
0137       }
0138     }
0139   }
0140
0141  private:
0142   /// @brief This function records the material effect
0143   ///
0144   /// @param [in] interaction Interaction cache container
0145   /// @param [in, out] result Result storage
0146   void recordResult(const detail::PointwiseMaterialInteraction& interaction,
0147                     result_type& result) const {
0148     result.materialInX0 += interaction.slab.thicknessInX0();
0149     result.materialInL0 += interaction.slab.thicknessInL0();
0150
0151     // Record the interaction if requested
0152     if (!recordInteractions) {
0153       return;
0154     }
0155
0156     MaterialInteraction mi;
0157     mi.position = interaction.pos;
0158     mi.time = interaction.time;
0159     mi.direction = interaction.dir;
0160     mi.deltaP = interaction.nextP - interaction.momentum;
0161     mi.sigmaPhi2 = interaction.variancePhi;
0162     mi.sigmaTheta2 = interaction.varianceTheta;
0163     mi.sigmaQoP2 = interaction.varianceQoverP;
0164     mi.surface = interaction.surface;
0165     mi.volume = InteractionVolume();
0166     mi.pathCorrection = interaction.pathCorrection;
0167     mi.materialSlab = interaction.slab;
0168     result.materialInteractions.push_back(std::move(mi));
0169   }
0170
0171   /// @brief This function records the material effect
0172   ///
0173   /// @param [in] interaction Interaction cache container
0174   /// @param [in, out] result Result storage
0175   void recordResult(const detail::VolumeMaterialInteraction& interaction,
0176                     result_type& result) const {
0177     // Record the interaction if requested
0178     if (!recordInteractions) {
0179       return;
0180     }
0181
0182     MaterialInteraction mi;
0183     mi.position = interaction.pos;
0184     mi.time = interaction.time;
0185     mi.direction = interaction.dir;
0186     mi.surface = nullptr;
0187     mi.volume = interaction.volume;
0188     mi.pathCorrection = interaction.pathCorrection;
0189     mi.materialSlab = interaction.slab;
0190     result.materialInteractions.push_back(std::move(mi));
0191   }
0192
0193   /// @brief This function update the previous material step
0194   ///
0195   /// @param [in,out] state The state object
0196   /// @param [in] stepper The stepper instance
0197   /// @param [in, out] result Result storage
0198   template <typename propagator_state_t, typename stepper_t>
0199   void updateResult(propagator_state_t& state, const stepper_t& stepper,
0200                     result_type& result) const {
0201     // Update the previous interaction if requested
0202     if (!recordInteractions) {
0203       return;
0204     }
0205
0206     Vector3 shift = stepper.position(state.stepping) -
0207                     result.materialInteractions.back().position;
0208     double momentum = stepper.direction(state.stepping).norm();
0209     result.materialInteractions.back().deltaP =
0210         momentum - result.materialInteractions.back().direction.norm();
0211     result.materialInteractions.back().materialSlab.scaleThickness(
0212         shift.norm());
0213     result.materialInteractions.back().updatedVolumeStep = true;
0214     result.materialInX0 +=
0215         result.materialInteractions.back().materialSlab.thicknessInX0();
0216     result.materialInL0 +=
0217         result.materialInteractions.back().materialSlab.thicknessInL0();
0218   }
0219 };
0220
0221 }  // namespace Acts