A measurement of the $K^{+} \to \pi^{+} \mu^{+} \mu^{-}$ decay

A sample of $2.8 \times 10^4$ $K^{+} \to \pi^{+} \mu^{+} \mu^{-}$ candidates with negligible background was collected by the NA62 experiment at the CERN SPS in 2017--2018. The model-independent branching fraction is measured to be $(9.15 \pm 0.08)\times 10^{-8}$, a factor three more precise than previous measurements. The decay form factor is presented as a function of the squared dimuon mass. A measurement of the form factor parameters and their uncertainties is performed using a description based on Chiral Perturbation Theory at $\mathcal{O}(p^6)$.

Studies of the K πee and K πµµ decay form factors contribute to experimental tests of lepton flavour universality [5,6]. The first lattice QCD calculation of the form factor value at a specific lepton pair mass (lying outside the K πµµ kinematic region) using physical light-quark masses is presented in [7]. Future methodology optimizations together with advances in computing technology are expected to provide competitive lattice QCD predictions of the form factor.
The E787 collaboration at the Brookhaven National Laboratory reported the first observation of the K πµµ decay in 1997 [8], which was followed by the E865 [9] and HyperCP [10] measurements. The E865 result established the vector nature of the decay form factor, while HyperCP studied both K + πµµ and K − πµµ decays and measured the CP violating decay rate asymmetry, found to be compatible with zero. The most precise study [11] of K πµµ was performed by the NA48/2 collaboration at the CERN SPS. The K πee decay was first observed at the CERN PS by the Geneva-Saclay collaboration in 1975 [12], and subsequently measured by the E777 [13], E865 [14] and NA48/2 [15] experiments. A summary of form factor measurements can be found in [16].
Improved measurements of the K πµµ model-independent branching fraction and form factor parameters, based on the dataset collected in 2017-2018 by the NA62 experiment at the CERN SPS, are presented in the following. The forward-backward asymmetry of the decay with respect to angle θ Kµ between the K + and the µ − three-momenta in the µ + µ − rest frame, is also measured.

