API Reference

Classes

class AbstractVolume : public Acts::Volume

AbstractVolume description inside the tracking realm.

This is the purely geometrical object volume.

The Acts::AbstractVolume is constructed by giving a pointer to a Transform3 and a pointer to Acts::VolumeBounds, this implies that the ownership of the objects pointed to is passed as well. For memory optimisation, the AbstractVolume can also be constructed with shared_ptr objects.

A Acts::AbstractVolume is at first a collection class of Acts::BoundarySurface, the vector of Acts::BoundarySurface is returned by the Acts::VolumeBounds that carry a decompose method.

Boundary surfaces can be shared between AbstractVolumes to enhance automatic navigation between AbstractVolumes, therefore they are reference counted by a std::shared_ptr holder class.

Public Functions

AbstractVolume() = delete
AbstractVolume(const AbstractVolume &vol) = default
AbstractVolume(const Transform3 &transform, VolumeBoundsPtr volbounds)

Constructor with shared Transform3*, VolumeBounds*.

Parameters:
  • transform – is the global 3d transformation into the volume frame

  • volbounds – is the boundary definition

~AbstractVolume() override = default
const std::vector<BoundarySurfacePtr> &boundarySurfaces() const

Method to return the BoundarySurfaces.

Returns:

the vector of boundary surfaces

AbstractVolume &operator=(const AbstractVolume &vol) = delete
class AnnulusBounds : public Acts::DiscBounds

Class that implements a (potentially asymmetric) bounds with difference between surface bound center and surface coordinate center.

These bounds combine two different systems:

  • module system : radial bounds centred on the moduleOrigin

  • strip system : phi bounds centred on the stripOrigin

The measurement will be done in the strip system, with r/phi local coordinates.

Public Types

enum BoundValues

Values:

enumerator eMinR
enumerator eMaxR
enumerator eMinPhiRel
enumerator eMaxPhiRel
enumerator eAveragePhi
enumerator eOriginX
enumerator eOriginY
enumerator eSize

Public Functions

AnnulusBounds() = delete
AnnulusBounds(const AnnulusBounds &source) = default
AnnulusBounds(const std::array<double, eSize> &values) noexcept(false)

Constructor - from parameters array.

Parameters:

values – The parameter array

inline AnnulusBounds(double minR, double maxR, double minPhiRel, double maxPhiRel, const Vector2 &moduleOrigin = {0, 0}, double avgPhi = 0) noexcept(false)

Default constructor from parameters.

Note

For morigin you need to actually calculate the cartesian offset

Parameters:
  • minR – inner radius, in module system

  • maxR – outer radius, in module system

  • minPhiRel – right angular edge, in strip system, rel to avgPhi

  • maxPhiRel – left angular edge, in strip system, rel to avgPhi

  • moduleOrigin – The origin offset between the two systems.

  • avgPhi – (Optional) internal rotation of this bounds object’s local frame

inline virtual double binningValuePhi() const final

Return a reference radius for binning.

inline virtual double binningValueR() const final

Return a reference radius for binning.

std::vector<Vector2> corners() const

This method returns the four corners of the bounds in polar coordinates Starting from the upper right (max R, pos locX) and proceeding clock-wise i.e.

(max R; pos locX), (min R; pos locX), (min R; neg loc X), (max R: neg locX)

inline virtual bool coversFullAzimuth() const final

Returns true for full phi coverage.

inline double get(BoundValues bValue) const

Access to the bound values.

Parameters:

bValue – the class nested enum for the array access

virtual bool inside(const Vector2 &lposition, const BoundaryCheck &bcheck) const final

Inside check for the bounds object driven by the boundary check directive Each Bounds has a method inside, which checks if a LocalPosition is inside the bounds Inside can be called without/with tolerances.

Parameters:
  • lposition – Local position (assumed to be in right surface frame)

  • bcheck – boundary check directive

Returns:

boolean indicator for the success of this operation

inline virtual bool insideRadialBounds(double R, double tolerance = 0.) const final

Checks if this is inside the radial coverage given the a tolerance.

Vector2 moduleOrigin() const

Returns moduleOrigin, but rotated out, so averagePhi is already considered.

The module origin needs to consider the rotation introduced by averagePhi

Returns:

The origin of the local frame

inline double phiMax() const

Returns the left angular edge of the module.

Returns:

The left side angle

inline double phiMin() const

Returns the right angular edge of the module.

Returns:

The right side angle

inline virtual double rMax() const final

This method returns outer radius.

inline virtual double rMin() const final

This method returns inner radius.

virtual std::ostream &toStream(std::ostream &sl) const final

Outstream operator.

Parameters:

sl – is the ostream to be dumped into

inline virtual SurfaceBounds::BoundsType type() const final

Return the bounds type - for persistency optimization.

Returns:

is a BoundsType enum

inline virtual std::vector<double> values() const final

Return the bound values as dynamically sized vector.

Returns:

this returns a copy of the internal values

virtual std::vector<Vector2> vertices(unsigned int lseg) const override

This method returns the xy coordinates of the four corners of the bounds in module coordinates (in x/y) Starting from the upper right (max R, pos locX) and proceeding clock-wise i.e.

(max R; pos locX), (min R; pos locX), (min R; neg loc X), (max R: neg locX)

Note

that that if lseg > 0, the extrema points are given, which may slightly alter the number of segments returned

Parameters:

lseg – the number of segments used to approximate and eventually curved line

Returns:

vector for vertices in 2D

template<int NComponents, int PolyDegree>
class AtlasBetheHeitlerApprox

This class approximates the Bethe-Heitler distribution as a gaussian mixture.

To enable an approximation for continuous input variables, the weights, means and variances are internally parametrized as a Nth order polynomial.

Public Types

using Data = std::array<PolyData, NComponents>

Public Functions

inline constexpr AtlasBetheHeitlerApprox(const Data &lowData, const Data &highData, bool lowTransform, bool highTransform, double lowLimit = 0.1, double highLimit = 0.2)

Construct the Bethe-Heitler approximation description with two parameterizations, one for lower ranges, one for higher ranges.

Is it assumed that the lower limit of the high-x/x0 data is equal to the upper limit of the low-x/x0 data.

Parameters:
  • lowData – data for the lower x/x0 range

  • highData – data for the higher x/x0 range

  • lowTransform – whether the low data need to be transformed

  • highTransform – whether the high data need to be transformed

  • lowLimit – the upper limit for the low data

  • highLimit – the upper limit for the high data

inline auto mixture(ActsScalar x) const

Generates the mixture from the polynomials and reweights them, so that the sum of all weights is 1.

Parameters:

x – pathlength in terms of the radiation length

inline constexpr auto numComponents() const

Returns the number of components the returned mixture will have.

inline constexpr bool validXOverX0(ActsScalar x) const

Checks if an input is valid for the parameterization.

Parameters:

x – pathlength in terms of the radiation length

Public Static Functions

static inline auto loadFromFiles(const std::string &low_parameters_path, const std::string &high_parameters_path, double lowLimit = 0.1, double highLimit = 0.2)

Loads a parameterization from a file according to the Atlas file description.

Parameters:
  • low_parameters_path – Path to the foo.par file that stores the parameterization for low x/x0

  • high_parameters_path – Path to the foo.par file that stores the parameterization for high x/x0

  • lowLimit – the upper limit for the low x/x0-data

  • highLimit – the upper limit for the high x/x0-data

struct PolyData

Public Members

std::array<ActsScalar, PolyDegree + 1> meanCoeffs
std::array<ActsScalar, PolyDegree + 1> varCoeffs
std::array<ActsScalar, PolyDegree + 1> weightCoeffs
class AtlasStepper

