The cross section of the process $e^+e^- \to p\bar{p}$ in the vicinity of charmonium $\psi(3770)$ including three-gluon and $D$-meson loop contributions

The total cross section of the process $e^+e^- \to p\bar{p}$ is calculated within the energy range close to the mass of $\psi(3770)$ charmonium state. It was shown that the main contribution to this cross section comes from three gluon mechanism which is also responsible of the large phase with respect to the Born amplitude. This phase provides the characteristic dip behaviour of the cross section in contrast to the usual Breit-Wigner peak shape. OZI-allowed mechanism with D-mesons in the intermediate state was also estimated and gives relatively small contribution to the cross section and to the phase.


Introduction
The electron-positron colliders allow to study different physical processes with final particles in pure J P C = 1 −− state. For example, single charmonium production gives the possibility to look at the relativistic bound state of cquarks, elaborate and test new ideas of confinement description and to search for possible exotic admixtures in the wave function of the state.
One of the intriguing states of charmonium is ψ(3770) which was studied by many collaborations (for example, by KEDR-VEPP-4M [1,2], CLEO [3] and more recently by BES III [4,5]). Especially the latter one [5] made precise measurement of the total cross section at the specific kinematics of ψ(3770) mass and showed that instead of a Breit-Wigner peak one can see a dip in the dependence of the total invariant mass squared s (see Fig. 2 in [5]). This observation immediately lead to the conclusion that there must be some mechanism which generates large relative phase φ between resonant and continuum terms in the amplitude. And indeed soon some possible explanation of this phase was suggested in Ref. [6] where is was attributed to the OZI-violated three gluon mechanism of charmonium ψ(3770) transition into final proton-antiproton state. This type of mechanism was already shown to give large phase for another charmonium χ 2 ( 3 P 2 ) in Ref. [7]: it was shown that two gluon mechanism can produce the phase of order of φ ∼ −90 • .
Here we want to continue the elaboration of three gluon model from Ref. [6] and to refine D-meson loop calculation. We compare our estimations of the total cross section with the precise experimental data of ψ(3770) production obtained at BES III [5].
The paper is organized in the following manner: in Section 2 the total cross section of the process e + e − → pp in Born approximation (that is so called continuum part) is obtained and some kinematical notations are introduced; Section 3 is the main calculation part of the paper, the extra resonant contribution to the amplitude is added and the formalism of how to take it into account in the cross section is developed. There are two subsections: Subsection 3.1 describes the details of calculation of possible OZI-allowed mechanism with Dmeson loop, while Subsection 3.2 briefly reminds some key steps of three gluon mechanism calculation from [6] (with some minor corrections of typos and errors in [6]); Section 4 gives some numerical estimations and comparison of our calculation with experimental data from BES III [5]; Section 5 summarizes our results and proposes some possible extension of this work in the future.

