Physics/ElectroMagnetic/BetheHeitler.hpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2018-2021 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 http://mozilla.org/MPL/2.0/.
0008
0009 #pragma once
0010
0011 #include "Acts/Material/MaterialSlab.hpp"
0012 #include "Acts/Utilities/UnitVectors.hpp"
0013 #include "ActsFatras/EventData/Particle.hpp"
0014
0015 #include <array>
0016 #include <cmath>
0017 #include <random>
0018
0019 namespace ActsFatras {
0020
0021 /// Simulate electron energy loss using the Bethe-Heitler description.
0022 ///
0023 /// Bethe-Heitler for electron bremsstrahlung description as described here:
0024 /// "A Gaussian-mixture approximation of the Bethe–Heitler model of electron
0025 /// energy loss by bremsstrahlung" R. Frühwirth
0026 struct BetheHeitler {
0027   using Scalar = Particle::Scalar;
0028   using Vector3 = Particle::Vector3;
0029
0030   /// A scaling factor to
0031   double scaleFactor = 1.;
0032
0033   // Simplified angle evaluation
0034   bool uniformHertzDipoleAngle = false;
0035
0036   /// Simulate the photon emission
0037   ///
0038   /// @param [in] particle The unmodified electron
0039   /// @param [in] gammaE Energy of the photon
0040   /// @param [in] rndPsi Random number for the azimuthal angle
0041   /// @param [in] rndTheta1 Random number for the polar angle
0042   /// @param [in] rndTheta2 Random number for the polar angle
0043   /// @param [in] rndTheta3 Random number for the polar angle
0044   Particle bremPhoton(const Particle &particle, Scalar gammaE, Scalar rndPsi,
0045                       Scalar rndTheta1, Scalar rndTheta2,
0046                       Scalar rndTheta3) const;
0047
0048   /// Simulate energy loss and update the particle parameters.
0049   ///
0050   /// @param[in]     generator is the random number generator
0051   /// @param[in]     slab      defines the passed material
0052   /// @param[in,out] particle  is the particle being updated
0053   /// @return Produced photon.
0054   ///
0055   /// @tparam generator_t is a RandomNumberEngine
0056   template <typename generator_t>
0057   std::array<Particle, 1> operator()(generator_t &generator,
0058                                      const Acts::MaterialSlab &slab,
0059                                      Particle &particle) const {
0060     // Take a random gamma-distributed value - depending on t/X0
0061     std::gamma_distribution<double> gDist(slab.thicknessInX0() / std::log(2.0),
0062                                           1.0);
0063
0064     const auto u = gDist(generator);
0065     const auto z = std::exp(-u);  // MARK: fpeMask(FLTUND, 1, #2346)
0066     const auto sampledEnergyLoss =
0067         std::abs(scaleFactor * particle.energy() * (z - 1.));
0068
0069     std::uniform_real_distribution<Scalar> uDist(0., 1.);
0070     // Build the produced photon
0071     Particle photon =
0072         bremPhoton(particle, sampledEnergyLoss, uDist(generator),
0073                    uDist(generator), uDist(generator), uDist(generator));
0074     // Recoil input momentum
0075     particle.setDirection(particle.direction() * particle.absoluteMomentum() -
0076                           photon.energy() * photon.direction());
0077
0078     // apply the energy loss
0079     particle.correctEnergy(-sampledEnergyLoss);
0080
0081     return {photon};
0082   }
0083 };
0084
0085 }  // namespace ActsFatras