The fragmentation functions for a quark into a spin-singlet quarkonium: Different flavor case

In the paper, we calculate the fragmentation functions for a quark to fragment into a spin-singlet quarkonium, where the flavor of the initial quark is different from that of the constituent quark in the quarkonium. The ultraviolet divergences in the phase space integral are removed through the operator renormalization under the modified minimal subtraction scheme. The fragmentation function $D_{q \to \eta_Q}(z,\mu_F)$ is expressed as a two-dimensional integral. Numerical results for the fragmentation functions of a light quark or a bottom quark to fragment into the $\eta_c$ are presented. As an application of those fragmentation functions, we study the processes $Z \to \eta_c+q\bar{q}g(q=u,d,s)$ and $Z \to \eta_c+b\bar{b}g$ under the fragmentation and the direct nonrelativistic QCD approaches.


I. INTRODUCTION
According to QCD factorization theorem, the cross section for the inclusive production of a hadron H with high transverse momentum (p T ) in a high-energy collision is dominated by the single parton fragmentation [1], i.e., where ⊗ denotes a convolution in the momentum fraction z, the sum extends over all species of partons. dσ A+B→i+X indicates the partonic cross section that can be calculated in perturbation theory, while D i→H indicates the fragmentation function for the parton i into a hardron H. µ F denotes the factorization scale which is introduced to separate the energy scales of the two parts. The factorization formula (1) was first derived by Collins and Soper for light hadron production [2]. This factorization formula can be equally applied to the heavy quarkonium production. The proof of the factorization formula (1) for the quarkonium production was presented by Nayak, Qiu, and Sterman [3]. The factorization formula (1) is called as leading power (LP) factorization because it gives the LP contribution in the expansion in powers of m H /p T . The factorization formula for the next-to-leading power (NLP) correction was derived in Refs. [4][5][6][7], and the NLP contribution comes from the double-parton fragmentation.
Fragmentation functions play an important role in the calculation of the cross sections under the LP factorization. Unlike the fragmentation functions for the production of the light hadrons which are nonperturbative in nature, the fragmentation functions for the heavy quarkonium production can be calculated through the nonrelativistic QCD (NRQCD) factorization [8]. Under NRQCD factorization, the fragmentation functions for a parton to fragment into a quarkonium can be written as (2) where d i→(QQ)[n] are short-distance coefficients (SDCs) which can be expanded as powers of α s (m Q ), and O H (n) are long-distance matrix elements (LDMEs).
The fragmentation functions for quarkonia have been studied extensively. Most of the fragmentation functions for the S-wave and P-wave quarkonia are known up to α 2 s order [9][10][11][12][13][14][15][16][17][18][19][20][21][22][23][24][25][26][27], and a few fragmentation functions for quarkonia were calculated up to α 3 s order [28][29][30][31][32][33][34][35]. Among these studies, the NLO corrections to the fragmentation functions for g → QQ[ 1 S [1,8] 0 ] have been calculated recently by three groups [28,[30][31][32]. The NLO fragmentation functions for g → QQ[ 1 S [1,8] 0 ] are important in prediction of the production of the η c,b and the h c,b at the LHC. However, the fragmentation functions for a quark into the η c,b are only available up to α 2 s order. Those fragmentation functions are also important to the precision prediction of the η c,b production at the LHC. Moreover, for the production of the η c,b in e + e − collisions, the quark fragmentation contribution is more important than the gluon fragmentation contribution due to that the cross section for a quark is α 0 s order but the cross section for a gluon is α s order. In this paper, we will calculate the fragmentation functions for a quark q into the η Q , where Q = c, b but q = Q.
In the early calculations of fragmentation functions for doubly heavy mesons [9,15], the fragmentation functions are determined through comparing the cross section calculated based on the NRQCD factorization with that calculated based on the factorization formula (1) for a process containing the doubly heavy meson being produced. In fact, the fragmentation functions can be defined through the matrix elements of non-local gaugeinvariant operators [2]. The operator definition for the fragmentation functions was first applied to calculations of the fragmentation functions for doubly heavy mesons by Ma [17]. The calculations based on the operator definition are particularly convenient to extend to higher orders. Therefore, we will calculate the fragmentation functions based on the operator definition suggested by Collins and Soper.
The paper is organized as follows. Following the Introduction, in Sec.II, we present the definition and the analytical calculation for the fragmentation functions. In Sec.III, we present the numerical results for the fragmentation functions D q→ηc (z, µ F )(q = u, d, s) and D b→ηc (z, µ F ), and apply the fragmentation functions to processes Z → η c + qqg(q = u, d, s) and Z → η c + bbg. Sec.IV is reserved for a summary.

II. THE ANALYTICAL CALCULATION FOR THE FRAGMENTATION FUNCTIONS
A. The definition of fragmentation function The fragmentation functions are usually defined in the light-cone coordinate system. In this coordinate sys-tem, a d-dimensional vector V is expressed as Then the product of two vectors is The gauge-invariant fragmentation function for a quark q to fragment into a hadron H is defined as [2] where Ψ is the field of initial quark, A µ a is the gluon field, and t a (a = 1 · · · 8) are SU (3)-color matrices. The longitudinal momentum fraction is defined as z ≡ P + /K + , where K is the momentum of the initial quark. The fragmentation function is defined in a reference frame in which the transverse momentum of the hadron H vanishes. It is convenient to introduce a light-like momentum whose expression is n µ = (0, 1, 0 ⊥ ) in the reference frame where the definition of the fragmentation function carried out. Then, z can be expressed as a Lorentz invariant, i.e., z = P · n/K · n. The Feynman rules can be derived from the definition (3) directly, and we have presented the Feynman rules in a previous paper [34].

B. The calculation of fragmentation function
The definition (3) is gauge invariant. However, for the practical calculation, the gauge should be specified. We adopt the usual Feynman gauge throughout the paper. There are ultraviolet (UV) divergences in the calculation. To deal with the UV divergences, we adopt dimensional regularization with d = 4 − 2ǫ, then the UV divergences appear as the pole terms in ǫ.
In the calculation, we first calculate the fragmentation function for an on-shell QQ pair in 1 S Then the fragmentation function D q→ηQ can be obtained from D q→(QQ)[ 1 S [1] 0 ] through replacing the LDME 0 ](p 1 ) + g(p 2 ) + q(p 3 ) under the Feynman gauge, which can be collectively represented by four diagrams in Fig.1. The squared amplitudes, corresponding to four diagrams in Fig.1, can be written as (1)

