Namespace ActsFatras

namespace ActsFatras

Typedefs

using Cell = std::pair<unsigned int, Acts::ActsScalar>

A single cell definition: index, cell central value.

using CombineAnd = detail::CombineSelectors<true, std::logical_and<bool>, selectors_t...>

Select objects that fullfil all selectors.

using CombineOr = detail::CombineSelectors<false, std::logical_or<bool>, selectors_t...>

Select objects that fullfil at least one selector.

using GaussianMixtureScattering = detail::ScatteringImpl<detail::GaussianMixture>
using GeneralMixtureScattering = detail::ScatteringImpl<detail::GeneralMixture>
using HighlandScattering = detail::ScatteringImpl<detail::Highland>
using SingleParameterSmearFunction = std::function<Acts::Result<std::pair<double, double>>(double, generator_t&)>

Smearing function definition for single track parameters.

The function takes the unsmeared parameter and returns the smeared value and a standard deviation.

Template Parameters

generator_t – The type of the random generator.

using StandardChargedElectroMagneticInteractions = InteractionList<detail::StandardScattering, detail::StandardBetheBloch, detail::StandardBetheHeitler>

Standard set of electro-magnetic interactions for charged particles.

Scattering must come first so it is computed with the unmodified initial energy before energy loss is applied.

Warning

The list has no cuts on input particle charge or kinematics, i.e. it relies on the simulator to preselect relevant charged particles before application.

Enums

enum DigitizationError

Values:

enumerator SmearingOutOfRange
enumerator SmearingError
enumerator UndefinedSurface
enumerator MaskingError
enum ProcessType

Process type identifier.

Encodes the type of process that generated a particle.

Values:

enumerator eUndefined
enumerator eDecay
enumerator ePhotonConversion
enumerator eBremsstrahlung
enumerator eNuclearInteraction

Functions

float findCharge(Acts::PdgParticle pdg)

Find the charge for a given PDG particle number.

Returns

Charge in native units or NaN if not available.

float findMass(Acts::PdgParticle pdg)

Find the mass for a given PDG particle number.

Returns

Mass in native units or zero if not available.

std::string_view findName(Acts::PdgParticle pdg)

Find a descriptive particle name for a given PDG particle number.

Returns

Particle name or empty if not available.

std::error_code make_error_code(DigitizationError e)
StandardChargedElectroMagneticInteractions makeStandardChargedElectroMagneticInteractions(double minimumAbsMomentum)

Construct the standard electro-magnetic interactions for charged particles.

Parameters

minimumAbsMomentum – lower p cut on output particles

template<typename signal_t, size_t kSize>
const std::vector<Channel<signal_t, kSize>> mergeChannels(const std::vector<Channel<signal_t, kSize>> &channels)

Generic implementation of a channel merger, currently only additive channel merging.

Template Parameters
  • signal_t – The type of signal, needs operator+= to be defined

  • kSize – the dimensonality of the object (cluster)

Parameters

channels – The channels from one cluster

Returns

A cluster containing the parameter set and cluster size

std::ostream &operator<<(std::ostream &os, ProcessType processType)
std::ostream &operator<<(std::ostream &os, const Particle &particle)
template<Acts::PdgParticle Pdg>
struct AbsPdgExcluder
#include <ActsFatras/Selectors/ParticleSelectors.hpp>

Select all particles except for (anti-)particles of one specific type.

Public Functions

inline bool operator()(const Particle &particle) const
template<Acts::PdgParticle Pdg>
struct AbsPdgSelector
#include <ActsFatras/Selectors/ParticleSelectors.hpp>

Select particles and antiparticles of one specific type.

Public Functions

inline bool operator()(const Particle &particle) const
class Barcode : public Acts::MultiIndex<uint64_t, 12, 12, 16, 8, 16>
#include <ActsFatras/EventData/Barcode.hpp>

Particle identifier that encodes additional event information.

The barcode has to fulfill two separate requirements: be able to act as unique identifier for particles within an event and to encode details on the event structure for fast lookup. Since we only care about tracking here, we need to support two scenarios:

  • Identify which primary/secondary vertex particles belong to. No information on intermediate/unstable/invisible particles needs to be retained. This information is already available in the underlying generator event and should not be duplicated.

  • If (visible) particles convert, decay, or interact with the detector, we need to be able to identify the initial (primary) particle. Typical examples are pion nuclear interactions or electron/gamma conversions.

The vertex information is encoded as two numbers that define the primary and secondary vertex. The primary vertex must be non-zero. Particles with a zero secondary vertex originate directly from the primary vertex.

Within one vertex (primary+secondary) each particle is identified by a particle, generation, and sub-particle number. Particles originating from the vertex must have zero generation and zero sub-particle number; a consequence is that only non-zero generation can have non-zero sub-particle numbers. A non-zero generation indicates that the particle is a descendant of the original particle, e.g. from interactions or decay, while the sub-particle number identifies the descendant particle.

With this encoding, non-primary particles and their primary parent can be easily identified at the expense of not storing the exact decay history.

A barcode with all elements set to zero (the default value) is an invalid value that can be used e.g. to mark missing or unknown particles.

Example

A particle generated in a primary interaction might have the barcode

2|0|14|0|0 -> vertex=2 (primary), particle=14, generation=0, sub=0

A simulation module might generate an interaction and create two new particles. These are descendants of the initial particle and the simulation module can generate the new barcodes directly by increasing the generation number and chosing sub-particle identifiers:

2|0|14|1|0 -> vertex=2 (primary), particle=14, generation=1, sub=0
2|0|14|1|1 -> vertex=2 (primary), particle=14, generation=1, sub=1

If these secondary particles generate further tertiary particles the barcode would be e.g.

2|0|14|2|0 -> vertex=2 (primary), particle=14, generation=2, sub=0

Possible issues

The hierachical nature of the barcode allows barcode creation without a central service. Since the full history is not stored, generated barcodes for higher-generation particles can overlap when generated by independent interactions. Assuming an initial primary particle with barcode

3|4|5|0|0 -> particle=5

a first interaction might create a secondary particle by increasing the generation number (without destroying the initial particle)

3|4|5|1|0 -> particle=5, generation+=1, first sub-particle

