Revisiting heavy quark radiative energy loss in nuclei within the high-twist approach

We revisit the calculation of multiple parton scattering of a heavy quark in nuclei within the framework of recently improved high-twist factorization formalism, in which gauge invariance is ensured by a delicate setup of the initial partons' transverse momenta. We derive a new result for medium modified heavy quark fragmentation functions in deeply inelastic scattering. It is consistent with the previous calculation of light quark energy loss in the massless limit, but leads to a new correction term in the heavy quark case, which vanishes in the soft gluon radiation limit. We show numerically the significance of the new correction term in the calculation of heavy quark energy loss as compared to previous studies and with soft gluon radiation approximation.


I. INTRODUCTION
In the past decade, tremendous progress has been made in understanding the jet quenching phenomena via various observables in high-energy nucleus-nucleus collisions [1][2][3]. Within perturbative-QCD-based theoretical frameworks for multiple parton scattering and parton energy loss in a nuclear medium, one is able to use the experimental data on jet quenching to probe the fundamental properties of the cold nuclei and the hot dense medium created in high-energy heavy-ion collisions [4][5][6][7]. One seminal study is the systematic extraction of jet quenching parameter,q, by the JET Collaboration [8], in which a global fitting to the experimental data on single inclusive hadron production at the Relativistic Heavy-Ion Collider and the Large Hadron Collider has been performed based on several jet quenching models for multiple parton scattering and parton energy loss.
Several of the successful jet quenching models in explaining the experimental data on hadron or jet production are based on the high-twist expansion approach [9][10][11][12], where contributions from multiple scattering between a propagating jet and medium partons can be effectively factorized as higher twist corrections to the vacuum fragmentation functions. These models utilize the generalized twist-4 factorization formalism developed by Qiu and Sterman [13,14]. The first calculation within this framework is performed in the process of semi-inclusive electron-nucleus deep inelastic scattering (SIDIS), where the light quark energy loss is encoded into the medium modified quark fragmentation functions [9,10]. This calculation has been further improved by going beyond the helicity amplitude approximation [11] and including multiple scattering [12] and double quark scattering [15]. The same techniques were also applied to evaluate heavy quark energy loss in SIDIS by considering charged-current interaction [16,17]. Similarly, one can also study the effect of initial state parton energy loss within the same framework, which has been evaluated in Drell-Yan dilepton production in proton-nucleus collisions. In analogy to the calculation in SIDIS, the effect of initial state parton energy loss in the Drell-Yan process is encoded in the medium modified beam quark distribution function [18]. In all these calculations, the main goal is to investigate the effect of parton energy loss in nuclear medium from either final state or initial state multiple scatterings. Therefore, only a subset of Feynamn diagrams at next-to-leading order (NLO) are considered to simplify the calculation by choosing appropriate physical gauge for the radiated gluon. However, without the inclusion of the complete NLO Feynman diagrams, a consistent check of gauge invariance is impossible.
The first complete NLO calculation at twist 4 has been carried out for the transverse momentum weighted differential cross section in SIDIS [19,20], and later extended to Drell-Yan lepton pair production in proton-nucleus collisions [21]. The calculated observable is directly related to the transverse momentum broadening in SIDIS off nuclear targets and heavy-ion collisions. In these two calculations, the authors have included all Feynman diagrams whose contributions are enhanced by the nuclear size. One technical aspect of these twist-4 NLO calculations is that appropriate initial partons' transverse momenta flow has to be assigned to ensure the gauge invariance during the course of collinear expansion. This is of particular importance for the subprocess of interference between soft and hard rescatterings. The calculation of symmetric subprocesses (soft-soft and hardhard double scattering) is less ambiguous; one can use either way to reach a gauge invariant result.
In this paper, we apply the newly developed twistexpansion technique to recalculate the effect of final state parton energy loss in SIDIS. In particular, we focus on the channel of charged-current interaction. This allows us to study light quark and heavy quark radiative energy loss on the same footing. We compare our results with the previous studies for both the light quark and heavy quark.
The rest of this paper is organized as follows. In Sec. II, we introduce our notations, and review the generalized factorization formalism at twist 4. In Sec. III, we present the details of our calculation at twist 4 by including both the quark-gluon double scattering and the interference between single and triple scatterings. We show how the principle of gauge invariance guides the transverse momentum flow of the four initial partons from the nuclear target. In Sec. IV, we give the final result of medium modified heavy quark fragmentation functions. In Sec. V, we illustrate numerically the importance of our gauge invariant result for heavy quark energy loss by comparing with the previous study and in the soft limit. A summary is given in Sec. VI.

