Evidence for $\eta_c$(2S)$\to\pipieta$ decay

The decay $\eta_c(2S)\to\pipieta$ is searched for through the radiative transition $\psi(3686) \to\gamma\eta_c(2S)$ using 448 million $\psi$(3686) events accumulated at the BESIII detector. The first evidence of $\eta_c(2S)\to\pi^+\pi^-\eta$ is found with a statistical significance of 3.5$\sigma$. The product of the branching fractions of $\psi(3686)\to\gamma\eta_c(2S)$ and $\eta_c(2S)\to\pipieta$ is measured to be $Br(\psi(3686)\to\gamma\eta_c(2S))\times Br(\eta_c(2S)\to\pipieta)=(2.97\pm0.81\pm0.26)\times10^{-6}$, where the first uncertainty is statistical and the second one is systematic. The branching fraction of the decay $\eta_c(2S)\to\pipieta$ is determined to be $Br(\eta_c(2S)\to\pipieta)=(42.4\pm11.6\pm3.8\pm30.3)\times10^{-4}$, where the third uncertainty is transferred from the uncertainty of the branching fraction of $\psi(3686)\to\gamma\eta_c(2S)$.


I. INTRODUCTION
Charmonium states play an important role in understanding the strong interaction, since their masses reside on the boundary between the perturbative and nonperturbative energy regions.The states below the opencharm threshold are better understood than the ones above, and their mass spectrum can be described by the quark potential model [1,2].Our knowledge about the spin singlets, however, including the P -wave state h c , the S-wave ground state η c and its first radial excitation η c (2S), is still limited [3].The η c (2S) was first observed by the Belle collaboration via the decay [4], more than two decades after its prediction [5].It was confirmed later by the CLEO [6] and BaBar [7] collaborations in the two-photon fusion process γγ → η c (2S) → K 0 S K ± π ∓ and in addition by the BaBar collaboration in the double charmonium production process e + e − → J/ψcc [8].The production of η c (2S) is also expected through the magnetic dipole (M1) transition of ψ(3686) [9].The radiative transition ψ(3686)→ γη c (2S), with η c (2S)→ K + K − π 0 and K 0 S K ± π ∓ , was reported by the BESIII collaboration in 2012 [10].
The decay of charmonium states into light hadrons is believed to be dominated by the annihilation of the cc pair into two or three gluons.The so-called "12% rule" states that the ratio of the inclusive branching fractions of light hadron final states between ψ(3686) and J/ψ is about 12% [11].Violations of this rule have been observed in various decay channels, especially in the ψ → ρπ process [12].That so-called "ρ − π puzzle" has not been solved yet.Similarly, one would expect a similar ratio of the hadronic branching fractions between η c (2S) and η c due to their analogous decaying dynamics in comparison to ψ(3686) and J/ψ.According to ref. [13], for any normal light hadronic channel h, while ref.[14] argues that this ratio should be close to one if no mixing with glueballs is considered.Recently, it was found that the experimental data significantly differ from both theoretical predictions [15], i.e., most of the ratios are obviously greater than 12% and less than one, except the ones with pp final states.Up to now, the total measured branching fraction of η c (2S) decays is small, i.e., less than 5% according to the report from the particle data group (PDG), and the uncertainties of all the available experimental measurements are greater than 50% [16].Significantly improving the measurement precision of any single channel is difficult because of the still limited statistics of η c (2S) samples and, e.g., the difficulty in tagging the very soft radiative photon from the ψ(3686) transition.Searching for more decay modes of the η c (2S) is therefore desirable to reduce the uncertainty of the averaged value of these ratios.

II. BESIII EXPERIMENT AND MONTE CARLO SIMULATION
The BESIII detector [18] records symmetric e + e − collisions provided by the BEPCII storage ring [19], which operates in the center-of-mass energy range from 2.0 to 4.95 GeV, with a peak luminosity of 1 × 10 33 cm −2 s −1 achieved at √ s = 3.77 GeV.BESIII has collected large data samples in this energy region [20].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 solenoid 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 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.
Simulation samples produced with a geant4based [21] Monte Carlo (MC) package, which includes the geometric description of the BESIII detector and the detector response, are utilized to determine reconstruction efficiencies and to estimate background contributions.The simulation models the beam-energy spread and initial state radiation (ISR) in the e + e − annihilation with the generator kkmc [22].The inclusive MC sample includes the production of the ψ(3686) resonance, the ISR production of the J/ψ, and the continuum processes incorporated in kkmc [22].The decay modes are modelled with evtgen [23] using the known averaged branching fractions [16], and the unknown charmonium decays are modelled with lundcharm [24].Final state radiation (FSR) from charged final state particles is incorporated using the photos package [25].Event type analysis of the inclusive MC samples with a generic tool, TopoAna [26], is used to study the potential background.The exclusive decays of ψ(3686) → γX(X = χ c1,2 and η c (2S)) are generated taking into the angular distribution, while the decay X → π + π − η is generated uniformly in phase space (PHSP).

