Observation of $CP$ violation in two-body $B_{(\mathrm{s})}^{\,0}$-meson decays to charged pions and kaons

The time-dependent $CP$ asymmetries of $B^0\to \pi^+\pi^-$ and $B^0_\mathrm{s}\!\to K^+K^-$ decays are measured using a data sample of $pp$ collisions corresponding to an integrated luminosity of 1.9 fb$^{-1}$, collected with the LHCb detector at a centre-of-mass energy of 13 TeV. The results are \begin{eqnarray*} C_{\pi\pi}&=&-0.311 \pm 0.045 \pm 0.015, \\ S_{\pi\pi}&=&-0.706 \pm 0.042 \pm 0.013, \\ C_{KK}&=&\phantom{-}0.164 \pm 0.034 \pm 0.014, \\ S_{KK}&=&\phantom{-}0.123 \pm 0.034 \pm 0.015, \\ \cal {A}^{\mathrm{\Delta\Gamma}}_{\mathit{KK}}&=&-0.83\phantom{0} \pm 0.05\phantom{0} \pm 0.09, \end{eqnarray*} where the first uncertainties are statistical and the second systematic. The same data sample is used to measure the time-integrated $CP$ asymmetries of $B^0\to K^+\pi^-$ and $B^0_\mathrm{s}\to K^-\pi^+$ decays and the results are \begin{eqnarray*} A_{CP}^{B^0}&=&-0.0824 \pm 0.0033 \pm 0.0033, \\ A_{CP}^{B^0_\mathrm{s}}&=&\phantom{-}0.236\phantom{0} \pm 0.013\phantom{0} \pm 0.011. \end{eqnarray*} All results are consistent with earlier measurements. A combination of LHCb measurements provides the first observation of time-dependent $CP$ violation in $B^0_\mathrm{s}$ decays.


Introduction
Charge-parity (CP ) asymmetries of charmless B 0 (s) -meson decays to two-body charged final states are important inputs to the validation of the Cabibbo-Kobayashi-Maskawa (CKM) mechanism [1,2], which models CP violation in charged-current quark transitions. Deviations from Standard Model (SM) predictions may reveal the presence of phenomena not included in the SM, manifested as modifications to the amplitudes of these decays. [3][4][5][6][7][8][9]. The CP asymmetry in the B 0 → π + π − decay is a fundamental input to the isospin analysis of B → ππ decays that allows the determination of the CKM angle α [10][11][12]. The analysis can be extended by exploiting the approximate U-spin symmetry [13] that relates the hadronic parameters entering the decay amplitudes of the B 0 → π + π − and B 0 s → K + K − decays. 1 It has been shown that, by incorporating the CP asymmetry and branching fraction of the B 0 s → K + K − decay into the standard isospin analysis, stringent constraints on the CKM angle γ and on the CP -violating phase −2β s can be set, even when allowing for U-spin breaking effects [14,15]. Furthermore, a substantial reduction of uncertainties on the determination of −2β s can be achieved by combining the CP asymmetries of the B 0 → π + π − and B 0 s → K + K − decays with information provided by the semileptonic decays B 0 → π − + ν and B 0 s → K − + ν [16,17]. The CP asymmetries and branching fractions of the B 0 → K + π − and B 0 s → K − π + provides the test of the SM, assuming U-spin symmetry, proposed in Ref. [7]. The CP asymmetry of the B 0 → K + π − decay is also a key input to the long-standing B → Kπ puzzle [18][19][20]. Strategies have been proposed to combine information from several decays of the B → ππ and B → Kπ systems in order to investigate the presence of physics beyond the SM [21][22][23].
This paper presents measurements of time-dependent CP asymmetries in B 0 → π + π − and B 0 s → K + K − decays and of time-integrated CP asymmetries in B 0 → K + π − and B 0 s → K − π + decays. The analysis is based on a data sample of pp collisions corresponding to an integrated luminosity of 1.9 fb −1 , collected with the LHCb detector at a centre-ofmass energy of 13 TeV during 2015 and 2016. These results are combined with previous LHCb results, published in Ref. [24], based on a sample corresponding to 3.0 fb −1 , collected at 7 and 8 TeV in the Run 1 data taking.
In decays of B 0 (s) mesons to a final state f , where f is a CP eigenstate (f = f ), CP violation originates from the interference between the decay and B 0 (s) -B 0 (s) mixing. The latter can be modelled by an effective Hamiltonian whose mass eigenstates are linear combinations of the two flavour eigenstates, p|B 0 (s) ± q|B 0 (s) , where p and q are complex parameters, normalised such that |p| 2 + |q| 2 = 1. The CP asymmetry as a function of decay time for B 0 (s) → f decays is given by where ∆m d(s) and ∆Γ d(s) are the mass and width differences of the mass eigenstates of the B 0 (s) system. In accordance with current experimental knowledge, the value of ∆Γ d is assumed to be negligible. The quantities C f , S f and A ∆Γ f are defined as where λ f is given by and A f (A f ) is the decay amplitude for the B 0 (s) (B 0 (s) ) → f transition. As current experimental determinations [25][26][27] confirm the SM expectation [28,29] of negligible CP violation in the B 0 (s) -B 0 (s) mixing (implying |q/p| = 1), a nonzero value of C f and S f indicates the presence of CP violation in the decay and in the interference between mixing and decay, respectively. The quantities C f , S f and A ∆Γ f are related through the unitary condition (C f ) 2 + (S f ) 2 + A ∆Γ f 2 = 1. This constraint is not imposed in this analysis and is instead used as a cross-check of the consistency of the results. Previous determinations of C ππ and S ππ were performed by BaBar [30], Belle [31] and LHCb [24] experiments, while only LHCb has measured C KK , S KK and A ∆Γ KK [24]. The time-integrated CP asymmetry for a B 0 (s) decay to a flavour-specific final state f , such as B 0 → K + π − and B 0 s → K − π + , is defined as Measurements of A CP for the B 0 → K + π − decay (A B 0 CP ) were carried out by BaBar [30], Belle [32], CDF [33] and LHCb [24], while A CP for the B 0 s → K − π + decay (A B 0 s CP ) was measured only by CDF [33] and LHCb [24]. This paper is organised as follows. The LHCb detector, its trigger system and the simulation process are briefly introduced in Sec. 2, while the sample selection is described in Sec. 3. The CP asymmetries are determined by means of unbinned maximum-likelihood fits to the invariant-mass and decay-time distributions of B 0 (s) candidates reconstructed in the π + π − , K + K − and K ± π ∓ final states. In order to measure the time-dependent CP asymmetries, it is necessary to determine the flavour of the B 0 (s) meson at its production. In addition, a precise determination of the B 0 (s) decay time is important, in particular for the B 0 s meson, due to its fast oscillation frequency. The flavour-tagging algorithms and their calibration are presented in Sec. 4, while the determination of the decay-time resolution is discussed in Sec. 5. The models used in the fits are described in Sec. 6. Two measurements of the CP -violating parameters for the B 0 → π + π − and B 0 s → K + K − decays are performed with different experimental techniques. The first method, referred to as the simultaneous method, fits all the signal decays simultaneously and uses a fit model similar to that described in Ref. [24]. The second method, referred to as the per-candidate method, describes the selection efficiency as a function of the decay time of the B 0 (s) meson on a per-candidate basis using the swimming technique [34][35][36][37]. The determination of the detection asymmetry between the B 0 → K + π − and B 0 s → K − π + decays and their charge-conjugate final states, necessary to measure A CP , is discussed in Sec. 7. The results are given in Sec. 8 and the assessment of systematic uncertainties is presented in Sec. 9. The statistical and systematic uncertainties on the simultaneous method are found to be, in general, smaller than those for the per-candidate method. The results from the simultaneous method are therefore given as the main results of this paper. The final results and their combination with previous LHCb measurements from Ref. [24] are presented in Sec. 10, while considerations on the combined measurements are reported in Sec. 11.

