Measurements of the absolute branching fractions and $CP$ asymmetries for $D^+\rightarrow K_{S,L}^0 K^+(\pi^0)$

Using $e^+e^-$ collision data corresponding to an integrated luminosity of 2.93 fb$^{-1}$ taken at a center-of-mass energy of 3.773 GeV with the BESIII detector, we determine the absolute branching fractions $\mathcal{B}(D^+\rightarrow K_S^0 K^+)$ = $(3.02 \pm 0.09 \pm 0.08)\times10^{-3}$, $\mathcal{B}(D^+\rightarrow K_S^0 K^+ \pi^0)$ = $(5.07 \pm 0.19 \pm 0.23)\times10^{-3}$, $\mathcal{B}(D^+\rightarrow K_L^0 K^+)$ = $(3.21 \pm 0.11 \pm 0.11)\times10^{-3}$, and $\mathcal{B}(D^+\rightarrow K_L^0 K^+ \pi^0)$ = $(5.24 \pm 0.22 \pm 0.22)\times10^{-3}$, where the first and second uncertainties are statistical and systematic, respectively. The branching fraction \wbrksk is consistent with the world average value and the other three branching fractions are measured for the first time. We also measure the $CP$ asymmetries for the four decays and do not find a significant deviation from zero.


I. INTRODUCTION
Experimental studies of hadronic decays of charm mesons shed light on the interplay between the strong and weak forces. In the standard model (SM), the singly Cabibbo-suppressed (SCS) D ± meson hadronic decays are predicted to exhibit CP asymmetries of the order of 10 −3 [1]. Direct CP violation in SCS D ± meson decays can arise from the interference between tree-level and penguin decay processes [2]. However, the doubly Cabibbo-suppressed and Cabibbo-favored D ± meson decays are expected to be CP invariant because they are dominated by a single weak amplitude. Consequently, any observation of CP asymmetry greater than O(10 −3 ) in the SCS D ± meson hadronic decays would be evidence for new physics beyond the SM [3]. In theory, the branching fractions of two-body hadronic decays of D mesons can be calculated within SU(3) flavor symmetry [4]. An improved measurement of the branching fraction of the SCS decay D + →K 0 K + will help to test the theoretical calculations and benefit the understanding of the violation of SU(3) flavor symmetry in D meson decays [4]. In this paper, we present measurements of the absolute branching fractions and the direct CP asymmetries of the SCS decays of D + → K 0 S K + , K 0 S K + π 0 , K 0 L K + and K 0 L K + π 0 . In this analysis, we employ the "double-tag" (DT) technique, which was first developed by the MARK-III Collaboration [5,6], to measure the absolute branching fractions. First, we select "single-tag" (ST) events in which either a D orD meson is fully reconstructed in one of several specific hadronic decays. Then we look for the D meson decays of interest in the presence of the ST D events; the so called the DT events in which both the D andD mesons are fully reconstructed. The ST and DT yields (N ST and N DT ) can be described by where N D + D − is the total number of D + D − pairs produced in data, ǫ ST and ǫ DT are the efficiencies of reconstructing the ST and DT candidate events, and B tag and B sig are the branching fractions for the tag mode and the signal mode, respectively. The absolute branching fraction for the signal decay can be determined by where ǫ = ǫ DT /ǫ ST is the efficiency of finding a signal candidate in the presence of a STD, which can be obtained from MC simulations. With the measured absolute branching fractions of D + and D − meson decays (B + sig and B − sig ), the CP asymmetry for the decay of interest can be determined by (3)

