University of Birmingham Search for CP violation in D 0-+ 0 decays with the energy test

Citation for published version (Harvard): LHCb Collaboration, Aaij, R, Adeva, B, Adinolfi, M, Affolder, A, Ajaltouni, Z, Akar, S, Albrecht, J, Alessio, F, Alexander, M, Ali, S, Alkhazov, G, Alvarez Cartelle, P, Alves, AA, Amato, S, Amerio, S, Amhis, Y, An, L, Anderlini, L, Anderson, J, Andreassen, R, Andreotti, M, Andrews, JE, Appleby, RB, Aquines Gutierrez, O, Archilli, F, Artamonov, A, Artuso, M, Aslanides, E, Auriemma, G, Baalouch, M, Bachmann, S, Back, JJ, Badalov, A, Baesso, C, Baldini, W, Barlow, RJ, Barschel, C, Barsuk, S, Barter, W, Batozskaya, V, Bifani, S, Farley, N, Griffith, P, Kenyon, IR, Lazzeroni, C, Mazurov, A, McCarthy, J, Parkinson, CJ, Pescatore, L & Watson, NK 2015, 'Search for CP violation in D 0-+0 decays with the energy test', Physics Letters B, vol. 740, pp. 158-167. https://doi.org/10.1016/j.physletb.2014.11.043


Introduction
The decay D 0 → π − π + π 0 (charge conjugate decays are implied unless stated otherwise) proceeds via a singly Cabibbosuppressed c → dud transition with a possible admixture from a penguin amplitude. The interference of these amplitudes may give rise to a violation of the charge-parity symmetry (CP violation), which may be observed as an asymmetry in the total rates, or in the distribution of events over the Dalitz plot. Contributions from particles that are not described in the Standard Model (SM) and participate in the loops of the penguin amplitude can enhance the O (10 −3 ) CP violation effects expected within the SM [1]. Therefore CP violation in D 0 → π − π + π 0 decays provides sensitivity to such non-SM physics.
In addition to this direct CP violation, time-integrated CP asymmetry in D 0 → π − π + π 0 decays can also receive an indirect contribution arising from either the D 0 -D 0 mixing or interference in decays following mixing. While direct CP asymmetry depends on the decay mode, indirect CP violation is expected to be the same for all CP eigenstates. Recent time-dependent measurements of D 0 → π − π + , K − K + decays constrain the indirect CP asymmetry to the O (10 −3 ) level [2].
Singly Cabibbo-suppressed charm decays have recently received a significant attention in the literature (see e.g. Ref. [1]). A particular interest in the decay D 0 → π − π + π 0 was pointed out in Ref. [4], where the authors derived isospin relations between the different D → ρπ amplitudes and discussed the possibility of identifying contributions from non-SM physics using CP violation measurements.
Previously, the most sensitive search for CP violation in this decay was performed by the BaBar Collaboration [5]. Their result excluded CP-violating effects larger than a few percent. The results presented here are based on a signal sample that is about eight times larger and have higher precision. This is the first CP violation analysis performed at the LHCb experiment with decays involving π 0 mesons. As π − π + π 0 is a self-conjugate final state and accessible to both D 0 and D 0 decays, flavour tagging of the D mesons is performed through the measurement of the soft pion (π s ) charge in the D * + → D 0 π + s decay. The method exploited in this Letter, called the energy test [6,7], verifies the compatibility of the observed data with CP symmetry. It is sensitive to local CP violation in the Dalitz plot and not to global asymmetries. The unbinned technique applied here is used for the first time. A visualisation method is also used that allows identification of regions of the Dalitz plot in which CP violation is observed. As this model-independent method cannot identify which amplitudes contribute to the observed asymmetry, a modeldependent analysis would be required in the case of a signal for a non-zero CP asymmetry.