Detector, trigger and simulation
The LHCb detector [38,39] is a single-arm forward spectrometer covering the pseudorapidity in the range between 2 and 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 silicon-strip detectors and straw drift tubes placed downstream of the magnet [40,41]. The tracking system provides a measurement of momentum, p, of charged particles with a relative uncertainty that varies from 0.5% at low momentum 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 (RICH) detectors [42]. Photons, electrons and hadrons are identified by a calorimeter system consisting of scintillating-pad and preshower detectors, an electromagnetic calorimeter 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 [43], which consists of a hardware stage, based on information from the calorimeter and muon systems, followed by a software stage, which applies a full event reconstruction.
At the hardware trigger stage, events are required to have a muon with high p T , or a hadron, photon or electron with high transverse energy in the calorimeters. For hadrons, the transverse energy threshold is 3.5 GeV. The software trigger requires the presence in the event of at least one charged particle with p T > 1.6 GeV/c and inconsistent with originating from any PV. The tracks identified at this stage are used by a trigger selection dedicated for two-body b-hadron decays. The selection algorithm imposes requirements on the quality of the reconstructed tracks, their p T and minimum χ 2 IP with respect to every PV in the event, where the χ 2 IP is defined as the difference in the vertex-fit χ 2 of a given PV reconstructed with and without the track under consideration. Pairs of oppositely charged tracks must have a small distance of closest approach and a large scalar sum of their p T in order to be eligible to form a B 0 (s) candidate. Finally, the B 0 (s) candidates are required to pass criteria based on their p T , χ 2 IP , flight distance with respect to their associated PV, and angle between the direction of the B 0 (s) candidate momentum and the direction defined by its decay vertex and associated PV. Candidates are associated with the PV that is most consistent with their flight direction.
Simulation is used to study the discrimination between signal and background candidates, and to assess differences between signal and calibration decays. The pp collisions are generated using Pythia [44,45] with a specific LHCb configuration [46]. Decays of hadronic particles are described by EvtGen [47], in which final-state radiation is generated using Photos [48]. The interaction of the generated particles with the detector, and its response, are implemented using the Geant4 toolkit [49] as described in Ref. [50].

Selection
The B 0 (s) candidates selected by the dedicated software trigger are further filtered, requiring that either the decay products or particles from the rest of the event are responsible for the positive decision of the hadronic hardware trigger. Candidates are then classified into mutually exclusive samples of different final states (π + π − , K + K − and K ± π ∓ ) using particle identification (PID) information. Finally, a boosted decision tree (BDT) algorithm [51,52] is used to separate signal candidates from combinatorial background candidates for each of the final states.
Four types of background contributions are considered: two-body b-hadron decays with misidentified pions, kaons or protons in the final state (cross-feed background); pairs of randomly associated and oppositely charged tracks (combinatorial background); pairs of oppositely charged tracks from partially reconstructed three-body decays of b hadrons (three-body background); B 0 (s) mesons produced in B + c decays rather than at a PV, whose measured decay time is biased due to the finite lifetime of the B + c meson. Given the small production rate of B + c mesons [53], this background contribution is neglected in the analysis and a systematic uncertainty is assessed in Sec. 9. Since the three-body background candidates give rise to B 0 (s) candidates with invariant-mass values well separated from the mass peak, the candidate selection is customised to reject mainly the cross-feed and combinatorial background candidates, as they affect the invariant-mass region around the B 0 and B 0 s nominal masses. The requirements imposed on the PID variables, used to identify the π + π − and K + K − samples, are optimised using pseudoexperiments that take into account the different background contributions. First the PID efficiencies and misidentification probabilities for kaons and pions are determined, for different requirements, using samples of D * + → D 0 (→ K − π + )π + decays [54] and are used to estimate the cross-feed background yields in each of the final states. The results of the PID calibration and the fitting model described in Sec. 6 are used to generate pseudoexperiments that are fitted with the same model. The results of the fits are used to find the configuration of PID requirements giving the best trade-off between the statistical sensitivity to the CP -violation parameters of the B 0 → π + π − and B 0 s → K + K − decays and the systematic effects due to large contributions of cross-feed background candidates. The PID selection used to identify the K ± π ∓ samples is, instead, optimised to reduce the amount of the B 0 → π + π − and B 0 s → K + K − cross-feed background yields to approximately 10% of the B 0 s → K − π + yield. The BDT algorithm exploits the following properties of the B 0 (s) decay products: the p T of the two tracks; the χ 2 IP of each track with respect to their associated PV; the distance of closest approach between the two tracks, and the quality of their common vertex. The BDT classifier also uses properties of the reconstructed B 0 (s) candidate, particularly the p T , the χ 2 IP and the χ 2 of the flight distance with respect to the associated PV. Separate BDT algorithms are trained and optimised for the selection of the B 0 → π + π − and the B 0 s → K + K − decays. Simulated events of the two decay modes are used to model the signal candidates, while data from their high-mass sidebands (from 5.6 GeV/c 2 to 6.2 GeV/c 2 ) are used to model the combinatorial background candidates. The optimal threshold on the response of the BDT algorithm is chosen to maximise S/ √ S + B, where S and B represent the estimated yield of signal and combinatorial background candidates within ±60 MeV/c 2 (corresponding to about ±3 times the invariant-mass resolution) around the known B 0 (s) mass. The K ± π ∓ samples are selected using the BDT classifier optimised for the B 0 → π + π − decay. 2 Multiple candidates are present in less than 0.06% of the events satisfying the offline selections. Only one candidate is accepted at random from each event.
The optimisation of the selection criteria preferentially rejects short-lived candidates over longer lived ones. This introduces a distorted decay-time efficiency that must be corrected for. The selection criteria present in the analysis that produce this efficiency are the requirements on the χ 2 IP of all particles, the χ 2 of the B 0 (s) flight distance, the direction defined by its decay vertex and associated PV, and the outputs of the BDT algorithm. In addition, there are also decay-time biasing selection criteria due to the geometry of the detector. These are the limit on the radial flight distance of the B 0 (s) , which is required to avoid secondary interactions with the vertex detector material, and the minimal number of the vertex-detector sensors required to have track hits, which is imposed by the software triggers. The bias introduced by the radial flight distance is only present in the per-candidate method.