III. EVENT SELECTION AND BACKGROUND ANALYSIS
The η c (2S) candidates studied in this analysis are produced via the M1 transition of ψ(3686) → γη c (2S)(η c (2S) → π + π − η), the η is reconstructed via its two-photon decay (η → γγ).Therefore, the final state is γγγπ + π − , i.e., the candidate events are required to have two charged tracks with a net charge of zero, and at least three photons.
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, the distance of the closest approach to the interaction point must be less than 10 cm along the z-axis, and less than 1 cm in the transverse plane.Particle identification (PID) for charged tracks combines measurements of the energy deposited in the MDC (dE/dx) and the flight time in the TOF to form likelihoods L(h) (h = p, K, π) for each hadron h hypothesis.Tracks are identified as pions when the pion hypothesis has the greatest likelihood (L(π) > L(K) and L(π) > L(p)).Each event is required to have one π + and one π − .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 40 MeV in the end cap region (0.86 < | cos θ| < 0.92).To exclude showers that originate from charged tracks, the angle subtended by the EMC shower and the position of the closest charged track at the EMC must be greater than 10 degrees as measured from the interaction point.To suppress electronic noise and showers unrelated to the event, the difference between the EMC time and the event start time is required to be within [0, 700] ns.
A kinematic fit with five constraints (5C) on each candidate is performed, where the total energy-momentum of final states is constrained to the initial four-momentum and the invariant mass of the two photons is constrained to the nominal η mass.The fit loops over all possible η → γγ and M1 γ combinations, and the combination with the minimal χ 2 5C (3γ) is selected.In order to suppress the background channels ψ(3686)→ π + π − η and ψ(3686)→ γγπ + π − η, with one less or one more photon in the final state than the signal, the value of χ 2 5C (3γ) is required to be less than χ 2 5C (2γ) and χ 2 5C (4γ).Also, a selection is set on the χ 2 5C (3γ) distribution.As can be seen in the invariant mass spectrum of π + π − η in Fig. 1, the background contribution of the channel ψ(3686) → π + π − η with a fake photon would appear as a peak close to the η c (2S) signal.This peak can be reduced, if the energy of the candidate M1 photon is not used in the kinematic fit, by adjusting the 5C fit to a modified 4C fit (m4C) [10].By this method, this background channel is reduced in the η c (2S) mass region as shown in Fig. 1.In this analysis, the invariant mass distribution after applying the m4C fit is used for further study.
The signal suffers significantly from background contributions associated with J/ψ decaying to µ + µ − (γ FSR ), such as ψ(3686) → ηJ/ψ and π 0 π 0 J/ψ.As shown in Fig. 2, there is a large enhancement around the J/ψ resonance in the invariant mass of γ(M1)π + π − .Events with M γ(M1)π + π − ≤ 3.00 GeV/c 2 are accepted for further  The decay of ψ(3686) → π + π − ηγ FSR with an FSR photon cannot be reduced because of the same final state as the signal decay mode.It could potentially contaminate the signal channel with a long tail in the invariant mass of π + π − η.The contribution of this channel depends on the FSR ratio (R FSR ), which is defined as R FSR = NFSR Nno FSR ignoring any kinematic dependence, where N FSR and N no FSR are the numbers of events with and without an FSR photon [27].A control sample of ψ(3686) → γχ c0 → γγ FSR 2(π + π − ) is selected to study the difference of R FSR between data and MC sample.From the analysis, we find f FSR = R Data FSR /R MC FSR = 1.70 ± 0.07, with a purely statistical uncertainty.The line-shape of ψ(3686) → (γ FSR )π + π − η is described by the sum of MC simulated shapes of ψ(3686) → π + π − η and ψ(3686) → γ FSR π + π − η with the FSR ratio corrected by the factor f FSR .It is then used to describe the corresponding background contribution later in the fit.
In addition, background channels such as ψ(3686) → γχ c1 (π + π − π 0 ) and ψ(3686)→ ω(γπ 0 )π + π − with a π 0 in the final state contribute only insignificantly.Furthermore, the line shapes of these two decay modes are smooth in the invariant mass of π + π − η according to the study of MC samples, so no π 0 veto is applied to ensure the significance of the signal.

