Search for the suppressed decays $B^{+}\rightarrow K^{+}K^{+}\pi^{-}$ and $B^{+}\rightarrow \pi^{+}\pi^{+}K^{-}$

A search is made for the highly-suppressed B meson decays $B^{+}\rightarrow K^{+}K^{+}\pi^{-}$ and $B^{+}\rightarrow \pi^{+}\pi^{+}K^{-}$ using a data sample corresponding to an integrated luminosity of 3.0 $fb^{-1}$ collected by the LHCb experiment in proton-proton collisions at centre-of-mass energies of 7 and 8 TeV. No evidence is found for the decays, and upper limits at 90\% confidence level are determined to be $\mathcal{B}(B^{+}\rightarrow K^{+}K^{+}\pi^{-})<1.1\times 10^{-8}$ and $\mathcal{B}(B^{+}\rightarrow \pi^{+}\pi^{+}K^{-})<4.6\times 10^{-8}$.


Introduction
Transitions of the type b → ssd and b → dds are rare in the Standard Model (SM) [1,2]. The calculation of the b → ssd amplitude results in branching fractions of at most O(10 −11 ), the exact value depending on the unknown relative phase between t and c quark contributions in the W ± -exchange box [3], as shown in Fig. 1. The magnitude of the b → dds amplitude is expected to be even smaller due to the relative |V td /V ts | factor, leading to predicted branching fractions of O (10 −14 ) [4].
Physics beyond the Standard Model (BSM) could result in enhanced branching fractions that can be detected at current experiments. Theoretical models that have been investigated include the Minimal Supersymmetric Standard Model with and without R-Parity Violation, variations of the two Higgs doublet model, and extensions of the SM where an additional flavour changing Z neutral boson is present [3,4]. In these SM extensions, for certain plausible values of the model parameters, b → ssd and b → dds transitions may lead to branching fractions of up to 10 −8 and 10 −7 , respectively. It has been suggested by these theoretical studies that the most suitable three-body decay modes to see the effects of BSM physics in such transitions are the B + → K + K + π − and B + → π + π + K − decays, where two particles in the final state have the same flavour and charge. 1 An upper limit of 1.29 × 10 −4 at 90% confidence level on the branching fraction for the B + → π + π + K − decay was first determined by OPAL [5]. Improvements in sensitivity were obtained by the e + e − B-factories, and currently the 90% confidence level upper limits for B + → K + K + π − and B + → π + π + K − decays are 1.6 × 10 −7 and 9.5 × 10 −7 , respectively [6][7][8].
This paper reports on the search for the suppressed decays B + → K + K + π − and B + → π + π + K − using data samples corresponding to 1.0 and 2.0 fb −1 collected by LHCb at √ s = 7 and 8 TeV, respectively. The corresponding unsuppressed decays B + → K + K − π + and B + → π + π − K + are used for normalisation. 1 The inclusion of charge-conjugate decays is implied throughout.

Detector and simulation
The LHCb detector [9,10] 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 silicon-strip detectors and straw drift tubes placed downstream of the magnet. The polarity of the dipole magnet is reversed periodically throughout data-taking. The configuration with the magnetic field vertically upwards (downwards) bends positively (negatively) charged particles in the horizontal plane towards the centre of the LHC. 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 impact parameter (IP) is the minimum distance of a track to a primary vertex (PV) and 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 [11]. 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 [12].
The online event selection is performed by a trigger [13], 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 stage, the candidates are triggered either one of the particles from the b candidate decay depositing a transverse energy of at least 3500 MeV in the calorimeter, or by other activity in the event, mainly associated with the decay products of the other b hadron produced in the pp primary interaction.
Simulated B + → K + K ± π ∓ and B + → π + π ± K ∓ decays, generated uniformly in phase space, are used to optimize the suppressed signal selections and to evaluate the ratios of the efficiencies for each suppressed decay mode relative to their corresponding unsuppressed decay modes. In the simulation, pp collisions are generated using Pythia 8 [14,15] with a specific LHCb configuration [16]. Decays of hadronic particles are described by EvtGen [17] in which final state radiation is generated by Photos [18]. The interaction of the generated particles with the detector and its response are implemented using the Geant4 toolkit [19] as described in Ref. [20].