the AtlasStepper implementation for the

Public Types

using BoundState = std::tuple<BoundTrackParameters, Jacobian, double>
using Covariance = BoundSquareMatrix
using CurvilinearState = std::tuple<CurvilinearTrackParameters, Jacobian, double>
using Jacobian = BoundMatrix

Public Functions

inline AtlasStepper(std::shared_ptr<const MagneticFieldProvider> bField)
inline double absoluteMomentum(const State &state) const

Absolute momentum accessor.

Parameters:

state – [in] The stepping state (thread-local cache)

inline Result<BoundState> boundState(State &state, const Surface &surface, bool transportCov = true, const FreeToBoundCorrection &freeToBoundCorrection = FreeToBoundCorrection(false)) const

Create and return the bound state at the current position.

Parameters:
  • state[in] State that will be presented as BoundState

  • surface[in] The surface to which we bind the state

  • transportCov[in] Flag steering covariance transport

  • freeToBoundCorrection[in] Correction for non-linearity effect during transform from free to bound

Returns:

A bound state:

  • the parameters at the surface

  • the stepwise jacobian towards it

  • and the path length (from start - for ordering)

inline double charge(const State &state) const

Charge access.

Parameters:

state – [in] The stepping state (thread-local cache)

inline CurvilinearState curvilinearState(State &state, bool transportCov = true) const

Create and return a curvilinear state at the current position.

Parameters:
  • state[in] State that will be presented as CurvilinearState

  • transportCov[in] Flag steering covariance transport

Returns:

A curvilinear state:

  • the curvilinear parameters at given position

  • the stepweise jacobian towards it

  • and the path length (from start - for ordering)

inline Vector3 direction(const State &state) const
inline Result<Vector3> getField(State &state, const Vector3 &pos) const

Get the field for the stepping It checks first if the access is still within the Cell, and updates the cell if necessary, then it takes the field from the cell.

Parameters:
  • state[inout] is the stepper state associated with the track the magnetic field cell is used (and potentially updated)

  • pos[in] is the field position

inline double getStepSize(const State &state, ConstrainedStep::Type stype) const

Get the step size.

Parameters:
  • state – [in] The stepping state (thread-local cache)

  • stype – [in] The step size type to be returned

inline State makeState(std::reference_wrapper<const GeometryContext> gctx, std::reference_wrapper<const MagneticFieldContext> mctx, const BoundTrackParameters &par, double ssize = std::numeric_limits<double>::max(), double stolerance = s_onSurfaceTolerance) const
inline Vector3 momentum(const State &state) const
inline std::string outputStepSize(const State &state) const

Output the Step Size - single component.

Parameters:

state[inout] The stepping state (thread-local cache)

inline double overstepLimit(const State &state) const

Overstep limit.

Parameters:

state – The stepping state (thread-local cache)

inline const ParticleHypothesis &particleHypothesis(const State &state) const

Particle hypothesis.

Parameters:

state – [in] The stepping state (thread-local cache)

inline Vector3 position(const State &state) const
inline double qOverP(const State &state) const
inline void releaseStepSize(State &state, ConstrainedStep::Type stype) const

Release the Step size.

Parameters:
  • state[inout] The stepping state (thread-local cache)

  • stype[in] The step size type to be released

inline void resetState(State &state, const BoundVector &boundParams, const BoundSquareMatrix &cov, const Surface &surface, const double stepSize = std::numeric_limits<double>::max()) const

Resets the state.

Parameters:
  • state[inout] State of the stepper

  • boundParams[in] Parameters in bound parametrisation

  • cov[in] Covariance matrix

  • surface[in] Reset state will be on this surface

  • stepSize[in] Step size

inline void setIdentityJacobian(State &state) const

Method that reset the Jacobian to the Identity for when no bound state are available.

Parameters:

state[inout] State of the stepper

template<typename propagator_state_t, typename navigator_t>
inline Result<double> step(propagator_state_t &state, const navigator_t&) const

Perform the actual step on the state.

Parameters:

state – is the provided stepper state (caller keeps thread locality)

inline double time(const State &state) const

Time access.

inline void transportCovarianceToBound(State &state, const Surface &surface, const FreeToBoundCorrection& = FreeToBoundCorrection(false)) const

Method for on-demand transport of the covariance to a new curvilinear frame at current position, or direction of the state.

Parameters:
  • state[inout] State of the stepper

  • surface[in] is the surface to which the covariance is forwarded to

inline void transportCovarianceToCurvilinear(State &state) const

Method for on-demand transport of the covariance to a new curvilinear frame at current position, or direction of the state.

Parameters:

state[inout] State of the stepper

inline void update(State &state, const FreeVector &parameters, const BoundVector &boundParams, const Covariance &covariance, const Surface &surface) const

The state update method.

Parameters:
  • state[inout] The stepper state for

  • parameters[in] The new free track parameters at start

  • boundParams[in] Corresponding bound parameters

  • covariance[in] The updated covariance matrix

  • surface[in] The surface used to update the pVector

inline void update(State &state, const Vector3 &uposition, const Vector3 &udirection, double qop, double time) const

Method to update momentum, direction and p.

Parameters:
  • state – The state object

  • uposition – the updated position

  • udirection – the updated direction

  • qop – the updated momentum value

  • time – the update time

template<typename object_intersection_t>
inline void updateStepSize(State &state, const object_intersection_t &oIntersection, Direction, bool release = true) const

Update step size.

It checks the status to the reference surface & updates the step size accordingly

Parameters:
  • state – [in,out] The stepping state (thread-local cache)

  • oIntersection – [in] The ObjectIntersection to layer, boundary, etc

  • release – [in] boolean to trigger step size release

inline void updateStepSize(State &state, double stepSize, ConstrainedStep::Type stype, bool release = true) const

Update step size - explicitly with a double.

Parameters:
  • state[inout] The stepping state (thread-local cache)

  • stepSize[in] The step size value

  • stype[in] The step size type to be set

  • release – [in] Do we release the step size?

inline Intersection3D::Status updateSurfaceStatus(State &state, const Surface &surface, std::uint8_t index, Direction navDir, const BoundaryCheck &bcheck, ActsScalar surfaceTolerance = s_onSurfaceTolerance, const Logger &logger = getDummyLogger()) const

Update surface status.

This method intersect the provided surface and update the navigation step estimation accordingly (hence it changes the state). It also returns the status of the intersection to trigger onSurface in case the surface is reached.

Parameters:
  • state[inout] The stepping state (thread-local cache)

  • surface[in] The surface provided

  • index[in] The surface intersection index

  • navDir[in] The navigation direction

  • bcheck[in] The boundary check for this status update

  • surfaceTolerance[in] Surface tolerance used for intersection

  • logger[in] Logger instance to use

struct State

Nested State struct for the local caching.

Public Functions

State() = delete

Default constructor - deleted.

template<typename Parameters>
inline State(const GeometryContext &gctx, MagneticFieldProvider::Cache fieldCacheIn, const Parameters &pars, double ssize = std::numeric_limits<double>::max(), double stolerance = s_onSurfaceTolerance)

Constructor.

Template Parameters:

Type – of TrackParameters

Parameters:
  • gctx[in] The geometry context tof this call

  • fieldCacheIn[in] The magnetic field cache for this call

  • pars[in] Input parameters

  • ssize[in] the steps size limitation

  • stolerance[in] is the stepping tolerance

Public Members

Covariance cov = Covariance::Zero()
const Covariance *covariance = nullptr
bool covTransport = false
bool debug = false

Debug output the string where debug messages are stored (optionally)

std::size_t debugMsgWidth = 50
std::size_t debugPfxWidth = 30

