Amplitude analysis and branching-fraction measurement of $D_{s}^{+} \to \pi^{+}\pi^{0}\eta^{\prime}$

Using data collected with the BESIII detector in $e^+e^-$ collisions at center-of-mass energies between 4.178 and 4.226 GeV and corresponding to 6.32~fb$^{-1}$ of integrated luminosity, we report the amplitude analysis and branching-fraction measurement of the $D^+_s \to \pi^+ \pi^0 \eta^{\prime}$ decay. We find that the dominant intermediate process is $D^+_s \to\rho^+ \eta^{\prime}$ and the significances of other resonant and nonresonant processes are all less than $3\sigma$. The upper limits on the branching fractions of $S$-wave and $P$-wave nonresonant components are set to $0.10\%$ and $0.74\%$ at the $90\%$ confidence level, respectively. In addition, the branching fraction of the $D^+_s \to \pi^+ \pi^0 \eta^{\prime}$ decay is measured to be $(6.15\pm0.25(\rm stat.)\pm0.18(\rm syst.))\%$, which receives significant contribution only from $D_s^+\to \rho^+\eta^{\prime}$ according to the amplitude analysis.


Introduction
Hadronic decays of the D ± s meson probe the interplay of short-distance weak-decay matrix elements and long distance QCD interactions. Measurements of the branching fractions (BFs) of these decays provide direct knowledge of the amplitudes and phases in the decay process [1][2][3]. In addition, an improved understanding of D ± s decays is particularly valuable for studies of the B 0 s meson, which mainly decays to final states involving D ± s mesons [4]. There are two kinds of topological diagrams for D + s → ρ + η , including tree (T )-and annihilation (A)-diagrams, as shown in Fig. 1 [5]. Based on reference [6], the topological amplitude (A) expressions of D + s → ρ + η, D + s → ρ + η and D + s → π + ω satisfy the sum rule: Here, φ is the mixing angle between η and η : where η q and η s are defined by η q = 1 √ 2 (uu + dd) and η s = ss. Considering the BFs of D + s → π + ω and D + s → ρ + η and noting a simple triangular inequality in Eq. (1.1), one obtains the bounds (2.19 ± 0.27)% < B(D + s → ρ + η ) < (4.51 ± 0.38)% [6]. The predictions of the BF of D + s → ρ + η from several theoretical approaches [7,8] and the corresponding BFs from experimental measurements are shown in Table 1. The theoretical predictions for B(D + s → ρ + η ) are lower than the experimental measurement by around 2σ as shown in Table 1. A possible way to reconcile the predictions with the measured values would be to take account of the QCD flavor-singlet hairpin contribution shown in Fig. 2 [5]. A more precise measurement of the BF of D + s → ρ + η will be very valuable in establishing whether indeed the existing predictions are incorrect.   Experiment D + s → π + π 0 η 5.6 ± 0.5 ± 0.6 CLEO [9] D + s → ρ + η 5.8 ± 1.4 ± 0.4 BESIII [10] D + s → π + π 0 η < 5.1 (nonresonant) (90% confidence level) Previously, BESIII reported the BF measurement of D + s → ρ + η performed through the process e + e − → D + s D − s , with a 482 pb −1 data sample collected at center-of-mass (C.M.) energy √ s = 4.009 GeV and CLEO measured the BF of D + s → π + π 0 η using 586 pb −1 of e + e − collisions recorded at C.M. energy √ s = 4.17 GeV. In this paper, we perform the first amplitude analysis of D + s → π + π 0 η and improve the BF measurement of this decay via the process e + e − → D * ± s D ∓ s by using data samples corresponding to an integrated luminosity of 6.32 fb −1 collected by the BESIII detector at C.M. energies √ s = 4.178 − 4.226 GeV. Charge-conjugate states are implied throughout this paper.