The initial particle gets simulated further and at another step a second interaction also creates a new particle. Since it knows nothing about the previously created particle (no central service), it will generate

3|4|5|1|0 -> particle=5, generation+=1, first sub-particle

which is identical to the previously create barcode. These cases can be easily solved by renumbering the sub-particle identifier within each generation to contain unique values. However, this can only be done when all particles are known.

Public Types

using Value = T

The type of ther underlying storage value.

Public Functions

inline Barcode()
Barcode(const Barcode&) = default
Barcode(Barcode&&) = default
inline constexpr Value generation() const

Return the generation identifier.

inline Barcode makeDescendant(Value sub = 0u) const

Construct a new barcode representing a descendant particle.

Parameters

sub – sub-particle index of the new barcode.

Barcode &operator=(const Barcode&) = default
Barcode &operator=(Barcode&&) = default
inline constexpr Value particle() const

Return the particle identifier.

inline constexpr Barcode &setGeneration(Value id)

Set the particle identifier.

inline constexpr Barcode &setParticle(Value id)

Set the parent particle identifier.

inline constexpr Barcode &setSubParticle(Value id)

Set the process identifier.

inline constexpr Barcode &setVertexPrimary(Value id)

Set the primary vertex identifier.

inline constexpr Barcode &setVertexSecondary(Value id)

Set the secondary vertex identifier.

inline constexpr Value subParticle() const

Return the sub-particle identifier.

inline constexpr Value vertexPrimary() const

Return the primary vertex identifier.

inline constexpr Value vertexSecondary() const

Return the secondary vertex identifier.

struct BetheBloch
#include <ActsFatras/Physics/ElectroMagnetic/BetheBloch.hpp>

Simulate energy loss using the Bethe-Bloch/Landau description.

Energy loss is computed using the most probable value and appropriate fluctuations from a Landau distribution. No secondaries are generated for the removed energy.

Public Functions

template<typename generator_t>
inline std::array<Particle, 0> operator()(generator_t &generator, const Acts::MaterialSlab &slab, Particle &particle) const

Simulate energy loss and update the particle parameters.

Parameters
  • generator[in] is the random number generator

  • slab[in] defines the passed material

  • particle[inout] is the particle being updated

Template Parameters

generator_t – is a RandomNumberEngine

Returns

Empty secondaries containers.

Public Members

double scaleFactorMPV = 1.

Scaling for most probable value.

double scaleFactorSigma = 1.

Scaling for Sigma.

struct BetheHeitler
#include <ActsFatras/Physics/ElectroMagnetic/BetheHeitler.hpp>

Simulate electron energy loss using the Bethe-Heitler description.

Bethe-Heitler for electron bremsstrahlung description as described here: “A Gaussian-mixture approximation of the Bethe–Heitler model of electron

energy loss by bremsstrahlung” R. Frühwirth

Public Types

using Scalar = Particle::Scalar
using Vector3 = Particle::Vector3

Public Functions

Particle bremPhoton(const Particle &particle, Scalar gammaE, Scalar rndPsi, Scalar rndTheta1, Scalar rndTheta2, Scalar rndTheta3) const

Simulate the photon emission.

Parameters
  • particle[in] The unmodified electron

  • gammaE[in] Energy of the photon

  • rndPsi[in] Random number for the azimuthal angle

  • rndTheta1[in] Random number for the polar angle

  • rndTheta2[in] Random number for the polar angle

  • rndTheta3[in] Random number for the polar angle

template<typename generator_t>
inline std::array<Particle, 1> operator()(generator_t &generator, const Acts::MaterialSlab &slab, Particle &particle) const

Simulate energy loss and update the particle parameters.

Parameters
  • generator[in] is the random number generator

  • slab[in] defines the passed material

  • particle[inout] is the particle being updated

Template Parameters

generator_t – is a RandomNumberEngine

Returns

Produced photon.

Public Members

double scaleFactor = 1.

A scaling factor to.

bool uniformHertzDipoleAngle = false
template<typename generator_t, size_t kSize>
struct BoundParametersSmearer
#include <ActsFatras/Digitization/UncorrelatedHitSmearer.hpp>

Uncorrelated smearing algorithm for fast digitisation of bound parameters.

The smearer takes a single simulated Hit and generates a smeared parameter vector and associated covariance matrix.

Template Parameters
  • generator_t – Random number generator type

  • kSize – Number of smeared parameters

Public Types

using CovarianceMatrix = Acts::ActsSymMatrix<kSize>
using ParametersVector = Acts::ActsVector<kSize>
using Result = Acts::Result<std::pair<ParametersVector, CovarianceMatrix>>
using Scalar = Acts::ActsScalar

Public Functions

inline Result operator()(generator_t &rng, const Hit &hit, const Acts::Surface &surface, const Acts::GeometryContext &geoCtx) const

Generate smeared measured for configured parameters.

Parameters
  • rng – Random number generator

  • hit – Simulated hit

  • surface – Local surface on which the hit is smeared

  • geoCtx – Geometry context

Return values
  • Smeared – parameters vector and associated covariance on success

  • Error – code for failure

Public Members

std::array<Acts::BoundIndices, kSize> indices = {}

Parameter indices that will be used to create the smeared measurements.

std::array<SingleParameterSmearFunction<generator_t>, kSize> smearFunctions = {}

Public Static Functions

static inline constexpr size_t size()
template<typename signal_t, size_t kSize>
struct Channel
#include <ActsFatras/Digitization/DigitizationData.hpp>

A channel definition: Cell identification, readout word, links.

Template Parameters
  • signal_t – Type of the signal, requires += operator

  • kSize – Number of channel coordinates

Public Functions

inline Channel(std::array<Cell, kSize> cellId_, signal_t value_, std::unordered_set<unsigned int> links_ = {})

Channel constructor.

Parameters
  • cellId_ – The Cell idenficiation and position

  • value_ – The Cell value

  • links_ – The (optional) links to e.g. truth indices

Channel() = delete

Public Members

std::array<Cell, kSize> cellId

The cell identification in sizeof..(kParameters) dimensions.

std::unordered_set<unsigned int> links = {}

The potential (truth) links.

