Searching for a $D \bar D$ bound state with the $\psi(3770) \to \gamma D^0 \bar D^0$ decay

We perform a calculation of the mass distribution in the $\psi(3770) \to \gamma D \bar D$ decay, studying both the $D^+ D^- $ and $D^0 \bar D^0 $ decays. The electromagnetic interaction is such that the tree level amplitude is null for the neutral channel, which leaves the $\psi(3770) \to \gamma D^0 \bar D^0$ proceeding through a loop involving the $D^+ D^- \to D^0 \bar D^0$ scattering amplitude. We use the results for this amplitude of a theoretical model that predicts a $D \bar D$ bound state and find a $D^0 \bar D^0 $ mass distribution in the decay drastically different than phase space. The rates obtained are relatively large and the experiment is easily feasible in the present BESIII facility. The performance of this experiment could provide an answer to the issue of this much searched for state, which is the analogue of the $f_0(980)$ resonance.


I. INTRODUCTION
Molecular states made from mesons or mesons and baryons have become some of the important objects in the present plethora of hadronic states. Reviews on this topic can be found in [1,2] and more recently in [3][4][5][6]. One of the interesting predicted states is a bound state of DD [7][8][9] for which there is no clear experimental evidence so far. The state is analogous to the KK bound state which was claimed in [10] to be a good representation of the f 0 (980) state. This early claim found support later with the results of the chiral unitary approach for the meson-meson interaction [11][12][13][14], where the meson meson interaction in coupled channels was studied, and other states, as the f 0 (500), K * 0 (700), a 0 (980), were also found dynamically generated from that interaction. There has been some search for this state and in [15] it was shown that an accumulation of strength close to the DD threshold in the e + e − → J/ψDD reaction [16] found a natural explanation in terms of the predicted bound state of [7] with a mass of 3730 MeV, with large uncertainties. Hopes where raised that an update of the experiment in [17] would constrain the predictions, but it was shown in [18] that this is not the case, and there is a large ambiguity in the conclusions.
In view of this, there have been some works proposing new reactions that would give evidence for this elusive state. In [19] the radiative decay of the ψ(3770) resonance into γ and the DD bound state was proposed and the feasibility of the reaction with present production rates of the ψ(3770) was assessed. There is of course the problem of which should be the ideal channel to observe the bound state. In [20] three different reactions were suggested to detect that state. In [21] the B 0 → D 0D0 K 0 , B + → D 0D0 K + reactions were suggested to find evidence for the state, looking into the mass distribution of DD production close to threshold. The analysis found a good agreement with experiment for the B + → D 0D0 K + reaction, but it was shown there that the B 0 → D 0D0 K 0 reaction was better suited to search for the DD state, because there is no tree level contribution in this later reaction, and hence the amplitude for that mechanism is proportional to the DD amplitude which contains the pole of the DD bound state. The mass distribution then was rather different from that of phase space, and its precise measurement close to threshold should give and answer to the question.
In the present work we wish to combine the lessons drawn from the works of [19] and [21] and study the DD mass distribution close to threshold from the ψ(3770) → γDD decay. Anticipating the results, we will find an interesting situation in which the D 0D0 production does not proceed at tree level, while D + D − has contribution from tree level. As a consequence, the D 0D0 production is directly influenced by the DD pole below threshold and exhibits a behavior close to threshold very different from phase space. For the D + D − production the tree level part is very important and the behavior is quite different, and in addition it shows the infrared divergence behavior when the photon energy goes to zero, which, again, is not the case for D 0D0 production. The reaction is, thus, suited for investigation of the DD bound state and the present rates of ψ(3770) production make the experimental investigation feasible.

II. FORMALISM
The state ψ(3770) decays into DD [22] with a width Γ = 27.2 MeV, 52% of which goes to D 0D0 and 41% to D + D − . Its shape in e + e − production and the decay width have been object of intense study [23][24][25][26][27][28][29]. In [29] the cc component is allowed to get hadronized into meson-meson components, and the strength of the hadronization is fitted to the ψ(3770) lineshape. One of the conclusions in [29] is that from the experimental data one can induce that the ψ(3770) is largely a cc state and the weight of the meson-meson components is only of the order of 15%. The ψ(3770) state bears some similarity to the φ(1020), which is assumed to be a ss state and decays into KK. The decay mode ψ(3770) → DDγ necessarily has much resemblance to the φ → γKK, γπ 0 π 0 decays, which have also been a subject of much study [30][31][32][33][34][35][36][37]. As in the φ → γπ 0 π 0 reaction, which does not proceed via tree level, we shall also see that ψ(3770) → γD 0D0 does not get contribution from tree level and both processes proceed via a similar loop mechanism.
A. Tree level for ψ(3770) → γD + D − The tree level mechanism in ψ(3770) → γD + D − is shown diagrammatically in Fig. 1. The ψ → D + D − elementary vertex is given by The ψ → D + D − decay width is given by where the sum and average of |t| 2 calculated from Eq. (1) gives and q is the D + momentum in the ψ decay at rest. Adjusting to the experimental D + D − decay width we find Considering also the γD with e the electron charge, e 2 /4π = α = 1/137, the ψ → γD + D − amplitude of the diagram of Fig. 1 is given by where the term g µν corresponds to the diagram of Fig. 1(c) and is introduced to respect gauge invariance. The photon has zero coupling to D 0 D 0 and hence there is no tree level for ψ → γD 0D0 .