Detector and data sets
The BESIII detector [11] records symmetric e + e − collisions provided by the BEPCII storage ring [12], which operates in the C.M. energy range from 2.00 to 4.95 GeV. BESIII has collected large data samples in this energy region [13]. The cylindrical core of the BESIII detector covers 93% of the full solid angle and consists of a helium-based multilayer drift chamber (MDC), a plastic scintillator time-of-flight system (TOF), and a CsI(Tl) electromagnetic calorimeter (EMC), which are all enclosed in a superconducting solenoidal magnet providing a 1.0 T magnetic field. The solenoid is supported by an octagonal flux-return yoke with resistive plate counter muon identification modules interleaved with steel. The charged-particle momentum resolution at 1 GeV/c is 0.5%, and the specific energy loss (dE/dx) resolution is 6% for electrons from Bhabha scattering. The EMC measures photon energies with a resolution of 2.5% (5%) at 1 GeV in the barrel (end-cap) region. The time resolution in the TOF barrel region is 68 ps, while that in the end-cap region is 110 ps. The end-cap TOF system was upgraded in 2015 using multi-gap resistive plate chamber technology, providing a time resolution of 60 ps [14][15][16].
The data samples used in this analysis are listed in Table 2 [17]. Since the cross section of D * ± s D ∓ s production in e + e − annihilation is about a factor of twenty larger than that of D + s D − s [18] at C.M. energies √ s = 4.178 − 4.226 GeV, and the D * ± s meson decays to γD ± s with a dominant BF of (93.5 ± 0.7)% [4], the signal events discussed in this paper are selected from the process Simulated data samples produced with a geant4-based [19] Monte Carlo (MC) package, which includes the geometric description of the BESIII detector and the detector response, are used to determine detection efficiencies and to estimate backgrounds. The simulation models the beam energy spread and initial-state radiation (ISR) in the e + e − annihilations with the generator kkmc [20,21]. The inclusive MC sample includes the production of open-charm processes, the ISR production of vector charmonium(-like) states, and the continuum processes incorporated in kkmc. The known decay modes are modelled with evtgen [22,23] using BFs taken from the Particle Data Group (PDG) [4], and the re-maining unknown charmonium decays are modeled with lundcharm [24][25][26][27][28]. Final-state radiation (FSR) from charged final state particles is incorporated using photos [29].

Event selection
The data samples were collected just above the D * ± s D ∓ s threshold. The tag method [30] allows clean signal samples to be selected, providing an opportunity to perform amplitude analyses and to measure the absolute BFs of the hadronic D + s meson decays. In the tag method, a single-tag (ST) candidate requires only one of the D ± s mesons to be reconstructed via a hadronic decay; a double-tag (DT) candidate has both D + s D − s mesons reconstructed via hadronic decays. The DT candidates are required to have the D + s meson decaying to the signal mode D + s → π + π 0 η and the D − s meson decaying to twelve tag modes listed in Table 3. Table 3. Requirements on the tagging D − s mass (M tag ) for various tag modes, where the η and η subscripts denote the decay modes used to reconstruct these particles.