signal_t value = 0

The signal value, as complex as possible, but need += operator and double() cast for the weight.

struct Channelizer
#include <ActsFatras/Digitization/Channelizer.hpp>

The Channelizer splits a surface segment, i.e.

after projection onto the readout surface into channel segments.

Public Types

using Bin2D = std::array<unsigned int, 2>

Shorthand for a 2D bin.

using BinDelta2D = std::array<int, 2>

shorthand for a 2D bin delta

using Segment2D = std::array<Acts::Vector2, 2>

Shorthand for a 2D segment.

Public Functions

std::vector<ChannelSegment> segments(const Acts::GeometryContext &geoCtx, const Acts::Surface &surface, const Acts::BinUtility &segmentation, const Segment2D &segment) const

Divide the surface segment into channel segments.

Note

Channelizing is done in cartesian coordinates (start/end)

Note

The start and end cartesian vector is supposed to be inside the surface bounds (pre-run through the SurfaceMasker)

Note

The segmentation has to be 2-dimensional, even if the actual readout is 1-dimensional, in latter case one bin in the second coordinate direction is required.

Parameters
  • geoCtx – The geometry context for the localToGlobal, etc.

  • surface – The surface for the channelizing

  • segmentation – The segmentation for the channelizing

  • segment – The surface segment (cartesian coordinates)

Returns

a vector of ChannelSegment objects

struct ChannelSegment
#include <ActsFatras/Digitization/Channelizer.hpp>

Nested struct for representing channel steps.

Public Functions

inline ChannelSegment(Bin2D bin_, Segment2D path2D_, double activation_)

Constructor with arguments.

Parameters
  • bin_ – The bin corresponding to this step

  • path2D_ – The start/end 2D position of the segement

  • activation_ – The segment activation (clean: length) for this bin

Public Members

double activation = 0.

The (clipped) value (uncorrected: path length)

Bin2D bin = {0, 0}

The bin of this segment.

Segment2D path2D

The segment start, end points.

struct ChannelStep
#include <ActsFatras/Digitization/Channelizer.hpp>

Nested struct for stepping from one channel to the next.

Public Functions

inline ChannelStep(BinDelta2D delta_, Acts::Vector2 intersect_, const Acts::Vector2 &start)

Constructor with arguments for a ChannelStep.

Parameters
  • delta_ – The bin delta for this step

  • intersect_ – The intersect with the channel boundary

  • start – The start of the surface segment, for path from origin

inline bool operator<(const ChannelStep &cstep) const

Smaller operator for sorting the ChannelStep objects.

The ChannelStep objects can be compared with its path distance from the start (surface segment origin)

Parameters

cstep – The other ChannelStep to be compared

Public Members

BinDelta2D delta = {0, 0}

This is the delta to the last step in bins.

Acts::Vector2 intersect

The intersection with the channel boundary.

double path = 0.

The patlength from the start.

struct ChargedSelector
#include <ActsFatras/Selectors/ParticleSelectors.hpp>

Select all charged particles.

Public Functions

inline bool operator()(const Particle &particle) const
template<typename signal_t, size_t kSize>
struct Cluster
#include <ActsFatras/Digitization/DigitizationData.hpp>

A (simulated) cluster with its constituents.

Template Parameters
  • signal_t – Type of the signal carried, see above

  • kSize – Number of cluster coordinates

Public Types

using CovarianceMatrix = Acts::ActsSymMatrix<kSize>
using ParametersVector = Acts::ActsVector<kSize>
using Scalar = Acts::ActsScalar

Public Functions

template<typename parameters_t, typename covariance_t>
inline Cluster(const Eigen::MatrixBase<parameters_t> &p, const Eigen::MatrixBase<covariance_t> &c, std::array<unsigned int, kSize> cSize, std::vector<Channel<signal_t, kSize>> cChannels)

Cluster constructor.

Parameters
  • p – Measured parameters

  • c – Measurement covariance

  • cSize – The cluster size definition

  • cChannels – The channel

Cluster() = delete

Public Members

std::vector<Channel<signal_t, kSize>> channels

The constituating signal channels.

std::array<unsigned int, kSize> clusterSize

The resulting cluster size along each channel dimension.

CovarianceMatrix covariance = CovarianceMatrix::Zero()

Measurement covariance.

ParametersVector parameters = ParametersVector::Zero()

Measured parameters.

template<typename physics_t, typename input_particle_selector_t, typename output_particle_selector_t, typename child_particle_selector_t = output_particle_selector_t>
struct ContinuousProcess
#include <ActsFatras/Kernel/ContinuousProcess.hpp>

A continuous simulation process based on a physics model plus selectors.

The physics model type must provide a call operator with the following signature

<Particle Container>
operator()(
    generator_t& generator,
    const Acts::MaterialSlab& slab,
    Particle& particle) const

The return type can be any Container with Particle elements.

The input selector defines whether the process is applied while the output selector defines a break condition, i.e. whether to continue simulating the particle propagation. The child selector is used to filter the generated child particles.

Note

The output and child particle selectors are identical unless the child particle selector is explicitely specified.

Template Parameters
  • physics_t – is the physics model type

  • input_particle_selector_t – is the input particle selector

  • output_particle_selector_t – is the output particle selector

  • child_particle_selector_t – is the child particle selector

Public Functions

template<typename generator_t>
inline bool operator()(generator_t &generator, const Acts::MaterialSlab &slab, Particle &particle, std::vector<Particle> &generated) const

Execute the physics process considering the configured selectors.

Parameters
  • generator[in] is the random number generator

  • slab[in] is the passed material

  • particle[inout] is the particle being updated

  • generated[out] is the container of generated particles

Template Parameters

generator_t – must be a RandomNumberEngine

Returns

Break condition, i.e. whether this process stops the propagation

Public Members

physics_t physics

The physics interactions implementation.

child_particle_selector_t selectChildParticle

Child selection: if a generated child particle should be kept.

input_particle_selector_t selectInputParticle

Input selection: if this process applies to this particle.

output_particle_selector_t selectOutputParticle

Output selection: if the particle is still valid after the interaction.

struct EveryParticle
#include <ActsFatras/Selectors/ParticleSelectors.hpp>

