Measurements of 
e+e−→KS0K±π∓π0
 and 
KS0K±π∓η
 at center-of-mass energies from 3.90 to 4.60 GeV

M. Ablikim, M. N. Achasov, S. Ahmed, M. Albrecht, M. Alekseev, A. Amoroso, F. F. An, Q. An, Y. Bai, O. Bakina, R. Baldini Ferroli, Y. Ban, K. Begzsuren, D.W. Bennett, J. V. Bennett, N. Berger, M. Bertani, D. Bettoni, F. Bianchi, I. Boyko, R. A. Briere, H. Cai, X. Cai, A. Calcaterra, G. F. Cao, S. A. Cetin, J. Chai, J. F. Chang, W. L. Chang, G. Chelkov, G. Chen, H. S. Chen, J. C. Chen, M. L. Chen, P. L. Chen, S. J. Chen, Y. B. Chen, W. Cheng, G. Cibinetto, F. Cossio, H. L. Dai, J. P. Dai, A. Dbeyssi, D. Dedovich, Z. Y. Deng, A. Denig, I. Denysenko, M. Destefanis, F. De Mori, Y. Ding, C. Dong, J. Dong, L. Y. Dong, M. Y. Dong, Z. L. Dou, S. X. Du, P. F. Duan, J. Z. Fan, J. Fang, S. S. Fang, Y. Fang, R. Farinelli, L. Fava, S. Fegan, F. Feldbauer, G. Felici, C. Q. Feng, M. Fritsch, C. D. Fu, Y. Fu, Q. Gao, X. L. Gao, Y. Gao, Y. G. Gao, Z. Gao, B. Garillon, I. Garzia, A. Gilman, K. Goetzen, L. Gong, W. X. Gong, W. Gradl, M. Greco, L. M. Gu, M. H. Gu, Y. T. Gu, A. Q. Guo, L. B. Guo, R. P. Guo, Y. P. Guo, A. Guskov, Z. Haddadi, S. Han, X. Q. Hao, F. A. Harris, K. L. He, F. H. Heinsius, T. Held, Y. K. Heng, Z. L. Hou, H. M. Hu, J. F. Hu, T. Hu, Y. Hu, G. S. Huang, J. S. Huang, X. T. Huang, X. Z. Huang, Z. L. Huang, T. Hussain, N. Hsken, W. Ikegami Andersson, M. Irshad, Q. Ji, Q. P. Ji, X. B. Ji, X. L. Ji, X. S. Jiang, X. Y. Jiang, J. B. Jiao, Z. Jiao, D. P. Jin, S. Jin, Y. Jin, T. Johansson, A. Julin, N. Kalantar-Nayestanaki, X. S. Kang, M. Kavatsyuk, B. C. Ke, I. K. Keshk, T. Khan, A. Khoukaz, P. Kiese, R. Kiuchi, R. Kliemt, L. Koch, O. B. Kolcu, B. Kopf, M. Kornicer, M. Kuemmel, M. Kuessner, A. Kupsc, M. Kurth, W. Kühn, J. S. Lange, P. Larin, L. Lavezzi, S. Leiber, H. Leithoff, C. Li, Cheng Li, D. M. Li, F. Li, F. Y. Li, G. Li, H. B. Li, H. J. Li, J. C. Li, J. W. Li, Ke Li, Lei Li, P. L. Li, P. R. Li, Q. Y. Li, T. Li, W. D. Li, W. G. Li, X. L. Li, X. N. Li, X. Q. Li, Z. B. Li, H. Liang, Y. F. Liang, Y. T. Liang, G. R. Liao, L. Z. Liao, J. Libby, C. X. Lin, D. X. Lin, B. Liu, B. J. Liu, C. X. Liu, D. Liu, D. Y. Liu, F. H. Liu, Fang Liu, Feng Liu, H. B. Liu, H. L. Liu, H. M. Liu, Huanhuan Liu, Huihui Liu, J. B. Liu, J. Y. Liu, K. Y. Liu, Ke Liu, Q. Liu, S. B. Liu, X. Liu, Y. B. Liu, Z. A. Liu, Zhiqing Liu, Y. F. Long, X. C. Lou, H. J. Lu, J. D. Lu, J. G. Lu, Y. Lu, Y. P. Lu, C. L. Luo, M. X. Luo, T. Luo, X. L. Luo, S. Lusso, X. R. Lyu, F. C. Ma, H. L. Ma, L. L. Ma, M.M. Ma, Q. M. Ma, X. N. Ma, X. X. Ma, X. Y. Ma, Y. M. Ma, F. E. Maas, M. Maggiora, S. Maldaner, Q. A. Malik, A. Mangoni, Y. J. Mao, Z. P. Mao, S. Marcello, Z. X. Meng, J. G. Messchendorp, G. Mezzadri, J. Min, T. J. Min, R. E. Mitchell, X. H. Mo, Y. J. Mo, C. Morales Morales, N. Yu. Muchnoi, H. Muramatsu, A. Mustafa, S. Nakhoul, Y. Nefedov, F. Nerling, I. B. Nikolaev, Z. Ning, S. Nisar, S. L. Niu, S. L. Olsen, Q. Ouyang, S. Pacetti, Y. Pan, M. Papenbrock, P. Patteri, M. Pelizaeus, J. Pellegrino, H. P. Peng, K. Peters, J. Pettersson, J. L. Ping, R. G. Ping, A. Pitka, R. Poling, V. Prasad, H. R. Qi, M. Qi, T. Y. Qi, S. Qian, C. F. Qiao, N. Qin, X. S. Qin, Z. H. Qin, J. F. Qiu, S. Q. Qu, K. H. Rashid, C. F. Redmer, M. Richter, M. Ripka, A. Rivetti, M. Rolo, G. Rong, Ch. Rosner, M. Rump, A. Sarantsev, M. Savrié, K. Schoenning, W. Shan, X. Y. Shan, M. Shao, C. P. Shen, P. X. Shen, X. Y. Shen, H. Y. Sheng, X. Shi, J. J. Song, W.M. Song, X. Y. Song, S. Sosio, C. Sowa, S. Spataro, G. X. Sun, J. F. Sun, L. Sun, S. S. Sun, X. H. Sun, Y. J. Sun, Y. K. Sun, Y. Z. Sun, Z. J. Sun, Z. T. Sun, Y. T. Tan, C. J. Tang, G. Y. Tang, X. Tang, M. Tiemens, B. Tsednee, I. Uman, B. Wang, B. L. Wang, C.W. Wang, D. Y. Wang, Dan Wang, K. Wang, L. L. Wang, L. S. Wang, M. Wang, Meng Wang, P. Wang, P. L. Wang, W. P. Wang, X. F. Wang, Y. Wang, Y. F. Wang, Z. Wang, Z. G. Wang, Z. Y. Wang, Zongyuan Wang, T. Weber, D. H. Wei, P. Weidenkaff, S. P. Wen, U. Wiedner, M. Wolke, L. H. Wu, L. J. Wu, Z. Wu, L. Xia, X. Xia, Y. Xia, D. Xiao, Y. J. Xiao, Z. J. Xiao, Y. G. Xie, Y. H. Xie, X. A. Xiong, Q. L. Xiu, G. F. Xu, J. J. Xu, L. Xu, Q. J. Xu, X. P. Xu, F. Yan, L. Yan, W. B. Yan, W. C. Yan, Y. H. Yan, H. J. Yang, H. X. Yang, L. Yang, R. X. Yang, S. L. Yang, Y. H. Yang, Y. X. Yang, Yifan Yang, Z. Q. Yang, M. Ye, M. H. Ye, J. H. Yin, Z. Y. You, B. X. Yu, C. X. Yu, J. S. Yu, C. Z. Yuan, Y. Yuan, A. Yuncu, A. A. Zafar, Y. Zeng, B. X. Zhang, B. Y. Zhang, C. C. Zhang, D. H. Zhang, H. H. Zhang, H. Y. Zhang, J. Zhang, J. L. Zhang, J. Q. Zhang, J. W. Zhang, J. Y. Zhang, J. Z. Zhang, K. Zhang, L. Zhang, S. F. Zhang, T. J. Zhang, X. Y. Zhang, Y. Zhang, Y. H. Zhang, Y. T. Zhang, Yang Zhang, Yao Zhang, Yu Zhang, Z. H. Zhang, Z. P. Zhang, Z. Y. Zhang, G. Zhao, J. W. Zhao, J. Y. Zhao, J. Z. Zhao, Lei Zhao, Ling Zhao, M. G. Zhao, Q. Zhao, S. J. Zhao, T. C. Zhao, Y. B. Zhao, Z. G. Zhao, A. Zhemchugov, B. Zheng, J. P. Zheng, W. J. Zheng, Y. H. Zheng, B. Zhong, L. Zhou, Q. Zhou, X. Zhou, X. K. Zhou, X. R. Zhou, PHYSICAL REVIEW D 99, 012003 (2019)


