Study of e + e − → D + D − π + π − at center-of-mass energies from 4.36 to 4.60 GeV

We report a study of the e + e − → D + D − π + π − process using e + e − collision data samples with an integrated luminosity of 2 . 5 fb − 1 at center-of-mass energies from 4.36 to 4.60GeV, collected with the

The numbers relevant to the Born cross section measurements, where the first uncertainties are statistical, the second are independent systematic uncertainties, and the third are common systematics. The index of 1 represents the process e + e − → D 1 (2420) + D − + c.c. → D + D − π + π − while the index of 2 represents the process e + e − → ψ(3770)π + π − → D + D − π + π − .

Introduction
Recent discoveries of charmonium-like states that do not fit naturally with the predicted charmonium states in the quark model have stirred up great experimental and theoretical interests [1][2][3][4][5]. Among these so-called X Y Z states, the observations of the Y (4260) [6] and Z c (4430) [7] states have drawn special attention, and stimulated extensive discussions on their structures. Some calculations indicate that the Y (4260) is possibly a D 1 (2420)D molecular state, while the Z c (4430) is possibly a D 1 (2420)D * molecular state [8][9][10][11]. Hence, more studies on the properties of the involved D 1 (2420), such as mass and width, are helpful to better understand the nature of these exotic candidate states.
The lightest charmonium state above the DD threshold is the ψ(3770) resonance, which is considered to have the quantum numbers of 1 3 D 1 [12,13]. Its spin-triplet partner 1 3 D 2 candidate, X(3823), has been observed in the process e + e − → X(3823)π + π − at BESIII [14]. Analogously, it is interesting to study the production of the ψ(3770) in the process e + e − → ψ(3770)π + π − [15], which is observed at √ s =4415.6 MeV at BE-SIII [16]. More precise measurements at different energy points are desired, as it provides an important way to investigate the intrinsic nature of the Y (4360) and ψ(4415) by studying the transitions between these charmonium(-like) states, such as In this analysis, we study the process e + e − → D + D − π + π − at the center-of-mass (c.m.) energies, E c.m. , from 4358.3 to 4599.5 MeV, as listed in Table 1. Compared to the process e + e − → D 0D0 π + π − , this final state has the advantage of being free from D * intermediate states, which greatly simplifies the analysis. We reconstruct the D + via its high branching fraction decay K − π + π + and adopt a recoil-mass technique to identify the D − and related resonant states. Unless explicitly mentioned otherwise, inclusion of charge conjugate mode is implied throughout the context. Clear signals of the D 1 (2420) + and ψ(3770) are extracted in this data 10  set via their decays to D + π + π − and D + D − , respectively. The resonance parameters of the D 1 (2420) + are measured. Additionally, the Born cross sections of e + e − → D 1 (2420)

The experiment and data sets
The BESIII detector is a magnetic spectrometer [17] located at the Beijing Electron Positron Collider (BEPCII) [18]. The cylindrical core of the BESIII detector 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 identifier modules interleaved with steel. The acceptance of charged particles and photons is 93% over 4π solid angle. The charged-particle momentum resolution at 1 GeV/c is 0.5%, and the dE/dx resolution is 6% for the 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 of the TOF barrel part is 68 ps, while that of the end cap part is 110 ps.
The E c.m. of the seven data sets are measured using di-muon events [19], and the corresponding luminosities are measured with large-angle Bhabha scattering events [21]. To optimize selection criteria, estimate the detection efficiency, and understand background contributions, we simulate the e + e − annihilation processes with the kkmc [22] generator, which takes into account continuum processes, initial state radiation (ISR), and inclusive D ( * ) (s) production. The known decay rates are taken from the Particle Data Group (PDG) [13], and the decays are modeled with evtgen [23]. The remaining decays are simulated with the lundcharm package [24].

Event selections
To reconstruct the D + meson, charged track candidates for one K − and two π + in the MDC are selected. For each track, the polar angle θ defined with respect to the e + beam is required to satisfy |cosθ| < 0.93. The closest approach to the e + e − interaction point is required to be within ±10 cm along the beam direction and within ±1 cm in the plane perpendicular to the beam direction. A track is identified as a π(K ) when the PID probabilities satisfy P(π ) > P(K ) (P(K ) > P(π )), according to the information of dE/dx and TOF. We reconstruct D + candidates by considering all possible combinations of the charged tracks which are required to originate from a common vertex. The quality of the vertex fit is required to satisfy χ 2 VF < 100. We constrain the reconstructed D + mass with a kinematic fit to the nominal D + mass [13], and require the fit quality χ 2 KF < 20. We then require the presence of one additional π + π − pair, with neither track used in the reconstructed D + . The identification of the signal process e + e − → D + D − π + π − is based on the recoil mass spectra of D + π + π − , R M(D + π + π − ), which are shown in Fig. 1. The rate of multiple candidates per event is about 10%, and is corrected for via the MC efficiency.
The peaks observed at 1.87 GeV/c 2 correspond to the D − meson signals. They are consistent with the MC simulations of the D + D − π + π − final state. The background contributions are due to random combinations of charged tracks. We further restrict the candidate events to the region 1.855 < R M(D + π + π − ) <