No-op particle selector that selects all particles.

Public Functions

inline bool operator()(const Particle&) const
struct EverySurface
#include <ActsFatras/Selectors/SurfaceSelectors.hpp>

Select every surface.

Public Functions

inline constexpr bool operator()(const Acts::Surface&) const
struct FailedParticle
#include <ActsFatras/Kernel/Simulation.hpp>

A particle that failed to simulate.

Public Members

std::error_code error

The associated error code for this particular failure case.

Particle particle

Initial particle state of the failed particle.

This must store the full particle state to be able to handle secondaries that are not in the input particle list. Otherwise they could not be referenced.

template<typename generator_t, size_t kSize>
struct FreeParametersSmearer
#include <ActsFatras/Digitization/UncorrelatedHitSmearer.hpp>

Uncorrelated smearing algorithm for fast digitisation of free parameters.

The smearer takes a single simulated Hit and generates a smeared parameter vector and associated covariance matrix.

Note

Uncorrelated smearing of the direction using each components individually is not recommended

Template Parameters
  • generator_t – Random number generator type

  • kSize – Number of smeared parameters

Public Types

using CovarianceMatrix = Acts::ActsSymMatrix<kSize>
using ParametersVector = Acts::ActsVector<kSize>
using Result = Acts::Result<std::pair<ParametersVector, CovarianceMatrix>>
using Scalar = Acts::ActsScalar

Public Functions

inline Result operator()(generator_t &rng, const Hit &hit) const

Generate smeared measured for configured parameters.

Parameters
  • rng – Random number generator

  • hit – Simulated hit

Return values
  • Smeared – parameters vector and associated covariance on success

  • Error – code for failure

Returns

Smeared free parameter set wrapped in a Result<…> object

Public Members

std::array<Acts::FreeIndices, kSize> indices = {}

Parameter indices that will be used to create the smeared measurements.

std::array<SingleParameterSmearFunction<generator_t>, kSize> smearFunctions

Public Static Functions

static inline constexpr size_t size()
class Hit
#include <ActsFatras/EventData/Hit.hpp>

A simulation hit on a surface.

This is the undigitized, truth hit, i.e. just a recording of the particle state at the surface intersection. Since Fatras is surface-based, the hit position is always constrained to a surface. Depending on the simulated interactions the momentum state before and after might differ and is thus stored as two separate four-vectors.

Public Types

using Scalar = Acts::ActsScalar
using Vector3 = Acts::ActsVector<3>
using Vector4 = Acts::ActsVector<4>

Public Functions

Hit() = default

Construct default hit with (mostly) invalid information.

inline Hit(Acts::GeometryIdentifier geometryId, Barcode particleId, const Vector4 &pos4, const Vector4 &before4, const Vector4 &after4, int32_t index_ = -1)

Construct from four-position and four-momenta.

All quantities are given in the global coordinate system. It is the users responsibility to ensure that the position correspond to a position on the given surface.

Parameters
  • geometryId – Geometry identifier of the surface

  • particleId – Particle identifier of the particle that created the hit

  • pos4 – Particle space-time four-vector on the surface

  • before4 – Particle four-momentum before the interaction

  • after4 – Particle four-momentum after the interaction

  • index_ – Hit index along the particle trajectory

Hit(const Hit&) = default
Hit(Hit&&) = default
inline Scalar depositedEnergy() const

Energy deposited by the hit.

Return values
  • positive – if the particle lost energy when it passed the surface

  • negative – if magic was involved

inline const Vector4 &fourPosition() const

Space-time position four-vector.

inline constexpr Acts::GeometryIdentifier geometryId() const

Geometry identifier of the hit surface.

inline constexpr int32_t index() const

Hit index along the particle trajectory.

Return values

negative – if the hit index is undefined.

inline const Vector4 &momentum4After() const

Particle four-momentum after the hit.

inline const Vector4 &momentum4Before() const

Particle four-momentum before the hit.

Hit &operator=(const Hit&) = default
Hit &operator=(Hit&&) = default
inline constexpr Barcode particleId() const

Particle identifier of the particle that generated the hit.

inline auto position() const

Three-position, i.e. spatial coordinates without the time.

inline Scalar time() const

Time coordinate.

inline Vector3 unitDirection() const

Average normalized particle direction vector through the surface.

inline Vector3 unitDirectionAfter() const

Normalized particle direction vector the hit.

inline Vector3 unitDirectionBefore() const

Normalized particle direction vector before the hit.

template<typename ...processes_t>
class InteractionList
#include <ActsFatras/Kernel/InteractionList.hpp>

Compile-time set of interaction processes for the simulation.

Two different type of interaction processes are supported: continuous and point-like interactions.

Continuous processes scale with the passed material. They tpyically describe effective results of a large number of small interactions such as multiple scattering or ionisation. Continous process types must provide a call operator with the following signature:

template <typename generator_t>
bool
operator()(
    generator_t& rng,
    const Acts::MaterialSlab& slab,
    Particle& particle,
    std::vector<Particle>& generatedParticles) const

If multiple continuous processes are defined, they are executed serially in the order in which they are given.

For point-like processes, the passed material only affects the probability with which they occur but not the interaction itself, e.g. photon conversion into electron pairs. They are simulated by first drawing a limit on the material paths and then executing the interaction with the shortest limit when the drawn amount of material has been passed. Point-like process types must provide the following two member functions:

// generate X0/L0 limits
template <typename generator_t>
std::pair<Scalar, Scalar>
generatePathLimits(
    generator& rng,
    const Particle& particle) const

// run the process simulation
template <typename generator_t>
bool
run(
    generator_t& rng,
    Particle& particle,
    std::vector<Particle>& generatedParticles) const

For both continuous and point-like interactions, the output particle is modified in-place (if needed) and the return value indicates a break condition in the simulation, i.e. the particle is dead (true) or alive (false) after the interaction.

The physics processes are extendable by the user to accomodate their specific requirements. While the set of available physics processes must be configured at compile-time, within that set, processes can again be selectively disabled at run-time. By default all processes are applied.

Note

If an interaction destroys the incoming particle, the process simulation should indicate this via the break condition only and not by reducing the particle momentum to zero. The incoming particle should retain its initial kinematic state; the final kinematic state before destruction is typically of more interest to the user and this simplifies validation.

