Observation of the B 0 → ρ 0 ρ 0 decay from an amplitude analysis of B 0 → ( π + π − )( π + π − ) decays

,


Introduction
The study of B meson decays to ρρ final states provides the most powerful constraint to date for the Cabibbo-Kobayashi-Maskawa (CKM) angle α ≡ arg [(V td V * tb )/(V ud V * ub )] [1][2][3].Most of the physics information is provided by the decay B 0 → ρ + ρ − as measured at the e + e − colliders at the Υ(4S) resonance [4,5], 1 for which the dominant decay amplitude, involving the emission of a W boson only (tree), exhibits a phase difference that can be interpreted as the sum of the CKM angles β + γ = π − α in the Standard Model.The subleading amplitude associated with the exchange of a W boson and a quark (penguin) must be determined in order to interpret the electroweak phase difference in terms of the angle α.This is realised by means of an isospin analysis involving the companion modes B + → ρ + ρ 0 [6,7] and B 0 → ρ 0 ρ 0 [8,9]. 2 In particular, the smallness of the amplitude of the latter leads to a better constraint on α.
The BaBar and Belle experiments reported evidence for the B 0 → ρ 0 ρ 0 decay [8,9] with an average branching fraction of B(B 0 → ρ 0 ρ 0 ) = (0.97 ± 0.24) × 10 −6 [8,9].Despite small observed signal yields, each experiment measured the fraction f L of decays yielding a longitudinally polarised final state through an angular analysis.The Belle collaboration did not find evidence for polarisation, f L = 0.21 +0.22  −0.26 [9], while the BaBar experiment measured a mostly longitudinally polarised decay, f L = 0.75 +0.12  −0.15 [8].These results differ at the level of 2.0 standard deviations.The large LHCb data set may shed light on this discrepancy.In addition, LHCb may confirm the hint of B 0 → ρ 0 f 0 (980) decays reported by Belle [9].Measurements of the B 0 → ρ 0 ρ 0 branching fraction and longitudinal polarisation fraction at LHCb can be used as inputs in the determination of α [2,3].
This work focuses on the search and study of the B 0 → (π + π − )(π + π − ) decay in which the two (π + π − ) pairs are selected in the low invariant mass range (< 1100 MeV/c 2 ).The B 0 → ρ 0 ρ 0 is expected to dominate the (π + π − ) mass spectrum.The (π + π − ) combinations can actually emerge from S-wave non-resonant and resonant contributions or other P-or D-wave resonances interfering with the signal.Hence, the determination of the B 0 → ρ 0 ρ 0 yields requires a two-body mass and angular analysis, from which the fraction of the longitudinally polarised final state can be measured.
The branching fraction is measured relative to the B 0 → φK * (892) 0 mode.The B 0 → φK * (892) 0 decay, which results in four light mesons in the final state, is similar to the signal, thus allowing for a cancellation of the uncertainties in the ratio of selection efficiencies.