Detector and reconstruction
The LHCb detector [8] 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 combined tracking system provides a measurement of momentum, p, with a relative uncertainty that varies from 0.4% at low momentum to 0.6% at 100 GeV/c. The minimum distance of a track to a primary vertex, the impact parameter (IP), is measured with a resolution of (15 + 29/p T ) μm, where p T is the component of p transverse to the beam, in GeV/c. Charged hadrons are identified using two ring-imaging Cherenkov detectors. Photon, electron and hadron candidates 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.
The trigger consists of a hardware stage, based on high-p T signatures from the calorimeter and muon systems, followed by a two-level software stage, which applies partial event reconstruction. The software trigger at its first level requires at least one good quality track associated with a particle having high p T and high χ 2 IP , defined as the difference in χ 2 of the primary pp interaction vertex (PV) reconstructed with and without this particle.
A dedicated second-level trigger partially reconstructs D 0 → π − π + π 0 candidates coming from D * + → D 0 π + s decays using only information from the charged particles. These requirements ensure the suppression of combinatorial background without distorting the acceptance in the decay phase space. The π − π + pair is combined with a pion to form a D * + candidate, which is accepted if it has p T greater than 2.5 GeV/c and a difference of invariant masses m(π − π + π + s ) − m(π − π + ) < 285 MeV/c 2 . All the chargedparticle tracks used must have good quality, p T > 0.3 GeV/c and p > 3.0 GeV/c, while the charged pions from the D 0 decays must also have high χ 2 IP . The π − π + combination is required to form a good quality secondary vertex significantly displaced from the PV, while the soft pion must originate from the PV. The inclusive D * + trigger was introduced at the start of 2012. The data collected by the LHCb experiment in 2012 corresponds to an integrated luminosity of 2 fb −1 of pp collisions collected at a centre-of-mass energy of 8 TeV. The magnetic field polarity is reversed regularly during the data taking with approximately half the data collected at each polarity to reduce the overall effect of any charge-dependent detection and reconstruction efficiency.
In this analysis, the two π 0 categories reconstructed in LHCb are exploited [9]. These are pions for which both final state photons are reconstructed separately (resolved pions), as well as pions that have higher momentum (typically p T > 2 GeV/c) and thus a smaller opening angle of the two photons (merged pions). These π 0 mesons are detected in the calorimeter as one merged cluster which is further split into two subclusters based on the expected shape of the photon shower. The merged pions make up about 30% of the reconstructed π 0 mesons. Among the resolved π 0 mesons there are also candidates made of photons which, after interacting with detector material, have converted into an e + e − pair. The two π 0 samples provide coverage of complementary regions of the D 0 → π − π + π 0 Dalitz plot and thus the use of both contributes significantly to the sensitivity of the analysis.

