Direct CP violation in charmless three-body decays of B ± mesons

Measurements of CP asymmetries in charmless three-body decays of B ± mesons are reported using proton-proton collision data collected by the LHCb detector, corresponding to an integrated luminosity of 5.9 fb − 1 . The previously observed CP asymmetry in B ± → π ± K + K − decays is confirmed, and CP asymmetries are observed with a significance of more than five standard deviations in the B ± → π ± π + π − and B ± → K ± K + K − decays, while the CP asymmetry of B ± → K ± π + π − decays is confirmed to be compatible with zero. The distributions of these asymmetries are also studied as a function of the three-body phase space and suggest contributions from rescattering and resonance interference processes. An indication of the presence of the decays B ± → π ± χ c 0 (1 P ) in both B ± → π ± π + π − and B ± → π ± K + K − decays is observed, as is CP violation involving these amplitudes.


I Introduction
The violation of charge-parity (CP ) symmetry has been observed in B + , B 0 and B 0 s decays, and yet our understanding of the dynamics involved in these asymmetries is incomplete. Some crucial questions remain to be understood: what is the dynamical origin of the weak phase associated to the decay amplitudes? Why has CP violation not been observed so far in baryon decays? What is the primary source of the strong phase difference in the observed CP asymmetries [1]?
The last question is directly related to the longstanding debate about the role of short-and long-distance contributions to the generation of the strong-phase differences, needed for direct CP violation to occur, and three-body decays of B mesons offer a way of finding the answer. A recent amplitude analysis found a large CP asymmetry related to the interference between the S -and P -wave contributions in B ± → π ± π + π − decays [2,3]. In a similar analysis, large CP violation involving ππ → KK rescattering was observed in B ± → π ± K + K − decays [4]. The patterns of localized CP asymmetries in the four charmless decay modes B ± → π ± π + π − , B ± → K ± K + K − , B ± → π ± K + K − and B ± → K ± π + π − , reported in Refs. [5,6], can be interpreted as originating from long-distance hadronic interactions [1].
Heavy meson three-body decays can be understood in terms of a superposition of intermediate states involving resonances. From this perspective, the phase-space-integrated CP asymmetries (A CP ) observed in charmless three-body decays of B ± decays are the sum of the asymmetries of each intermediate state.
In the B ± → π ± K + K − decay [6], the observed CP asymmetry would result predominantly from the ππ → KK rescattering amplitude, weighted by its fit fraction [4]. The same can be seen from the recent amplitude analysis of the B ± → π ± π + π − decay [2,3]. The sum of the CP asymmetries for the B ± → σπ ± and B ± → f 2 (1270)π ± decays weighted by their respective fractions is almost equal to the inclusive CP asymmetry observed for the B ± → π ± π + π − decay [6].
Another approach considers the global properties of the three-body amplitudes.
These two approaches are complementary and the role of each one must be determined from experimental results, such as the CP asymmetries and the branching fractions of charmless three-body B ± decays, recently presented by the LHCb collaboration [12].
In this article, CP violation is studied in four charmless B ± meson decays into three charged pseudoscalar particles: B ± → K ± π + π − , B ± → K ± K + K − , B ± → π ± K + K − and B ± → π ± π + π − . The analysis is based on proton-proton (pp) collision data at a center-ofmass energy of 13 TeV, corresponding to an integrated luminosity of 5.9 fb −1 , recorded with the LHCb detector between 2015 and 2018. Phase-space integrated CP asymmetries are measured, as well as CP asymmetries in specific regions of the Dalitz plots. The analysis uses the decay B ± → J/ψ(→ µ + µ − )K ± as a control channel, in order to take into account the production asymmetry of the B ± mesons.