Event selection
The candidate B + → K + K ± π ∓ and B + → π + π ± K ∓ decays are reconstructed using three charged tracks with mass hypotheses and total charge consistent with the decay. The final state particles are required to have a good track fit with a reduced chi-square χ 2 /ndf < 3. All three tracks must have momentum p > 1500 MeV/c, transverse momentum p T > 100 MeV/c, and χ 2 IP > 1 with respect to all PVs in the event. The quantity χ 2 IP is defined as the difference between the vertex-fit χ 2 of the PV reconstructed with and without the considered track. Combinatorial backgrounds are suppressed by requiring that the scalar sum of the p T of the tracks is greater than 4500 MeV/c and the sum of the tracks' χ 2 IP > 200. The track with the highest p T must have IP > 0.05 mm. The second highest track p T is required to be greater than 900 MeV/c. The maximum distance of closest approach between tracks has to be less than 0.2 mm.
The information from the RICH, calorimeter and muon systems is used for particle identification (PID). Muons are rejected by a veto applied to each track. Loose kaon and pion PID is required for the remaining charged tracks to reduce the number of wrong combinations before forming a B + candidate.
The reconstructed B + candidates are required to have an invariant mass in the range 5080 − 5480 MeV/c 2 , p T > 1700 MeV/c, χ 2 IP < 10, vertex fit χ 2 /ndf < 12, distance between the PV and the decay point (or secondary vertex, SV) greater than 3 mm, and a significant displacement between primary and secondary vertex, with the three dimensional χ 2 -distance between the two greater than 700. When more than one PV is reconstructed, the one with the minimum χ 2 IP for the B + candidate is chosen. The cosine of the angle θ f between the reconstructed momentum of the B + candidate and the B + candidate flight direction is required to be cos θ f > 0.99998. The pointing variable is required to be less than 0.05, where P B is the total momentum of the three-particle final state and 3 i q T i is the sum of the transverse momenta of the three tracks with respect to the momentum direction of the B + candidate. These requirements remove additional combinatorial background and partially reconstructed b hadron decays.
Depending on its charge and the polarity of the dipole magnet, a charged particle traversing the magnetic field can be bent horizontally into or out of the detector acceptance. To minimise charge-dependent differences in the reconstruction efficiencies for the signal or normalisation channels caused by the magnetic field, an additional criterion is placed on the x and z components 2 of the B + candidate momentum such that |p x | ≤ 0.317 × (p z − 2400 MeV/c) [21], which ensures that tracks of both charges are well contained in the detector acceptance.
A B + candidate is rejected if the invariant mass formed from two of the charged tracks with opposite charge is within 25 MeV/c 2 of the D 0 mass. The D 0 veto suppresses possible background from B + → D 0 h + decays. Backgrounds from B 0 → D − h + decays are found to be negligible. For B + → π + π − K + decays only, an additional invariant mass criterion of |3104 MeV/c 2 − m π + π − | > 20 MeV/c 2 is required to eliminate contamination J/ψ → µ + µ − , where the muons are misidentified as pions, and from J/ψ → π + π − .
The reconstructed candidates that meet the above criteria are filtered using a boosted decision tree (BDT) algorithm [22,23]. The BDT is trained with a sample of simulated signal candidates and a background sample of data candidates taken from the B + candidate invariant mass sideband above 5400 MeV/c 2 , which is dominated by combinatorial background. The training is performed separately for events that have been triggered by information from the signal decay only and for events that have been triggered by other particles. This is required as the two trigger scenarios have different sensitivities to the signal and background. The input variables are chosen to produce the best discrimination and to minimise the dependence of the BDT on the mass of the B + candidate, the PID variables and on the kinematic configuration of the B + candidate parametrised in the Dalitz plane. The twelve variables used by the BDT are: the B + candidate p T ; the flight distance between the PV and B + candidate SV; the angle θ f ; the B + candidate pointing angle θ p ; the p T of the track with the lowest p T ; the sum of the tracks' p T ; the sum of the χ 2 IP of the three tracks; the IP of the track with the highest p T , with respect to the PV; the momentum p of each of the three tracks; and the χ 2 IP of the track reconstructed with the π hypothesis for B + → K + K ± π ∓ decays or the K hypothesis for B + → π + π ± K ∓ decays. The same set of variables is found to result in a robust multivariate classifier for the four decay channels under consideration.
, is used to identify the optimal BDT selection criteria. Here is the simulated signal selection efficiency, N B is the number of background events that pass the selection and have a mass in a ±50 MeV/c 2 window around the B ± mass [24], and N σ is the required significance, expressed in terms of standard deviations from the hypothesis of a null signal [25]. The quantity FOM is optimized with N σ = 5, but the final result is robust against changes of one or two units in N σ . For events that pass the selection criteria described above, the optimal BDT selection criterion for B + → K + K + π − results in 85% of the simulated signal events being accepted and 71% of the background events rejected. For B + → π + π + K − , 67% of the simulated signal events are accepted and 91% of the background events are rejected.
After the BDT selection criterion has been applied, each track is required to pass PID criteria. Each track has a probability to be a kaon and a probability to be a pion, leading in total to six possible PID assignments for a B + candidate. The same FOM optimisation described above is performed for each of the six cases in turn, starting with the track with the highest p T . After the application of these criteria, less than 2-4% have more than one candidate, depending on the decay mode. For these multiple candidate events, one candidate is chosen at random and the others discarded.
The efficiencies of all the selection requirements are calculated with simulated events. The PID efficiency for hadrons is determined from data using large calibration samples of D * + → D 0 (→ K − π + )π + decays. The PID efficiencies in the simulated sample are corrected by reweighting the calibration sample to match the momentum and the pseudorapidity distributions of the final state particles in the signal decay and the multiplicity of tracks in the event. The effective kaon and pion PID efficiencies in this analysis are approximately 90% and 80%, respectively. The rates for misidentifying a kaon as a pion or a pion as a kaon are less than 0.1%.
After all selection criteria have been applied, the ratio of the B + → K + K − π + to B + → K + K + π − selection efficiencies, weighted by the integrated luminosity taken with different magnet polarities and beam energy, is 1.00 ± 0.02, while the ratio of the B + → π + π − K + to B + → π + π + K − selection efficiencies is 0.97 ± 0.01. The quoted uncertainties are statistical only and are due to the limited simulation sample size. For each of the four modes, the selection efficiency across the Dalitz plane is constant with a relative variation of less than 9%.