Event selection
The offline selection is split into a pre-selection, which follows the trigger selection, and a selection based on a boosted decision tree (BDT) [10,11]. All the D * + candidates are required to pass both levels of the software trigger. In addition, the pre-selection requires more stringent kinematic criteria than those applied in the trigger; in particular p T > 0.5 GeV/c is required for all the D 0 decay products to reduce the combinatorial background. For resolved π 0 mesons, the di-photon invariant mass has to be within 15 MeV/c 2 of the known π 0 mass, this corresponds to about three times the m(γ γ ) resolution. The invariant mass of merged photons, due to its lower resolution, is required to be within the range of 75-195 MeV/c 2 . The purity of the merged π 0 sample is nevertheless significantly higher with respect to the resolved π 0 , as it benefits from large transverse energy and much lower combinatorial background.
Cross-feed from D 0 → K − π + π 0 decays, with a kaon misidentified as a pion, is reduced with requirements on the π ± particle identification based on the ring-imaging Cherenkov detectors.
The D 0 candidates satisfying the above criteria and having invariant mass m(π − π + π 0 ) within 40 (60) MeV/c 2 of the known D 0 mass are accepted in the resolved (merged) sample; this range corresponds to approximately four times the m(π − π + π 0 ) resolution. The D * + candidates, formed with the D 0 and π + s candidate, have their entire decay chain refitted requiring that the D * + candidate originates from the corresponding PV, and π 0 and D 0 candidates have their nominal masses. This improves the D * + mass resolution and the resolution in the D 0 → π − π + π 0 Dalitz plot, while a requirement put on the fit quality efficiently suppresses the background. This requirement also suppresses the contribution from the D * + mesons originating from long-lived b-hadrons. The remainder of this component is not affected by CP asymmetries in b-hadrons since the flavour tag is obtained from the D * + meson. This preliminary selection is followed by a multivariate analysis based on a BDT. Signal and background samples used to train the BDT are obtained by applying the sPlot technique [12] to a quarter of the real data. The sWeights for signal and background separation are determined from a fit to the distribution of the mass differ- The optimal value of the BDT discriminant is determined by estimating the D * + signal significance for various requirements on the BDT output. It retains approximately 75% (90%) of the resolved (merged) signal events while removing 90% (55%) of the background. Fig. 1 shows the m distributions for the selected data set for events with resolved and merged π 0 candidates. The signal shapes, fitted in Fig. 1 with a sum of three Gaussian functions, significantly differ between both samples reflecting the different π 0 momentum resolutions. The lower momentum resolution of the merged π 0 mesons relative to the resolved π 0 mesons makes the core part of the merged signal distribution wider, while the low-p T π 0 mesons contributing to the resolved signal enlarge its tail component. The background shape is fitted using a secondorder polynomial multiplied by √ 1 − m π + / m. The final signal sample is selected requiring | m − 145.4| < 1.8 MeV/c 2 , which corresponds to roughly four times the effective m resolution. The effective resolution is similar for both resolved and merged π 0 samples when averaging the narrow and broad components of the peak. This gives 416 × 10 3 resolved and 247 × 10 3 merged signal candidates with a purity of 82% and 91%, respectively. The Dalitz plot of the final signal sample is shown in Fig. 2. The smaller number of candidates in the low m 2 (π + π − ) region compared with the high m 2 (π + π − ) region is due to acceptance effects related to the π 0 reconstruction as discussed in Section 5.

