(seeding_core)=
To reduce the time needed to reconstruct particle tracks, a track seed (henceforth: seed) is created which serves as the initial direction for the track reconstruction algorithm (henceforth: the tracking). The tracking then tries to find all measurements belonging to a single particle in this direction in order to reconstruct the track. This means, if no seed exists for a particle, this particle will not be reconstructed. On the other hand, finding too many seeds which either correspond to particles with already existing seeds or which do not correspond to particles at all increases the time needed for tracking.
A good seeding algorithm, therefore, has the following properties:
- It finds at least one seed for each particle that should be found
- It doesn't find many seeds which do NOT correspond to particles
- It doesn't find many seeds per particle
The most typical way to create seeds is to combine measurements. In a homogeneous magnetic field, 3 measurements perfectly describe the helical path of a charged particle. One such triplet of measurements would then constitute a seed and defines, in close bounds, where the tracking needs to look for additional measurements to create a track spanning the whole detector. The difficulty is in choosing the correct measurements, as a helix can be fitted through any 3 measurements in a collision event with potentially tens of thousands of measurements. Therefore, many constraints or “cuts” are defined to reduce the number of candidates. Cuts may define where particles originate or the range of energy of particles to be found or otherwise restrict the combination of measurements for seed creation.
The seeding implementation in Core/include/Acts/Seeding/
was written with a focus on parallelism and maintainability and as detector agnostic as possible, only assuming a (near) homogeneous magnetic field with particles originating from the central detector region. Cuts are configurable and can be plugged in as an algorithm which is called by the seeding. The seeding works on measurements or “SpacePoints” (SP), which need to provide
:::{figure} figures/seeding/3Dcoordinates.svg
:width: 550px
:align: center
Sketch of the detector with 3 layers. The interaction region is supposed to be located along the
:::{attention} Note that the seeding algorithm breaks down for particles with a particle track whose helix diameter is smaller than the detector radius until which seeds are to be created. This is due to ordering assumptions of SP locations as well as due to approximations which become inaccurate for lower energy particles. :::
The SPs in each detector layer are projected on a rectangular grid of
configurable granularity. The search for seed starts by selecting SP in the
middle detector layer. Then matching SPs are searched in the inner and outer
layers. Grouping of the SPs in the aforementioned grid allows to limit the
search to neighbouring grid cells thus improving significantly algorithm
performance (see {numref}tripletsFormation
). The number of neighboring bins
used in the SP search can be defined separately for the bottom and top layer
SPs in the
(tripletsFormation)=
:::{figure} figures/seeding/tripletsFormation.svg
:width: 450px
:align: center
Representation of the search for triplet combinations in the
The {func}Acts::SeedFinder::createSeedsForGroup
function receives three iterators
over SPs constructed from detector layers of increasing radii. The seedfinder will
then attempt to create seeds, with each seed containing exactly one SP returned by
each of the three iterators. It starts by iterating over SPs in the middle layer
(2nd iterator), and within this loop separately iterates once over the bottom SP
and once over the top SP. Within each of the nested loops, SP pairs are tested for
compatibility by applying a set of configurable cuts that can be tested with
two SP only (pseudorapidity, origin along
:::{doxygenfunction} Acts::SeedFinder::createSeedsForGroup(const Acts::SeedFinderOptions &options, SeedingState &state, const grid_t &grid, std::back_insert_iterator<container_t<Seed<external_spacepoint_t>>> outIt, const sp_range_t &bottomSPs, const std::size_t middleSPs, const sp_range_t &topSPs, const Acts::Range1D &rMiddleSPRange) const :::
:::{doxygenfunction} Acts::SeedFinder::createSeedsForGroup(const Acts::SeedFinderOptions &options, const grid_t &grid, const sp_range_t &bottomSPs, const std::size_t middleSPs, const sp_range_t &topSPs) const :::
For all pairs passing the selection the triplets of bottom-middle-top SPs are formed. Each triplet is then confronted with the helix hypothesis. In order to perform calculations only once, the circle calculation is spread out over the three loops.
(x-yCoordinates)= :::{figure} figures/seeding/x-yCoordinates.svg :width: 400px :align: center The x-y projection of the detector with the charged particle helical track originating from the centre of the detector. Signals left by the passage of the track through the detector layers are marked with green crosses. :::
From the helix circle (see {numref}x-yCoordinates
), particle energy and impact parameters can be estimated.
To calculate the helix circle in the
\begin{equation*} u = \frac{x}{x^2+y^2}, \quad \quad v = \frac{y}{x^2+y^2} , \end{equation*}
where the circle containing the three SPs are transformed into a line with equation
\begin{equation*} A = \frac{v_t-v_b}{u_t-u_b} , \end{equation*}
Then,
Inserting the coefficients in the circle equation and assuming that the circle goes through the origin we obtain:
\begin{equation*} (2R)^2 = \frac{A^2+1}{B^2} . \end{equation*}
Now we can apply a cut on the estimate of the minimum helix diameter minHelixDiameter2
without the extra overhead of conversions or computationally complex calculations. The seed is accepted if
\begin{equation*} \frac{A^2+1}{B^2} > (2 R^{min})^2 = \left ( \frac{2 \cdot p_T^{min}}{300 \cdot B_z} \right)^2 \equiv \textnormal{minHelixDiameter2} , \end{equation*}
where
(r-zCoordinates)= :::{figure} figures/seeding/r-zCoordinates.svg :width: 500px :align: center The r-z projection of the detector with the same charged particle track. The track is depicted with the same colours as in the previous figure. :::
The track is not an ideal helix. At each detector layer (or any other material)
scattering may occur making the helix approximate.
The algorithm will check if the triplet forms a nearly straight line
in the r-zCoordinates
) as the particle path in the
\begin{equation*} \left (\cot \theta_b - \cot \theta_t \right )^2 < \sigma^2_{p_T^{min}} + \sigma_f^2, \end{equation*}
This check takes into account the squared uncertainty in the difference between slopes (scatteringInRegion2
(sigmaScattering
, the configurable number of sigmas of scattering angle
to be considered, and maxScatteringAngle2
, which is evaluated from the
Lynch & Dahl correction1
\begin{equation*} \Theta_0^{min} = \frac{13.6 \text{MeV}}{p_T^{min}} \cdot \sqrt{\frac{z_q^2 L}{\beta^2 L_0}} \left(1+0.038 \ln \left(\frac{z_q^2 L}{\beta^2 L_0}\right) \right), \end{equation*}
The second part check against the multiple scattering term p2scatterSigma
(
\begin{equation*} \left (\cot \theta_b - \cot \theta_t \right ) ^2 < \sigma^2_{p_T^{estimated}} + \sigma_f^2, \end{equation*}
The last cut applied in this function is on the transverse impact parameter (or DCA -
distance of closest approach), which is the distance of the perigee of a track from
the interaction region in
(impactParameter)=
:::{figure} figures/seeding/impactParameter.svg
:width: 400px
:align: center
Helix representation in
Assuming the middle layer SP is in the origin of the impactParameter
.
The distance between the centre of the helix and the interaction point (IP) is given by
\begin{equation*}
(x_0 + r_m)^2 + y_0^2 = (R + d_0)^2 \quad \xrightarrow{R^2 = x_0^2 + y_0^2} \quad \frac{d_0^2}{R^2} + 2 \frac{d_0}{R} = \frac{2 x_0 r_m + r_m^2}{R^2} .
\end{equation*}
Considering that
\begin{equation*} d_0 \leq \left| \left( A - B \cdot r_M \right) \cdot r_M \right| \end{equation*}
After creating the potential seeds we apply a seed filter procedure that compares the seeds with other SPs compatible with the seed curvature.
This process ranks the potential seeds based on certain quality criteria and selects the ones that are more likely to produce high-quality tracks
The filter is divided into two functions {func}Acts::SeedFilter::filterSeeds_2SpFixed
and {func}Acts::SeedFilter::filterSeeds_1SpFixed
.
The first function compares the middle and bottom layer SPs of the seeds to other top layer SPs; seeds only differing in top SP are compatible if they have similar helix radius with the same sign (i.e. the same charge). The SPs must have a minimum distance in detector radius, such that SPs from the same layer cannot be considered compatible. The second function iterates over the seeds with only a common middle layer SP and selects the higher quality combinations.
:::{doxygenfunction} Acts::SeedFilter::filterSeeds_2SpFixed :outline: :::
This function assigns a weight (which should correspond to the likelihood that a seed is good) to all seeds and applies detector-specific selection of seeds based on weights. The weight is a “soft cut”, which means that it is only used to discard tracks if many seeds are created for the same middle SP. This process is important to improving computational performance and the quality of the final track collections by rejecting lower-quality seeds.
The weight is calculated by: \begin{equation*} w = (c_1 \cdot N_{t} - c_2 \cdot d_0 - c_3 |z_0| ) + \textnormal{detector specific cuts}. \end{equation*}
The transverse (
The {func}Acts::SeedFilter::filterSeeds_2SpFixed
function also includes a configurable {struct}Acts::SeedConfirmationRangeConfig
seed confirmation step that, when enabled,
classifies higher quality seeds as "quality confined" seeds if they fall within a predefined range of parameters (
The seed confirmation also sets a limit on the number of seeds produced for each middle SP, which retains only the higher quality seeds. If this limit is exceeded, the algorithm checks if there is any low-quality seed in the seed container of this middle SP that can be removed.
:::{doxygenfunction} Acts::SeedFilter::filterSeeds_1SpFixed (Acts::SpacePointData &spacePointData, CandidatesForMiddleSp<const InternalSpacePoint<external_spacepoint_t>>&,const std::size_t,std::back_insert_iterator<std::vector<Seed<external_spacepoint_t>>>) const :outline: :::
This function allows the detector-specific cuts to filter based on all seeds with a common middle SP and limits the number of seeds per middle SP to the configured limit. It sorts the seeds by weight and, to achieve a well-defined ordering in the rare case weights are equal, sorts them by location. The ordering by location is only done to make sure reimplementations (such as the GPU code) are comparable and return the bitwise exactly the same result.
Footnotes
-
M. Tanabashi et al. (Particle Data Group), Passage of Particles Through Matter, Phys. Rev. D 98 0300001 (2018) 2. ↩
-
G.R. Lynch and O.I Dahl, Nucl. Instrum. Methods B58, Phys. Rev. D 98 6 (1991) 2. ↩