II LHCb detector, dataset and simulation
The LHCb detector [13,14] is a single-arm forward spectrometer covering the pseudorapidity range 2 < η < 5, designed for the study of particles containing b or c quarks. The detector includes a high-precision tracking system consisting of a silicon-strip vertex detector surrounding the pp interaction region, a large-area silicon-strip detector located upstream of a dipole magnet with a bending power of about 4 Tm, and three stations of siliconstrip detectors and straw drift tubes placed downstream of the magnet. The tracking system provides a measurement of the momentum, p, of charged particles with a relative uncertainty that varies from 0.5% at low momentum (∼ 2 GeV/c) to 1.0% at 200 GeV/c. The minimum distance of a track to a primary pp collision vertex (PV), the impact parameter (IP), is measured with a resolution of (15 + 29/p T ) µm, where p T is the component of the momentum transverse to the beam, in GeV/c. Different types of charged hadrons are distinguished using information from two ring-imaging Cherenkov detectors. The corresponding efficiencies and misidentification rates are given in Ref. [15]. Photons, electrons and hadrons are identified by a calorimeter system consisting of scintillating-pad and preshower detectors, an electromagnetic and a hadronic calorimeter. Muons are identified by a system composed of alternating layers of iron and multiwire proportional chambers.
The online event selection is performed by a trigger, which consists of a hardware stage, based on information from the calorimeter system, followed by a software stage, which applies a full event reconstruction. At the hardware trigger stage, events are selected if a transverse energy, typically larger than 3.5 GeV, is deposited in the calorimeters by at least one hadron from the B ± → h ± h ′+ h ′− candidates. For the control channel, it is required that the hardware trigger decision must also be due to at least one muon with high p T . In all final states, events are also selected if the trigger decision is independent of the signal.
The software trigger requires a two-, three-or four-track vertex with a significant displacement from any PV. At least one charged particle must have a large p T and be inconsistent with originating from any PV. A multivariate algorithm is used for the identification of displaced vertices consistent with the decay of a b-hadron [16].
Samples of simulated events are used to model the effects of the detector efficiency and the selection requirements, to determine the fit models and to evaluate efficiencies. In the simulation, pp collisions are generated using Pythia 8 [17] with a specific LHCb configuration [18]. Decays of unstable particles are described by EvtGen [19], in which final-state radiation is generated using Photos [20]. The interaction of the generated particles with the detector, and its response, are implemented using the Geant4 toolkit [21,22] as described in Ref. [23].

III Candidate selection
The selection of the four B ± → h ± h ′+ h ′− samples is based on simulations and data-driven methods. The analyzed decay modes share the same topology and similar kinematics, leading to a common strategy for the selection. Requirements are made on the quality of the three-track vertices: the decay vertex must be detached from the PV, and the momentum vector of the reconstructed B ± candidate must point back to the PV. The final-state particles must have a significant IP with respect to the PV.
The combinatorial background is reduced by a multivariate analysis of a set of ten topological and kinematic variables, using boosted decision trees (BDT) [24] as classifiers. Samples of simulated events are used to describe the signal, whereas candidates satisfying m(h ± h ′ + h ′ − ) > 5.4 GeV/c 2 , are used to represent the combinatorial background. The BDT classifier is trained separately for each decay mode. A requirement on the BDT response for each mode is chosen to maximise the ratio ε sim / (S + B) data , where ε sim is the signal efficiency from simulation and (S + B) data is obtained by counting the candidates from data with invariant mass in the interval of ±40 MeV/c 2 centred on the known B + mass [25]. The most discriminating BDT variables are related to the B ± candidate distance of closest approach to the PV, pointing angle, decay vertex quality, and IP of each final-state particle.
Particle identification (PID) requirements are made on all three final-state particles, reducing the contamination from other B ± → h ± h ′+ h ′− decays to the per cent level. This background is due mainly to K − π misidentification. The analysis is restricted to particles with momentum below 100 GeV/c, since above this limit the discrimination power of the PID system is limited. In order to avoid the contamination from semileptonic decays, a muon veto is placed on all decay products by applying requirements on track reconstruction and trigger decisions. Finally, a method based on Artificial Neural Network is employed to remove electrons [26].
More stringent PID requirements are applied to the B ± → π ± π + π − and B ± → π ± K + K − decays. These decays have lower branching fractions compared to the B ± → K ± π + π − and B ± → K ± K + K − decays and are therefore more severely affected by the cross-feed from other There are significant backgrounds from two-body decays of B ± mesons into a D 0 or D 0 meson and a pseudoscalar particle, followed by the two-body decay of the D 0 meson into pions and kaons. This background is removed by excluding candidates for which the two-body invariant masses m(π + π − ), m(K ± π ∓ ) or m(K + K − ) are within ±30 MeV/c 2 of the known D 0 mass [25].
After the final selection, approximately 0.5% of the events contain more than one signal candidate. One candidate is chosen randomly in these cases.
The decay B ± → J/ψK ± , with subsequent decay J/ψ → µ + µ − is used to measure the B ± production asymmetry. The sample is selected with the same requirements as the B ± → K ± π + π − decay, with three exceptions: the particle identification, in which the pion PID is substituted by that of muon; the use of events for which the hardware trigger decision is due to a muon with high p T ; the requirement of the J/ψ candidate mass within the interval 3.05 < m(µ + µ − ) < 3.15 GeV/c 2 .

