Charmonium decays: J/psi -->e+e- and eta_c -->gamma gamma

We compute the O(alpha_s^2) correction to the decay rate eta_c -->gamma gamma and discuss its implications for precision quarkonium physics. We study the suitability of the ratio Gamma(J/psi -->e+e-)/Gamma(eta_c -->gamma gamma), in which the non-perturbative or soft effects cancel at O(alpha_s^{0,1}), for extracting fundamental parameters of QCD at low energies. We show that the QCD-based theory of charmonia is not capable of predicting this ratio with any degree of confidence.

The physics of cc mesons (charmonium) is a mature field with a long history.Discovery of the J/ψ resonance at SLAC and Brookhaven in the autumn of 1974 is often called "November revolution", to emphasize its importance for the development of QCD and of the Standard Model.Already the early theoretical papers on the subject interpreted the observed narrow resonance as the non-relativistic bound state of cc quarks and thus initiated studies of heavy quarks in QCD [1].
During almost thirty years since the discovery of the first charmonium, both experimental and theoretical studies of these mesons have been successfully pushed forward.The spectrum, lifetimes and branching ratios have been precisely measured.Further progress is expected at the planned dedicated facility CLEO-c [2].On the theoretical side, various attempts have been made to improve the description of these hadrons.In particular a lot of effort went into determining how well the cc bound states can be described if one starts directly from the QCD Lagrangian.
An important recent development has been the introduction of effective field theory techniques for describing hadrons consisting of two non-relativistic heavy quarks [3,4].This effective field theory, the Nonrelativistic Quantum Chromodynamics (NRQCD), connects the original QCD and the new NRQCD Lagrangian which takes full advantage of the fact that the quarks in the cc bound state are non-relativistic.As usual in effective field theories, the two Lagrangians are matched perturbatively at energy scales around the charm quark mass m c ≡ m ∼ 1.7 GeV.
In recent years the effective field theory approach to non-relativistic bound states has been extended further.
It has been noticed that two additional scales, the heavy quark momentum mv and the heavy quark binding energy mv 2 , exist in quarkonia and, for sufficiently heavy quarks, permit perturbative treatment.Unfortunately, this is not quite possible for charmonium because the c quark mass is too small and therefore mv ∼ mv 2 ∼ Λ QCD .However, we still have Λ QCD ≪ m and therefore there is a chance that integrating out hard (k ∼ m) modes and matching QCD at NRQCD perturbatively is a sensible thing to do.If the soft effects are universal, they cancel in ratios of various observables and a clean perturbative QCD prediction emerges.
For various observables, this approach has been taken at order O(α s ) and the common perception is that it works rather well.Let us consider the simplest decays of the ground state charmonia, J/ψ → e + e − and η c → γγ.To order O(α s ) the decay rates can be written as where ψ(r) is the charmonium wave function and x ψ,η are numbers which can be determined by perturbative matching of the QCD and NRQCD Lagrangians.Taking the ratio of these two decay widths, one obtains a prediction that is free from any non-perturbative uncertainties, and can be either compared to the data provided α s (m) is known or used to extract the value of α s .For example, the CLEO collaboration has recently [5] determined the strong coupling constant α s (m c ) from the branching ratio of η c → γγ assuming that the η c total decay width is saturated by the two-gluon channel.Since the ratio (2) is independent of the wave functions, its comparison with measurements has often been considered a solid test of perturbative QCD.
In the past few years the development of techniques for non-relativistic effective field theories has undergone an important transition and we can now study the next order in the strong coupling constant expansion.Interestingly, in doing so one encounters new conceptual difficulties.The most important problem is that the hard Wilson coefficient of the operator responsible for the decay in the leading order becomes infrared divergent at two loops, which implies that the wave function at the origin becomes scale dependent.If this scale dependence were the same for spin triplet J/ψ and singlet η c states, this would not pose a difficulty since it would cancel in the ratio.However, this is not the case and the divergent parts of the Wilson coefficients are spin-dependent.This immediately implies that with O(α 2 s ) accuracy the wave functions at the origin of J/ψ and η c become different and therefore the ratio of the corresponding decay widths is sensitive to some soft-scale effects.Since it is rather difficult to compute these effects accurately, the QCD-based prediction for the ratio of decay widths Γ(J/ψ → e + e − )/Γ(η c → γγ) becomes much less precise if one computes higher order corrections in O(α s ) -a somewhat paradoxical situation.This is the principal message we would like to get across in this Letter.
The remaining part of this Letter is organized as follows.We first consider hard renormalization factors of the non-relativistic operators responsible for the decays η c → γγ and J/ψ → e + e − .We then show how the soft effects are taken into account in our calculation and derive our final result for the ratio of the decay rates of J/ψ → e + e − and η c → γγ.
The hard renormalization factor for the spin-singlet decay operator has been evaluated very recently [6], extending an earlier QED result obtained for the parapositronium decay [7,8].We did not include the diagrams shown in Fig. 1 in our final result.We have checked with a rough approximation that this finite and gauge-invariant subset contributes only insignificantly.
Our result for the hard renormalization factor for the singlet decay operator is (we use dimensional regulariza- where α s = α s (m) and In the above equation µ is the factorization scale that separates relativistic and non-relativistic momenta in the NRQCD framework.
The matching coefficient for the vector current, relevant for the decay J/ψ → e + e − , can be found in Eqs.(10)(11)(12)(13)(14)(15) of Ref. [9] (see also [10]).Divergences in the matching coefficients of the two currents are different.In the MS renormalization scheme the ratio of the decay widths is (we use Obviously, the wave functions at the origin become renormalization scheme and factorization scale dependent and eq.( 5) indicates that they must be different for J/ψ and η c .Since they can no longer be eliminated in the ratio of the decay widths, one looses predictive power.Note also that eq. ( 5) involves the MS renormalized wave functions, which cannot be directly determined even on the lattice.Let us stress that to reach these conclusions one only has to integrate out relativistic degrees of freedom.
No assumption about the dynamics of the bound state (besides its quantum numbers) is necessary.
To cancel the divergences of the Wilson coefficients we have to calculate certain soft effects.For charmonium, which is not a Coulombic bound state, a completely "honest" calculation is not possible.However, we will show that with plausible assumptions one can determine the soft contribution to the ratio of two decay widths using experimental data on e + e − → cc.Two types of relativistic corrections have to be considered: to the amplitude and to the wave function at the origin.It turns out that the former are the same for J/ψ → e + e − and η c → γγ and cancel in the ratio (for this reason they were not included in eq. ( 5)).As for the latter, only spin-dependent effects can survive in the ratio.We assume that the only operator in the non-relativistic Hamiltonian responsible for the hyperfine splitting is where the Pauli matrices σ and Σ act respectively on spins of the charm quark and of the antiquark.This operator follows from a one-gluon exchange diagram and its QED analog is the hyperfine splitting operator in the Breit Hamiltonian.With this operator we compute the wave functions at the origin, where G(0, 0) is the reduced Green's function of the nonrelativistic cc state, computed at the ground state energy E 1 (at this level we neglect the difference between J/ψ and η c masses), and the prime means that the sum does not include the ground state J/ψ (n = 1).The dots in eq. ( 7) indicate that spin-independent corrections are also present but they drop out in the ratio of the decay widths.The Green's function at the origin is divergent and it is precisely the divergence needed to cancel that in the ratio of the two hard Wilson coefficients.To illustrate this, we first consider an academic example of ultra-heavy quarkonia, where binding effects can be computed in the Coulomb approximation.In this case, the Green's function at the origin G(0, 0) can be extracted from Ref. [11] and reads We now substitute this result in to eq. ( 7) and obtain the ratio of the decay widths of the spin-triplet and spinsinglet ultra-heavy quarkonia, where a s = α s (m Q ).The coefficient of the a 2 s ln a s term, 4C 2 F /3 ≃ 2.37, agrees with [12].We now turn to the more difficult case of charmonium, where one can not use the Coulomb approximation for the low energy dynamics.We will still use the hyperfine splitting operator, eq. ( 6), to describe the hyperfine splitting; for this reason eqs.( 7) are still valid.The challenge is to compute the Green's function G(0, 0) in dimensional regularization without invoking the Coulomb approximation.
We will perform such computation in two different ways dealing either with the observed cross section e + e − → cc or with a simple potential model constructed to describe charmonium.It will be clear that both calculations could be improved.Here, our primary goals are to demonstrate how they can be carried out using simple approximations and to give reasonable estimates of the second order corrections to the ratio of the two decay widths.
We first describe the calculation which utilizes experimental data on e + e − → cc.To this end, we separate the Green's function into a contribution of charmonium resonances and that of the continuum, The former is finite and can be computed using available data on the spectrum and e + e − decay widths of the spin 1 resonances, where we have employed the mass and width information on the first six ψ resonances [13].
The continuum contribution is, on the other hand, divergent.To determine this divergence using as little input information as possible, we proceed in the following way.We consider continuum contribution as a function of E, and take the derivative with respect to E. We then solve the resulting differential equation and obtain We can further use the relation between R c = σ(e + e − → cc)/σ(e + e − → µ + µ − ) and the Green's functions G(0, 0) to relate eq. ( 14) to experimental data.We obtain where the function K(E) is given by The divergence in G(0, 0) now resides in G cont Ei and the integral in eq. ( 15) is finite.For the initial condition G cont Ei we can choose a "deep Euclidean" point E i → −∞, where perturbative calculations are justified and where G cont Ei can be determined with, in principle, arbitrary precision.We therefore see that the divergent part of the Green's function can be extracted from perturbative calculations and the finite part can be obtained from R c (other options are to compute the finite part using a potential model for quarkonium or NRQCD on the lattice).This separation solves the problem in principle and provides a way to determine the MS charmonium wave functions without assuming that the bound state is Coulombic; the only true model-dependence remaining in this calculation is the form of the operator responsible for the hyperfine splitting, eq. ( 6).
The data on R c is not quite precise yet.Its high energy asymptotics (in the non-relativistic sense) is fixed since the dependence on the initial energy E i in eq. ( 15) should cancel.We therefore write For R npt c (E) we choose: with E 0 = √ s 0 − 2m and √ s 0 = 4 GeV.This ansatz is motivated by the data on R c in [14] where one sees that there is no need for a large second term in eq. ( 17) below 4 GeV.On the other hand we do need this term at higher energies, since otherwise "perturbative" and "non-perturbative" expressions do not match.
We obtain The soft contributions to the decay width ratio are obtained by employing eqs.( 12) and (19) in eq.(7).Before combining them with the hard contributions in (5), we note that the BLM effects [15] were computed in Ref. [16] for the rates of J/ψ → e + e − and η c → γγ.These corrections turn out to have different signs and are enhanced in the ratio.For this reason, we decided to eliminate them by choosing different scales for the strong coupling constant in the Wilson coefficients for J/ψ → e + e − and η c → γγ.We then find where a S=0 = α s (1.95 m) and a S=1 = α s (0.63 m) and the charm mass should be expressed in GeV.The scale of the coupling constant in the second order correction is not specified; we will use α s = 0.3 for the estimates.The inverse power of α s in the square brackets arises because we have used experimental data to compute the contribution of the resonances to G(0, 0) and also used the energy of the ground state E 1 = M J/ψ − 2m to estimate the continuum contribution.This spoils the homogeneity in α s .
Employing, for the sake of illustration, α s = 0.3, α s (0.63 m) = 0.35 and α s (1.95 m) = 0.26 in eq. ( 20), one finds where the three terms in the brackets are the tree level, the O(α s ) and the O(α 2 s ) corrections, respectively.As shown in Fig. 2, the second order correction depends strongly on the value of the charm quark mass.It is large and very sensitive to the form of R npt c assumed in eq.(17).For example, if we use the Coulomb approximation to estimate the value of the wave function in the continuum, we obtain a negative result for the second order correction.Within our model for the continuum, we estimate As we mentioned earlier, another possibility to obtain the Green's function necessary to compute the ratio of two decay widths is to use either QCD on the lattice or potential models.Here we would like to illustrate this point by considering a simple potential model.This will allow us to check that the estimate, eq. ( 22), is reasonable.Although the potential model below is really simple, the calculation can be repeated with more sophisticated potentials, provided that, at short distances, these potentials match the QCD analog of the Coulomb potential.
To describe charmonium, we will use the Schrödinger equation with the potential where b = 0.18 GeV 2 and V 0 is adjusted to give the correct mass of the 1S state (J/ψ) for given values of m and α s .
To compute the Green's function G(0, 0) we solve the Schrödinger equation following the treatment in [17] and obtain the following representation for the full (including the ground state) Green's function: where the function B(E) is derived from the large r limit of the ratio of two solutions of the Schrödinger equation with prescribed behavior at the origin (see [17]).In order to obtain the Green's function in dimensional regularization, which is needed for our purposes, we write We then take the limit E 1 → −∞, perturbatively compute the Green's function G E1 and derive lim where c is and γ E is the Euler constant.
In order to compute the reduced Green's function G(0, 0) eq. ( 8), we find the first eigenvalue numerically and remove the pole from B(E).
We have computed the ratio of two decay widths for different values of α s and the mass of the charm quark.We obtain We see that the result of the potential model calculation is relatively close to the result of eq. ( 22).The advantage of the potential model calculation is its fair stability against variations in α s and m.We take the result in eq. ( 28) as our final estimate.
In spite of the fact that the model leading to eq. ( 28) is quite simple, we believe that eq. ( 28) is important in that it clearly shows the magnitude of second order QCD corrections one might expect for such observables.
We conclude that eqs.(20,21,28), the principal results of this Letter, illustrate an unexpected problem in the theory of heavy quarkonia at the two-loop level.In recent years we have learned how to integrate out relativistic degrees of freedom efficiently and it seemed as if we could improve the accuracy of our predictions.This turns out not to be the case.The reason is that at O(α 2 s ) the soft and relativistic effects do not decouple completely, as it happens at O(α s ), and therefore, in general, one cannot avoid non-perturbative effects by taking ratios of different observables.We have shown how, in principle, the soft contribution can be estimated using experimental data or potential models.
With the QCD corrections as big as in eq. ( 28), the determination of α s (m c ) from charmonia decay rates, as e.g. in Ref. [5], does not look trustworthy, regardless of the fact that the numerical values of α s turn out to be in a theoretically sensible range.On the other hand, it is interesting to point out that the ratio of the decay rates of η c → γγ to η c → gg, actually used in Ref. [5] for determination of α s , is free from the soft effects we discussed in this Letter, since it refers to the same initial state.It would therefore be interesting to compute second order QCD corrections to this ratio since in this case the hard corrections alone might provide an unambiguous answer.
Among various charmonium decays, only J/ψ → e + e − and η c → γγ have now been studied to O(α 2 s ).Clearly these are the two simplest channels since they do not involve any complications related to the dynamics of hadrons in the final state.If the understanding of even those simplest decays encounters such difficulties, one should exercise great care when extracting physical information from more complicated charmonium decays.The fact that perturbative QCD appears to work well in the one-loop order is certainly insufficient to ensure that heavy quarkonia are well understood.

4 FIG. 2 :
FIG. 2:The second order correction in eq.(21) as a function of the charm quark mass.