Born approximation
We consider the process of proton-antiproton pair creation in an electronpositron collisions: where quantities in parenthesis are the 4-momenta of the corresponding particles. In Born approximation (see Fig. 1(a)) the electron-positron pair annihilates into virtual photon, which then produces proton-antiproton pair. The amplitude M B corresponding to this process has the form: where s = q 2 = (q + + q − ) 2 = (p + + p − ) 2 is the total invariant mass squared of the lepton pair (q is the momenta of intermediate photon). The quantities J eē→γ µ and J γ→pp µ from (2) are lepton and proton electromagnetic currents: where e is the modulus of electron charge e = √ 4πα and α ≈ 1/137 is the fine structure constant [8]. The proton current is described in terms of the proton electromagnetic form factors: where we use the notationâ ≡ a µ γ µ . Here M p is the mass of proton and functions F 1,2 (q 2 ) are the proton electromagnetic form factors normalized as F 1 (0) = 1 and F 2 (0) = µ p − 1, where µ p is the proton anomalous magnetic moment. In paper [6] we used the result of the paper [9] and assumed that energy range under our consideration, i.e. √ s ∼ 3 − 4 GeV, is not far from pp production threshold and thus we can neglect by F 2 and use the point-like proton approximation, putting F 1 = 1. However recent analysis [10] showed that the situation is not so evident and point-like proton approximation is not valid at this energy range. That is why we use the effective form factor G(s): (as it was done in [5] following to the results of [11]) which is obtained in the assumption that electric G E and magnetic G M form factors of the proton are equal: |G E | = |G M |. That assumption leads to F 1 (s) = G(s) and F 2 (s) = 0. The formfactor (6) is the pQCD inspired form factor [12,13], which is in a fair agreement with the cross section of the process (1) in the energy range √ s from 2 GeV to 3.07 GeV measured by BES [14]. In equation (6) the quantity Λ = 300 MeV is the QCD scale and C is a free parameter fitted in [5] to be equal to C = (62.6 ± 4.1) GeV 4 . This fit also agrees with more recent result [15] where this constant was fitted to be C = 72 GeV 4 while Λ = 520 MeV. During our numerical estimations we use the BES values of C and Λ.
Using the amplitude M B from (2) one can write down the cross section in the standard way: where the summation of the amplitude square |M B | 2 runs over all possible initial and final particles spin states. We also systematically neglect the mass of the electron m e in this paper. The phase volume of final particles dΦ 2 has the form: where φ p and θ p are the azimuthal and polar angles of the final proton in the center-of-mass reference frame (center of mass system, c.m.s), i.e. θ p is the angle between 3-momenta of the initial electron q − and the final proton p + (see Fig. 2) and |p| ≡ |p + | = |p − | = √ sβ/2 is the modulus of 3-momenta of the final proton or antiproton. The quantity β = 1 − 4M 2 p /s is the final proton velocity. Using the explicit form of M B from (2) and integrating over the final particles phase space (8) from (7) one obtains: The total cross section in Born approximation then reads as: which is in a good agreement with the experimental data on a wide energy range far from the resonance ψ(3770), see Fig. 3.

The quarkonium ψ(3770) intermediate state
As one can see in Fig. 3 the total cross section including only the electromagnetic mechanism (10) fails to describe this delicate behaviour in the vicinity of the charmonium resonance ψ(3770). Obviously in this region one should take into account the additional contribution to the amplitude which appears from the diagram with ψ(3770) in the intermediate state (see Fig. 1(b)) and is enhanced by Breit-Wigner propagator. This leads the total amplitude M of the process (1) to become the sum of two terms: BES III [5] BaBar [16] BaBar [17] Born (10) where M B is the Born amplitude from (2) (see Fig. 1(a)) and M ψ takes into account this second mechanism (see Fig. 1(b)): where M ψ and Γ ψ are the mass and the total decay width of ψ(3770) resonance and J eē→ψ µ and J ψ→pp µ are the currents which describe the transition of lepton pair into ψ(3770) resonance and the transition of the ψ(3770) resonance into proton-antiproton pair correspondingly. We notice that the second term in parenthesis in (12) does not contribute since the currents have to be conserved: q µ J eē→ψ µ = q µ J ψ→pp µ = 0. Next step is to assume that J eē→ψ µ has the same structure as J eē→γ µ from (3), i.e.
where the constant g e = F ψ→pp is the value of the form factor of the vertex ψ → pp at the ψ(3770) mass-shell (here we follow the same approximation as in the Born case and assume that F ψ→pp 2 (M 2 ψ ) = 0). This constant is defined via ψ → e + e − decay width Γ ψ→e + e − = 261 eV [8] which gives the following value: We neglect a possible imaginary part of vertex eē → ψ since it was shown in [7] that it is small, less then 10 % of the real part.
The strategy for calculating the contribution M ψ to the cross section is the following. If we introduce the relative phase φ between Born contribution M B and the additional contribution M ψ then we can write the cross section as: Knowing the Born cross section σ B from (10) and the interference contribution σ int with the phase φ one can calculate the total cross section including both contributions using (15) and evaluating σ ψ in the following manner: Thus we need to evaluate only the interference of the resonant amplitude M ψ with the Born amplitude M B which has the standard form: Our goal is to get the contribution to the total cross section thus we integrate over the final particles phase space: where s is the summation over the spin states of initial particles and s ′ is the summation over the final particles spin states. Next we use invariant integration trick: Applying the conservation of the currents J eē→γ µ and J eē→ψ ν , recalling that and using the phase volume from (8) we get the following simplified expression of the interference contribution to the total cross section: Thus we can present the interference contribution to the total cross section in the form: where S i (s) contains all the dynamics of the transformation of charmonium into proton-antiproton pair and has the following explicit form: The subscript index i above denotes the type of mechanism of this transformation. Since the mass of ψ(3770) is higher than the threshold of D-meson pair production it is natural to expect that the D-meson loop will be the main mechanism in this reaction (see Fig. 4). However this mechanism appears to be small since the mass of ψ(3770) exceeds the D-meson pair production threshold slightly (one can see that (M ψ − 2M D ) /M ψ ≈ 1%). Thus below, in Section 3.2, we consider also the OZI-violated three gluon mechanism (see Fig. 5) which appears to give the dominant contribution.