II. THE BESIII DETECTOR AND DATA SAMPLE
The analysis presented in this paper is based on a data sample with an integrated luminosity of 2.93 fb −1 [7] collected with the BESIII detector [8] at the center-of-mass (c.m.) energy of √ s = 3.773 GeV. The BESIII detector is a general-purpose detector at the BEPCII [9] with double storage rings. The detector has a geometrical acceptance of 93% of the full solid angle. We briefly describe the components of BESIII from the interaction point (IP) outward. A small-cell multi-layer drift chamber (MDC), using a helium-based gas to measure momenta and specific ionization of charged particles, is surrounded by a time-of-flight (TOF) system based on plastic scintillators which determines the time of flight of charged particles. A CsI(Tl) electromagnetic calorimeter (EMC) detects electromagnetic showers. These components are all situated inside a superconducting solenoid magnet, which provides a 1.0 T magnetic field parallel to the beam direction. Finally, a multilayer resistive plate counter system installed in the iron flux return yoke of the magnet is used to track muons. The momentum resolution for charged tracks in the MDC is 0.5% for a transverse momentum of 1 GeV/c. The specific energy loss (dE/dx) measured in the MDC has a resolution better than 6%. The TOF can measure the flight time of charged particles with a time resolution of 80 ps in the barrel and 110 ps in the end caps. The energy resolution for the EMC is 2.5% in the barrel and 5.0% in the end caps for photons and electrons with an energy of 1 GeV. The position resolution of the EMC is 6 mm in the barrel and 9 mm in the end caps. More details on the features and capabilities of BESIII can be found elsewhere [8].
A geant4-based [10] Monte Carlo (MC) simulation software package, which includes the geometric description of the detector and its response, is used to determine the detector efficiency and to estimate potential backgrounds. An inclusive MC sample, which includes the D 0D0 , D + D − , and non-DD decays of ψ(3770), the initial state radiation (ISR) production of ψ(3686) and J/ψ, the qq (q = u, d, s) continuum process, Bhabha scattering events, and di-muon and di-tau events, is produced at √ s = 3.773 GeV. The kkmc [11] package, which incorporates the beam energy spread and the ISR effects (radiative corrections up to next to leading order), is used to generate the ψ(3770) meson. Final state radiation of charged tracks is simulated with the photos package [12]. ψ(3770) → DD events are generated using evtgen [13,14], and each D meson is allowed to decay according to the branching fractions in the Particle Data Group (PDG) [15]. This sample is referred as the "generic" MC sample. Another MC sample of ψ(3770) → DD events, in which one D meson decays to the signal mode and the other one decays to any of the ST modes, is referred as the "signal" MC sample. In both the generic and signal MC samples, the two-body decays D + → K 0 S,L K + are generated with a phase space model, while the threebody decays D + → K 0 S,L K + π 0 are generated as a mixture of known intermediate decays with fractions taken from the Dalitz plot analysis of their charge conjugated decay D + → K + K − π + [16].

III. DATA ANALYSIS
The ST D ∓ mesons are reconstructed using six hadronic final states: S is reconstructed by its π + π − decay mode and π 0 with the γγ final state. The event selection criteria are described below.
Charged tracks are reconstructed within the MDC coverage | cos θ| < 0.93, where θ is the polar angle with respect to the positron beam direction. Tracks (except for those from K 0 S decays) are required to have a point of closest approach to the IP satisfying |V z | < 10 cm in the beam direction and |V r | < 1 cm in the plane perpendicular to the beam direction. Particle identification (PID) is performed by combining the information of dE/dx in the MDC and the flight time obtained from the TOF. For a charged π(K) candidate, the probability of the π(K) hypothesis is required to be larger than that of the K(π) hypothesis.
The K 0 S candidates are reconstructed from combinations of two tracks with opposite charges which satisfy |V z | < 20 cm, but without requirement on |V r |. The two charged tracks are assumed to be π + π − without PID and are constrained to originate from a common decay vertex. The π + π − invariant masses M π + π − are required to satisfy |M π + π − − M K 0 [15]. Finally, the K 0 S candidates are required to have a decay length significance L/σ L of more than two standard deviations, as obtained from the vertex fit.
Photon candidates are selected from isolated showers in the EMC with minimum energy larger than 25 MeV in the barrel region (| cos θ| < 0.80) or 50 MeV in the endcap region (0.86 < | cos θ| < 0.92). The shower timing is required to be no later than 700 ns after the event start time to suppress electronic noise and energy deposits unrelated to the event.
The π 0 candidates are reconstructed from pairs of photon candidates with invariant mass within 0.110 < M γγ < 0.155 GeV/c 2 . The γγ invariant mass is then constrained to the nominal π 0 mass [15] by a kinematic fit, and the corresponding χ 2 is required to be less than 20.