Flavour tagging
Tagging of the initial flavour of the B 0 (s) meson plays a crucial role in measuring the timedependent CP asymmetries of decays to CP eigenstates, since the sensitivity to the C f and S f coefficients, defined in Eq. (1), is related to the tagging performance. The flavour of the B 0 (s) candidates is inferred by two classes of the flavour-tagging algorithms called opposite-side (OS) and same-side (SS) taggers. The OS taggers [55] exploit the fact that in pp collisions beauty quarks are almost exclusively produced in bb pairs. Thus the flavour of the decaying signal B 0 (s) meson can be determined by looking at the decay products of the other b hadron in the event, for example, the charge of the lepton originating from semileptonic decays, the charge of the kaon from the b → c → s transition, or the charge of a charm hadron. An additional OS tagger is based on the inclusive reconstruction of the opposite b-hadron decay vertex by computing the p T -weighted average of the charges of all tracks associated to that vertex. The SS taggers are based on the identification of the particles produced in the hadronisation of the signal beauty quarks. In contrast to the OS taggers, which to a very good approximation act equally on B 0 and B 0 s mesons, SS taggers are specific to the light quark of the B 0 (s) meson under study. Additional d (d) or s (s) quarks produced in association with a B 0 (B 0 ) or a B 0 s (B 0 s ) meson, respectively, can form charged pions and protons, in the down-quark case, or charged kaons, in the strange-quark case. The so-called SSπ and SSp taggers [56] are used to determine the initial flavour of B 0 mesons, while the SSK tagger [57] is used for B 0 s mesons. For each tagger, the probability of misidentifying the flavour of the B 0 (s) meson at production, the mistag probability, η, is estimated by means of a multivariate classifier, and is defined in the range 0 ≤ η ≤ 0.5. The flavour-tagging performance of each tagger can be quantified by means of the tagging power, defined as where ξ i and η i are the tagging decision and the probability of misidentifying the flavour of the i-th out of N B 0 (s) candidates, respectively. The tagging decision ξ i takes the value of +1 when the candidate is tagged as B 0 (s) , −1 when the candidate is tagged as B 0 (s) , and zero for untagged candidates. Multivariate algorithms are used to determine the values of η for the OS and SS taggers, denoted as η OS and η SS . These are trained using specific B-meson decay channels and selections. The differences between the training samples and the selected signal B 0 (s) candidates can lead to an imperfect determination of the mistag probability. Hence, a more accurate estimate, denoted as ω hereafter, is obtained by means of a calibration procedure that takes into account the specific kinematics of selected signal B 0 (s) mesons. The relation between η and ω is calibrated using B + → D 0 π + , B 0 s → D − s π + and B 0 → D − π + decays for the OS, SSK, and SSπ and SSp taggers, respectively. The flavour for the B + meson is tagged by the charge of the pion in the final state. For the B 0 and B 0 s modes, which decay into flavour-specific final states, the amplitude of the tagged time-dependent asymmetry is proportional to 1 − 2ω. When the response of more than one OS tagger is available per candidate, the different decisions and associated calibrated mistag probabilities are combined into a unique decision ξ OS and a single η OS . A similar combination is also performed between the SSπ and SSp taggers to create a combined same-side tagger, SS c , where a combined tagging decision ξ SSc and mistag probability η SSc is evaluated, as discussed in App. A.2.
In the simultaneous method, the OS and SSc combinations are recalibrated in the final fit, discussed in Sec. 6, using the B 0 → K + π − decays in order to correct for possible correlations between the individual algorithms not taken into account in the combination procedure. For the SSK case, since the small yield of B 0 s → K − π + decays is insufficient for a reliable recalibration, the original calibration is kept and a systematic uncertainty is assigned. In the per-candidate method, the OS and SS combinations are further combined into a unique tagging decision and mistag probability using the calibrations determined by the simultaneous method. This combination is again recalibrated with the calibration samples. The description of the implementation of the flavour tagging into the fit models is presented in Sec. 6.

Decay-time resolution
The decay-time resolution is modelled with a Gaussian function, whose mean and width are calibrated with a sample of J/ψ → µ + µ − decays produced directly in pp collisions. The background contribution in the J/ψ sample is subtracted using the sPlot technique [58] with the dimuon invariant mass acting as a discriminating variable. The background-subtracted sample is separated in intervals of decay-time uncertainty, δ t , which is determined for each candidate from the kinematic fit used to measure the decay time. The decay-time distribution in each bin of δ t is fitted with a model comprising three Gaussian functions with shared mean and independent widths. According to Ref. [59] the parameters obtained from the fits are combined into an effective resolution, σ eff , such that a single-Gaussian resolution model with width σ eff gives the same dilution effect on the amplitude of the time-dependent asymmetry as the triple-Gaussian model. The value of σ eff is calibrated assuming all the signal decays have the same mixing frequency as the B 0 s meson. This assumption does not impact the analysis for B 0 mesons, since for them the effect of the decay-time resolution is negligible. Figure 1 shows the dependence of σ eff on δ t and is found to be well modelled with a linear function with an intercept q 0 and slope q 1 . The fit is repeated for different numbers of bins of δ t , and the obtained mean values of the slope and intercept are found to be 0.94 ± 0.02 and 1.64 ± 1.09 fs, respectively. Differences in the decay-time resolution between J/ψ → µ + µ − and two-body b-hadron decays are studied using samples of fully simulated where σ K + K − eff (δ t ) and σ µ + µ − eff (δ t ) are the effective resolution widths for the simulated B 0 s → K + K − and J/ψ → µ + µ − decays, respectively. For the per-candidate method, the calibrated resolution in Eq. 6 is applied to each candidate in the fit to the B 0 s →K + K − decay-time spectrum. 3 For the simultaneous method, the decay-time resolution is not used on a per-candidate basis, but an average model is used instead. The consequence of using the average model is a small loss in the statistical precision for C KK and S KK , corresponding to a relative 1% difference on the final uncertainties, while the effect on the other CP -violation parameters is negligible. The loss is compensated by a significant simplification of the fit model, as will be discussed in detail in Sec. 6.3. To obtain the average resolution, σ t (δ t ) in Eq. (6) is integrated over the distribution of δ t from background-subtracted B 0 s → K + K − decays, and an averaged resolution ofσ t = 42.9 ± 0.1 fs is obtained. A dependence of the resolution on the decaying particle mass is found when repeating the procedure using a sample of Υ (1S) → µ + µ − decays instead of the J/ψ → µ + µ − sample, which yieldsσ t = 44.1 ± 0.1 fs. The average between the two calibrations,σ t = 43.5 fs, is used in the fit to data with the simultaneous method, and the difference between them is considered in the determination of the related systematic uncertainty.
In the fit to the J/ψ → µ + µ − data sample, an offset of the mean of the triple-Gaussian model is observed and attributed to a misalignment in the vertex detector. The size of the bias, µ t = −6.5 fs, is used as mean value in the resolution model in both fit methods.

Fitting methods
Two independent methods, called simultaneous and per-candidate, are used to determine the CP -violation parameters in the B 0 → π + π − and B 0 s → K + K − decays, while the simultaneous method also determines the direct CP -asymmetries in B 0 → K + π − and B 0 s → K − π + . A comparison of their respective results serves as validation of the measurements. The common aspects of the two methods are described in Sec. 6.1 and 6.2, while the specific details of each one are discussed in Sec. 6.3 and 6.4.