D-meson loop mechanism
The D-meson loop mechanism can be illustrated by the diagram in Fig. 4. It is known that the ψ(3770) structure is somewhat complicated. For example it leads to a specific peak shape of the cross section of the process e + e − → DD in the vicinity of ψ(3770) state mass. One of the strong conjecture to interpret this fact is to suppose the presence of a ψ(2S) resonance with mass close to the ψ(3770) [18]. We however do not dive into this complicated consideration since D-meson loop contribution appears to be small in our case and one can limit oneself with single two quark state with spin 1. We introduce the following parameterizations of vertices in Fig. 4: where e µ is the polarization vector of charmonium ψ(3770). In the calculation we do not need to know the complete functional dependence of functions G ψDD (q 2 , k 2 1 , k 2 2 ) and G ΛDP k 2 , p 2 . We discuss this dependence in this Section below.
Let us write down the D-meson loop contribution to the amplitude from the diagram in Fig. 4: where M D and M Λ are masses of D-meson and Λ c -hyperon correspondingly. Comparing this amplitude with the general form (12) one can extract the current: and insert it into (23). This gives the following contribution of the D-meson loop to the cross section which is expressed is terms of the quantity S D from (22): where SpD(s, k 2 ) is the trace of γ-matrices over the baryon line: Now we need to evaluate the quantity Z D (s) from (29). In order to do this we use Cutkosky rule [19] for D-meson propagators: and obtain the imaginary part of this quantity: We notice that one can utilize these two δ-functions in (32) to significantly simplify the evaluation of Im Z D . For example, one can see that usage of δfunctions leads to the replacements: which we already used in (30). Performing the loop integrations one gets: with the integration over cosine of polar angle C k = cos θ k is evaluated numerically below. Appendix A details the derivation of eq. (35) and the definitions of the quantities C k , k (i) and D 1 . Now we consider the explicit expression of form factors which we need to evaluate Im Z D from (35). First we see that for ψ → DD vertex we need only dependence over charmonium virtuality q 2 = s, since D-meson legs are on-massshell. We can start with the normalization of function G ψDD (s, M 2 D , M 2 D ) to the decay of charmonium ψ(3770) into DD final state. Calculating this decay width one gets: where β D = 1 − 4M 2 D /M 2 ψ is the D-meson velocity in this decay. The numerical value of g ψDD in (36) is obtained by using the experimental value of the width of charmonium decay into charged (D + D − ) and neutral (D 0D0 ) mesons Γ ψ→DD = 25 MeV [8]. This let us to take into account all these types of D-mesons in the loop in Fig. 4.
The dynamical dependence of function G ψDD (s, M 2 D , M 2 D ) we propose following to [20] in the form: where constant C ψDD can be found from the normalization (36) and we get Finally we use the following explicit form of function: where scale Λ D we fix on the characteristic value of the reaction Λ D = 2M D . Next we consider function G ΛDP k 2 , p 2 from (35). And again the only dependence left after the application of Cutkosky rule is the off-mass-shellness of Λ c -hyperon in the scattering regime, since here k 2 = −k 2 (i) < 0. The dependence over the virtuality of the D-meson has disappeared, D-mesons are on mass shell. But still we need to take into account the remnants of this dependence. Following to [6] we use this form for ΛDP -vertex with the off-mass-shell Dmeson which was established in [21,22]: where f D ≈ 180 MeV. For quark masses the following values are used: m u ≈ 280 MeV and m c = 1.27 GeV [8]. The constant g DN Λ ≈ 6.74 was estimated in Ψ (cc) Figure 5: Three gluon mechanism.
[22] in scattering regime, i.e. for p 2 < 0. Thus in our calculation in eq. (35) we use the following expression: which do not take the effects of the Λ c -hyperon off-mass-shellness. We expect that these effects are not very important. Now we're able to calculate the imaginary part of Z D following to eq. (35). The real part of this quantity is restored by using the dispersion relation (see Appendix B).