IV Invariant mass fits
The yield and raw asymmetry of each decay channel are extracted from simultaneous unbinned extended maximum-likelihood fits to the B − and B + invariant mass distributions, as shown in Fig. 1. The signal components in each of the four channels are parameterized by a sum of a Gaussian and two Crystal Ball functions [27] with a common mean and different widths and tails on either side of the peak, accounting for asymmetric effects such as quark final-state radiation. The combinatorial backgrounds is described by an exponential 5.4 5.6 3 10 × 5.4 5.6 3 10 ×   Figure 1: function. The slopes and charge asymmetry of the combinatorial background are free parameters in the fit, the latter being compatible with zero when systematic uncertainties, production and detection asymmetries are taken into account. The background due to partially reconstructed four-body B decays, in which one particle is not reconstructed, is parameterized by an ARGUS function [28] convolved with a Gaussian resolution function. A contribution from four-body B 0 s decays is also present in the B ± → π ± K + K − sample. The dominant peaking backgrounds in the signal regions arise from other three-body B decays with one or more misidentified final-state particles. For B ± → K ± π + π − decays, the largest contribution comes from B ± → π ± π + π − with the pion misidentified as kaon, followed by the partially reconstructed decay B ± → η ′ (ρ 0 γ)K ± , in which the photon is undetected. Conversely, the only peaking contribution to the B ± → π ± π + π − spectrum comes from B ± → K ± π + π − decays. The peaking background to B ± → K ± K + K − is the result of single or double misidentification of pions as kaons from B ± → π ± K + K − and B ± → K ± π + π − decays. The B ± → π ± K + K − decay has the largest number of significant peaking background contributions due to its lower branching fraction and both K → π and π → K misidentification possibilities, which is why the three other channels can contribute to its spectrum.
The shapes and fractions of peaking backgrounds relative to the signal yield are obtained from simulation of the relevant decay modes and fixed in the fits, whereas the yields of the partially reconstructed background components vary freely in the fit. The resulting signal yields and raw asymmetries (A raw ) are listed in Table I. The data samples are larger than those presented in Ref. [6] due to both an increase in the integrated luminosity and the production cross-section, and the use of a more efficient selection. Table I: Total signal yields, raw asymmetries and average efficiency ratio R of charmless threebody B ± decays for the inclusive data set.