I. INTRODUCTION
With the experimental progress in the past decade, many charmoniumlike ðXYZÞ states were observed, which can not be accommodated within the naive quark model and are proposed as the candidate of the hidden-charm exotic mesons [1,2]. In this paper, we focus on the exotic states Yð4260Þ and Z c ð3900Þ.
The Z c ð3900Þ was observed in the J=ψπ AE invariant mass distribution of the e þ e − → π þ π − J=ψ process by the BESIII Collaboration [44]. Subsequently, additional Z c ð3900Þ decay channels were observed in DD Ã þ c:c: [45][46][47]. The J P of Z c ð3900Þ was determined to be 1 þ with a partial wave analysis of the π þ π − J=ψ final state [48]. Since its discovery, many interpretations on the nature of the Z c ð3900Þ have been proposed, such as a DD Ã molecule [49], a tetraquark state [50], a cusp effect [51], and dynamical generation through threshold effects [52,53]. Despite many interpretations, the natures of the Yð4260Þ and Z c ð3900Þ are still unclear. To comprehend these states, it is necessary to study more decay modes. All of the observed decay modes of the Yð4260Þ and Z c ð3900Þ are associated with the charm sector and no light hadron decay modes have been found yet [54][55][56]. A search for light hadron decay modes of the Yð4260Þ and Z c ð3900Þ is complementary to previous studies and may help to distinguish between different theoretical models and to understand strong interaction effects in this energy region. Among the large number of potential light hadron decay modes, the branching fractions (BFs) of charmonium states decaying into K S K AE π ∓ π 0 and K S K AE π ∓ η are usually large [57]. In four-body final state, there should be abundant intermediate states, which may supply more possible decay channels for searching Yð4260Þ and Z c ð3900Þ. Furthermore, the existence of charged and neutral pions in the final states enables a study of isospin multiplets. In this paper, we present a measurement of the Born cross sections (σ B ) for e þ e − → K 0 S K AE π AE π 0 and K 0 S K AE π ∓ η. We also report upper limits of e þ e − → Yð4260Þ