II. GENERAL FRAMEWORK
For simplicity, we consider the following process of heavy quark production via the charged-current interaction in deep inelastic scattering (DIS) off a large nucleus A, where l 1 and l 2 are the momenta of the incoming lepton and the outgoing neutrino, p ¼ ½p þ ; 0; 0 ⊥ is the momentum per nucleon of the target nucleus with atomic number A, and l H is the observed final state heavy meson (H) momentum. In the channel of charged-current interaction, the momentum transfer via the exchange of a W AE boson is given by q ¼ l 1 − l 2 with the invariant mass q 2 ¼ −Q 2 and Q 2 ≪ M 2 W is assumed. Notice that both the heavy quark flavor and the momentum scale of the exchanged vector boson are labeled as Q in this work.
The differential cross section for single inclusive heavy meson production can be written as where s ¼ ðp þ l 1 Þ 2 is the lepton-nucleon collision energy and G F stands for the four-fermion coupling constant. The charged-current leptonic tensor reads The semi-inclusive hadronic tensor is defined by where P X sums over all possible final states, J μ ¼ P fψ f γ μ Vψ f is the hadronic charged current, and V ¼ ð1 − γ 5 ÞV ij with V ij stands for the Cabibbo-Kobayashi-Maskawa flavor mixing matrix [22].
The above mechanism of single inclusive heavy meson production in DIS can be illustrated in Fig. 1(a): a W AE boson is radiated from the projectile lepton and collides with a light quark from the nuclear target, a heavy quark Q with mass M and momentum l Q is produced and then fragments into a heavy meson in the final state. In fixed order calculation, there is no collinear divergence for final state gluon radiation, because of the finite heavy quark mass that naturally serves as a regulator. However, in the case when the momentum scale Q is much larger than the heavy quark mass M, Q 2 ≫ M 2 , one encounters large logarithms, lnðQ 2 =M 2 Þ, which spoil the perturbative convergence. In this case, an all order resummation of such large logarithms has to be performed in the calculation. Such a resummation is normally done by solving the renormalization equation of final state heavy quark fragmentation functions [23]. Therefore, the hadronic tensor at leading twist can be factorized into the convolution of nuclear quark distributions f A q , heavy flavor fragmentation functions D Q→H , and the partonic cross section H μν W AE þq→QþX , where the superscript S stands for single scattering, the sum over initial state quark flavors is suppressed for simplicity, and μ I and μ F represent initial and final state factorization scale, respectively. At leading twist, the hard partonic part at lowest order is H ð0Þ μν ðx; p; q; MÞ where the Lorentz invariant variables are defined as Inside a large nucleus, the propagating heavy quark encounters additional scatterings with the nuclear target remnants as shown in Fig. 1(b), which lead to nontrivial medium modifications to heavy quark production. In this paper, we focus on the radiative energy loss due to the medium induced gluon radiations with 4-momenta l. This effect in general is a nuclear enhanced power correction to final state heavy quark fragmentation functions [16,17]. Such a power correction can be computed within the hightwist expansion formalism developed by Qiu and Sterman [13,14]. Within this framework, the collinear QCD dynamics of multiple parton interaction are contained in the medium modified splitting functions that are perturbatively calculable, while the medium property is contained in the high-twist nonperturbative multiparton correlation functions. Recently, the QCD evolution equation for the renormalized twist-4 quark-gluon correlation function was derived [19]. In this paper, we aim to derive the perturbative medium modified splitting functions through NLO computations.
We apply the recently improved twist-4 collinear factorization technique [20] to calculate the medium induced gluon radiation of the heavy quark in DIS. The leading contribution from double scattering processes can be obtained by taking collinear expansion of the hard partonic cross section with respect to the transverse momenta of initial partons where the superscript D stands for double scattering, and k 2T and k 3T are the relative transverse momenta carried by gluons from the nucleus in the double scattering (see the next section for their definitions). The hard partonic part of central-cut diagrams can be written in the following general form: We apply collinear approximation to simplify the evaluation of the trace of γ matrices, After integrating out x; x 1 ; x 2 , and l AE with the help of contour integration and δ-functions from final state phase space, the partonic hard part at twist-4 can be factorized into the product of leading-order hard part for V þ quark H ð0Þ μν ðx; p; q; MÞ at leading twist as shown in Eq. (6) and the heavy-quark-gluon rescattering partH D , Hence in the following we only show the rescattering partH D . Notice that we have verified that the same result can be obtained by going beyond the collinear approximation, i.e., through exact calculation by contracting the hadronic tensor W μν with the so-called "metric" contribution proportional to g μν in the charged-current leptonic tensor L μν in Eq. (3) [20].