Decay mode
Total yield Neither the selection efficiencies nor the observed raw asymmetries are uniform across the Dalitz plot, and efficiencies may differ for B + and B − . Therefore, the efficiency of selecting a signal decay is parameterized in the two-dimensional square Dalitz plot separately for B + and B − . Such a transformation improves the resolution in the regions with high sensitivity by expanding the corners and the borders of the nominal Dalitz plot relative to the sparsely populated central region. These maps are primarily obtained using simulated decays; however, effects arising from the hardware trigger and PID efficiency are determined using data calibration samples. The efficiency shows a variation from 0.9% to 1.7% for The efficiency correction applied to the integrated raw asymmetries is determined from the ratio between the B − and B + average efficiencies, as listed in Table I. These averages reflect the data distributions across the phase space by weighting the efficiency maps by the data events. The phase-space integrated CP asymmetries A CP are obtained by correcting the raw asymmetry for selection efficiency effects to get the efficiency-corrected raw asymmetries A corr raw , and for the production asymmetry A P . The efficiency-corrected raw asymmetry is defined as where N corr B − and N corr B + are the number of B − and B + signal events obtained from the invariant mass fit, N B − and N B + , corrected for the efficiency. They are related to the asymmetries by Substituting Eqs. (2) and (3) into Eq. (1), the physical CP asymmetry can be expressed in terms of A corr raw and A P : To obtain A P , the decay B ± → J/ψK ± is used as a control channel and its raw asymmetry is obtained from a fit to the invariant mass distribution. This decay has larger signal yield and lower background contributions compared to B ± → h ± h ′+ h ′− decays. Its background is parameterized as the sum of a combinatorial component, represented by a single exponential, and partially-reconstructed B → J/ψK * (892) ±,0 decay, described by an ARGUS function. The B ± → J/ψK ± model is parameterized with the same function used for the B ± → h ± h ′+ h ′− channels. Systematic uncertainties of the raw asymmetry are obtained by varying the signal fit model, leaving the background asymmetry to vary in the fit, and looking at variations from different trigger samples of the data. The total systematic uncertainty is taken as the sum in quadrature of the individual uncertainties. The raw asymmetry of the control channel is measured to be A raw (B ± → J/ψK ± ) = −0.0118 ± 0.0008 +0.0007 −0.0008 , where the first uncertainty is statistical and the second systematic.
To obtain the B ± production asymmetry, this raw asymmetry is corrected by its efficiency ratio, calculated using a sample of simulated events produced without any B ± production asymmetry, to obtain A corr raw as before, and the world average value of the B ± → J/ψK ± CP asymmetry, 0.0018 ± 0.0030 [25], is subtracted The measured B meson production asymmetry is A P = −0.0070 ± 0.0008 +0.0007 −0.0008 ± 0.0030, where the last uncertainty is due to the CP asymmetry of B ± → J/ψK ± decays [25].
Finally, the CP asymmetries of the four B ± → h ± h ′+ h ′− modes are measured to be where the statistical uncertainties are obtained from propagation of Eq. (4), assuming no correlation term.

VI Systematic uncertainties and results
Several sources of systematic uncertainties are considered and can be broadly divided into three groups: potential mismodelling of the invariant mass distributions, phase-space efficiency corrections and knowledge of the B ± production asymmetry. The systematic uncertainties due to the mass fit models are quantified by taking the difference in the CP asymmetry resulting from variations of the model. The alternative fits have good quality and describe the data accurately. To estimate the uncertainty due to the choice of the signal mass function, the initial model is replaced by an alternative empirical distribution [27]. The contribution associated with the peaking background fractions reflects the uncertainties in the expected yields determined from simulation and it is evaluated by varying the fractions within their statistical uncertainties. In addition, the systematic uncertainty associated to the fact that the peaking background asymmetry is fixed to zero is estimated by setting it to the value obtained in the previous analysis [6], within the corresponding uncertainties. The uncertainty due to the choice of an exponential function to model the combinatorial component is estimated by repeating the fit using a second order polynomial function. Table II: Individual components and total systematic uncertainties related to the mass fit model, efficiency corrections for the integrated phase space and B ± production asymmetry. The total is the sum in quadrature of these components. The systematic uncertainty related to the efficiency correction procedure consists of two parts: the statistical uncertainty in the detection efficiency due to the finite size of the simulated samples, and the uncertainty due to the finite size of the bins, which is evaluated by varying the binning used in the efficiency correction. The former term dominates for all decays.

Source of uncertainty
The systematic uncertainty due to the B ± production asymmetry measurement using the B ± → J/ψK ± control channel is obtained by summing in quadrature the experimental statistical and systematic uncertainties of A P .
The systematic uncertainties are summarized in Table II, where the total is the sum in quadrature of the individual contributions. The results for the integrated CP asymmetries are A CP (B ± → K ± π + π − ) = +0.011 ± 0.002 ± 0.003 ± 0.003, A CP (B ± → K ± K + K − ) = −0.037 ± 0.002 ± 0.002 ± 0.003, A CP (B ± → π ± π + π − ) = +0.080 ± 0.004 ± 0.003 ± 0.003, A CP (B ± → π ± K + K − ) = −0.114 ± 0.007 ± 0.003 ± 0.003, where the first uncertainty is statistical, the second is systematic and the third is due to the limited knowledge of the CP asymmetry of the B ± → J/ψK ± control channel [25]. The significance of the CP asymmetries is computed by dividing the central values by the sum in quadrature of the uncertainties, yielding 8.5 standard deviations (σ) for B ± → K ± K + K − decays, 14.1σ for B ± → π ± π + π − decays and 13.6σ for B ± → π ± K + K − decays. This is the first observation of CP asymmetries in these channels. For the B ± → K ± π + π − decays, the significance of the CP asymmetry is 2.4σ, consistent with CP conservation.