Public Functions

template<typename generator_t>
inline Selection armPointLike(generator_t &rng, const Particle &particle) const

Arm the point-like interactions by generating limits and select processes.

Template Parameters

generator_t – must be a RandomNumberEngine

Parameters
  • rng[in] is the random number generator

  • particle[in] is the initial particle state

Returns

X0/L0 limits for the particle and the process index that should be executed once the limit has been reached.

inline void disable(size_t process)

Disable a specific process identified by index.

template<typename process_t>
inline void disable()

Disable a specific process identified by type.

Note

Disables only the first element, if multiple elements of the same type exist.

template<size_t kProcess>
inline std::tuple_element_t<kProcess, Processes> &get()

Access a specific process identified by index.

template<typename process_t>
inline process_t &get()

Access a specific process identified by type.

Warning

This function only works if all configured processes have different types.

template<typename generator_t>
inline bool runContinuous(generator_t &rng, const Acts::MaterialSlab &slab, Particle &particle, std::vector<Particle> &generated) const

Simulate the combined effects from all continuous interactions.

Template Parameters

generator_t – must be a RandomNumberEngine

Parameters
  • rng[in] is the random number generator

  • slab[in] is the passed material

  • particle[inout] is the particle being updated

  • generated[out] is the container of generated particles

Returns

Break condition, i.e. whether a process stoped the propagation

template<typename generator_t>
inline bool runPointLike(generator_t &rng, size_t processIndex, Particle &particle, std::vector<Particle> &generated) const

Simulate the effects from a single point-like interaction.

The process index is expected to originate from a previous armPointLike(...) call, but this is not enforced. How to select the correct process requires more information that is not available here.

Template Parameters

generator_t – must be a RandomNumberEngine

Parameters
  • rng[in] is the random number generator

  • processIndex[in] is the index of the process to be executed

  • particle[inout] is the particle being updated

  • generated[out] is the container of generated particles

Returns

Break condition, i.e. whether a process killed the particle

struct Selection
#include <ActsFatras/Kernel/InteractionList.hpp>

Point-like interaction selection.

Public Members

Particle::Scalar l0Limit = std::numeric_limits<Particle::Scalar>::infinity()
size_t l0Process = SIZE_MAX
Particle::Scalar x0Limit = std::numeric_limits<Particle::Scalar>::infinity()
size_t x0Process = SIZE_MAX
class LandauDistribution
#include <ActsFatras/Utilities/LandauDistribution.hpp>

Draw random numbers from a Landau distribution.

Implements the same interface as the standard library distributions.

Public Types

using result_type = double

The type of the generated values.

Public Functions

inline LandauDistribution(double location, double scale)

Construct directly from the distribution parameters.

inline LandauDistribution(const param_type &cfg)

Construct from a parameter object.

LandauDistribution() = default
LandauDistribution(const LandauDistribution&) = default
LandauDistribution(LandauDistribution&&) = default
inline result_type max() const

The maximum value the distribution generates.

inline result_type min() const

The minimum value the distribution generates.

template<typename Generator>
inline result_type operator()(Generator &generator)

Generate a random number from the configured Landau distribution.

template<typename Generator>
inline result_type operator()(Generator &generator, const param_type &params)

Generate a random number from the given Landau distribution.

LandauDistribution &operator=(const LandauDistribution&) = default
LandauDistribution &operator=(LandauDistribution&&) = default
inline param_type param() const

Return the currently configured distribution parameters.

inline void param(const param_type &cfg)

Set the distribution parameters.

inline void reset()

Reset any possible internal state. Noop, since there is no internal state.

Friends

inline friend bool operator!=(const LandauDistribution &lhs, const LandauDistribution &rhs)
inline friend bool operator==(const LandauDistribution &lhs, const LandauDistribution &rhs)

Provide standard comparison operators.

struct param_type
#include <ActsFatras/Utilities/LandauDistribution.hpp>

Parameter struct that contains all distribution parameters.

Public Types

using distribution_type = LandauDistribution

Parameters must link back to the host distribution.

Public Functions

inline param_type(double location_, double scale_)

Construct from parameters.

param_type() = default
param_type(const param_type&) = default
param_type(param_type&&) = default
param_type &operator=(const param_type&) = default
param_type &operator=(param_type&&) = default

Public Members

double location = 0.0

Location parameter.

Warning

This is neither the mean nor the most probable value.

double scale = 1.0

Scale parameter.

Friends

inline friend bool operator!=(const param_type &lhs, const param_type &rhs)
inline friend bool operator==(const param_type &lhs, const param_type &rhs)

Parameters should be EqualityComparable.

template<typename cast_t>
struct Max
#include <ActsFatras/Selectors/SelectorHelpers.hpp>

Select all objects with an extracted value below the cut.

Public Functions

template<typename T>
inline bool operator()(const T &thing) const

Public Members

double valMax = std::numeric_limits<double>::max()
template<typename cast_t>
struct Min
#include <ActsFatras/Selectors/SelectorHelpers.hpp>

Select all objects with an extracted value equal or larger than the cut.

Public Functions

template<typename T>
inline bool operator()(const T &thing) const

Public Members

double valMin = 0.
struct NegativeSelector
#include <ActsFatras/Selectors/ParticleSelectors.hpp>

Select negatively charged particles.

Public Functions

inline bool operator()(const Particle &particle) const
struct NeutralSelector
#include <ActsFatras/Selectors/ParticleSelectors.hpp>

Select neutral particles.

Public Functions

inline bool operator()(const Particle &particle) const
struct NoDecay
#include <ActsFatras/Physics/Decay/NoDecay.hpp>

Decay module that treats all particles as stable.

Public Functions

template<typename generator_t>
inline constexpr Particle::Scalar generateProperTimeLimit(generator_t&, const Particle&) const

Generate the proper time limit for a particle.

Returns

Always returns infinity as limit.

template<typename generator_t>
inline constexpr std::array<Particle, 0> run(generator_t&, const Particle&) const

Decay the particle without generating any descendant particles.

struct NoSurface
#include <ActsFatras/Selectors/SurfaceSelectors.hpp>

