Measurement of the branching fraction and $\it CP$ asymmetry of $B^{0} \rightarrow \pi^{0} \pi^{0}$ decays using $198 \times 10^6$ $B\overline{B}$ pairs in Belle II data

We report measurements of the branching fraction and $\it CP$ asymmetry in $B^{0} \to \pi^{0} \pi^{0}$ decays reconstructed at Belle II in an electron-positron collision sample containing $198 \times 10^{6}$ $B\overline{B}$ pairs. We measure a branching fraction $\mathcal{B}(\Bpipi) = (1.38 \pm 0.27 \pm 0.22) \times 10^{-6}$ and a $\it CP$ asymmetry $\Acp(\Bpipi) = 0.14 \pm 0.46 \pm 0.07$, where the first uncertainty is statistical and the second is systematic.

The study of decay-time-dependent CP asymmetries in decays dominated by the b → u transition, specifically hadronic decays of bottom mesons into charmless two-body final states, is currently the most precise way to measure the least known angle of the unitarity triangle, ϕ 2 (or α) ≡ arg −V td V * tb /V ud V * ub .Here, V ij are elements of the Cabibbo-Kobayashi-Maskawa (CKM) quark-mixing matrix [1].Improved measurements of ϕ 2 will test the unitarity of the CKM matrix and constrain possible flavor-structure extensions of the standard model (SM).One approach is to measure the timedependent decay-rate asymmetry between B 0 and B 0 mesons that decay to π + π − final states.This asymmetry would be proportional to sin(2ϕ 2 ) if the decay involved only tree-level b → u processes.However, the asymmetry is affected by an unknown and difficult-to-predict shift with respect to the desired ϕ 2 angle due to the presence of b → d loop ('penguin') contributions.The tree-level and penguin amplitudes have similar magnitudes, so the shift is sizable and complicates the determination of ϕ 2 .The penguin and tree contributions can be disentangled using the B → ππ isospin relations [2, 3] where A ij and Āij are amplitudes for the decays B → π i π j and B → π i π j , respectively.Here, B and π indicate charged or neutral bottom-mesons and pions, respectively, while i and j refer to electric charge.Taking advantage of these relations requires precise measurements of the branching fraction B and CP asymmetries of each B → ππ decay mode.The greatest limitation to exploiting the isospin relations lies in the uncertainty of the B 0 → π 0 π 0 inputs, B and the time-integrated CP asymmetry, where Γ is the decay width.The world-average values B(B 0 → π 0 π 0 ) = (1.59 ± 0.26) × 10 −6 and A CP (B 0 → π 0 π 0 ) = 0.33 ± 0.22 [4] combine measurements reported by the BaBar [5] and Belle [6] collaborations.Additional measurements would improve our knowledge of ϕ 2 .Theoretical predictions for B(B 0 → π 0 π 0 ) based on QCD factorization [7][8][9][10] and perturbative QCD [11,12] are approximately five times smaller than the world average value.Furthermore, the ratio of color-suppressed to color-allowed tree amplitudes, as inferred from other charmless two-body decay modes, does not agree well with expectations [13].
This might indicate large electroweak-penguin contributions, which are difficult to explain in the SM [14,15].Various approaches, which predict a wide range of values for B and A CP , have been proposed as possible solutions to this disagreement [16][17][18][19].More precise measurements of these quantities would help in discriminating among the various solutions proposed to address this discrepancy.In addition, a better understanding of the color-suppressed tree amplitude could help resolve the so-called B → Kπ puzzle [20][21][22] In this paper, we present a measurement of B and A CP for the B 0 → π 0 π 0 decay using a data sample consisting of (198.0 ± 3.0) × 10 6 BB pairs [23] collected from 2019 through 2021 [24].The sample is collected with the Belle II detector, located at the SuperKEKB asymmetricenergy e + e − collider [25].A full description of the Belle II detector is given in Ref. [26].The detector consists of several subdetectors arranged in a cylindrical structure around the beam pipe.The z axis of the lab frame is defined as the symmetry axis of a superconducting solenoid, which generates a 1.5 T uniform field along the beam di-rection.The positive direction is given by the electronbeam direction, and the polar angle, θ, is defined with respect to the +z axis.The detector is divided into three regions, and in increasing order of θ, they are referred to as the forward endcap, barrel, and backward endcap.The inner subdetectors are a silicon pixel detector surrounded by a four-layer double-sided silicon strip detector and a central drift chamber (CDC).These subdetectors are used to reconstruct charged particles and measure their momentum.A time-of-propagation counter [27] and an aerogel ring-imaging Cherenkov detector cover the barrel and forward endcap regions, respectively, and are used for charged particle identification.The electromagnetic calorimeter (ECL) is a segmented array of 8736 thallium-doped cesium iodide [CsI(Tl)] crystals arranged in a projective geometry toward the interaction point and covering about 90% of the solid angle in the center-ofmass (c.m.) frame.The ECL identifies electrons and photons in an energy range of 20 MeV to 4 GeV and occupies the remaining volume inside the superconducting solenoid.Resistive plate chambers and scintillating fibers to identify muons and K 0 L mesons are installed in the flux return of the magnet.
We use GEANT4-based [28] simulated samples to optimize event selection criteria, compare distributions observed in data with expectations, determine fit models, calculate signal efficiencies, and study sources of background.To study the signal, we use 10 7 Υ (4S) → B 0 B 0 decays generated with EvtGen [29], where one B meson decays as B 0 → π 0 π 0 .To study backgrounds, we use a simulated sample approximately five times larger than the data sample.This sample consists of e + e − → Υ (4S) → BB processes and continuum e + e − → qq background, generated with EvtGen and PYTHIA [30] where q denotes a u, d, s, or c quark.To account for a large observed τ + τ − background, we use a sample of e + e − → τ + τ − events generated with KKMC [31] and TAUOLA [32] that is the same size as the continuum sample.To validate our analysis, we use the B 0 → D 0 (→ K + π − π 0 )π 0 decay as a control mode, as it contains two π 0 particles in the final state and has an order of magnitude more yield.We use a simulated sample of 5 × 10 6 control-mode events generated with EvtGen.To calibrate and validate our photon reconstruction, we use the The data are processed with the Belle II analysis software framework [33,34].This is the first measurement of this channel at Belle II.
Measuring B 0 → π 0 π 0 decay properties is challenging, as the decay is both CKM-suppressed and colorsuppressed.As the final state consists of photons with no tracks, it is difficult to reconstruct.In addition, the large number of neutral pions produced in e + e − → q q continuum events can be combined to mimic the B 0 → π 0 π 0 signal.Finally, the reconstruction is susceptible to extraneous photons arising from beam interactions with the beam pipe and residual gas; this is referred to as beaminduced background.Hence, conventional and machine-learning based approaches, validated on data, are employed to achieve optimized selections.The signal yield and A CP are determined by performing a maximum likelihood fit to the data.
The online-event selection requires all events to pass criteria based on total energy and neutral-particle multiplicity.In the offline analysis we identify photon candidates by requiring that the number of crystals in an ECL energy deposition (cluster), which can be fractional as a result of energy splitting with nearby clusters, be greater than 1.5.The cluster timing is required to be within of approximately +2.0σ and −2.5σ about the known π 0 mass.The mass requirement is asymmetric as the reconstructed π 0 mass has a slight negative skew due to energy leakage from the ECL calorimeter.We improve the momentum resolution of the π 0 candidates by performing a kinematic fit that constrains their mass to the known value [4].Signal B 0 candidates are reconstructed by combining two π 0 candidates.To select signal B 0 candidates, two kinematic variables are defined, where E beam is the beam energy and (E B , ⃗ p B ) is the reconstructed four-momentum of the B candidate.All quantities are calculated in the c.m. frame of the Υ (4S) resonance.The M bc and ∆E distributions of signal decays peak at the B mass and zero, respectively.Candidate B mesons are required to have 5.26 < M bc < 5.29 GeV/c 2 and −0.3 < ∆E < 0.2 GeV.The ∆E requirement is not centered around zero because of energy leakage from the ECL cluster.
In data, photon-energy corrections are applied to correct for ECL miscalibration.
Studies of the mode are used to validate the corrections.The π 0 momentum is predicted using the momenta of the charged pions and energy-momentum conservation.In simulation, this predicted momentum is typically closer to the true momentum than the measured momentum as the momentum of charged pions is measured more precisely in the CDC than photon energies are measured in the ECL.When the corrections are applied, the data-simulation difference between the predicted and measured π 0 momentum decreases, and the difference between the ∆E peak position in data and simulation is approximately 1 MeV.
The sample includes a large continuum background.To reduce this background, we use topological variables that take advantage of the jet-like nature of qq events and the spherical distribution of BB events.We train a BDT classifier to analyze 28 variables comprising modified Fox-Wolfram moments [36], sphericity-related quantities [37], thrust-related quantities [38], and sets of concentric cones with various opening angles centered around the thrust axis.Variables that show correlations with ∆E and M bc greater than 5% are excluded.We train the continuum classifier to identify statistically significant signal and background features using simulated signal samples and sideband data; the latter consists of events that satisfy all selection criteria but are in a signal-depleted region 5.22 < M bc < 5.27 GeV/c 2 and 0.1 < ∆E < 0.5 GeV.We use these simulated samples to determine a minimum threshold C min of the continuum classifier output C that minimizes the expected statistical uncertainty of the A CP measurement.This selection rejects 93% of the background while retaining 76.5% of the signal.The continuum classifier output is transformed into a Gaussianlike shape according to where C max is the maximum value of the continuum clas-sifier output.Candidate B mesons are required to have |T c | < 3. The T c distributions of signal candidates and continuum are expected to peak at one and zero, respectively.
After suppression of continuum background, 1.3% of events have more than one B 0 candidate.For such events, the average multiplicity is 2.03 candidates per event.We choose the candidate with the minimum sum of the absolute deviations of the reconstructed π 0 masses from the known value [4].This requirement is 56% efficient in selecting the correct B 0 candidate.Following all selections, 35.5% of signal events remain, of which 99.0% are correctly reconstructed.The high fraction of correctly reconstructed events is due to the low percentage of B B events in which there are three high momentum π 0 candidates.
The resulting event sample consists of four main components: signal, continuum, background from non-signal B decays (BB background), and τ + τ − events.From studies of simulated samples, we find that 90% of BB background is from B + B − events in which the B + meson decays into a ρ + π 0 final state and the charged pion from the subsequent ρ + → π + π 0 decay is not reconstructed.The remaining 10% is dominated by B 0 B 0 events in which the B 0 meson decays into a K 0 S (→ π 0 π 0 )π 0 final state, and one π 0 from the K 0 S decay is not reconstructed.The BB background peaks at similar values of M bc and T c as the signal, but its ∆E distribution is shifted to negative values due to the energy lost from the signal-decay daughter that is not reconstructed.In addition, 2.9% of signal candidates arise from τ + τ − events; these candidates are treated as part of the continuum background because their respective M bc , ∆E, and T c distributions are nearly identical.
To measure A CP , the flavor of the signal B 0 is determined by reconstructing the accompanying (tag-side) B meson in each event using the category-based algorithm described in Ref. [39].The tagging information is encoded in two parameters: the b or b flavor of the tag-side B (q) and the purity (r).The value q = +1 tags a B 0 , whereas q = −1 tags a B 0 .The value of r is the algorithm's confidence for an assigned q value.It is defined as r = 1 − 2w, where w is the fraction of wrongly tagged events and ranges from zero, for no flavour distinction between B 0 and B 0 , to one for an unambiguous flavor assignment.For example, for a sample of events having r = 0, there would be an equal number of correctly and incorrectly tagged events; for a sample having r = 1, there would be no incorrectly tagged events.We divide signal candidates into seven intervals of r, with intervals having approximately equal numbers of events.The signal yield and A CP values are determined by performing a three-dimensional (M bc , ∆E, T c ) unbinned extended maximum likelihood fit simultaneously to events in the , for all seven r bins combined.The result of the fit to the data is shown as a solid blue curve.The fit components are shown as a red dashed curve (signal), blue dotted curve (continuum background), green dash-dotted curve (BB background), and magenta soliddotted curve (crossfeed).The plots are signal-enhanced, which correspond to candidates with 5.275 < M bc < 5.285 GeV/c 2 , −0.10 < ∆E < 0.05 GeV, and 0 < T c < 3. When the respective variable is displayed, the selections on that variable are not applied.The difference between observed and fit value divided by the uncertainty from the fit (pulls) are shown below each distribution.
seven intervals of r.The likelihood function is given by where i is the number of candidates, j is the sample component in terms of signal (s), continuum (c), and BB, and k indicates the r interval.Here, N j denotes the yield for component j, N k denotes the number of candidates in the kth bin, f j k is the fraction of candidates in the kth bin for the jth component, and P j k (M i bc , ∆E i , T i c , q i ) is the probability density function (PDF) to have the ith event of the jth component in the kth bin.The values of f c k implicitly include a factor of one-half due to the division of the data into positive and negative q values for each r intervals.Sideband data are used to determine f c k , while f BB k is obtained from large simulated samples.The PDF for the signal component is where w k is the fraction of signal events incorrectly tagged (wrong-tag), ∆w k is the difference in the wrongtag fraction between positive and negative tags, and ∆ϵ k = ∆ϵ k /2ϵ k is the asymmetry of the tagging efficiency.Here, ϵ k is the tagging efficiency and ∆ϵ k is the difference in the tagging efficiency between positive and negative tags.The fraction of signal events in each r interval (f s k ), along with w k , ∆w k , and ∆ϵ k , are fixed to values obtained from a fit to B 0 → D ( * )− h + decays, where h + stands for a π + or K + , following Ref.[39].The CP asymmetry in data is diluted by a factor (1−2w) due to incorrect tagging, and by a factor of (1 − 2χ d ) due to B 0 B 0 mixing, where χ d = 0.1875 ± 0.0017 is the timeintegrated B 0 B 0 -mixing probability [4].
The T c PDFs of the signal, continuum, and BB components are each modeled using the sum of a Gaussian and a bifurcated Gaussian function with independent mean and width parameters.The T c PDFs for the signal and BB are modeled independently for each r bin using simulated data; this accounts for an observed dependence of T c on r.The T c PDF for continuum events is the same for all r bins and is taken from the data sideband.
For the (M bc , ∆E) modeling, a small correlation between M bc and ∆E for the signal is taken into account by using a two-dimensional kernel-density shape.
To simplify the M bc and ∆E modeling of the B + B − and B 0 B 0 backgrounds, we assume that they follow the same distributions as the dominant B + → ρ + π 0 and B 0 → K 0 S (→ π 0 π 0 )π 0 decays, respectively.The PDFs for the BB backgrounds are the sum of two ARGUS functions in M bc and a kernel-density shape in ∆E.All signal and BB PDF parameters are fixed to those obtained from fits to large samples of simulated events.
The upper endpoint of the M bc distribution depends on the beam energy, which varied throughout the course of data taking.To account for this, the continuum is modeled with eight ARGUS functions that have endpoints evenly spaced from 5.287 to 5.290 GeV/c 2 .The contribution of each ARGUS function is fixed to the fraction of events reconstructed at each of the corresponding c.m. energies.Using eight ARGUS functions models well the variation of the M bc endpoint and provides a good fit to the data.The ∆E distribution of the continuum is modeled with a straight line.We determine the parameters of the continuum PDF for all r bins by fitting to the data sideband region.A small dependence of ∆E on q • r found in simulated samples is neglected, as there are insufficient events in the higher r bins of the data sideband for a reliable fit.The q • r distribution for the continuum events shows an asymmetry that could bias the A CP results.This asymmetry, defined similarly to Eq.( 2), is determined to be 0.033 ± 0.002.To account for this, we include a q • r asymmetry term in the continuum PDF that is equal to the q • r asymmetry term extracted from the sideband data.From simulated experiments, small biases of 1% in the branching fraction and 0.02 in A CP are found; we treat these biases as systematic uncertainties.
The reconstruction and fitting procedure is further validated using B 0 → D 0 (→ K + π − π 0 )π 0 decays.This control sample includes a small crossfeed component in which a particle from the tag side is mistakenly included in the signal reconstruction.All photon and π 0 selections are the same, with the exception of the 1.5 GeV/c threshold on π 0 momentum, which is removed since the π 0 from the D 0 has significantly lower momentum than the π 0 from a signal decay.We determine the branching fraction to be B(B 0 → D 0 (→ K + π − π 0 )π 0 ) = (3.66 ± 0.21) × 10 −5 , and the direct CP asymmetry to be A CP (B 0 → D 0 π 0 ) = 0.01 ± 0.16.The uncertainties for the control mode measurements are statistical only.These values agree with previous measurements [40].Figure 1 shows signal-enhanced projections of the fits to data.The signal-enhanced region is defined as 5.275 < M bc < 5.285 GeV/c 2 , −0.10 < ∆E < 0.05 GeV, and 0 < T c < 3; for each plot, the selection on the plotted variable is not applied.On average, these signalenhanced regions contain 47% of signal decays but only 11% of background.The control mode is also used to calibrate the ∆E width of the signal mode, which is taken from simulation.
We apply the fit described above to the 3177 selected B 0 → π 0 π 0 candidate events.The signal yield, A CP , and continuum yield are free to vary, while the BB yield is fixed to the expectation from simulations.We obtain a signal yield of 93 ± 18 events.Figure 2 shows the signalenhanced projections of the fits to data, separately for positive and negative q tags.The signal-enhanced region for the B 0 → π 0 π 0 signal decay is the same as that for the B 0 → D 0 (→ K + π − π 0 )π 0 control mode and rejects approximately 96% of the continuum background.To determine the signal significance, we convolve the statistical and additive systematic uncertainties and calculate the test statistic 2(log L m − log L 0 ) = 32.0 with two degrees of freedom, where log L m is the log-likelihood of the measured signal yield and log L 0 is determined by fixing the signal yield to zero.The second degree of freedom is lost due to A CP = 0 when there is no signal.A total significance of 5.2 standard deviations is obtained.
The branching fraction is calculated using where N s is the signal yield, ε is the signal reconstruction and selection efficiency, N BB is the number of BB pairs produced, B(π 0 → γγ) [4] is the π 0 → γγ branching fraction, and f +− /f 00 is the ratio of the branching fractions for the decay of Υ (4S) to B + B − and B 0 B 0 .The ratio f +− /f 00 is determined to be 1.065 ± 0.012 ± 0.019 ± 0.047 [41], where the first and second uncertainties are statistical and systematic, respectively, and the third uncertainty is due to the assumption of isospin symmetry in B → J/ψ (→ ℓℓ)K, where ℓ = e or µ.Inserting the values N s = 93 ± 18, ε = (35.5 ± 4.7)%, N BB = (198.0± 3.0) × 10 6 , and B(π 0 → γγ) = (98.823± 0.034)% [4], we obtain where the first and second uncertainties are statistical and systematic (discussed below), respectively.The uncertainty in ε is due to the systematic uncertainty associated with π 0 reconstruction and continuum classifier efficiency.The main sources of systematic uncertainties are listed in Table I and are evaluated as follows.A 3.4% systematic uncertainty associated with the π 0 reconstruction efficiency is determined from data using the decays where the π 0 selection is identical to that of the signal.The π 0 reconstruction efficiency as a function of momentum is also measured using τ − → 3ππ 0 ν and τ − → 3πν decays.A difference of 4.7% in efficiency is observed between the measurement based on D decays and the measurement based on τ leptons.This difference increases the systematic uncertainty for a total of 5.8% per pion.The total systematic uncertainty associated with the π 0 reconstruction efficiency is then 11.6%, as there are two pions and their errors are fully correlated.
The systematic uncertainty associated with the continuum parametrization accounts for the uncertainty in each 0 → π 0 π 0 candidates, for all seven r bins combined, with positive (top) and negative (bottom) q tags.The result of the fit to the data is shown as a solid blue curve.The fit components are shown as a red dashed curve (signal), blue dotted curve (continuum background), and green dash-dotted curve (BB background).The plots are signal-enhanced, which correspond to candidates with 5.275 < M bc < 5.285 GeV/c 2 , −0.10 < ∆E < 0.05 GeV, and 0 < T c < 3. When the respective variable is displayed, the selections on that variable are not applied.The difference between observed and fit value divided by the uncertainty from the fit (pulls) are shown below each distribution.
of the data-driven continuum PDF parameters.The contribution of each parameter is determined by refitting on simulated data with the parameter used in the continuum PDF fluctuated by its one-standard-deviation uncertainties.All other continuum PDF parameters are correspondingly shifted according to their correlation with the fluctuated parameter.The systematic uncertainty is the sum in quadrature of the change in signal yield for each parameter.
The systematic uncertainty associated with the efficiency of the continuum classifier is determined using B 0 → D 0 (→ K + π − π 0 )π 0 decays.The efficiencies of the classifier selection in data and simulation are consistent within the statistical uncertainties.The overall statistical uncertainty is assigned as systematic uncertainty.The f +− /f 00 systematic uncertainty combines the original systematic uncertainty and the uncertainty due to the assumption of isospin symmetry.The systematic uncertainty associated with the fixed BB background yield is determined by refitting on simulated data with the generated BB yield fluctuated by its one-standard-deviation uncertainties.The systematic uncertainty associated with the fixed signal fractions for r bins is determined by refitting simulated data with the signal fractions fluctuated by their one-standard-deviation uncertainties.The systematic uncertainty is the sum in quadrature of the change in signal yield for each bin.A similar procedure to determine the systematic uncertainty associated with fixing the continuum fractions in the r bins is also performed.The systematic uncertainty associated with the photon-energy corrections is determined by refitting on data with the values of the corrections fluctuated by their uncertainties.The largest change in yield is taken as the systematic uncertainty.The systematic uncertainty associated with the assumption of independence of ∆E from r is determined by refitting on simulated data with the ∆E slope for each r bin separately estimated using large simulated samples.The procedure for estimating the uncertainty in the number of BB meson pairs is described in Ref. [24].The systematic uncertainty associated with the choice of (M bc , ∆E) signal models is determined by refitting on simulated data with two uncorrelated Crystal Ball functions [42].A small bias in the calculated branching fraction due to the limitations of the PDFs used to model the data is included as a systematic uncertainty.The systematic uncertainty associated with a possible bias due to best candidate selection is determined by refitting on data with the best candidate randomly selected.The total systematic uncertainty is taken to be the sum in quadrature of the individual contributions (16.2%).