Three gluon mechanism
The three gluon mechanism was considered in details in [6] and thus we just briefly recall a few steps of this calculation in order to correct some misprints and minor mistakes in formulas of [6] which do not affect the conclusions of the paper.
The three gluon mechanism presented in Fig. 5 gives the following contribution to the quantity S 3g from (22) (which coincide with equations (16) and (17) from [6]): where the quantity Sp3g is the product of traces over proton and c-quark lines: witĥ

× 2 variants
where the permutations over gluon vertices are performed in the gray block in Fig. 5. First we notice that there is a misprint in second proton propagator denominator in eq. (18) in [6] which is corrected here in (42). This misprint leads to the following mistakes in angular integrals in eq. (21) and in the Appendix B of [6]. The angular integration must be performed over dc 1 dc 3 instead of dc 1 dc 2 .
Next we also note that the quantity φ in α 3g (s) in (41) is similar to the quantity R from eq. (17) of [6]. The difference comes from the calculations of ψ → 3g decay which we use for the normalization of the decay constant φ. In paper [6] (see Appendix A there) the factor 2/3 was missed in the amplitude of ψ → 3g decay (see eq. (A.7) there). This factor comes from missed 2 in the right-hand-side of the sum over spins in the expression below eq. (A.6) in [6] and from the complete missing of symmetrization factor 1/ √ 3 in the color wave function of charmonium (q 1 q 1 +q 2 q 2 +q 3 q 3 ) / √ 3 which is necessary for correct normalization of this wave function (see, for example, eq. (5) in [23]). Thus, after correction of this missing factor one gets: A correction should be also applied to the color factor g col in (41). This factor was evaluated in eq. (11) of Ref. [6], taking into account the summation over gluon and quark colors. But the probability of gluon connection with the proton is evaluated by simple multiplication of a factor of 3 which comes from 3 quarks in the proton. Taking into account all the possibilities (see Fig. 6) we get a factor of 27: The most important correction concerns the final proton-antiproton state. The three gluons obtained from ψ(3770) decay produce three quark-antiquark pairs and it is implicitly assumed that they form the proton-antiproton final state. The details of this process is not touched in Ref. [6] and it is mostly not important for the relative phase. However here we want to reproduce the absolute value of the cross section and thus we need to implement this mechanism somehow. In general it is the mechanism of transition of three gluons (with total angular momentum equal to 1) into final proton-antiproton pair. We suggest that this mechanism has much in common with proton-antiproton pair production from the photon, i.e. the electromagnetic vertex γ * → pp in time-like region. Thus we insert into (41) the extra form factor, similar to (6), but with some other value of the parameter C ψ : We find the value of this parameter C ψ in Section 4. The parameter Λ is still the QCD scale parameter since it takes into account the running of α s coupling at microscopic scale. The strategy of the calculation of the quantity Z 3g (s) is the same as for the D-meson loop contribution in the previous section: first we calculate its imaginary part Im Z 3g (s) by the use of Cutkosky rule for gluon propagators: then we restore the real part Re Z 3g (s) by using dispersion relation technic. This strategy is the same as in Ref. [6]. See all the details of these calculation there.