IV. BRANCHING FRACTION DETERMINATION
The branching fraction is calculated by (2) where N sig is the observed number of signal events, N ψ(3686) is the total number of ψ(3686) decays [17], is the detection efficiency, BR 1 and BR 2 are the branching fractions of ψ(3686) → γη c (2S) and η → γγ [16], respectively.These values are listed in Table I.
Nsig N ψ(3686) (×10 6 ) BR1(×10 −4 ) BR2(×10 −2 ) 106 ± 29 448.1 ± 2.9 0.202 7 ± 5 39.41 ± 0.20 The signal yield N sig is extracted by an unbinned maximum likelihood fit to the invariant mass distribution of π + π − η (see Fig. 3).The fit range is from 3.35 to 3.70 GeV/c 2 , which includes the χ c1 and χ c2 signals.The line shapes of η c (2S) and χ c1,2 are described by the following formulas, and respectively, Here, m is the invariant mass of π + π − η, is the energy of the transition photon in the rest frame of ψ(3686), BW 0 , BW 1 and BW 2 denote the Breit-Wigner functions for η c (2S), χ c1 and χ c2 with the mass and width fixed to the reported averaged values of these three resonances [16], DGaus(δm, σ) is a double-Gaussian function describing the mass shift and resolution whose parameters are shared by χ cJ and η c (2S), f d (E γ ) is a damping function to suppress the di- the mean energy of the transition photon chosen for this work.The contributions of χ c1,2 backgrounds are determined by the fit, and the results are consistent with the world averaged values [16].The remaining background due to the FSR process ψ(3686) → (γ FSR )π + π − η is described by the MC shape.All smooth background contributions are described by an Argus function with the threshold fixed and the other parameters floated.
The fit result is shown in Fig. 3 and the goodness of fit is χ 2 /ndf = 67/59 = 1.1, where ndf denotes the number of degrees of freedom.The yield of the η c (2S) signal is N sig = 106 ± 29 with a statistical significance of 3.5σ that is obtained from the difference of the logarithmic likelihoods [29], taking into account the difference of ndf .We obtain a product branching fraction of Br(ψ(3686 FIG.3: The result of a fit to the invariant mass distribution of π + π − η.The black dots with error bars are BESIII data, the red and black solid curves denote the total fit curve and the shape of the smooth background contributions, respectively.The green dash-dot-dotted curve, the blue dotted, and the cyan dash-dotted denote the decay modes of χc1, χc2 and ηc(2S), while the pink dashed curve denotes the contributions of the FSR process ψ(3686) → (γFSR)π + π − η.

V. SYSTEMATIC UNCERTAINTIES
The sources of systematic uncertainties for the branching fraction measurement include the number of ψ(3686) events, tracking, PID, photon detection, the kinematic fit, the branching fractions of the intermediate decays, and the fit to the invariant mass distribution of π + π − η.
• Tracking and PID efficiency.The pion tracking efficiency was studied by using a control sample of J/ψ → ppπ + π − .The difference in the tracking efficiencies between MC simulation and data is 1% per pion [30].The pion PID efficiency was investigated with a control sample of J/ψ → ρπ [31].The difference in the PID efficiencies between data and MC simulation is 1% per pion.
• Photon detection.The uncertainty due to photon reconstruction is determined to be 1% per photon by the study of the process J/ψ → π + π − π 0 [32].
• Kinematic fit.A control sample of ψ(3686) → γχ c1 , χ c1 → π + π − η is selected to estimate the uncertainty associated with the kinematic fit.The difference of the efficiencies with and without the kinematic fit is taken as the corresponding uncertainty, determined to be 1.6%.
Seven potential sources of systematic uncertainties are considered.First, the uncertainty from the background shape is estimated by using an alternative Chebychev function.The change of the fitted signal yield, 3.8%, is assigned as the uncertainty.Second, the uncertainty of the signal shape is estimated by changing the width of η c (2S) with ±1 standard deviation away from the nominal value, and is determined to be 4.7%.Third, the uncertainty caused by the shapes of χ c1 and χ c2 is estimated by changing the widths of χ c1 and χ c2 with ±1 standard deviation away from their nominal values, it is determined to be 1.9%.Fourth, the damping function is changed to an alternative form used by CLEO [33], ) with β = (65.0±2.5)MeV, the resulting difference in the fit is 3.8%, and is assigned as the systematic uncertainty.Fifth, the systematic uncertainty from the double-Gaussian function is estimated by changing it to a triple-Gaussian function, and the resulting difference is 0.9%.Sixth, the systematic uncertainty associated with the f FSR is estimated by varying the f FSR with ±1 standard deviation that is mainly caused by the statistics of the control samples, and the difference is 0.9%.Seventh, the fit range is varied, and the maximum differences in the fitted yields are considered as the associated systematic uncertainties.It is determined to be 1.9%.Among all sources of systematic uncertainties, by far the largest one comes from the branching fraction of ψ(3686) → γη c (2S), and is therefore treated separately.All other sources of systematic uncertainties are assumed to be independent of each other and combined in quadrature to obtain the overall systematic uncertainty as listed in Table II.With the (448.1±2.9)×10 6ψ(3686) data sample, the process of η c (2S) → π + π − η is searched for by utilizing the M1 transition ψ(3686) → γη c (2S).Evidence for the decay η c (2S) → π + π − η is found for the first time, with a statistical significance of 3.5σ.The product of the branching fractions is measured to be Br(ψ(3686 With the branching fraction Br(η c → π + π − η) =(1.7±0.5)%[16], the ratio of the branching fractions of η c and η c (2S) decaying into π + π − η is calculated to be Br(ηc(2S)→π + π − η) Br(ηc→π + π − η) = 0.25 ± 0.20.Combining the ratios of other hadronic decay modes of η c (2S) to η c [16,34], the averaged value of all these ratios including this measurement is determined to be 0.30±0.10(see Fig. 4).This ratio agrees neither with the prediction in Ref. [13] nor in Ref. [14].The observed discrepancy reflects our limited knowledge of the decay mechanisms of the spin singlet charmonium states.More searches on new decay modes and more precise measurements of the η c (2S) decays are required to shed light on this puzzle.With about 3 billion ψ(3686) events to be accumulated, BESIII will make a further substantial contribution to this field [35].A better precision of the branching fraction of ψ(3686) → γη c (2S) is hereby crucial.
FIG. 4: Estimation of the averaged value of the ratio of Br(ηc(2S) → h) to Br(ηc → h).Here, h means various hadronic final states, as shown on this figure.Except the branching fraction of ηc(2S) → π + π − η of this work, the other results are quoted from Refs.[16,34].The shade is the averaged value of these five decay modes with one standard deviation.Both statistical and systematic uncertainties have been included.

FIG. 1 :
FIG.1:The comparison of M π + π − η obtained from the m4C and 5C kinematic fits of MC samples (Color Online).The red solid curve and the red dashed curve denote the m4C and 5C kinematic results of the signal channel, respectively; the black dotted curve and the black dash-dotted curve denote the m4C and 5C kinematic results of ψ(3686)→ π + π − η, respectively.

FIG. 2 :
FIG. 2: The distributions of M γ M1 π + π − of data and signal MC samples.The black dots with error bars denote the data sample, the red histogram denotes the signal MC sample, and the blue vertical line indicates the corresponding selection criterion.

e
Also at Key Laboratory for Particle Physics, Astrophysics and Cosmology, Ministry of Education; Shanghai Key Laboratory for Particle Physics and Cosmology; Institute of Nuclear and Particle Physics, Shanghai 200240, People's Republic of China f Also at Key Laboratory of Nuclear Physics and Ion-beam Application (MOE) and Institute of Modern Physics, Fudan University, Shanghai 200443, People's Republic of China g Also at State Key Laboratory of Nuclear Physics and Technology, Peking University, Beijing 100871, People's Republic of China h Also at School of Physics and Electronics, Hunan University, Changsha 410082, China i Also at Guangdong Provincial Key Laboratory of Nuclear Science, Institute of Quantum Matter, South China Normal University, Guangzhou 510006, China j Also at Frontiers Science Center for Rare Isotopes, Lanzhou University, Lanzhou 730000, People's Republic of China k Also at Lanzhou Center for Theoretical Physics, Lanzhou University, Lanzhou 730000, People's Republic of China l Also at the Department of Mathematical Sciences, IBA, Karachi , Pakistan

TABLE I :
Input values for calculating the branching fraction Br

TABLE II :
Summary of systematic uncertainties.