Search for $CP$ violation in $D_{(s)}^{+}\rightarrow h^{+}\pi^{0}$ and $D_{(s)}^{+}\rightarrow h^{+}\eta$ decays

Searches for $CP$ violation in the two-body decays $D_{(s)}^{+}\rightarrow h^{+}\pi^{0}$ and $D_{(s)}^{+} \rightarrow h^{+}\eta$ (where $h^{+}$ denotes a $\pi^{+}$ or $K^{+}$ meson) are performed using $pp$ collision data collected by the LHCb experiment corresponding to either 9 fb$^{-1}$ or 6 fb$^{-1}$ of integrated luminosity. The $\pi^{0}$ and $\eta$ mesons are reconstructed using the $e^{+}e^{-}\gamma$ final state, which can proceed as three-body decays $\pi^{0}\rightarrow e^{+}e^{-}\gamma$ and $\eta \rightarrow e^{+}e^{-}\gamma$, or via the two-body decays $\pi^{0}\rightarrow \gamma\gamma$ and $\eta\rightarrow\gamma\gamma$ followed by a photon conversion. The measurements are made relative to the control modes $D_{(s)}^{+}\rightarrow K_{S}^{0} h^{+}$ to cancel the production and detection asymmetries. The $CP$ asymmetries are measured to be \begin{align} \mathcal{A}_{CP}(D^{+}\rightarrow \pi^{+}\pi^{0})&= (-1.3 \pm 0.9 \pm 0.6 )\%, \end{align} \begin{align} \mathcal{A}_{CP}(D^{+}\rightarrow K^{+}\pi^{0})&= (-3.2 \pm 4.7 \pm 2.1 )\%, \end{align} \begin{align} \mathcal{A}_{CP}(D^{+}\rightarrow \pi^{+}\eta)&= (-0.2 \pm 0.8 \pm 0.4 )\%, \end{align} \begin{align} \mathcal{A}_{CP}(D^{+}\rightarrow K^{+}\eta)&= (-6 \pm 10 \pm 4 )\%, \end{align} \begin{align} \mathcal{A}_{CP}(D_{s}^{+}\rightarrow K^{+}\pi^{0})&= (-0.8 \pm 3.9 \pm 1.2 )\%, \end{align} \begin{align} \mathcal{A}_{CP}(D_{s}^{+}\rightarrow \pi^{+}\eta)&= (0.8 \pm 0.7 \pm 0.5 )\%, \end{align} \begin{align} \mathcal{A}_{CP}(D_{s}^{+}\rightarrow K^{+}\eta)&= (0.9 \pm 3.7 \pm 1.1 )\%, \end{align} where the first uncertainties are statistical and the second systematic. These results are consistent with no $CP$ violation and mostly constitute the most precise measurements of $\mathcal{A}_{CP}$ in these decay modes to date.