Signal extraction
The 2-dimensional distributions of M(D + π + π − ) versus R M(D + ) for the D 1 (2420) + D − are shown in Fig. 3. The vertical band corresponds to the D 1 (2420) − signal and the horizontal band corresponds to the D 1 (2420) + . The projection to the R M(D + ) axis (Fig. 2) consists of a prominent D 1 (2420) − peak and a corresponding broad bump. The contributions of D 1 (2420) + D − and ψ(3770)π + π − in the selected data are determined using fits to the R M(D + ) one-dimensional distribution. The shape of this distribution is described using templates obtained from the signal MC simulation. In order to perform a likelihood scan of the resonance parameters, we generate a series of D 1 (2420) + signal MC with different values of mass and width, and smear these template shapes with a Gaussian function to take into account the resolution difference between data and MC simulations. The width of the Gaussian function is fixed to the difference of resolution in R M(D + ) for the control sample of e + e − → D + D − . The signal shape for the mode ψ(3770)π + π − is obtained from the MC simulation, where the resonance parameters of the ψ(3770) are taken from the PDG [13].
A simultaneous unbinned maximum likelihood fit to the data samples is performed at three high luminosity energy points of (black) solid lines are the sum of statistical uncertainties and independent systematic uncertainties in quadrature, the (red) dot lines are total uncertainties. ble 1. Here, the contribution of the non-resonant four-body process e + e − → D + D − π + π − is neglected in the fit, as an alternative fit including this process gives its size consistent with zero.
In addition, we analyze the data samples at E c.m. = 4487.4, 4467.1, 4527.1 and 4574.5 MeV with relatively low luminosities. We apply the same strategy to extract the signal yields of the D 1 (2420) + D − and ψ(3770)π + π − , except that we fix the resonance parameters for the D 1 (2420) + according to the aforementioned fit results.

Cross section measurement
The Born cross section is calculated with where index i denotes the respective signal process, n sig i is the observed signal yield, L is the integrated luminosity, B is the branching fraction B(D + → K − π + π + ) = (9.38 ± 0.16)% [13], ε i is the detection efficiency, (1 + δ rad i ) is the radiative correction factor which is obtained from a QED calculation using the line shape of the data cross section of signal process as input in an iterative procedure, and 1 |1− | 2 is the vacuum polarization factor [26]. The trigger efficiencies for the two processes are 100%, as there are at least 5 charged tracks detected [27]. The processes e + e − → D 1 (2420) + D − + c.c. → D + D − π + π − and e + e − → ψ(3770)π + π − → D + D − π + π − are denoted with index i = 1 and i = 2, respectively. The calculated Born cross sections are given in Table 1 and plotted in Fig. 4. We evaluate the statistical significance by the ratio of the maximum likelihood value and the likelihood value for a fit with a null-signal hypothesis. For the energy points with low statistical significances, we determine the upper limits for the cross sections which are calculated by using the signal yield upper limits n UL in Eq. (1). The upper limit n UL at 90% confidence level is obtained with a Bayesian approach scanning the expected signal yield. The probability is calculated from the Gaussian-smeared likelihood to take into account the systematic uncertainty.

Systematic uncertainties
The systematic uncertainties of the measurement of the D 1 (2420) + resonance parameters and the Born cross sections listed in Tables 2 and 3 include correlated (common) contributions, from tracking, PID, luminosity measurements, vacuum polarization factors, interference effect and the input branching fraction, as well as uncorrelated (independent) contributions from background shapes, mass scaling, detector resolution, signal shape due to the angular distributions, and radiative corrections.
• Uncertainties of tracking and PID are each 1% per track [28].   Summary of systematic uncertainties on the D 1 (2420) + resonance parameters and the Born cross sections for the high luminosity energy points.   • The uncertainty of modeling the angular distributions of the signal processes are studied by repeating the analysis procedure on the basis of new signal model. For e + e − → D 1 (2420) + D − , we considered two extreme cases of 1 + cos 2 θ D 1 and 1 − cos 2 θ D 1 , where θ D 1 is the helicity angle of the D 1 (2420) + in the rest frame of the initial e + e − system.
For e + e − → ψ(3770)π + π − , a model, named as JPIPI [23] in evtgen, is considered. The maximum changes on the results are taken as systematic uncertainties.
• The uncertainty of luminosity measurement is 1%, as given in Ref. [21].
• The uncertainty of radiative correction is calculated by using the generator kkmc. Initially, the observed signal events are assumed to originate from the Y (4260) resonance to obtain the efficiency and ISR correction factor. Then, the measured line shape is used as input to calculate the efficiency and ISR correction factor again. This procedure is repeated until the difference between the subsequent iterations is comparable with the statistical uncertainty. We take the difference of the radiative correction factors between the last two iterations as the systematic uncertainty.
• We take 0.1% as the uncertainty of the vacuum polarization factor, which is calculated in Ref. [26].
• The input branching fraction of D + → K − π + π + in PDG has the relative uncertainty of 1.7%, which is taken into account.
The systematic uncertainties are summarized in Tables 2 and 3; the sum of different uncertainties are obtained by adding up all the relevant contributions in quadrature.
The Born cross sections of e + e − → D 1 (2420) + D − + c.c. → D + D − π + π − and e + e − → ψ(3770)π + π − → D + D − π + π − are measured as functions of the center-of-mass energy. The cross section line shape is consistent with previous BESIII measurement based on full reconstruction method [16]. There are some indications of enhanced cross sections for both processes between 4.36 and 4.42 GeV, where the reported states Y (4360) and ψ(4415) locate. Hence, the measured cross sections can be useful inputs to the properties of these states.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.