Energy test method
Model-independent searches for local CP violation are typically carried out using a binned χ 2 approach to compare the relative density in the Dalitz plot of a decay and its CP-conjugate sample (see for example [5,13]). A model-independent unbinned statistical method called the energy test was introduced in Refs. [6,7]. Ref. [14] suggests applying this method to Dalitz plot analyses and demonstrates the potential to obtain improved sensitivity to CP violation over the standard binned approach. This Letter describes the first application of this technique to experimental data.
In this method a test statistic, T , is used to compare average distances in phase space, based on a metric function, ψ ij , of pairs of events ij belonging to two samples of opposite flavour. It is defined as where the first and second terms correspond to a metric-weighted average distance of events within n events of one flavour and n events of the opposite flavour, respectively. The third term measures the weighted average distance of events in one flavour sample to events of the opposite flavour sample. The normalisation factors in the denominator remove the impact of global asymmetries. If the distributions of events in both flavour samples are identical, T will fluctuate around a value close to zero. The metric function should be falling with increasing distance d ij between events i and j, in order to increase the sensitivity to local asymmetries. A Gaussian metric is chosen, defined as ij /2σ 2 with a tunable parameter σ , which describes the effective radius in phase space within which a local asymmetry is measured. Thus, this parameter should be larger than the resolution of d ij and small enough not to dilute locally varying asymmetries. The distance between two points in phase space is usually measured as the distance in the Dalitz plot. However, this distance depends on the choice of the axes of the Dalitz plot. This dependence is removed by using all three invariant masses to determine the distance, d ij , calculated as the length of the displacement vec- subscripts indicate the final-state particles. Using all three invariant masses does not add information, but it avoids an arbitrary choice that could impact the sensitivity of the method to different CP violation scenarios.
In the case of CP violation, the average distances entering in the third term of Eq. (1) are larger, which, because of the characteristics of the metric function, leads to a reduced magnitude of this term. Therefore larger CP asymmetries lead to larger values of T . This is translated into a p-value under the hypothesis of CP symmetry by comparing the nominal T value observed in data to a distribution of T values obtained from permutation samples, where the flavour of each candidate is randomly reassigned to simulate samples without CP violation. The p-value for the no CP violation hypothesis is obtained as the fraction of permutation T values greater than the nominal T value.
A statistical uncertainty of the p-value is obtained as a binomial standard deviation. If large CP violation is observed, the observed T value is likely to lie outside the range of permutation T values. In this case the permutation T distribution can be fitted with a generalised extreme value (GEV) function, as demonstrated in Refs. [6,7] and verified in large simulation samples for this analysis. The GEV function is defined as with normalisation N, location parameter μ, scale parameter δ, and shape parameter ξ . This function is set to zero for T > μ − δ/ξ for ξ > 0, and for T < μ − δ/ξ for ξ < 0. Fig. 3 shows an example T value distribution with a GEV function fit for a simulated data set including CP violation (see Section 5).
The p-value from the fitted T distribution can be calculated as the fraction of the integral of the function above the nominal T value. The uncertainty on the p-value is obtained by randomly resampling the fit parameters within their uncertainties, taking into account their correlations, and by extracting a p-value for each of these generated T distributions. The spread of the resulting p-value distribution is used to set 68% confidence uncertainties. A 90% confidence upper limit is quoted where no significantly nonzero p-value can be obtained from the fit.
The number of available permutations is constrained by the available computing time. The default p-value extraction uses the counting method as long as at least three permutation T values are found to be larger than the observed T value. Beyond that, the p-value is extracted by integrating the fitted GEV function. The p-values presented here are based on 1000 permutations for the default data results and on 100 permutations for the sensitivity studies.
A visualisation of regions of significant asymmetry is obtained by assigning an asymmetry significance to each event. The contributions of a single event of one flavour, T i , and a single event of the opposite flavour, T i , to the total T value are given by Example T i distributions for the simulated CP violating data sets are shown in Fig. 4(a), (c); events contributing with T i of the largest magnitude point to CP violation regions. However, the CP asymmetry arising from the ρ + amplitude difference (Fig. 4(a) ) distribution that corresponds to the significance level. The same procedure is applied to the T i distribution, leading to a Dalitz plot with an inverted asymmetry pattern.
The asymmetry significances for each simulated event are plotted on a Dalitz plot (see Fig. 4(b), (d)). If an amplitude difference exists between CP-conjugate states of a resonance, the region of significant asymmetry appears as a band around the mass of the resonance on the plot. If a phase difference is present instead, regions of positive and negative asymmetry appear around the resonance on the plot, indicating the phase shift.
The practical limitation of this method is that the number of mathematical operations scales quadratically with the sample size. Furthermore, a significant number of permutations is required to get a sufficient precision on the p-value. In this analysis the method is implemented using parallelisation on graphics processing units (GPUs) [15].