Introduction
The observation of Charge-Parity (CP ) violation in two-body decays of neutral D mesons [1] motivates searches for similar effects in charged D meson decays. The twobody D + (s) → h + π 0 and D + (s) → h + η decays, where h + denotes a π + or K + meson, 1 are mediated by Cabibbo favoured (CF), singly Cabibbo suppressed (SCS) or doubly Cabibbo suppressed (DCS) processes. The contributing decay topologies are shown in Fig. 1. The SCS modes D + s → K + π 0 , D + → π + η and D + s → K + η receive contributions from two different weak phases, proportional to the products of the CKM matrix elements V cd V * ud and V cs V * us , allowing CP violation at tree-level. In the Standard Model (SM), the CP asymmetries are expected to be of the order 10 −4 -10 −3 [2][3][4][5][6][7]. The CF mode D + s → π + η and the DCS modes D + → K + π 0 and D + → K + η receive contributions from only one weak phase at tree-level. The D + s → π + π 0 mode proceeds via an annihilation topology decay and is therefore highly suppressed.  The SCS D + → π + π 0 mode is of particular interest as the CP asymmetry in the SM is expected to be zero as a result of isospin constraints [3][4][5][6]. The CP asymmetries of the signal decays are defined to be where Γ is the partial decay rate and h 0 denotes either a π 0 or an η meson. A non-zero value of A CP (D + → π + π 0 ), coupled with a verification that the isospin sum rule B 00 1 Inclusion of charge conjugated processes is implied throughout, except when discussing asymmetry definitions.
In this article measurements of CP asymmetries of seven D + (s) → h + π 0 and D + (s) → h + η modes are performed, using samples corresponding to either 9 fb −1 or 6 fb −1 of integrated luminosity, respectively, collected by the LHCb experiment in proton-proton (pp) collisions at the LHC. The 6 fb −1 data set comprises data collected during 2015-2018 (Run 2) at a centre-of-mass energy of 13 TeV, whilst the 9 fb −1 data set additionally includes data collected during 2011-2012 (Run 1) at centre-of-mass energies of 7 TeV and 8 TeV. The neutral π 0 and η mesons are reconstructed via decays to the e + e − γ final state. The reconstruction of electron and positron tracks, in addition to the charged hadron track from the D + (s) meson decay, enables the determination of the displaced D + (s) meson decay vertex and suppresses background from particles originating from the primary pp interaction. The signal receives contributions from the suppressed three-body Dalitz decays π 0 → e + e − γ and η → e + e − γ with branching fractions (1.174 ± 0.035)% and (6.9 ± 0.4) × 10 −3 , respectively [11], as well as the more common π 0 → γγ and η → γγ decays with branching fractions (98.823 ± 0.034)% and (39.41 ± 0.20)% [11], where one of the photons subsequently interacts with the detector material and is converted to an e + e − pair. Converted photons have been previously exploited at LHCb [12][13][14][15][16], but this is the first measurement to use converted photons to reconstruct π 0 and η mesons.
Experimentally, the raw asymmetry of each signal mode is measured, which is defined to be where N is the signal yield. This can be approximated by where A Prod (D + (s) ) and A Det (h + ) represent the production and detection asymmetries of the corresponding hadrons. In order to cancel the production and detection asymmetries, the raw asymmetry of D + (s) → K 0 S h + control decays is subtracted, approximated by where the extra term A Mix (K 0 ) arises due to the CP asymmetry induced by mixing and decay of the neutral K 0 S meson [17]. As the nuisance asymmetries are known to be kinematically dependent, the D + (s) → K 0 S h + samples are weighted to match the kinematic distributions of the signal candidates to optimally reduce the impact of the production and detection asymmetries. The CP asymmetry for the signal modes can then be determined as where A w Raw represents the raw asymmetry determined from weighted samples, the values of A CP (D + (s) → K 0 S h + ) are accounted for using external inputs with sub-percent precision [18], and A Mix (K 0 ) is calculated using a description of the detector material and the distribution of K 0 S decay times and momentum in the selected data, as detailed in Refs. [18,19]. This article is structured as follows: the LHCb experiment is described in Section 2; the requirements used to reconstruct the signal samples are given in Section 3; a description of the fits to the invariant mass distributions can be found in Section 4; the treatment of the D + (s) → K 0 S h + control modes is given in Section 5; the sources of systematic uncertainty are detailed in Section 6; and finally the results and conclusions are summarised in Section 7.

Detector
The LHCb detector [20,21] is a single-arm forward spectrometer covering the pseudorapidity range between 2 and 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 (VELO), 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 tracking system provides a measurement of the 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 pp collision vertex (PV), the impact parameter (IP), 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 detectors. Photons, electrons and hadrons are identified by a calorimeter system consisting of scintillating-pad and preshower detectors, an electromagnetic and a hadronic calorimeter. Muons are identified by a system composed of alternating layers of iron and multiwire proportional chambers. The online event selection is performed by a trigger, 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.
Simulation is required to determine the invariant-mass distributions of the signal decays, develop the selection and constrain the yields of background from other particles misidentified as the signal-decay products. In the simulation, pp collisions are generated using Pythia [22] with a specific LHCb configuration [23]. Decays of unstable particles are described by EvtGen [24], in which final-state radiation is generated using Photos [25]. The interaction of the generated particles with the detector, and its response, are implemented using the Geant4 toolkit [26] as described in Ref. [27]. The underlying pp interaction is reused multiple times, with an independently generated signal decay for each [28].