Data sets and selection requirements
The analysed data correspond to an integrated luminosity of 1.0 fb −1 and 2.0 fb −1 from pp collisions recorded at a centre-of-mass energy of 7 TeV, collected in 2011, and 8 TeV, collected in 2012, by the LHCb experiment at CERN.The LHCb detector [10,11] is a single-arm forward spectrometer covering the pseudorapidity range 2 < η < 5, designed for the study of particles containing b or c quarks.It includes a high-precision tracking system consisting of a silicon-strip vertex detector surrounding the pp interaction region [12], 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 [13] placed downstream of the magnet.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 vertex, the impact parameter, 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 [14].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 [15].The online event selection is performed by a trigger [16], 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.
In this analysis two categories of events that pass the hardware trigger stage are considered: those where the trigger decision is satisfied by the signal b-hadron decay products (TOS) and those where only the other activity in the event determines the trigger decision (TIS).The software trigger requires a two-, three-or four-track secondary vertex with large transverse momenta of charged particles and a significant displacement from the primary pp interaction vertices (PVs).At least one charged particle should have p T > 1.7 GeV/c and is required to be inconsistent with originating from any primary interaction.A multivariate algorithm [17] is used for the identification of secondary vertices consistent with the decay of a b hadron.
Further selection criteria are applied offline to reduce the number of background events with respect to the signal.The (π + π − ) candidates must have transverse momentum larger than 600 MeV/c, with at least one charged decay product with p T > 1000 MeV/c.The two (π + π − ) pairs are then combined to form a B 0 candidate with a good vertex quality and transverse momentum larger than 2500 MeV/c.The invariant mass of each pair of oppositecharge pions forming the B 0 candidate is required to be in the range 300-1100 MeV/c 2 .The identification of the final-state particles (PID) is performed with dedicated neuralnetworks-based discriminating variables that combine information from the RICH detectors and other properties of the event [14].The combinatorial background is further suppressed with multivariate discriminators based on a boosted decision tree algorithm (BDT) [18,19].The BDT is trained with simulated B 0 → ρ 0 ρ 0 (where ρ 0 → π + π − ) events as signal sample and candidates reconstructed with four-body mass in excess of 5420 MeV/c 2 as background sample.The discriminating variables are based on the kinematics of the B decay candidate (B p T and the minimum p T of the two ρ 0 candidates) and on geometrical vertex measurements (quality of the B candidate vertex, impact parameter significances of the daughters, B flight distance significance and B pointing to the primary vertex).
The optimal thresholds for the BDT and PID discriminating variables are determined simultaneously by means of a frequentist estimator for which no hypothesis on the signal yield is assumed [20].The B 0 meson candidates are accepted in the mass range 5050-5500 MeV/c 2 .
The normalisation mode B 0 → φK * (892) 0 is selected with similar criteria, requiring in addition that the invariant mass of the (K + π − ) candidate is found in a range of ±150 MeV/c 2 around the known value of the K * (892) 0 meson mass [21] and the invariant mass of the (K + K − ) pair is in a range of ±15 MeV/c 2 centred at the known value of the φ meson mass [21].A sample enriched in B 0 → (K + π − )(π + π − ) events is selected using the same ranges in (π + π − ) and (K + π − ) masses to estimate the background with one misidentified kaon.
The presence of (π + π − ) pairs originating from J/ψ , χ c0 and χ c2 charmonia decays is vetoed by requiring the invariant masses M of all possible (π + π − ) pairs to be |M − M 0 | > 30 MeV/c 2 , where M 0 stands for the corresponding known values of the J/ψ , χ c0 and χ c2 meson masses [21].Similarly, the decays D 0 → K − π + and D 0 → π + π − are vetoed by requiring the corresponding invariant masses to differ by 25 MeV/c 2 or more from the known D 0 meson mass [21].To reduce contamination from other charm backgrounds and from the B 0 → a + 1 (→ ρ 0 π + )π − decay, the invariant mass of any three-body combination in the event is required to be larger than 2100 MeV/c 2 .Simulated B 0 → ρ 0 ρ 0 and B 0 → φK * (892) 0 decays are also used for determining the relative reconstruction efficiencies.The pp collisions are generated using Pythia [22] with a specific LHCb configuration [23].Decays of hadronic particles are described by EvtGen [24].The interaction of the generated particles with the detector and its response are implemented using the Geant4 toolkit [25] as described in Ref. [26].