buffer & formatting for consistent output

std::string debugString = ""
Vector3 field
MagneticFieldProvider::Cache fieldCache

It caches the current magnetic field cell and stays (and interpolates) within as long as this is valid.

See step() code for details.

std::reference_wrapper<const GeometryContext> geoContext

Cache the geometry context.

double jacobian[eBoundSize * eBoundSize] = {}
double maxPathLength = 0
bool mcondition = false
bool needgradient = false
bool newfield = true
double parameters[eBoundSize] = {0., 0., 0., 0., 0., 0.}

Storage pattern of pVector /dL0 /dL1 /dPhi /the /dCM /dT X ->P[0] dX / P[ 8] P[16] P[24] P[32] P[40] P[48] Y ->P[1] dY / P[ 9] P[17] P[25] P[33] P[41] P[49] Z ->P[2] dZ / P[10] P[18] P[26] P[34] P[42] P[50] T ->P[3] dT/ P[11] P[19] P[27] P[35] P[43] P[51] Ax ->P[4] dAx/ P[12] P[20] P[28] P[36] P[44] P[52] Ay ->P[5] dAy/ P[13] P[21] P[29] P[37] P[45] P[53] Az ->P[6] dAz/ P[14] P[22] P[30] P[38] P[46] P[54] CM ->P[7] dCM/ P[15] P[23] P[31] P[39] P[47] P[55] Cache: P[56] - P[59].

ParticleHypothesis particleHypothesis
double pathAccumulated = 0.
double previousStepSize = 0.
std::array<double, 60> pVector = {}
bool state_ready = false
double step = 0
ConstrainedStep stepSize
double tolerance = s_onSurfaceTolerance

The tolerance for the stepping.

bool useJacobian = false
class BinUtility

The BinUtility class that translated global and local position into a bins of a BinnedArray, most performant is equidistant binning without a transform, however, optionally a transform can be provided, e.g.

for binning on shifted object, the transform is usually shared with the geometric object the Array is defined on, for performance reasons, also the inverse transform is stored.

Public Functions

inline BinUtility()

Constructor for equidistant.

BinUtility(BinUtility &&sbu) = default
inline BinUtility(const BinningData &bData, const Transform3 &tForm = Transform3::Identity())

Constructor from BinningData directly.

Parameters:
  • bData – is the provided binning data

  • tForm – is the (optional) transform

BinUtility(const BinUtility &sbu) = default

Copy constructor.

Parameters:

sbu – is the source bin utility

inline BinUtility(const Transform3 &tForm)

Constructor with only a Transform3.

Parameters:

tForm – is the local to global transform

inline BinUtility(std::size_t bins, float min, float max, BinningOption opt = open, BinningValue value = binX, const Transform3 &tForm = Transform3::Identity())

Constructor for equidistant.

Parameters:
  • bins – is the number of bins

  • min – in the minimal value

  • max – is the maximal value

  • opt – is the binning option : open, closed

  • value – is the binninb value : binX, binY, binZ, etc.

  • tForm – is the (optional) transform

inline BinUtility(std::vector<float> &bValues, BinningOption opt = open, BinningValue value = binPhi, const Transform3 &tForm = Transform3::Identity())

Constructor for arbitrary.

Parameters:
  • bValues – is the boundary values of the binning

  • opt – is the binning option : open, closed

  • value – is the binninb value : binX, binY, binZ, etc.

  • tForm – is the (optional) transform

~BinUtility() = default

Virtual Destructor.

inline std::size_t bin(const Vector2 &lposition, std::size_t ba = 0) const

Bin from a 2D vector (following local parameters defintitions)

  • no optional transform applied

  • USE WITH CARE !!

You need to make sure that the local position is actually in the binning frame of the BinUtility

Parameters:
  • lposition – is the local position to be set

  • ba – is the bin dimension

Returns:

bin calculated from local

inline std::size_t bin(const Vector3 &position, std::size_t ba = 0) const

Bin from a 3D vector (already in binning frame)

Parameters:
  • position – is the 3D position to be evaluated

  • ba – is the bin dimension

Returns:

is the bin value

inline const std::vector<BinningData> &binningData() const

Return the binning data vector.

inline BinningValue binningValue(std::size_t ba = 0) const

The type/value of the binning.

Parameters:

ba – is the binaccessor

Returns:

the binning value of the accessor entry

inline std::size_t bins() const

Return the total number of bins.

inline std::size_t bins(std::size_t ba) const

Number of bins.

Parameters:

ba – is the binaccessor

Returns:

std::size_t is the bins of the accessor entry

inline std::array<std::size_t, 3> binTriple(const Vector3 &position) const

Bin-triple fast access.

  • calculate the bin triple with one transform

Parameters:

position – is the 3D position to be evaluated

Returns:

is the bin value in 3D

inline std::size_t dimensions() const

First bin maximal value.

Returns:

the dimension of the binning data

inline bool inside(const Vector3 &position) const

Check if bin is inside from Vector2 - optional transform applied.

Parameters:

position – is the global position to be evaluated

Returns:

is a boolean check

inline std::size_t max(std::size_t ba = 0) const

First bin maximal value.

Parameters:

ba – is the binaccessor

Returns:

std::size_t is the maximal bin of the accessor entry

inline int nextDirection(const Vector3 &position, const Vector3 &direction, std::size_t ba = 0) const

Return the other direction for fast interlinking.

Parameters:
  • position – is the global position for the next search

  • direction – is the global position for the next search

  • ba – is the bin accessor

Returns:

the next bin

inline BinUtility &operator+=(const BinUtility &gbu)

Operator+= to make multidimensional BinUtility.

Parameters:

gbu – is the additional BinUtility to be chosen

BinUtility &operator=(BinUtility&&) = default
inline BinUtility &operator=(const BinUtility &sbu)

Assignment operator.

Parameters:

sbu – is the source bin utility

inline bool operator==(const BinUtility &other) const

Equality operator.

inline std::size_t serialize(const std::array<std::size_t, 3> &bin) const

Serialize the bin triple.

  • this creates a simple std::size_t from a triple object

Parameters:

bin – is the bin to be serialized

inline std::ostream &toStream(std::ostream &sl, const std::string &indent = "") const

Output Method for std::ostream, to be overloaded by child classes.

Parameters:
  • sl – is the ostream to be dumped into

  • indent – the current indentation

Returns:

the input stream

inline std::string toString(const std::string &indent = "") const

Output into a string.

Parameters:

indent – the current indentation

Returns:

a string with the stream information

inline const Transform3 &transform() const

Transform applied to global positions before lookup.

Returns:

Shared pointer to transform

class BinnedSurfaceMaterial : public Acts::ISurfaceMaterial

It extends the ISurfaceMaterial base class and is an array pf MaterialSlab.

This is not memory optimised as every bin holds one material property object.

Public Functions

BinnedSurfaceMaterial() = delete

Default Constructor - deleted.

BinnedSurfaceMaterial(BinnedSurfaceMaterial &&bsm) = default

Copy Move Constructor.

Parameters:

bsm – is the source object to be copied

BinnedSurfaceMaterial(const BinnedSurfaceMaterial &bsm) = default

Copy Constructor.

Parameters:

bsm – is the source object to be copied

BinnedSurfaceMaterial(const BinUtility &binUtility, MaterialSlabMatrix fullProperties, double splitFactor = 0., MappingType mappingType = MappingType::Default)

Explicit constructor with only full MaterialSlab, for two-dimensional binning.

The split factors:

  • 1. : oppositePre

  • 0. : alongPre ===> 1 Dimensional array