VII Localised CP asymmetries
The two-dimensional phase space referred to as the Dalitz plot, formed by the invariant masses squared of two of the three possible particle pairs, reflects directly the dynamics of the decay and allows the study of the different resonant and nonresonant components to inspect CP asymmetry effects. The Dalitz plot distributions of signal yield, obtained with the sPlot technique [29], for B ± → π ± π + π − , B ± → K ± π + π − , B ± → π ± K + K − and The colour scale indicates the number of events. Fig. 2. In the symmetric channels, the phase space distribution and its projections are presented with the two axes being the squares of the low-mass m low and high-mass m high combinations of the opposite-sign particle pairs, for visualization purposes.
Most of the candidates are concentrated in the low-mass regions, as expected for charmless decays dominated by resonant contributions. The gap from the vetoed potential J/ψ contributions is visible in the B ± → K ± π + π − channel, as well as the gaps from D 0 regions excluded in all modes.
In order to visualise localized asymmetries, the A CP in bins of the phase space [30] is constructed. Adaptive binning is employed, such that the efficiency-corrected signal yield, also obtained with the sPlot technique, is approximately equal in all bins. There is no specific rule for choosing the best binning except for requiring a minimum bin occupancy. Figure 3 reveals a rich pattern of large and localized asymmetries which result from interference between the contributions, as well as possible ππ → KK rescattering, as was observed in the amplitude analyses of Refs. [2][3][4] and elastic scattering experiments [31,32].
Different regions of the Dalitz plots in Fig. 3 are defined to study the localized CP asymmetries. The rescattering region [31] is defined in the Dalitz plot in the two-kaon or two-pion invariant mass range 1.1-2.25 GeV 2 /c 4 for B ± → K ± K + K − due to the presence of the ϕ(1020) resonance and 1.0-2.25 GeV 2 /c 4 for the other three channels, as listed in Ta-  Figure 3: Asymmetry distribution in bins of the Dalitz plot for (a) B ± → π ± π + π − , using 400 bins with approximately 229 events/bin; (b) B ± → K ± π + π − , using 1728 bins with approximately 276 events/bin; (c) B ± → π ± K + K − , using 256 bins with approximately 127 events/bin; and (d) B ± → K ± K + K − , using 729 bins with approximately 461 events/bin. ble III. The data are studied in intervals of the other two-body invariant mass. The vertical axis projection depends on the decay channel, with two regions defined for B ± → π ± π + π − , B ± → K ± π + π − and B ± → K ± K + K − , and one for B ± → π ± K + K − . There are also two other regions covering higher masses for B ± → π ± π + π − and B ± → π ± K + K − decays. All regions for the different decay channels are defined in Table III. The invariant-mass fits used to compute A CP in each region follow the same fit procedure as described in Sec. IV, where the efficiency ratios are computed separately per region. Most of the rescattering regions give CP asymmetry values in excess of at least five Gaussian standard deviations. For B ± → π ± π + π − decays, the rescattering region when analysed in the m 2 (π + π − ) high projection presents a CP asymmetry that is positive in region 1 and negative in region 2, as shown in Fig. 4.
For B ± → K ± π + π − decays, the two regions in the m 2 (K + π − ) axis reveal a flip of the asymmetry sign as shown in Fig. 5. Figure 6 illustrates the m 2 (K + π − ) axis projection for the rescattering region for B ± → π ± K + K − decays, where a nearly constant CP asymmetry is observed. The B ± → K ± K + K − decays also reveal CP asymmetries in the m 2 (K + K − ) high projections in the rescattering region, as shown in Fig. 7. Here, another change in the asymmetry sign is observed, in an opposite direction with respect to B ± → K ± π + π − decays. Region 2 of the B ± → K ± K + K − decay is the only rescattering projection that shows small asymmetry.
The high-mass regions of B ± → π ± π + π − and B ± → π ± K + K − as defined in Table III  Table III: Definitions of phase-space regions in units of GeV 2 /c 4 . In the B ± → K ± K + K − channel, the ϕ(1020) vector meson is excluded by requiring m 2 (K + K − ) low > 1.1 GeV 2 /c 4 .