Numerical results
First we present the main ingredient of our calculation, the quantities Z D (s) from (29) and Z 3g (s) from (42) as a function of total energy √ s in the range starting from the threshold of the reaction s = 4M 2 p up to 4.5 GeV, see Fig. 7(a) and Fig. 7(b). The position of the ψ(3770) resonance is marked by a vertical dashed line. One can see that the quantity Z 3g is much bigger than Z D and both have large real and imaginary parts in that region. The imaginary part of Z D below DD threshold (i.e. at s < 4M 2 D ) is zero and grows slightly up from this point with energy giving small contribution to the observables.
Next we plot the cross section itself and compare it with the data from the BESIII collaboration presented in [5]. In Fig. 8 we present experimental data of : The comparison with the fit of BESIII [5]. The solid points correspond to the Solution 1 in Table 2 in [5], while the hollow points correspond to the Solution 2.
BESIII scan of total cross section of the reaction e + e − → pp around the mass of ψ(3770) resonance and with the specific precise measurement of the point at exactly the mass of ψ(3770). We plot the Born cross section according to (10) with the proton electromagnetic form factor G(s) from (6) as a gray line in Fig. 8. We see that it overestimates the precise point at ψ(3770) resonance mass. Adding the D-meson loop contribution (dashed line) does not reproduce the desired ψ(3770) resonance point. But if one replaces the D-meson loop contribution with the three-gluon mechanism (dotted line) then the tendency of the curve becomes rather similar to the data. And if we take into account all three contributions (solid line) then we see that our total curve reproduce the precise ψ(3770) resonance point and follows the side shoulders, i.e. the experimental points that seem to be higher below the ψ(3770) and lower above the ψ(3770) resonance. In fact this coincidence of our total curve with the precise data point is the result of fitting of the constant C ψ from (46) to the best measured experimental point at ψ(3770) mass. Since the dependence of the total cross section of the constant C ψ is quadratic we have two solutions for it. Using the data point errorbars we can estimate that the range of the values of this constant should be 33.7 GeV 4 < C ψ < 59.1 GeV 4 with the best fit to the central point at two values: Fit 2: We find this values to be rather close to the one from the electromagnetic vertex case, see the value of parameter C below eq. (6). In all our numerical estimations we use the value of Fit 1 from (48). The usage of Fit 2 gives the slight change in the plot shoulders but keeps the central point untouched. Now we can compare our calculation with the fit of the expression for the total cross section (1) in [5]: where the quantities σ ψ and φ ψ are considered as parameters to be fitted on the data. The fitting procedure gave two solutions (see Table 2 in [5]): We plot these two solutions in Fig. 9 at the ψ(3770) resonance mass (in Fig. 9(b) we shifted the points slightly aside from ψ(3770) resonance mass to make the points seen separately). Assuming the form (50) we extract the quantities σ ψ and φ ψ from our result for total cross section (15) and plot them in Fig. 9. One can see that both phases are consistent with our curve (see Fig. 9(b)) while only one value for σ ψ agrees with our result (see Fig. 9(a)). Namely the solution 1 from Table 2 in [5] agrees with our numbers: Finally we note that our calculation of three gluon mechanism contains the QCD coupling constant α s at charmonium scale (i.e. at s ∼ M 2 c ) in rather high degree (see eq. (41)) and thus it is very sensitive to its value. We use the value α s (M c ) = 0.28 which is expected by the QCD evolution of α s from the b-quark scale to the c-quark scale. We should note that this value differs from the one for J/ψ charmonium for which one must use much smaller value of parameter α s (M c ) = 0.19 [23].