Event selection
To reconstruct the D + (s) meson candidate a well-identified kaon or pion track is combined with a neutral meson to form a secondary decay vertex displaced from any PV. The neutral π 0 and η candidates are formed from two oppositely charged electron tracks that are combined with a photon candidate to create a neutral-meson decay vertex. A bremsstrahlung-recovery algorithm associates additional deposits from soft photons to those produced by the electrons in the electromagnetic calorimeter. To improve the resolution, the electron tracks must include a track segment within the VELO.
At the hardware trigger level, candidates are selected by either directly identifying high transverse-momentum deposits from the signal in the electromagnetic or hadronic calorimeters, or by independently identifying another energetic particle produced in the pp collision. Inclusive multivariate (MVA) software triggers ensure the presence of well-reconstructed tracks that are inconsistent with originating from any PV. A second high-level software trigger performs a full event reconstruction to form the D + (s) candidates. In Run 1, no dedicated exclusive triggers for the signal modes were implemented, but small samples of D + (s) → h + π 0 candidates are reconstructable as a result of the overlap with existing exclusive two-and three-body D-meson-decay triggers. No attempt is made to reconstruct D + (s) → h + η candidates using the Run 1 data set. In Run 2, dedicated exclusive software triggers were added to form both D + (s) → h + π 0 and D + (s) → h + η signal candidates. These require the presence of a photon and three well-reconstructed tracks, inconsistent with originating from any PV. The invariant masses of the π 0 (η) meson candidates are required to be in the range 70 < m(e + e − γ) < 210 MeV/c 2 (450 < m(e + e − γ) < 650 MeV/c 2 ) with p T > 200 MeV/c (500 MeV/c). The D + (s) candidate is required to have a transverse momentum p T > 3000 MeV/c and a good quality vertex with an associated p-value of greater than 0.0018, created by first combining the e + e − candidates to form the photon conversion or h 0 Dalitz decay vertex, which is then further combined with a photon and charged hadron to create the D + (s) decay vertex. Offline, the D + (s) candidate selection is refined by requiring that the momentum of the tracks is in the range 3 < p < 100 GeV/c and their pseudorapidity is between 1.5 and 5.0. The D + (s) candidates are required to have a mass in the range 1600 < m(D + (s) ) < 2200 MeV/c 2 , be consistent with originating at a primary interaction and have a proper decay time of t > 0.15 ps (0.25 ps) for D + (s) → h + π 0 (D + (s) → h + η) candidates. Additionally, the angle between the momentum direction and the vector joining the PV and D + (s) decay vertex, referred to as the direction angle, must be smaller than 10 mrad.
Fiducial requirements are placed on the charged-hadron tracks to remove regions of large detection asymmetries, for example regions where a track of one charge would be bent out of the acceptance by the magnetic field whilst the opposite charge would be detected; the same criteria are used as in the previous measurements of the control modes [18].
Particle identification (PID) requirements are applied using MVA-based PID variables for the charged particles and the photon to reduce the amount of combinatorial and misidentification background [29,30]. Loose PID requirements are applied to the pion and electron tracks. Tighter requirements are applied to kaon candidates to reduce the rate of π + → K + misidentification from the more abundant pion modes into the suppressed kaon modes. When reconstructing π 0 mesons, a loose requirement is placed on an MVA-based photon-quality variable [31], whilst for η mesons, a tighter condition is required to reduce the level of combinatorial background. Requirements are placed on electron bremsstrahlung PID variables that match the bremsstrahlung calorimeter deposit to the electron track before passing through the magnetic field to ensure that the correct photon deposits are recovered. Decays with a total of either zero or one bremsstrahlung photon per e + e − pair are used in this analysis. For D + → h + π 0 (D + → h + η) decays this corresponds to 62% and 38% (31% and 48%) of the reconstructed candidates, respectively. Decays with two or more bremsstrahlung photons per e + e − pair are removed as they result in a poor D + (s) invariant mass resolution and high background level. The same offline selection requirements are used for candidates selected with different numbers of bremsstrahlung photons, and also between candidates decaying via photon conversions or three-body h 0 → e + e − γ decays. The requirements give a reasonable compromise between the efficiency of each type of decay with efficiencies of the order O(10 −6 ) in Run 1 and O(10 −5 ) in Run 2.
After the full selection has been applied, approximately 3% (2%) of events are found to have multiple D + (s) → h + π 0 (D + (s) → h + η) candidates predominately due to combinations with alternative photon candidates, of which all are retained. The signal decays are found to be dominated by π 0 → γγ and η → γγ decays followed by a photon conversion, rather than the three-body Dalitz decays π 0 → e + e − γ and η → e + e − γ, with approximately 86% of the candidates resulting from photon conversions.