Parameters:
  • binUtility – defines the binning structure on the surface (copied)

  • fullProperties – is the vector of properties as recorded (moved)

  • splitFactor – is the pre/post splitting directive

  • mappingType – is the type of surface mapping associated to the surface

BinnedSurfaceMaterial(const BinUtility &binUtility, MaterialSlabVector fullProperties, double splitFactor = 0., MappingType mappingType = MappingType::Default)

Explicit constructor with only full MaterialSlab, for one-dimensional binning.

The split factors:

  • 1. : oppositePre

  • 0. : alongPre ===> 1 Dimensional array

Parameters:
  • binUtility – defines the binning structure on the surface (copied)

  • fullProperties – is the vector of properties as recorded (moved)

  • splitFactor – is the pre/post splitting directive

  • mappingType – is the type of surface mapping associated to the surface

~BinnedSurfaceMaterial() override = default

Destructor.

inline const BinUtility &binUtility() const

Return the BinUtility.

inline const MaterialSlabMatrix &fullMaterial() const

Retrieve the entire material slab matrix.

virtual const MaterialSlab &materialSlab(const Vector2 &lp) const final

Return method for full material description of the Surface.

  • from local coordinate on the surface

Parameters:

lp – is the local position used for the (eventual) lookup

Returns:

const MaterialSlab

virtual const MaterialSlab &materialSlab(const Vector3 &gp) const final

Return method for full material description of the Surface.

  • from the global coordinates

Parameters:

gp – is the global position used for the (eventual) lookup

Returns:

const MaterialSlab

inline virtual const MaterialSlab &materialSlab(std::size_t bin0, std::size_t bin1) const final

Direct access via bins to the MaterialSlab.

Parameters:
  • bin0 – is the material bin in dimension 0

  • bin1 – is the material bin in dimension 1

virtual BinnedSurfaceMaterial &operator*=(double scale) final

Scale operator.

Parameters:

scale – is the scale factor for the full material

BinnedSurfaceMaterial &operator=(BinnedSurfaceMaterial &&bsm) = default

Assignment Move operator.

BinnedSurfaceMaterial &operator=(const BinnedSurfaceMaterial &bsm) = default

Assignment operator.

virtual std::ostream &toStream(std::ostream &sl) const final

Output Method for std::ostream, to be overloaded by child classes.

class BinningData

This class holds all the data necessary for the bin calculation.

phi has a very particular behaviour:

  • there’s the change around +/- PI

  • it can be multiplicative or additive multiplicative : each major bin has the same sub structure i.e. first binnning

structure is equidistant additive : sub structure replaces one bin (and one bin only)

Public Functions

BinningData() = default
inline BinningData(BinningOption bOption, BinningValue bValue, const std::vector<float> &bBoundaries, std::unique_ptr<const BinningData> sBinData = nullptr)

Constructor for non-equidistant binning.

Parameters:
  • bOption – is the binning option : open / closed

  • bValue – is the binning value : binX, binY, etc.

  • bBoundaries – are the bin boundaries

  • sBinData – is (optional) sub structure

inline BinningData(BinningOption bOption, BinningValue bValue, std::size_t bBins, float bMin, float bMax, std::unique_ptr<const BinningData> sBinData = nullptr, bool sBinAdditive = false)

Constructor for equidistant binning and optional sub structure can be multiplicative or additive.

Parameters:
  • bOption – is the binning option : open, closed

  • bValue – is the binning value: binX, binY, etc.

  • bBins – is number of equidistant bins

  • bMin – is the minimum value

  • bMax – is the maximum value

  • sBinData – is (optional) sub structure

  • sBinAdditive – is the prescription for the sub structure

inline BinningData(BinningValue bValue, float bMin, float bMax)

Constructor for 0D binning.

Parameters:
  • bValue – is the binning value: binX, binY, etc.

  • bMin – is the minimum value

  • bMax – is the maximum value

inline BinningData(const BinningData &bdata)

Copy constructor.

Parameters:

bdata – is the source object

~BinningData() = default
inline std::size_t bins() const

Return the number of bins - including sub bins.

inline const std::vector<float> &boundaries() const

Return the boundaries - including sub boundaries.

Returns:

vector of floats indicating the boundary values

inline float center(std::size_t bin) const

Get the center value of a bin.

Parameters:

bin – is the bin for which the center value is requested

Returns:

float value according to the bin center

inline float centerValue(std::size_t bin) const

access to the center value this uses the bin boundary vector, it also works with sub structure

Parameters:

bin – is the bin for which the value is requested, if bin > nbins it is set to max

Returns:

the center value of the bin is given

inline bool inside(const Vector2 &lposition) const

Check if bin is inside from Vector2.

Parameters:

lposition – is the search position in global coordinated

Returns:

boolean if this is inside() method is true

inline bool inside(const Vector3 &position) const

Check if bin is inside from Vector3.

Parameters:

position – is the search position in global coordinated

Returns:

boolean if this is inside() method is true

inline int nextDirection(const Vector3 &position, const Vector3 &dir) const

Layer next direction is needed.

Parameters:
  • position – is the start search position

  • dir – is the direction

Returns:

integer that indicates which direction to move

inline BinningData &operator=(const BinningData &bdata)

Assignment operator.

Parameters:

bdata – is the source object

inline bool operator==(const BinningData &bData) const

Equality operator.

Parameters:

bData – is the binning data to be checked against

Returns:

a boolean indicating if they are the same

inline std::size_t search(float value) const

Generic search - forwards to correct function pointer.

Parameters:

value – is the searchvalue as float

Returns:

bin according tot this

inline std::size_t searchGlobal(const Vector3 &position) const

Generic search from a 3D position &#8212; corresponds to global coordinate schema.

Parameters:

position – is the search position in global coordinated

Returns:

bin according tot this

inline std::size_t searchLocal(const Vector2 &lposition) const

Generic search from a 2D position &#8212; corresponds to local coordinate schema.

Parameters:

lposition – is the search position in local coordinated

Returns:

bin according tot this

inline std::size_t searchWithSubStructure(float value) const

Generic search with sub structure.

  • forwards to correct function pointer

Parameters:

value – is the searchvalue as float

Returns:

bin according tot this

inline std::string toString(const std::string &indent) const

String screen output method.

Parameters:

indent – the current indentation

Returns:

a string containing the screen information

inline float value(const Vector2 &lposition) const

Take the right float value.

Parameters:

lposition – assumes the correct local position expression

Returns:

float value according to the binning setup

inline float value(const Vector3 &position) const

Take the right float value.

Parameters:

position – is the global position

Returns:

float value according to the binning setup

inline float width(std::size_t bin) const

Get the width of a bin.

Parameters:

bin – is the bin for which the width is requested

Returns:

float value of width

Public Members

BinningValue binvalue = {}

binning value: binX, binY, binZ, binR …

float max = {}

maximum value

float min = {}

minimum value

BinningOption option = {}

binning option: open, closed

float step = {}

binning step

bool subBinningAdditive = {}

sub structure: additive or multipicative

std::unique_ptr<const BinningData> subBinningData

sub structure: describe some sub binning

BinningType type = {}

binning type: equidistant, arbitrary

bool zdim = {}

zero dimensional binning : direct access

class BoundaryCheck

The BoundaryCheck class provides boundary checks and distance calculations for aligned box-like and polygonal boundaries on local surfaces.

Different types of boundary checks are supported and are transparently selected when calling the isInside(...) and distance(...) methods:

  • Hard checks w/o any tolerances

  • Tolerance-based checks in one or in both local coordinates

  • Chi2-based checks based on a covariance matrix. Non-vanishing correlations are correctly taken into account.

With a defined covariance matrix, the closest point and the distance are not defined along the usual Euclidean metric, but by the Mahalanobis distance induced by the covariance.