A. ST yields
The ST D ∓ candidates are formed by the combinations of K ± π ∓ π ∓ , K ± π ∓ π ∓ π 0 , K 0 S π ∓ , K 0 S π ∓ π 0 , K 0 S π ± π ∓ π ∓ and K ± K ∓ π ∓ . Two variables are used to identify ST D mesons: the energy difference ∆E and the beam-energy constrained mass M BC , which are defined as Here, p D and E D are the reconstructed momentum and energy of the D candidate in the e + e − c.m. system, and E beam is the beam energy. Signal events are expected to peak around zero in the ∆E distribution and around the nominal D mass in the M BC distribution. In the case of multiple candidates in one event, the one with the smallest |∆E| is chosen. Tag mode-dependent ∆E requirements as used in Ref. [17] are imposed on the accepted ST candidate events, as summarized in Table I. To obtain the ST yield for each tag mode in data, a binned maximum likelihood fit is performed on the M BC distribution, where the signal of D meson is described by a MC-simulated shape and the background is modeled by an ARGUS function [18]. The MC-simulated shape is convolved with a Gaussian function with free parameters to take into account the resolution difference between data and MC simulation. Figures 1 and 2 illustrate the resulting fits to the M BC distributions for ST D + and D − candidate events in data, respectively. The fitted ST yields of data are presented in Table I, too.