Tag mode
Mass window (GeV/c 2 ) Charged tracks detected in the MDC are required to be within a polar angle (θ) range of |cosθ| < 0.93, where θ is defined with respect to the z-axis which is the symmetry axis of the MDC. For charged tracks not originating from K 0 S decays, the distance of closest approach to the interaction point is required to be less than 10 cm along the beam direction and less than 1 cm in the plane perpendicular to the beam. Particle identification (PID) for charged tracks combines measurements of the dE/dx in the MDC and the flight time in the TOF to form a probability L(h) (h = p, K, π) for each hadron h hypothesis. Charged kaons and pions are identified by comparing the probability for the two hypotheses, L(K) > L(π) and L(π) > L(K), respectively.
The K 0 S candidates are selected by looping over all pairs of tracks with opposite charges, whose distances to the interaction point along the beam direction are within 20 cm. These two tracks are assumed to be pions without PID applied. A primary vertex and a secondary vertex are reconstructed and the decay length between the two vertexes is required to be greater than twice its uncertainty. This requirement is not applied for the D − s → K 0 S K − decay due to the low combinatorial background. Candidate K 0 S particles are required to have the vertex fit and an invariant mass of the π + π − pair (M π + π − ) in the range [0.487, 0.511] GeV/c 2 . To prevent an event being doubly counted in the D − s → K 0 S K − and D − s → K − π + π − selections, the value of M π + π − is required to be outside of the mass range [0.487, 0.511] GeV/c 2 for D − s → K − π + π − decay. Photon candidates are identified using showers in the EMC. The deposited energy of each shower must be more than 25 MeV in the barrel region (|cosθ| < 0.80) and more than 50 MeV in the end cap region (0.86 < |cosθ| < 0.92). The opening angle between the position of each shower in the EMC and the closest extrapolated charged track must be greater than 10 degrees to exclude showers that originate from charged tracks. The difference between the EMC time and the event start time is required to be within [0, 700] ns to suppress electronic noise and showers unrelated to the event.
The π 0 (η) candidates are reconstructed through π 0 → γγ (η → γγ) decays, with at least one photon falling in the barrel region. The invariant mass of the photon pair for π 0 and η candidates must be in the ranges [0.115, 0.150] GeV/c 2 and [0.500, 0.570] GeV/c 2 , respectively, which are about three times larger than the detector resolution. A kinematic fit that constrains the γγ invariant mass to the π 0 or η known mass [4] is performed to improve the mass resolution. The η candidates are also reconstructed through η → π + π − π 0 and the invariant mass of π + π − π 0 are required to satisfy the range of [0.530, 0.560] GeV/c 2 . The ρ 0 candidates are selected via the decay ρ 0 → π + π − with an invariant mass window [0.620, 0.920] GeV/c 2 . The η candidates are formed from the π + π − η and γρ 0 combinations with an invariant mass within a range of [0.946, 0.970] GeV/c 2 .
Twelve tag modes are reconstructed and the corresponding mass windows on the tagging D − s mass (M tag ) are listed in Table 3. The D ± s candidates with M rec lying within the mass windows listed in Table 2 are retained for further study. The quantity M rec is the recoil mass of D ± s and is defined as where E cm is the initial energy of the e + e − C.M. system, p Ds is the three-momentum of the D ± s candidate in the e + e − C.M. frame, and m Ds is the D ± s known mass [4].

Further event selection
The ST D − s mesons are reconstructed using the first eight hadronic decays as shown in Table 3 and the following selection criteria are further applied in order to obtain data samples with high purities for the amplitude analysis. The selection criteria discussed in this section are not used in the BF measurement since the BF measurement is dominated by statistical uncertainty.
After a tag D − s is identified, the signal candidate is selected by requiring one η candidate, one track identified as a charged pion and one π 0 candidate, where π 0 and η candidates are selected by the same requirements described in section 3, but include the decay η → π + π − η; η → γγ only. Then, an nine-constraint (9C) kinematic fit is performed to the process e + e − → D * ± s D ∓ s → γD + s D − s assuming the D − s decays to one of the tag modes and the D + s decays to the signal mode. Two hypotheses are considered: that the signal D + s comes from a D * + s meson or the tag D − s comes from a D * − s meson. The invariant masses of (γγ) π 0 , (γγ) η , (π + π 0 η) η , tag D − s , and D * ± s candidates are constrained to the corresponding known masses [4] and the constraints of four-momentum conservation in the e + e − C.M. system are also applied. The D * ± s D ∓ s combination with the minimum χ 2 9C is chosen. In order to ensure that all candidates fall within the phase-space boundary, the constraint of the signal D + s mass is added to the 9C kinematic fit and the updated four-momenta are used for the amplitude analysis.
To suppress background from fake η candidates, we check the invariant-mass distributions of the γγ combination (M recombined ) which can be with one photon from the signal η (η sig ) and the other photon from the signal π 0 , D * s or the η/π 0 on the tag side. Events with where M η sig is the invariant mass of the photon pair for η sig candidates, while M π 0 and M η are the π 0 and η known masses [4], respectively. GeV. The signal is described by a MC-simulated shape convolved with a Gaussian resolution function, and the background is described by a linear function. Finally, a mass window, [1.92, 2.00] GeV/c 2 , is applied on the signal D + s candidates. A total of 411 events are retained for the amplitude analysis with a purity, w sig , of (96.1 ± 0.9)%. GeV. This plot is obtained using first eight tag modes in Table 3. The black points with error bars are data. The blue line is the total fit. The red and black dashed lines are the fitted signal and background, respectively. The pair of red arrows indicate the selected signal region.

