$\chi_{cJ}\rightarrow e^{+}e^{-}$ decays revisited

We present a calculation of the width for $\chi_{cJ}\rightarrow e^{+}e^{-}$ decay. The amplitude of the process is computed within the NRQCD framework. The leading-order contribution is described by two terms associated with the two different integration domains in the electromagnetic loop describing two-photon annihilation of the heavy quark-antiquark pair. The corresponding operators are defined in the framework of NRQCD. The matrix element of one of these operators describes a configuration with an ultrasoft photon and can be associated with the higher Fock state of the heavy meson. In order to compute this contribution we use the heavy hadron chiral perturbation theory. We obtain that this contribution is numerically dominant. The obtained estimates for the decay widths of the $\chi_{c1}$ and $\chi_{c2}$ states are $0.09$ eV and $0.06$ eV, respectively.


Introduction
The leptonic decays of C-even charmonium states into a lepton pair have a very small branching ratio because they can only occur via a two-photon intermediate state χ cJ → γ * γ * → e + e − .However with the high-luminosity BEPC-II e + e − -accelerator operating on the charmonium energy region such measurements of direct production cross section e + e − → χ cJ become feasible.The study of the mechanism of the production of C-even quarkonium states is especially interesting in view of the production of higher resonances such as the exotic charmonium-like state X(3872).
The decays χ cJ → l + l − have already been studied long time ago in Ref. [1].The authors considered the quarkonium description and a vector dominance model (VDM) in order to describe the decay amplitudes of charmonium states with J = 1, 2. It was found that a naive quarkonium description is problematic because of infrared logarithmic divergencies arising in the integrals describing the quark-photon annihilation loop.Such divergencies indicate that the corresponding contribution is also sensitive to long distance physics.The corresponding integrals in Ref. [1] have been regularized by introducing the binding energy M χ − 2m Q 500 MeV.The numerical estimates obtained in this way give very small value of the widths, smaller thn bounds derived from analyticity and unitarity in the same work.Probably, a more realistic estimate was obtained using generalized VDM which yields a larger numerical value for the widths, consistent with the unitarity constrains: Γ[χ c1 → e + e − ] 0.46 eV and Γ[χ c2 → e + e − ] 0.014 eV [1].
Recently the decay rate of the χ c1 state was again estimated in [2] using the VDM approach, with result Γ[χ c1 → e + e − ] 0.1 eV.The short distance contributions describing a configuration with highly virtual photons were considered in this framework as unknown contact vertices giving rise to a theoretical uncertainty.
The aim of the present work is to provide a more systematic description of the decay amplitudes for χ cJ → γ * γ * → e + e − process using the NRQCD factorization framework [4,5], see also review [6] and references theirin.This technique allows one to perform a systematic description of heavy quarkonium states using the small relative velocity of heavy quarks and the small QCD running coupling at short distances.Within this framework we associate the IR-divergencies found in [1] with a specific quarkphoton operator which describes a configuration with an ultrasoft photon.The matrix element of this operator describes the overlapping with the higher Fock state which consists of heavy quark-antiquark and photon.This allows us to perform a systematic separation and description of the different contributions relevant in a leading-order expansion in small velocity v.
Our paper is organized as follows.In Sec. 2 we briefly describe our notation and provide definitions of various quantities used in the following.Sec. 3 is devoted to the investigation of the one-loop integral which describes the leading-order contribution.In this section we establish the dominant regions and provide a description of the amplitude within the NRQCD factorization framework.Sec. 4 is devoted to the calculation of the ultrasoft photon matrix element in the heavy hadron chiral effective theory (HHχPT).Furthermore, a estimate of the decay rates is given.In Sec. 5 we briefly summarize our results.