B. DT yields
On the recoiling side against the ST D ∓ mesons, the hadronic decays of D ± → K 0 S,L K ± (π 0 ) are selected using the remaining tracks and neutral clusters. The charged kaon is required to have the same charge as the signal D meson candidate. To suppress backgrounds, no extra good charged track is allowed in the DT candidate events. The signal D ± candidates are also identified with the energy difference and the beam energy constrained mass. In the following, the energy difference and the beam-energy constrained mass of the particle combination for the ST/signal side are denoted as ∆E tag/sig and M tag/sig BC , respectively. In each event, if there are multiple signal candidates for D ± → K 0 S K ± (π 0 ), the one with the smallest |∆E sig | is selected. The ∆E sig is required to be within (−0.031, 0.031) GeV and (−0.057, 0.040) GeV for D ± → K 0 S K ± and D ± → K 0 S K ± π 0 , respectively. Due to its long lifetime, very few K 0 L decay in the MDC. However, most K 0 L will interact in the material of the EMC, which gives their position but no reliable measurement of their energy. Thus, to select the candidates of D ± → K 0 L K ± (π 0 ), the momentum direction of the K 0 L particle is inferred by the position of a shower in the EMC, and a kinematic fit imposing momentum and energy conversation for the observed particles and a missing K 0 L particle is performed to select the signal, where the K 0 L particle is of known mass and momentum direction, but of unknown momentum magnitude. We perform the kinematic fit individually for all shower candidates in the EMC that are not used in the ST side and do not form a π 0 candidate with any other shower candidate with invariant mass within (0.110, 0.155) GeV/c 2 [17]. The candidate with the minimal chi-square of the kinematic fit (χ 2 K 0 L ) is selected. To minimize the correlation between M tag BC and M sig BC , the momentum of the K 0 L candidate is not taken from the kinematic fit, but inferred by constraining ∆E sig to be zero. In order to suppress backgrounds due to cluster candidates produced mainly from electronics noise, the energy of the K 0 L shower in the EMC is required to be greater than 0.1 GeV. Finally, DT candidate events are imposed with the optimized, and ST and signal mode dependent χ 2 K 0 L requirements, as summarized in Table II.  Figure 3 illustrates the distribution of M tag BC versus M sig BC for the DT candidate events of D + → K 0 S π + , summed over the six ST modes. Candidate signal events concentrate around the intersection of M tag [15]. Candidate events with one correctly reconstructed and one incorrectly reconstructed D meson are spread along the vertical band with M sig BC = M D + or horizontal band with M tag BC = M D + , respectively (named BKGI thereafter). Other candidate events, smeared along the diagonal, are mainly from the continuum process e + e − → qq (named BKGII thereafter). To determine the DT signal yield, we perform an unbinned two-dimensional (2D) maximum likelihood fit on the distribution of M tag BC versus M sig BC of the selected events. In the fit, the probability density functions for the signal, BKGI and BKGII components are constructed as follows: • Signal: a(x, y) ⊗ g(x; x 0 , σ x0 ) ⊗ g(y; y 0 , σ y0 ), . The signal is described with a MC-simulated shape a(x, y) convolved with two independent Gaussian functions g(x; x 0 , σ x0 ) and g(y; y 0 , σ y0 ) representing the resolution difference between data and MC simulation in the variables M sig BC and M tag BC , respectively. The parameters of the Gaussian functions x 0 , σ x0 , y 0 , σ y0 are determined by performing 1D fits on the M sig BC and M tag BC distributions of data, individually. The shape of BKGI b(x, y) is determined from the generic DD MC sample. In particular, in the studies of D + → K 0 L K + (π 0 ), irreducible and peaking backgrounds come mainly from D + → K 0 S K + (π 0 ) with K 0 S → π 0 π 0 . Since their shape is too similar to be separated from the signal in the fit, their size and shape are fixed. To take into account possible differences between data and MC simulation, both shapes and magnitudes of the D + → K 0 S K + (π 0 ) background events are re-estimated as follows. The background shapes are determined by imposing the same selection criteria as for  data on the MC samples of D + → K 0 S K + (π 0 ) with K 0 S decaying inclusively. The background magnitudes are estimated by using the samples of D + → K 0 S K + (π 0 ) with K 0 S → π 0 π 0 selected from data and MC samples, from which the event yields N DT The 2D fit is performed on the M sig BC versus M tag BC distribution for each ST mode individually. Figure 4 shows the projections on the M sig BC and M tag BC distributions of the 2D fits summed over all six ST modes. The detection efficiencies of D ± → K 0 S,L K ± (π 0 ) are determined by MC simulation. In our previous work [17], differences of the K 0 S,L reconstruction efficiencies between data and MC simulation (called data-MC difference) were found, due to differences in nuclear interactions of K 0 andK 0 mesons. The detection efficiencies were investigated for K 0 → K 0 S,L andK 0 → K 0 S,L separately. To compen-sate for these differences, the signal efficiencies are corrected by the K 0 S,L momentum-weighted data-MC differences of the K 0 S,L reconstruction efficiencies. The efficiency correction factors are about 2% and 10% for D ± → K 0 S K ± (π 0 ) and D ± → K 0 L K ± (π 0 ), respectively. The DT signal yields in data (N DT ) and the corrected detection efficiencies (ǫ) of D ± → K 0 S,L K ± (π 0 ) are presented in Table III. C. Branching fraction and CP asymmetry According to Eq. (2) and taking into account the numbers of N ST , N DT , and ǫ listed in Tables I and III, the branching fractions of D + and D − decays for the individual ST modes are calculated. The average branching fractions of D + and D − decays as well as combination of charged conjugation modes are obtained by using the standard weighted least-squares method [15], and are summarized in Table IV. We also determine the CP asymmetries with Eq. (3) based on the average branching fractions of D + and D − decays, and the results are listed in Table IV, too.