Public Types

enum class Type

Values:

enumerator eNone

disable boundary check

enumerator eAbsolute

absolute cut

enumerator eChi2

chi2-based cut with full correlations

Public Functions

inline explicit BoundaryCheck(bool check)

Construct either hard cut in both dimensions or no cut at all.

inline BoundaryCheck(bool checkLocal0, bool checkLocal1, double tolerance0 = 0, double tolerance1 = 0)

Construct a tolerance based check.

Parameters:
  • checkLocal0 – Boolean directive to check coordinate 0

  • checkLocal1 – Boolean directive to check coordinate 1

  • tolerance0 – Tolerance along coordinate 0

  • tolerance1 – Tolerance along coordinate 1

inline BoundaryCheck(const SquareMatrix2 &localCovariance, double sigmaMax = 1)

Construct a chi2-based check.

Parameters:
  • localCovariance – Coverance matrix in local coordinates

  • sigmaMax – Significance for the compatibility test

template<typename Vector2Container>
inline Acts::Vector2 computeClosestPointOnPolygon(const Acts::Vector2 &point, const Vector2Container &vertices) const
inline SquareMatrix2 covariance() const
inline double distance(const Vector2 &point, const Vector2 &lowerLeft, const Vector2 &upperRight) const

Calculate the signed, weighted, closest distance to an aligned box.

If a covariance is defined, the distance is the corresponding Mahalanobis distance. Otherwise, it is the Eucleadian distance.

Parameters:
  • point – Test point

  • lowerLeft – Minimal vertex of the box

  • upperRight – Maximal vertex of the box

Returns:

Negative value if inside, positive if outside

template<typename Vector2Container>
inline double distance(const Vector2 &point, const Vector2Container &vertices) const

Calculate the signed, weighted, closest distance to a polygonal boundary.

If a covariance is defined, the distance is the corresponding Mahalanobis distance. Otherwise, it is the Eucleadian distance.

Parameters:
  • point – Test point

  • vertices – Forward iterable container of convex polygon vertices. Calling std::begin/ std::end on the container must return an iterator where *it must be convertible to an Acts::Vector2.

Returns:

Negative value if inside, positive if outside

inline bool isInside(const Vector2 &point, const Vector2 &lowerLeft, const Vector2 &upperRight) const

Check if the point is inside a box aligned with the local axes.

The check takes into account whether tolerances or covariances are defined for the boundary check.

Parameters:
  • point – Test point

  • lowerLeft – Minimal vertex of the box

  • upperRight – Maximal vertex of the box

template<typename Vector2Container>
inline bool isInside(const Vector2 &point, const Vector2Container &vertices) const

Check if the point is inside a polygon.

The check takes into account whether tolerances or covariances are defined for the boundary check.

Parameters:
  • point – Test point

  • vertices – Forward iterable container of convex polygon vertices. Calling std::begin/ std::end on the container must return an iterator where *it must be convertible to an Acts::Vector2.

inline operator bool() const
inline bool operator!() const
inline const Vector2 &tolerance() const
inline Type type() const

Broadcast the type.

class ConeBounds : public Acts::SurfaceBounds

Bounds for a conical surface, the opening angle is stored in \( \tan(\alpha) \) and always positively defined.

The cone can open to both sides, steered by \( z_min \) and \( z_max \).

../_images/ConeBounds.gif

Public Types

enum BoundValues

Values:

enumerator eAlpha
enumerator eMinZ
enumerator eMaxZ
enumerator eHalfPhiSector
enumerator eAveragePhi
enumerator eSize

Public Functions

ConeBounds() = delete
ConeBounds(const std::array<double, eSize> &values) noexcept(false)

Constructor - from parameters array.

Parameters:

values – The parameter array

ConeBounds(double alpha, bool symm, double halfphi = M_PI, double avphi = 0.) noexcept(false)

Constructor - open cone with alpha, by default a full cone but optionally can make a conical section.

Parameters:
  • alpha – is the opening angle of the cone

  • symm – is the boolean indicating if the cone is symmetric in +/- z

  • halfphi – is the half opening angle (default is pi)

  • avphi – is the phi value around which the bounds are opened (default=0)

ConeBounds(double alpha, double minz, double maxz, double halfphi = M_PI, double avphi = 0.) noexcept(false)

Constructor - open cone with alpha, minz and maxz, by default a full cone but can optionally make it a conical section.

Parameters:
  • alpha – is the opening angle of the cone

  • minz – cone expanding from minimal z

  • maxz – cone expanding to maximal z

  • halfphi – is the half opening angle (default is pi)

  • avphi – is the phi value around which the bounds are opened (default=0)

~ConeBounds() override = default
inline double get(BoundValues bValue) const

Access to the bound values.

Parameters:

bValue – the class nested enum for the array access

virtual bool inside(const Vector2 &lposition, const BoundaryCheck &bcheck = BoundaryCheck(true)) const final

inside method for local position

Parameters:
  • lposition – is the local position to be checked

  • bcheck – is the boundary check directive

Returns:

is a boolean indicating if the position is inside

inline double r(double z) const

Return the radius at a specific z values.

Parameters:

z – is the z value for which r is requested

Returns:

is the r value associated with z

inline double tanAlpha() const

Return tangent of alpha (pre-computed)

virtual std::ostream &toStream(std::ostream &sl) const final

Output Method for std::ostream.

Parameters:

sl – is the ostrea into which the dump is done

Returns:

is the input object

virtual BoundsType type() const final

Return the bounds type - for persistency optimization.

Returns:

is a BoundsType enum

inline virtual std::vector<double> values() const final

Return the bound values as dynamically sized vector.

Returns:

this returns a copy of the internal values

class ConeLayer : public virtual Acts::ConeSurface, public Acts::Layer

Class to describe a conical detector layer for tracking, it inhertis from both, Layer base class and ConeSurface class.

Public Functions

ConeLayer() = delete
ConeLayer(const ConeLayer &cla) = delete
~ConeLayer() override = default
ConeLayer &operator=(const ConeLayer&) = delete
virtual const ConeSurface &surfaceRepresentation() const override

Transforms the layer into a Surface representation for extrapolation.

virtual ConeSurface &surfaceRepresentation() override

Public Static Functions

static inline MutableLayerPtr create(const Transform3 &transform, std::shared_ptr<const ConeBounds> cbounds, std::unique_ptr<SurfaceArray> surfaceArray, double thickness = 0., std::unique_ptr<ApproachDescriptor> ad = nullptr, LayerType laytyp = Acts::active)

Factory for shared layer.

Parameters:
  • transform – is the 3D transform that positions the layer in 3D frame

  • cbounds – is the conical bound description

  • surfaceArray – is the array of sensitive surfaces

  • thickness – is the layer thickness along the normal axis

  • ad – is the approach descriptor for navigation towards the layer

  • laytyp – is the layer type

Returns:

is a shared pointer to a layer

class ConeSurface : public Acts::RegularSurface

Class for a conical surface in the Tracking geometry.

It inherits from Surface.

The ConeSurface is special since no corresponding Track parameters exist since they’re numerical instable at the tip of the cone. Propagations to a cone surface will be returned in curvilinear coordinates.

Subclassed by Acts::ConeLayer

Public Functions

ConeSurface() = delete
~ConeSurface() override = default
virtual AlignmentToPathMatrix alignmentToPathDerivative(const GeometryContext &gctx, const FreeVector &parameters) const final

Calculate the derivative of path length at the geometry constraint or point-of-closest-approach w.r.t.

alignment parameters of the surface (i.e. local frame origin in global 3D Cartesian coordinates and its rotation represented with extrinsic Euler angles)

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • parameters – is the free parameters