Four-body mass fit
The four-body mass spectrum M (π + π − )(π + π − ) is fit with an unbinned extended likelihood.The fit is performed simultaneously for the two data taking periods together with the normalisation channel M (K + K − )(K + π − ) and PID misidentification control channel M (K + π − )(π + π − ) mass spectra.The four-body invariant mass models account for B 0 and possible B 0 s signals, combinatorial backgrounds, signal cross-feeds and background contributions arising from partially reconstructed b-hadron decays in which one or more particles are not reconstructed.
The B 0 and B 0 s meson shapes are modelled with a modified Crystal Ball distribution [27].A second power-law tail is added on the high-mass side of the signal shape to account for imperfections of the tracking system.The model parameters are determined from a simultaneous fit of simulated signal events that fulfill the trigger, reconstruction and selection chain, for each data taking period.The values of the tail parameters are identical for the B 0 and B 0 s mesons.Their mass difference is constrained to the value from Ref. [21].The mean and width of the modified Crystal Ball function are free parameters of the fit to the data.
The data are represented by the black dots.The fit is represented by the solid blue line, the B 0 signal by the solid red line and the B 0 s by the solid green line.The combinatorial background is represented by the pink dotted line, the partially reconstructed background by the cyan dotted line and the cross-feed by the dark blue dashed line.
The combinatorial background in each four-body spectrum is described by an exponential function where the slope is allowed to vary in the fit.
The misidentification of one or more final-state hadrons may result in a fully reconstructed background contribution to the corresponding signal spectrum, denoted signal cross-feed.The magnitude of the branching fractions of the signal and control modes as well as the two-body mass selection criteria make these signal cross-feeds negligible, with one exception: the misidentification of the kaon of the decay B 0 → (K + π − )(π + π − ) as a pion yields a significant contribution in the M (π + π − )(π + π − ) mass spectrum.The mass shape of B 0 → (K + π − )(π + π − ) decays reconstructed as B 0 → (π + π − )(π + π − ) is modelled by a Crystal Ball function, whose parameters are determined from simulated events.The yield of this signal cross-feed is allowed to vary in the fit.The measurement of the actual number of reconstructed B 0 → (K + π − )(π + π − ) events multiplied by the data-driven estimate of the misidentification efficiency is consistent with the measured yield.
The partially reconstructed background is modelled by an ARGUS function [28] convolved with a Gaussian function accounting for resolution effects.Various mass shape parameterisations are examined.The best fit is obtained when the endpoint of the ARGUS function is fixed to the value expected when one pion is not attributed to the decay.The other shape parameters of the ARGUS function are free parameters of the fit, common to the two data taking periods.The floating width parameter of the signal mass shape is constrained to be equal to the width of the Gaussian function used in the convolution.
Figure 1 displays the M (π + π − )(π + π − ) and M (K + K − )(K + π − ) spectra with the fit results overlaid.The signal event yields are shown in Table 1.Aside from the prominent signal of the ) is observed with a statistical significance of more than 10 standard deviations.The statistical significance is evaluated by taking the ratio of the likelihood of the nominal fit and of the fit with the signal yield fixed to zero.A systematic uncertainty due to the fit model is associated to the measured yields.The dominant uncertainties arise from the knowledge of the signal and signal crossfeed shape parameters determined from simulated events.Several pseudoexperiments are generated while varying the shape parameters within their uncertainties, and the systematic uncertainties on the yields are estimated from the differences in results with respect to the nominal fit.