Determination of the signal yields
Signal yields are determined from simultaneous, unbinned and extended maximum likelihood fits to the invariant mass distributions m h + h ± h ∓ of the suppressed and unsuppressed decays, where h and h denote π or K. Separate fits are made for B + → K + K + π − and B + → K + K − π + decays, and for B + → π + π + K − and B + → π + π − K + decays.
The signal component for both the suppressed and unsuppressed decays is parameterised as the sum of a Gaussian function and a Crystal Ball function [26]. The values of the signal function parameters are constrained to be the same for both the suppressed and unsuppressed signal components.
For all decay modes, the combinatorial background is parameterised with an exponential function. Partially reconstructed backgrounds, largely due to B decays with four or more particles in the final state, where one or more particles are not reconstructed, appear predominantly at m h + h ± h ∓ masses below 5150 MeV/c 2 and are modelled with an ARGUS function [27] convolved with a Gaussian resolution function. In the case of the suppressed B + → π + π + K − mode, an additional component, modelled in the same way, is used to account for B 0 s four-body decays such as B 0 s → D − s π + , with D − s → K − π + π − , where one of the decay particles is not reconstructed. These appear at m h + h ± h ∓ masses below 5250 MeV/c 2 . The slope parameter of this ARGUS function is fixed from a fit to simulated decays. The signal yields, background yields, and all signal and background parameters (apart from the fixed slope parameter of the B 0 s → D − s π + component) are allowed to float in the fit.
To investigate the presence of any peaking background in the signal mass region, a total of 350 million events from 200 B, B 0 s , Λ 0 b , and Ξ b decays are simulated and reconstructed using the same selection criteria used for the signal decay modes. In addition, the data are divided into two samples, above and below the signal mass region 5230 < m h + h ± h ∓ < 5320 MeV/c 2 . The Dalitz plot distributions in each sample are studied to identify any resonances that might be present in the signal region. There is no evidence for any peaking background in the signal mass region.
The performance of the fit procedure is tested by creating ensembles of simulated datasets generated from the functions fitted to the data. A large number of datasets is generated and fits are performed with the number of signal and background events left free to fluctuate according to a Poisson distribution. The fit biases on the signal yields extracted from the pseudoexperiments are 0.05 ± 0.10 and −0.48 ± 0.15 events for the B + → K + K + π − and B + → π + π + K − decays, respectively, where the uncertainty is statistical only. Figure 2 shows the fit to the B + → K + K + π − and B + → K + K − π + decay candidates. The signal yields are −7.2 ± 4.6 and 955 ± 75 for B + → K + K + π − and B + → K + K − π + respectively. Figure 3 shows the fit to the B + → π + π + K − and B + → π + π − K + decay candidates. The signal yields are 2.7 ± 10.0 and 24 044 ± 193 for B + → π + π + K − and B + → π + π − K + , respectively. The B + → K + K + π − and B + → π + π + K − signal yields have been corrected for the fit bias introduced by the fitting procedure.
The branching fractions of the suppressed decays are calculated using where N unsup sig and N sup sig are the numbers of fitted signal events for the unsuppressed and suppressed modes (corrected for fit bias), while unsup and sup are the selection efficiencies calculated from simulated events and corrected for differences in selection efficiency between simulation and data [28]. Finally, B unsup is the known branching fraction for the unsuppressed reference mode [24].