Returns:

Derivative of path length w.r.t. the alignment parameters

virtual Vector3 binningPosition(const GeometryContext &gctx, BinningValue bValue) const final

The binning position method - is overloaded for r-type binning.

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • bValue – defines the type of binning applied in the global frame

Returns:

The return type is a vector for positioning in the global frame

virtual const ConeBounds &bounds() const final

This method returns the ConeBounds by reference.

Result<Vector2> globalToLocal(const GeometryContext &gctx, const Vector3 &position, const Vector3 &direction, double tolerance = s_onSurfaceTolerance) const final

Convert a global position to a local one this is the most generic interface, which is implemented by all surfaces.

Note

The position is required to be on-surface, which is indicated by the Result return value.

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • position – is the global position to be converted

  • direction – is the direction of the local position (ignored for RegularSurface)

  • tolerance – is the tolerance for the on-surface check

Returns:

Result type containing local position by value

virtual Result<Vector2> globalToLocal(const GeometryContext &gctx, const Vector3 &position, double tolerance = s_onSurfaceTolerance) const final

Global to local transformation.

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • position – is the global position to be transformed

  • tolerance – optional tolerance within which a point is considered valid on surface

Returns:

a Result<Vector2> which can be !ok() if the operation fails

Result<Vector2> globalToLocal(const GeometryContext &gctx, const Vector3 &position, double tolerance = s_onSurfaceTolerance) const = 0

Convert a global position to a local one.

Note

The position is required to be on-surface, which is indicated by the Result return value.

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • position – is the global position to be converted

  • tolerance – is the tolerance for the on-surface check

Returns:

Result type containing local position by value

virtual SurfaceMultiIntersection intersect(const GeometryContext &gctx, const Vector3 &position, const Vector3 &direction, const BoundaryCheck &bcheck = BoundaryCheck(false), double tolerance = s_onSurfaceTolerance) const final

Straight line intersection schema from position/direction.

If possible returns both solutions for the cylinder

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • position – The position to start from

  • direction – The direction at start

  • bcheck – the Boundary Check

  • tolerance – the tolerance used for the intersection

Returns:

SurfaceMultiIntersection object (contains intersection & surface)

virtual ActsMatrix<2, 3> localCartesianToBoundLocalDerivative(const GeometryContext &gctx, const Vector3 &position) const final

Calculate the derivative of bound track parameters local position w.r.t.

position in local 3D Cartesian coordinates

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • position – The position of the parameters in global

Returns:

Derivative of bound local position w.r.t. position in local 3D cartesian coordinates

virtual Vector3 localToGlobal(const GeometryContext &gctx, const Vector2 &lposition) const final

Local to global transformation.

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • lposition – is the local position to be transformed

Returns:

The global position by value

Vector3 localToGlobal(const GeometryContext &gctx, const Vector2 &lposition) const = 0

Local to global transformation.

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • lposition – local 2D position in specialized surface frame

Returns:

The global position by value

Vector3 localToGlobal(const GeometryContext &gctx, const Vector2 &lposition, const Vector3 &direction) const final

Local to global transformation.

This is the most generic interface, which is implemented by all surfaces.

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • lposition – local 2D position in specialized surface frame

  • direction – global 3D momentum direction (ignored for RegularSurface)

Returns:

The global position by value

virtual std::string name() const override

Return properly formatted class name for screen output.

virtual Vector3 normal(const GeometryContext &gctx, const Vector2 &lposition) const final

Return method for surface normal information.

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • lposition – is the local position at normal vector request

Returns:

Vector3 normal vector in global frame

Vector3 normal(const GeometryContext &gctx, const Vector2 &lposition) const = 0

Calculate the normal vector of the surface This overload requires an on-surface local position.

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • lposition – is the local position where the normal vector is constructed

Returns:

normal vector by value

Vector3 normal(const GeometryContext &gctx, const Vector3 &pos, const Vector3 &direction) const final

Calculate the normal vector of the surface This overload is fully generic, fulfills the Surface interface and accepts a global position and a direction.

For RegularSurface this is equivalent to the normal overload, ignoring the direction

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • pos – is the global position where the normal vector is constructed

  • direction – is the direction of the normal vector (ignored for RegularSurface)

virtual Vector3 normal(const GeometryContext &gctx, const Vector3 &position) const final

Return method for surface normal information.

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • position – is the global position as normal vector base

Returns:

Vector3 normal vector in global frame

Vector3 normal(const GeometryContext &gctx, const Vector3 &position) const = 0

Calculate the normal vector of the surface This overload accepts a global position.

Parameters:
  • position – is the global position where the normal vector is constructed

  • gctx – The current geometry context object, e.g. alignment

Returns:

normal vector by value

ConeSurface &operator=(const ConeSurface &other)

Assignment operator.

Parameters:

other – is the source surface for the assignment

virtual double pathCorrection(const GeometryContext &gctx, const Vector3 &position, const Vector3 &direction) const final

The pathCorrection for derived classes with thickness.

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • position – is the global potion at the correction point

  • direction – is the momentum direction at the correction point

Returns:

is the path correction due to incident angle

virtual Polyhedron polyhedronRepresentation(const GeometryContext &gctx, std::size_t lseg) const override

Return a Polyhedron for the surfaces.

Note

that a surface transform can invalidate the extrema in the transformed space

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • lseg – Number of segments along curved lines, it represents the full 2*M_PI coverange, if lseg is set to 1 only the extrema are given

Returns:

A list of vertices and a face/facett description of it

virtual RotationMatrix3 referenceFrame(const GeometryContext &gctx, const Vector3 &position, const Vector3 &direction) const final

Return the measurement frame - this is needed for alignment, in particular for StraightLine and Perigee Surface.

  • the default implementation is the RotationMatrix3 of the transform

Parameters:
  • gctx – The current geometry context object, e.g. alignment

  • position – is the global position where the measurement frame is constructed

  • direction – is the momentum direction used for the measurement frame construction

Returns:

matrix that indicates the measurement frame

virtual Vector3 rotSymmetryAxis(const GeometryContext &gctx) const
Parameters:

gctx – The current geometry context object, e.g. alignment

virtual SurfaceType type() const override

Return the surface type.

class ConeVolumeBounds : public Acts::VolumeBounds

Volume bound class for describing conical volumes either with cylindrical inlay or outer boundary, it also allows for a sectoral description.

Public Types

enum BoundValues

for readability

Values:

enumerator eInnerAlpha
enumerator eInnerOffsetZ
enumerator eOuterAlpha
enumerator eOuterOffsetZ
enumerator eHalfLengthZ
enumerator eAveragePhi
enumerator eHalfPhiSector
enumerator eSize

Public Functions

ConeVolumeBounds() = delete
ConeVolumeBounds(const ConeVolumeBounds &cobo) = default
inline ConeVolumeBounds(const std::array<double, eSize> &values) noexcept(false)

Constructor - from a fixed size array.

Parameters:

values – The bound values

ConeVolumeBounds(double cylinderR, double alpha, double offsetZ, double halflengthZ, double averagePhi, double halfPhiSector) noexcept(false)

Constructor - for general cylidner-cone setups.

Note

depending on cylinderR > coneR it is constructing a cone with cylindrical cutout or a cylinder with conical cutout

Parameters:
  • cylinderR – The inner radius of the cylinder

  • alpha – The opening angle of the cone (0 if no cone)

  • offsetZ – The tip z position in of the cone, w.r.t center

  • halflengthZ – The minimum z value of the inner and outer cones

  • averagePhi – The phi orientation of the sector (defaulted to 0)

  • halfPhiSector – The opening angle phi sector

