Evidence for the decay $B^0\to p\bar{p}\pi^0$

We report a search for the charmless baryonic decay $B^0\to p\bar{p}\pi^0$ with a data sample corresponding to an integrated luminosity of 711~$\rm fb^{-1}$ containing $(772\pm 10)\times 10^6$ $B^0\bar{B}^0$ pairs. The data was collected by the Belle experiment running on the $\Upsilon(4S)$ resonance at the KEKB $e^+e^-$ collider. We measure a branching fraction $\mathcal{B}(B^0\to p\bar{p}\pi^0)= (5.0\pm1.8\pm0.6 )\times 10^{-7}$, where the first uncertainty is statistical and the second is systematic. The signal has a significance of 3.1 standard deviations and constitutes the first evidence for this decay mode. We also search for the intermediate two-body decays $B^{0}\to\Delta^+\bar{p}$ and $B^0\to\bar{\Delta}^-p$, and set an upper limit on the branching fraction: $\mathcal{B}(B^0\to \Delta^+\bar{p})+\mathcal{B}(B^0\to\bar{\Delta}^-p)<1.6\times10^{-6}$ at 90% confidence level.

We report a search for the charmless baryonic decay B 0 → ppπ 0 with a data sample corresponding to an integrated luminosity of 711 fb −1 containing (772±10)×10 6 BB pairs. The data was collected by the Belle experiment running on the Υ (4S) resonance at the KEKB e + e − collider. We measure a branching fraction B(B 0 → ppπ 0 ) = (5.0 ± 1.8 ± 0.6) × 10 −7 , where the first uncertainty is statistical and the second is systematic. The signal has a significance of 3.1 standard deviations and constitutes the first evidence for this decay mode. We also search for the intermediate two-body decays B 0 → ∆ +p and B 0 →∆ − p, and set an upper limit on the branching fraction: B(B 0 → ∆ +p ) + B(B 0 →∆ − p) < 1.6 × 10 −6 at 90% confidence level.
PACS numbers: 13.25.Hw, 11.30.Er The first observed charmless baryonic B decay was B + → ppK + [1]. Following this first observation, many other charmless baryonic B decays have been found [2]. Except for B + → pΛπ 0 and pΛγ decays, all the channels reported to date are entirely reconstructed from charged particles in the final state. A noticeable hierarchy is also observed in the branching fractions of these decays: three-body decays are usually more frequent than their two-body counterparts but less frequent than four-body decays [3,4]. This phenomenon can be understood in terms of the so-called "threshold effect," which refers to the fact that the B meson prefers to decay into a dibaryon pair with low invariant mass accompanied by a fast recoil meson [3,5,6]. This peaking behavior was unexpected, and has led to various speculations about possible mechanisms [3]. Studying additional three-body baryonic decays might provide a better understanding of the dynamics of B decays and the aforementioned threshold-effect. These decays are also useful for CP violation studies.
This paper reports a search for a three-body charmless baryonic B 0 decays to the ppπ 0 final state [7] using a data set corresponding to an integrated luminosity of 711 fb −1 collected with the Belle detector [8] at the Υ (4S) resonance at the KEKB asymmetric-energy e + e − (3.5 on 8.0 GeV) collider [9]. So far, the decay B 0 → ppπ 0 has not been studied by any experiment. No theoretical prediction for the branching fraction of this process is yet available. A glance at the known branching fractions for B decays [2] shows the three-body charmless baryonic decays to occur in the several times 10 −6 range, indicating that the discovery of the mode B 0 → ppπ 0 might be possible with the currently available data set.
The Belle detector is a large-solid-angle magnetic spectrometer consisting of a silicon vertex detector (SVD), a 50-layer central drift chamber (CDC), an array of aerogel threshold Cherenkov counters (ACC), a barrel-like arrangement of time-of-flight scintillation counters (TOF), and an electromagnetic calorimeter (ECL) comprising CsI(Tl) crystals. These detector components are located inside a superconducting solenoid coil that provides a 1.5 T magnetic field. An iron flux-return located outside the coil is instrumented (KLM) to detect K 0 L mesons and to identify muons. Two inner detector configurations were used: a 2.0 cm radius beampipe and a three-layer SVD were used for the first 152 × 10 6 BB pairs of data, while a 1.5 cm radius beampipe, a four-layer SVD, and a small-cell inner drift chamber were used for the remaining 620×10 6 BB pairs of data. The detector is described in detail in Ref. [8]. Event selection requirements are optimized using Monte Carlo (MC) simulations. MC events are generated using EvtGen [10], and the detector response is modeled using Geant3 [11]. Final-state radiation is taken into account using the Photos package [12].
The reconstruction of B 0 → ppπ 0 proceeds by first reconstructing π 0 → γγ candidates. An ECL cluster not matched to any track in the CDC is identified as a photon candidate. Such candidates are required to have an energy greater than 50 MeV in the barrel region and greater than 100 MeV in the end-cap regions, where the barrel region covers the polar angle 32 • < θ < 130 • and the end-cap regions cover the ranges 12 • < θ < 32 • and 130 • < θ < 157 • . To reject showers produced by neutral hadrons, the energy deposited in the 3 × 3 array of ECL crystals centered on the crystal with the highest energy must exceed 80% of the energy deposited in the corresponding 5 × 5 array of crystals. We require that the γγ invariant mass be within 20 MeV/c 2 (about 3.5σ in resolution) of the π 0 mass [2]. To improve the π 0 momentum resolution, we perform a mass-constrained fit and require that the resulting χ 2 be less than 30. This requirement is relatively loose, retaining more than 99% of candidates.
We subsequently combine π 0 candidates with two oppositely charged tracks, identified as a proton-antiproton pair. Such tracks are identified using requirements on the distance of closest approach with respect to the interaction point along the z axis (antiparallel to the e + beam) of |dz| < 3.0 cm, and in the transverse plane of dr < 0.3 cm. In addition, charged tracks are required to have a minimum number of SVD hits (> 2 along the z axis and > 1 in the transverse direction). Particle identification is achieved using information from the CDC, the TOF, and the ACC subdetectors. This information is combined to form a hadron likelihood L h ; a charged track with likelihood ratios of L p /(L p + L K ) > 0.9 and L p /(L p +L π ) > 0.9 is regarded as a proton or antiproton. Furthermore, we reject tracks consistent with either the electron or muon hypothesis. The proton identification efficiency is 75% and the probability for a kaon (pion) to be misidentified as a proton is 6% (2%).
Candidate B 0 mesons are identified using the beamenergy-constrained mass, M bc = E 2 beam − | p B c| 2 /c 2 , and the energy difference ∆E = E B − E beam , where E beam is the beam energy, and E B and p B are the reconstructed energy and momentum, respectively, of the B 0 candidate. All quantities are evaluated in the center-of-mass (CM) frame. To improve the M bc resolution, the momentum p B is calculated as where m π 0 is the nominal π 0 mass [2], E h and p h are the energy and momentum of the hadron h (h = p,p, π 0 ). In addition, a vertex fit is performed to the charged tracks to form a B 0 vertex. We require that the χ 2 from the fit be less than 200. Events with M bc > 5.25 GeV/c 2 and −0.20 GeV < ∆E < 0.15 GeV are retained for further analysis. The signal yield is calculated in a smaller region M bc ∈ (5.272, 5.286) GeV/c 2 and ∆E ∈ (−0.12, +0.06) GeV. In order to reject contributions from charmonium states (e.g., η c , J/ψ , ψ(2S), χ c0 , χ c1 and χ c2 ), we apply a "charmonium veto" and exclude the regions of 2.850 GeV/c 2 < m(pp) < 3.128 GeV/c 2 and 3.315 GeV/c 2 < m(pp) < 3.735 GeV/c 2 from the event sample.
Charmless hadronic decays suffer from large amount of continuum background, arising from light quark production (e + e − → qq, q = u, d, s, c). To suppress this background, we use a multivariate analyzer based on a neural network (NN) [13] that distinguishes jet-like continuum events from more spherical BB events. The NN uses the following input variables: the cosine of the angle between the thrust axis [14] of the B 0 candidate and the thrust axis of the rest of the event; the cosine of the angle between the B 0 thrust axis and the +z axis; the cosine of the angle between the +z axis and the B 0 candidate flight direction; a set of 18 modified Fox-Wolfram moments [15]; the ratio of the second to zeroth (unmodified) Fox-Wolfram moments; the separation along the z axis between the two B vertices; and the B-flavor tagging information [16]. All but for the last two quantities are evaluated in the CM frame. The NN is trained using MC simulated signal events and qq background events. The NN generates a single output variable (C NN ) that ranges from −1 for background-like events to +1 for signal-like events. We require C NN > −0.5, which rejects approximately 86% of the qq background while retaining 94% of the signal. We then translate C NN to a new variable where C min NN = −0.5 and C max NN = 1.0. This translation is advantageous as the C NN distribution for both signal and background is well described by a sum of Gaussian functions.
After applying all selection criteria, approximately 7% of the events have multiple B 0 candidates. For these events, we retain the candidate having the smallest sum of χ 2 values obtained from the π 0 → γγ mass-constrained fit and the B 0 vertex-constrained fit. According to MC simulation, this criterion selects the correct B 0 candidate in 83% of multiple-candidate events.
We measure the signal yield by performing an unbinned extended maximum likelihood fit to the variables M bc , ∆E, and C NN . The likelihood function is defined where Y j is the yield of component j; P j (M i bc , ∆E i , C i NN ) is the probability density function (PDF) of component j for event i; j runs over all signal and background components; and i runs over all events in the sample (N ). The background components consist of continuum events, b → c (generic B) processes, and rare charmless processes. The latter two backgrounds are small compared to the continuum events and are studied using MC simulations. The rare charmless background shows a peaking structure in the M bc distribution, most of which arises from B + → ppρ + decays. As correlations among the variables M bc , ∆E, and C NN are found to be small, the three-dimensional PDFs P j (M i bc , ∆E i , C i NN ) are factorized into the product of separate one-dimensional PDFs.
The PDF of signal events consists of two parts: one for candidates that are correctly reconstructed, and one for those incorrectly reconstructed, i.e., at least one daughter originates from the other (tag-side) B. For the former case, the M bc and ∆E distributions are modeled with Gaussian and Crystal Ball (CB) [17] functions, respectively, while the C NN distribution is modeled with a sum of Gaussian and bifurcated Gaussian functions having a common mean. The peak positions and resolutions of the M bc , ∆E, and C NN PDFs are adjusted to account for data-MC differences observed in a high-statistics control sample of B 0 → D 0 (→ K + π − )π 0 decays. For the latter case, the correlated two-dimensional M bc -∆E distribution is modeled with a non-parametric PDF [18], and the C NN component is modeled with a Gaussian function. The fraction of incorrectly reconstructed decays (∼ 4% in the signal region) is taken from MC simulation. For the rare charmless background, the C NN component is modeled with a bifurcated Gaussian function. The M bc and ∆E components are modeled by a joint two-dimensional non-parametric PDF. We model the M bc , ∆E, and C NN distributions of continuum background with an ARGUS [19] function having its endpoint fixed to 5.29 GeV/c 2 , a first-order polynomial, and a sum of two Gaussians having a common mean, respectively. For the generic B background, we use a bifurcated Gaussian function to model the C NN shape, while the similar shapes as of continuum background are used to model the M bc and ∆E distributions. In addition to the fitted yields Y j , all shape parameters for continuum background are also floated. All other parameters are fixed to the corresponding MC values.
The projections of the fit are shown in Fig. 1. From the fit, we extract 40.5 ± 14.2 signal events, 1490.3 ± 34.5 continuum, 100.6 ± 35.0 generic B, and 6.5 ± 10.1 rare charmless background events in the M bc − ∆E signal region. The resulting branching fraction is calculated as where Y sig represents the extracted signal yield, N BB = (772 ± 11) × 10 6 is the total number of BB events, ε = (10.53 ± 0.04)% is the reconstruction efficiency. The efficiency is corrected to account for possible differences in particle identification (PID) and π 0 detection efficiencies between data and simulations. In Eq. (3) we assume equal production of B 0 B 0 and B + B − pairs at the Υ (4S) resonance. The result is where the first uncertainty is statistical and the second is systematic. This is the first measurement of this branching fraction. The signal significance is calculated as −2 ln(L 0 /L max ), where L 0 is the likelihood value when the signal yield is fixed to zero, and L max is the likelihood value of the nominal fit. To include systematic uncertainties in the significance, we convolve the likelihood distribution with a Gaussian function whose width is set to the total systematic uncertainty that affects the signal yield. The resulting significance is 3.1 standard deviations. Thus, our measurement constitutes the first evidence for this decay mode.
The systematic uncertainty in B(B 0 → ppπ 0 ) arises from several sources, as listed in Table I. The uncertainty due to the fixed parameters in the PDF is estimated by varying them individually according to their statistical uncertainties. For each variation, the branching fraction is recalculated, and the difference with the nominal value is taken as the systematic uncertainty associated with that parameter. The smoothing parameters of the non-parametric functions are also varied. The differences in the fit results are included as systematic uncertainties. We add all uncertainties in quadrature to obtain the overall uncertainty due to PDF parametrization. The uncertainties due to errors in the calibration factors used to account for data-MC differences in the signal PDF are evaluated separately but in a similar manner. To test the stability of our fitting procedure, we generate and fit a large ensemble of pseudoexperiments. We find a potential fit bias of +2.1%. We attribute this bias to neglecting small correlations among the fitted observables. We assign a 1.5% systematic uncertainty due to π 0 reconstruction; this is determined from a study of τ − → π − π 0 ν τ decays [20]. The systematic uncertainty due to the track reconstruction efficiency is 0.35% per track, as determined from a study of partially reconstructed D * + → D 0 π + , D 0 → K 0 S π + π − decays. A 0.6% systematic uncertainty is assigned due to the particle identification efficiency of the proton-antiproton pair; this is determined from a study of Λ → pπ − decays. We determine the systematic uncertainty due to the C NN selection by applying different C NN criteria and comparing the results with that of the C NN nominal selection. The uncertainty due to the estimated fraction of incorrectly reconstructed signal events is obtained by varying this fraction by ±50%. The systematic uncertainty due to the total number of BB pairs is 1.4%, and the uncertainty due to MC used to evaluate the reconstructed efficiency is 0.4%. The total systematic uncertainty is obtained by adding each source in quadrature, as they are assumed to be uncorrelated.  Figure 2 shows the background-subtracted and efficiency-corrected distribution of m(pp), where the charmonium veto is removed. For the background subtraction, we use the sPlot technique [21], with M bc , ∆E and C NN as the discriminating variables. As expected, an