Kinematics and notation
Let us start from the description of the decay kinematics χ cJ (P ) → e + (l 1 )e − (l 2 ).The initial state momentum can be written as where ω denotes the charmonium velocity.In the following, we consider the charmonium rest frame which implies ω = (1, 0). ( The small relative velocity of heavy quarks in the bound state is denoted as v. Neglecting lepton masses, the lepton momenta can be written as where n and n denote the light-like vectors which satisfy (n • n) = 2. Any 4-vector V µ can be expanded as where V ⊥ denotes the components transverse to the light-like vectors : Similarly, one can also write a decomposition where V denotes the component which is orthogonal to the velocity ω: (ω • V ) = 0.In the following we assume that in the rest frame The momenta of the heavy quark and antiquark with mass m which form a quarkonium state can be written as where the relative momentum ∆ satisfies The heavy quarks which create a bound state are non-relativistic, implying that the relative velocity v ∼ ∆/m is quite small: v 1.The power counting rules of NRQCD has been established in [3,4].Following these arguments we assume that the mass m is large enough and that the most important scales such as mass m, typical three-momentum of the heavy quark ∼ mv and its typical kinetic energy ∼ mv 2 satisfy Integrating out the modes with hard momenta p h ∼ m c one passes onto the effective theory NRQCD which describes the modes with the soft momenta p s ∼ mv.If the scale mv Λ QCD then one can integrate over the soft region together with potential gluons with momenta [5,[7][8][9][10] After this one obtains a new effective theory which is known as potential NRQCD (pNRQCD).For a more detailed information about these effective theories see Ref. [6] and references therein.
The charm quark mass m c 1.5 GeV is not large enough compared to Λ QCD therefore in this case one can only factorize the effects at momentum scales of order m c .However, in the QED sector one can also consider the possibility to integrate over the soft region too.As we will show further on, such situation is relevant for the factorization of the electromagnetic loop describing the transition cc → γ * γ * → e + e − .
After factorization of the hard contribution, the nonpertubative QCD dynamics is described by the matrix elements of appropriate operators defined in an effective theory.In the following, we will need the following set of NRQCD operators which describe the matrix elements between the charmonium and the vacuum states: where the covariant derivative Furthermore, we use the covariant four-component fields ψ ω , χ ω to describe the soft quark and antiquarks within the NRQCD framework.These fields satisfy Using Eq.( 2), one can show that the operators in Eqs.( 11)-( 14) can be reduced to the set of well-known non-relativistic operators constructed from two-component Pauli spinors.
The matrix elements of these operators are well known in the literature and can be written as 0| O( The constants O( 3 S 1 ) and O( 3 P 0 ) are related to the value of the charmonium wave functions at the origin where R nl (r) is the radial part of the Schrödinger wave function and R nl (r) denotes its derivative.The rhs of Eqs.( 18)-20) depends on of the same constant O( 3 P 0 ) due to the spin symmetry of the leading non-relativistic action [4].The polarization vectors σ ψ (λ), β χ (λ) and αβ χ (λ) correspond to spin-1 and spin-2 charmonium states, respectively.They are normalized to satisfy with X = {J/ψ, χ c1 } and with M αβ = −g αβ + P α P β /M 2 χc2 .The decay amplitudes χ cJ → e + e − are defined as with and where ūn and vn denotes the spinors of the massles lepton and antilepton, respectively and The leading-order contribution to these amplitudes arises from the annihilation of heavy quarks into two photons which create the outgoing lepton pair, see Fig. 1.If the one-loop integral is dominated by the hard region where both photons and heavy quark are highly off-shell then one can expect that such process can be described within the NRQCD approach and the amplitude can be factorized into hard and soft parts.In the next section we consider this possibility in more detail.

Factorization of decay amplitudes in NRQCD
The leading-order in α s diagrams describing the e + e − decay of C-even charmonia are shown in Fig. 1.These one-loop diagrams are constructed from the photon, lepton and heavy quark (double lines).The diagrams in Fig. 1 can be computed in the heavy quark mass limit m → ∞, performing an expansion in the small parameter ∆/m ∼ v. Let us to start from a naive guess that the dominant contribution is only provided by the hard region where the loop momentum k µ ∼ m, and therefore all propagators are far of off-shell.The leading-order contribution in 1/m is provided by projections onto the leading-order operators O( 3 P J ) described in Eqs.( 18)- (20).The technical details are well known in the literature, see e.g.[1].The resulting expressions can be presented as where the square brackets for the propagators denote the standard Feynman prescription The corresponding contribution to the amplitude χ c0 → e + e − vanishes and therefore is suppressed by a power of v and will not be considered it in this work.
The total structure of the integrands in expressions (30) and (31) can be divided into the lepton and heavy quark parts.The lepton part has ūn . . .v n in the numerator and includes the photon and lepton propagators in the denominators.The heavy quark part is given by Tr P i Γ αβµ (k) .We introduced the projections P J onto charmonium states The expression for Γ αβµ (k) reads where e Q is the charge of the heavy quark (e c = 2/3).The small squared relative momentum ∆ 2 ∼ (mv) 2 which appears in the heavy quark propagator provides an IR-regularization and can be neglected if it is not required.With this regularization the traces and loop integrals are computed in four dimensions with dk ≡ d 4 k/(2π) 4 .The expressions (30) and (31) have been obtained by expanding the heavy quark fields in position space c(y) e −im(ωy) 1 + y and projecting the soft quark fields χ † ω and ψ ω onto leading-order operators ( 12)-( 14).The terms ∼ y•∂ψ ω (arising from the multipole expansion of the soft quark field arguments ) lead to the expansion of the integrand with respect to small relative momentum ∆ giving the contribution D µ Q .The terms proportional to ∼ 1 2m / D give the contribution with DQ .The evaluation of the integrals in Eqs. ( 30) and (31) gives These expressions are in agreement with the results obtained in Ref. [1].We obtain that both amplitudes depend on the large logarithm ∼ ln ∆ 2 /m 2 which is sensitive to the soft scale ∆ 2 .This shows that the starting assumption about one dominant region k ∼ m is incorrect.There must be at least one more domain where some propagators in the loop integral are soft.One can expect that the additional region is associated with the configuration when one of the photons is soft.In this case the propagator of the heavy quark is also soft and the hard configuration is described by the tree level subdiagram describing the annihilation cc → e + e − through one photon.
In order to get an idea about the explicit definition of this region it is useful to investigate the integrals of diagrams in Fig. 1 within the threshold expansion technique worked out in Ref. [5].According to this analysis the threshold kinematics is described by the following regions hard : soft : The same regions can also be considered for the photon with momentum P − k.These regions can be associated with the fields appearing in the effective Lagrangians, see e.g.Ref. [6].
According to the threshold expansion prescription an integrand is expanded in each domain to a required accuracy and the resulting integral is computed in dimensional regularization.A detailed analysis of the full expressions in Eqs. ( 30) and ( 31) is quite similar.To be definite let us consider the integral which enters in Eq.( 30) Keeping the denominators of the heavy quark propagators in Γ αβµ unexpanded where P 0 ∼ m, P 2 0 /4 − m 2 ∼ (mv) 2 , ∆ ∼ mv.In the hard region, the small scalar products with ∆ and the term P 2 0 /4 − m 2 m 2 can be neglected resulting in which appear in the expressions (35) and (36) (up to small regularization term ∆ 2 ).From dimensional counting one immediately finds where Γ denotes the Dirac structure.Computing the hard integral J h in dimensional regularization one finds the IR poles 1/ε.These singularities must cancel in the sum with other contribution.
Expanding the integrand (44) in the soft region (41) yields where Γ αβµ s (k) is given by (34) with Calculating the trace and performing the contractions in the numerator results in As the integral in (51) is scaleless it therefore vanishes in the dimensional regularization, i.e.J s = 0.
In the potential region (42), the expansion of the heavy quark propagator reads The computation of the corresponding integral then yields However the poles in k 0 in the integrand of Eq.( 53) lie in the same imaginary half-plane and therefore the integral over k 0 vanishes.This observation is also true for the higher order contributions in v appearing from this domain.We can therefore conclude that the potential region cannot contribute in this case.
In the ultasoft domain (43), the heavy quark propagators are expanded as Performing the expansion of the integrand one gets This integral has the same scaling behavior ∼ m −3 as the hard integral J h in Eq.( 48).One can also see that the integral in Eq.( 55) is UV divergent.The similar analysis can also be carried out for the second photon with momentum k − P .Therefore we conclude that the exact integral must be given by sum where J us denotes the contributions from the both ultrasoft domains.This conclusion can be checked by explicit calculations.A similar conclusion for the two-photon diagrams in Fig. 1 has also been obtained in Ref. [11].
Guided by this consideration we suggest that the additional relevant domain is described by the ultrasoft region.In order to find the description of the appropriate operator in the effective theory one has to integrate out hard and soft photons and leptons.After that the description of QED sector includes only collinear leptons and ultrasoft photons.The integration of the soft photon with the lepton and quark must be described in the framework of the effective theory.
Within the above picture the factorization of the decay amplitudes can be described as a sum of two contributions The first term on rhs of this equation corresponds to the hard domain with the hard photons, C (J) γγ denotes the corresponding hard coefficient function.
The second term on rhs of Eq.( 57) corresponds to the domain with the ultrasoft photon.The operator O σ ( 3 S 1 ) is defined in Eq. (11).The outgoing collinear leptons are described by fields ξn and ξ n which defined as ξn The photon Wilson lines Y † n and Y n describe the interaction of the ultrasoft longitudinal photons with the energetic lepton and antilepton and read where B us µ denotes the ultrasoft photon field.The appearance of these Wilson lines is related with the fact that in a general gauge the tree level diagram with attachments of n • B us photon to the collinear field ξn describing the outgoing lepton1 are resummed to the P-ordered exponents The leading-order hard coefficient function C γ is defined by the diagram in Fig. 2 and reads The soft and collinear modes in the effective action describing the QED sector are decoupled.This property is well known in the soft-collinear effective theory, see e.g.Refs.[12][13][14].This allows us to contract the lepton fields in the second matrix element in Eq.(57) and rewrite it as with The presence of the soft scale ∆ 2 in Eqs. ( 38) and (39) can be explained by the contribution with ultrasoft photon.Therefore in order to find the hard coefficient functions C  matching onto the configuration described by Eq.(57).For that purpose we need to compute the ultrasoft matrix element (62) in the effective theory.
The interaction of ultrasoft photons with quarks are described within the pNRQED.The ultrasoft photons have momentum p ∼ mv 2 so that photon field scales as The scaling of the quark fields reads Using this counting one finds At the same time However the pure quark operator O σ ( 3 S 1 ) is C-odd and therefore it cannot contribute to the matrix element with a C-even charmonium state In order to obtain a nontrivial contribution one needs to consider at least one interaction of an ultrasoft photon with the quark in pNRQED.We only need the two-particle sector describing the electromagnetic interactions of quarks (in rest frame ω = (1, 0)) and analogous contributions with antiquark fields.The arguments of the ultrasoft photon field are expanded because the space components of the quark fields varies at x ∼ 1/mv, the measure scales as dx 0 ∼ 1/mv 2 , d 3 x ∼ (mv) −3 .With these rules one finds that L em 0 ∼ v 0 and L em 1 ∼ v 1 .The leading-order term (69) provides the soft quark propagator A nontrivial contribution to the matrix element 0| ... |χ cJ can be obtained from T -product which is of the same order as the hard contribution in Eq.(66).Calculation of this T -product gives diagrams shown in Fig. 3.The dashed lines can be associated with the collinear leptons or equivalently with the ultrasoft Wilson lines (59).These diagrams induce a mixing of the operators O σ γ ( 3 S 1 ) and O( 3 P J ) due to electromagnetic interaction in the framework of pNRQED.
In order to perform the matching onto operators according to formula (57) one has also to compute the contribution of the diagrams in Fig. 3.The simplest way to proceed is to follow the same technique as we used above for diagrams in Fig. 1.Let us consider χ c2 as initial state.In this calculation we set P 0 = 2m and only keep the relative momentum ∆.Then the sum of all four diagrams gives Computing these diagrams we project the soft quarks fields on the operator O( 3 P 2 ) and substitute the corresponding matrix element which gives the factor σρ χ i O( 3 P 0 ) , the coefficient 4 arises from the sum of the four diagrams shown in Fig. 3, the ultrasoft loop integral reads This integral coincides with the ultrasoft integral of Eq.( 55) obtained within the threshold expansion approach up to term P 2 0 /4 − m 2 which vanishes because we set P 0 = 2m.The integral in Eq.( 74) is UV-divergent and we use dimension regularization D = 4 − 2ε in order to compute it.The result reads where µ F is the factorization scale.The 1/ε pole is the UV-pole which describes UV-mixing of the operators O σ γ ( 3 S 1 ) and O σ ( 3 P J ), schematically where [O] R on the lhs of Eq.(76) denotes the renormalized operator.Furthermore, Z J ∼ e 2 /ε is the corresponding renormalization constant.Assuming M S-subtraction scheme one finds Hence we obtain The soft matrix element for the χ c1 can be computed in the same way, resulting in The hard coefficients C (J) γγ are given by The doted lines denote potential gluons with momenta given by Eq.( 42).
where the expressions for A J are given by Eqs.( 38) and (39).The important check of the factorization formula (57) is the cancellation of the ultrasoft scale ∆ 2 in the expressions for C (J) γγ obtained from Eq.( 80).Substituting the computed expressions in Eq.( 80) we obtain These expressions are the main result of this section.We observe that the soft scale cancel in Eqs.( 81) and ( 82) as it must be.Hence the factorization formula described by Eq.( 57) describes properly the ultrasoft region of the one-loop diagram.
The coefficient function C γγ has an imaginary part which originates from the two-photon cut.Such mechanism can not work for χ c1 state, therefore C (1) γγ is real.The hard coefficient functions depend on the factorization scale µ F .Therefore and the independence of the amplitude A J on µ F yields the evolution equation The solution of this equation depends on the initial condition defined at some scale µ 0 .Performing numerical estimates one has to fix a value of this scale.By derivation this scale separates the hard region (two hard photons) from the ultrasoft region (hard and ultrasoft photons).Therefore it is natural to associate this scale with the virtuality of the ultrasoft photon and to set µ 0 to be of order 300 − 500 MeV.Then the matrix element of the operator O σ γ ( 3 S 1 ) on the rhs of Eq.( 83) describes only the ultrasoft nonperturbative contribution which can be only estimated within some low-energy effective theory or model.Similar to the well known color octet mechanism, the operator O σ γ ( 3 S 1 ) can also be associated with the electromagnetic mechanism.The corresponding matrix element can be interpreted as an overlap with the higher Fock state |Q Qγ which includes a dynamical photon while the matrix elements of the operators O( 3 P J ) describe the coupling to the dominant quark-antiquark state.Therefore the full description of the leptonic decay requires a knowledge on the subleading structure of the quarkonium state.
In the large mass limit m → ∞ one can consider a specific situation known as the Coulomb limit when the binding energy is larger then the typical hadronic scale E ∼ mv 2 Λ QCD .In this case the strong coupling is quite small α s (mv) ∼ v and ultrasoft contribution can be estimated within the pNRQCD.Then one has to compute the diagram as in Fig. 4 resumming the interactions with Coulomb gluons.Such calculation has been carried out for the radiation function in Ref. [15].Perhaps, such calculation might also be interesting here in order to get an idea about the relative value of this matrix element in the Coulomb limit.In present paper we will obtain an estimate of the ultrasoft matrix element using the so-called heavy hadron chiral perturbation theory (HHχPT) framework in the next section.

Calculation of the ultrasoft matrix element in the heavy hadron chiral perturbation theory
In order to provide a numerical estimate of the decay rate we need to estimate the ultrasoft matrix elements which describes an overlap with the higher Fock component of the charmonium state χ cJ in which a dynamical photon is present.One can expect that the soft photon has already quite large wavelength and therefore it interacts with heavy charmonium as with a point-like source.Then it is natural to expect that the relevant dynamical degrees of freedom in this case are associated with mesonic fields and the corresponding low energy dynamics is described by the most generic effective action compatible with the symmetries of NRQCD.Such an approach is known as heavy hadron chiral perturbation theory in Refs.[16,17] for the heavy-light mesons and then generalized on quarkonia in Refs.[18][19][20].This framework can also be used for the calculation of the matrix element in Eq.( 85).
For our purpose we need only the electromagnetic sector of the HHχPT described by the effective action which includes the kinetic terms for J/ψ and ψ states and the vertices describing the electromagnetic vertices χ cJ J/ψ γ and χ cJ ψ γ2 .As before we assume the rest frame for the initial state χ cJ .The kinetic Lagrangian reads with the residual masses ∆M = (M describes the residual motion of the heavy J/ψ and ψ particles and satisfy ω µ ψ (ω) = 0.The leading-order in 1/m effective Lagrangian describing the radiative decays χ cJ → J/ψ + γ and ψ → χ cJ + γ reads [20] with and The currents J S , J S and J µ P describe particles from S-and P -wave multiplets, respectively.In Eq.(87) we introduced the dimensionless couplings f γ and f γ .The fields χ J describe charmonium states χ cJ .Computing the trace in Eq.(87) one finds where we show only the relevant terms.Our calculations involve operators O(3 P J ) and O( 3 S 1 ) which have also to be matched onto physical quarkonium fields.The spin symmetry in the heavy quark limit yields where αβ are spinor indices.Taking the matrix element and computing the traces one can see that Eqs.( 93) and (94) reproduce correctly the matrix elements ( 16)- (20).Using these results one finds Computing the T -product in Eq.(96) gives the diagrams in Fig. 5.These diagrams are UV-divergent and we use in our calculation the dimensional regularization and M S subtraction scheme.The results read where iε ⊥ [σα] ≡ iε σαβρ n β nρ /2 and From these results one sees that the spin symmetry of NRQCD relates the soft photon matrix elements for J = 1 and J = 2 which are defined by the same nonperturbative couplings f γ O( 3 S 1 ) and f γ O ( 3 S 1 ) .The imaginary part in Eq.( 99) corresponds to the photon-quarkonium (J/ψ) cut in the diagrams in Fig. 5. Hence we conclude that the two-photon cut appears only in the hard photon contribution.The contribution with ψ has no physical cut and in the diagram this is provided by the negative value of ∆ M .The UV-poles in the HHχPT diagrams appear due to the mixing of the operators O σ γ ( 3 S 1 ) and χ cJ .Therefore this UV-pole can be absorbed into renormalization of the chiral constant in front of the operators χ cJ 3 .The expression for the total amplitude now reads We set the value of the chiral scale µ χ = µ 0 , defining the chiral couplings CJ (µ 0 ) as product of the two-photon hard coefficient functions C (J) γγ and constant O( This defines the expressions for the amplitudes which will be used for our numerical estimates.

Numerical estimates
In order to perform numerical estimates we need the values of the nonperturbative parameters O( 3 S 1 ) , O( 3 P 0 ) and f γ .Two of them are related to the values of quarkonium wave function at the origin, see Eqs.( 21) and (23).Their absolute values have been estimated in Ref. [21] using different models for the potential.In our numerical calculations we use the values obtained for Buchmüller-Tye potential [23] |R 21 (0)| 2 0.075GeV 5 , (103) We also assume that they correspond to positive values: The absolute values of the electromagnetic couplings f γ and f γ can be estimated from the decays χ cJ → J/ψγ and ψ → χ cJ γ.Using for widths Γ[χ c1 ] = 0.84 × 10 −3 GeV, Γ[χ c2 ] = 1.93 × 10 −3 GeV and branching fractions Br[χ c1 → J/ψγ] = 0.340 , Br[χ c2 → J/ψγ] = 0.192 from [22] we obtain where k 0 = (M 2 χ cJ − M 2 ψ )/2M χ cJ is the photon energy.Similarly, using widths Γ[ψ ] = 0.299 × 10 −3 GeV and branching fractions Br[ψ → χ c1 γ] = 0.096 , Br[ψ → χ c2 γ] = 0.091 from [22] we obtain We also need to know the sign of this coupling which can only be defined by a specific nonperturbative calculation.It turns out that this coupling can be represented as an overlap integral of the radial wave functions.Comparing our results for the decay amplitudes χ cJ → J/ψγ with the ones computed in Ref. [24] we find where factors √ 2M appear due to relativistic normalizations of the hadronic states.The analogous expression also holds for the coupling f γ .The overlap integral has been computed in the framework of potential models, see e.g.Refs.[25,26].Its value is found to be positive for f γ and negative for f γ .Therefore we assume in the following that f γ > 0 and f γ < 0.
The expressions for the decay width read where C γ is given by Eq.( 61).We use m c = 1.5 GeV for the mass of the charm quark and compute h(µ 0 ) by substituting M χ = (M χc1 + M χc2 )/2 in the expression (99).
Our numerical results are presented in Table 1 for different values of µ 0 .The subscripts s and h denote contributions from the soft and hard photon terms and hs corresponds to the interference of these contributions.In all cases the largest numerical contribution is provided by the ultrasoft matrix element.This contribution is relatively large and it weakly depends on the factorization scale µ 0 .Our estimates for Γ[χ c1 → e − e + ] is approximately factor 5 smaller then the estimate in Ref. [1] and in a good agreement with the estimate in Ref. [2].For Γ[χ c2 → e − e + ] our result is five times larger than one obtained in Ref. [1].
From Table 1 one can observe that the interference of the hard and ultrasoft contributions is numerically large for χ c2 width and relatively small for χ c1 .This can be explained as following.The imaginary part of h(µ 0 ) is numerically much larger than the real one, see Eq.(99).Further, the hard coefficient function C (1) γγ is real and therefore corresponding interference in the width depends only from the real part of h(µ 0 ).The imaginary part of C (2) γγ is not zero and therefore in this case the interference depends on the large imaginary part h(µ 0 ) and turns out numerically large.This observation allows one to conclude that the decay width χ c2 is quite sensitive to the relative sign of parameters R 10 and R 21 .In our estimate we used that these parameters has the same sign, see Eq.(105).However if they have opposite sign then the interference contribution is negative and this reduces the numerical value of the Γ[χ c2 → e − e + ] by factor 2.
In Ref. [1] it was shown that unitarity and analyticity allows one to constrain the minimal values of decay widths In the presented formalism these constrains are always satisfied because the soft contribution has a cut which yields the imaginary part required for the saturation of the bounds in Eqs.( 110) and ( 111).Therefore all our estimates shown in Table I are in agreement with these inequalities.As one can see from Table I the hard two-photon contribution is always smaller than the limiting value in both cases.
The same observation was also made in Ref. [1].This clearly indicates that the soft photon configuration provides a critically important contribution to these decay amplitudes.The derived approach can also be used for a description of leptonic decays of bottonium states χ bJ .These particles have almost the same branching fractions for χ bJ → Υ(1S)γ decay, see e.g.[22], but at present the widths of these states are not yet measured.Therefore, we cannot extract the decay coupling f (b) γ using experimental data.Instead, we use the estimates for the corresponding widths obtained in the model with a Cornell potential [25].The corresponding values can be found in Ref. [26] 115) is again defined with the help of of the overlap representation as in Eq.( 108) and corresponding estimates given in Ref. [26].The corresponding radial wave functions at the origin read [21] |R 21 (0)| 2 2.067GeV We observe that in this case the contribution of the ultrasoft configuration also remains larger than the hard one.

Conclusions
The decay width Γ[χ cJ → e − e + ] was computed using a factorization NRQCD approach.The dominant partonic subprocess was described by the annihilation of the heavy quark-antiquark pair into two photons: cc → γ * γ * → e + e − .The corresponding contribution is given by the one-loop diagram with two photons in the intermediate state.The dominant regions in the loop integral are associated with two configurations: hard photons and one ultrasoft and hard photons.The soft part of the contribution with ultrasoft photon overlaps with the higher Fock state |Q Qγ of the heavy meson, while the hard contribution overlaps with the leading two quark state.We have demonstrated that these contributions can be factorized and described by two different operators in NRQCD effective theory.The ultrasoft photon contribution is estimated using framework of the heavy hadron chiral perturbation theory.This allows us to obtain numerical estimates using a minimal set of the known nonperturbative parameters.Our estimates for charmonia χ c1 and χ c2 show that the ultrasoft photon configurations provide the numerically dominant contribution.This explains why the obtained numerical results for Γ[χ c1 → e − e + ] are in good agreement with the estimates obtained in Ref. [2] where only the usoft contribution was considered.We also expect that the developed formalism can be helpful to perform a more systemic analysis of decays if charmoniumlike state such as X(3872) → e − e + .

Aknowlegements
We are grateful to Achim Denig for useful discussions.This work is supported by the Helmholtz Institute Mainz.

Figure 1 :
Figure 1: One-loop diagrams describing the annihilation of χ cJ into an e + e − pair.

Figure 2 :
Figure 2: The hard one-photon exchange diagram.

Figure 3 :
Figure 3: The diagrams generated by the T -product (72) in pNRQED.The crossed circle denotes the vertex of the operator O γ ( 3 S 1 ), the dashed lines represent the Wilson lines associated with light-like directions n and n.The small crosses on the dashed lines show all possible attachments of the photon.

Figure 4 :
Figure 4: An example of the diagram in pNRQCD describing the ultrasoft matrix element of Eq. (62).The doted lines denote potential gluons with momenta given by Eq.(42).

Figure 5 :
Figure 5: The diagrams which describe the matrix element (96) in HHχPT.The crossed box denotes the operator (95), dashes lines describe the Wilson lines Y † n and Y n, black circle corresponds to the interaction vertices generated by L em SP (87).