Fit method
The composition of intermediate resonances in the decay D + s → π + π 0 η is determined by an unbinned maximum-likelihood fit to data. The likelihood function is constructed with a probability density function (PDF), which depends on the momenta of the three daughter particles. The amplitude of the n th intermediate state (A n ) is given by where S n and F r(Ds) n are the spin factor and the Blatt-Weisskopf barriers of the intermediate state (the D ± s meson), respectively, and P n is the propagator of the intermediate resonance. The nonresonant amplitudes N L where L denotes the orbital angular momentum between the π + π 0 system with L = 0 (S wave), 1 (P wave), 2 (D wave) are similar to A n in Eq. (4.1) but do not contain resonant propagator terms P n : The total amplitude M is then the coherent sum of the amplitudes of intermediate processes, M = ρ n e iφn A n + ρ L e iφ L N L , where the parameters ρ n and φ n are the magnitudes and phases of the n th resonance, while ρ L and φ L correspond to the magnitudes and phases of the nonresonant contribution with angular momentum L.
The signal PDF f S (p j ) is written as where (p j ) is the detection efficiency parameterized in terms of the final four-momenta p j . The index j refers to the different particles in the final states, and R 3 (p j ) is the standard element of three-body phase space. The normalization integral is determined by a MC integration, where k is the index of the k th event and N MC is the number of the selected MC events. Here M g (p j ) is the PDF used to generate the MC samples in MC integration. To account for any bias caused by differences in tracking and PID efficiencies, and π 0 and η reconstruction efficiencies between data and MC simulation, each MC event is weighted with a ratio, γ (p), between the efficiency of data and MC simulation. Then the MC integral becomes (4.5) A signal-background combined PDF is introduced to account for the background in this analysis. The background PDF is given by The background events in the signal region from the inclusive MC sample are used to model the corresponding background in data. This background description is validated by comparing the M π + η , M π + π 0 and M π 0 η distributions of events outside the M sig signal region between the data and the inclusive MC samples. The distributions of background events from the inclusive MC sample within and outside the M sig signal region are also examined.
They are found to be compatible within statistical uncertainties. The background shape B(p j ) is a probability density function sampled from a multidimensional histogram by using RooHistPdf implemented in RooFit [31]. This background PDF is then added to the signal PDF incoherently and the combined PDF is written as A efficiency-corrected background shape, B (p j ) ≡ B(p j )/ (p j ) is introduced in order to factorize the (p j ) term out from the combined PDF. In this way, the (p j ) term, which is independent of the fitted variables, is regarded as a constant and can be dropped during the log-likelihood fit. As a consequence, the combined PDF becomes Next, the integration in the denominator of the background term can also be handled by the MC integration method in the same way as for the signal only sample: The final log-likelihood function is written as 10) where N D is the number of candidate events in data.

Blatt-Weisskopf barrier factors
For the process a → bc, the Blatt-Weisskopf barrier F L (p j ) [32] is parameterized as a function of the angular momenta L and the momenta q of the daughter b or c in the rest system of a, F L=2 (q) = z 4 0 + 3z 2 0 + 9 z 4 + 3z 2 + 9 ,

Spin factors
The spin-projection operators are defined as [34] P    The quantities p a , p b , and p c are the momenta of particles a, b, and c, respectively, and r a = p b − p c . The covariant tensors are given bỹ The spin factors for S, P , and D wave decays are where theT (l) factors have the same definitions ast (l) . The tensor describing the D + s decay is denoted byT and that of the a decay is denoted byt.