II. DETECTORS AND DATA SAMPLES
The BESIII detector [58], operating at the BEPCII collider [59], is a general purpose spectrometer with a geometrical acceptance of 93% of 4π solid angle. It has four main components: (1) a small-cell, helium-based (60% He, 40% C 3 H 8 ) multilayer drift chamber with 43 layers providing an average single-hit resolution of 135 μm, a momentum resolution of 0.5% at 1.0 GeV=c in a 1.0 T magnetic field, and a specific ionization energy loss (dE=dx) resolution better than 6%, (2) a time-of-flight (TOF) detector constructed of 5 cm thick plastic scintillators, with 176 strips of 2.4 m length in two layers in the barrel and 96 fans of the end caps with time resolutions of 80 and 110 ps, respectively, which provide a 2σ K=π separation for momenta up to ∼1.0 GeV=c, (3) an electromagnetic calorimeter (EMC) consisting of 6240 CsI crystals in a cylindrical barrel structure and two end caps with an energy resolution of 2.5% (5%) at 1.0 GeV and a position resolution of 6 mm (9 mm) in the barrel (end caps), and (4) a muon counter consisting of resistive plate chambers in nine barrel and eight end-cap layers, which provide a 2 cm position resolution. More details of the BESIII detector can be found in Ref. [58].
This analysis is based on 5.2 fb −1 e þ e − annihilation data samples [60] collected with the BESIII detector at centerof-mass energies ( ffiffi ffi s p ) from 3.90 to 4.60 GeV [61], which are listed in Table I. Monte Carlo (MC) simulations are used to optimize the event selection criteria, to study the detector acceptance and to understand the potential backgrounds. The GEANT4-based [62] simulation software BOOST [63] is implemented to simulate the detector response, describe geometry and material, realize digitization, and incorporate time-dependent beam backgrounds. Six generic MC samples, equivalent to the integrated luminosity of the data at the energy points 4.009, 4.230, 4.260, 4.360, 4.420 and 4.600 GeV are generated to study the backgrounds. The primary known decay channels are TABLE I. Data sets and results of the Born cross section measurement for e þ e − → K 0 S K AE π ∓ π 0 . The table includes the integrated luminosity L, the number of observed signals events N sig , the total efficiency ϵ, the ISR correction factor ð1 þ δ ISR Þ, the vacuum polarization correction factor 1 j1−ΠðsÞj 2 , and the Born cross section σ B . The first errors are statistical and the second ones are systematic. The details of systematic uncertainties are described in Sec. III D.
ffiffi ffi  [64] with the BFs set to the world average values [57] while the unknown decay modes are generated with LUNDCHARM [65]. Continuum hadronic events are generated with KKMC [66] and QED processes such as Bhabha scattering, dimuon, and digamma events are generated with KKMC and BABAYAGA [67]. To study the efficiency of each final state, a sample of 1 × 10 5 signal events is generated at each energy point using KKMC, which simulates e þ e − annihilation, including beam energy spread and ISR effects.
and π 0 =η → γγ are selected according to the following steps. First, K 0 S candidates are selected by looping over all pairs of oppositely charged tracks, which are assumed to be pions. Next, primary and secondary vertex fits [68] are performed and the decay length of the secondary vertex fit is required to be greater than twice its uncertainty. Furthermore, the invariant mass of the pion pair is required denotes the nominal mass of the K 0 S [57]. If there are multiple K 0 S candidates in one event, the one with the smallest χ 2 from the secondary vertex fit is selected.
In addition to the two charged tracks that make up the K 0 S , two oppositely charged tracks are required. For the latter two charged tracks, the polar angle θ must satisfy j cos θj < 0.93 and the distance of closest approach to the interaction point must be less than 10.0 and 1.0 cm along the beam direction and in the plane perpendicular to the beam direction, respectively. The particle type for each charged track is determined by selecting the hypotheses with the highest probability, which is calculated with the combined information from TOF and dE=dx measurements for different particle hypotheses. One charged track must be identified as a kaon and the other as a pion.
Photons are reconstructed from clusters deposited in the EMC, with the energy measured in the TOF included to improve reconstruction efficiency and energy resolution. At least two photons are required per event. The energy of a photon candidate is required to be larger than 25 MeV in the barrel region (j cos θj < 0.80) or 50 MeV in the end-cap region (0.86 < j cos θj < 0.92). The cluster timing is required to be between 0 and 700 ns to suppress electronic noise and energy depositions unrelated to the event of interest. To eliminate showers associated with charged particles, the opening angle between a photon candidate and the extrapolated position of the closest charged track should be larger than 20 degrees.
Finally, a four-constraint (4C) kinematic fit imposing energy-momentum conservation is performed to the final states. Only events with χ 2 4C < 60 are accepted. For events with more than two photon candidates, the photon pair with the smallest χ 2 4C from the kinematic fit is accepted. After the 4C kinematic fit, no peaking background is observed in the generic MC samples. The invariant mass distributions of π þ π − versus γγ and the projections onto Mðπ þ π − Þ and MðγγÞ after 4C fit are shown in Fig. 1 (at 4.258 GeV as an example), in which obvious K 0 S and π 0 =η peaks are observed. The signal regions are defined as Mðπ þ π − Þ ∈ ð0.488; 0.508Þ GeV=c 2 , MðγγÞ ∈ ð0.12;0.15Þ GeV=c 2 (for the π 0 mode) and MðγγÞ ∈ ð0.52; 0.58Þ GeV=c 2 (for the η mode). The sideband regions are defined as Mðπ þ π − Þ ∈ ð0.463; 0.483Þ ∪ ð0.513; 0.533Þ GeV=c 2 , MðγγÞ ∈ ð0.08; 0.11Þ ∪ ð0.16; 0.19Þ GeV=c 2 (for the π 0 mode) and MðγγÞ ∈ ð0.44; 0.50Þ ∪ ð0.60; 0.66Þ GeV=c 2 (for the η mode). The signal yields at each energy, presented in Tables I and II, are obtained according where N is the number of events and the subscript A denotes the signal region, and the subscripts B and C denote the sideband regions.
The Born cross section is calculated from where L is the integrated luminosity, ϵ is the detection efficiency, B is the product of the BF of K 0 S → π þ π − and that of π 0 =η → γγ [57], 1 j1−ΠðsÞj 2 is the vacuum polarization correction factor [69], and ð1 þ δ ISR Þ is the ISR correction factor [70] which is determined by the MC simulation programmer KKMC. The ISR factors are set to 1.0 to get the initial cross section line shape as input to KKMC. From KKMC, the updated ISR factors are obtained, then the cross section line shape is updated too. We repeat this process until both ISR factors and cross section converge.
The  Tables I  and II, and illustrated in Fig. 4.
B. Upper limits of e + e − → Yð4260Þ → K 0 S K AE π ∓ π 0 and e + e − → Yð4260Þ → K 0 S K AE π ∓ η Since there is no obvious structure in the line shapes of the Born cross sections for e þ e − → K 0 S K AE π ∓ π 0 and e þ e − → K 0 S K AE π ∓ η, as shown in Fig. 4, the upper limits of Yð4260Þ → K 0 S K AE π ∓ π 0 and Yð4260Þ → K 0 S K AE π ∓ η are determined by fitting the line shapes with the func- describes the continuum process e þ e − → K 0 S K AE π ∓ π 0 =η, the parameters p 0 is left free in the fit, while p 1 is fixed to the result from a fit that only uses fð ffiffi ffi s p Þ to fit the line shapes. BWð ffiffi s p Þ given in Eq. (2) is a Breit-Wigner function describing the resonance Yð4260Þ, where M, Γ tot , and Γ e þ e − are the mass, full width, and electronic width of Yð4260Þ, respectively; B is the branching fraction of the decay Yð4260Þ → K 0 S K AE π ∓ π 0 =η. The mass and the full width of Y(4260) are set to the world average values 4230 AE 8 MeV=c 2 and 55 AE 19 MeV=c 2 [57]. The product Γ e þ e − B increases from 0 to 0.5 eV in step length of 0.001 eV. For each value of it, a fitting estimator Q 2 defined by Eq. (3) is obtained. Here σ B and σ fit B are the measured and fitted Born cross sections, δ i is the energy dependent part of the total uncertainty, which includes the statistical uncertainty and the energy dependent part of systematic uncertainty,   the δ c is the energy independent part of the systematic uncertainty (the systematic uncertainties are described in detail in Sec. III D), h is a free parameter introduced to take into account the correlation of different energy points, and the subscript i indicates the index of each energy point [71]. The Q 2 is used to calculate the likelihood L ¼ e −0.5Q 2 , whose normalized distribution is used to get the upper limits of Γ e þ e − B at the 90% confidence level (C.L.), which is determined to be 0.050 and 0.19 eV for the π 0 mode and the η mode, respectively.
C. Upper limits on σ B ðe + e − → π 0;∓ Z 0;AE c ð3900Þ; Z 0;AE c ð3900Þ → K 0 S K AE π ∓;0 =ηÞ Since there is no obvious Z c ð3900Þ signal in the invariant mass distributions of K 0 S K AE π ∓;0 (π 0 mode) and K 0 S K AE η (η mode), as shown in Figs The normalized likelihood distribution of the Born cross section Lðσ B Þ is determined by changing the number of signal events from 0 to 150 with a step size of 1. The upper limit (UL) at the 90% C.L. is calculated by solving the equation The final upper limits are shown in Table III, where all of the systematic uncertainties have been considered, the details of which are explained in Sec. III D. The ratio is also given in Table III, where the cross sections for e þ e − → π ∓ Z c ð3900Þ AE → π þ π − J=ψ and e þ e − → π 0 Z c ð3900Þ 0 → π 0 π 0 J=ψ are from Refs. [48] and [72], respectively.
GeV as examples). The black error bars are the data, the black solid lines are the fit curve, the green dashed lines denote the background shape, and the red dashed lines denote the signal shape.

