$Z_c(3900)$: what has been really seen?

The $Z^\pm_c(3900)/Z^\pm_c(3885)$ resonant structure has been experimentally observed in the $Y(4260) \to J/\psi \pi\pi$ and $Y(4260) \to \bar{D}^\ast D \pi$ decays. This structure is intriguing since it is a prominent candidate of an exotic hadron. Yet, its nature is unclear so far. In this work, we simultaneously describe the $\bar{D}^\ast D$ and $J/\psi \pi$ invariant mass distributions in which the $Z_c$ peak is seen using amplitudes with exact unitarity. Two different scenarios are statistically acceptable, where the origin of the $Z_c$ state is different. They correspond to using energy dependent or independent $\bar D^* D$ $S$-wave interaction. In the first one, the $Z_c$ peak is due to a resonance with a mass around the $D\bar D^*$ threshold. In the second one, the $Z_c$ peak is produced by a virtual state which must have a hadronic molecular nature. In both cases the two observations, $Z^\pm_c(3900)$ and $Z^\pm_c(3885)$, are shown to have the same common origin, and a $\bar D^* D$ bound state solution is not allowed. Precise measurements of the line shapes around the $D\bar D^*$ threshold are called for in order to understand the nature of this state.

The resonant-like structure Z c (3900) ± was first seen simultaneously by the BESIII and Belle collaborations [1,2] in the J/ψπ spectrum produced in the e + e − → Y(4260) → J/ψπ + π − reaction. An analysis [3] based on CLEO-c data for the e + e − → ψ(4160) → J/ψπ + π − reaction confirmed the presence of this structure as well, although with a somewhat lower mass. Under a different name, Z c (3885) ± , a similar structure, with quantum numbers favored to be J P = 1 + , has also been reported by the BESIII collaboration [4,5] in theD * D spectrum of e + e − →D * Dπ at different e + e − center-of-mass (c.m.) energies [including the production of Y(4260)]. Because there is a little difference in the central values of the masses and in particular the widths of these two structures, whether they correspond to the same state is still unknown. As will be shown in this Letter, the two structures have indeed the same common origin. We generically denote it here as Z c . Evidence for a neutral partner of this structure was first reported in Ref. [3], and more recently in Ref. [6].
If this resonant structure happens to be a real state as argued in Ref. [7], it is one of the most interesting hadron resonances, since it couples strongly to charmonium and yet it is charged, thus it is something clearly distinct of a conventional cc state -its minimal constituent quark content should be four quarks, ccud (for Z + c ). A discussion of possible internal structures is given in Ref. [8]. It has been interpreted as a molecularD * D state [9][10][11], as a tetraquark of various configurations [12] or as a simple kinematical effect [13], although this possibility has been ruled out in Ref. [7]. Distinct consequences of some of these different models have been discussed in Ref. [14]. It has been also searched for in lattice QCD though with negative results so far [15].
Being a candidate for an explicitly exotic hadron, the Z c (3900) definitely deserves a detailed and careful study. Indeed, the last years have witnessed an intense theoretical activity aiming at understanding the actual nature of this state. What is still missing, however, is a simultaneous study of the two reactions analysed by BESIII and mentioned above in which the Z c structure has been seen. 1 The goal of this work is to perform such a study, and, from it, to extract information about this seemingly resonant intriguing structure. We will first settle aD * D, J/ψπ coupled channel formalism, considering that the Z c emerges from theD * D interaction, and that its coupling to J/ψπ proceeds through the former intermediate state. The resulting T -matrix will enter the calculation of the amplitudes for the reactions Y(4260) → J/ψππ,D * Dπ. We will assume that the Y(4260) state is dominantly a D 1 (2420)D+c.c. bound state [9,16] and use the ideas of Ref. [9] to compute the relevant amplitudes.
Let us denote with 1 and 2 the J/ψπ andD * D channels, respectively, with I = 1 and J PC = 1 +− (here and below, the C-parity refers to the neutral member of the isospin triplet). The coupled-channel T -matrix can be written as where G is the loop function diagonal matrix, and the matrix elements of the potential read where m i n is the mass of the the nth particle in the channel i, and the mass factors are included to account for the non-relativistic normalization of the heavy meson fields. The J/ψπ → J/ψπ interaction strength is known to be tiny [17,18], and we neglect the direct coupling of this channel, C 11 = 0. Such a treatment was also done in Ref. [19] in a coupled-channel analysis of the Z b states. For the inelas-ticD * D → J/ψπ S -wave interaction, we make the simplest possible assumption, that amounts to take it to be a constant, C 12 ≡ C. In a momentum expansion, the lowest order contact potential for theD * D →D * D transition is simply a constant as well, denoted by C 22 ≡ C 1Z [20]. However, it can be shown that even with two coupled channels, no resonance can be generated in the complex plane above threshold with only constant potentials. To that end, we will also allow some energy dependence for the V 22 term, introducing a new parameter b, and writing with E the total c.m. energy. The new term is of higher order in low-momentum expansion in comparison with C 1Z . The interactions considered here need to be regularized in some way, and hence we employ a standard gaussian regulator [21], where the c.m. momentum squared of the channel i is denoted by q i 2 . We adopt a relativistic (non-relativistic) definition of the latter for the i = 1 (i = 2) channel, i.e., being µ the reduced mass of theD * D system. Since the interaction for this channel is derived from a non-relativistic field theory, we take cutoff values Λ 2 = 0.5 − 1 GeV [20]. At the Z c energy, the c.m. momentum of the J/ψπ channel is q 1 ≃ 0.7 GeV, and hence we use a different cutoff for it. For definiteness, we set Λ 1 = 1.5 GeV, although the specific value is not very relevant as we have checked since changes in the cut-off can be reabsorbed in the strength of the transition potential controlled by the undertermined C 12 low energy constant. With this convention for the regulator, the loop functions in the matrix G read with ω n = l 2 + m 2 1n . The DD * channel loop function G 2 is computed in the non-relativistic approximation.
For the e + e − annihilations at the Y(4260) mass, both BE-SIII and Belle have reported the Z c structure in the J/ψπ final state [1,2], but only BESIII provides data for theD * D channel [4,5]. Hence, for consistency, we will only study the BESIII data. In particular, we will consider the most recent double-D-tag data of Ref. [5], in which the D * is reconstructed from several decay modes, whereas in Ref. [4] the presence of the D * is only inferred from energy conservation. Hence, in the former data the background in the higher energy D * D invariant mass regions is much reduced. For definiteness, we will consider the reported spectra of the D * − D 0 and J/ψπ − final states, and set m D * = m D * − , m D = m D 0 , and m π = m π ± . This implicitly assumes that isospin breaking effects are neglected. These data are taken at a c.m. energy equal to the nominal Y(4260) mass, so the decays to π(J/ψπ,D * D) proceed mainly through the formation of this resonance. The mechanisms for the Y(4260) decays are shown in Fig. 1. The coupling YD 1 D, whose value is not important here to describe the lines shapes, is taken from Ref. [9], where the Y(4260) is considered to be dominantly aDD 1 + c.c. bound state. The subsequent D 1 D * π coupling can also be found there.
We denote M 1 (M 2 ) to the amplitude for the Y → J/ψπ + π − (Y → D * − D 0 π + ) decay, and s and t, respectively, to the invariant masses squared of J/ψπ − and J/ψπ + (D * − D 0 and D * − π + ) in the first (second) decay. Up to some common irrelevant constant, both amplitudes can be written (after the appropriate sum and average over polarizations) as: where , and θ denotes the relative angle between the two pions in the Y(4260) rest frame. Further, I 3 (s) is the scalar three-meson non-relativistic loop function, for which details can be found in Ref. [22]. One first notes that M 1 (s, t) is symmetric under s ↔ t. The term with α represents diagram (1a), and it acts as a non-resonant background amplitude, added coherently to the rest of the diagrams. It has the same dependence on the external momenta and polarization vectors as that of diagrams (1b)-(1e). The first term in M 1 (s, t) is the amplitude of diagrams (1b)+(1c), the second term is the one from diagrams (1d)+(1e), and the last one is their interference. In M 2 , the first summand of the first term corresponds to diagram (2a) in Fig. 1, whereas the second one, which includes theD * D final state interaction (FSI), is the contribution from diagrams (2b)+(2c). Diagrams (2a)-(2c) proceed through the formation of D 1 , but we also consider some non-resonantD * Dπ production by means of diagram (2d). TheD * D rescattering effects in this last diagram give rise, in turn, to diagrams (2e) and (2f). The term with β in Eq. (8) represents these latter three diagrams. The parameters α and β in Eqs. (7) and (8) are unknown. Note that the effect of D 1 width, Γ D 1 = (25 ± 6) MeV, is negligible here since m D 1 + m D − Γ D 1 /2 is well above 4.26 GeV. 2 The spectrum for both reactions can be obtained as a contribution from the amplitudes (A i ) plus a background (B i ): where t i,± (s) are the limits of the t Mandelstam variables for the decay mode i. The two global constants K i could be related if the event selection efficiencies of the two spectra analyzed in this work were known. If the latter were roughly the same, then one would have K 1 ≃ 5K 2 (due to the different bin sizes). If both parameters are considered free, a large correlation arises between K 1 and C, since K 1 | C| 2 basically determines the total strength of the event distribution N 1 . This is due to the fact that the influence of C in the shape of the Tmatrix elements, and thus of the signal of Z c in the spectrum, is small. To obtain a reasonable estimate of this coupling constant, we consider a further experimental input from Ref. [4], and estimate this ratio as that is, as the ratio of the background subtracted areas of each physical spectrum around the Z c mass, namely in the range √ s = (3900 ± 35) MeV. In principle, the double-D-tag technique ensures that all thē D * D spectrum events in Ref. [5] contain aD * D pair, so there is no background due to wrong identification of the final state. There could be, however, contributions to the spectrum from higher waves other than the S -wave. In any case, an inspection of Fig. 2 shows that the tail of the spectrum is small, and we set B 2 = 0. We shall come to this point later on. For the J/ψπ spectrum, B 1 is parameterized with a symmetric smooth threshold function as used in the experimental work of Ref. [1]: with m 1− = m J/ψ + m π and m + = m Y − m π , i.e., the limits of the available phase space for the reaction. The parameters B 1 and d 1 are free.
We have three free parameters directly related to our Tmatrix (C 1Z , C, and b), and six (B 1 , d 1 , α, β and K 1,2 ) related to the background and the overall normalization. These nine free parameters are adjusted to reproduce the data of Refs. [1,5] (a total of 104 data points). In this work, two errors are given. The first error is statistical and it is computed from the hessian matrix of the χ 2 merit function. The second error is systematic, and to estimate it we have considered two different uncertainty sources. First, we have varied the J/ψπ background function [Eq. (13)] and used other smooth functions. The second source of uncertainties is related to the tail of thē D * D spectrum, and it is estimated as follows. The central value of the parameters is computed by fitting this spectrum up to √ s = 4025 MeV. Then, we vary this limit between √ s = 3975 MeV and m + (the maximum allowed invariant mass), and repeat the fit. In all cases, we find statistically acceptable fits and the difference between the new fitted parameters and the central ones is used to determine the systematic error. The same method is applied to estimate the systematic error of our predictions for the spectra and the mass and width of the Z c state, to be presented below.
We perform four different fits, corresponding to the two cases of keeping the parameter b, which controls the energy dependence of theD * D potential, free or set to zero, and for each of these, we choose Λ 2 to be 0.5 or 1 GeV [20]. Results from the four fits are compiled in Table I, where only the parameters that are directly related to our T -matrix are shown. One first notes that the reduced χ 2 is very close to unity in all four cases. Indeed, the description of the experimental spectra is very good in all cases, as can be seen in the top panels of Fig. 2, where the results from one of the fits (b free and Λ 2 = 0.5 GeV) are shown and confronted with the data. In particular, the effect of the Z c is nicely reproduced in theD * D spectrum above threshold and in the J/ψπ spectrum around theD * D threshold. Its reflection can also be appreciated in the J/ψπ distribution around 3.5 GeV. The other fits lead to results similar to those shown in Fig. 2. The largest differences can be found in theD * D spectrum between the b 0 and b = 0 cases, which are compared for Λ 2 = 0.5 GeV in the bottom right panel of the same figure. In any case, we see that we are able to simultaneously reproduce the two available BESIII data sets related to the Z ± c (3900)/Z ± c (3885) state with a single structure for the very first time.
The numerical values are shown in Table II. threshold. For the case b = 0, however, the situation is quite different. While the description of the experimental data is still quite good with χ 2 /d.o.f. ∈ [1.3, 1.4], the pole in this case is located below threshold, with a small imaginary part (around 8 MeV), and in the (01) Riemann sheet. If the J/ψπ channel were now switched off ( C = 0), this pole would move into the real axis in the unphysical Riemann sheet of the elastic amplitude T 22 . In this sense, the obtained pole does not qualify as a resonance, and we see it as a virtual or anti-bound DD * state. It does not correspond to a particle in the sense that its wave function, unlike that of a bound state, is not localized. However, it produces observable effects at the DD * threshold similar to those produced by a near threshold resonance or bound state. 3 Indeed, scattering experiments alone, in principle, cannot distinguish between virtual and bound states, but the difference is not a purely academic one since they can produce different line shapes in inelastic open channels [24]. The line shapes of a virtual state and a near-threshold resonance are different since the former is peaked exactly at the threshold while the latter, in principle, is above. This can be seen in the left bottom panel of Fig. 2 where the Jψπ − spectrum for the two fits b = 0 and b 0 are shown (for the case Λ 2 = 0.5 GeV). Although the two curves are different, each one would approximately lie within the error band of the other. Clearly, very precise data with a good energy resolution and small bin size are necessary to distinguish among them.
Without taking sides, and given that both natures for the Z c structure (resonance or virtual state) arise in fits of good quality, it must be stated that the experimental information available at this time cannot fully discriminate between both scenarios and, hence, claims about the Z c structure should be made with caution. Nevertheless, the resonance scenario seems to be statistically slightly preferred. It is also clear that more experimental information is needed to elaborate on the nature of Z c . In particular, the spectrum of J/ψπ with narrower bins would be highly desirable to have a good resolution on its line shape. If it is finally shown to be a virtual state, then it cannot be a tetraquark, since it does not correspond to a normal particle, and it can only have a hadronic molecular nature, in the sense that it appears only because of the DD * interaction.
Summarizing, we have studied the two decays (Y(4260) → J/ψπ + π − , D * − D 0 π + ) in which the Z ± c resonant-like structure is seen. We have presented the first simultaneous study of the invariant mass distributions of the J/ψπ andD * D channels with fully unitarized amplitudes. We find that these data sets are well reproduced in two different scenarios. In the first one, in which there is an energy dependence in theD * D →D * D potential, the Z c appears as a dynamically generatedD * D resonance. In the second one, however, when the aforementioned energy dependence is not allowed, it appears as a virtual state, with the pole located below the DD * threshold. In any case, it is demonstrated that both data sets can be reproduced with only one Z c state, so that the two experimentally observed structures Z ± c (3900) and Z ± c (3885), in different channels, are proven to correspond to the same state. Moreover, both fits do not allow aD * D bound state solution. 4 Since the virtual state can only be of hadronic-molecule type, it is really important to discriminate between these two scenarios. For that purpose, one needs a very precise measurement of the line shapes around, in particular slightly above, the DD * threshold. Such a measurement is foreseen when more data are collected at BESIII.