Amplitude analysis
An amplitude analysis is used to determine the vector-vector (VV) contribution B 0 → ρ 0 ρ 0 by using two-body mass spectra and angular variables.The four-body mass spectrum is first analysed with the sPlot technique [29] to subtract statistically the background under the B 0 → (π + π − )(π + π − ) signal.
For the two-body mass spectra, contributions from resonant and non-resonant scalar (S), resonant vector (V ) and tensor (T ) components are considered in the amplitude fit model through complex mass propagators, M (m i ), where the label i = 1, 2 are the first and second pion pairs, which are assigned randomly in every decay since they are indistinguishable.The P-wave lineshape model comprises the ρ 0 meson, described using the Gounaris-Sakurai parameterisation M ρ (m i ) [30], and the ω meson, parameterised with a relativistic spin-1 Breit-Wigner M ω (m i ).The D-wave lineshape M f 2 (m i ) accounts for the f 2 (1270), modelled with a relativistic spin-2 Breit-Wigner.The S-wave model includes the f 0 (980) propagator M f (980) (m i ), described using a Flatté parameterisation [31,32], and a low-mass component.The latter includes the broad low-mass resonance f 0 (500) and a non-resonant contributions, which are jointly modelled in the framework of the K-matrix formalism [33] and referred as M (ππ) 0 (m i ).Following the K-matrix formalism, the amplitude for the low-mass π + π − S-wave can be written as where κ is measured to be −0.07 ± 0.24 from a fit to the inclusive π + π − mass distribution and m 0 and Γ are the nominal mass and mass-dependent width of the f 0 (500), as determined in Ref. [34].The functions ρ(m) and q(m), defined in Ref. [33], are the phase space factor and the relative momentum of a pion in the ρ 0 centre-of-mass system.By convention, the phase of the M (ππ) 0 (m i ) mass propagator is set to zero at the ρ 0 nominal mass.The signal sample is described by considering the dominant amplitudes of the signal decay.The B → V V component contains the B → ρ 0 ρ 0 and B 0 → ρ 0 ω amplitudes.The B → V S component accounts for B 0 → ρ 0 (π + π − ) 0 and B 0 → ρ 0 f 0 (980) amplitudes and the B → V T contribution is limited to the purely longitudinal amplitude of the B 0 → ρ 0 f 2 (1270) transition.Because of the broad natural width of the a ± 1 particle, a small contamination from the decays B 0 → a ± 1 π ∓ remains in the sample.This contribution with a ± 1 → ρ 0 π ± in S-wave is considered along with its interference with the other amplitudes.This is done by introducing the CP -even eigenstate from the linear combination of individual amplitudes of the decays B 0 → a + 1 π − and B 0 → a − 1 π + , as defined in Ref. [35].The contribution of the decays ) ,⊥ are assumed to be negligible, where the and ⊥ subindices indicate the parallel and perpendicular amplitudes of the decay.The choice of the baseline model was made prior to the measurement of the physical parameters of interest after comparing a set of alternative parameterisations according to a dissimilarity statistical test [36].
The differential decay rate for B 0 → (π + π − )(π + π − ) decays at the B 0 production time t = 0 is given by where the variables θ 1 , θ 2 and ϕ are the helicity angles, described in Fig. 2, and Φ 4 is the four-body phase space factor.The notations of the complex amplitudes, A i , and the expressions of their related angular distributions, f i , are displayed in Table 2.The mass propagators included in the f i functions are normalised to unity in the fit range.
For the CP conjugated mode, B 0 → (π + π − )(π + π − ), the decay rate is obtained under the transformation A i → η i A i , where η i is the CP eigenvalue of the CP eigenstate i, shown in Table 2.
The untagged time-integrated decay rate of B 0 and B 0 to four pions, assuming no CP violation, can be written as Table 2: Amplitudes, A i , CP eigenvalues, η i , and mass-angle distributions, f i , of the B 0 → (π + π − )(π + π − ) model.The indices ijkl indicate the eight possible combinations of pairs of opposite-charge pions.The angles α kl , β ij and Φ kl are defined in Ref. [37].where δ ij = 1 when i = j and δ ij = 0 otherwise.The efficiency of the selection of the final state B 0 → (π + π − )(π + π − ) varies as a function of the helicity angles and the two-body invariant masses.To take into account variations in the efficiencies, four event categories k are defined according to their hardware trigger decisions (TIS or TOS) and data taking period (2011 and 2012).
The acceptance is accounted for through the complex integrals where f i are the functions given in Table 2 and the overall efficiency.The integrals are computed with simulated events of each of the four considered categories, selected with Parameter Definition Fit result the same criteria as those applied to data, following the method described in Ref. [38].
The coefficients ω k ij are used to determine the efficiency and to build a probability density function for each category, which is defined as The four event categories are used in the simultaneous unbinned maximum likelihood fit which depends on the 19 free parameters indicated in Table 3.
Systematic effects are estimated by fitting with the angular model an ensemble of 1000 pseudoexperiments generated with the same number of events as observed in data.The biases are for the parameters of interest consistent with zero.A systematic uncertainty is assigned by taking 50% of the fit bias or the uncertainty on the rms when the latter is bigger in order to account for possible statistical fluctuations.
Several model related uncertainties are envisaged.The B 0 → a ± 1 π ∓ angular model requires knowledge of the lineshape of the a ± 1 meson.The a ± 1 natural width is chosen to be 400 MeV/c 2 .The difference to the fit results obtained by varying the width from 250 to 600 MeV/c 2 is taken as the corresponding systematic uncertainty.In addition, a systematic uncertainty is obtained by introducing the CP -odd component in the fit model of the decay amplitude B 0 → a ± 1 π ∓ by fixing the relative amplitudes of B 0 → a + 1 π − and B 0 → a − 1 π + components to the values measured in Ref. [39].Another source of uncertainty originates in the modelling of the low mass (π + π − ) S-wave lineshape.The f 0 (500) mass and natural width uncertainties from Ref. [34] and the uncertainty on the parameter that quantifies the non-resonant contribution are propagated to the angular analysis parameters by generating and fitting 1000 pseudoexperiments in which these input values are varied according to a Gaussian distribution having their uncertainties as widths.The root mean square of the distribution of the results is assigned as a systematic uncertainty.The same strategy is followed to estimate the systematic uncertainties originating from the ρ 0 , f 0 (500) and ω lineshape parameters.
The uncertainty related to the background subtraction method is estimated by varying within their uncertainties the fixed parameters of the mass fit model and studying the resulting angular distributions and two-body mass spectra.The difference to the fit results is taken as a systematic uncertainty.An alternative subtraction of the background estimated from the high-mass sideband is performed, yielding compatible results.
The knowledge of the acceptance model described in Eq.( 6) comes from a finite sample of simulated events.An ensemble of pseudoexperiments is generated by varying the acceptance weights according to their covariance matrix.The root mean square of the distribution of the results is assigned as a systematic uncertainty.
The resolution on the helicity angles is evaluated with pseudoexperiments resulting in a negligible systematic uncertainty.The systematic uncertainty related to the (π + π − ) mass resolution is estimated with pseudoexperiments by introducing a smearing of the (π + π − ) mass.Differences in the parameters between the fit with and without smearing are taken as a systematic uncertainty.
Table 4 details the contributions to the systematic uncertainty in the measurement of the fraction of B 0 → ρ 0 ρ 0 signal decays in the B 0 → (π + π − )(π + π − ) and its longitudinal polarisation fraction.
The B 0 → ρ 0 ρ 0 signal significance is measured to be 7.1 standard deviations.The significance is obtained by dividing the value of the purity by the quadrature of the statistical and systematic uncertainties.No evidence for the B 0 → ρ 0 f 0 (980) decay mode  , cos θ 1,2 and ϕ distributions.The black dots correspond to the four-body background-subtracted data and the black line is the projection of the fit model.The specific decays B 0 → ρ 0 ρ 0 (brown), B 0 → ωρ 0 (dashed brown), B 0 → V S (dashed blue), B 0 → SS (long dashed green), B 0 → V T (orange) and B 0 → a ± 1 π ∓ (light blue) are also displayed.The B 0 → ρ 0 ρ 0 contribution is split into longitudinal (dashed red) and transverse (dotted red) components.Interference contributions are only plotted for the total (black) model.The efficiency for longitudinally polarized B 0 → ρ 0 ρ 0 events is ∼5 times smaller than for the transverse component.