Beam, detector and data sample
The layout of the NA62 beamline and detector [17] is shown schematically in Figure 1. An unseparated secondary beam of π + (70%), protons (23%) and K + (6%) is created by directing 400 GeV/c protons extracted from the CERN SPS onto a beryllium target in spills of 3 s effective duration. The target position defines the origin of the NA62 reference system: the beam travels along the Z axis in the positive direction (downstream), the Y axis points vertically up, and the X axis is horizontal and directed to form a right-handed coordinate system. The central beam momentum is 75 GeV/c, with a momentum spread of 1% (rms).
Beam kaons are tagged with a time resolution of 70 ps by a differential Cherenkov counter (KTAG), which uses nitrogen gas at 1.75 bar pressure contained in a 5 m long vessel as radiator. Beam particle positions, momenta and times (to better than 100 ps resolution) are measured by a silicon pixel spectrometer consisting of three stations (GTK1,2,3) and four dipole magnets. A toroidal muon sweeper, called scraper (SCR), is installed between GTK1 and GTK2. A 1.2 m thick steel collimator (COL) with a 76 × 40 mm 2 central aperture and 1.7 × 1.8 m 2 outer dimensions is placed upstream of GTK3 to absorb hadrons from upstream K + decays; a variable aperture collimator of 0.15 × 0.15 m 2 outer dimensions was used up to early 2018. Inelastic interactions of beam particles in GTK3 are detected by an array of scintillator hodoscopes (CHANTI). A dipole magnet (TRIM5) providing a 90 MeV/c horizontal momentum kick is located in front of GTK3. The beam is delivered into a vacuum tank evacuated to a pressure of 10 −6 mbar, which contains a 75 m long fiducial volume (FV) starting 2.6 m downstream of GTK3. The beam angular spread at the FV entrance is 0.11 mrad (rms) in both horizontal and vertical planes. Downstream of the FV, undecayed beam particles continue their path in vacuum. Three-momenta of charged particles produced in K + decays are measured by a magnetic spectrometer (STRAW) located in the vacuum tank downstream of the FV. The spectrometer consists of four tracking chambers made of straw tubes, and a large aperture dipole magnet (M), located between the second and third chamber, that provides a horizontal momentum kick of 270 MeV/c. The momentum resolution is σ p /p = (0.30 ⊕ 0.005 · p)%, with the momentum p expressed in GeV/c.
A ring-imaging Cherenkov detector (RICH) consisting of a 17.5 m long vessel filled with neon at atmospheric pressure (with a Cherenkov threshold of 12.5 GeV/c for pions) provides particle identification, charged particle time measurements (to a 70 ps accuracy for particles well above the Cherenkov threshold), and the trigger time. Two scintillator hodoscopes (CHOD), which include a matrix of tiles and two planes of slabs arranged in four quadrants located downstream of the RICH, provide trigger signals and time measurements. The tile matrix hodoscope has a time resolution of 1 ns, while the slab hodoscope measures time with 200 ps precision.
A 27X 0 thick quasi-homogeneous liquid krypton (LKr) electromagnetic calorimeter is used for particle identification and photon detection. The calorimeter has an active volume of 7 m 3 segmented in the transverse direction into 13248 projective cells of 2 × 2 cm 2 size, and provides an energy resolution σ E /E = (4.8/ √ E ⊕ 11/E ⊕ 0.9)%, with E expressed in GeV. To achieve hermetic acceptance for photons emitted in K + decays in the FV at angles up to 50 mrad from the beam axis, the LKr calorimeter is supplemented by annular lead glass detectors (LAV) installed in 12 positions inside and downstream of the vacuum tank, and two lead/scintillator sampling calorimeters (IRC, SAC) located close to the beam axis. An iron/scintillator sampling hadronic calorimeter formed of two modules (MUV1,2) and a muon detector consisting of 148 scintillator tiles located behind an 80 cm thick iron wall (MUV3) are used for particle identification. The eight smaller tiles of MUV3 adjacent to the beam pipe are referred to as the inner tiles, while the remaining 140 regular tiles are called the outer tiles.
The data sample used for this analysis is obtained from 0.84 × 10 6 SPS spills collected in 2017-2018, with the typical beam intensity increasing over time from 1.5 × 10 12 to 2.2 × 10 12 protons per spill. The latter value corresponds on average to a 500 MHz instantaneous beam particle rate at the FV entrance, and a 3.7 MHz K + decay rate in the FV. The main trigger of NA62 is dedicated to the collection of the very rare K + → π + νν decays [18]. Multi-track (MT) and di-muon multi-track (2µMT) triggers considered in this analysis operate concurrently, downscaled by typical factors of 100 and 2, respectively. The downscaling factors of both triggers were varied throughout the data taking. The MT trigger line selects the K + → π + π + π − (K 3π ) decays, used for normalization, while the 2µMT line selects the K πµµ signal decays. The low-level hardware (L0) trigger [19] for both lines is based on RICH signal multiplicity and coincidence of signals in two opposite CHOD quadrants. The 2µMT line additionally involves a requirement of signal coincidence in two outer MUV3 tiles. The high-level software (L1) trigger requires K + identification by KTAG, and reconstruction of a negatively charged STRAW track for both MT and 2µMT trigger lines. A detailed description of the NA62 trigger system and its performance is given in [20].
Monte Carlo (MC) simulations of particle interactions with the detector and its response are performed using a software package based on the Geant4 toolkit [21]. In addition, the accidental activity is simulated, and the response of both trigger lines is emulated.