Components of the fit models
For each component, the distributions of the final-state invariant mass, decay time and flavour-tagging assignment with the associated mistag probability are modelled for B 0 (s) candidates. Signal components are B 0 → K + π − and B 0 s → K − π + decays in the K ± π ∓ samples, the B 0 → π + π − decay in the π + π − sample, and the B 0 s → K + K − decay in the K + K − sample. In the π + π − and K + K − samples, a small contribution from B 0 s → π + π − and B 0 → K + K − decays is present and must be taken into account. Cross-feed, combinatorial and three-body background contributions are described by the model. Apart from B-meson decays, the only relevant source of cross-feed background is the Λ 0 b → pK − decay with the proton misidentified as a kaon in the K + K − sample. Considering the PID efficiencies, the branching fractions and the relative hadronisation probabilities [25], the contribution of this background component is expected to be about 2.5% relative to the B 0 s → K + K − decay and is included in the fit. Components describing partially reconstructed three-body B 0 (s) -meson decays and combinatorial background candidates are necessary in all of the three final states.

Decay-time model for two-body B 0 (s) decays
The time-dependent decay rate of a flavour-specific B → f decay and of its CP conjugate B → f , as for the B 0 → K + π − and B 0 s → K − π + decays, is given by the probability density function (PDF) where K FS is a normalisation factor and the discrete variable ψ assumes the value +1 for the final state f and −1 for the final state f . The functions H ± , Ω sig and Ω sig are defined below. The direct CP asymmetry, A CP , is defined in Eq. (4), while the final-state detection asymmetry, A D , and the B 0 (s) -meson production asymmetry, A P , are defined as where ε tot is the time-integrated efficiency in reconstructing and selecting the final state f or f , and ) is the production cross-section of the given B 0 (s) (B 0 (s) ) meson. The asymmetry A P arises because production rates of B 0 (s) and B 0 (s) mesons are not identical in pp collisions. It is measured to be of the order of one percent at LHC energies [60]. From the time-dependent fit it is possible to determine simultaneously A P and the sum A CP + A D . The contribution of A D is subtracted a posteriori as described in Sec. 7.
The variable ξ = (ξ OS , ξ SS ) is the pair of flavour-tagging assignments of the OS and SS algorithms used to identify the B 0 (s) -meson flavour at production, and η = (η OS , η SS ) is the pair of associated mistag probabilities defined in Sec. 4. The functions Ω sig (t, ξ, η) and Ω sig (t, ξ, η) describe how the flavour tagging modifies the time-dependent decay rate. The functions H + (t) and H − (t) are defined as where Γ d and Γ s are the B 0 and B 0 s mean decay widths, respectively, R (t − t ) is the decay-time resolution model described in Sec. 5 and ⊗ denotes the convolution product.
In the case of a decay to a CP eigenstate f , as it is for the B 0 → π + π − and B 0 s → K + K − decays, the decay-time PDF is given by where K CP is a normalisation factor and the functions I + (t) and I − (t) are In this case f is equal to f , hence the final-state detection asymmetry A D is zero. The parameters ∆m d(s) , Γ d(s) , and ∆Γ d(s) are fixed in the fit to data to the values reported in Table 1.

Simultaneous fit method
The simultaneous method relies on a concurrent fit to all the final-state samples (π + π − , K + K − and K ± π ∓ ), modelling the multidimensional space defined by the final-state invariant mass, B 0 (s) decay time, flavour-tagging decision and associated mistag probability for the signal and background components. The models used in the fit are a modification of those described in Ref. [24].
The model describing the invariant-mass shape of the signal components comprises a sum of two Gaussian functions and a Johnson function [61], while the model for cross-feed background is based on a kernel estimation (KDE) method [62] and tuned with simulated decays. The normalisation of each cross-feed background component is determined by rescaling the yields of the decay reconstructed with the correct mass hypothesis by the ratio between the misidentification probability and the PID efficiency for the wrong and correct mass hypotheses.  [59] used in the two methods. For Γ s and ∆Γ s the correlation factor, ρ, between the two quantities is also reported. The decay width difference ∆Γ d is fixed to zero.

Parameter
Value The decay-time model of the signal components is also used for the cross-feed background components originating from the signal decays reconstructed with the wrong mass hypothesis. This is valid under the assumption that the decay-time calculated under the wrong mass hypothesis is equal to that calculated using the correct hypothesis, and is verified using samples of simulated decays. The flavour-tagging assignments and related mistag probabilities for OS and SS taggers enter the time-dependent decay rates of Eqs. (7) and (10) through the functions Ω sig (t, ξ, η) and Ω sig (t, ξ, η). These functions are the same as already used in Ref. [24] with the only difference being that they now depend on the decay time, as do the efficiencies of the SS taggers. This dependence is accommodated using separate efficiencies: one independent of the SS-tagger decision and one specific for the candidates tagged by the SS taggers. More details are reported in App. A.
The decay-time efficiency, ε sig (t), is sculpted by the selection criteria presented in Sec. 3. It is parameterised using an empirical function determined using the B 0 → K + π − calibration decay, whose time-dependent decay rate is independent of the flavour-tagging decision and described by an exponential distribution with Γ d = 0.6588 ± 0.0017 ps −1 [25]. A sample of background-subtracted B 0 → K + π − candidates is obtained from the K ± π ∓ sample in the invariant-mass window 5.23 < m(K ± π ∓ ) < 5.32 GeV/c 2 . The contributions of the combinatorial background, the only non-negligible background in this region, is subtracted by injecting, with negative weights, candidates from the sideband m(K ± π ∓ ) > 5.6 GeV/c 2 . As explained above, the procedure is repeated for the subsample with ξ SS = 0, in order to model the time dependence of the SS-tagging efficiency. For the B 0 → π + π − and B 0 s → K + K − decays, a small correction is applied to the efficiency in order to take into account the differences between signal and calibration modes. The correction for a given mode is a product of the efficiency determined from the B 0 → K + π − data and the ratio between the efficiencies of this mode and of the B 0 → K + π − decay, as determined from simulation.
The final difference with respect to the model used in Ref. [24] is that the decay-time resolution is no longer modelled on a per-candidate basis. This change is made since a correlation between the distributions of the decay-time and decay-time error is observed for the combinatorial background candidates. A full description of this correlation would imply a considerable complication of the fitting model that outweighs the small loss in statistical power that the use of an average decay-time resolution implies. A systematic uncertainty is established in order to cover for possible biases coming from using an average rather than per-candidate decay-time resolution.
The invariant-mass model for the combinatorial background components for each decays is an exponential function, with its slope depending on the decay time, in order to take into account a slight correlation between invariant mass and decay time observed in the high-mass sideband. The time dependence of the slope is studied using a two-dimensional unbinned maximum-likelihood fit to the invariant mass and decay time of the sample in the high-mass sideband above 5.6 GeV/c 2 , where only combinatorial background candidates contribute. The obtained time-dependent mass slope is used for the combinatorial background model in the entire invariant-mass window, going from 5.0 to 6.2 GeV/c 2 . The relative normalisation of each candidate in the sideband is scaled to reproduce that in the total invariant-mass window. A KDE method is applied to the weighted candidates and the output is used to model the decay-time shape of the combinatorial-background component. A dependence of the decay-time shape of combinatorial background candidates on the tagging assignment of the OS-and SS-taggers is also observed. Hence the time dependence of the mass slope is studied separately for the subsamples corresponding to the tagging decision (|ξ OS |, |ξ SS |) = {(1, 1), (1, 0), (0, 1), (0, 0)}. Different weights are determined for each subsample, and also the KDE method is applied separately to each of them. The weighting procedure is the same as employed for the background subtraction used to study the decay-time efficiency for B 0 (s) decays. The functions taking into account the flavour-tagging assignment and mistag probabilities are the same used in Ref [24], but are generalised to consider all the possible combinations of (|ξ OS |, |ξ SS |). Finally, in the case of the K ± π ∓ samples, possible asymmetries in the flavour-tagging or reconstruction efficiencies for the two charge-conjugate final states are taken into account.
The invariant-mass model of partially reconstructed B 0 (s) decays is the same as that used in Ref. [24], comprising the sum of two Gaussian functions, which are defined using the same parameters as in the signal model and are convolved with ARGUS functions [63]. For the K ± π ∓ sample two three-body background components are used: one describing three-body B 0 and B + decays and another describing three-body B 0 s decays. For the π + π − and K + K − samples a single ARGUS component is found to be sufficient to describe the invariant-mass shape in the low-mass region. The shape of the decay-time distribution is obtained by applying a KDE method to the candidates in the low-mass sideband below 5.2 GeV/c 2 , after subtracting the combinatorial background contribution, as explained above. This is repeated separately for the candidates with |ξ SS | = 0 and |ξ SS | = 0, since a difference in the decay-time shape is observed in data for the two subsamples. The functions used to take into account the flavour-tagging information are the same as used for the combinatorial background model, but with independent parameters. Also for this component possible differences in flavour-tagging and reconstruction efficiencies between the K + π − and π + K − final states are taken into account in the same way as used for the combinatorial background model.