III. MEDIUM INDUCED GLUON RADIATION AT TWIST 4
In this section we present the details of the calculation of the twist-4 contribution to heavy meson production in DIS. A complete twist-4 calculation at NLO contains contributions involving both quark-gluon and gluon-gluon twist-4 matrix elements. In this paper, we neglect contribution from gluon-gluon double scattering as it is irrelevant to the power correction to final state heavy quark fragmentation, which is the main focus of this paper. There are four different kinds of subprocesses in quark-gluon double scattering with real gluon radiation in the final state: soft-soft double scattering, hard-hard double scattering, and the interferences between them, referred as soft-hard and hard-soft double scattering. 1 This classification is based on whether the momenta of the exchanged gluons (i.e., k g and k 0 g on the amplitude and complex conjugate of the amplitude, respectively) have 0 or finite momenta in the collinear limit k 2T ¼ k 3T ¼ 0.
The computation of heavy quark energy loss in DIS has already been performed in Ref. [17] by setting the same transverse momentum for the rescattered gluons from the nucleus. This setting has been shown to be valid, in the case of light hadron production, for symmetric double scattering subprocesses, such as hard-hard and soft-soft double scatterings [20]. We have verified explicitly that this statement also holds for heavy flavor meson production. Therefore we do not present exhaustive calculations of the subprocesses for hard-hard and soft-soft double scattering, but instead concentrate on the asymmetric subprocesses, i.e., the interference between soft and hard rescatterings. In order to ensure gauge invariant results at twist-4, on-shell conditions for the initial partons (associated with the 2 → 2 hard scattering) have to be satisfied up to Oðk 2T Þ and Oðk 3T Þ as specified in Ref. [20] for light hadron production, which is also true for heavy meson production as explained below. This requirement is fulfilled, as proposed in Ref. [20], by a delicate assignment of the transverse momenta for the initial partons from the nuclear target. We apply the same assignment in the calculation of interference subprocesses for heavy meson production in DIS.
To simplify the calculation by means of less Feynman diagram evaluation, we use the following polarization tensor for the final state radiated gluon in light-cone gauge (n · A ¼ 0), where a particular reference vector n ¼ ½1; 0 − ; ⃗ 0 ⊥ has to be implemented. In principle, the final result is gauge independent if one considers all Feynman diagrams. In this paper, we aim to extract the large logarithmic term, lnðQ 2 =M 2 Þ, which is only related to final state gluon radiation. Thus we only consider Feynman diagrams with gluon radiated from the final state heavy quark. This simplification is realized by choosing an appropriate light-cone reference vector as specified above. With this particular gauge choice, all Feynman diagrams with gluon radiated from the initial state light quark, as well as the interference between initial and final state gluon radiations, contribute to collinear divergence (in q 2 T ≪ Q 2 limit) related to the renormalization of initial state multiparton correlation function, which is not the focus of this paper.