Event selection
Kinematic similarities of the signal (K πµµ ) and normalization (K 3π ) decays allow for substantial overlap between the signal and normalization event selections, which results in first-order cancellation of most detector and trigger inefficiencies, thus reducing the systematic uncertainties in the measurement.
The following selection criteria are common to the K πµµ and K 3π event selections.
• Each STRAW track is assigned a time computed as a weighted average of the associated CHOD hodoscope signals. The weights are obtained from the time resolutions of the CHOD hodoscopes. Triplets of STRAW tracks compatible with a common origin in the FV are combined into three-track vertices. Vertex time is defined as the weighted average of the times of CHOD signals associated with the vertex tracks.
• Exactly one three-track vertex with the following properties is required to be present: total charge q = 1, time within 6 ns of the trigger time, Z position between 110 m and 180 m from the target, total momentum compatible with the mean beam momentum within 2.5 GeV/c, total transverse momentum with respect to the beam axis below 30 MeV/c, and vertex distance from the beam axis below 5 cm. The beam axis and momentum are monitored throughout the data taking with fully reconstructed K 3π decays. Only the three tracks forming the chosen vertex are considered in the following.
• All track times must be within 12 ns of the vertex time, and the vertex time is required to be within 6 ns of a KTAG kaon signal.
• The tracks must be within the geometrical acceptance of all STRAW chambers, and extrapolate to lie within the CHOD, LKr, and MUV3 acceptances.
• The track momenta should exceed 10 GeV/c to ensure track reconstruction efficiency above 90%. The angles between each track and the beam axis must be smaller than 9 mrad to reduce background to the K πµµ sample from K 3π decays followed by π ± → µ ± ν decays.
• The spatial separation between each pair of vertex tracks must be at least 15 mm in the plane of the first STRAW chamber and 200 mm in the LKr front plane to suppress photon conversions and the overlap of energy deposits.
The following particle identification criteria are employed.
• A track is identified as a charged pion if it has no spatially associated MUV3 signals within 10 ns of the vertex time, and the ratio of the associated LKr cluster energy to the track momentum is E/p < 0.9.
• A track with E/p < 0.2 is identified as a muon if it has a spatially associated MUV3 signal in an outer tile within 6 ns of both the vertex and the trigger times.
The following criteria are specific to the K πµµ event selection.
• To reduce the background from K 3π decays occurring upstream of the FV, the track identified as a π + is extrapolated backward to the COL plane, taking into account the TRIM5 magnetic field. The extrapolated position is required to lie outside a rectangle defined by |X| < 40 mm and |Y | < 25 mm.
• Further K 3π background suppression is achieved by requiring the momenta of both muon tracks to be below 45 GeV/c.
• The invariant mass m(πµµ) of the three selected tracks is reconstructed with a 1.1 MeV/c 2 resolution and must be within 8 MeV/c 2 of the nominal K + mass m K [22].
The following criteria are specific to the K 3π event selection.
• In order to minimize differences between the signal and normalization selections, only one positive track, chosen at random, is required to be identified as a π + .
• The identified π + track extrapolated to the COL plane must satisfy the same requirements as the π + in the K πµµ selection.
• The invariant mass m(3π) of the three selected tracks is reconstructed with a 0.8 MeV/c 2 resolution and must be within 8 MeV/c 2 of m K .
For both selections, simulated events are required to be accepted by a set of software algorithms emulating the conditions employed in the online trigger system.

Signal and normalization samples
The reconstructed mass spectra of the data and simulated events passing the signal and normalization event selections are shown in Figure 2. The selected K 3π data sample, contaminated by background decays to a negligible level of 10 −6 , is used together with the simulated K 3π events with inner bremsstrahlung included [23], to obtain the effective number of kaon decays in the FV where the index i runs over data taking periods defined by constant trigger downscaling factors, N i 3π are the numbers of K 3π events selected with the MT trigger with downscaling factor D i MT , D i 2µMT are the downscaling factors of the 2µMT trigger, and A 3π = (6.58 ± 0.16)% and B 3π = (5.583 ± 0.024)% are the acceptance (obtained from simulation) and the branching fraction [22] of the K 3π decay, respectively. The statistical errors in A 3π and N K are negligible, while the systematic uncertainties are dominated by the accuracy of the CHOD detector efficiency in the simulation. The external error on N K stems from the uncertainty on the K 3π branching fraction.
The m(πµµ) signal region contains 27679 data events with a background contamination of about 8 events, estimated from simulation. The arrows indicate the selected mass regions. The contribution from the simulated K πµµ decays is scaled according to the PDG branching fraction [22]. Bottom: ratios of data and simulated spectra for signal (left) and normalization (right).

Interpretation of the data Decay width and form factor parameterization
The one-photon-inclusive K πµµ differential decay width expressed in terms of the normalized µ + µ − invariant mass squared z = m 2 (µ + µ − )/m 2 K reads [2,3,4,24] dΓ(z) is the form factor of the K + → π + µ + µ − transition, and g(z) is a function describing the decay kinematics [3] and including next-to-leading order electromagnetic effects in terms of radiative corrections. While the µ + µ − interactions are fully taken into account by virtual and bremsstrahlung corrections for the lepton and meson contributions, discussed in [4] and extended beyond the soft-photon approximation, the semiclassical Coulomb corrections, summarized for example in [25], are applied to the π + µ + and π + µ − pairs. These last corrections have opposite sign and the same average magnitude; their combined effect on the results of the present analysis is found to be negligible. The hard-photon 4-body (K + → π + µ + µ − γ) part of the phase-space is separated from the soft-photon 3-body (K + → π + µ + µ − ) part by the condition (P π + P γ ) 2 − m 2 π > 100 MeV 2 , where P π and P γ are 4-momenta of the π + and γ, respectively. The cutoff value is optimized with respect to the experimental resolution. The resulting ratio of the 4-body to 3-body integrated decay widths is (1.64±0.02)%, where the uncertainty comes mainly from the accuracy of the theoretical description dΓ 4-body (z)/dz of the 4-body decay [24]. In the present analysis, the 4-body decay width, depending non-trivially on the form factor, is approximated by a unique function displayed in Figure 3-left. Effects of this approximation are treated as systematic uncertainties.