Signal modes and fit model
The raw asymmetries of the signal modes are measured using two-dimensional extended simultaneous unbinned maximum-likelihood fits to the invariant mass m(e + e − γ) and corresponds to the known π 0 and η masses [11]. The quantity m(h + h 0 ) is constructed to reduce the correlations between the two dimensions, and is referred to as the D + The fits are performed simultaneously on candidates in categories that depend on the running period, the presence of bremsstrahlung photons, charged-hadron type (pion or kaon) and the candidate charge. All D + (s) → h + η candidates were collected during Run 2. The D + (s) → h + π 0 candidates are split into three running period categories, 2011, 2012 and Run 2, where the centre-of-mass energies were 7, 8, and 13 TeV, respectively. Candidates with either zero or one bremsstrahlung photon per e + e − pair are split into two categories as they have different mass resolutions. The fits are performed on candidates with π + and K + mesons simultaneously to allow the signal yields in either category to determine the misidentification-background yields in the corresponding category.
Two-dimensional probability density functions (PDFs) are used to model different contributions within the mass windows. These contributions can be categorised as signal decays, misidentification background, partially reconstructed low-mass background and combinatorial background. The sum of positively-and negatively-charged candidate yields and raw asymmetry of all signal and background components are free to vary in the fits. A component for D + s → π + π 0 signal is included in the fit, but due to the insignificant yield no corresponding raw asymmetry is measured. The PDFs are assumed to be the same for positively and negatively charged candidates, but otherwise allowed to differ for the other categories of the simultaneous fit. In the fit to D + (s) → h + π 0 candidates the same raw asymmetries are shared between different running periods. The signal modes are modelled by the sum of a two-dimensional Gaussian function and two two-dimensional Crystal Ball functions [32]. The shape parameters and fraction of each function are determined from fits to simulated decays passing the full selection. To account for residual correlations between m(h + h 0 ) and m(e + e − γ) resulting in part from radiative tails, the mean h 0 (D + (s) ) mass is allowed to vary quadratically as a function of the D + (s) (h 0 ) mass in the fits to D + (s) → h + π 0 (D + (s) → h + η) candidates. When performing fits to data, freely varying scaling factors are applied to the widths of the PDFs, and freely varying offsets are added to the mean positions and quadratic correlation coefficients to account for differences between data and simulation. Different parameters are introduced for each running period and bremsstrahlung category. When determining PDFs from simulated decays, the candidates are weighted to account for the PID requirements using input from calibration samples [30].
The fit model accounts for misidentified signal decays, where a π + track has been incorrectly assigned the K + mass hypothesis, or vice versa, using the same two-dimensional parameterisation as the signal shapes. The PDF parameters are determined from fits to the corresponding simulated signal decays passing the full selection for the charged hadron with the wrong mass hypothesis, including weights to account for the misidentifi- cation probabilities. When performing the fits to data, the yield of the misidentification background is constrained to the yield of signal in the other charged-hadron category multiplied by the relevant ratio of efficiencies determined from simulated decays and PID calibration samples. The yields of misidentification background contributions are below approximately 3% of the corresponding signal yields.
Combinatorial background resulting from random combinations of tracks and photons is modelled with an exponential function in the m(h + h 0 ) dimension and a second-order Chebychev polynomial function in the m(e + e − γ) dimension. The exponential coefficient and Chebychev polynomial coefficients freely vary in the fit. In the fit to D + (s) → h + π 0 candidates, it is found necessary to include a combinatorial component comprising a real π 0 meson combined with an unrelated track. The PDF is constructed from a peaking distribution in the m(e + e − γ) dimension and an exponential function in the m(h + h 0 ) dimension. The peaking distribution is constructed from the sum of two Crystal Ball functions, whose shape is determined from one-dimensional fits to the simulated signal decays. However, when fitting data a freely varying mass offset and resolution scaling factor are included to allow the π 0 mass distribution to differ from that of the signal decays. No significant contribution from combinatorial decays with a real η meson and D + → π + π 0 740 ± 60 2 240 ± 120 25 750 ± 430 −1.64 ± 0.93 D + s → π + π 0 20 ± 30 −50 ± 50 450 ± 120 -D + → K + π 0 10 ± 13 90 ± 30 2 440 ± 110 −2.53 ± 4.75 D + s → K + π 0 54 ± 13 150 ± 30 2 580 ± 90 −0.25 ± 3.87 an unrelated track is found when fitting D + (s) → h + η decays, therefore no corresponding component is included.
Decays of charm mesons to h + h 0 X final states, where X is at least one unreconstructed particle, appear as partially reconstructed background below the D + (s) meson masses. Using external input on branching fractions and charm-meson production crosssection ratios [11,33] it is determined that only the decay D + s → π + ηπ 0 has a significant contribution in the fit to D + (s) → π + η candidates. To account for this component, a shape comprising an exponential function in the m(h + η) dimension with a freely varying coefficient and a peaking m(e + e − γ) distribution constructed from two Crystal Ball functions is added.
The fit to D + (s) → h + π 0 (D + (s) → h + η) candidates includes 91 (54) freely varying parameters. The models are validated using pseudo-experiments and no significant biases in the values or statistical uncertainties of the raw asymmetries are observed. The projections of the fits to D + (s) → h + π 0 and D + (s) → h + η candidates, summed over all relevant categories, are shown in Figs. 2 and 3, respectively. The pull distributions are examined for each category of the fit in both projections and in two dimensions, and no significant biases are seen. The goodness-of-fit is quantified by calculating the χ 2 value for each projection and category separately, and combining to determine χ 2 /N dof = 0.90 and χ 2 /N dof = 1.06 for the fits to D + (s) → h + π 0 and D + (s) → h + η candidates, where N dof is the total number of degrees of freedom. The corresponding signal yields and raw asymmetries are listed in Table 1. The D + and D + s signal distributions overlap, leading to small correlations between the measured raw asymmetries. The correlation coefficients are listed in Table 2 and the largest correlation is 10%.