ConeVolumeBounds(double innerAlpha, double innerOffsetZ, double outerAlpha, double outerOffsetZ, double halflengthZ, double averagePhi, double halfPhiSector) noexcept(false)

Constructor - for general cone-cone setups.

Parameters:
  • innerAlpha – The opening angle of the inner cone (0 if no cone)

  • innerOffsetZ – The tip z position in of the inner cone, w.r.t center

  • outerAlpha – The opening angle of the outer cone (0 if no cone)

  • outerOffsetZ – The tip z position in of the outer cone, w.r.t center

  • halflengthZ – The minimum z value of the inner and outer cones

  • averagePhi – The phi orientation of the sector

  • halfPhiSector – The opening angle phi sector

~ConeVolumeBounds() override = default
virtual Volume::BoundingBox boundingBox(const Transform3 *trf = nullptr, const Vector3 &envelope = {0, 0, 0}, const Volume *entity = nullptr) const final

Construct bounding box for this shape.

Parameters:
  • trf – Optional transform

  • envelope – Optional envelope to add / subtract from min/max

  • entity – Entity to associate this bounding box with

Returns:

Constructed bounding box

inline double get(BoundValues bValue) const

Access to the bound values.

Parameters:

bValue – the class nested enum for the array access

inline double innerRmax() const
inline double innerRmin() const
inline double innerTanAlpha() const
virtual bool inside(const Vector3 &pos, double tol = 0.) const final

This method checks if position in the 3D volume frame is inside the cylinder.

Parameters:
  • pos – is the position in volume frame to be checked

  • tol – is the absolute tolerance to be applied

ConeVolumeBounds &operator=(const ConeVolumeBounds &cobo) = default
virtual OrientedSurfaces orientedSurfaces(const Transform3 &transform = Transform3::Identity()) const final

Oriented surfaces, i.e.

the decomposed boundary surfaces and the according navigation direction into the volume given the normal vector on the surface

It will throw an exception if the orientation prescription is not adequate

Parameters:

transform – is the 3D transform to be applied to the boundary surfaces to position them in 3D space

Returns:

a vector of surfaces bounding this volume

inline double outerRmax() const
inline double outerRmin() const
inline double outerTanAlpha() const
virtual std::ostream &toStream(std::ostream &sl) const final

Output Method for std::ostream.

Parameters:

sl – is ostream operator to be dumped into

inline virtual VolumeBounds::BoundsType type() const final

Return the bounds type - for persistency optimization.

Returns:

is a BoundsType enum

inline virtual std::vector<double> values() const final

Return the bound values as dynamically sized vector.

Returns:

this returns a copy of the internal values

class ConstantBField : public Acts::MagneticFieldProvider

This class implements a simple constant magnetic field.

The magnetic field value has to be set at creation time, but can be updated later on.

Public Functions

inline explicit ConstantBField(Vector3 B)

Construct constant magnetic field from field vector.

Parameters:

B[in] magnetic field vector in global coordinate system

inline Vector3 getField() const

Get the B field at a position.

inline virtual Result<Vector3> getField(const Vector3 &position, MagneticFieldProvider::Cache &cache) const override

Retrieve magnetic field value at a given location.

Requires a cache object created through makeCache().

Note

The position is ignored and only kept as argument to provide a consistent interface with other magnetic field services.

Parameters:
  • position[in] global 3D position for the lookup

  • cache[inout] Field provider specific cache object

Returns:

magnetic field vector at given position

inline virtual Result<Vector3> getFieldGradient(const Vector3 &position, ActsMatrix<3, 3> &derivative, MagneticFieldProvider::Cache &cache) const override

Retrieve magnetic field value its its gradient.

Requires a cache object created through makeCache().

Note

The position is ignored and only kept as argument to provide a consistent interface with other magnetic field services.

Note

currently the derivative is not calculated

Parameters:
  • position[in] global 3D position

  • derivative[out] gradient of magnetic field vector as (3x3) matrix

  • cache[inout] Field provider specific cache object

Returns:

magnetic field vector

inline bool isInside(const Vector3 &position) const

check whether given 3D position is inside look-up domain

Parameters:

position[in] global 3D position

Returns:

Always true for constant magnetic field

inline virtual Acts::MagneticFieldProvider::Cache makeCache(const Acts::MagneticFieldContext &mctx) const override

Make an opaque cache for the magnetic field.

Instructs the specific implementation to generate a Cache instance for magnetic field lookup.

Parameters:

mctx – The magnetic field context to generate cache for

Returns:

Cache The opaque cache object

inline void setField(const Vector3 &B)

update magnetic field vector

Parameters:

B[in] magnetic field vector in global coordinate system

struct Cache

Public Functions

inline Cache(const MagneticFieldContext&)

constructor with context

class ConstrainedStep

A constrained step class for the steppers.

This class is symmetrical for forward and backward propagation. The sign of the propagation direction should not enter here but rather be applied the step is actually taken.

As simple as this class looks it hides a few very important details:

  • Overstepping handling. The step size sign will flip if we happened to pass our target.

  • Convergence handling. Smaller and smaller step sizes have to be used in order to converge on a target.

Because of the points mentioned above, the update function will always prefer negative step sizes. A side effect of this is that we will propagate in the opposite direction if the target is “behind us”.

The hierarchy is:

  • Overstepping resolution / backpropagation

  • Convergence

  • Step into the void with std::numeric_limits<Scalar>max()

Public Types

using Scalar = ActsScalar
enum Type

the types of constraints from actor - this would be a typical navigation step from aborter - this would be a target condition from user - this is user given for what reason ever

Values:

enumerator actor
enumerator aborter
enumerator user

Public Functions

constexpr ConstrainedStep() = default
inline explicit constexpr ConstrainedStep(Scalar value)

constructor from Scalar

Parameters:

value – is the user given initial value

inline constexpr Scalar accuracy() const

Access the accuracy value.

inline constexpr void release(Type type)

release a certain constraint value

Parameters:

type – is the constraint type to be released

inline constexpr void releaseAccuracy()

release accuracy

inline constexpr void setAccuracy(Scalar value)

set accuracy by one Scalar

this will set only the accuracy, as this is the most exposed to the Propagator

Parameters:

value – is the new accuracy value

inline constexpr void setUser(Scalar value)

set user by one Scalar

Parameters:

value – is the new user value

inline std::ostream &toStream(std::ostream &os) const
inline std::string toString() const
inline constexpr void update(Scalar value, Type type, bool releaseStep = false)

Update the step size of a certain type.

Only navigation and target abortion step size updates may change the sign due to overstepping

Parameters:
  • value – is the new value to be updated

  • type – is the constraint type

  • releaseStep – Allow step size to increase again

inline constexpr Scalar value() const

returns the min step size

inline constexpr Scalar value(Type type) const

Access a specific value.

Parameters:

type – is the requested parameter type

Public Members

std::size_t nStepTrials = std::numeric_limits<std::size_t>::max()

Number of iterations needed by the stepsize finder (e.g.

Runge-Kutta) of the stepper.

template<int N>
class ConvexPolygonBounds : public Acts::ConvexPolygonBoundsBase

This is the actual implementation of the bounds.

It is templated on the number of vertices, but there is a specialization for dynamic number of vertices, where the underlying storage is then a vector.

Template Parameters:

N – Number of vertices

Public Types

using value_array = std::array<double, eSize>

Type that’s used to store the vertices, in this case a fixed size array.

using vertex_array = std::array<Vector2, num_vertices>

Type that’s used to store the vertices, in this case a fixed size array.

Public Functions