Systematic uncertainties
The measurements of the branching fractions of the suppressed modes depend on the ratios of suppressed to unsuppressed signal yields and selection efficiencies. Since the final state is the same, apart from the charge assignment, the ratio of unsuppressed and suppressed selection efficiencies is close to unity, and most of the systematic uncertainties cancel.
The main systematic uncertainties on the branching fractions are due to the unsuppressed branching fraction uncertainties, PID charge dependence, discrepancies in PID between data and simulated events, fit biases, alternative mass fit models, simulated event statistics, and assumptions concerning the Dalitz plot distributions.
Studies have shown a small PID dependence on the track charge with differences in efficiency below 0.2% for π ± and below 0.4% for K ± . Systematic uncertainties of 0.04×10 −8 are assigned to both the B + → K + K − π + and the B + → π + π − K + branching fractions, derived from a 0.2% systematic uncertainty per pion and 0.4% per kaon, added linearly. The PID efficiency is extracted from data and systematic uncertainties of 0.09 × 10 −8 and 0.11 × 10 −9 are applied to B + → K + K + π − and B + → π + π + K − , respectively, to account for different data-taking conditions.
The process used to determine the selection criteria for the BDT output and the six PID probabilities is validated by changing the order in which the PID criteria are optimized, adjusting the FOM for the predicted suppressed signal yield (using published branching fraction upper limits [6][7][8]) rather than the signal reconstruction efficiency , and reweighting the simulation samples to match PID distributions in data [28]. The process is shown to be robust and no systematic uncertainty is applied.
To allow for possible differences in reconstruction efficiency as a function of position in the Dalitz plane due to resonances and interference, simulated events are generated with a distribution of resonances taken from previous published results, which are however only available for the B + → π + π − K + decay [29,30]. The average values of the reconstruction efficiencies for the phase-space and resonance-allowed Dalitz plots differ by (8 ± 2)%, which is taken as the difference between a phase-space and a resonance-allowed distribution for all modes. The resulting uncertainties are 0.15 × 10 −8 and 0.22 × 10 −9 for B + → K + K + π − and B + → π + π + K − , respectively.
The fit bias from the ensemble of simulated pseudoexperiments generated from the fit to data is used to correct the signal yield observed in data. The systematic uncertainty on this correction is defined as half the fit bias added in quadrature with the fit bias statistical uncertainty. The systematic uncertainties introduced by this procedure are 0.05 × 10 −8 and 0.58 × 10 −9 for the B + → K + K + π − and B + → π + π + K − branching fractions, respectively.
To understand the impact of the fit model on the results, the components of the default models are changed independently. The signal component is replaced with a Crystal Ball function. The slope parameter of the B 0 s → D − s π + ARGUS function is allowed to float and the B 0 s → D − s π + component is replaced with a bifurcated Gaussian function. The combinatorial background distribution is modelled with a second order polynomial instead of an exponential. Studies of cross-feed from simulated unsuppressed B + → h + h − h + decays, where the flavour of one or more particles is misidentified, indicate 5.7 ± 2.7 events in the B + → K + K + π − decay mode and 3.1 ± 1.8 events in B + → π + π + K − decay mode. The uncertainty is dominated by limited simulation sample size. Cross-feed events are shifted by a minimum of ∼ 40 MeV/c 2 above and below the B mass and the majority of the events do not appear in the signal region. To confirm this, the fit is repeated with two additional Gaussian functions centred around 5240 MeV/c 2 and 5320 MeV/c 2 , respectively. As a cross-check, an additional fit is performed with the means and widths of the Gaussian component allowed to vary. The fitted yields are compatible with zero. Systematic uncertainties of 0.45 × 10 −8 and 0.29 × 10 −9 are assigned for the B + → K + K + π − and B + → π + π + K − decays, respectively.

Results and conclusions
Including all statistical and systematic uncertainties, the ratios of branching fractions are calculated to be where the first uncertainties are statistical and the second are systematic. Using the above and the world average of the unsuppressed branching fractions [24] and using Eq. 1, we calculate the branching fractions where the first uncertainties are statistical, the second are systematic without the unsuppressed decay branching fraction uncertainties, and the third are associated with the present knowledge of the B + → K + K − π + and B + → π + π − K + branching fractions.
To obtain upper limits on the branching fractions, the frequentist approach of Feldman and Cousins [32] is used to determine 90% and 95% confidence region bands that relate the true values of the branching fractions to the measured numbers of signal events. These bands are constructed using the results of simulation studies that account for relevant biases in the fit procedure and include statistical and systematic uncertainties. The construction of the confidence region bands is shown in Fig. 4. The 90% (95%) confidence level (CL) upper limits are found to be In summary, searches are presented for the highly-suppressed decays B + → K + K + π − and B + → π + π + K − using a data sample of 3.0 fb −1 collected by the LHCb experiment in proton-proton collisions at the centre-of-mass energies of 7 and 8 TeV. No evidence is found for these decays and upper limits are placed on B + → K + K + π − and B + → π + π + K − branching fractions. The results are approximately fourteen and twenty times more stringent, than previous measurements and constrain various extensions of the SM.