Control modes
The impact of production and detection asymmetries of the signal modes is accounted for using large samples of D + (s) → K 0 S h + decays. The samples are selected using similar requirements to the signal modes, where possible. Candidates are built at the high-level software trigger stage by first combining two well-reconstructed hadronic tracks that are D + → π + π 0 D + → K + π 0 D + s → K + π 0 D + → π + π 0 1.00 D + → K + π 0 −0.01 1.00 D + s → K + π 0 −0.09 0.10 1.00 inconsistent with originating from any PV to create the K 0 S decay vertex. Similar to the electrons, these tracks must also have track segments within the VELO. The K 0 S candidate is combined with a hadronic track with either the pion or kaon mass hypothesis to form the D + (s) decay vertex. The same momentum, pseudorapidity and fiducial requirements are placed on the tracks as used for the signal. The candidates are required to have 482 < m(π + π − ) < 512 MeV/c 2 and 1800 < m(K 0 S h + ) < 2050 MeV/c 2 , a proper decay time of t > 0.25 ps, and the same direction angle and p T requirements as the signal. Tighter PID requirements are placed on the control mode candidates than the signal to remove larger contamination from misidentification background.
The kinematic distributions of the signal and control candidates are determined using the sPlot technique [34] with m(K 0 S h + ) as the discriminating variable for the latter. Binned maximum-likelihood fits are performed on the control-mode candidates using signal models comprising a Gaussian function and Johnson S U function [35] as described in Ref. [18]. The results are shown in Fig. 4. The weighting procedure is performed separately for Run 1 and Run 2 to allow for differences in the signal selection during these periods. To ensure the cancellation of the production and detection asymmetries, the relevant D + (s) and h + kinematics (p, azimuthal angle and pseudorapidity) are weighted to match those of the signal. Due to the large correlation between the D + (s) and h + kinematics the weights for each variable are determined using a two-dimensional binning of the D + (s) and h + distributions. In addition to the kinematics, weights are determined for the trigger category and IP distributions for the D + (s) candidates. At the hardware trigger stage the candidates can be split into exclusive categories according to the origin of the positive trigger decision: the first category contains any candidate with a calorimeter deposit associated to the h 0 or K 0 S decay; the second category contains any remaining candidate with a deposit not associated to any of the signal particles; and the third category contains candidates still remaining with a high p T deposit associated to the charged pion or kaon. The control-mode candidates are weighted to reproduce the populations of signal candidates in each of these three categories.
The IP of the D + (s) candidate is indicative of whether the meson was produced in the primary interaction, or as a product of a b-hadron decay, and therefore with a significant IP with respect to the PV, referred to as a secondary decay. In the latter case the production asymmetry of the parent b-hadron could differ from that of the D + or D + s meson. The signal and control mode selections require that the D + (s) candidates are consistent with originating at a PV, suppressing the fraction of candidates from secondary decays to less than 10%. If the fraction of D + (s) candidates from the primary interaction and secondary decays varies between the signal and control mode then the production asymmetries may not exactly cancel, therefore the control sample is weighted to match the IP distribution of the signal.
Binned maximum-likelihood fits are performed to the charge-split samples to determine the raw asymmetries separately for Run 1 and Run 2. The signals are described using the sum of a Gaussian function and Johnson S U function, using the same model as described in Ref. [18]. The fits are performed after the samples have been weighted to match the kinematics of the signal modes, and the statistical uncertainty is calculated using the weights to account for the loss of precision resulting from the weighting procedure.