Fit results
The Dalitz plots of M 2 π + π 0 versus M 2 η π + for the data samples and the signal MC samples generated based on the results of the amplitude analysis are shown in Fig. 4 (a) and Fig. 4 (b), respectively. One can see a clear ρ + resonance. Therefore we choose the D + s → ρ + η amplitude as a reference, and fix the magnitude and the phase of its amplitude to 1.0 and 0.0, respectively, while those of other amplitudes are floated. The masses and widths of all resonances are fixed to the corresponding PDG averages [4], and w sig are fixed to the purities discussed in Sec. 4.1. Then we test other possible intermediate resonances, such as ρ(1450), a 0 (1450), π 1 (1600), a 2 (1320), etc., by adding them one by one. We also examined the possible combinations of these intermediate resonances to check their significances, correlations and interferences. We use the difference of log-likelihoods of fits with and without these amplitudes to calculate the significance and find that in all cases these significances are less than three standard deviations. The significance of each intermediate resonance tested is listed in Table 4. Hence the final model consists only of the mode D + s → ρ + η . The mass projections of the fit results are shown in Fig. 5.

Systematic uncertainties for amplitude analysis
The following four sources of potential bias are considered when assigning systematic uncertainties.
i Resonance parameters. The uncertainties related to the fixed parameters in the amplitudes are estimated by varying the masse and width of the ρ + resonance by ±1σ [4].
ii The ρ + lineshape. The uncertainties related to the lineshape of the ρ + are estimated by using a Breit-Wigner function instead of the Gounaris-Sakurai description.
iii The uncertainties associated with γ (p) are obtained by performing alternative amplitude analyses varying PID and tracking efficiencies according to their uncertainties. The systematic uncertainty from this source is found to be negligible.
The assigned systematic uncertainties on the fit fractions (FF) for the S-wave and P -wave nonresonant components are summarized in Table 5. The FF for the L-wave non-resonant amplitude is defined as where N gen is the number of phase-space MC events at generator level. It involves the phase-space MC truth information without detector acceptance or resolution effects. Table 5. Systematic uncertainties on FFs for S-wave and P -wave nonresonant components.