A. Heavy-quark-gluon double scattering
With the light-cone gauge choice for the radiated gluon as shown in Eq. (12), half of the diagrams that we need to consider in the process of interference between soft and hard gluon rescattering are shown in Fig. 2. All the relevant double scattering diagrams can be found in Refs. [10,17]. The soft rescattering subprocess at amplitude level can be viewed as two factorized scatterings: the first V þ q → Q þ g hard scattering, and the second soft rescattering Q þ g → Q by exchanging a soft gluon. To ensure gauge invariance, the initial quark associated with the first scattering is required to be on its mass shell, i.e., k 2 1 ¼ 0 [24], while no on-shell requirement is needed for the initial gluon in the second soft scattering since it is essentially a leading-order QCD interaction. On the other hand, the hard rescattering subprocess can be factorized as the first vector boson-quark hard scattering V þ q → Q and the second hard rescattering Q þ g → Q þ g. In this case, to ensure gauge invariance for the second Q þ g → Q þ g process, the exchanged gluon is required to be on its mass shell (k 2 g ¼ 0) up to the order in which we perform collinear expansion, and no requirement is needed for the initial quark associated with the first scattering. In addition, to avoid complications from the "minus" components of the momenta for both exchanged gluons in the amplitude and its complex conjugate, we could simply set the momenta of the two exchanged gluons as k g ¼ x 3 p þ k 3T and k 0 g ¼ x 2 p þ k 2T . To preserve momentum conservation on both sides of the cut line, we have to assign a transverse momentum of k 3T − k 2T to the initial quark associated with the leading-order V þ q → Q hard scattering.
We take soft-hard double scattering as an example to outline the essential steps in the calculation of twist-4 contributions, and all the other subprocesses could be evaluated in the same manner. In this subprocess as shown in Fig. 2, based on the analysis specified above, we can assign the momenta for the initial partons as Fig. 1(b)]. The two on-shell propagators labeled by the short bar and cross mark that lead to soft and hard initial gluons can be expressed as follows: On the other hand, the on-shell condition for the final state quark gives where x M and x B are defined in Eq. (7), and One can then integrate over x; x 1 , and x 2 , where two of the integrations are carried out by contour integrations, and the last one is fixed by the δ function from final state on-shell condition. Then the initial state partons' momentum fractions are fixed as follows, x ¼ x L þ x B þ x M =z; It is clear that in the collinear limit k 2T ¼ k 3T ¼ 0, the momentum fraction for the initial gluon on the left-and right-hand side of the cut line becomes 0 and remains finite, respectively. This explains why we refer to this process as soft-hard double scattering.
The key point in high-twist calculation is to perform collinear expansion. With the fixed momentum fractions for softhard double scattering as shown in Eq. (18), we have Notice that we have neglected contributions that are power suppressed, such as the derivative terms on the twist-4 quarkgluon correlation functions, and that all the momentum fractions (x's) appearing in the phase factor here and below have been shifted to x B þ x M as Eq. (6) requires. Now, let us compare our result to the one in Ref. [17], which is derived from the naive setup k 2T ¼ k 3T . 2 The difference from soft-hard double scattering is As one can see the difference is proportional to ð1 − zÞ 2 M 2 , which vanishes in the soft limit z → 1 or massless limit M ¼ 0.
Following the same logic, one can assign initial state parton momenta in the process of hard-soft double scattering. The final result can be obtained from the soft-hard double scattering via the replacement e −iðx L þð1−zÞx M =zÞp þ ðy − −y − 1 Þ → e −iðx L þð1−zÞx M =zÞp þ y − 2 in Eq. (20). Thus the difference to Ref. [17] in hard-soft double scattering is which, again, vanishes in the soft or massless limit. For soft-soft and hard-hard double scattering subprocesses, one can use the same assignments of the initial partons' transverse momentum flow as in Ref. [17], and we obtain the same results, namely, For hard-hard double scattering, we have checked that three different settings of the transverse momenta for initial state partons reach exactly the same result. By combining all contributions from quark-gluon double scattering with a central cut together, we obtain the final result where 2 There are some typos in Ref. [17]; the correct final results can be obtained from Eqs. (20)- (25) and (28)-(31).
The first two terms on the right-hand side of Eq. (25) represent the contributions from the soft-soft and hard-hard subprocesses, respectively, which remain the same as in Ref. [17]. The third term represents the contributions from their interferences, which give different quark mass dependence as compared to previous calculations in Ref. [17]. It is worth noting that a factor 1=2 is needed when they are compared with Eq. (26) in Ref. [17].

B. Interference from single and triple scatterings
To complete the calculation, we also need to consider the asymmetric-cut diagrams (left cut and right cut), which represent interferences between single and triple scatterings. All the possible interference diagrams can be found in Ref. [10]. We can obtain the rescattering partH D RðLÞ of all those asymmetric-cut diagrams in the guidance of the gauge invariance as demonstrated above. The calculation techniques follow exactly the same as those in the double scattering process presented in the previous subsection. Thus we neglect the details and list the final results below, Comparing to the results in Ref. [17], one can see that the gauge invariant collinear expansion also leads to additional terms, Again, the new correction terms vanish in the soft or massless limit.
− zÞ þx M if we choose the factorization scale μ 2 ¼ Q 2 , and the constraint ofx L integration comes from the requirement of x L < 1 as in Ref. [17].
To evaluate heavy quark energy loss numerically, we choose charm quark mass M ¼ 1.5 GeV and x A ¼ 0.04 for a nucleus with a radius R A ¼ 5 fm. In the left panel of Fig. 3, we show in the red solid curve the relative difference between the new result in this study and that in Ref. [17] as a function of Q 2 with fixed x B ¼ 0.1. One can see that the new correction term due to consideration of gauge invariance leads to significant additional contribution to the heavy quark energy loss in the small Q 2 region. However, the difference becomes negligible in the large Q 2 region. This is understandable from Eq. (44). The difference is proportional to M 2 , and therefore suppressed in the large Q 2 region like any other higher-twist effect. The difference between our new result and the commonly used soft gluon radiation limit is shown in the blue dashed curve, which is as large as 16% in the large Q 2 region. It is therefore important to consider contributions beyond the soft limit in more precise phenomenological studies of heavy quark energy loss in heavy-ion collisions. Shown in the right panel of Fig. 3 is the relative difference between our new results and the previous calculation in Ref. [17] as a function of x B for fixed Q 2 ¼ 10 GeV 2 . The difference becomes significant for large values of x B (small initial heavy quark energy) as shown in the red solid curve. The contribution beyond the soft gluon limit is also appreciable at large x B (small heavy quark energy) as shown in the blue dashed curve.
In order to discuss the difference between heavy and light quark energy loss, we show the ratio of charm quark and light quark energy loss R ¼ hΔz Q g iðx B ; μ 2 Þ= hΔz q g iðx B ; μ 2 Þ in Fig. 4, where the light quark energy loss hΔz q g iðx B ; μ 2 Þ can be obtained by setting M ¼ 0 in Eq. (45). In Fig. 4, we show the dependence of R with Q 2 for fixed x B ¼ 0.1 (left panel) and with x B for fixed Q 2 ¼ 10 GeV 2 (right panel). We observe the reduction of heavy quark energy loss due to the effect of the dead cone in our new result (red solid curve) as compared to that in Ref. [17] (blue dashed curve). Such a reduction is due to the The Q 2 and x B dependence of the relative correction of charm quark energy loss δhΔz Q g i=hΔz Q g i as compared to that from old work [17] (red solid) and with soft gluon approximation (blue dashed). heavy quark mass and therefore should disappear at high Q 2 and large initial quark energy (small x B ). Again, contributions from beyond the soft gluon limit are shown to be important by comparing the red solid curve (our new result) and the green dashed-dotted curve [soft gluon radiation limit, which is also needed to be taken in light quark energy loss hΔz q g iðx B ; μ 2 Þ at the same time].

VI. SUMMARY
In this paper, we revisited a series of studies [9][10][11]16,17] on quark energy loss induced by multiple parton scattering in DIS off a nuclear target, using the recently improved framework of the generalized factorization formalism for twist-4 processes. By performing the gauge invariant collinear expansion, we found that the light quark energy loss is not affected, but new correction terms arise for heavy quark energy loss beyond the soft gluon limit. The correction terms come from the interference between soft and hard rescatterings. In the soft gluon limit, the new correction terms can be safely neglected since they are proportional to ð1 − zÞ 2 . This validates the phenomenological implementations of heavy quark energy loss from high-twist calculations in heavy-ion collisions as have been presented in Refs. [25][26][27].
To demonstrate the significance of the correction terms, we evaluated numerically the heavy quark energy loss and compared with that in Ref. [17]. We found significant correction in the small Q 2 and large x B (small heavy quark energy) regions. Our new result was also compared with that with soft gluon approximation. The noticeable difference between these two suggests the importance of implementing the complete result (beyond soft limit) for more precise calculation of heavy flavor jet quenching in heavy-ion collisions. This also has phenomenological consequences in precise extraction of the jet transport coefficient from light and heavy flavor data in heavy-ion collisions.