Systematic uncertainties
The systematic uncertainty on the CP asymmetries receives contributions from a number of sources, including the signal and background parameterisations, the control modes and selection requirements. The assumptions used when creating the signal and background parameterisations are varied and the corresponding systematic uncertainty is quantified using the resulting difference in the raw asymmetries in the fits to data. This includes using: alternative signal parameterisation comprising Johnson S U functions instead of Crystal Ball functions; different pure-combinatorial m(h + h 0 ) parameterisations of a constant plus exponential function; alternative pure-combinatorial m(e + e − γ) parameterisations of a third-order Chebychev polynomial function; alternative real-π 0 combinatorial parameterisation of a double Johnson S U ; and different misidentification-background parameterisations using Johnson S U functions instead of Crystal Ball functions. The efficiencies used to constrain the level of misidentification background are varied within the corresponding uncertainties in 100 fits and the spread in the raw asymmetries is used to estimate the systematic uncertainty. The impact on the raw asymmetries is quantified when various neglected background components are included in the model, including semileptonic D + (s) → h 0 e + ν e and D + (s) → h 0 µ + ν µ decays, partially reconstructed D 0 → K − π + π 0 decays and a combinatorial component with a real-η distribution. Additionally, the assumption that the pure-combinatorial m(h + h 0 ) exponential slope is independent of m(e + e − γ) is relaxed by allowing a linear dependence. The signal tail parameters that are fixed to values obtained from simulation are allowed to vary with an overall scaling factor and the impact on the raw asymmetries is quantified. The assumption that the mean D + (s) mass positions are the same for D + (s) and D − (s) candidates is tested by allowing different values. The systematic uncertainty from the fit model is dominated by the fixed tail parameters for D + → π + π 0 , the fixed misidentification efficiency ratio for D + → K + π 0 and D + s → K + π 0 decays, the signal parameterisation for D + → π + η decays and the lack of real-η combinatorial contribution for D + s → π + η, D + → K + η and D + s → K + η decays. The selection of the signal and control modes uses different requirements for the PID variables. Tighter conditions are needed for the control modes to reduce misidentification background such as Λ + c → pK 0 S decays. The size of a possible charge asymmetry induced by these different requirements is quantified by first computing the asymmetry of the PID efficiencies, PID , when determined separately for positively and negatively charged hadrons, . Then, the difference in PID asymmetry when calculated using signal and control mode PID requirements, ∆A PID = A signal PID −A control PID , is used to quantify the corresponding systematic uncertainty. Additionally, the difference in the raw asymmetries when not performing the IP weighting is used to quantify the systematic uncertainty arising from the secondary decays.
The asymmetries for D + (s) → h + π 0 decays are determined from simultaneous fits to data sets taken during Run 1 and Run 2, with a single CP asymmetry shared between the categories for each mode. In contrast, the control-mode fits are performed separately for Run 1 and Run 2 and then a weighted average is performed to combine the measurements, where the weighting is determined from the yields of signal mode decays. The systematic uncertainty arising from this method is quantified by performing the signal fits separately for Run 1 and Run 2, taking the appropriate difference with the control-mode asymmetries and then combining the Run 1 and Run 2 results to get an alternative estimate.
The control-mode weighting is performed in nearly equally populated bins. The binning scheme is varied to determine the associated systematic uncertainty. After performing the weighting procedure, the remaining discrepancies in the kinematic distributions are quantified by summing the difference in the normalised distributions of signal and control modes, multiplied by the local asymmetry minus the average asymmetry. The fit model used to measure the control-mode raw asymmetries is varied from the sum of a Johnson S U function and a Gaussian function to the sum of a Crystal Ball function and a Gaussian function. The contribution to the control mode raw asymmetry from the neutral-kaon mixing and decay asymmetry is calculated and the corresponding uncertainty of this calculation is dominated by the knowledge of the detector material. The uncertainties of the external values of the control mode A CP are included as systematic uncertainties.
The systematic uncertainties are listed for the D + (s) → h + π 0 modes in Table 3 and for the D + (s) → h + η modes in Table 4. These are dominated by the fit-model uncertainty in most cases, except for the mode D + s → π + η which is dominated by the uncertainty arising from the control mode D + s → K 0 S π + , the smallest of the control samples. As a crosscheck, the fits are performed in various subsamples: split by year of data taking; magnet polarity; trigger category; bremsstrahlung category; D + (s) kinematics and h + kinematics. No significant biases are found with respect to the nominal results.  [18,[36][37][38][39][40]. For comparison the unweighted control asymmetries are A Raw (D + → K 0 S π + ) = −0.45 ± 0.02, A Raw (D + → K 0 S K + ) = 0.47 ± 0.05 and A Raw (D + s → K 0 S K + ) = 0.51 ± 0.04.

Results and conclusions
The CP asymmetries are calculated using Eq. 6, where for each mode the corresponding control channel A w Raw , independently measured A CP (D + (s) → K 0 S h + ) and calculated A Mix (K 0 ) are taken. The final results are listed in Tables 5 and 6. The results are shown with the corresponding statistical uncertainty from the fits and the total systematic uncertainty as listed in Tables 3 and 4. The systematic uncertainties attributed to and A Mix (K 0 ) are listed separately. In summary, measurements of CP asymmetries in D + (s) → h + π 0 and D + (s) → h + η decays are performed using pp collision data corresponding to 9 fb −1 and 6 fb −1 of integrated luminosity collected at the LHCb experiment, respectively. The neutral mesons are reconstructed using the e + e − γ final state, allowing the D + (s) decay vertex to be reconstructed. The production and detection asymmetries are cancelled using large samples of D + (s) → K 0 S h + decays, weighted to match the kinematics of the signal modes. The CP asymmetries are determined to be A CP (D + → π + π 0 ) = (−1.3 ± 0.9 ± 0.6)%, A CP (D + → K + π 0 ) = (−3.2 ± 4.7 ± 2.1)%, where the first uncertainty is statistical and the second systematic. All of the results are consistent with no CP asymmetry and the first five constitute the most precise measurements to date. Very recently the Belle collaboration has also reported precise measurements of A CP (D + s → K + π 0 ), A CP (D + s → π + η) and A CP (D + s → K + η) [41]. The result for A CP (D + → π + π 0 ) is consistent with the SM expectation and the previous measurement by the Belle collaboration [10]. Using the relevant lifetimes, branching  [18,[36][37][38][39][40]. For comparison the unweighted control asymmetries are A Raw (D + → K 0 S π + ) = −0.45 ± 0.02, A Raw (D + s → K 0 S π + ) = −0.13 ± 0.17, A Raw (D + → K 0 S K + ) = 0.47 ± 0.05 and A Raw (D + s → K 0 S K + ) = 0.51 ± 0.04.