Search for the decay mode $B^0\rightarrow p p \bar{p} \bar{p}$

A search is presented for the four-body decay $B^0\rightarrow p p \bar{p} \bar{p}$ in a sample of $471$ million $B\bar{B}$ pairs collected with the BaBar detector, operated at the SLAC PEP-II asymmetric-energy $e^+e^-$ collider. The center-of-mass energy is $10.58$ GeV. From a fit to the distribution of the energy-substituted mass $m_{\mathrm{ES}}$, the result $\mathcal{B}(B^0\rightarrow p p \bar{p} \bar{p})=(1.1 \pm 0.5 \pm 0.2) \times 10^{-7}$ is extracted, where the first uncertainty is statistical and the second is systematic. The significance of the signal is $2.9$ standard deviations. The upper limit on the branching fraction is determined to be $2.0 \times 10^{-7}$ at $90\%$ confidence level.

collected with the BABAR detector, operated at the SLAC PEP-II asymmetric-energy e + e − collider.The center-of-mass energy is 10.58 GeV.From a fit to the distribution of the energy-substituted mass mES, the result B(B 0 → pppp) = (1.1±0.5±0.2)×10−7 is extracted, where the first uncertainty is statistical and the second is systematic.The significance of the signal is 2.9 standard deviations.The upper limit on the branching fraction is determined to be 2.0 × 10 −7 at 90% confidence level.
PACS numbers: 13.25.Hw, 13.60.Rj, 14.40.Nd The inclusive branching fraction of B mesons decaying into final states with at least one baryon-antibaryon pair is approximately 7% [1], while the sum of all measurements of exclusive baryonic channels is less than 1% [2].Recent measurements from the LHCb experiment [3][4][5] have raised new interest in this field.Studying exclusive baryonic decays of B mesons provides a deeper insight into the mechanism of hadronization into baryons and may allow a better understanding of the threshold enhancement effect, which is a dynamical enhancement, relative to the pure phase space expectation, of the production rate of baryon-antibaryon pairs at their invariant mass threshold.So far, this process is only qualitatively understood.Theoretical models, e.g., the quantum chromodynamics (QCD) sum rule [6] and the perturbative QCD approach [7], need validation and input from experimental data.Although the threshold effect is also observed in B decays to charmed baryons [8], its effect is not as pronounced as in charmless three-body baryonic decays, where the peak at the threshold of the invariant baryon-antibaryon mass distribution was first observed [9,10].This enhancement could explain the hierarchy trend of the branching fractions for baryonic B decays.It has been observed that many three-body final states have larger rates than their two-body counterparts, and also that some three-body decays are suppressed compared to the four-body case [11][12][13].The phenomenological approaches describe these observations in terms of gluonic and fragmentation mechanisms [14] and pole models [15].For final states with a pp pair, a threshold enhancement could possibly arise from an intermediate X(1835) baryonium resonance, as proposed in Ref. [16].
In this paper, we report on a search for B 0 → pppp decays (the inclusion of charge conjugate processes is implied throughout).The data were collected with the BABAR detector [17,18] at the SLAC PEP-II asymmetricenergy e + e − collider.The decay of a B meson to two baryon-antibaryon pairs has not yet been observed.No quantitative predictions for this process are yet available.The measurement of a four baryon decay mode would provide useful information to help discriminate between existing models and aid in the development of phenomenological models for four-baryon production.Previously, we performed a search for B 0 → Λ + c ppp decays, setting a 90% confidence level (C.L.) upper limit on the decay branching fraction of 2.8×10 −6 [19].Based on this result, using a scaling factor to account for the Cabibbo suppression for the b → u decay, and also taking into account the larger phase space of the final state, we estimate the branching fraction for the B 0 → pppp decay mode to be on the order of 10 −7 .We use this assumption to optimize the selection criteria.The threshold effect has been found to be enhanced for large values of q, which is the available momentum in the rest frame of the decaying B, with no visible enhancement for q values below about 200 MeV/c [20].This feature could explain the absence of an observed signal in B 0 → Λ + c ppp decays and at the same time could enhance the branching fraction for B 0 → pppp.Moreover, the B 0 → pppp decay rate may benefit from the low-invariant-mass enhancement of the double pp system and from presence of nontrivial intermediate bound states [21].
The analysis is based on the full data set collected with BABAR at center-of-mass energy 10.58 GeV, corresponding to the peak of the Υ (4S) resonance.The event sample contains N BB = 471 × 10 6 BB pairs, corresponding to integrated luminosity of 424 fb −1 [22].Chargedparticle momenta are measured by means of a five-layer double-sided silicon vertex tracker and a 40-layer multiwire drift chamber, both operating in the 1.5 T magnetic field of a superconducting solenoid.The particle identification (PID) for protons, kaons, and pions uses the specific energy loss measured in the tracking devices and the measurement of the Cherenkov angle provided by the internally reflecting, ring-imaging Cherenkov detector.We use Monte Carlo (MC) simulated events of the processes e + e − → B B, where the B mesons decay generically according to known branching fractions and decay amplitudes [2], and e + e − → q q (with q = u, d, s, c) to model the background.These samples correspond to at least three times the integrated luminosity of the data.In addition, we generate a sample of 687 000 signal decays e + e − → B 0 B 0 , where one of the B mesons decays into pppp (referred to as the signal MC sample).Monte Carlo events are simulated with the EvtGen and Jetset [23,24] event generators, with the response of the detector simulated using the Geant4 suite of programs [25].Signal and background MC samples are used for the signal efficiency determination and for the modeling of the signal and background distributions.
A B meson candidate is reconstructed by combining four charged tracks, two identified as protons and two as antiprotons, kinematically fitting them to a common vertex and requiring the fit probability to exceed 0.1%.The direction of the reconstructed B meson is required to originate from the interaction region, which is constrained to the beam-spot size in the laboratory frame.Tracks are rejected if the combination of two oppositely charged tracks is found to be consistent with K 0 S or Λ hypotheses.
Loose preselection requirements are applied to the kinematic variables [20] 2 GeV, where p * B and E * B are, respectively, the momentum and energy of the reconstructed B candidate in the CM frame, and E * beam is half the CM energy.The study is performed as a blind analysis, which means that the selection is optimized without examining the data in the signal region, 5.27 < m ES < 5.29 GeV/c 2 .
The PID efficiency for protons is larger than 99% and the misidentification of kaons and pions as protons is less than 1%.The difference between the PID performance in data and simulation is evaluated using events from high-purity channels, which form the control samples (CS) for a given particle type.For example, events with Λ 0 → pπ − decays form the CS for the validation of proton PID, K 0 S → π + π − for pion PID, and D * + → π + D 0 (D 0 → K − π + ) for kaon PID.The PID efficiency of the MC-simulated events is corrected to match that observed in data by applying the weight CS,Data / CS,MC , where CS,Data and CS,MC are the PID efficiencies evaluated from the CS in data and simulation, respectively.
After applying the particle identification and preselection requirements, the fraction of misidentified signal candidates in simulation is found to be negligible (< 0.2%).The main background is combinatorial, from genuine protons in continuum (e + e − → q q) events.The continuum background is further reduced by imposing a signal-like selection on the output of a multivariate boosted decision tree (BDT) algorithm.The BDT classifier uses the following input variables: ∆E, cos θ * B , with θ * B the polar angle of the B meson candidate with respect to the beam axis in the CM frame, and the event shape variables R 2 and | cos θ TH |, where R 2 is the ratio of the second to the zeroth Fox Wolfram moments [26] and θ TH is the angle between the thrust axis [27] of the B candidate and that of the rest of the event in the Υ (4S) rest frame.These kinematic and topological variables are effective in discriminating between spherically shaped events from B B decays and jet-like q q events.In the BDT output, signal (background) events peak at positive (negative) values (Figure 1).The optimal selection on the BDT output is determined by maximizing the figure of merit S/ √ S + B, where S and B are the number of expected signal and background events, respectively.The number of signal events is estimated assuming the signal branching fraction of 10 −7 mentioned above.The distributions of the input variables, before applying the BDT selection, are shown in Figure 2, where the signal and background MC samples have been nor- events.The two background components, from continuum and generic BB events, are stacked, with the total background prediction scaled to correspond to the number of selected events in the data.For purposes of visibility, the signal distribution has been multiplied by a factor of 100.The selection on the BDT output is indicated by the black vertical line.
malized to match the number of selected events in data.The selection is optimized using the MC samples and is validated by comparing the distributions for the background MC samples to the data in the control region m ES < 5.27 GeV/c 2 .Good MC-data agreement is observed.
The total number of selected data events in the interval 5.2 < m ES < 5.3 GeV/c 2 is 117.The signal efficiency, evaluated from simulation, is found to be = 0.2068 ± 0.0004 (stat).
We investigate the potential presence of peaking backgrounds from the baryonic modes B → pph + h − , recently measured by the LHCb Collaboration [4], where h is a hadron other than a proton.These decays can potentially enter the background if the h + h − pair is erroneously identified as a pp pair.This background is evaluated by applying the event selection to the simulated MC samples for the modes reported in Table I and determing the selection efficiencies p ph + h − .The number of expected background events in the data for each channel is estimated as N p ph + h − = p ph + h − B N B B (Table I), where B is the branching fraction measured in Ref. [4] and N B B is the total number of BB pairs in the initial data sample.The expected contamination from these decays is found to be negligible.To describe the m ES distribution in data, we use a probability density function (PDF) corresponding to the sum of the signal and the background components.The signal component is described by a Gaussian function, whose mean and width are fixed to values determined from a fit to simulated signal events.The combinatorial background component is described FIG. 2. Comparison between data (dots) and MC (stacked histograms) for the BDT input variable distributions, before the BDT selection is applied.The total MC prediction for backgrounds has been scaled to correspond to the number of selected events in the data before the BDT selection, while the signal MC (shaded histogram) is multiplied by a factor of 100 for better visibility.by the empirical ARGUS function [28], which depends on two parameters: a shape parameter and a cutoff parameter.The cutoff parameter is set equal to the endpoint in the m ES spectrum, 5.289 GeV/c 2 .The shape parameter is determined in the fit, along with the signal and background event yields, N sig and N bkg , respectively.
The signal yield is extracted by performing an unbinned extended maximum likelihood fit to the m ES distribution in the range 5.2 < m ES < 5.3 GeV/c 2 (Figure 3).The logarithm of the extended likelihood is written as: where x corresponds to the measured m ES distribution and the f j (x) are the corresponding PDFs for the signal and background components.The sum of the signal and background yields is constrained to the total number of observed events n.The result is N sig = 11.1 ± 4.6 (stat) events, from which the corresponding branching fraction is calculated as: where the uncertainty is statistical only and we assume equal production of B 0 B 0 and B + B − in Υ (4S) decays.Therefore, 2 N B 0 B 0 = N BB , where N BB is the total number of BB pairs in the initial data sample.The value of 2 takes into account that the charge-conjugate decay is also reconstructed.The experimental values for N sig and N BB are listed in Table II.
To evaluate the statistical significance of the branching fraction result, we fit the data under a backgroundonly hypothesis and determine the corresponding change −2(∆ ln L) with respect to the standard fit, where L is the likelihood function.The statistical significance is found to be 2.9 standard deviations.The systematic uncertainty in the ARGUS cutoff is taken into account and is found to not affect the signal significance.
Systematic uncertainties in the branching fraction arise from the uncertainty in N BB , from the fit procedure, and from the uncertainty in the signal efficiency.The relative systematic uncertainties for the considered sources are listed in Table III.The systematic uncertainty in N BB is estimated to be 0.6% [29].
Potential systematic uncertainties associated with the fit procedure arise from the choice made for the signal PDF shape and from variation of the parameters held constant in the fit.Variations of the form chosen to model the signal PDF are found to have a negligible impact on the result, while the uncertainty associated with the ARGUS cutoff value, evaluated by varying the cutoff value within its uncertainty of 0.5 MeV/c 2 , is 0.9% To determine the systematic uncertainty in the signal efficiency, several sources are taken into account: the statistical uncertainty from the MC samples, the PID performance, the track finding efficiency, the BDT method, and the decay model used for the generation of the signal MC sample.The finite size of the signal MC sample results in a relative systematic uncertainty of 0.2%.The PID performance contribution is taken as the effect of the full data-to-MC correction mentioned above and corresponds to a relative uncertainty of 0.9%.The systematic uncertainty related to the track finding efficiency is a function of the particle momentum [30] and amounts to 0.9% for protons of approximately 1 GeV/c momentum [30].The systematic uncertainty in the signal efficiency introduced by the BDT method is evaluated by reweighting, separately for each of the four input variables of the BDT classifier, the shape of the MC distribution to match that observed in data.The weights are calculated in the control region m ES < 5.27 GeV/c 2 , before the BDT selection, as the bin-by-bin ratio between data and MC events, and are applied to the corresponding distribution of the signal MC sample.The difference in the efficiency computed with and without the weights applied provides a systematic uncertainty of 2.2%.The systematic uncertainty related to the unknown dynamics of B 0 → pppp decays is evaluated by comparing the pure phase space decay MC sample to a model in which the decay proceeds through an intermediate spinless resonance, B → X(→ pp)X(→ pp).Weights, binned in the 4-dimensional space of the magnitudes of the momenta of the four tracks, are obtained by dividing the momentum distribution resulting from the resonant model by that from the phase space model, and are applied to the proton momentum distribution of the signal MC.The systematic uncertainty is obtained from the difference in the efficiency computed from the weighted and unweighted samples.The largest relative difference is obtained for an X mass of 2.6 GeV/c 2 .It amounts to 14% and is the largest contribution to the total systematic uncertainty.
The upper limit at 90% C.L. on the signal yield is computed by integrating the likelihood as a function of N sig , up to the value N UL sig such that the equality L(N sig )dN sig is satisfied.This calculation is based on the Bayesian approach, assuming a uniform prior for N sig > 0 and 0 otherwise.The 90% C.L. upper limit on N sig is N UL sig = 19 events, corresponding to a 90% C.L. upper limit on the signal branching fraction 2.0 × 10 −7 .
We use pseudoexperiments to establish the robustness of our result for the upper limit against systematic variation.The signal yield N sig is varied according to the signal PDF computed from the fitted likelihood function and the signal efficiency is randomly smeared ac-cording to a Gaussian distribution with mean 0.2068 and with width 0.03, corresponding to its absolute systematic uncertainty.For each pseudoexperiment, the branching fraction is calculated from the input values for N sig and randomly selected from the above defined PDFs and its distribution is integrated up to 90% C.L., which shows a smeared upper limit consistent within its uncertainty with the unsmeared result.
In summary, we have performed a search for B meson decays to the pppp final state, obtaining 11.1 signal events.The statistical significance of the result is 2.9 standard deviations.The branching fraction is measured to be B = (1.1 ± 0.5(stat) ± 0.2(syst)) × 10 −7 .The corresponding 90% C.L. upper limit is B(B 0 → pppp) < 2.0 × 10 −7 .Our result can provide important input for QCD models of hadronization and improve understanding of the threshold enhancement effect.
We are grateful for the extraordinary contributions of our PEP-II colleagues in achieving the excellent luminosity and machine conditions that have made this work possible.The success of this project also relies critically on the expertise and dedication of the computing organizations that support BABAR.The collaborating institutions wish to thank SLAC for its support and the kind hospitality extended to them.This work is supported by the US Department of Energy and National Science Foundation, the Natural Sciences and Engineering Research Council (Canada), the Commissariat à l'Energie Atomique and Institut National de Physique Nucléaire et de Physique des Particules (France), the Bundesministerium für Bildung und Forschung and Deutsche Forschungsgemeinschaft (Germany), the Istituto Nazionale di Fisica Nucleare (Italy), the Foundation for Fundamental Research on Matter (The Netherlands), the Research Council of Norway, the Ministry of Education and Science of the Russian Federation, Ministerio de Economía y Competitividad (Spain), the Science and Technology Facilities Council (United Kingdom), and the Binational Science Foundation (U.S.-Israel).Individuals have received support from the Marie-Curie IEF program (European Union) and the A. P. Sloan Foundation (USA).

FIG. 1 .
FIG.1.The BDT output distribution for simulated signal (shaded histogram) and background (filled histogram) events.The two background components, from continuum and generic BB events, are stacked, with the total background prediction scaled to correspond to the number of selected events in the data.For purposes of visibility, the signal distribution has been multiplied by a factor of 100.The selection on the BDT output is indicated by the black vertical line.

FIG. 3 .
FIG.3.Fit to the data mES distribution (dots) in the interval 5.2 < mES < 5.3 GeV/c 2 .The bottom plot shows the pull distribution, which is the bin-by-bin difference between the data and fitted distribution normalized by the corresponding statistical uncertainty from the fit.

TABLE II .
Experimental inputs used for the branching fraction calculation.

TABLE III .
Relative systematic uncertainties in the signal branching fraction.The total systematic uncertainty is determined by summing the individual contributions in quadrature.