Branching fraction determination
The branching fraction of the decay mode B 0 → ρ 0 ρ 0 relative to the decay B 0 → φK * (892) 0 can be expressed as where the factor λ f L corrects for differences in detection efficiencies between experimental and simulated data due to the polarisation hypothesis of the B 0 → ρ 0 ρ 0 sample, P (B 0 → ρ 0 ρ 0 ) and P (B 0 → φK * (892) 0 ) are the fractions of B 0 → ρ 0 ρ 0 and B 0 → φK * (892) 0 signals in the samples of ) decays as determined from a fit to the four-body mass distributions, weighted for each data-taking period by the efficiencies of the signal and normalisation channels obtained from their respective simulated data.Finally, B(φ → K + K − ), B(K * (892) 0 → K + π − ) and B(ρ 0 → π + π − ) denote known branching fractions [21].
The product λ f L • P (B 0 → ρ 0 ρ 0 ) is determined from the amplitude analysis to be 1.13 ± 0.19 (stat) ± 0.10 (syst).This quantity is mainly related to the modelling of the S-wave component, and dominates the systematic uncertainty of the parameters of interest.
The fraction of B 0 → φK * (892) 0 present in the B 0 → (K + K − )(K + π − ) sample is taken from Ref. [40].A 1% systematic uncertainty is added, accounting for differences in the selection acceptance for P-and S-wave contributions.
The amounts of B 0 → (π + π − )(π + π − ) and B 0 → (K + K − )(K + π − ) candidates are determined from the four-body mass spectra analysis and their associated statistical and systematical uncertainties are propagated quadratically to the branching fraction uncertainty estimate.
The limited size of the simulated events samples that meet all selection criteria result in a systematic uncertainty of 1.7% (2.6%) on the measurement of the relative branching fraction for the 2011 (2012) data-taking period.The impact of the discrepancies between experimental and simulated data related to the B 0 meson kinematical properties is 0.6% (1.2%).The efficiencies of the particle-identification requirements are determined from control samples of data with a systematic uncertainty of 0.5%, mostly originating from the limited size of the calibration samples.An additional 1% systematic uncertainty on the tracking efficiency is added accounting for different interaction lenghts between π and K.