VIII Partial width difference
A relation between CP violation in B ± → K ± π + π − and B ± → π ± K + K − , as well as that between B ± → K ± K + K − and B ± → π ± π + π − , based on U-spin symmetry [7][8][9], is investigated. This symmetry is an SU(2) subgroup of the flavour SU(3) group, under which the (d,s) pairs of quarks form a doublet, similar to (u,d) in isospin symmetry. U-spin symmetry predicts the relationship between the differences of the partial decay widths for the decays B ± → K ± π + π − and B ± → π ± K + K − to be A similar prediction is made between the symmetric decay channels B ± → π ± π + π − and These relations state that the difference in the number of events associated to CP violation must be the same, independently of the branching fraction of each decay. The above relations can be rewritten in terms of the measured integrated asymmetries, branching fractions and the B meson lifetime. For the B ± → π ± K + K − decay the result is where τ (B ± ) is the B ± lifetime and B is the average branching fraction. Using the measured values of the integrated asymmetry and branching fractions [12], the U-spin relationship between B ± → K ± π + π − and B ± → π ± K + K − is A CP (B ± → π ± K + K − )B(B ± → π ± K + K − ) A CP (B ± → K ± π + π − )B(B ± → K ± π + π − ) = −0.92 ± 0.18 .
The results presented here reveal a central value with an opposite sign with respect to the U-spin prediction of −1. However, one may not forget that the U-spin symmetry prediction occurs in a scenario with momentum-independent amplitudes, in which no hadronic final-state interaction where two pairs of channels are coupled through the ππ ↔ KK interaction is considered. When those contributions are included, one may find ratios that differ from −1 without necessarily claiming the observation of U-spin symmetry breaking. All these implications are further discussed in Refs. [10,11].

IX Summary and conclusions
In summary, the inclusive CP asymmetries of the B ± → K ± π + π − , B ± → π ± K + K − , B ± → K ± K + K − and B ± → π ± π + π − charmless three-body decays are measured with data corresponding to 5.9 fb −1 of integrated luminosity, collected by the LHCb detector between 2015 and 2018. Significant inclusive CP asymmetries are found for the latter three B decay channels, the last two observed for the first time.
The asymmetry in localized regions of phase space shows that the CP asymmetries are not uniformly distributed in the phase space, with positive and negative A CP appearing in the same charged B decay channel. The results of the previous LHCb analysis [6] are confirmed with high localized CP asymmetries in several regions of the Dalitz plot. It is also confirmed that, in the ππ → KK rescattering region, between 1 and 2.25 GeV 2 /c 4 in both π + π − and K + K − invariant masses, significant CP violation is present in the four analysed channels. For the B ± → π ± π + π − , B ± → K ± π + π − and B ± → K ± K + K − decays, an A CP sign change is observed across the phase space.
Additionally, an indication of the χ c0 (1P ) resonance is seen in the high π + π − and K + K − invariant mass regions of B ± → π ± π + π − and B ± → π ± K + K − decays for the first time. There is an indication of CP violation involving this resonance, as predicted several years ago in Refs. [33,34]. The combination of the present A CP results with the recent LHCb relative branching fraction results for these decays, shows good agreement with the U-spin symmetry proposed in Refs. [7][8][9].
Research Funds for the Central Universities, and Sci. & Tech. Program of Guangzhou (China); GVA, XuntaGal and GENCAT (Spain); SRC (Sweden); the Leverhulme Trust, the Royal Society and UKRI (United Kingdom).