[1,2]Jenny Taylor

1]GSI-FFN, GSI Helmholtzzentrum fur Schwerionenforschung GmbH, Planckstraße 1,
Darmstadt, 642 91, Hesse, Germany

2]Department of Physics and Astronomy, Uppsala University, Lägerhyddsvagen 1,
Uppsala, 752 37, Sweden

#  KinFit
A Kinematic Fitting Package for Hadron Physics Experiments

Waleed Esmail<EMAIL_ADDRESS>Jana Rieger<EMAIL_ADDRESS><EMAIL_ADDRESS>Malin Bohman  Karin Schönning  [  [

###### Abstract

A kinematic fitting package, KinFit, based on the Lagrange multiplier
technique has been implemented for generic hadron physics experiments. It is
particularly suitable for experiments where the interaction point is unknown,
such as experiments with extended target volumes. The KinFit package includes
vertex finding tools and fitting with kinematic constraints, such as mass
hypothesis and four-momentum conservation, as well as combinations of these
constraints. The new package is distributed as an open source software via
GitHub.

This paper presents a comprehensive description of the KinFit package and its
features, as well as a benchmark study using Monte Carlo simulations of the
$pp\rightarrow pK^{+}\Lambda\rightarrow pK^{+}p\pi^{-}$ reaction. The results
show that KinFit improves the parameter resolution and provides an excellent
basis for event selection.

###### keywords:

Kinematic Fitting, Lagrange Multiplier Technique, Hyperons, Hadrons

##  1 Introduction

Kinematic fitting is a powerful tool widely used in particle and nuclear
physics analyses, recognised for its ability to improve the resolution of
measured particle track parameters, suppress background, reconstruct
undetected particles and to determine the position of vertices. Available
information from measurements, such as momenta, angles and energy, combined
with physics constraints, such as four-momentum conservation in a production
or displaced decay vertex, or the mass of an undetected unstable particle
through its decay products, are exploited. With this information, a
mathematical minimization problem can be formulated and solved using e.g., the
Lagrange multiplier technique.

Kinematic fitting techniques have been available in particle physics since the
1960’s [1] and gained momentum during the bubble chamber days. However,
existing packages like RAVE [2] from ILD or the full decay chain fitters from
Belle II [3] or PANDA [4] are often embedded in the respective experiment
software and specialized for the detector setup. In the measurement of complex
decays of heavy particles, it is beneficial apply kinematic tree fits [3] to
complex decay chains or jets [2, 5].

In many hadron physics experiments like the HADES experiment [6] at GSI, the
interaction point is not fixed to a point but can be located within a target
volume that covers several centimeters of the path along the beamline. Hence,
the interaction vertex position needs to be determined from the measured track
parameters of the outgoing particles. In addition, weakly decaying particles
such as hyperons produce daughter particles in a secondary vertex, located a
measurable distance away from the interaction point. These are crucial to
reconstruct since many hyperons, e.g. $\Lambda$ and $\Sigma^{0}$, are neutral
and could otherwise escape detection. Figure 1 illustrates an example where a
$\Lambda$ hyperon is produced in the target volume and subsequently decays in
a secondary vertex after travelling some distance. Tracks from secondary
particles are more challenging to reconstruct and hence, the reconstructed
track parameters have larger measurement uncertainties compared to those from
primary particles. All this calls for kinematic fitting. A new package,
KinFit, has been developed to facilitate the hyperon program at HADES [7] but
is provided in an experiment-agnostic way and is therefore applicable also to
other experiments. Hadron physics experiments at J-PARC [8], JLab [9] or
COMPASS [10] and AMBER at CERN, may profit from this package.

In addition to kinematic fitting, tools are provided to construct a particle
candidate from its daughter tracks. A functionality for running a kinematic
fit for event selection on a ROOT [11] file in an automated way is included as
well.

This paper is outlined as follows: The methodology including the fitting
procedure and available constraint equations are described in Section 2. In
Section 3, the classes contained in the package are described and in Section
4, a benchmark study is provided for $\Lambda$ hyperon reconstruction to
demonstrate the performance of the fitter. In addition, a comprehensive user’s
guide is provided in Appendix A.

Figure 1: Illustration of a $\Lambda$ hyperon produced by a proton beam
hitting a proton target. The position of the interaction point (IP, marked in
blue) and the $\Lambda$ decay vertex (marked in green) are determined by the
point of closest approach of the $p_{1}$ and $K$ and $p_{2}$ and $\pi^{-}$
tracks, respectively.

##  2 Methodology

The goal of a kinematic fitting algorithm is to find an improved set of track
parameters as close as possible to the true values, such that they fulfill a
set of kinematic and/or geometric constraints. This concept can be translated
quantitatively into a $\chi^{2}$ minimization expressed as follows:

$\chi^{2}=(\vec{y}-\vec{\eta})^{T}\>V^{-1}(\vec{y})\>(\vec{y}-\vec{\eta})=\mathrm{minimum}\>,$
(1)

where $\vec{y}$ is a vector of the $N\in\mathbb{N}$ measured track parameters
provided by the tracking algorithm, $V(\vec{y})$ is the corresponding
covariance matrix and $\vec{\eta}$ are the estimated track parameters.

In addition, the kinematic and geometric constraints are implemented in the
procedure by constraint equations $\vec{f}$. These are in general represented
by a set of $K\in\mathbb{N}$ continuously differentiable functions of the
estimated parameters $\vec{\eta}$ and a set of $J\in\mathbb{N}_{0}$ unmeasured
parameters combined in $\vec{\xi}$:

$\vec{f}=\vec{f}(\vec{\eta},\,\vec{\xi})=0$

The objective is to minimize the $\chi^{2}$ from Equation 1 subject to these
constraints.

###  2.1 The Iterative Lagrange Multiplier Method

The method of Lagrange multipliers [12] is well-suited for addressing such
problems. This technique enables the transformation of the constrained
minimization into the minimization of a single Lagrange function,
$\mathcal{L}$:

$\displaystyle\mathcal{L}$
$\displaystyle=(\vec{y}-\vec{\eta})^{T}V^{-1}(\vec{y}-\vec{\eta})+2\vec{\lambda}^{T}f(\vec{\eta},\vec{\xi})\>,$
(2)

where the $K$ additional variables summarized in $\vec{\lambda}$ are
introduced, referred to as Lagrange multipliers.

Minimizing the Lagrange function (Equation 2) involves finding the derivatives
of $\mathcal{L}$ with respect to all unknowns $\vec{\eta}$, $\vec{\xi}$,
$\vec{\lambda}$. By setting these derivatives to zero and subsequently solving
for $\vec{\eta}$ and $\vec{\xi}$, the minimization can be achieved. As this
problem is generally non-linear, an iterative procedure is employed, with each
iteration yielding improved approximations for $\vec{\eta}$ and $\vec{\xi}$.
Assuming the values of all quantities of the iteration $\nu$ have already been
extracted, we want to express the quantities of the subsequent iteration
$\nu+1$ in terms of the values of iteration $\nu$. We then proceed as follows
[12]:

  1. 1.

First, the following notations are introduced

$\displaystyle\vec{r}$
$\displaystyle=\vec{f}^{\nu}+F^{\nu}_{\eta}(\vec{y}-\vec{\eta}^{\nu})$
$\displaystyle S$ $\displaystyle=F^{\nu}_{\eta}V(F^{\nu}_{\eta})^{T}\>,$

where $F_{\eta}$ is a $K\times N$ Jacobian matrix
($F_{\eta}=D_{\eta}\vec{f}$), i.e. the derivative of the constraint equations
with respect to $\vec{\eta}$.

  2. 2.

The updated, unmeasured variables, $\vec{\xi}^{\nu+1}$, are obtained by

$\displaystyle\vec{\xi}^{\nu+1}=\vec{\xi}^{\nu}-(F_{\xi}^{T}S^{-1}F_{\xi})^{-1}F_{\xi}^{T}S^{-1}\vec{r}\>,$

where $F_{\xi}$ is a $K\times J$ Jacobian matrix ($F_{\xi}=D_{\xi}\vec{f}$),
i.e. the derivative of the constraint equations with respect to $\vec{\xi}$.

  3. 3.

The updated Lagrange multipliers, $\vec{\lambda}^{\nu+1}$, are obtained from

$\displaystyle\vec{\lambda}^{\nu+1}=S^{-1}\Big{(}\vec{r}+F_{\xi}(\vec{\xi}^{\nu+1}-\vec{\xi}^{\nu})\Big{)}\>.$

  4. 4.

The updated measured parameters, $\vec{\eta}^{\nu+1}$, are calculated as

$\displaystyle\vec{\eta}^{\nu+1}=\vec{y}-VF_{\eta}^{T}\vec{\lambda}^{\nu+1}\>.$

  5. 5.

Finally, the new $\mathcal{L}$ is calculated and the results compared with the
previous iteration. To decide when the solution is sufficiently close to the
minimum, convergence criteria are defined, see Section 2.4.

In the end, the new covariance matrix, $V^{\nu+1}$, is calculated

$V^{\nu+1}\>=\>V-V[F^{T}_{\eta}S^{-1}F_{\eta}\\\
\>-\>(F^{T}_{\eta}S^{-1}F_{\xi})(F^{T}_{\xi}S^{-1}F_{\xi})^{-1}(F^{T}_{\eta}S^{-1}F_{\xi})^{T}]V\>.$
(3)

This ansatz assumes a quadratic minimum of the $\chi^{2}$ function. The
function will look different if it deviates too much from its minimum.
Currently there is no guard against this but the user should be mindful of
potential deviations and assess the validity of the results accordingly.

###  2.2 Track Parametrization

The input objects for KinFit include track parameters such as momentum and
angles, defined in polar coordinates, as well as the point of closest approach
to the beam axis. Functions for converting from Cartesian coordinates are
provided to ensure compatibility with other experiments (see the User’s Guide
in Appendix A). In addition, the covariance matrix is required as input. The
particles are assumed to move in a region free of magnetic fields, meaning
they are propagated as straight tracks. The parameters employed here to
uniquely describe a track include:

  * •

Inverse momentum $1/p$ [MeV-1c].

  * •

Polar angle $\theta$ [radians] defined from $0$ to $\pi$.

  * •

Azimuthal angle $\phi$ [radians] defined from $-\pi$ to $\pi$ relative to the
beam ($z$) axis.

  * •

$R$ [mm], the distance between the beam axis and the point of closest approach
of the track to the beam axis.

  * •

$Z$ [mm], the $z$ coordinate of the point of closest approach of the track to
the beam axis.
and

  * •

Covariance matrix, where the diagonal entries correspond to the uncertainties
in the parameters $1/p$, $\theta$, $\phi$, $R$ and $Z$.

For the purely kinematic fits, only the track parameters $1/p$, $\theta$ and
$\phi$ are required.

###  2.3 Constraints

KinFit offers a variety of constraints that can be selected for kinematic
fitting.

####  2.3.1 1C: Vertex Constraint

This purely geometrical vertex constraint minimizes the distance between two
tracks, ensuring they originate from the same vertex. A straight line in 3D is
uniquely defined by a base vector that points from the origin of a chosen
coordinate system to a coordinate on the line and the direction of the line.
The components of these vectors are:

$\begin{cases}b_{x}=R\cdot\cos(\phi+\pi/2),\\\
b_{y}=R\cdot\sin(\phi+\pi/2),\\\ b_{z}=Z,\end{cases}$ (4)

and

$\begin{cases}d_{x}=\sin(\theta)\cdot\cos(\phi),\\\
d_{y}=\sin(\theta)\cdot\sin(\phi),\\\ d_{z}=\cos(\theta).\end{cases}$ (5)

In the base vector, $\pi/2$ is added to $\phi$ to ensure the base vector is
constructed in the HADES coordinate system. A suitable constraint equation can
then be formulated as follows:

$f=(\vec{d_{1}}\times\vec{d_{2}})\cdot(\vec{b_{1}}-\vec{b_{2}})=0.$ (6)

This expression is proportional to the minimum distance between the straight
lines parameterized by the respective base and direction vectors
$\vec{b_{1}}$, $\vec{b_{2}}$, $\vec{d_{1}}$ and $\vec{d_{2}}$.

The vertex fit can be performed either separately or as part of a series of
consecutive fitting procedures. In addition, there is the possibility to
perform a mass and a geometrical vertex fit simultaneously.

####  2.3.2 4C: Four-momentum conservation in the beam-target interaction

The 4C constraint demands that the sum of the four-momenta of the final state
particles equals that of the initial beam-target system (in the following
denoted with the suffix ini). Since all parameters are treated as measured
parameters with uncertainties, there are four over-constraints. For $N$
particles, the constraint equations are

$f=\begin{cases}\sum\limits^{N}_{i=1}p_{i}\sin\theta_{i}\cos\phi_{i}-p_{\mathrm{ini},x}=0\text{
}\text{ }(p_{x}),\\\
\sum\limits^{N}_{i=1}p_{i}\sin\theta_{i}\sin\phi_{i}-p_{\mathrm{ini},y}=0\text{
}\text{ }(p_{y}),\\\
\sum\limits^{N}_{i=1}p_{i}\cos\theta_{i}-p_{\mathrm{ini},z}=0\text{ }\text{
}(p_{z}),\\\
\sum\limits^{N}_{i=1}\sqrt{p_{i}^{2}+m_{i}^{2}}-E_{\mathrm{ini}}=0\text{
}\text{ }(E).\end{cases}$ (7)

####  2.3.3 3C: Four-momentum conservation in a displaced vertex

The 3C constraint utilizes four-momentum conservation at a given decay vertex,
where a mother particle (M) decays weakly into $N$ particles. This procedure
relies on measured information about the primary vertex and the decay vertex,
and on a mass hypothesis of the mother particle as e.g. $\Lambda$ or $K_{s}$.
The angles of the mother particle, $\theta_{M}$ and $\phi_{M}$, are determined
by the direction of the vector pointing from the primary to the decay vertex
and are therefore treated as measured parameters, whereas the momentum $p_{M}$
is unmeasured and hence obtained from the fit. This results in three over-
constraints (3C). The initial value is estimated from energy conservation of
the decay products using

$p_{M}=\sqrt{\quantity(\sum_{i}\sqrt{p_{i}^{2}+m_{i}^{2}})^{2}-m^{2}_{M}}\;.$
(8)

In this expression, the subscript $M$ represents the mother particle, while
$i$ refers to the decay products. The user must provide the mass of the mother
particle as a mass hypothesis. The uncertainties in the vertex positions must
estimated in the previous step and are propagated to the uncertainties in the
parameters. Since the momentum of the mother particle $p_{M}$ is an unmeasured
quantity, its uncertainties are not known a priori but are, as the momentum
itself, obtained from the fit. The constraint equations ensure four-momentum
conservation in the decay of the mother:

$f=\begin{cases}\sum\limits^{N}_{i=1}p_{i}\sin\theta_{i}\cos\phi_{i}-p_{M}\sin\theta_{M}\cos\phi_{M}=0\text{
}\text{ }(p_{x}),\\\
\sum\limits\limits^{N}_{i=1}p_{i}\sin\theta_{i}\sin\phi_{i}-p_{M}\sin\theta_{M}\sin\varphi_{M}=0\text{
}\text{ }(p_{y}),\\\
\sum\limits^{N}_{i=1}p_{i}\cos\theta_{i}-p_{M}\cos\theta_{M}=0\text{ }\text{
}(p_{z}),\\\
\sum\limits^{N}_{i=1}\sqrt{p_{i}^{2}+m_{i}^{2}}-\sqrt{p_{M}^{2}+m_{M}^{2}}=0\text{
}\text{ }(E).\par\end{cases}$ (9)

####  2.3.4 1C: Four-momentum conservation with a missing particle

Four-momentum conservation at the interaction point enables the reconstruction
of one undetected or missing (indicated by the suffix miss) particle with a
mass hypothesis $m_{\mathrm{miss}}$. I all other inital and final state
particles are known or measured, the four-momentum conservation results in the
following constraint equations:

$f=\begin{cases}\sum\limits^{N}_{i=1}p_{i}\sin\theta_{i}\cos\phi_{i}+p_{\mathrm{miss},x}-p_{\mathrm{ini},x}=0\text{
}\text{ }(p_{x}),\\\
\sum\limits^{N}_{i=1}p_{i}\sin\theta_{i}\sin\phi_{i}+p_{\mathrm{miss},y}-p_{\mathrm{ini},y}=0\text{
}\text{ }(p_{y}),\\\
\sum\limits^{N}_{i=1}p_{i}\cos\theta_{i}+p_{\mathrm{miss},z}-p_{\mathrm{ini},z}=0\text{
}\text{ }(p_{z}),\\\
\sum\limits^{N}_{i=1}\sqrt{p_{i}^{2}+m_{i}^{2}}+\sqrt{p_{\mathrm{miss}}^{2}+m_{\mathrm{miss}}^{2}}-E_{\mathrm{ini}}\\\
\hskip 20.0pt=0\text{ }\text{ }(E).\end{cases}$ (10)

Here one has measured particles with indices $i$, a missing particle with mass
$m_{\mathrm{miss}}$ and the initial beam-target four vector
$p^{\mu}_{\mathrm{ini}}$. The mass, $m_{\mathrm{miss}}$ needs to be provided
by the user as a hypothesis as well as $p^{\mu}_{\mathrm{ini}}$.

The four-momentum of the missing particle can be extracted after the fit.
Since there are three unmeasured variables, i.e. the three-momentum of the
missing particle, only one over-constraint (1C) remains.

####  2.3.5 1C: Missing Mass Constraint

The missing mass constraint is suitable in the case when there is one
undetected final state particle whose momentum is not of interest. Instead, a
missing mass constraint can be formulated from the mass hypothesis of the
missing particle:

$f\>=\>\sqrt{\quantity(E_{\mathrm{ini}}-\sum^{N}_{i=1}{E_{i}})^{2}-\quantity(\vec{p}_{\mathrm{ini}}-\sum^{N}_{i=1}{\vec{p}_{i}})^{2}}\\\
\hskip 5.97226pt\>-m_{\mathrm{miss}}=0.$ (11)

####  2.3.6 1C: Invariant Mass Constraint

The mass constraint can be used to constrain a set of particles to originate
from a common mother particle, whose mass is known. The invariant mass of this
set of particles is then constrained to the mass of the hypothetical mother
particle:

$f=\sqrt{\quantity(\sum^{N}_{i=1}{E_{i}})^{2}-\quantity(\sum^{N}_{i=1}{\vec{p}_{i}})^{2}}-m_{\mathrm{mother}}=0\;.$
(12)

###  2.4 Convergence

The iterative fitting procedure is terminated when convergence is achieved.
There are three different convergence criteria: The difference in $\chi^{2}$
of the Lagrange function between consecutive iterations, the Euclidean norm of
the sum of all constraint equations and the Euclidean norm of the difference
of all track parameters, normalized to the initial measurement, between two
consecutive iterations. The default value for the convergence criteria is
$10^{-4}$ but can be customized. If convergence is not reached before a
specified number of iterations (default 20), the fitting procedure is exited.

###  2.5 Goodness of fit

The performance of the fit is assessed by the final $\chi^{2}_{final}$ of the
fit and the so-called pull distributions for all fitted variables. The value
of $\chi^{2}_{final}$ should be small, on the order of the number of over-
constraints, $N_{C}$. For an ensemble of events, $N_{C}$ defines the shape of
the $\chi^{2}_{final}$ distribution. There is a one-to-one relation between
$\chi^{2}_{final}$ for a given $N_{C}$ and the corresponding probability, as
defined by the probability density function [13]. Ideally, the probability
distribution should be uniform if the correct particle hypotheses have been
used in the fit and if the covariance matrices accurately describe the
measurement precision. However, if the particle hypotheses are incorrect, then
a peak towards zero probability should be discernible. If elements of the
covariance matrix are over- or underestimated, the distribution will decrease
or increase towards larger values of probability.

The pull distributions are defined as:

$pull=\frac{\eta-y}{\sqrt{\sigma^{2}(y)-\sigma^{2}(\eta)}}\;,$ (13)

for each variable, where $\eta$/$y$ are the fitted/measured variables,
respectively, and $\sigma(\eta)$/$\sigma(y)$ are the corresponding
uncertainties. The pull distribution should follow a normal distribution with
a mean value of 0 and a standard deviation of 1.

##  3 Class Descriptions

The KinFit package is written in C++ and based on ROOT [11] (version 6) and
uses CMake [14] (version 3.0 or newer) for the installation. It is available
online at:

⬇

https://github.com/KinFit/KinFit.git

Documentation about the usage of the provided tools can be found in the README
of the git repository and in the User’s Guide in Appendix A.

The properties of the particle candidates needed as input are the particle’s
track parameters, mass, and covariance matrix. Functions are provided to
transform Cartesian track parameters to the $R,Z$ track parameters.

### KFitParticle

The track parameters of particle candidates are organized in KFitParticle
objects that provide all information needed by the fitter class. It inherits
from the ROOT class TLorentzVector. For each candidate, the values for the
track parameters described in Section 2.2 have to be set. Optionally, an
arbitrary particle ID and track ID can be chosen by the user for later
reference.

### KFitDecayCandFinder

KFitDecayCandFinder calculates properties of an unmeasured candidate that
decays in a displaced vertex to be used later in the fit. The angles, $\theta$
and $\phi$ are calculated from the line segment connecting the primary to the
displaced decay vertex. The uncertainties for these angles are propagated from
the uncertainties in the vertex positions in X, Y and Z using the matrix
formalism for error propagation. The user needs to provide the uncertainties
for the vertex positions.

### KFitVertexFinder

KFitVertexFinder finds the vertex by calculating the point of closest approach
between at least two tracks by a matrix formalism. For more than two tracks,
the vertex is taken as the center of gravity, i.e. the point that is
simultaneously closest to all tracks.. This is calculated from a least square
method. The vertex finding code is based on a procedure in HYDRA
111https://subversion.gsi.de/hades/hydra2/, the HADES software.

### KinFitter

KinFitter is the main class containing the fitting functions. It can perform a
vertex, 4C, 3C, missing particle, missing mass, mass and a combined
vertex+mass fit (see Section 2.3). The constraint equations and Jacobi
matrices are implemented here. This is where the iterative fitting procedure
is carried out. The maximum number of iterations or convergence criteria can
be changed to custom values.

### KFitAnalyzer

KFitAnalyzer is a user interface class. It contains user settings and performs
the event loop with a fit of choice. The input particles need to be stored as
KFitParticles in a TClonesArray. For each event, particles are selected based
on their particle ID (PID), which is also provided to the KFitAnalyzer. A
KFitDecayBuilder object is created which takes the selected particles and the
desired constraint as input. In the end the fitted particles are retrieved
from the KFitDecayBuilder and stored in an output file together with the fit
probability. The KFitAnalyzer currently performs all fits except the 3C fit.

### KFitDecayBuilder

This class is responsible for building all possible combinations of the
particles within one event. Each combination is passed to the KinFitter which
performs the selected fit. The KFitDecayBuilder selects the combination with
the best fit probability.

##  4 Performance Study

The performance of KinFit is benchmarked using Monte Carlo (MC) simulations.
These represent an ideal scenario where a finite detector resolution is added
by smearing, and thus known, whereas in real data, it might be difficult to
estimate the covariance matrix exactly. This approach allows for a focused
investigation of the quality of the KinFit tools. The reaction investigated is
$pp\rightarrow pK^{+}\Lambda$, $\Lambda\rightarrow p\pi^{-}$, as illustrated
in Figure 1. This reaction provides an opportunity to examine all tools
contained in the package. All fits performed in this section use the default
settings of KinFitter, described in Section 2.4.

The event generator Pluto [15] was used to generate $100\,000$ events of the
$pp\rightarrow pK^{+}\Lambda$ $\Lambda\rightarrow p\pi^{-}$ reaction at a beam
kinetic energy of T=$4.5\,\mathrm{GeV}$. In our simulations, the produced
particles and their decay products are distributed isotropically across the
available phase space. The primary vertex is generated at the origin, i.e.
(0,0,0). A full 4$\pi$ acceptance is assumed and no material effects are
included. The $\Lambda$ hyperons decay at a displaced decay vertex according
to the mean $\Lambda$ hyperon life time $c_{\tau}(\Lambda)=7.98\,\mathrm{cm}$
[13]. The track parameters of the final state particles are smeared according
to a Gaussian distribution, simulating uncertainties from e.g., detector
resolution in a controlled way. The uncertainty for each track parameter is
listed in Table 1. Since the results presented here do not depend on a
specific experiment, the uncertainties are set to be constant, except the
momentum uncertainty, which is momentum dependent.

Track parameter | $1/p$ | $\theta$ | $\phi$ | $R$ | $Z$
---|---|---|---|---|---
$\sigma$ | $0.025\cdot 1/p$ | $0.0009\,\mathrm{rad}$ | $0.0009\,\mathrm{rad}$ | $0.5\,\mathrm{mm}$ | $1\,\mathrm{mm}$
Table 1: Gaussian standard deviations of the track parameters in the MC sample
after smearing.

###  4.1 Mass Fit using KFitAnalyzer

In each event, there are four measured particles: a proton ($p_{1}$) and a
kaon from the interaction point and a pion and a proton ($p_{2}$) from the
$\Lambda$ hyperon decay. This means there are two possibilities to combine a
pion with a proton, either $p_{1}\pi^{-}$ and $p_{2}\pi^{-}$, of which the
latter is the correct $\Lambda$ decay candidate. The KFitAnalyzer is used to
find the pion-proton combination that comes from the $\Lambda$ hyperon decay
through the application of a mass fit.
The pion-proton pair giving the largest mass fit probability (see Fig. 2) are
selected as $\Lambda$ daughters. In $99\,\%$ of the events, the correct proton
$p_{2}$ is chosen. The probability distribution is uniform, as expected.

Figure 2: Probability distribution of a $\Lambda$ hyperon mass fit in the
reaction $pp\rightarrow pK^{+}\Lambda$ $\Lambda\rightarrow p\pi^{-}$ using
KFitAnalyzer (top) and the proton selected to originate from the $\Lambda$
hyperon decay (bottom)

###  4.2 Reconstruction of $\Lambda$ Hyperon from Vertex Positions and a 3C
Fit

The $\Lambda$ hyperon reconstruction procedure includes three steps:

  1. 1.

Vertex finding
As a first step the position of the primary vertex, where the $\Lambda$
hyperon was produced is determined from the point of closest approach between
the other tracks coming from the production vertex. In this case, these are
the proton ($p_{1}$) and kaon tracks. This task is performed by
KFitVertexFinder. In a similar way, the $\Lambda$ hyperon decay vertex is
obtained from the point of closest approach between the proton ($p_{2}$) and
pion tracks originating from the $\Lambda$ hyperon decay. The vertex positions
are shown in Fig. 3. All primary particles are produced at the origin which
means that the distribution in the top panel reflects the finite resolution of
the vertex estimation. The dominating effect in the $\Lambda$ decay vertex
distribution is the distance travelled by $\Lambda$ before decaying.

Figure 3: Reconstructed production vertex from proton and $K^{+}$ (top) and
location of reconstructed decay vertices of the $\Lambda$ hyperon from the
secondary proton and $\pi^{-}$ (bottom) as reconstructed by the
KFitVertexFinder. The positions are projected onto the x-z plane.

  2. 2.

Reconstruction of the neutral candidate
The direction of the $\Lambda$ candidate is given by the vector pointing from
the primary to the decay vertex. The magnitude of the $\Lambda$ hyperon
momentum is estimated by the of the sum of the three-momenta of its daughter
particles, see Equation 8. The KFitDecayCandFinder carries out this task and
calculates the uncertainties of the track parameters of the $\Lambda$
candidate. The latter requires information on the vertex uncertainties which
were estimated by Gaussian fits to the vertex resolution in $x$, $y$ and $z$
for each vertex.

  3. 3.

3C fit
A kinematic fit ensuring four-momentum conservation at the $\Lambda$ hyperon
decay vertex is the final step and is executed by KinFitter. This improves the
resolution of the $\Lambda$ hyperon track parameters significantly (Figure 4)
and can be used to reject false combinations of protons in the vertices.

Figure 4: Resolution of the $\Lambda$ hyperon track parameters before and
after the 3C fit.

Only the correct combination of particles was used in this example, i.e.
correct assignment of the protons to the vertices. The pull distributions for
the pion track are shown as an example in Fig. 5. These nearly follow a normal
distribution, as expected. Fig. 6 shows the probability distribution of the 3C
fit. The probability deviates slightly from a uniform distribution and there
are some events for which the fit converged, but gives a low probability. This
effect and the slightly larger deviation of the pull distributions from a
normal distribution are only seen in the 3C fit and originate from the
estimation of the vertex resolutions, which are not exactly described by a
Gaussian.

Figure 5: Pull distributions for the $\pi^{-}$ track parameters after the 3C
fit with respective mean and standard deviation. Figure 6: Probability
distribution from the 3C fit.

###  4.3 4C Fit of the Reaction $pp\rightarrow pK^{+}\Lambda$

The 4C fit is used to constrain the four-momenta of all final state particles
to that of the initial beam-target system. In this example, the final state
particles are a proton and a kaon from the primary vertex and a proton and a
pion from the $\Lambda$ hyperon decay vertex. Fig. 7 shows the momentum
resolution for the kaon and the pion before and after the fit. The improvement
is more substantial for the kaon momentum resolution compared to the other
particles. This is due to the larger uncertainty of the kaon associated with
its larger total momentum. The probability distribution is uniform and the
pull distributions follow a normal distribution, as shown in Figure 8.

Figure 7: Momentum resolution for the $K^{+}$ (top) and $\pi^{-}$ (bottom)
before and after the 4C fit.

Figure 8: 4C fit. Probability distribution (top) and example pull distribution
of the $\pi^{-}$ momentum with mean and standard deviation (bottom).

###  4.4 Missing Particle Fit of $K^{+}$ in the Reaction $pp\rightarrow
pK^{+}\Lambda$.

In this scenario, it is assumed that the kaon is not detected. The missing
particle fit is employed to constrain the four-momenta of the detected
particles along with the undetected particle to match the beam-target system,
as well as to estimate the momentum of the missing particle. Figure 9 presents
the probability distribution of the fit and compares the kaon momentum
resolution from the initial guess and after the fit. The initial guess for the
kaon momentum is calculated from three-momentum conservation in the
interaction point. The momentum resolution before the fit is worse than in the
example in Section 4.3, where the momentum was measured directly. A noticeable
improvement in the momentum resolution can be observed after the fit. However,
the resolution is still worse than in the case of the 4C fit, as the missing
particle fit has fewer over-constraints. The resolutions of the other track
parameters and the pull distributions appear quite similar to those of the 4C
fit.

Figure 9: Missing particle fit. Probability distribution (top) and resolution
of the estimated and fitted $K^{+}$ momentum (bottom).

###  4.5 Vertex Fit

The vertex fit aims to constrain the track parameters of final state particles
that originate from the same single point in space. To perform the vertex fit,
at least two outgoing particles from the same vertex must be measured. In this
example, these particles are the kaon and proton ($p_{1}$) produced at the
beam-target interaction point. Fig. 10 illustrates the resolution of the
estimated R-parameter, defined in Section 2.2, for the proton. The probability
distribution for this vertex fit, along with an example pull distribution of
the R-parameter of $p_{1}$, is depicted in Fig. 11. As anticipated, the
probability is uniformly distributed across its range, while the pull follows
a normal distribution. The KinFitAnalyzer was used to choose the proton track
that has the largest probability to come from the same vertex as the kaon
track. The correct proton ($p_{1}$) is chosen in $85\,\%$ of the events, which
is lower than the fraction of correctly identified combinations by the mass
fit presented in Section 4.1. This shows that the mass constraint is more
powerful than the vertex fit.

Figure 10: Vertex fit. Resolution of the R-parameter of the primary proton
($p_{1}$) before and after the fit.

Figure 11: Vertex fit. Probability distribution (top) and example pull
distribution of the R-parameter for the primary proton ($p1$) (bottom).

###  4.6 Runtime Performance

The runtime performance of various components of KinFit was evaluated on a
laptop with the following hardware specifications:

  * •

Processor: Intel Core i7-1185G7, 3.00 GHz

  * •

Memory: 32 GB RAM

The results are displayed in Table 2. The average runtime per 1000 events was
calculated for those in which the fit converged. A low overhead of
approximately a microsecond was observed.

Stage | Fitting | Initialization | Vertex Finding | Mother Candidate Finding
---|---|---|---|---
Time / event [$\mu$s] | 10 | 1 | 1 | 1

Table 2: Runtime per event for different parts of KinFit.

##  5 Conclusions

The KinFit package is a versatile kinematic fitting package that has been
developed to provide tools for reconstruction of vertices, momenta and masses
of particles in a generic hadron physics experiment. In particular, it is
capable of fitting an unknown production vertex from an extended beam-target
interaction volume, and of combining kinematic and geometric constraints. The
tools provided by KinFit have been tested with toy MC simulations,
demonstrating outstanding performance in terms of improved track parameter
resolutions and the ability to select correct particle hypotheses and
-combinations.

While the focus has been on providing suitable tools for hadronic interactions
and specifically hyperon decays, the package can be applied for other types of
reactions as well. Though KinFit uses the same track parametrization as HADES,
but it is generally considered experiment-independent.

The minimal overhead introduced by adding KinFit to an analysis suggests that
the package is suitable for running on a local machine.

##  6 Outlook

Although the majority of particle physics analyses are based on C++, the
demand for similar Python tools is currently growing due to the rapid
expansion of the Python ecosystem. Therefore, the KinFit package is planned to
be integrated into the SciKit-HEP project [16], which is a community-driven
Python ecosystem for data analysis.

## Author Contributions

This work is conceptualized by Waleed Esmail and Karin Schönning. The
development and implementation of the C++ KinFit package was done by Esmail,
Jana Rieger and Jenny Taylor. The first paper draft preparation, was done by
Esmail, Rieger and Taylor and the drafting was coordinated by Rieger. The
project was supervised by Schönning. Benchmark studies and bug fixes were
performed by Rieger, Taylor and Malin Bohman, while all authors contributed to
the review of the paper’s intellectual content and the editing.

## Acknowledgements

We would like to thank Tord Johansson for useful discussions on kinematic
fitting of $\Lambda$ hyperons. Furthermore, we are grateful to Anar Rustamonov
for the suggestion to make KinFit available to a larger community as an
independent package. Finally, we would like to express our appreciation of the
fruitful discussions and the encouragements we received from James Ritman and
Piotr Salabura.

The publication is funded by the Deutsche Forschungsgemeinschaft (DFG, German
Research Foundation) – 491382106 , and by the Open Access Publishing Fund of
GSI Helmholtzzentrum fuer Schwerionenforschung. This project has received
funding from the Knut and Alice Wallenberg Foundation, Contract No. 2016.0157
and 2021.0299 (Sweden), and the Swedish Research Council, Contract No.
2019-04594 (Sweden).

## References

  *   * Bock [1960] Bock, R.K.: Application of a Generalized Method of Least Squares for Kinematical Analysis of Tracks in Bubble Chambers. CERN Yellow Reports: Monographs. CERN, Geneva (1960). https://doi.org/10.5170/CERN-1960-030
  * Moser et al. [2009] Moser, F., Waltenberger, W., Regler, M., Mitaroff, W.: Implementation and application of kinematic vertex fitting in the software environment of ILD. arXiv (2009). https://doi.org/10.48550/ARXIV.0901.4020
  * Krohn et al. [2020] Krohn, J.-F., et al.: Global decay chain vertex fitting at belle II. Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment 976, 164269 (2020) https://doi.org/10.1016/j.nima.2020.164269
  * PANDA Collaboration [2021] PANDA Collaboration: Panda phase one. Eur. Phys. J. A 184(57) (2021) https://doi.org/10.1140/epja/s10050-021-00475-y
  * Radkhorrami and List [2021] Radkhorrami, Y., List, J.: Kinematic Fitting for ParticleFlow Detectors at Future Higgs Factories. arXiv (2021). https://doi.org/10.48550/ARXIV.2110.13731
  * HADES Collaboration [2009] HADES Collaboration: The high-acceptance dielectron spectrometer HADES. The European Physical Journal A 41(2), 243–277 (2009) https://doi.org/10.1140/epja/i2009-10807-5
  * The HADES Collaboration and PANDA Collaboration [2021] The HADES Collaboration and PANDA Collaboration: Production and electromagnetic decay of hyperons: a feasibility study with HADES as a phase-0 experiment at FAIR. Eur. Phys. J. A 57(4), 138 (2021) https://doi.org/10.1140/epja/s10050-021-00388-w arXiv:2010.06961 [nucl-ex]
  * Ohnishi et al. [2020] Ohnishi, H., Sakuma, F., Takahashi, T.: Hadron physics at j-parc. Progress in Particle and Nuclear Physics 113, 103773 (2020) https://doi.org/10.1016/j.ppnp.2020.103773
  * [9] JLab Experiments. https://www.jlab.org/physics/experiments. Accessed: 2023-05-05
  * Abbon et al. [2015] Abbon, P., et al.: The COMPASS setup for physics with hadron beams. Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment 779, 69–115 (2015) https://doi.org/10.1016/j.nima.2015.01.035
  * Brun and Rademakers [1997] Brun, R., Rademakers, F.: ROOT: An object oriented data analysis framework. Nucl. Instrum. Meth. A 389, 81–86 (1997) https://doi.org/10.1016/S0168-9002(97)00048-X
  * Frodesen and Skjeggestad [1979] Frodesen, A.G., Skjeggestad, O.: Probability and Statistics in Particle Physics. Universitetsforlaget, Bergen, Norway (1979)
  * Zyla et al. [2020] Zyla, P.A., et al.: Review of Particle Physics. PTEP 2020(8), 083–01 (2020) https://doi.org/10.1093/ptep/ptaa104
  * [14] CMake. https://cmake.org/. Accessed: 2022-12-15
  * Fröhlich et al. [2010] Fröhlich, I., et al.: Design of the pluto event generator. Journal of Physics: Conference Series 219(3), 032039 (2010) https://doi.org/10.1088/1742-6596/219/3/032039
  * Rodrigues et al. [2020] Rodrigues, E., et al.: The Scikit HEP Project – overview and prospects. EPJ Web Conf. 245, 06028 (2020) https://doi.org/10.1051/epjconf/202024506028 arXiv:2007.03577 [physics.comp-ph]

##  Appendix A User’s Guide

The source code is available online in the KinFit git repository. It can be
downloaded via

⬇

git clone https://github.com/KinFit/KinFit.git

and installed using CMake. A ROOT 6 installation is required.

There are two options how the tools provided by this package may be applied by
the user. One of them is to use the provided user functions through
KFitAnalyzer. These make it straight-forward to apply the fitter using basic
constraints. However since the type of analysis where kinematic fitting can
yield large improvements often requires a careful and customized event
selection, the provided classes can also be applied directly by the user,
which makes a fine tuning of the fitting process possible but requires a
deeper understanding of the procedure.

###  A.1 Automated application of a kinematic fit

In order to use KinFit in an automated way, the user needs to provide a root
file with a TTree called ”data” that contains a TClonesArray ”KFitParticle” of
KFitParticles as input. The macro analysis_user.C illustrates how to set up
the KFitAnalyzer.

⬇

#include "KFitAnalyzer.h"

Int_t analysis_user(TString infile, TString outfile, Int_t evts){

KFitAnalyzer RootAnalyzer(infile, outfile, evts);

std::vector<int> pids;

pids.push_back(14); pids.push_back(9);

Double_t mass = 1.11568;

RootAnalyzer.doFitterTask("Mass", pids, mass);

return 0;

}

A KFitAnalyzer object is created. It takes the input file(s), the output file
name and the number of events to be processed as input. The only command that
needs to be executed to start the analysis is doFitterTask. Here the type of
fit, the PIDs of the particles to be fit and additional information that is
required for the specific fit is requested. In this example, an invariant mass
fit of a proton and a pion is performed, assuming that they originate from a
$\Lambda$ decay. After the analysis procedure the fitted tracks of the best
particle combination and the fit probability are written to the output file.
This requires that the kinematic fit converges and that that combination yield
the highest fit probability of the event.

Other available fits are

⬇

doFitterTask("4C", pids, -1, ppSystem);

doFitterTask("MassVtx", pids, mLambda); doFitterTask("MissingMass", pids,
mass, ppSystem);

doFitterTask("MissingParticle", pids, mass, ppSystem);

doFitterTask("Vertex", pids, -1);

###  A.2 Using the individual classes

In an event-based analysis, the kinematic fitting tools are typically applied
as one of the first analysis steps inside the ”event loop”, i.e. while
iterating through all events, selecting the most suitable particle candidates
from each.

Creating a KFitParticle
After a set of particle candidates are selected whose track parameters are to
be fit, a KFitParticle object is created for each candidate. In addition, the
covariance matrix has to be assigned to each candidate. This is done either by
calling the constructor KFitParticle(TLorentzVector cand, double R, double Z)
or KFitParticle(TLorentzVector cand, double X, double Y, double Z) and setting
the covariance matrix using the setCovariance() function directly, or by using
a suitable FillData function, useful when the covariance was estimated and is
not known for each particle individually. An example of how a FillData
function can look like to set the attributes of the KFitParticle is shown
below.

⬇

void FillData(TLorentzVector *cand, KFitParticle *outcand, double R, double Z,
double arr[], double mass, int trackID, int pid)

/** Cov(5,5)Covariance matrix of the KFitParticle

* Diagonal entries correspond to the covariances

* in the parameters in the following order

*

* \-----------------------------

* | 1/p |

* | theta |

* | phi |

* | R |

* | Z |

* \-----------------------------

*

* Off diagonal elements corresponds to the

* correlations between the parameters

* arr[] can be adjusted acordingly to set the

* off diagonal elements

*/

TMatrixD cov(5, 5);

cov(0, 0) = std::pow(arr[0], 2);

cov(1, 1) = std::pow(arr[1], 2);

cov(2, 2) = std::pow(arr[2], 2);

cov(3, 3) = std::pow(arr[3], 2);

cov(4, 4) = std::pow(arr[4], 2);

outcand->SetXYZM(

cand->P() * std::sin(cand->Theta()) * std::cos(cand->Phi()),

cand->P() * std::sin(cand->Theta()) * std::sin(cand->Phi()),

cand->P() * std::cos(cand->Theta()), mass);

outcand->setThetaRad(cand->Theta());

outcand->setPhiRad(cand->Phi());

// outcand->setThetaDeg(cand->Theta()); // Depending on the unit of the angles
from cand the functions setThetaRad(double val) and setThetaDeg(double val)
can be used iterchangably

// outcand->setPhiDeg(cand->Phi()); // Depending on the unit of the angles
from cand the functions setPhiRad(double val) and setPhiDeg(double val) can be
used iterchangably

outcand->setR(R);

outcand->setZ(Z);

outcand->setCovariance(cov);

//optional

outcand->setTrackId(trackID);

outcand->setPid(pid);

In this example, R and Z are given as an input by the user. Alternatively the
creation point of the particle candidate can be given in cartesian coordinates
which are internally converted to $R$ and $Z$. SetXYZM() sets the attributes
of the TLorentzVector that the KFitParticle inherits from whereas
setThetaRad(), setPhiRad(), setR(), setZ() and setCovariance() operate
directly on the KFitParticle object. The covariance matrix needs to be
provided by the user as well as a mass hypothesis. A PID and track-ID can be
set optionally.

Handling the KinFitter
To use fitting functions, the KFitParticle objects that will be used need to
be placed in a vector of the type std::vector from the C++ STD library. This
vector needs to be passed to the constructor of KinFitter. How this is done is
illustrated in the examples below. For most fit options, an arbitrary number
of particles can be added. However, for the vertex fit, exactly two particles
need to be added. Following this, the user must choose one of the constraints,
either; the 3C, the 4C, mass, missing mass, the vertex or the missing
particle. Adjusting the number of iterations as well as the convergence
criteria is possible but optional; if not set, the default criteria are
applied. Then the fit() function is called. This performs the actual fitting
and returns true if the fit converged. Some information can be obtained
independently of the individual fit function that is used, e.g. getChi2()
returns the $\chi^{2}$ of the fit and getProb() returns the corresponding
probability. getPull(int val) returns the pull of variable $v$ for a daughter
particle $p$, $val=5\cdot d+v$. The function isConverged() returns a boolean;
true if the fit has converged and false otherwise, this can be used for event
or sample selections. The initialization of KinFitter settings that can be
used before calling fit() are summarized below.

⬇

KinFitter(const std::vector<KFitParticle> &cands);

setNumberOfIterations(int val); // default 20

setConvergenceCriteria(double val, double val2, double val3);

// default 1e-4

setVerbosity(int val); //between 0 and 2, default 0

After the fitting procedure called by fit(), the fit result can be assessed by
the following functions. Updated track parameters are obtained by the
KFitParticle objects returned by the first function

⬇

KFitParticle getDaughter(int val); // Fitted daughter particle #val

double getChi2(); // Returns chi2

double getProb(); // Returns fit probability

double getPull(int val); // Returns pull of variable val

//* The pulls are returned in the following order:

* val = i + 0 : 1/ p of particle i

* val = i + 1 : theta of particle i

* val = i + 2 : phi of particle i

* val = i + 3 : R of particle i

* val = i + 4 : Z of particle i

* Where the particles appear in the same order they

* were entered into the fit

*/

bool isConverged(); // Returns true if fit converged

int getIteration(); // Returns number of iterations until convergence or
maximal number of iterations

An example of how to perform the vertex fit is performed is shown below.

⬇

#include "KinFitter.h"

KinFitter fitter(cand_vector);

fitter.addVertexConstraint(); // Choose vertex constraint

fitter.fit();

To perform the missing particle fit, the mass of the missing particle and a
TLorentzVector corresponding to the initial beam-target system must be given
as input to the fitter in addition to a mass hypothesis. An example of how the
TLorentzVector can be constructed is in the following way:

TLorentzVector ppSystem(p1,p2,p3,E);

where the first three entries correspond to the momentum in each Cartesian
direction and the last entry is the total energy of the system. After the fit
has been performed, the missing daughter can be retrieved as a TLorentzVector

⬇

#include "TLorentzVector.h"

#include "KinFitter.h"

Double_t mass;

TLorentzVector ppSystem(p1,p2,p3,E);

KinFitter fitter(cand_vector);

void addMissingParticleConstraint(ppSystem, mass); // Choose momentum
constraint for missing particle

fitter.fit();

TLorentzVector getMissingDaughter(); // Retrieve the missing daughter after
the fit

To perform the 4C fit, in addition to the vector of final state particles, a
TLorentzVector corresponding to the initial beam-target system must be passed
to the fitter.

⬇

#include "TLorentzVector.h"

#include "KinFitter.h"

TLorentzVector ppSystem(p1,p2,p3,E);

KinFitter fitter(cand_vector);

fitter.add4Constraint(ppSystem); // Choose 4C fit

fitter.fit();

To perform a mass fit, in addition to the vector of final state particles, the
mass of the particle must be passed to the fitter.

⬇

#include "KinFitter.h"

Double_t mass;

KinFitter fitter(cand_vector);

fitter.addMassConstraint(mass); // Choose mass fit

fitter.fit();

To perform a missing mass fit, in addition to the vector of final state
particles, the mass of the missing particle and a TLorentzVector corresponding
to the initial beam-target system must be passed to the fitter.

⬇

#include "TLorentzVector.h"

#include "KinFitter.h"

Double_t mass;

TLorentzVector ppSystem(p1,p2,p3,E);

KinFitter fitter(cand_vector);

fitter.addMissingMassConstraint(ppSystem, mass); // Choose missing mass fit

fitter.fit();

Since a mother particle needs to be constructed in order to perform the 3C
fit, this is a bit more involved and an example of how to perform this fit and
run the 3C fit is given below

⬇

#include "KinFitter.h"

#include "KFitVertexFinder.h"

#include "KFitDecayCandFinder.h"

#include "TLorentzVector.h"

#include <vector>

#include <algorithm>

#include <map>

#include <iostream>

#include <iomanip>

#include <math.h>

// Create a vector for particle pairs from the interaction point and decay
vertex, respectively

std::vector<KFitParticle> cands1, cands2;

cands1.push_back(proton1_fit);

cands1.push_back(kaon_fit);

cands2.push_back(proton2_fit);

cands2.push_back(pion_fit);

// Initialize a KFitVertexFinder for each vertex

KFitVertexFinder vtx1finder(cands1);

KFitVertexFinder vtx2finder(cands2);

// Find vertex and retrieve it

TVector3 vtx1 = vtx1finder.getVertex();

TVector3 vtx2 = vtx2finder.getVertex();

// Vertex resolutions

Double_t vtx1_xres, vtx1_yres, vtx1_zres, vtx2_xres, vtx2_yres, vtx2_zres;

// Find the decaying candidate

KFitDecayCandFinder lambdafinder(cands2, 1.115683, vtx1, vtx2, vtx1_xres,
vtx1_yres, vtx1_zres, vtx2_xres, vtx2_yres, vtx2_zres);

KFitParticle lambda_cand =

lambdafinder.getDecayCand();

// Do 3C fit in decay vertex

KinFitter fitter(cands2);

fitter.add3Constraint(lambda_cand);

fitter.fit();

// Get fit result

KFitParticle fcand1 = fitter.getDaughter(0);

// proton

KFitParticle fcand2 = fitter.getDaughter(1); // pion

KFitParticle lambda_fit = fitter.getMother();

// lambda

TMatrixD test = lambda_fit.getCovariance();

// lambda covariance

##  Appendix B QA Plots

###  B.1 Pull Distributions of 4C Fit

Figure 12: Pull distributions for the track parameters of all particles after
the 4C fit with respective mean and standard deviation.

###  B.2 Pull Distributions of Missing Particle Fit

Figure 13: Pull distributions for the track parameters of all particles after
the missing $K^{+}$ fit with respective mean and standard deviation.

###  B.3 Pull Distributions of 3C Fit in $\Lambda$ Decay Vertex

Figure 14: Pull distributions for the track parameters of all particles after
the 3C fit with respective mean and standard deviation.

###  B.4 Pull Distributions of the Vertex Fit in the Interaction Point

Figure 15: Pull distributions for the track parameters of the primary proton
tracks after the vertex fit with respective mean and standard deviation.


# Counting Schur Rings over Cyclic Groups of Semiprime Order

Joseph Keller, Andrew Misseldine111Andrew Misseldine, Southern Utah
University<EMAIL_ADDRESS>(telephone) 435-865-8228, (fax)
435-865-8666, Max Sullivan

###### Abstract

In this paper, we continue the enumeration of Schur rings over cyclic groups.
Cyclic groups of semiprime order $pq$, where $p$ and $q$ are distinct primes,
are considered. Additionally, cyclic groups of order $4p$ are considered.

Keywords: Schur ring, cyclic group, association scheme, lattice of subgroups

AMS Classification: 20c05, 05c25, 05e30, 05a15 05e16 20k27

##  1 Introduction

Let $G$ be a finite group, and let $\mathbb{Q}[G]$ denote the rational group
algebra. Let $\mathcal{L}(G)$ denote the lattice of subgroups of $G$. For any
subset $C\subseteq G$, we may associate to it an element of the group algebra,
namely $\sum_{g\in C}g\in\mathbb{Q}[G]$. Such an element is called a _simple
quantity_. When the context is clear, we will likewise denote the simple
quantity associated to $C$ as $C$ itself. Define $C^{*}:=\\{x^{-1}\mid x\in
C\\}$ for all $C\subseteq G$. Let $\\{C_{1},C_{2},\ldots,C_{r}\\}$ be a
partition of $G$, and let $\mathfrak{S}$ be the subspace of $\mathbb{Q}[G]$
spanned by $\\{C_{1},C_{2},\ldots,C_{r}\\}$. We say that $\mathfrak{S}$ is a
_Schur ring_ over $G$ if

  1. 1.

$C_{1}=\\{1\\}$;

  2. 2.

For each $i$, there is a $j$ such that $C_{i}^{*}=C_{j}$;

  3. 3.

For each $i$ and $j$,
$C_{i}C_{j}=\displaystyle\sum_{k=1}^{r}\lambda_{ijk}C_{k}$, for
$\lambda_{ijk}\in\mathbb{N}$.

The sets $C_{1},C_{2},\ldots,C_{r}$ are called the _$\mathfrak{S}$ -classes_
(or _primitive sets_ of $\mathfrak{S}$). Note that a Schur ring is uniquely
determined by its associated partition of $G$. We will denote this partition
as $\mathcal{D}(\mathfrak{S})$.

Schur rings over cyclic groups have been of great interest for the last few
decades because of their connection to algebraic graph theory and association
schemes (see [7]). In [4], the second author develops a general technique of
enumerating cyclic Schur rings that we will employ in the case of $n=pq,4p$,
where $p$ and $q$ are distinct primes. This paper is then a continuation of
work begun by the second author in [4, 5, 6].

Let $\mathcal{Z}_{n}=\langle z\rangle$ denote the cyclic group of order $n$.
Let $\Omega(n)$ denote the number of Schur rings over $\mathcal{Z}_{n}$. In
this paper, we provide formulas for $\Omega(pq)$ and $\Omega(4p)$. We present
first the semiprime case $pq$:

######  Theorem 1.1.

Let $p$ and $q$ be distinct primes such that
$p=\prod_{i=1}^{n}r_{i}^{k_{i}}+1$ and $q=\prod_{i=1}^{n}r_{i}^{\ell_{i}}+1$,
where $\\{r_{1},r_{2},\ldots,r_{n}\\}$ is a list of distinct primes. Then

$\Omega(pq)=\prod_{i=1}^{n}\sum_{j=0}^{\min(k_{i},\ell_{i})}\phi(r_{i}^{j})(k_{i}-j+1)(\ell_{i}-j+1)+2\prod_{i=1}^{n}(k_{i}+1)(\ell_{i}+1)+1,$

where $\phi$ denotes Euler’s totient function.

Be aware that in the above decompositions of $p$ and $q$, the exponents
$k_{i}$ and $\ell_{i}$ may possibly be zero, allowing for a common family of
primes $\\{r_{1},r_{2},\ldots,r_{n}\\}$ between $p$ and $q$.

We list next some useful simplifications of Theorem 1.1 when special
conditions are placed on the primes $p$ and $q$. The proofs of the following
corollaries are immediate from Theorem 1.1.

######  Corollary 1.2.

Let $p$ and $q$ be distinct primes such that $p=2^{k}a+1$ and $q=2^{\ell}b+1$,
where $a$ and $b$ are both odd integers and $\gcd(a,b)=1$. Let $x$ and $y$ be
the number of divisors of $p-1$ and $q-1$, respectively. Then

$\Omega(pq)=\left[3(k+1)(\ell+1)+\sum_{j=1}^{\min(k,\ell)}2^{i}(k-j+1)(\ell-j+1)\right]\left(\dfrac{xy}{(k+1)(\ell+1)}\right)+1.$

######  Corollary 1.3.

Let $p\neq 2$ be a prime, and let $x$ be the number of divisors of $p-1$. Then

$\Omega(2p)=3x+1.$

Corollary 1.2 is particular useful when $p$ is a _Fermat prime_ , that is,
$p=2^{k}+1$. There are only five know Fermat primes: $3$, $5$, $17$, $257$,
and $65537$. It is widely conjectured that these are the only Fermat primes.
We illustrate the simplification of Corollary 1.2 for the Fermat primes $3$
and $5$.

######  Corollary 1.4.

Let $p\neq 3$ be a prime such that $p=2^{k}a+1$ where $a$ is odd, and let $x$
be the number of divisors of $p-1$. Then

$\Omega(3p)=\left(\dfrac{7k+6}{k+1}\right)x+1$

When $p\equiv 3\pmod{4}$, $\Omega(3p)=\dfrac{13}{2}x+1$.

######  Corollary 1.5.

Let $p\neq 5$ be a prime such that $p=2^{k}a+1$ where $a$ is odd, and let $x$
be the number of divisors of $p-1$. Then

$\Omega(5p)=\left(\dfrac{13k+7}{k+1}\right)x+1$

When $p\equiv 3\pmod{4}$, $\Omega(5p)=10x+1$.

We mention that Corollary 1.2 is also applicable when $p$ is a _safe prime_ ,
that is, $p=2r+1$, where $r$ is itself a prime222In this case, $r$ is
necessarily a Sophie Germain prime.. It is widely conjectured that there are
infinitely many safe primes, the first few being:333We have intentionally
omitted $5$ from the list of safe primes as it is the only safe prime which is
Fermat. As a consequence, it is the only safe prime $p$ for which the number
of divisors of $p-1$ is 3 instead of 4. $7$, $11$, $23$, $47$, $59$, $83$, and
$107$. Note that by Corollary 1.2, if $p$ and $q$ are both safe primes, then
$\Omega(pq)=53$. Likewise, if $p$ is a safe prime, then $\Omega(2p)=13$,
$\Omega(3p)=27$ and $\Omega(5p)=41$, by Corollaries 1.3, 1.4, and 1.5,
respectively.

Using Corollaries 1.3, 1.4, and 1.5 and the above discussion of safe primes,
one can easily compute the number of Schur rings over $\mathcal{Z}_{pq}$ for
all semiprimes under 100. These are listed in Table 1.1. The one exception
here is $n=91=7\cdot 13=(2\cdot 3+1)(2^{2}\cdot 3+1)$. In this case,
$\Omega(91)$ can be computed directly using Theorem 1.1.

Table 1.1: Number of Schur Rings over $Z_{pq}$ _n_ | $\boldsymbol{\Omega(n)}$ | _n_ | $\boldsymbol{\Omega(n)}$ | _n_ | $\boldsymbol{\Omega(n)}$ | _n_ | $\boldsymbol{\Omega(n)}$ | _n_ | $\boldsymbol{\Omega(n)}$
---|---|---|---|---|---|---|---|---|---
6 | 7 | 26 | 19 | 46 | 13 | 65 | 67 | 86 | 25
10 | 10 | 33 | 27 | 51 | 35 | 69 | 27 | 87 | 41
14 | 13 | 34 | 16 | 55 | 41 | 74 | 28 | 91 | 97
15 | 21 | 35 | 41 | 57 | 40 | 77 | 53 | 93 | 53
21 | 27 | 38 | 19 | 58 | 19 | 82 | 25 | 94 | 13
22 | 13 | 39 | 41 | 62 | 25 | 85 | 60 | 95 | 61

We next present the counting formula for $4p$:

######  Theorem 1.6.

Let $p$ be an odd prime such that $p=2^{k}a+1$, where a is an odd integer and
$x$ the number of divisors of $p-1$. Then

$\Omega(4p)=\frac{15k+14}{k+1}x+3.$

Table 1.2: Number of Schur Rings over $Z_{4p}$ _n_ | $\boldsymbol{\Omega(n)}$ | _n_ | $\boldsymbol{\Omega(n)}$ | _n_ | $\boldsymbol{\Omega(n)}$ | _n_ | $\boldsymbol{\Omega(n)}$
---|---|---|---|---|---|---|---
12 | 32 | 28 | 61 | 52 | 91 | 76 | 90
20 | 47 | 44 | 61 | 68 | 77 | 92 | 61

In the special case that $p=2^{k}+1$ is a Fermat prime, the above formula
simplifies to $\Omega(4p)=15k+17$. For safe primes, we always have
$\Omega(4p)=61$. In Table 1.2 we list all integers of the form $4p$ less than
100.

The proof of Theorem 1.1 and Theorem 1.6 (which proofs can be found in Section
3 and Section 4, respectively) can be summarized as following. By the
Fundamental Theorem of Schur Rings over Cyclic Groups (due to Leung and Man
[2, 3]), all Schur rings over cyclic groups belong to one of four families,
which we call the _traditional Schur rings_ : namely, the trivial Schur ring,
automorphic Schur rings, direct products of Schur rings, and wedge products of
Schur rings (see the next section for definitions). Because these four
families often overlap, special care is taken to ensure that an exact count is
made. The trivial case, as the name suggests, is easy to consider. The
families of wedge and direct products are considered recursively. The
automorphic Schur rings are in one-to-one correspondence with the subgroups of
$\operatorname{Aut}(\mathcal{Z}_{n})\cong\prod_{i=1}^{k}\operatorname{Aut}(\mathcal{Z}_{p_{i}^{e_{i}}})\cong\prod_{i=1}^{k}\mathcal{Z}_{(p_{i}-1)p_{i}^{e_{i}-1}}$,
where $n=\prod_{i=1}^{r}p_{i}^{e_{i}}$ is the prime factorization of $n$.

In [8], Ziv-Av enumerates all Schur rings over small finite groups up to order
63. In [4], the second author enumerates all Schur rings over cyclic groups up
to order 400. In both cases, this was accomplished by computer software. In
all instances, the two enumerations agree with the formulas found herein.

##  2 Traditional Schur Rings

In this section, we remind the reader of important counting techniques
introduced in [4].

For a Schur ring $\mathfrak{S}$ over $G$ and a subgroup $H\leq G$, we say $H$
is a an $\mathfrak{S}$-_subgroup_ if $H$ can be partitioned using the
primitive sets of $\mathfrak{S}$. Then
$\mathfrak{S}_{H}:=\mathfrak{S}\cap\mathbb{Q}[H]$ is a Schur ring over
$\mathcal{H}$ and is called a _Schur subring_ of $\mathfrak{S}$. We say a
Schur ring $\mathfrak{S}$ is _primitive_ if the only $\mathfrak{S}$-subgroups
are $1$ and $G$.

For any subgroup $\mathcal{H}\leq\operatorname{Aut}(G)$, let $G^{\mathcal{H}}$
denote the _automorphic Schur ring_ associated to $\mathcal{H}$. In the case
$\mathcal{H}=1$, $G^{1}=\mathbb{Q}[G]$, the group algebra itself. For
simplicity of notation, this Schur ring, called the _discrete Schur ring_ , is
simply denoted as $G$ as the coefficient ring will provide little consequence.
In the case that $G$ is abelian and $\mathcal{H}=\langle*\rangle$, we denote
$G^{\mathcal{H}}$ as $G^{\pm}$.

For any $n$, there is exactly one _trivial Schur ring_ over $G$, namely
$G^{0}:=\operatorname{Span}\\{1,G\smallsetminus 1\\}$. The only primitive
Schur ring over $\mathcal{Z}_{n}$ when $n$ is composite is this trivial ring.
Thus, the trivial family is disjoint from the other three traditional families
and contains exactly one ring. On the other hand, all Schur rings over
$\mathcal{Z}_{p}$, where $p$ is prime, are primitive. As observed in [4, 6],
all Schur rings over $\mathcal{Z}_{p}$ are automorphic
$\left(\mathcal{Z}_{p}^{0}=\mathcal{Z}_{p}^{\operatorname{Aut}(\mathcal{Z}_{p})}\right)$
and $\Omega(p)=x$, where $x$ is the number of divisors of $p-1$.

If $G=H\times K$, $\mathfrak{S}$ is a Schur ring over $H$, and $\mathfrak{T}$
is a Schur ring over $K$, then
$\mathfrak{S}\times\mathfrak{T}=\mathfrak{S}\otimes_{\mathbb{Q}}\mathfrak{T}$
denotes the _direct product_ of $\mathfrak{S}$ and $\mathfrak{T}$. Note that
both $H$ and $K$ are necessarily subgroups of
$\mathfrak{S}\times\mathfrak{T}$. In fact,
$(\mathfrak{S}\times\mathfrak{T})_{H}=\mathfrak{S}$ and
$(\mathfrak{S}\times\mathfrak{T})_{K}=\mathfrak{T}$. Let $H$ and $K$ be groups
and let $\mathcal{H}\leq\operatorname{Aut}(H)$ and
$\mathcal{K}\leq\operatorname{Aut}(K)$. If $G=H\times K$, then we may
naturally view $\mathcal{H}\times\mathcal{K}$ as a subgroup of
$\operatorname{Aut}(G)$. As observed in [4], if $G=H\times K$ then
$\mathfrak{S}=\mathfrak{S}_{H}\times\mathfrak{S}_{K}$ is automorphic if and
only if $\mathfrak{S}_{H}$ and $\mathfrak{S}_{K}$ are automorphic.

A _section_ $U$ of $G$ is a pair of subgroups $U=[K,H]$ such that $1\leq K\leq
H\leq G$ and $K\trianglelefteq G$. We say that a section $U=[K,H]$ is _proper_
if $1<K\leq H<G$. We say that a section is _trivial_ if $K=H$. As all
subgroups of $\mathcal{Z}_{n}$ are normal and uniquely determined by their
orders, we shall denote the section $[\mathcal{Z}_{d},\mathcal{Z}_{e}]$ simply
as $[d,e]$.

Given any proper section $U=[K,H]$, a Schur ring $\mathfrak{S}$ over $H$, and
a Schur ring $\mathfrak{T}$ over $K$, we form the _wedge product_
$\mathfrak{S}\wedge_{U}\mathfrak{T}$ by constructing the common refinement of
$\mathcal{D}(\mathfrak{S})$ and $\mathcal{D}(\pi^{-1}(\mathfrak{T}))$, where
$\pi:G\to G/K$ is the natural map. To guarantee that
$(\mathfrak{S}\wedge_{U}\mathfrak{T})_{H}=\mathfrak{S}$ and
$\pi(\mathfrak{S}\wedge_{U}\mathfrak{T})=\mathfrak{T}$, the extra
compatibility condition $\pi(\mathfrak{S})=\mathfrak{T}_{H/K}$ is required in
this construction. When $U=[H,H]$ is trivial, compatibility is automatic. We
say a Schur ring is _wedge-decomposable_ when there exists a proper section
such that the Schur ring can be expressed as a wedge product of two other
Schur rings. In particular, $H$ and $K$ are
$\mathfrak{S}\wedge_{U}\mathfrak{T}$-subgroups. Otherwise, we say the Schur
ring is _wedge-indecomposable_.

As explained in [4], every Schur ring has a _wedge core_ , which is a maximal
wedge-indecomposable Schur subring. Let $\Omega(n,\mathfrak{S})$ be the number
of Schur rings over $\mathcal{Z}_{n}$ for which $\mathfrak{S}$ is its wedge
core. Let $\Omega(n,d)=\Omega(n,\mathcal{Z}_{d})$. Clearly,
$\Omega(n)=\sum_{\mathfrak{S}}\Omega(n,\mathfrak{S})$, where the sum ranges
over the indecomposable Schur subrings of $\mathcal{Z}_{n}$. We note that the
indecomposable Schur rings are necessarily trivial, direct products of
indecomposable Schur rings, or automorphic (of course, not every automorphic
ring is indecomposable). The first case is trivial to count and the second is
also enumerated recursively.

Lastly, to count automorphic Schur rings over $\mathcal{Z}_{n}$, we note that
there is a one-to-one correspondence between subgroups of
$\operatorname{Aut}(\mathcal{Z}_{n})$ and automorphic Schur rings. So, it
suffices to count the number of subgroups of
$\operatorname{Aut}(\mathcal{Z}_{n})$, that is, compute
$|\mathcal{L}(\operatorname{Aut}(\mathcal{Z}_{n}))|$. The problem of counting
the number of subgroups of an abelian group is a well studied problem in the
literature (for example, see the references in [4]). The following lemma due
to Călugăreanu [1] will be useful in our enumeration of automorphic Schur
rings.

######  Lemma 2.1.

The number of subgroups of $\mathcal{Z}_{p^{k}}\times\mathcal{Z}_{p^{\ell}}$
is given as

$\left|\mathcal{L}\middle(\mathcal{Z}_{p^{k}}\times\mathcal{Z}_{p^{\ell}}\middle)\right|=\sum_{j=0}^{\min(k,\ell)}\phi(p^{j})(k-j+1)(\ell-j+1),$

where $\phi$ denotes Euler’s totient function.

##  3 Proof of Theorem 1.1

We proceed to the proof of Theorem 1.1. Consider $\mathcal{Z}_{pq}$, where
$p=\prod_{i=1}^{n}r_{i}^{k_{i}}+1$ and $q=\prod_{i=1}^{n}r_{i}^{\ell_{i}}+1$
are primes. By the Fundamental Theorem, all the Schur rings over this group
are traditional. So, we now consider each of these families. As observed
above, since $pq$ is composite, the trivial Schur ring contributes one to the
count of $\Omega(pq)$ and does not intersect the other three families.
Additionally, as all Schur rings over $\mathcal{Z}_{p}$ and $\mathcal{Z}_{q}$
are necessarily automorphic, all direct product rings over $\mathcal{Z}_{pq}$
are automorphic and will be counted under the automorphic family.

We consider next the wedge-decomposable Schur rings over $\mathcal{Z}_{pq}$.
Over $\mathcal{Z}_{pq}$, the only possible proper sections are $[p,p]$ and
$[q,q]$, which are trivial. In the first case, any possible Schur ring over
$\mathcal{Z}_{p}$ could be wedged with any possible Schur ring over
$\mathcal{Z}_{q}$. This produces $\Omega(p)\Omega(q)=xy$ distinct Schur rings.
The subgroups of each of these Schur rings will be exactly $1$,
$\mathcal{Z}_{p}$, and $\mathcal{Z}_{pq}$. Note that $\mathcal{Z}_{q}$ is
missing since elements of order $q$ and $pq$ are fused together in cosets of
$\mathcal{Z}_{p}$. Thus, these Schur rings are not automorphic. The second
case is similar and produces $xy$ distinct Schur rings from the automorphic
ring and the wedge products already accounted for. Therefore, there are $2xy$
wedge products over $\mathcal{Z}_{pq}$.

Note that by the above decomposition of $p$, $p-1$ has
$\prod_{i=1}^{n}(k_{i}+1)$ divisors. Likewise, $q-1$ has
$\prod_{i=1}^{n}(\ell_{1}+1)$ divisors. Considering the wedge-decomposable and
trivial families, we have already accounted for
$2\prod_{i=1}^{n}(k_{i}+1)(\ell_{i}+1)+1$ distinct, non-automorphic Schur
rings. Thus, to prove Theorem 1.1 it suffices to count the number of distinct
automorphic Schur rings over $\mathcal{Z}_{pq}$. As these Schur rings are in
direct one-to-one correspondence with the subgroups of
$\operatorname{Aut}(\mathcal{Z}_{pq})$, we see that

$\Omega(pq)=|\mathcal{L}(\operatorname{Aut}(\mathcal{Z}_{pq}))|+2\prod_{i=1}^{n}(k_{i}+1)(\ell_{i}+1)+1.$
(3.1)

For a finite group $G$, let $G=\prod_{i=1}^{k}P_{i}$ be its primary
decomposition. Then it is well-known that
$\mathcal{L}(G)\cong\prod_{i=1}^{k}\mathcal{L}(P_{i})$. In the case of
$\mathcal{Z}_{pq}$, we see that
$\operatorname{Aut}(\mathcal{Z}_{pq})\cong\mathcal{Z}_{p-1}\times\mathcal{Z}_{q-1}$.
Hence, the primary decomposition for $\operatorname{Aut}(\mathcal{Z}_{pq})$ is
given as

$\operatorname{Aut}(\mathcal{Z}_{pq})\cong\prod_{i=1}^{n}(\mathcal{Z}_{r_{i}^{k_{i}}}\times\mathcal{Z}_{r_{i}^{\ell_{i}}}).$

By Lemma 2.1,

$\left|\mathcal{L}\middle(\mathcal{Z}_{r_{i}^{k_{i}}}\times\mathcal{Z}_{r_{i}^{\ell_{i}}}\middle)\right|=\sum_{j=0}^{\min(k_{i},\ell_{i})}\phi(r_{i}^{j})(k_{i}-j+1)(\ell_{i}-j+1).$

Therefore,

$|\mathcal{L}(\operatorname{Aut}(\mathcal{Z}_{pq}))|=\prod_{i=1}^{n}\sum_{j=0}^{\min(k_{i},\ell_{i})}\phi(r_{i}^{j})(k_{i}-j+1)(\ell_{i}-j+1).$
(3.2)

Finally, Theorem 1.1 follows immediately from (3.1) and (3.2), which finishes
the proof.

######  Example 3.1.

We present a complete enumeration of the Schur rings over $\mathcal{Z}_{21}$
as an example to illustrate the previous proof. There are $\Omega(3)=2$ Schur
rings over $\mathcal{Z}_{3}$, namely $\mathcal{Z}_{3}^{0}$ and
$\mathcal{Z}_{3}$. There are $\Omega(7)=4$ Schur rings over $\mathcal{Z}_{7}$,
namely, $\mathcal{Z}_{7}^{0}$, $\mathcal{Z}_{7}^{\langle 2\rangle}$,
$\mathcal{Z}_{7}^{\pm}$, and $\mathcal{Z}_{7}$.

Below we list the $\Omega(21)=27$ Schur rings over $\mathcal{Z}_{21}$:

$\mathcal{Z}_{3}^{0}\wedge\mathcal{Z}_{7}^{0},\mathcal{Z}_{3}^{0}\wedge\mathcal{Z}_{7}^{\langle
2\rangle},\mathcal{Z}_{3}^{0}\wedge\mathcal{Z}_{7}^{\pm},\mathcal{Z}_{3}^{0}\wedge\mathcal{Z}_{7},\mathcal{Z}_{3}\wedge\mathcal{Z}_{7}^{0},\mathcal{Z}_{3}\wedge\mathcal{Z}_{7}^{\langle
2\rangle},\mathcal{Z}_{3}\wedge\mathcal{Z}_{7}^{\pm},\mathcal{Z}_{3}\wedge\mathcal{Z}_{7},$
$\mathcal{Z}_{7}^{0}\wedge\mathcal{Z}_{3}^{0},\mathcal{Z}_{7}^{\langle
2\rangle}\wedge\mathcal{Z}_{3}^{0},\mathcal{Z}_{7}^{\pm}\wedge\mathcal{Z}_{3}^{0},\mathcal{Z}_{7}\wedge\mathcal{Z}_{3}^{0},\mathcal{Z}_{7}^{0}\wedge\mathcal{Z}_{3},\mathcal{Z}_{7}^{\langle
2\rangle}\wedge\mathcal{Z}_{3},\mathcal{Z}_{7}^{\pm}\wedge\mathcal{Z}_{3},\mathcal{Z}_{7}\wedge\mathcal{Z}_{3},$
$\mathcal{Z}_{3}^{0}\times\mathcal{Z}_{7}^{0},\mathcal{Z}_{3}^{0}\times\mathcal{Z}_{7}^{\langle
2\rangle},\mathcal{Z}_{3}^{0}\times\mathcal{Z}_{7}^{\pm},\mathcal{Z}_{3}^{0}\times\mathcal{Z}_{7},\mathcal{Z}_{3}\times\mathcal{Z}_{7}^{0},\mathcal{Z}_{3}\times\mathcal{Z}_{7}^{\langle
2\rangle},\mathcal{Z}_{3}\times\mathcal{Z}_{7}^{\pm},\mathcal{Z}_{3}\times\mathcal{Z}_{7}\
(\cong\mathcal{Z}_{21}),$ $\mathcal{Z}_{21}^{0},\mathcal{Z}_{21}^{\langle
5\rangle},\mathcal{Z}_{21}^{\pm}.$

##  4 Proof of Theorem 1.6

Let $p=2^{k}a+1$, where $a$ is odd. Let $x$ denote the number of divisors of
$p-1$. Then $x/(k+1)$ is the number of divisors of $a$.

The construction of wedge products will be more complicated in the case $n=4p$
as wedge-compatibility is not necessarily trivial. There is also the concern
that different sections can potentially create the same wedge product, e.g. if
$U=[2,4]$ and $U^{\prime}=[4,4]$, then
$\mathcal{Z}_{4}\wedge_{U}\mathcal{Z}_{6}=\mathcal{Z}_{4}\wedge_{U^{\prime}}\mathcal{Z}_{3}$.
Additionally, the wedge product construction is associative, e.g.
$(\mathcal{Z}_{2}\wedge\mathcal{Z}_{2})\wedge\mathcal{Z}_{3}=\mathcal{Z}_{2}\wedge(\mathcal{Z}_{2}\wedge\mathcal{Z}_{3})$.
To avoid these issues, we will restrict our attention to wedge-decomposable
Schur rings $\mathfrak{S}\wedge_{U}\mathfrak{T}$ where $U=[K,H]$, $K$ is a
minimal $\mathfrak{S}$-subgroup, and $\mathfrak{S}$ is wedge-indecomposable,
that is, $\mathfrak{S}$ is the wedge-core of
$\mathfrak{S}\wedge_{U}\mathfrak{T}$.

We will organize Schur rings over $\mathcal{Z}_{4p}$ according to the order of
its wedge-core. We need to consider $\Omega(4p,d)$ for the divisors
$d=2,p,4,2p$ ($d=1,4p$ are omitted as these are not proper sections). If
$\mathfrak{S}$ is primitive, the only possible section would be trivial. As
such, $\Omega(n,\mathfrak{S})=\Omega(n/d)$. For example,
$\Omega(4p,2)=\Omega(2p)=3x+1$ and $\Omega(4p,p)=\Omega(4)=3$. In fact, as all
the $x$ many Schur rings $\mathfrak{S}$ over $\mathcal{Z}_{p}$ are primitive,
$\Omega(4p,\mathfrak{S})=\Omega(4)$. Hence, there are $\Omega(p)\Omega(4)=3x$
many Schur rings over $\mathcal{Z}_{4p}$ with a core of order $p$.

The three Schur rings over $\mathcal{Z}_{4}$ are $\mathcal{Z}_{4}^{0}$,
$\mathcal{Z}_{2}\wedge\mathcal{Z}_{2}$, and $\mathcal{Z}_{4}$, where the first
and last ones are indecomposable. As $\mathcal{Z}_{4}^{0}$ is primitive,
$\Omega(4p,\mathcal{Z}_{4}^{0})=\Omega(p)=x$. For
$\mathfrak{S}=\mathcal{Z}_{4}$ we select the section $U=[2,4]$. Hence,
$\Omega(4p,4)$ counts the number of Schur rings of the form
$\mathcal{Z}_{4}\wedge_{U}\mathfrak{T}$ where $\mathfrak{T}$ is a Schur ring
over $\mathcal{Z}_{2p}$ and $\mathcal{Z}_{2}$ is a subring. To count the
number of such rings $\mathfrak{T}$ we revisit the proof of Theorem 1.1 with
the simplification $q=2$, which shows that there are exactly $x$ many Schur
rings of the form $\mathcal{Z}_{2}\wedge\mathfrak{S}$, $x$ many Schur rings of
the form $\mathfrak{S}\wedge\mathcal{Z}_{2}$, $x$ many Schur rings of the form
$\mathcal{Z}_{2}\times\mathfrak{S}$, and one trivial ring
$\mathcal{Z}_{2p}^{0}$. Hence, there are $2x$ many Schur rings over
$\mathcal{Z}_{2p}$ that contain the subring $\mathcal{Z}_{2}$, that is,
$\Omega(4p,4)=2x$. Therefore, there are $3x$ many Schur rings over
$\mathcal{Z}_{4p}$ with a core of order 4.

It is also observed that there are $x+1$ indecomposable Schur rings over
$\mathcal{Z}_{2p}$, which $x$ have the form
$\mathcal{Z}_{2}\times\mathfrak{S}^{\prime}$ and one is trivial. Note
$\Omega(4p,\mathcal{Z}_{2p}^{0})=\Omega(2)=1$, that is,
$\mathcal{Z}_{2p}^{0}\wedge\mathcal{Z}_{2}$ is the only option. If
$\mathfrak{S}=\mathcal{Z}_{2}\times\mathfrak{S}^{\prime}$ is the core, then
$\mathfrak{S}$ has two distinct minimal subgroups, namely, $\mathcal{Z}_{2}$
and $\mathcal{Z}_{p}$. As such, two proper sections need to be considered,
$[2,2p]$ and $[p,2p]$. For fixed $\mathfrak{S}^{\prime}$, there are $2$ Schur
rings over $\mathcal{Z}_{2p}$ which contain $\mathfrak{S}^{\prime}$ as a
subring, namely $\mathfrak{S}^{\prime}\wedge\mathcal{Z}_{2}$ and
$\mathfrak{S}^{\prime}\times\mathcal{Z}_{2}$. Hence, there are two
possibilities for
$(\mathcal{Z}_{2}\times\mathfrak{S}^{\prime})\wedge_{[2,2p]}\mathfrak{T}$. For
$(\mathcal{Z}_{2}\times\mathfrak{S}^{\prime})\wedge_{[p,2p]}\mathfrak{T}$, we
count the number of Schur rings over $\mathcal{Z}_{4}$ which contain
$\mathcal{Z}_{2}$ as a subgroup. There are two such rings,
$\mathcal{Z}_{2}\wedge\mathcal{Z}_{2}$ and $\mathcal{Z}_{4}$. Lastly, there is
$\Omega(2)=1$ ring of the form
$(\mathcal{Z}_{2}\times\mathfrak{S}^{\prime})\wedge_{[2p,2p]}\mathfrak{T}$,
namely $(\mathcal{Z}_{2}\times\mathfrak{S}^{\prime})\wedge\mathcal{Z}_{2}$.
Therefore, there are $2+2-1=3$ many Schur rings over $\mathcal{Z}_{4p}$ with
$\mathcal{Z}_{2}\times\mathfrak{S}^{\prime}$ as its core. This accounts for
$3x+1$ many Schur rings over $\mathcal{Z}_{4p}$ with a core of order $2p$.

In summary, this accounts for $(3x+1)+3x+3x+(3x+1)=12x+2$ many wedge-
decomposable Schur rings over $\mathcal{Z}_{4p}$. It remains to consider the
indecomposable Schur rings over $\mathcal{Z}_{4p}$. There is, of course, the
trivial Schur ring $\mathcal{Z}_{4p}^{0}$, as well as the indecomposable
automorphic Schur rings and those direct products of the form
$\mathcal{Z}_{4}^{0}\times\mathfrak{S}$ for some Schur rings $\mathfrak{S}$ of
order $p$. In regard to the automorphic Schur rings, we know the total count
is equal to $|\mathcal{L}(\operatorname{Aut}(\mathcal{Z}_{4p}))|$. The
sublattice
$\mathcal{L}(\operatorname{Aut}(\mathcal{Z}_{4}))\times\mathcal{L}(\operatorname{Aut}(\mathcal{Z}_{p}))$
will consist of two lattice-isomorphic copies of
$\mathcal{L}(\mathcal{Z}_{p-1})$. The full lattice
$\mathcal{L}(\operatorname{Aut}(\mathcal{Z}_{4p}))$ contains these two layers
and all the diagonal entries that sit between the top and bottom layers. Those
non-diagonal automorphic Schur rings in the top layer have the form
$\mathcal{Z}_{4}\times\mathfrak{S}$, for some Schur ring $\mathfrak{S}$ of
order $p$, and are indecomposable. Those non-diagonal automorphic Schur rings
in the bottom layer have the form
$(\mathcal{Z}_{2}\wedge\mathcal{Z}_{2})\times\mathfrak{S}$ and are wedge-
decomposable. Notice that these decomposable automorphic Schur rings are in
one-to-one correspondence with those indecomposable Schur rings of the form
$\mathcal{Z}_{4}^{0}\times\mathfrak{S}$. Thus, if every diagonal automorphic
Schur ring over $\mathcal{Z}_{4p}$ is indecomposable, which we claim, then the
number of indecomposable Schur rings over $\mathcal{Z}_{4p}$ is
$1+|\mathcal{L}(\operatorname{Aut}(\mathcal{Z}_{4p}))|$. Using Lemma 2.1, we
see that
$\operatorname{Aut}(\mathcal{Z}_{4p})\cong\mathcal{Z}_{2}\times\mathcal{Z}_{p-1}\cong(\mathcal{Z}_{2}\times\mathcal{Z}_{2^{k}})\times\mathcal{Z}_{a}$
and

$|\mathcal{L}(\operatorname{Aut}(\mathcal{Z}_{4p}))|=|\mathcal{L}(\mathcal{Z}_{2}\times\mathcal{Z}_{2^{k}})||\mathcal{L}(\mathcal{Z}_{a})|=(2(k+1)+k)\left(\dfrac{x}{k+1}\right)=\dfrac{3k+2}{k+1}x.$

To prove the claim, we introduce the representation
$\omega:\mathbb{Q}[\mathcal{Z}_{n}]\to\mathbb{Q}(\zeta_{n})$ which maps a
generator of $\mathcal{Z}_{n}$ to $\zeta_{n}:=e^{2\pi i/n}$. We remind the
reader that in [4] we saw that an automorphic Schur ring $\mathfrak{S}$ is
wedge-decomposable if and only if
$\omega(\mathfrak{S})\leq\mathbb{Q}(\zeta_{d})$ for some proper divisor $d\mid
n$ (excluding, of course, the case when $n$ is prime). By definition, the
diagonal automorphic Schur rings are not contained in
$\mathcal{L}(\operatorname{Aut}(\mathcal{Z}_{4}))\times\mathcal{L}(\operatorname{Aut}(\mathcal{Z}_{p}))$,
which implies that their image is not contained in
$\mathcal{L}(\mathbb{Q}(\zeta_{4}))$ or $\mathcal{L}(\mathbb{Q}(\zeta_{p}))$.
This implies they are all indecomposable, as claimed.

As we have now exhausted all possibilities, we see that

$\Omega(4p)=(12x+2)+1+\dfrac{3k+2}{k+1}x=\dfrac{15k+14}{k+1}x+3,$

which finishes the proof of Theorem 1.6.

######  Example 4.1.

We present a complete enumeration of the Schur rings over $\mathcal{Z}_{12}$
as an example to illustrate the previous proof. Note $12=4(3)$ and
$3=2^{1}\cdot 1+1$. There are $\Omega(3)=2$ Schur rings over
$\mathcal{Z}_{3}$, namely $\mathcal{Z}_{3}^{0}$ and $\mathcal{Z}_{3}$. There
are $\Omega(6)=7$ Schur rings over $\mathcal{Z}_{6}$, namely

$\mathcal{Z}_{2}\wedge\mathcal{Z}_{3}^{0},\mathcal{Z}_{2}\wedge\mathcal{Z}_{3},\mathcal{Z}_{3}^{0}\wedge\mathcal{Z}_{2},\mathcal{Z}_{3}\wedge\mathcal{Z}_{2},\mathcal{Z}_{6}^{0},\mathcal{Z}_{2}\times\mathcal{Z}_{3}^{0}\
(=\mathcal{Z}_{6}^{\pm}),\mathcal{Z}_{2}\times\mathcal{Z}_{3}\
(=\mathcal{Z}_{6}).$

Below we list the $\Omega(12)=32$ Schur rings over $\mathcal{Z}_{12}$:

$\mathcal{Z}_{2}\wedge\mathcal{Z}_{2}\wedge\mathcal{Z}_{3}^{0},\mathcal{Z}_{2}\wedge\mathcal{Z}_{2}\wedge\mathcal{Z}_{3},\mathcal{Z}_{2}\wedge\mathcal{Z}_{3}^{0}\wedge\mathcal{Z}_{2},\mathcal{Z}_{2}\wedge\mathcal{Z}_{3}\wedge\mathcal{Z}_{2},\mathcal{Z}_{2}\wedge\mathcal{Z}_{6}^{0},\mathcal{Z}_{2}\wedge\mathcal{Z}_{6}^{\pm},\mathcal{Z}_{2}\wedge\mathcal{Z}_{6},$
$\mathcal{Z}_{3}^{0}\wedge\mathcal{Z}_{4}^{0},\mathcal{Z}_{3}^{0}\wedge\mathcal{Z}_{2}\wedge\mathcal{Z}_{2},\mathcal{Z}_{3}^{0}\wedge\mathcal{Z}_{4},\mathcal{Z}_{3}\wedge\mathcal{Z}_{4}^{0},\mathcal{Z}_{3}\wedge\mathcal{Z}_{2}\wedge\mathcal{Z}_{2},\mathcal{Z}_{3}\wedge\mathcal{Z}_{4},$
$\mathcal{Z}_{4}^{0}\wedge\mathcal{Z}_{3}^{0},\mathcal{Z}_{4}^{0}\wedge\mathcal{Z}_{3},\mathcal{Z}_{4}\wedge_{[2,4]}(\mathcal{Z}_{2}\wedge\mathcal{Z}_{3}^{0}),\mathcal{Z}_{4}\wedge_{[2,4]}(\mathcal{Z}_{2}\wedge\mathcal{Z}_{3}),\mathcal{Z}_{4}\wedge_{[2,4]}\mathcal{Z}_{6}^{\pm},\mathcal{Z}_{4}\wedge_{[2,4]}\mathcal{Z}_{6},$
$\mathcal{Z}_{6}^{0}\wedge\mathcal{Z}_{2},\mathcal{Z}_{6}^{\pm}\wedge\mathcal{Z}_{2},\mathcal{Z}_{6}^{\pm}\wedge_{[2,6]}\mathcal{Z}_{6}^{\pm},\mathcal{Z}_{6}^{\pm}\wedge_{[3,6]}\mathcal{Z}_{4},\mathcal{Z}_{6}\wedge\mathcal{Z}_{2},\mathcal{Z}_{6}\wedge_{[2,6]}\mathcal{Z}_{6},\mathcal{Z}_{6}\wedge_{[3,6]}\mathcal{Z}_{4},$
$\mathcal{Z}_{12}^{0},\mathcal{Z}_{4}^{0}\times\mathcal{Z}_{3}^{0},\mathcal{Z}_{4}^{0}\times\mathcal{Z}_{3},\mathcal{Z}_{4}\times\mathcal{Z}_{3}^{0},\mathcal{Z}_{4}\times\mathcal{Z}_{3}\
(=\mathcal{Z}_{12}),\mathcal{Z}_{12}^{\pm}.$

## References

  * [1] Grigore Călugăreanu.  The total number of subgroups of a finite Abelian group.  Sci. Math. Jpn., 60(1):157–167, 2004.
  * [2] Ka Hin Leung and Shin Hing Man.  On Schur rings over cyclic groups II.  Journal of Algebra, 183:273–285, 1996.
  * [3] Ka Hin Leung and Shin Hing Man.  On Schur rings over cyclic groups.  Israel Journal of Mathematics, 106:251–267, 1998.
  * [4] Andrew Misseldine.  On counting Schur rings over cyclic groups.  pre-print (submitted).
  * [5] Andrew Misseldine.  Algebraic and Combinatorial Properties of Schur Rings over Cyclic Groups.  PhD thesis, Brigham Young University, 2014.
  * [6] Andrew Misseldine.  Counting Schur rings over cyclic groups.  Journal of Algebraic Combinatorics, 51:155–169, Feb 2020.
  * [7] Mikhail Muzychuk and Ilia Ponomarenko.  Schur rings.  European Journal of Combinatorics, 30:1526–1539, 2009.
  * [8] Matan Ziv-Av.  Enumeration of schur rings over small groups.  In Computer algebra in scientific computing. 16th international workshop, CASC 2014, Warsaw, Poland, September 8–12, 2014. Proceedings, pages 491–500. Berlin: Springer, 09 2014.