D. Systematic uncertainties
Various sources of systematic uncertainty are investigated in the e þ e − → K 0 S K AE π ∓ π 0 =η cross section measurements. We assume that the systematic uncertainties associated with the physics model used in the MC simulation, the luminosity, tracking, PID, γ reconstruction efficiency, K 0 S reconstruction efficiency, ISR correction factor, vacuum polarization factor and quoted BFs are energy independent, while the other systematic effects are energy dependent.
For the π 0 mode, a data-driven MC method is developed to obtain the efficiency. To estimate the uncertainty of this method, one thousand testing samples of e þ e − → K 0 S K AE π ∓ π 0 are generated with eighteen different physics processes with random ratios, the ratio of each process is generated using uniform distribution between 0 to 1 and then normalized by the summation of these eighteen ratios. The difference between the estimated and the real efficiencies is fitted with a Gaussian function. The fit results give a mean of 0.4% which is neglected, and a width of 0.9% which is taken as the systematic uncertainty from the datadriven MC method. For the η mode, which has much lower statistics than the π 0 mode, alternative mixing ratios are used to generate a new MC sample and the efficiency difference between the two MC samples is adopted as the systematic uncertainty.
The uncertainty on the integrated luminosity is estimated to be 1.0% using Bhabha events [60].
The uncertainty associated with the K 0 S reconstruction is studied using J=ψ → K Ã ð892Þ AE K ∓ ; K Ã ð892Þ AE → K 0 S π AE and J=ψ → ϕK 0 S K AE π ∓ control samples and is estimated to be 1.2% [75].
The ISR correction factor introduces a 1.0% uncertainty since the termination condition of the recursion method used to get the correction factor is 1.0% between the last two iterations.
The uncertainty due to the vacuum polarization factor is found to be negligible [69]. The uncertainties of the quoted BFs are also considered.
The energy dependent ones include the systematic uncertainties from the choosing about mass window and sideband regions of K 0 S , π 0 , and η and the kinematic fit. The uncertainties associated with the K 0 S , π 0 , and η invariant mass regions are determined by changing them from (0.488, 0.508) to ð0.483; 0.513Þ GeV=c 2 , (0.12, 0.15) to ð0.115; 0.155Þ GeV=c 2 and (0.52, 0.58) to ð0.51; 0.59Þ GeV=c 2 for the K 0 S , π 0 and η, respectively. The differences in the efficiencies are taken as the corresponding systematic uncertainties.
The uncertainty associated with the kinematic fit is determined by comparing the efficiencies with and without corrections to the track helix parameters [76].
Assuming all sources of systematic uncertainties are independent, the total uncertainties are the sums of the individual values in quadrature (Table IV).
The systematic uncertainties that affect the upper limits on σ B ðe þ e − → πZ c ð3900Þ; Z c ð3900Þ → K 0 S Kπ=ηÞ are considered in two categories: multiplicative and nonmultiplicative. The nonmultiplicative systematic uncertainties on the signal shape and the background shape are considered by changing the signal shape to a Breit-Wigner function and varying the fit range, the parameters of the Z c ð3900Þ, and the order of the polynomial functions in the fit. The maximum upper limits are adopted for all combinations of these variations. The intermediate states in the Z c ð3900Þ decay are considered by generating signal MC samples with alternative processes Z c ð3900Þ → K Ã ð892ÞK, K Ã ð892Þ → KðK 0 S Þπ (π 0 mode), and Z c ð3900Þ → a 0 ð980Þη, a 0 ð980Þ → K 0 S K (η mode). The efficiency difference is considered as a multiplicative systematic uncertainty. All of the systematic uncertainties, which are listed in Table IV, excluding the side-band item and mixing MC item, are considered as the multiplicative systematic uncertainties. The effects of multiplicative systematic uncertainties are taken into account by convolving the distribution of Lðσ B Þ with a probability distribution function of sensitivity (S), which is assumed to be a Gaussian function with central valueŜ and standard deviation δ S [77]: Here S is the sensitivity that refers to the denominator of Eq. (1) and δ s is the total multiplicative systematic uncertainty. L 0 ðσ B Þ is the likelihood distribution of the Born cross section after the multiplicative systematic uncertainties are incorporated.