Conclusions
The full data set collected by the LHCb experiment in 2011 and 2012, corresponding to an integrated luminosity of 3.0 fb −1 , is analysed to search for the B 0 → ρ 0 ρ 0 decay.A yield of 634 ± 28 ± 8 B 0 → (π + π − )(π + π − ) signal decays with π + π − pairs in the 300-1100 MeV/c 2 mass range is obtained.An amplitude analysis is conducted to determine the contribution from B 0 → ρ 0 ρ 0 decays.This decay mode is observed for the first time with a significance of 7.1 standard deviations.In the same π + π − pairs mass range, B 0 s → (π + π − )(π + π − ) decays are also observed with a statistical significance of more than 10 standard deviations.
The measured longitudinal polarisation fraction is consistent with the measured value from BaBar [8] while it differs by 2.3 standard deviations from the value obtained by Belle [9].The branching fraction measurement is in agreement with the values measured by both BaBar [8] and Belle [9] collaborations.

Table 1 :
Yields from the simultaneous fit for the 2011 and 2012 data sets.The first and second uncertainties are the statistical and systematic contributions, respectively.

Table 3 :
Results of the unbinned maximum likelihood fit to the angular and two-body invariant mass distributions.The first uncertainty is statistical, the second systematic.
Università di Bologna, Bologna, Italy e Università di Cagliari, Cagliari, Italy f Università di Ferrara, Ferrara, Italy g Università di Firenze, Firenze, Italy h Università di Urbino, Urbino, Italy i Università di Modena e Reggio Emilia, Modena, Italy j Università di Genova, Genova, Italy k Università di Milano Bicocca, Milano, Italy l Università di Roma Tor Vergata, Roma, Italy m Università di Roma La Sapienza, Roma, Italy n Università della Basilicata, Potenza, Italy o AGH -University of Science and Technology, Faculty of Computer Science, Electronics and Telecommunications, Kraków, Poland p LIFAELS, La Salle, Universitat Ramon Llull, Barcelona, Spain q Hanoi University of Science, Hanoi, Viet Nam r Università di Padova, Padova, Italy s Università di Pisa, Pisa, Italy t Scuola Normale Superiore, Pisa, Italy u Università degli Studi di Milano, Milano, Italy v Politecnico di Milano, Milano, Italy d