B. Loop mechanism
There is, however, a loop mechanism that allows the ψ → γD 0D0 decay which is depicted in Fig. 2.
FIG. 2: Loop mechanisms for ψ → γDD production. In parenthesis the momenta of the particles.
This follows exactly the same trend as in [30][31][32][33][34][35][36][37] for φ → γπ 0 π 0 , where the intermediate state is K + K − and the final DD are replaced by π 0 π 0 . The diagram of Fig. 2(d) is also demanded by gauge invariance of the loops. Gauge invariance plays an important role in this process and thanks to it there is an efficient computational scheme which requires only the evaluation of diagrams (a) and (b) of Fig. 2, which give the same contribution, and shows that the result of the loop integral is convergent [34,38,39]. The derivation goes as follows: The full amplitude for the diagrams of Fig. 2 has the structure and T µν must be a tensor that can be written in terms of the two independent momenta P and k, the momentum of the ψ and γ respectively. The most general form for T µν is given by Gauge invariance, substituting ν (γ) by k ν and demanding leads to which implies two independent equations The b term does not contribute because of the Lorentz condition µ (ψ)P µ = 0 and the c and e terms do not contribute because of the Lorentz condition on the photon ν (γ)k ν = 0. Hence, only the a and d terms of Eq. (8) contribute to the amplitude and it is enough to calculate only the a or d coefficient. It is easy to see that only the diagrams (a), (b) of Fig. 2 contribute to the d coefficient and since two external momenta P µ k ν are factorized out of the integral, for dimensional reasons this means two powers of q less in the integral, which renders it convergent. In addition, if we work at the end, as we do, in the Coulomb gauge, 0 (γ) = 0, (γ) · k = 0, then the term dP i k j j (ψ) i (γ) = 0 in the ψ rest frame, and the whole amplitude is given by It is customary to perform the integration of the loop integral using Feynman parametrization, but here we must divert from this formalism because the DD → DD scattering matrix regularized with a cut off, q max , transfers a structure Θ(q max − |q|)Θ(q max − |p|) to the T matrix [40] and we must implement a cut off in the loop integral. On the other hand, we can benefit from the fact that the D mesons are heavy particles, they are close to on-shell in the loops and we can just keep the positive energy part of their propagators with ω(q) = q 2 + m 2 D .
The contribution of the two diagrams of Fig. 2(a), 2(b) with D 0D0 in the final state is given by where t D + D − →D 0D0 is a function of the D 0D0 invariant mass. We can now take the propagators of Eq. (14) and perform the q 0 integration analytically using Cauchy's integration, and we find, keeping only the i (γ) transverse components that we shall have in the Coulomb gauge (i, j=1, 2, 3) and in order to get the d coefficient we must look at the k i P j component of this integral. The first term in Eq. (18), with q i P j , provides this structure immediately since when P → 0 for the rest frame of the ψ, the integral only depends on k, hence d 3 q q i f (q, k) = k i d 3 q (k · q)/k 2 f (q, k). The second integral for q i q j is a bit more involved, and since we will have P → 0 at the end, we can make an expansion for P small of all the terms. Then we have an integral at the end of the type and only the δ jl k i will contribute to the term k i P j that we look for. Finally we find with with ω 1 = q 2 + m 2 D and ω 2 = (q + k) 2 + m 2 D . The loop t matrix for ψ(3770) → γD 0D0 can then be put as For the case of ψ(3770) → γD + D − we have the same formalism for the loop substituting t D + D − ,D 0D0 by t D + D − ,D + D − , but we have to add the tree level contribution of Eq. (6). Given the structure of the amplitude it is convenient to evaluate the phase space in terms of the energy of the photon and the D 0 both in the ψ rest frame and we obtain at the end, where E 1 is the energy of the D 0 and A is the cosine of the angle between the photon and D 0 given by The sum and average over spins of |t| 2 is given for the case of D + D − production by For the case of D 0D0 production only t A L + t B L has to be kept and the |t| 2 gives us 2 3 |dM ψ k| 2 , which is what we directly obtain from Eq. (23). All terms appearing in Eq. (26) can be calculated in terms of M inv (D 0D0 ), E 1 and A of Eq. (25).