Do not select any surface, ever.

Public Functions

inline constexpr bool operator()(const Acts::Surface&) const
struct NuclearInteraction
#include <ActsFatras/Physics/NuclearInteraction/NuclearInteraction.hpp>

This class provides a parametrised nuclear interaction.

The thereby required parametrisation needs to be set and is not provided by default.

Note

This class differs between two different processes labelled as nuclear interaction. Either the initial particle survives (soft) or it gets destroyed (hard) by this process.

Public Types

using Scalar = Particle::Scalar

Public Functions

template<typename generator_t>
inline std::pair<Scalar, Scalar> generatePathLimits(generator_t &generator, const Particle &particle) const

This method evaluates the nuclear interaction length L0.

Template Parameters

generator_t – The random number generator type

Parameters
  • generator[inout] The random number generator

  • particle[in] The ingoing particle

Returns

valid X0 limit and no limit on L0

template<typename generator_t>
inline bool run(generator_t &generator, Particle &particle, std::vector<Particle> &generated) const

This method performs a nuclear interaction.

Template Parameters

generator_t – The random number generator type

Parameters
  • generator[inout] The random number generator

  • particle[inout] The ingoing particle

  • generated[out] Additional generated particles

Returns

True if the particle was killed, false otherwise

Public Members

detail::MultiParticleNuclearInteractionParametrisation multiParticleParameterisation

The storage of the parameterisation.

unsigned int nMatchingTrials = 100

The number of trials to match momenta and inveriant masses.

unsigned int nMatchingTrialsTotal = 1000
class Particle
#include <ActsFatras/EventData/Particle.hpp>

Particle identity information and kinematic state.

Also stores some simulation-specific properties.

Public Types

using Scalar = Acts::ActsScalar
using Vector3 = Acts::ActsVector<3>
using Vector4 = Acts::ActsVector<4>

Public Functions

Particle() = default

Construct a default particle with invalid identity.

inline Particle(Barcode particleId, Acts::PdgParticle pdg, Scalar charge, Scalar mass)

Construct a particle at rest with explicit mass and charge.

Warning

It is the users responsibility that charge and mass match the PDG particle number.

Parameters
  • particleId – Particle identifier within an event

  • pdg – PDG id

  • charge – Particle charge in native units

  • mass – Particle mass in native units

Particle(Barcode particleId, Acts::PdgParticle pdg)

Construct a particle at rest from a PDG particle number.

Charge and mass are retrieved from the particle data table.

Parameters
  • particleId – Particle identifier within an event

  • pdg – PDG particle number

Particle(const Particle&) = default
Particle(Particle&&) = default
inline constexpr Scalar absoluteMomentum() const

Absolute momentum.

inline constexpr Scalar charge() const

Particle charge.

inline Particle &correctEnergy(Scalar delta)

Change the energy by the given amount.

Energy loss corresponds to a negative change. If the updated energy would result in an unphysical value, the particle is put to rest, i.e. its absolute momentum is set to zero.

inline Scalar energy() const

Total energy, i.e. norm of the four-momentum.

inline Vector4 fourMomentum() const

Energy-momentum four-vector.

inline constexpr const Vector4 &fourPosition() const

Space-time position four-vector.

inline constexpr Scalar mass() const

Particle mass.

inline constexpr operator bool() const

Check if the particle is alive, i.e. is not at rest.

inline constexpr bool operator!() const

Check if the particle is dead, i.e is at rest.

Particle &operator=(const Particle&) = default
Particle &operator=(Particle&&) = default
inline constexpr Barcode particleId() const

Particle identifier within an event.

inline constexpr Scalar pathInL0() const

Accumulated path within material measured in interaction lengths.

inline constexpr Scalar pathInX0() const

Accumulated path within material measured in radiation lengths.

inline constexpr Acts::PdgParticle pdg() const

PDG particle number that identifies the type.

inline auto position() const

Three-position, i.e. spatial coordinates without the time.

inline constexpr ProcessType process() const

Which type of process generated this particle.

inline constexpr Scalar properTime() const

Proper time in the particle rest frame.

inline Particle &setAbsoluteMomentum(Scalar absMomentum)

Set the absolute momentum.

inline Particle &setCharge(Scalar charge)

Set the particle charge.

inline Particle &setDirection(const Vector3 &direction)

Set the direction three-vector.

inline Particle &setDirection(Scalar dx, Scalar dy, Scalar dz)

Set the direction three-vector from scalar components.

inline constexpr Particle &setMaterialPassed(Scalar pathInX0, Scalar pathInL0)

Set the accumulated material measured in radiation/interaction lengths.

Parameters
  • pathInX0 – accumulated material measured in radiation lengths

  • pathInL0 – accumulated material measured in interaction lengths

inline Particle &setParticleId(Barcode barcode)

Set the particle ID.

inline Particle &setPosition4(const Vector4 &pos4)

Set the space-time position four-vector.

inline Particle &setPosition4(const Vector3 &position, Scalar time)

Set the space-time position four-vector from three-position and time.

inline Particle &setPosition4(Scalar x, Scalar y, Scalar z, Scalar time)

Set the space-time position four-vector from scalar components.

inline Particle &setProcess(ProcessType proc)

Set the process type that generated this particle.

inline constexpr Particle &setProperTime(Scalar properTime)

Set the proper time in the particle rest frame.

Parameters

properTime – passed proper time in the rest frame

inline Scalar time() const

Time coordinate.

inline Scalar transverseMomentum() const

Absolute momentum in the x-y plane.

inline const Vector3 &unitDirection() const

Unit three-direction, i.e. the normalized momentum three-vector.

inline Particle withParticleId(Barcode particleId) const

Construct a new particle with a new identifier but same kinematics.

Note

This is intentionally not a regular setter. The particle id is used to identify the whole particle. Setting it on an existing particle is usually a mistake.

template<Acts::PdgParticle Pdg>
struct PdgExcluder
#include <ActsFatras/Selectors/ParticleSelectors.hpp>

Select all particles except one specific type.

Particle and Antiparticle are treated as two separate types.

Public Functions