FIG. 1. The cut diagrams for the fragmentation function
where Π 1 is the spin-singlet projector where 1 is the unit matrix of the SU (3) c group. There is a common factor X µν ab which arises from the annihilation of a virtual gluon into a (QQ)[ 1 S [1] 0 ] pair and a real gluon. The factor X µν ab can be expressed as where We employ the package FeynCalc [36,37] to carry out the Dirac and color traces, and then the total squared where .
The coefficients c i (s 1 , y, z) can be easily extracted, and we do not list their expressions here. The differential phase space for the fragmentation The contributions from the cut diagrams shown in Fig.1 can be calculated through is a factor from the definition of the fragmentation function. The integral on the right hand side of Eq. (15) is UV divergent with d = 4. This UV divergence is regularized by dimensional regularization with d = 4 − 2ǫ, and the integral generates 1/ǫ terms. To perform the integration in Eq. (15), it is important to choose proper parametrization for the phase space. We present a parametrization for the phase space in Appendix A.
The differential phase space given in Eq.(A8) can be expressed as follows where N g (p 1 , p 2 ) is defined as and dφ 2 (p 1 , p 2 ) is defined as 1 Here, we associate the scale factor µ 4−d with each dimensionally regulated integration in d space-time dimensions. In our previous papers [33,34], this scale factor was put in the squared amplitudes.
where the range of s 1 is from (4m 2 Q y/z) to ∞. dφ 2 (p 1 , p 2 ) stands for the differential phase space for a gluon with longitudinal momentum yK · n to fragment into a (QQ)[ 1 S [1] 0 ]-pair with longitudinal momentum zK · n at LO. According to Eqs. (16), (17), (18) and (A8), the expression of dφ (3) (p 1 , p 2 , p 3 ) can be derived The range of y is from z to 1, and the range of s is from [s 1 /y + m 2 q /(1 − y)] to ∞. The integrations over Ω 3⊥ and s of A (1−4) can be performed using the method introduced in Ref. [34]. Then, we obtain where N g dφ 2 (p 1 , p 2 ) ≡ N g (p 1 , p 2 )dφ 2 (p 1 , p 2 ). The expression of f (s 1 , y, z) are given in Appendix B. The contribution D (z) contains a UV pole, it should be removed through the operator renormalization [38]. We carry out the renormalization using the MS procedure. Then, the fragmentation function under the MS scheme can be obtained through where µ F is the factorization scale, the expression of the splitting function P gq (y) is and , i.e., where the expression of A g→(QQ) The integral over s 1 in Eq. (23)  is consistent with that obtained in Refs. [12,28,31,32]. Applying Eqs. (20) and (23) to Eq. (21), we can obtain the fragmentation function under the MS scheme. It is found that the UV pole of D (z) is exactly canceled by the UV pole of the counter term from the operator renormalization. The remaining integrals no longer generate divergence, we can set ǫ = 0 before carrying out the integrations. Multiplying the fragmentation function D q→(QQ)[ 1 S [1] 0 ] (z, µ F ) for the (QQ)[ 1 S [1] 0 ] pair by a factor O ηQ ( 1 S (0)| 2 /(4π), we obtain the fragmentation function D q→ηQ (z, µ F ), i.e., where Here, R (QQ) S (0) is the radial wave function at the origin for the (QQ) bound state.

III. NUMERICAL RESULTS AND DISCUSSION
In this section, we will present the numerical results for the fragmentation functions and apply the fragmentation functions to the cross sections for the η c production.
The value of |R cc S (0)| 2 is taken from the potential model calculation [39]. For the strong coupling constant, we adopt two-loop formula as used in our previous paper [33], where α s (2m c ) = 0.259.

A. The fragmentation functions
The fragmentation function for a light quark into η c , where µ F = 2m c , 4m c and 6m c , is shown in Fig.2. In the numerical calculation, the mass of the light quark is neglected, and the strong coupling is taken as α s (2m c ). From Fig.2, we can see that the fragmentation function is sensitive to the factorization scale. When µ F = 2m c , the fragmentation function increases first (z < 0.96) and then decreases (z > 0.96) with the increase of z, and the fragmentation function is less than 0 in the most z region; When µ F = 4m c , the fragmentation function also increases first (z < 0.12) and then decreases (z > 0.12) with the increase of z, but the fragmentation function is greater than 0 in the most z region; When µ F = 6m c , the fragmentation function decreases monotonically with the increase of z, and the fragmentation function is greater than 0 for z ∈ (0, 1). The fragmentation function has a singularity at z = 0.
In order to understand the dependence of the fragmentation function on the factorization scale, we have calculated the coefficient of ln(µ 2 F /m 2 c ) in the fragmentation function D q→ηc (z, µ F ). The coefficient is shown as a function of z in Fig.3 with a linear scale (upper one) and with a logarithmic scale (lower one). We can see that, the coefficient decreases monotonically with the increase of z and is positive for z ∈ (0, 1). Like the fragmentation function, the coefficient of ln(µ 2 F /m 2 c ) also has a singularity at z = 0. When µ F is very large, the fragmentation function is dominated by the ln(µ 2 F /m 2 c ) term. Therefore, the behavior of the fragmentation function is similar to that of the coefficient of ln(µ 2 F /m 2 c ) when µ F is large enough.
In Fig.4, the fragmentation function D b→ηc (z, µ F ) for µ F = m b , m b + 2m c and m b + 4m c is presented. We can see that the fragmentation function D b→ηc (z, µ F ) is also sensitive to the factorization scale. Actually, from Eqs. (25) and (26), we can see that the coefficient of ln(µ 2 F /m 2 Q ) in D b→ηc (z, µ F ) is the same as that in D q→ηc (z, µ F ).

B. Application to the decay widths
In this subsection, we will apply the obtained fragmentation functions to the decay widths for the processes Z → η c + qqg and Z → η c + bbg.
Here, we only present the calculation formulas for Z → η c +qqg, the formulas for Z → η c +bbg are similar. Under the fragmentation approach, the differential decay width dΓ/dz for Z → η c + qqg at LO can be written as where the energy fraction is defined as z ≡ E ηc /E max ηc , and E ηc and E max ηc are the energy and the maximum energy of the η c in the rest frame of the Z boson. Γ Z→q+q denotes the LO decay width for the process Z → q +q, and the factor of 2 is due to that the contributions from the q fragmentation andq fragmentation are the same. In Eq.(28), we have used the fact dΓ Z→q+q /dy = Γ Z→q+q δ(1 − y) at LO. dΓ Z→g+qq /dy is the differential decay width for the inclusive production of a gluon associated with a light-quark pair. Neglecting the light quark mass, we obtain the differential decay width under the MS factorization scheme D g→ηc (z) and D q→ηc (z, µ F ) are the LO fragmentation functions. The expression of D g→ηc (z) can be found in Ref. [12], and the expression of D q→ηc (z, µ F ) has been given in Eq. (25). It is easy to check that the logarithm terms of µ 2 F in Eq. (28) are canceled by each other which results in that dΓ Frag Z→ηc+qqg /dz is independent of µ F .
The physical picture of Eq.(28) is as follows: The first term gives the contribution from that the Z boson decays into a light quark and a light antiquark, with energies m Z /2 on a distance scale of order 1/m Z . Subsequently one of the light quark and the light antiquark decays into an η c on a distance scale of order 1/m c . The second term gives the contribution from that the Z boson decays into a light quark-antiquark pair and a gluon on a distance scale of order 1/m Z , and the gluon decays into an η c on a distance scale of order 1/m c . The two terms of Eq.(28) share the same Feynman diagrams which are shown in Fig.5, but they come from different regions of the phase space. When the invariant mass of the virtual light quark (antiquark) is very small compared to µ F , the contribution is given by the first term of Eq. (28). When the invariant mass of the virtual light quark (antiquark) is very large compared to µ F , the contribution is given by the second term of Eq. (28) (3) (2) (1) In addition to the fragmentation approach, we can also calculate the decay width directly based on the NRQCD factorization, i.e., where we use "Direct" to denote the results from the direct calculation based on the NRQCD factorization.
The dominant contributions to the decay widths for Z → η c + qqg and Z → η c + bbg come from four fragmentation diagrams shown in Fig.5. The contributions from the non-fragmentation diagrams are suppressed by powers of m c /m Z compared to the fragmentation contributions. For simplicity, under the direct NRQCD calculation, we only consider the contributions from the four fragmentation diagrams shown in Fig.5.
The differential decay widths for Z → η c + qqg and Z → η c + bbg under the fragmentation approach and the direct NRQCD approach are presented in Fig.6. From the figure, we can see that the curves from the fragmentation calculation are very close to those from the direct NRQCD calculation, especially for z ≥ 0.2. The total decay widths can be obtained through integrating the differential decay widths dΓ/dz over z. The total decay widths under the fragmentation approach and the direct NRQCD approach are given in Table I. We can see that the total decay widths obtained from two approaches are also very close. Therefore, the fragmentation approach provides a good approximation to the direct NRQCD calculation.  3   TABLE I. The decay widths of Z → ηc +qqg and Z → ηc +bbg under the fragmentation approach and the direct NRQCD approach. For the Z → ηc + qqg case, the contributions for q = u, d, s are summed.

IV. SUMMARY
In the present paper, we have calculated the fragmentation functions for a (heavy or light) quark into a spinsinglet quarkonium, where the flavor of the initial quark is different from that of the constituent quark in the quarkonium. There are UV divergences in the phasespace integral, which are removed through the operator renormalization of the fragmentation function. We have carried out the renormalization under the MS scheme. The fragmentation function D q→ηQ (z, µ F ) is given as a two-dimensional integral, and this two-dimensional integral can be calculated easily through numerical integration. Numerical results for a light quark or a bottom quark into the η c with several factorization scales are analyzed. The results show that these fragmentation func-tions are sensitive to the factorization scale. Especially, when µ F is small, the fragmentation functions are negative at small z values. There is a singularity at z = 0 for these fragmentation functions.
We have applied the obtained fragmentation functions to the decay widths for the processes Z → η c + qqg(q = u, d, s) and Z → η c + bbg. The differential decay widths and total widths are calculated under the fragmentation and direct NRQCD approaches. It is found that the results under two approaches are close to each other. Therefore, the fragmentation approach provides a good approximation to the direct NRQCD calculation. Moreover, the fragmentation functions obtained in this paper can be used to the studies on the production of η c and η b at high-energy colliders.