IV. SUMMARY
The Born cross sections for e þ e − → K 0 S K AE π ∓ π 0 and K 0 S K AE π ∓ η are measured with data samples collected at center-of-mass energies from 3.90 to 4.60 GeV. Since no clear structure is observed, the upper limits of the product Γ e þ e − B ðYð4260Þ → K 0 S K AE π ∓ π 0 Þ at 90% C.L. are estimated to be less than 0.05 eV and that of Γ e þ e − BðYð4260Þ → K 0 S K AE π ∓ ηÞ is estimated to be smaller than 0.19 eV. Reference [9] reported four solutions of the product Γ e þ e − BðYð4260Þ → π þ π − J=ψÞ, in which the maximum is 13.3 AE 1.4 eV and the minimum is 1.5 AE 0.3 eV. Comparing them with our results, the branching fraction of the Yð4260Þ decaying into K 0 S K AE π ∓ π 0 and K 0 S K AE π ∓ η is much smaller, which indicates a much smaller coupling of the Yð4260Þ to the light hadrons K 0 S K AE π ∓ π 0 and K 0 S K AE π ∓ η. We also search for e þ e − → πZ c ð3900Þ; Z c ð3900Þ → K 0 S Kπ=η and no obvious Z c ð3900Þ signal is observed in the charged or neutral mode.  IV. Summary of systematic uncertainties (in %). The systematic uncertainties listed in this table are related to the born cross section measurements of processes e þ e − → K 0 S K AE π ∓ π 0 =η and the upper limits extracted work of e þ e − → Yð4260Þ → K 0 S K AE π ∓ π 0 =η and e þ e − → πZ c ð3900Þ → K 0 S K AE π ∓ π 0 =η.