Per-candidate fit method
The per-candidate method relies on independent fits to the π + π − and K + K − samples with all background components statistically subtracted using the sFit technique [58,64] with the π + π − and K + K − invariant mass as the discriminating variable. Hence only the decay-time distributions are modelled for the signal modes B 0 → π + π − and B 0 s → K + K − . The invariant-mass distributions of the B 0 s → K + K − and B 0 → π + π − signal components are modelled with the sum of two Crystal Ball functions [65] where the tail parameters are fixed to the values obtained from the simulation. The mean and width of the Gaussian core are allowed to vary in the fit for the B 0 s → K + K − and B 0 → π + π − signal modes, while these parameters are constrained for the B 0 → K + K − and B 0 s → π + π − signal components using the known mass difference between B 0 and B 0 s and the ratio of resolutions obtained from simulations, respectively. The decay-time model for the signal components is the function described in Sect. 6.2, multiplied on a per-candidate basis with the acceptance functions described below.
The invariant-mass distributions of the misidentified background candidates from other two-body B 0 (s) decays are modelled with templates from simulations and their yields are constrained using efficiencies measured in data calibration samples. The three-body background components, which are the same as in the simultaneous method, are modelled using an exponentially modified Gaussian PDF.
The decay-time resolution consists of a single Gaussian function with its width varying candidate by candidate, depending on the decay-time error δ t for each candidate and calibrated as presented in Sec. 5. The per-candidate acceptance function is determined with the swimming method [34][35][36][37] by artificially changing the decay time of the B 0 (s) meson and re-evaluating whether the candidate would have been accepted by the selection requirements that are known to bias the decay-time measurement. The decay time is changed by moving the position of every PV in the event along the direction of the B 0 (s) momentum vector. For decay times for which the candidate is accepted the efficiency is 1, otherwise the efficiency is 0. By scanning a range of hypothetical decay times, a series of top-hat functions are constructed for each candidate as it changes from being rejected, to being accepted, finally to being rejected again. 4 The procedure is re-evaluated in steps of 50 µm along the B 0 (s) momentum vector and, when the selection decision changes, the position at which this change occurs is determined with a finer granularity, giving an overall resolution of 0.5 µm on the decay-time efficiency. The effective lifetime measured on the fully simulated B 0 s → K + K − events, assuming an exponential decay-time model and using only the swimming-based efficiency for this simulation, is found to be 1.416 ps. Compared to a generated effective lifetime of 1.394 ps it exhibits a bias of 1.5%. This arises from effects that are not fully modelled in the swimming method and can result in an incorrect measurement of the parameter A ∆Γ KK , for which high precision is expected. To correct for this, an additional decay-time efficiency weight is applied by comparing the decay-time efficiency extracted using the swimming method for the B 0 → K + π − data with the decay-time efficiency determined from the ratio of background-subtracted B 0 → K + π − events and the unbiased decay-time PDF. The unbiased decay-time PDF consists of an exponential function, whose decay time is fixed to the known B 0 lifetime, convolved with a Gaussian function to account for the intrinsic decay-time resolution. The width of the Gaussian is fixed to the effective decay-time resolution as detailed in Sec. 5. The ratio of these efficiencies is modelled with an empirical function where t is the decay time of the candidate and p {0,1,2,3} are free parameters measured in the fit to the ratio. Applying this weight to the swimming-based efficiency allows to correctly recover the effective lifetime of the simulated B 0 s → K + K − decays and the mean lifetime of B 0 → K + π − decays extracted from the K ± π ∓ samples. 7 Detection asymmetry between K − π + and K + π − final states In order to extract the CP asymmetries A B 0 CP and A B 0 s CP from the asymmetries measured through the simultaneous fit, an estimation of the nuisance experimental detection asymmetry is required as indicated in Eq. 7. This asymmetry is a consequence of the different efficiency for selecting the B 0 → K + π − and B 0 s → K − π + decays and their charge-conjugate final states. To an excellent approximation, it can be expressed as the sum of two contributions where A Kπ det is the asymmetry between the selection efficiencies without the application of the PID requirements and A Kπ PID is the asymmetry between the efficiencies of the PID requirements selecting the two final states. The convention used in the following to determine A Kπ det and A Kπ PID is such that a positive value of the asymmetry means a larger efficiency for the K − π + pair with respect to the K + π − pair. As a consequence of this convention, the values reported below for the B 0 and B 0 s asymmetries must be used with an inverted sign for the B 0 s → K − π + decay. The final-state detection asymmetry is determined using D + → K − π + π + and D + → K 0 π + control modes, with the neutral kaon decaying to π + π − , following the strategy used in Ref. [24]. Assuming negligible CP violation in these Cabibbo-favoured D-meson decays, the raw asymmetries between the measured yields of D + and D − decays can be written as where A D + P is the asymmetry between the production cross-sections of D + and D − mesons, and A π det (A K 0 det ) is the asymmetry between the detection efficiencies of π + (K 0 ) and π − (K 0 ) mesons. The difference between Eqs. (14) and (15) leads to The asymmetry A K 0 det includes the effects from the kaon mixing and CP violation, and was estimated to be (0.054 ± 0.014) % [66]. The asymmetries A D + P and A π det can depend on the kinematics of the D + and π + mesons. To obtain a better cancellation of these nuisance asymmetries in Eq. (16), the momentum and the transverse momentum of the D + and π + mesons from the D + → K − π + π + sample are simultaneously weighted to match the corresponding distributions in the D + → K 0 π + sample. The A Kπ det is determined in intervals of the kaon momentum, to account for the kinematic-dependent variation of the interaction cross-sections of positive and negative kaons with the detector material. This binned asymmetry is averaged over the momentum distribution of the kaon in the B 0 → K + π − and B 0 s → K − π + decays, giving no difference between the absolute values of the corrections for the two modes. The final-state detection asymmetry values for the  The asymmetry between the PID efficiencies is computed in intervals of momentum, pseudorapidity and azimuthal angle of the two final-state particles, using the D * + → D 0 (K − π + )π + calibration samples, as discussed in Sec. 3. The computation is repeated using several binning schemes, and then the average and standard deviation of the PID asymmetries determined in each scheme are used as the central value and associated uncertainty for A Kπ PID , respectively. The PID asymmetry is calculated taking into account the differences in the running conditions of the two years of data taking and the numerical  results are:

Fit results
The results obtained from unbinned maximum likelihood fits to data of the models described in Secs. 6 are presented in the following. Their comparison is also discussed.

Simultaneous method
The simultaneous fit to the final-state invariant mass, the B 0 (s) decay time, and the tagging decisions and their associated mistag probabilities of the π + π − , K + K − and K ± π ∓ samples determines the coefficients C ππ , S ππ , C KK , S KK , A ∆Γ KK and the CP asymmetries A B 0 CP and A B 0 s CP . The signal yields are N (B 0 → π + π − ) = 45 620 ± 260, N (B 0 s → K + K − ) = 70 310 ± 320, N (B 0 → K + π − ) = 140 340 ± 420 and N (B 0 s → K − π + ) = 10 580 ± 150, where uncertainties are statistical only. The distributions of the mass and decay time of the selected candidates are shown in Figs. 2, 3 and 4, for the K ± π ∓ , π + π − and K + K − samples, respectively. The time-dependent asymmetries, obtained separately by using the OS or the SS tagging decisions, for the B 0 (s) candidates in the region 5.20 < m(K ± π ∓ ) < 5.32 GeV/c 2 , dominated by the B 0 → K + π − decay, are shown in (−0.60 ± 0.49)% and (−1.2 ± 1.5)%, respectively, where uncertainties are statistical only. They are consistent with the expectations from Ref. [60]. The time-dependent asymmetries for the π + π − candidates with 5.20 < m(π + π − ) < 5.35 GeV/c 2 , and for the K + K − candidates with 5.30 < m(K + K − ) < 5.45 GeV/c 2 , dominated by the corresponding B 0 → π + π − and B 0 s → K + K − signal components, are shown in Fig. 6, again separately for the OS and SS tagging decision. The effective tagging powers for the B 0 → π + π − and B 0 s → K + K − decays are

Per-candidate method
The signal yields in the B 0 → π + π − and B 0 s → K + K − decays, used to determine the CP -violating parameters with the per-candidate method, are in agreement with those of the simultaneous method. The parameters ∆m d(s) , Γ d(s) , and ∆Γ d(s) are fixed to the values reported in Table 1. The value of the production asymmetry is fixed to that measured by 1.9 fb sample − K + K tagging K SS Figure 6: Time-dependent asymmetries for (top) π + π − and (bottom) K + K − candidates with 5.20 < m(π + π − ) < 5.35 GeV/c 2 and 5.30 < m(K + K − ) < 5.44 GeV/c 2 , respectively: (left) using the OS-tagging decision and (right) using either the SSc-tagging decision (for the π + π − candidates) or the SSK-tagging decision (for the K + K − candidates). The result of the simultaneous fit is overlaid. The asymmetry for the K + K − candidates is folded into one mixing period 2π/∆m s and the parameter t 0 = 0.2 ps corresponds to the minimum value of the decay-time used in the fit.
the simultaneous method. The fits to the π + π − and K + K − invariant-mass spectra are shown in Figs. 7 and 8 along with the decay-time fits to the B 0 (s) mesons having their flavours tagged.
The results for the CP -violating parameters using the per-candidate method are C ππ = −0.338 ± 0.048, S ππ = −0.673 ± 0.043, where the uncertainties are statistical only.

Comparison
To evaluate the compatibility of the results from the two methods, their statistical correlation is determined from 500 simulated pseudoexperiments. The correlation is found to be approximately 84% for all CP -violating parameters. This is used to determine the uncorrelated statistical uncertainty on the difference between the results of the two methods. The pseudoexperiments also confirm the smaller total uncertainty observed by the simultaneous method. A sizeable difference between the two results is observed for A ∆Γ KK . This difference is reduced to approximately 1.5 standard deviations when taking into account the systematic uncertainties due to the determination of the decay-time efficiency (see Tab. 2), which are completely uncorrelated between the two methods. Adding in quadrature the uncorrelated statistical and systematic uncertainties, the results are found to be compatible within one standard deviation. The resulting contour plots from measuring C ππ , S ππ , C KK and S KK are given in Fig. 9. Given the large correlation between the two determinations, the values obtained from the simultaneous method are quoted as the LHCb results. They are chosen due to the slightly smaller total uncertainty and the fact that the simultaneous method gives also the direct CP asymmetries allowing for a complete combination with the results published in Ref. [24].  Figure 9: Two-dimensional 68% and 95% confidence-level regions for the measured CP -violating parameters of the B 0 → π + π − (left) and B 0 s → K + K − (right) decays from the two methods. The simultaneous method is shown in purple while the per-candidate method in green. Previous measurements of these parameters are also shown, with the LHCb Run 1 result in blue, the Belle result in yellow and the BaBar result in red. The confidence-level regions are calculated using only the statistical uncertainties of all the measurements. The correlation is found to be approximately 84% for all CP -violating parameters between the simultaneous and per-event methods.