Conclusion
We considered the process of electron-positron annihilation into protonantiproton pair in the vicinity of charmonium ψ(3770) resonance. Besides the Born mechanism, which is the pure QED, there are two contributions related with intermediate charmonium ψ(3770) state. One of them is the D-meson loop and the other is three gluon mechanism.
We showed that D-meson loop mechanism, being the most probable candidate to describe this process, fails to reproduce the value and the shape of the cross section.
It has been shown that the main contribution comes from three gluon mechanism which can reproduce the position of the precise experimental point at the ψ(3770) resonance mass and also grasps the shape of the curve shoulders around it.
Having all the calculation in the hands we are able to select one of the two fit solutions obtained in [5] and to provide a solid basis for the consideration of similar processes with binary final states and with charmonium in the intermediate state. This will be the subject of our future works. Appendix A. Integration over dk in (32) using δ-functions In this section we show how one can obtain (35) from (32) utilizing the benefits of δ-functions in the integrand. In order to do this we consider the following expression: where C k ≡ cos θ k , θ k and φ k are the polar and azimuthal angles of loop momentum k which are measured from the direction of final proton momentum p + (see Fig. A.10). In this section we use the notation k ≡ |k|. The function F (k 0 , k) contains all non-trivial part of the integrand from (32): First we transform the second δ-function in (A.1) to the form: where the quantity k (1) 0 is one of the poles of the argument of δ-function with respect to the variable k 0 : The remaining δ-function can be transformed into the following form: where k (i) are the poles of the argument of δ-function with respect to the variable k: We note that if k = k (1,2) then D = s and thus k 0 = k (1) 0 = 0. The requirement that D 1 > 0 leads to the restriction that C k > C (1) k or C k < C Finally, we notice that integration over d cos θ p and dφ k is trivial and gives 4π: where the threshold condition (s > 4M 2 D ) comes from the fact that the integral is not zero only if C (1) k < 1. Comparing (A.1) and (A.9) we prove that (35) follows from (32).

Appendix B. Dispersion relations
In this section we derive dispersion relations to restore the real parts of the quantities Z D from (29) and Z 3g from (42). We apply dispersion relations over the variable s with subtraction at the point s = 0: Re Z i (s) = Re Z i (0) + s π P ∞ smin ds 1 s 1 (s 1 − s) Im Z i (s 1 ), (B.1) where s min is the minimal threshold at which Im Z i (s 1 ) becomes non-zero. The substraction constant Re Z i (0) vanishes (i.e. Re Z D (0) = 0) since there are no open charm in the proton and thus in the Compton limit the vertex ψ → pp must vanish. Next we substitute the variable β: Substituting these replacements into (B.1) one gets: where β min = 1 − 4M 2 p /s min . In order to improve the numerical stability of this integral we transform it: we add and subtract the regular part of the numerator in the point β: Re Z i (β) = P π 1 βmin 2β 1 dβ 1 β 2 1 − β 2 (Im Z i (β 1 ) − Im Z i (β) + Im Z i (β)) ,

(B.3)
and evaluate the part with Im Z i (β) in the numerator: Re Z i (β) = 1 π Im Z i (β) log 1 − β 2 β 2 min − β 2 + The integral in braces is already regular at the point β 1 = β and can be evaluated with ease by any standard numerical method of integration. We note that for the D-meson loop contribution the imaginary part Im Z D (β) from (35) is non-zero above threshold (s > 4M 2 D ), thus the lower limit of integration in (B.4) is β min = 1 − M 2 p /M 2 D . For the three gluon contribution the threshold for the imaginary part Im Z 3g (β) coincides with the threshold of the reaction, i.e. s min = 4M 2 p and thus lower limit of integration is β min = 0.