tions.
We also search for the intermediate two-body decay B 0 → ∆ + (→ pπ 0 )p. Events with m(pπ 0 ) < 1.4 GeV/c 2 are selected for this search. No significant signal is observed in this mass range. We set an upper limit on the branching fraction of B 0 → ∆ +p at 90% confidence level (C.L.) using a Bayesian approach. The limit is obtained by integrating the likelihood function from zero to infinity; the value that corresponds to 90% of this total area is taken as the 90% C.L. upper limit. We include the systematic uncertainty in the calculation by convolving the likelihood distribution with a Gaussian function whose width is set equal to the total systematic uncertainty of B(B 0 → ppπ 0 ). As we do not know the flavor of the B meson at decay, we express our result as a sum of final states containing either a ∆ + or a∆ − . The result is This is the first such limit and is in agreement with the theoretical predictions [22,23].
In summary, using the full set of Belle data, we report a measurement of the branching fraction for B 0 → ppπ 0 decays. We obtain B(B 0 → ppπ 0 ) = (5.0±1.8±0.6)×10 −7 , where the first uncertainty is statistical and the second is systematic. The significance of this result is 3.1 standard deviations, and thus this measurement constitutes the first evidence for this decay. We also search for the intermediate two-body decays B 0 → ∆ +p and B 0 →∆ − p, and set an upper limit on the branching fraction, B(B 0 → ∆ +p ) + B(B 0 →∆ − p) < 1.6 × 10 −6 at 90% C.L.
We thank the KEKB group for the excellent operation of the accelerator; the KEK cryogenics group for the efficient operation of the solenoid; and the KEK computer group, and the Pacific Northwest National Laboratory (PNNL) Environmental Molecular Sciences Laboratory (EMSL) computing group for strong computing support; and the National Institute of Informatics, and Science Information NETwork 5 (SINET5) for valuable network support. We acknowledge support from the Ministry of Education, Culture, Sports, Science, and Technology (MEXT) of Japan, the Japan Society for the Promotion of Science (JSPS), and the Tau