ConvexPolygonBounds() = delete
ConvexPolygonBounds(const std::vector<Vector2> &vertices) noexcept(false)

Constructor from a vector of vertices, to facilitate construction.

This will throw if the vector size does not match num_vertices. This will throw if the vertices do not form a convex polygon.

Parameters:

vertices – The list of vertices.

ConvexPolygonBounds(const value_array &values) noexcept(false)

Constructor from a fixed size array of parameters This will throw if the vertices do not form a convex polygon.

Parameters:

values – The values to build up the vertices

ConvexPolygonBounds(const vertex_array &vertices) noexcept(false)

Constructor from a fixed size array of vertices.

This will throw if the vertices do not form a convex polygon.

Parameters:

vertices – The vertices

~ConvexPolygonBounds() override = default
virtual const RectangleBounds &boundingBox() const final

Return a rectangle bounds object that encloses this polygon.

Returns:

The rectangular bounds

virtual bool inside(const Vector2 &lposition, const BoundaryCheck &bcheck) const final

Return whether a local 2D point lies inside of the bounds defined by this object.

Parameters:
  • lposition – The local position to check

  • bcheck – The BoundaryCheck object handling tolerances.

Returns:

Whether the points is inside

virtual BoundsType type() const final

Return the bounds type - for persistency optimization.

Returns:

is a BoundsType enum

virtual std::vector<Vector2> vertices(unsigned int lseg = 1) const final

Return the vertices.

Note

the number of segments is ignored in this representation

Parameters:

lseg – the number of segments used to approximate and eventually curved line

Returns:

vector for vertices in 2D

Public Static Attributes

static constexpr std::size_t eSize = 2 * N

Expose number of parameters as a template parameter.

static constexpr std::size_t num_vertices = N

Expose number of vertices given as template parameter.

class ConvexPolygonBoundsBase : public Acts::PlanarBounds

base class for convex polygon bounds

This class serves as a base class for the actual bounds class. The only deriving type is the templated ConvexPolygonBounds.

Subclassed by Acts::ConvexPolygonBounds< N >, Acts::ConvexPolygonBounds< PolygonDynamic >

Public Functions

virtual std::ostream &toStream(std::ostream &sl) const final

Output Method for std::ostream.

Parameters:

sl – is the ostream to be written into

virtual std::vector<double> values() const final

Return the bound values as dynamically sized vector.

Returns:

this returns a copy of the internal values

class CuboidVolumeBounds : public Acts::VolumeBounds

Bounds for a cubical Volume, the orientedSurfaces(…) method creates a vector of 6 surfaces:

BoundarySurfaceFace [index]:

  • negativeFaceXY [0] : Rectangular Acts::PlaneSurface, parallel to \( xy \) plane at negative \( z \)

  • positiveFaceXY [1] : Rectangular Acts::PlaneSurface, parallel to \( xy \) plane at positive \( z \)

  • negativeFaceXY [2] : Rectangular Acts::PlaneSurface, attached to \( yz \) plane at negative \( x \)

  • positiveFaceXY [3] : Rectangular Acts::PlaneSurface, attached to \( yz \) plane at negative \( x \)

  • negativeFaceXY [4] : Rectangular Acts::PlaneSurface, parallel to \( zx \) plane at negative \( y \)

  • positiveFaceXY [5] : Rectangular Acts::PlaneSurface, parallel to \( zx \) plane at positive \( y \)

Public Types

enum BoundValues

for streaming and access

Values:

enumerator eHalfLengthX
enumerator eHalfLengthY
enumerator eHalfLengthZ
enumerator eSize

Public Functions

CuboidVolumeBounds() = delete
CuboidVolumeBounds(const CuboidVolumeBounds &bobo)

Copy Constructor.

Parameters:

bobo – is the source volume bounds to be copied

inline CuboidVolumeBounds(const std::array<double, eSize> &values) noexcept(false)

Constructor - from a fixed size array.

Parameters:

values – iw the bound values

CuboidVolumeBounds(double halex, double haley, double halez) noexcept(false)

Constructor - the box boundaries.

Parameters:
  • halex – is the half length of the cube in x

  • haley – is the half length of the cube in y

  • halez – is the half length of the cube in z

~CuboidVolumeBounds() override = default
virtual double binningBorder(BinningValue bValue) const final

Binning borders in double.

Parameters:

bValue – is the binning schema used

Returns:

float offset to be used for the binning

virtual Volume::BoundingBox boundingBox(const Transform3 *trf = nullptr, const Vector3 &envelope = {0, 0, 0}, const Volume *entity = nullptr) const final

Construct bounding box for this shape.

Parameters:
  • trf – Optional transform

  • envelope – Optional envelope to add / subtract from min/max

  • entity – Entity to associate this bounding box with

Returns:

Constructed bounding box

inline double get(BoundValues bValue) const

Access to the bound values.

Parameters:

bValue – the class nested enum for the array access

inline virtual bool inside(const Vector3 &pos, double tol = 0.) const override

This method checks if position in the 3D volume frame is inside the cylinder.

Parameters:
  • pos – is the position in volume frame to be checked

  • tol – is the absolute tolerance to be applied

CuboidVolumeBounds &operator=(const CuboidVolumeBounds &bobo)

Assignment operator.

Parameters:

bobo – is the source volume bounds to be assigned

virtual OrientedSurfaces orientedSurfaces(const Transform3 &transform = Transform3::Identity()) const override

Oriented surfaces, i.e.

the decomposed boundary surfaces and the according navigation direction into the volume given the normal vector on the surface

It will throw an exception if the orientation prescription is not adequate

Parameters:

transform – is the 3D transform to be applied to the boundary surfaces to position them in 3D space

Returns:

a vector of surfaces bounding this volume

virtual std::ostream &toStream(std::ostream &sl) const override

Output Method for std::ostream.

Parameters:

sl – is ostream operator to be dumped into

inline virtual VolumeBounds::BoundsType type() const final

Return the bounds type - for persistency optimization.

Returns:

is a BoundsType enum

inline virtual std::vector<double> values() const final

Return the bound values as dynamically sized vector.

Returns:

this returns a copy of the internal values

class CutoutCylinderVolumeBounds : public Acts::VolumeBounds

Class which implements a cutout cylinder.

This shape is basically a cylinder, with another, smaller cylinder subtracted from the center. ——————&#8212; rmax | | | |——&#8212;| | rmed | | | | —&#8212; —&#8212; rmin &#8212; hlZc &#8212; ——&#8212; hlZ —-&#8212;

Public Types

enum BoundValues

for streaming and access

Values:

enumerator eMinR
enumerator eMedR
enumerator eMaxR
enumerator eHalfLengthZ
enumerator eHalfLengthZcutout
enumerator eSize

Public Functions

CutoutCylinderVolumeBounds() = delete
inline CutoutCylinderVolumeBounds(const std::array<double, eSize> &values) noexcept(false)

Constructor - from a fixed size array.

Parameters:

values – The bound values

inline CutoutCylinderVolumeBounds(double rmin, double rmed, double rmax, double hlZ, double hlZc) noexcept(false)

Constructor from defining parameters.

Parameters:
  • rmin – Minimum radius at the “choke points”

  • rmed – The medium radius (outer radius of the cutout)

  • rmax – The outer radius of the overall shape

  • hlZ – The longer halflength of the shape

  • hlZc – The cutout halflength of the shape

~CutoutCylinderVolumeBounds() override = default
virtual Volume::BoundingBox boundingBox(const Transform3 *trf = nullptr, const Vector3 &envelope = {0, 0, 0}, const Volume *entity = nullptr) const final

Construct bounding box for this shape.

Parameters:
  • trf – Optional transform

  • envelope – Opt