IV. SYSTEMATIC UNCERTAINTY
Due to the use of the DT method, those uncertainties associated with the ST selection are cancelled. The relative systematic uncertainties in the measurements of absolute branching fractions and the CP asymmetries of the decay D ± → K 0 S,L K ± (π 0 ) are summarized in Table V and are discussed in detail below.
The efficiencies of K ± tracking and PID in various K ± momentum ranges are investigated with K ± samples selected from DT hadronic DD decays. In each momentum range, the data-MC difference of efficiencies ǫ data /ǫ MC −1 is calculated. The data-MC differences weighted by the K ± momentum in the decays D ± → K 0 S,L K ± (π 0 ) are assigned as the associated systematic uncertainties.
The π 0 reconstruction efficiency is studied by the DT control sample D 0 → K − π + π 0 versusD 0 → K + π − or D 0 → K + π − π − π + using the partial reconstruction technique. The data-MC difference of the π 0 reconstruction efficiencies weighted according to the π 0 momentum distribution in D + → K 0 S,L K + π 0 is assigned as the systematic uncertainty in π 0 reconstruction. The branching fractions of K 0 S → π + π − and π 0 → γγ are taken from the Particle Data Group [15]. Their uncertainties are 0.07% and 0.03%, respectively, which are negligible in these measurements.
As described in Ref. [17], the correction factors of K 0 S,L reconstruction efficiencies are determined with the two control samples of J/ψ → K * (892) ∓ K ± with K * (892) ∓ → K 0 S,L π ∓ and J/ψ → φK 0 S,L K ± π ∓ decays. Since the efficiency corrections are imposed in this analysis, the corresponding statistical uncertainties of the cor-TABLE III. DT yields in data (NDT) and efficiencies (ǫ) of reconstructing the signal decays, where the uncertainties are statistical only. The efficiencies include the branching fractions for K 0 S → π + π − and π 0 → γγ.
As described in Ref. [17], in the determination of the correction factor of the K 0 L efficiency, we perform a kinematic fit to select the K 0 L candidate with the minimal χ 2  Table II brings an uncertainty. The momentum-weighted uncertainty of the χ 2 K 0 L selection according to the K 0 L momentum distribution of signal events is assigned as the associated systematic uncertainties.
In the analysis of multi-body decays, the detection efficiency may depend on the kinematic variables of the final-state particles. The possible difference of the kinematic variable distribution between data and MC simulation causes an uncertainty on detection efficiency. For the three-body decays D + → K 0 S,L K + π 0 , the nominal efficiencies are estimated by analyzing an MC sample composed of the decays D + → K * (892) +K 0 , D + →K * (892) 0 K + , D + →K * (1430) 0 K + , and D + → K * 2 (1430) 0 K + . The fractions of these components are taken from the Dalitz plot analysis of the charge conjugated decay D + → K + K − π + [16]. The differences of the nominal efficiencies to those estimated with an MC sample of their dominant decays of D + → K * (892) + K 0 S,L [15] are taken as the systematic uncertainties due to the MC model.
To evaluate the systematic uncertainty associated with the ST yields, we repeat the fit on the M BC distribution of ST candidate events by varying the resolution of the Gaussian function by one standard deviation. The resulting change on the ST yields is found to be negligible.
The systematic uncertainties in the 2D fit on the M tag BC versus M sig BC distribution are evaluated by repeating the fit with an alternative fit range (1.8400, 1.8865) GeV/c 2 , varying the resolution of the smearing Gaussian function by one standard deviation, and varying the endpoint of the ARGUS function by ±0.2 MeV/c 2 , individually, and the sum in quadrature of the changes in DT yields are taken as the systematic uncertainties.
As described in Sec. III B, the dominant peaking backgrounds for D ± → K 0 L K ± (π 0 ) are found to be from D ± → K 0 S K ± (π 0 ) with K 0 S → π 0 π 0 , whose contributions are about 3% (5%). Their sizes are estimated based on MC simulation after considering the branching fraction of the background channel and are fixed in the fits. Other peaking backgrounds like D ± → K 0 L π ± (π 0 ) are found to have contributions of less than 0.5%. The uncertainties due to these peaking backgrounds are estimated by varying the branching fractions of the peaking background channels by ±1σ, and the changes of the DT signal yields are assigned as the associated systematic uncertainties.
In the studies of D ± → K 0 S K ± (π 0 ), a ∆E requirement in the signal side is applied to suppress the background. The corresponding uncertainty is studied by comparing the DT yields with and without the ∆E requirement for an ST mode with low background, i.e. D ± → K ∓ π ± π ± . The resulting difference of relative change of DT yields between data and MC simulation is assigned as the systematic uncertainty.
For each signal mode, the total systematic uncertainty of the measured branching fraction is obtained by adding all above individual uncertainties in quadrature, as summarized in Table V. In the determination of the CP asymmetries, the uncertainties arising from π 0 reconstruction, χ 2 K 0 L requirement of the K 0 L selection, MC model of D ± → K 0 S,L K ± π 0 , M BC fit for ST events and ∆E requirement in signal side are canceled. The total systematic uncertainties in the measured CP asymmetries are also listed in Table V.
We also examine the CP asymmetries for the threebody decay D ± → K 0 S,L K ± π 0 in different regions across the Dalitz plot. For this study, a further kinematic fit constraining the masses of K 0 S and D + candidates to their nominal masses [15] is performed in the selection of D ± → K 0 S K ± π 0 . To select signal D ± → K 0 L K ± π 0 events, a kinematic fit constraining the D + to its nominal mass is performed in addition to the kinematic fit to select the K 0 L shower as described in Sec. III B. The recoiling mass of the K 0 S,L K ± π 0 system, which should equal the mass of the ST D meson in correctly reconstructed signal events, is used to identify the signal, where q 0 and q D are the four-momentum of the e + e − system and the selected D + candidate, respectively. This procedure ensures that D candidates have the same phase space (PHSP), regardless of whether M rec is in the signal or sideband region. Figure 5 shows the fits to the M rec distributions and the Dalitz plot of event candidates in the M rec signal region defined as (1.864, 1.877) GeV/c 2 . In the M rec distribution of D ± → K 0 S K ± π 0 , there is a significant tail above the D + mass due to ISR effects. For ISR events in D ± → K 0 L K ± π 0 , the momentum of the K 0 L becomes larger than what it should be due to the constraint of ∆E = 0, which leads to a significant tail below the D ± mass in the M rec distribution. The M rec distributions are fitted with a MC-derived signal shape convolved with a Gaussian function for the signal, together with an AR-GUS function for the combinatorial background.
The Dalitz plot of D ± → K 0 S,L K ± π 0 is further divided into three regions to examine the CP asymmetries. The DT yields in data are obtained by counting the numbers of events in the individual Dalitz plot regions in the M rec signal region, and then subtract the numbers of background events in the M rec sideband regions (shown in Fig. 5). MC studies show that the peaking backgrounds in the study of D ± → K 0 S K ± π 0 are negligible. For the study of D ± → K 0 L K ± π 0 , however, there are non-negligible peaking backgrounds dominated by D ± → K 0 S K ± π 0 with K 0 S → π 0 π 0 . These peaking backgrounds are estimated by MC simulations as described previously and are also subtracted from the data DT yields.
The background-subtracted DT yields in data N DT , the signal efficiencies ǫ, the calculated branching fractions B and the CP asymmetries A CP in the different Dalitz plot regions are summarized in Tables VI and VII. Here, the branching fractions and the CP asymmetries are calculated by Eq. (2) and Eq. (3), respectively. The corresponding systematic uncertainties are assigned after considering the different behaviors of K ± and K 0 S,L reconstruction in the detector. We use the same method as described in Sec. IV to estimate the systematic uncertainties on the CP asymmetries in the individual Dalitz plot regions, all of which are listed in Table VIII. No evidence for CP asymmetry is found in individual regions.

VI. SUMMARY
Using an e + e − collision data sample of 2.93 fb −1 taken at √ s = 3.773 GeV with the BESIII detector, we present the measurements of the absolute branching fraction TABLE VI. Background-subtracted DT yields in data (NDT) and detection efficiencies (ǫ) in different Dalitz plot regions for D ± → K 0 S,L K ± π 0 , where the uncertainties are statistical only.

VII. ACKNOWLEDGMENTS
The BESIII collaboration thanks the staff of BEPCII and the IHEP computing center for their strong support. This work is supported in part by National