FIG. 1 .
FIG.1.Distributions of M bc (left), ∆E (middle), and T c (right) for the B 0 → D 0 (→ K + π − π 0 )π 0 candidates, for all seven r bins combined.The result of the fit to the data is shown as a solid blue curve.The fit components are shown as a red dashed curve (signal), blue dotted curve (continuum background), green dash-dotted curve (BB background), and magenta soliddotted curve (crossfeed).The plots are signal-enhanced, which correspond to candidates with 5.275 < M bc < 5.285 GeV/c 2 , −0.10 < ∆E < 0.05 GeV, and 0 < T c < 3. When the respective variable is displayed, the selections on that variable are not applied.The difference between observed and fit value divided by the uncertainty from the fit (pulls) are shown below each distribution.

FIG. 2 .
FIG.2.Distributions of M bc (left), ∆E (middle), and T c (right) for the B 0 → π 0 π 0 candidates, for all seven r bins combined, with positive (top) and negative (bottom) q tags.The result of the fit to the data is shown as a solid blue curve.The fit components are shown as a red dashed curve (signal), blue dotted curve (continuum background), and green dash-dotted curve (BB background).The plots are signal-enhanced, which correspond to candidates with 5.275 < M bc < 5.285 GeV/c 2 , −0.10 < ∆E < 0.05 GeV, and 0 < T c < 3. When the respective variable is displayed, the selections on that variable are not applied.The difference between observed and fit value divided by the uncertainty from the fit (pulls) are shown below each distribution.