Sensitivity studies
The interpretation of the results requires a study of the sensitivity of the present data sample to different types of CP violation. The sensitivity is examined based on simplified Monte Carlo samples generated according to the model described in Ref. [3] using the generator package Laura++ [16].
The selection efficiency has to be taken into account in these studies as it varies strongly across phase space. This efficiency is measured using a sample of events based on the full LHCb detector simulation. In the simulation, pp collisions are generated using Pythia [17,18] with a specific LHCb configuration [19]. Decays of hadronic particles are described by EvtGen [20], in which final-state radiation is generated using Photos [21]. The interaction of the generated particles with the detector and its response are implemented using the Geant4 toolkit [22,23] as described in Ref. [24]. The efficiency is shown to vary as a function of the π + π − invariant mass, while it is found to be constant with respect to other variables. The resulting efficiency curves for merged  Table 1 Overview of sensitivities to various CP violation scenarios. A and φ denote, respectively, change in amplitude and phase of the resonance R.
and resolved π 0 mesons as well as for the combined sample are shown in Fig. 5. The small m(π + π − ) range corresponds to π 0 candidates with high momentum and is primarily covered by the merged π 0 sample. The overall efficiency in this region is low as the integrated merged π 0 identification efficiency is lower than for the resolved sample and decreases after its turn-on as the π 0 momentum increases. This affects in particular the ρ 0 resonance, as it lies entirely within the low acceptance region. For further studies, the efficiency, based on the fitted curve, is then applied to simplified Monte Carlo data sets by randomly discarding events based on the candidate's position in phase space. Background events are simulated by resampling phase-space distributions extracted from m sideband regions. Inclusion of the background does not significantly reduce the sensitivity to the CP violation scenarios discussed below. This is due to the low level of background and due to it being CP symmetric within the present sensitivity.
Various CP asymmetries are introduced by modifying, for a chosen D 0 flavour, either the amplitude or the phase of one of the three intermediate ρ resonances dominating the π − π + π 0 phase space. The resulting sensitivities are shown in Table 1. The p-values, including their statistical uncertainties, are obtained from fits of GEV functions to the T value distributions and 90% confidence limits are given in addition.
The sensitivity is comparable to that of the BaBar analysis [5] for the ρ 0 resonance and significantly better for the ρ + and ρ − resonances. This is expected due to the variation of the selection efficiency across phase space, which disfavours the ρ 0 region. The optimal σ value may vary with different CP violation scenarios. Hence the final results are also quoted for several values of σ .
The standard binned method [13] is also applied to the simulated data sets. This study shows that the energy test provides results compatible with, and equally or more precise than the binned method.
There are two main sources of asymmetry that may degrade or bias the results. One is an asymmetry that may arise from background events and the other is due to particle detection asymmetries that could vary across phase space.
Background asymmetries are tested by applying the energy test to events in the upper m sideband, m > 150 MeV/c 2 . No significant asymmetry is found. In addition, simplified simulation data sets are produced by generating signal candidates without CP violation and background candidates according to background distributions in data, separately for D 0 and D 0 candidates and thus allowing for a background-induced asymmetry. These samples show a distribution of p-values consistent with the absence of any asymmetry. Further tests using a binned approach [13] confirm this conclusion. These are carried out on the m sideband data sample as well as on background samples obtained using the sPlot technique based on the m fits in Fig. 1. Both approaches show no indication of a background asymmetry. As the background present in the signal region is found to be CP symmetric, it is simply included in the T value calculation discussed in Section 4. Local asymmetries are expected to arise, at a level below the current sensitivity, due to the momentum dependence of π + /π − detection asymmetries in combination with the different kinematic distributions of π + and π − in certain regions of phase space.
These effects are tested using the Cabibbo-favoured decay D 0 → K − π + π 0 as a control mode. This channel is affected by kaon detection asymmetries, which are known to be larger than pion detection asymmetries and thus should serve as a conservative test. The data sample is split into eight subsets, each of which contains approximately the same amount of data as the signal sample. The energy test yields p-values between 3% and 74%, which is consistent with the assumption that detection asymmetries are below the current level of sensitivity. A further test is conducted by splitting the control mode data sample by the polarity of the spectrometer dipole magnet, which yields two large approximately equal-sized samples. The resulting p-values of 8% and 15% show no evidence of sizable biases due to detector asymmetries.

Results and conclusions
The application of the energy test to all selected D 0 → π − π + π 0 candidates using a metric parameter of σ = 0.3 GeV 2 /c 4 yields T = 1.84 × 10 −6 . The permutation T value distribution is shown in Fig. 6(a). By counting the fraction of permutations with a T value above the nominal T value in the data, a p-value of (2.6 ± 0.5) × 10 −2 is extracted. Alternatively, extrapolation from a fit to the GEV function gives a p-value of (2.1 ± 0.3) × 10 −2 . The significance levels of the T i values are shown in Fig. 6(b). A small phase-space region dominated by the ρ + resonance contains candidates with a local positive asymmetry exceeding 1σ significance.
Varying the metric parameter results in the p-values listed in Table 2; all the p-values are at the 10 −2 level. The data sample has been split according to various criteria to test the stability of the results. Analyses of sub-samples with opposite magnet polarity, with different trigger configurations, and with fiducial selection requirements removing areas of high local asymmetry of the tagging soft pion from the D * + decay all show good consistency of the results.  In summary, a search for time-integrated CP violation in the Cabibbo-suppressed decay D 0 → π − π + π 0 is performed using a novel unbinned model-independent technique. The analysis has the best sensitivity from a single experiment to CP violation in this decay. The data are found to be consistent with the hypothesis of CP symmetry with a p-value of (2.6 ± 0.5)%.