Systematic uncertainties
The systematic uncertainties are evaluated for both the simultaneous and the per-candidate methods, and the total systematic uncertainties for both results are given in Table 2. A full description of the systematic uncertainties is only given for the simultaneous method since it is used as the LHCb result and for combination with the Run 1 measurement. Hence the description given in this section and the breakdown of the individual components in Table 2 refers to that method. The main differences in systematic uncertainties between the two methods are briefly discussed at the end of this section.
The systematic uncertainties on the CP -violating parameters are determined following two approaches. In the first case the fit to data is repeated a large number of times, each time modifying the values of the input parameters. This approach is used to account for the knowledge of external inputs whose values are fixed in the fit. In the second case, pseudoexperiments are performed according to the default model and both the default model and modified models are used to fit the generated data. This strategy is used to account for the systematic uncertainties due the the assumptions on the fitting model. In both cases the difference between the default and alternative results for the CP asymmetries is measured, and the mean and width of the obtained distribution is used to assign a systematic uncertainty.
Three sources of systematic uncertainty are considered on the invariant-mass model. First, the systematic uncertainty due to a possibly imperfect description of the massresolution function, used for both signal and cross-feed background components, is determined by replacing the double Gaussian function with a single Gaussian model. Second, the systematic uncertainty associated to the combinatorial background model is assessed using an alternative model with no correlation between decay time and invariant mass. Finally, a systematic uncertainty associated with the model adopted for the three-body background components is determined by fitting a set of pseudoexperiments, after removing the candidates with an invariant mass below 5.2 GeV/c 2 and ignoring the components describing this background contributions in the model.
The PID efficiencies and misidentification probabilities govern the amount of crossfeed background components. A systematic uncertainty related to their calibration is determined by repeating the fit to data changing those values according to their uncertainties estimated from the calibration samples.
The effect of ignoring the small fraction of B 0 s candidates originating from decays of the B + c meson is studied by injecting simulated B + c → B 0 s X decays (where X stands for any additional particle in the final state) into the pseudoexperiments, where the relative B + c yield is determined from Ref. [53]. No systematic uncertainty is assigned for the B 0 CP -violating parameters since the B + c → B 0 X decay is Cabibbo suppressed. Systematic uncertainties associated with the calibration of the OS and SSc flavourtagging responses are determined using an alternative relation between η OS(SS) and the calibrated mistag probability ω OS(SS) . The linear relation connecting the two quantities in Eq. (23) is replaced with a second-order polynomial. A similar approach is also used for the SSK tagger, but the values of the parameters of the alternative relations are first determined from the B 0 s → D − s π + sample and then used in the fit to data. For the SSK tagger an additional systematic uncertainty associated with the calibration of the flavour-tagging response is determined by varying the calibration parameters according to their uncertainties and correlations.
Regarding the decay-time model, a systematic uncertainty associated with the uncertainties on the parameters reported in Table 1 is determined by repeating the simultaneous fit using different fixed values, generated according to their uncertainties and correlations. The systematic effect due to the decay-time resolution can be decomposed into three contributions: one due to the calibration of the resolution width, another one due to the calibration of the bias in the determination of the decay time, and the last due to the usage of an average decay-time resolution instead of a per-candidate value. The first effect is estimated varying the value of the averaged decay-time resolution width according to a Gaussian distribution with mean equal to the default value, reported in Section 5, and with a width equal to the difference between the decay-time resolution for the fully simulated J/ψ → µ + µ − and B 0 s → K + K − decays. The second effect is determined varying the mean of the decay-time resolution model according to a Gaussian centered at the default value and with the width of 2 fs. The last contribution, due to the usage of an average decay-time resolution instead of a per-candidate value, is evaluated by fitting a set of pseudoexperiments with both the decay-time resolution models.
Three sources of systematic uncertainty related to the knowledge of the decay-time efficiency are identified. A systematic uncertainty on the chosen model is assessed by replacing the effective function with a cubic-spline polynomial in an alternative model. Second, a systematic uncertainty arising from the limited calibration-sample size is computed by varying the parameters governing the decay-time acceptance according to their uncertainties and correlations. An additional systematic uncertainty due to the imperfect description of the ratios between the decay-time efficiency of the various signal and calibration modes, determined from fully simulated samples, is estimated. In this case, the alternative model is created assuming that all the decay-time efficiencies are equal to that of the B 0 → K + π − decay.
To determine the systematic uncertainty associated with the choice of the decaytime model for the cross-feed background component, an alternative model is created by disabling the oscillating component of the cross-feed background model. This means assuming no CP violation for both, the B 0 → K + π − component in the π + π − and K + K − samples, and for the B 0 → π + π − and B 0 s → K + K − components in the K ± π ∓ sample. Finally, the uncertainties associated with the detection and PID asymmetries reported in Eqs. (17) and (18) are accounted for as systematic uncertainties on A B 0 CP and A B 0 s CP . The total systematic uncertainties are obtained as the quadratic sum of the individual contributions, and they are smaller than the corresponding statistical uncertainties for all the CP -violating parameters apart from A ∆Γ KK . The dominating systematic uncertainty for this parameter is related to the knowledge of how the efficiency varies with the decay time, whose knowledge is limited by the size of the calibration sample of B 0 → K + π − decays.
Most of the sources of systematic uncertainties related to the per-candidate method are the same as those on the simultaneous method. The systematic uncertainties are also similar in size. The main difference is a smaller uncertainty related to the decaytime acceptance in the per-candidate method. This uncertainty, which is uncorrelated between the two methods due to their different strategies, mainly affects A ∆Γ KK and largely cancels in the other parameters. The second most important difference is due to systematic uncertainties related to flavour tagging, where the uncertainties are larger in the per-candidate method, which arises from the different approaches of incorporating this information in the two fits. These differences in systematic uncertainties illustrate the strength of validating the result with two different methods. Table 2: Systematic uncertainties on the CP -violating parameters. The values given for each individual contribution to the systematic uncertainty are those for the simultaneous method. The total systematic uncertainties are given both for the simultaneous and the per-candidate methods. The dash indicates that the uncertainty is not applicable.

Results
The final results for the time-dependent CP violation in B 0 → π + π − and B 0 s → K + K − decays, and of the CP asymmetries in B 0 → K + π − and B 0 s → K − π + decays are     Tables 3 and 4, respectively.
The results are compatible with the previous LHCb measurement in Ref. [24]. A combination is performed between the results in this paper and those based on the Run 1 data sample reported in Ref. [24]. Since the values of Γ s and ∆Γ s used as input to the fit have changed with respect to Ref. [24], the Run 1 analysis is updated to account for the new values. The main variation is observed for the central value of A ∆Γ KK that changes from −0.79 ± 0.07 to −0.97 ± 0.07. The large variation of A ∆Γ KK is expected, given its correlation of 0.91 with Γ s , and the significant change in the value of Γ s from 0.6654 ± 0.0022 ps −1 to 0.6563 ± 0.0021 ps −1 [59]. The only other variation is for C KK , moving from 0.20 ± 0.06 to 0.19 ± 0.06. The compatibility between the updated Run 1 result and the numbers reported in Eq. (19) is computed by means of χ 2 test statistic, finding the two sets of results in agreement with a p-value of 0.68.
The full statistical and systematic covariance matrices of the two results are taken into account in the combination. The only relevant correlation between the two results is related to the values of the input parameters in Table 1, hence the corresponding systematic uncertainties are removed from the covariance matrices of the two results, before combining them. The systematic uncertainty due to these input parameters is included again by summing the corresponding covariance matrix to the covariance matrix  and their correlation matrix is reported in Table 5.

Concluding remarks
The time-dependent CP asymmetries of B 0 → π + π − and B 0 s → K + K − decays and the time-integrated CP asymmetries in B 0 → K + π − and B 0 s → K − π + decays are measured using a data sample of pp collisions corresponding to an integrated luminosity of 1.9 fb −1 , collected with the LHCb detector at a centre-of-mass energy of 13 TeV. The measurements are compatible with previous LHCb determinations of the same quantities obtained with Run 1 data [24] and are combined with them. The measurements of C ππ , S ππ , A B 0 CP and A B 0 s CP are in good agreement with previous results from other experiments [30][31][32][33] and are the most precise from a single experiment to date.
A χ 2 test statistic is used to determine the significance for (C KK , S KK , A ∆Γ KK ) to differ from (0, 0, −1) and for (C KK , S KK ) to differ from (0, 0). The significance for the combined LHCb results is found to be of 6.5 and 6.7 standard deviations, respectively. According to the test of the SM proposed in Ref. [7], the following sum must be satisfied where B (B 0 → K + π − ) and B (B 0 s → K − π + ) are CP -averaged branching fractions. The LHCb measurements of the relative fragmentation-fraction ratio between B 0 s and B 0 mesons f s /f d = 0.259 ± 0.015 [67], [68] and Γ s /Γ d [59] are used in this test along with the measurements of A B 0 CP and A B 0 s CP . The value ∆ = −0.085 ± 0.025 ± 0.035 is obtained, where the first uncertainty is from the measurements of the CP asymmetries and the second is from the other inputs in Eq. (19). With the present experimental precision, ∆ is in agreement with zero within two standard deviations.
Owing to the measurements reported in this paper, improved constraints on the CKM angles and B 0 s mixing phase can be obtained, as outlined in Refs. [9,14,15]. The comparison of these precises determinations, based on decays receiving sizeable loop-level contributions, with those provided by the study of the decays dominated by tree-level amplitudes, will constitute a stringent test of the SM hypothesis.