inline bool operator()(const Particle &particle) const
template<Acts::PdgParticle Pdg>
struct PdgSelector
#include <ActsFatras/Selectors/ParticleSelectors.hpp>

Select particles of one specific type.

Particle and Antiparticle are treated as two separate types.

Public Functions

inline bool operator()(const Particle &particle) const
class PhotonConversion
#include <ActsFatras/Physics/ElectroMagnetic/PhotonConversion.hpp>

This class handles the photon conversion.

It evaluates the distance after which the interaction will occur and the final state due the interaction itself.

Public Types

using Scalar = ActsFatras::Particle::Scalar

Public Functions

template<typename generator_t>
Particle::Scalar generateFirstChildEnergyFraction(generator_t &generator, Scalar gammaMom) const
template<typename generator_t>
std::pair<Scalar, Scalar> generatePathLimits(generator_t &generator, const Particle &particle) const

Method for evaluating the distance after which the photon conversion will occur.

Template Parameters

generator_t – Type of the random number generator

Parameters
  • generator[inout] The random number generator

  • particle[in] The particle

Returns

valid X0 limit and no limit on L0

template<typename generator_t>
std::pair<Particle::Scalar, Particle::Scalar> generatePathLimits(generator_t &generator, const Particle &particle) const
template<typename generator_t>
bool run(generator_t &generator, Particle &particle, std::vector<Particle> &generated) const

This method evaluates the final state due to the photon conversion.

Template Parameters

generator_t – Type of the random number generator

Parameters
  • generator[inout] The random number generator

  • particle[inout] The interacting photon

  • generated[out] List of generated particles

Returns

True if the conversion occured, else false

Public Members

Scalar childEnergyScaleFactor = 2.

Scaling factor of children energy.

Scalar conversionProbScaleFactor = 0.98

Scaling factor for photon conversion probability.

struct PlanarSurfaceDrift
#include <ActsFatras/Digitization/PlanarSurfaceDrift.hpp>

The PlanarSurfaceDrift takes an intersection in the nominal surface and projects the ends into the readout surface, which can be at : -1, 0, 1.

A Lorentz drift angle can be applied.

Public Types

using Segment2D = std::array<Acts::Vector2, 2>

Shorthand for a 2D segment.

Public Functions

Segment2D toReadout(const Acts::GeometryContext &gctx, const Acts::Surface &surface, double thickness, const Acts::Vector3 &pos, const Acts::Vector3 &dir, const Acts::Vector3 &driftdir) const

Drift the full 3D segment onto a surface 2D readout plane.

Note

a drift direction of (0,0,0) is drift to central plane any other a drift direction with driftDir.z() != 0. will result on a readout on either + 0.5*depletion or -0.5*depletion

Note

without masking

Parameters
  • gctx – The current Geometry context

  • surface – The nominal intersection surface

  • thickness – The emulated module/depletion thickness

  • pos – The position in global coordinates

  • dir – The direction in global coordinates

  • driftdir – The drift direction in local (surface) coordinates

Returns

a Segment on the readout surface

struct PlanarSurfaceMask
#include <ActsFatras/Digitization/PlanarSurfaceMask.hpp>

A brief struct that allows to apply a surface bound mask.

Public Types

using Segment2D = std::array<Acts::Vector2, 2>

Shorthand for a 2-d segment;.

Public Functions

Acts::Result<Segment2D> annulusMask(const Acts::AnnulusBounds &aBounds, const Segment2D &segment, bool firstInside) const

Apply the mask of an annulus disk.

Parameters
  • aBounds – The annulus disc for the masking

  • segment – The track segment (on surface)

  • firstInside – The indicator if the first is inside

Returns

a result wrapping a segment

Acts::Result<Segment2D> apply(const Acts::Surface &surface, const Segment2D &segment) const

Apply the mask on the segment.

  • If the semgent is full inside the surface, return unchanged

  • Otherwise mask/clip the segment to fit into the bounds

Note

Only PlaneSurface/DiscSurface are supported

Note

If both end points of the segment are inside, the segment is not clipped/masked, even if it would cross a surface boundary. Examples for those would be non-covex polygons or segments on a radial bound, where the radial boundary is crossed. Such segments do not occur in Digitization, as the hit has to be inside the surface bounds to start with.

Parameters
  • surface – The surface in question

  • segment – The track segment (on surface)

Returns

a result wrapping a segment

Acts::Result<Segment2D> polygonMask(const std::vector<Acts::Vector2> &vertices, const Segment2D &segment, bool firstInside) const

Apply the mask of a polygon.

Parameters
  • vertices – The vertices of the polygon

  • segment – The track segment (on surface)

  • firstInside – The indicator if the first is inside

Returns

a result wrapping a segment

Acts::Result<Segment2D> radialMask(const Acts::RadialBounds &rBounds, const Segment2D &segment, const Segment2D &polarSegment, bool firstInside) const

Apply the mask of a Radial disk.

Parameters
  • rBounds – The radial disc for the masking

  • segment – The track segment (on surface)

  • polarSegment – The track segmetn (on surface, in polar)

  • firstInside – The indicator if the first is inside

Returns

a result wrapping a segment

Public Members

Acts::detail::IntersectionHelper2D intersector = {}
struct PositiveSelector
#include <ActsFatras/Selectors/ParticleSelectors.hpp>

Select positively charged particles.

Public Functions

inline bool operator()(const Particle &particle) const
template<typename cast_t>
struct Range
#include <ActsFatras/Selectors/SelectorHelpers.hpp>

Select all objects with an extracted value within the range.

The range is defined as the left, half-open interval within the cuts.

Public Functions

template<typename T>
inline bool operator()(const T &thing) const

Public Members

double valMax = std::numeric_limits<double>::max()
double valMin = std::numeric_limits<double>::lowest()
template<typename charged_selector_t, typename charged_simulator_t, typename neutral_selector_t, typename neutral_simulator_t>
struct Simulation
#include <ActsFatras/Kernel/Simulation.hpp>

Multi-particle/event simulation.

The selector types for charged and neutral particles do not need to check for the particle charge. This is done automatically by the simulator to ensure consistency.