C. DD → DD transition amplitude
We use the Bethe-Salpeter equation in coupled channels where the channels used are D + D − , D 0D0 , D sDs , plus the ηη channel to account for the decay of the DD bound state into light meson-meson channels, as done in [21]. The V ij coefficients between D + D − , D 0D0 , D sDs are taken from [7] and we take V D + D − ,ηη = a, V D 0D0 ,ηη = a and all the other matrix elements zero [21]. The value of a is chosen at a = 42 in order to obtain a width Γ 36 MeV as found in [20]. The G function is a diagonal function of the meson meson loops for which either dimensional regularization or cut off regularization can be used. If the cut off regularization is used, q max = 830 MeV. As discussed above, the q max has to be also implemented in the q integral of the loop function. It is also interesting to note that since the ηη channel is only introduced to produce the width, it is sufficient to take The contribution of the real part of G ηη only introduces negligible changes in the results [18].

III. RESULTS
First we look at the t D + D − ,D + D − and t D + D − ,D 0D0 amplitudes. In Fig. 3 we plot |t D + D − ,D + D − | 2 and |t D + D − ,D 0D0 | 2 as a function of the DD invariant mass. We can see that the amplitudes are practically identical and have a peak around 3770 MeV corresponding to a DD bound state. This is due to the dominance of the isospin I = 0 contribution. The amplitudes also exhibit a fast fall around threshold corresponding to the opening of the DD decay channel (Flatté effect).
Next we plot the coefficients d 1 , d 2 in Fig. 4 given by Eqs. (26), (27). We see that the coefficient d 1 is bigger than d 2 by about a factor of two. At low invariant mass close to threshold they increase, reflecting the amplitude t DD,DD which is contained in the coefficients.
The most important result is shown in Fig. 5 where we show the results of dΓ/dM inv for the D 0D0 distribution. We see a concentration of the strength around threshold with a peak around 3735 MeV. In order to see that this structure is tied to the resonance below threshold we plot in Fig. 6 the phase space for ψ(3770) → γD 0D0 substituting |t| 2 of Eq. (25) by a constant. We also show the phase space for ψ(3770) → γD + D − keeping the physical masses for the D mesons. We can observe that the shape of the phase space distribution is drastically different from that predicted in the presence of a DD bound state. The phase space peaks around 3742 MeV instead of 3735 MeV for the distribution with the DD bound state. The shapes of the fall down of the two distributions are also different, the one with the DD bound state falling as a concave curve and the phase space as a convex one. Finally we show in Fig. 7 the mass distribution for ψ(3770) → γD + D − . The shape is quite different than the one for ψ(3770) → γD 0D0 and the reason is the contribution of the tree level, which is absent for D 0D0 production. It is clear that ψ(3770) → γD 0D0 is, thus, the reaction that better shows the presence of the DD bound state. Although not relevant for the present discussion, but we can see the infrared divergence associated to the limit k → 0 for the photon momentum, when the D propagator in the tree level amplitude becomes on shell.  The strengths of the dΓ/dM inv distribution obtained fall well within present capacity at BESIII. Indeed the integrated luminosity in the 3770 MeV region is about 5 fb −1 per year, which, together with the e + e − production cross section of the resonance of 7.2 nb [41], leads to about 3.5 × 10 7 accumulated events per year. The expected final data will have an accumulated 20 fb −1 luminosity [42] and with the planned BES upgrade 15 fb −1 per year (1.1×10 8 events per year) 1 . Our differential widths dΓ/dM inv of the order of 10 −6 are well measurable with these ψ(3770) production rates 2 .

IV. CONCLUSIONS
We have made a study of the ψ(3770) → γDD decay, looking at the D +D− and D 0D0 mass distributions in dΓ/dM inv close to theshold. We saw that the production of D 0D0 is particularly suited to study the dynamics of the DD interaction because the tree level contribution is zero and the process goes with a loop mechanism that involves the D + D − → D 0D0 scattering amplitude. We have used the results of a theory that predicts a DD bound state and this has as a consequence that the D 0D0 mass distribution accumulates close to threshold and diverts drastically from a phase space distribution. The rates that we obtain for the mass distribution are perfectly reachable with the present BESIII facility and we encourage the performance of the experiment that could shed light on the issue of this possible DD bound state, and in any case would provide information on the DD → DD interaction.