Appendix A Additional information on flavour-tagging A.1 Formalism
The functions Ω sig (t, ξ, η) and Ω sig (t, ξ, η) in Eqs. (7) and (10) are where η is the mistag probability computed by the flavour-tagging algorithms and discussed in Sec. 4, Ω sig OS and Ω sig OS are the same functions used in Ref. [24], i.e.
while Ω sig SS and Ω sig SS are Here, ε sig tag (ε sig tag ) is the probability that the flavour of a B 0 (s) (B 0 (s) ) meson is tagged, which in the case of the SS tagger depends on the decay time 5 ; ε(t) is the decay-time efficiency independent from the decision of the SS-tagger, such that ε(t) − ε sig SS (t) is the decay-time efficiency for candidates that have ξ SS = 0; ω sig tag (η tag ) and Ω sig tag (η tag ) are the calibrated mistag probabilities as a function of η tag for B 0 (s) and B 0 (s) mesons; h tag sig (η tag ) is the PDF describing the distribution of η tag for tagged candidates, and U (η tag ) is a uniform distribution of η tag . It is empirically observed that, to a good approximation, η tag and ω tag are related by a linear function, i.e.
whereη tag is a fixed value, chosen to be equal to the mean value of the η tag distribution to minimise the correlation among the parameters. To reduce the correlation among ε sig tag andε sig tag , and p tag 0 ,p tag 0 , p tag 1 , andp tag 1 , these variables are conveniently parameterised as ε sig =ε sig tag (1 + ∆ε sig tag ), (24) ε sig =ε sig tag (1 − ∆ε sig tag ), wherep tag 0,1 and ∆p tag 0,1 are the average and the asymmetry between p tag 0,1 andp tag 0,1 , andε sig tag and ∆ε sig tag are the average and the asymmetry between ε sig tag andε sig tag . The dependence on the decay-time is considered only for the averaged efficiencyε sig SS and not for the asymmetry ∆ε sig tag . The strategy used to determine the decay-time efficiencies ε(t) and ε sig SS (t) is reported in Sec. 6.3. The description of h SS sig (η) for the SS taggers is presented in Secs. A.2 and A.3, respectively.
The PDF Ω comb ( ξ, η) for the combinatorial background is empirically parameterised by where the functions Ω OS,j comb (ξ OS , η OS ) (j = 0, 1) are the PDF for the ξ OS and η OS , defined as The variables ε tag comb andε tag comb are the probabilities to tag a combinatorial background candidate as B 0 (s) or B 0 (s) and h tag comb (η tag ) is the distribution of η tag . The distribution is described using histograms taken from the right-hand sideband with invariant-mass range between 5.6 GeV/c 2 and 6.2 GeV/c 2 . The tagging efficiencies are parameterised by ε tag comb =ε tag comb 2 (1 + ∆ε tag comb ), ε tag comb =ε tag comb 2 (1 − ∆ε tag comb ), such that the fit to data determines the average probability to tag combinatorial background as B 0 (s) or B 0 (s) ,ε tag comb , and the asymmetry between the two probabilities, ∆ε tag comb . For the OS tagger, the distinction labelled by the index j = 0, 1 is used to differentiate the OS-tagger probability between cases that have ξ SS = 0 (j = 0) and ξ SS = 0 (j = 1). In the case of the K + π − and K − π + samples, Eq. × Ω comb ( ξ, η), where A comb raw is the total asymmetry between the combinatorial-background yields in the K + π − and K − π + samples, A comb OS and A comb SS are additional parameters that take into account the possibility that the flavour-tagging probabilities ε tag comb andε tag comb may depend on the final state.
The PDFs in Eqs. (20), (25) and (29) are valid if η OS and η SS are uncorrelated. This assumption is verified by means of background-subtracted [58] signal candidates, and of candidates from the high-and low-mass sidebands for the combinatorial and three-body backgrounds components, respectively.

A.2 Combination of the single SS and OS taggers
The calibration parameters governing the relations in Eqs. (7) are determined separately for the individual SS and OS taggers by means of a binomial regression to the tagged decay-time distribution of background-subtracted B 0 → D − π + decays. Then the extracted calibration parameters of the SSπ and SSp taggers are used to combine the two taggers into a unique one (SSc) with decision ξ SSc and mistag probability η SSc . The assumption of a linear relation between η tag and ω tag for each tagger is validated splitting the sample in bins of η SS(π,p) , estimating the average mistag fraction in each bin by means of the binomial regression. Similarly, the various OS taggers are combined together into a unique OS tagger with decision ξ OS and mistag probability η OS , and the same linearity check is performed.
The PDFs h SS sig (η SSc ) describing the η SSc distributions for the signal B 0 mesons are determined using background-subtracted distributions of B 0 → D − π + decays. As the pion and proton kinematics are correlated with those of the B 0 meson, the performance of the SSπ and SSp taggers also depend on the latter. The differences between the B 0 -meson kinematics and other relevant distributions in B 0 → D − π + and B 0 → π + π − decays, due to the different topologies and selection requirements, are taken into account by means of a weighting procedure to the B 0 → D − π + sample. It is empirically observed that the B meson transverse momentum and the number of hits in the SPD detector distributions need to be equalised.

A.3 Calibration of the SSK tagger
The natural control mode to calibrate the response of the SSK tagger would be the B 0 s → K − π + decay. However, since the signal yield of this decay is approximately 8% of that of the B 0 → K + π − decay and 20% of that of the B 0 s → K + K − decay, it would not be possible obtaining a reliable calibration. Furthermore, the calibration parameters of the SSK tagger would be affected by large uncertainties, limiting the precision on C KK and S KK . Therefore, the calibration is performed with a large sample of B 0 s → D − s π + decays. Analogously to the SSπ and SSp case, the SSK-calibration parameters are determined using a binomial regression to the tagged decay-time distribution of the B 0 s → D − s π + decays. Also in this case the regression is performed using the flavour-tagging information on a per-candidate basis, determining the calibration parameters directly, and a check of the linear relation between η SSK and ω SSK is performed. As described for the B 0 → D − π + sample in the previous section, a weighting procedure is applied to the B 0 s → D − s π + sample in order to equalise the signal distribution of the B meson transverse momentum and the number of hits in the SPD detector. The PDF h SSK sig (η SSK ) for B 0 s → K + K − decays is determined using a background-subtracted histogram of the same weighted sample of B 0 s → D − s π + decays used for the calibration.