Branching fraction measurement
The ST D − s mesons are reconstructed through all twelve hadronic decays as shown in Table 3 and the selection criteria are the same as those described in Sec. 3 for the branching fraction measurement. In addition, all pions are required to have momenta greater than 100 MeV/c to remove soft pions from D * + decays. The best tag candidate with M rec closest to the D * ± s known mass [4] is chosen if there are multiple ST candidates. The data sets are organized into three sample groups, 4.178 GeV, 4.189-4.219 GeV, and 4.226 GeV, that were acquired during the same year under consistent running condition. The yields for various tag modes are obtained by fitting the corresponding M tag distributions. As an example, the fits to the M tag distributions of the accepted ST candidates from the data sample at √ s = 4.178 GeV are shown in Fig. 7. In the fits, the signal is modeled by a MC-simulated shape convolved with a Gaussian function to account for differences in resolution between data and MC simulation. The background is described by a second-order Chebyshev polynomial. Inclusive MC studies show that there is no peaking background in any tag mode, except for D − → K 0 S π − and D − s → ηπ + π − π − faking the D − s → K 0 S K − and D − s → π − η tags, respectively. Therefore, the MC-simulated shapes of these two peaking background sources are added to the background polynomial functions.
Once a tag mode is reconstructed, we select the signal decay D + s → π + π 0 η . In the case of multiple candidates, the DT candidate with the average mass, (M sig + M tag )/2, closest to the D ± s nominal mass listed in the PDG [4] is retained. To measure the BF, we employ the following equations: where N ST tag is the ST yield for the tag mode; N DT tag,sig is the DT yield; N D + s D − s is the total number of D * ± s D ∓ s pairs produced from the e + e − collisions; B tag and B sig are the BFs of the tag and signal modes, respectively; ST tag is the ST efficiency to reconstruct the tag mode; and DT tag,sig is the DT efficiency to reconstruct both the tag and the signal decay modes. In the case of more than one tag mode and sample group, where α represents the tag mode in the i th sample group. By isolating B sig and replacing N D + s D − s shown in Eq. (5.1) , we find 4) where N ST α,i is obtained from the data sample, while N DT peaking , DT α,sig,i and ST α,i are obtained from the inclusive MC sample. The D + s → π + π 0 η simulated sample is generated according to the results of the amplitude analysis. The BFs B η→γγ , B π 0 →γγ and B η →π + π − η have been introduced to consider these sub-channels.
The N DT fitted is obtained from the fit to the M sig distribution of the selected D + s → π + π 0 η candidates. The fit result is shown in Fig. 8, where the signal shape is described by a MCsimulated shape convolved with a Gaussian function to account for differences in resolution between data and MC. The background shape is described by a MC-simulated shape which excludes peaking background from D + s → π + ηω π + π − π 0 . The number of peaking background events, N DT peaking is estimated from the inclusive MC sample. Thus, N DT fitted and N DT peaking are determined to be 837±35 and 5±1, respectively. Tables 6 -8 summarize the ST efficiencies, DT efficiencies, and ST yields in data samples at the C.M. energies √ s = 4.178−4.226 GeV. Taking into account the differences in π ± tracking/PID efficiencies, π 0 and η reconstruction efficiencies between data and MC simulation, we determine the BF B(D + s → π + π 0 η ) = (6.15 ± 0.25(stat.) ± 0.18(syst.))% according to Eq. (5.4).  Table 3. The data are represented by points with error bars, the total fit by the blue solid line, and the fitted signal and the fitted background by the red dotted and the black dashed lines, respectively. The pair of red arrows indicate the signal region.
The following sources of the systematic uncertainties are taken into account for the BF measurement.
• The number of ST D + s mesons. The systematic uncertainty due to the total yield of the ST D − s mesons is assigned to be 0.9% by taking into account the background fluctuation in the fit, and examining the changes of the fit yields by using alternative signal and background shapes.
• Background shape. The systematic uncertainty due to the MC-simulated background shape is studied by varying the relative fractions of the background from qq or non-D * + s D − s open charm by the statistical uncertainties of their related cross sections. It is found that the uncertainty arising from this source is 0.1% which is small enough to be neglected.
• π 0 , η reconstruction. The systematic uncertainty associated with the π 0 reconstruction efficiency is investigated by using a control sample of the process e + e − → K + K − π + π − π 0 . The same selection criteria described in Sec. 3 are used to reconstruct the two kaons and the two pions. The recoiling mass distribution of K + K − π + π − is fitted to obtain the total number of π 0 s and the π 0 selection is applied to determine the number of reconstructed π 0 mesons. The average ratio between data and MC efficiencies of π 0 reconstruction, weighted by the corresponding momentum spectra, is estimated to be 1.006 ± 0.009. Similarly, the average ratio between data and MC efficiencies of η reconstruction is estimated to be 1.011 ± 0.010. After correcting the efficiencies, the systematic uncertainties associated with reconstruction efficiencies are 0.9% for π 0 and 1.0% for η mesons.
• MC sample size. The uncertainty arising from the finite MC sample size is obtained by α (f α δ α α ) 2 , where f α is the tag-yield fraction, and α and δ α are the signal efficiency and the corresponding uncertainty of tag mode α, respectively.
• Amplitude model. The uncertainty from the amplitude model is estimated by varying the amplitude-model parameters. For the mass and width of ρ + resonance, we sample them with a Gaussian distribution in which the mean and width are set to the corresponding known value and uncertainty from PDG [4]. Meanwhile, we uniformly vary the effective radii of Blatt-Weisskopf Barrier within the range • Peaking background. The uncertainties caused by peaking background is studied by varying the BF of D + s → π + ηω π + π − π 0 from 0.85% to 1.39% based on the precision of the measured branching ratio [4]. The shift in DT yield is 0.2%, which is taken as the corresponding uncertainty.
All of the systematic uncertainties are summarized in Table 9. Adding them in quadrature gives a total systematic uncertainty in the BF measurement of 2.9%.

Summary
This paper presents the amplitude analysis of the decay D + s → π + π 0 η with 6.32 fb −1 of e + e − collision data samples at √ s = 4.178 − 4.226 GeV. The mode D + s → ρη is found to be the main intermediate process contributing to this final state. In addition, we also report the upper limits of the BFs of S−wave and P −wave nonresonant components of D + s → π + π 0 η to be B(D + s → (π + π 0 ) S η ) < 0.10% and B(D + s → (π + π 0 ) P η ) < 0.74% at the 90% confidence level, respectively.