Measurement of the model-independent branching fraction and form factor
The selected K πµµ signal sample with negligible background contamination is distributed in 50 equipopulated bins in z with widths ranging from 0.004 for z ≈ 0.25 to 0.066 for the last bin. The resolution in z increases linearly from zero to 0.0035 within the allowed kinematic range, and is always several times smaller than the corresponding bin width. The reconstructed differential decay width, shown in Figure 3-left, is given by where for each bin i: N πµµ,i is the number of K πµµ signal candidates, ∆z i is the bin width, A πµµ,i is the signal selection acceptance of the K πµµ decay (obtained from simulation, and equal to zero at both kinematic bounds of z while reaching the maximum of 12.5% around z = 0.2, see also , N K is the effective number of kaon decays in the FV collected by the 2µMT trigger (Equation 1), is the reduced Planck constant, and τ K = (1.238 ± 0.002) × 10 −8 s is the mean charged kaon lifetime [22]. The model-independent K πµµ branching fraction B πµµ = (9.15 ± 0.06 stat ) × 10 −8 is obtained from the reconstructed binned differential decay width (Equation 4, Figure 3-left) by integrating the spectrum over z and multiplying by τ K / .
The K πµµ data sample is also used to extract the |W (z)| 2 form factor (Figure 3-right). The values of the |W (z)| 2 function are reconstructed from the differential decay spectrum (Figure 3left) under the assumption that |W (z)| 2 is linear in each bin of z. This assumption defines the horizontal positions of the data points in Figure 3-right, which are different from the positions in Figure 3-left.

Forward-backward asymmetry measurement
The forward-backward asymmetry A FB of the K πµµ decay is defined in terms of the angle θ Kµ between the K + and the µ − three-momenta in the µ + µ − rest frame, as where the numbers of events N are obtained after correction for the non-uniform acceptance in the (cos θ Kµ , z) plane (Figure 4-left). The resulting cos θ Kµ spectrum of the data events and the distribution expected from the Standard Model (SM) are displayed in Figure 4-right.
The asymmetry is measured to be and shows no significant dependence on z. The statistical precision is at the level of the upper limits on A FB predicted by the Minimal Supersymmetric Standard Model [28] and by the calculation of the two-photon intermediate state K + → π + γ * γ * → π + µ + µ − [29].

Systematic and external uncertainties
The individual contributions to the total uncertainties are discussed in the following and listed in Table 1.

Trigger efficiency
The trigger behaviour is emulated with a set of software algorithms applied to simulated events. The algorithms are tuned and validated on K 3π events. The L0 RICH, L0 CHOD, L0 MUV3 and L1 KTAG trigger efficiencies (equal to 99.8%, 98.2%, 98.9% and 99.8%, respectively) are found to be independent of the decay kinematics. Data and simulation efficiencies agree within 0.3%. The L1 STRAW trigger efficiency is 94.7% and varies as a function of decay kinematics within O(1%). Data and simulation efficiencies agree within 0.5%. The similarity of the MT and 2µMT trigger lines results in substantial cancellation of triggerrelated systematic effects. The residual systematic uncertainties are estimated by either disabling the software trigger emulators in simulation (in the case of the L0 RICH and L0 CHOD conditions), or replacing them with simplified emulators (L0 MUV3, L1 KTAG, L1 STRAW).

Reconstruction and particle identification
The similarity of the signal and normalization selections allows for significant cancellation of most systematic effects coming from reconstruction and particle identification efficiencies.
Systematic uncertainties arising from differences between event reconstruction efficiencies in data and simulation are dominated by the three-track event reconstruction in the STRAW spectrometer. A dedicated K 3π event selection, relying on a reconstructed kaon track in the GTK and two pion tracks in the STRAW, is used to measure the efficiency of reconstructing the third pion track. The average measured efficiency is 84% and depends on the decay kinematics. The observed differences of up to 2% between the efficiencies in data and simulation are considered in evaluating the systematic effects resulting from the STRAW track reconstruction efficiency.
The CHOD and MUV3 reconstruction efficiencies are above 99%, with no more than 0.6% difference between data and simulation.
The differences between data and simulation in the hadronic shower development and energy reconstruction in the LKr are another source of systematic uncertainty. No significant difference is observed in the efficiency of the muon identification. The efficiency of the pion identification measured on data is 99%. The agreement between data and simulation varies with pion momentum within 1%. Residual effects due to different K πµµ and K 3π pion kinematics are treated as systematic uncertainties.

Beam and accidental activity simulation
Systematic uncertainties stemming from the quality of the simulation of the beam momentum spectrum and intensity profile, and from the accuracy of the simulation of the halo muons accompanying the beam, are combined into a single systematic uncertainty. The selected normalization sample of K 3π events is used for the beam momentum and intensity studies. The halo muons are selected from out-of-time STRAW tracks that have associated signals in MUV3 and are not compatible with originating from decays in the FV.

Background
The number of background events is estimated using simulation to be 7.8 ± 5.6, where the error comes from the limited statistics of simulated background decays. The background arises mainly from the K 3π contribution with two π ± → µ ± ν decays in flight. More details on the methods employed in the K 3π background estimation can be found in [30].
Systematic uncertainties from the background contamination are estimated conservatively as differences between the results obtained with background neglected and background subtracted.

External uncertainties
External uncertainties in the measured quantities originate from the K 3π branching fraction [22], from the accuracy of the radiative corrections to the K πµµ decay, including the numerical approximation of dΓ 4-body /dz, and from the pion loop term parameters α + and β + [26].

Comparison with earlier measurements
A comparison of the present results with those from previous measurements by E787, E865, HyperCP and NA48/2 is shown in Figure 5, Table 2, and Table 3. Note that the NA48/2 measurement [11], until now the most precise, used a different K 3π branching fraction [31], and did not simulate the inner bremsstrahlung radiation of K 3π decays. Implementing these conditions in the NA62 analysis has minor impact on the results, which would change by δa + = −0.001, δb + = −0.002, and δB πµµ = +0.03 × 10 −8 .
Furthermore, the analysis by NA48/2 did not simulate inclusive radiative corrections and the 4-body radiative decay K + → π + µ + µ − γ in the K πµµ sample, but implemented only the softphoton Coulomb corrections for all pairs of the K πµµ decay products. Adopting this approach Form factor parameter ) µ µ π NA62 ( ) µ µ π NA62 ( ) µ µ π NA48/2 ( ee) π NA48/2 ( ee) (stat. only) π E865 ( 68% CL contours: Figure 5: Comparison with earlier measurements. Left: the K πµµ branching fraction, with the PDG [22] average shown as a shaded band. Right: combined statistical and systematic 68% CL contours in the (a + , b + ) plane for the muon and electron modes. The NA48/2 value of ρ(a + , b + ) for the electron mode is used to construct the E865 contour, as E865 [14] does not quote the correlation. The systematic uncertainties in a + and b + are not provided by E865 [14]. changes the NA62 results by δa + = −0.006, δb + = +0.034, and δB πµµ = −0.06×10 −8 , where the 0.7% relative change in the branching fraction comes from the increase of the signal acceptance measured with the 3-body simulated K πµµ sample including Coulomb corrections.

Summary
A sample of 27679 K πµµ candidates with negligible background contamination was collected by the NA62 experiment in 2017-2018. The size of the K πµµ data sample is the main factor limiting the precision of the present analysis. The K πµµ model-independent branching fraction is measured to be (9.15 ± 0.08) × 10 −8 , consistent with previous measurements and at least a factor of three more precise.
The form factor parameters in the framework of the Chiral Perturbation Theory at O(p 6 ) are measured as a + = −0.575 ± 0.013, b + = −0.722 ± 0.043. Values and statistical errors of parameters in any other form factor model can be obtained from the reconstructed values of the |W (z)| 2 function. The present measurement is the first to employ inclusive radiative corrections in the simulation of the signal channel. The form factor parameters are consistent with those measured by NA48/2, as well as with the results obtained in the electron mode by other experiments, suggesting agreement with lepton flavour universality in the K π decays.
The forward-backward asymmetry of the K πµµ decay is measured to be A FB = (0.0 ± 0.7) × 10 −2 , a factor of 2.6 improvement in the precision with respect to NA48/2. The experimental precision reaches the level of the upper limits on A FB predicted by the Minimal Supersymmetric Standard Model and by the calculation of the two-photon intermediate state K + → π + γ * γ * → π + µ + µ − .