Template Parameters
  • charged_selector_t – Callable selector type for charged particles

  • charged_simulator_t – Single particle simulator for charged particles

  • neutral_selector_t – Callable selector type for neutral particles

  • neutral_simulator_t – Single particle simulator for neutral particles

Public Functions

inline Simulation(charged_simulator_t &&charged_, neutral_simulator_t &&neutral_)

Construct from the single charged/neutral particle simulators.

template<typename generator_t, typename input_particles_t, typename output_particles_t, typename hits_t>
inline Acts::Result<std::vector<FailedParticle>> simulate(const Acts::GeometryContext &geoCtx, const Acts::MagneticFieldContext &magCtx, generator_t &generator, const input_particles_t &inputParticles, output_particles_t &simulatedParticlesInitial, output_particles_t &simulatedParticlesFinal, hits_t &hits) const

Simulate multiple particles and generated secondaries.

This takes all input particles and simulates those passing the selection using the appropriate simulator. All selected particle states including additional ones generated from interactions are stored in separate output containers; both the initial state at the production vertex and the final state after propagation are stored. Hits generated from selected input and generated particles are stored in the hit container.

Note

Parameter edge-cases can lead to errors in the underlying propagator and thus to particles that fail to simulate. Here, full events are simulated and the failure to simulate one particle should not be considered a general failure of the simulator. Instead, a list of particles that fail to simulate is provided to the user. It is the users responsibility to handle them.

Note

Failed particles are removed from the regular output, i.e. they do not appear in the simulated particles containers nor do they generate hits.

Warning

Particle-hit association is based on particle ids generated during the simulation. This requires that all input particles must have generation and sub-particle number set to zero.

Parameters
  • geoCtx – is the geometry context to access surface geometries

  • magCtx – is the magnetic field context to access field values

  • generator – is the random number generator

  • inputParticles – contains all particles that should be simulated

  • simulatedParticlesInitial – contains initial particle states

  • simulatedParticlesFinal – contains final particle states

  • hits – contains all generated hits

Return values
  • Acts::Result::Error – if there is a fundamental issue

  • Acts::Result::Success – with all particles that failed to simulate

Template Parameters
  • generator_t – is the type of the random number generator

  • input_particles_t – is a Container for particles

  • output_particles_t – is a SequenceContainer for particles

  • hits_t – is a SequenceContainer for hits

Public Members

charged_simulator_t charged
neutral_simulator_t neutral
charged_selector_t selectCharged
neutral_selector_t selectNeutral
struct SimulationResult
#include <ActsFatras/Kernel/SimulationResult.hpp>

Single particle simulation result (and intermediate state).

This result struct is used by multiple components and is thus defined separately from its usage.

Public Members

std::vector<Particle> generatedParticles

Additional particles generated by interactions or decay.

std::vector<Hit> hits

Hits created by the particle.

bool isAlive = true
Particle::Scalar l0Limit = std::numeric_limits<Particle::Scalar>::quiet_NaN()
size_t l0Process = SIZE_MAX
Particle particle

Current/ final particle state.

Particle::Scalar properTimeLimit = std::numeric_limits<Particle::Scalar>::quiet_NaN()
Particle::Scalar x0Limit = std::numeric_limits<Particle::Scalar>::quiet_NaN()
size_t x0Process = SIZE_MAX
template<typename propagator_t, typename interactions_t, typename hit_surface_selector_t, typename decay_t>
struct SingleParticleSimulation
#include <ActsFatras/Kernel/Simulation.hpp>

Single particle simulation with fixed propagator, interactions, and decay.

Template Parameters
  • generator_t – random number generator

  • interactions_t – interaction list

  • hit_surface_selector_t – selector for hit surfaces

  • decay_t – decay module

Public Functions

inline SingleParticleSimulation(propagator_t &&propagator_, std::unique_ptr<const Acts::Logger> _logger)

Alternatively construct the simulator with an external logger.

template<typename generator_t>
inline Acts::Result<SimulationResult> simulate(const Acts::GeometryContext &geoCtx, const Acts::MagneticFieldContext &magCtx, generator_t &generator, const Particle &particle) const

Simulate a single particle without secondaries.

Template Parameters

generator_t – is the type of the random number generator

Parameters
  • geoCtx – is the geometry context to access surface geometries

  • magCtx – is the magnetic field context to access field values

  • generator – is the random number generator

  • particle – is the initial particle state

Returns

Simulated particle state, hits, and generated particles.

Public Members

decay_t decay

Decay module.

interactions_t interactions

Interaction list containing the simulated interactions.

std::unique_ptr<const Acts::Logger> logger

Logger for debug output.

propagator_t propagator

How and within which geometry to propagate the particle.

hit_surface_selector_t selectHitSurface

Selector for surfaces that should generate hits.

namespace Casts
struct AbsEta
#include <ActsFatras/Selectors/KinematicCasts.hpp>

Retrieve the direction absolute pseudo-rapidity.

Public Functions

inline double operator()(const Particle &particle) const
struct AbsVz
#include <ActsFatras/Selectors/KinematicCasts.hpp>

Retrieve the longitudinal absolute distance of the position to the origin.

Public Functions

inline double operator()(const Particle &particle) const
struct E
#include <ActsFatras/Selectors/KinematicCasts.hpp>

Retrieve the total energy.

Public Functions

inline double operator()(const Particle &particle) const
struct Eta
#include <ActsFatras/Selectors/KinematicCasts.hpp>

Retrieve the direction pseudo-rapidity.

Public Functions

inline double operator()(const Particle &particle) const
struct P
#include <ActsFatras/Selectors/KinematicCasts.hpp>

Retrieve the absolute momentum.

Public Functions

inline double operator()(const Particle &particle) const
struct Pt
#include <ActsFatras/Selectors/KinematicCasts.hpp>

Retrieve the transverse momentum.

Public Functions

inline double operator()(const Particle &particle) const
struct Vrho
#include <ActsFatras/Selectors/KinematicCasts.hpp>

Retrieve the transverse absolute distance of the position to the origin.

Public Functions

inline double operator()(const Particle &particle) const
struct Vz
#include <ActsFatras/Selectors/KinematicCasts.hpp>

Retrieve the longitudinal distance of the position to the origin.

Public Functions

inline double operator()(const Particle &particle) const