Higher order QCD corrections to charged-lepton deep-inelastic scattering and global ﬁts of parton distributions

We study the perturbative QCD corrections to structure functions of charged-lepton deep-inelastic scattering and their impact on global ﬁts of parton distributions. We include the logarithmically enhanced terms near threshold due to soft gluon resummation in the QCD corrections at next-to-next-to-leading order. We demonstrate that this approximation is su ﬃ cient to describe the available HERA data in most parts of the kinematic region. The threshold-enhanced next-to-next-to-leading order corrections improve the agreement between predictions based on global ﬁts of the parton distribution functions and the HERA collider data even in the small- x region.

In this letter we study structure functions in deep-inelastic scattering (DIS) of charged leptons off nucleons. We focus on the production of heavy quarks, like e.g., charm, which proceeds within perturbative Quantum Chromodynamics (QCD) predominantly through boson-gluon fusion [1,2]. This is a reaction of great interest because at moderate momentum transfer Q, it provides a direct probe of the gluon content of the nucleon over a wide range of Bjorken x. In the present work, we specifically include higher order QCD corrections to DIS heavy-quark production and we wish to determine their impact on our information about the parton distribution functions (PDFs) and, especially, the gluon PDF.
Our motivation stems from the high statistics data for the charm structure function F c 2 provided by the HERA experiments [3,4], where F c 2 accounts for a large fraction (up to 30%) of the total structure function F 2 . The presently available DIS data allows for high precision extractions of PDFs in global fits [5,6]. The treatment of the charm contribution in these fits is an important issue as it can induce potentially large effects also in the PDFs of light quarks and the gluon obtained from these global fits (see e.g. the recent review [7]).
In the standard factorization approach the heavy-quark contribution to the DIS structure functions can be written as a convolution of PDFs and coefficient functions, where a = 1 + 4m 2 /Q 2 and m, e q are the heavy-quark mass and charge. The strong coupling constant at the renormalization scale µ r is denoted α s = α s (µ r ) and the PDFs for the parton of flavor i are f i (x, µ 2 f ) at the factorization scale µ f . Both scales µ f and µ r are assumed equal throughout this work, i.e. µ = µ f = µ r . The kinematical variables η and ξ in Eq. (1) are given as and the partonic center-of-mass energy s = Q 2 (z/x − 1) in Eq. (1). The coefficient functions of the hard partonic scattering process enjoy an expansion in α s as The perturbative QCD predictions for the coefficient functions to the leading order (LO) are long known [1,2]. The next-to-leading order (NLO) radiative corrections are available since more than 15 years [8]. Likewise, the massless coefficient functions for the light-quark content of the DIS structure functions have also been calculated to the next-to-next-to-leading order (NNLO) some time ago [9][10][11][12][13]. More recently, the scale evolution of the PDFs has been matching the NNLO accuracy [14,15].
Although the full heavy-quark coefficient functions at two loops are the only unknown for a complete NNLO analysis, much can already be said about these terms at present. Because the structure functions for massive quarks contain two hard scales, that is the momentum transfer Q and the heavy-quark mass m, the study of particular kinematical limits yields valuable information.
For instance, at asymptotic values Q 2 , µ 2 ≫ m 2 one may treat the heavy quark as effectively massless. As an upshot, large logarithms ln(µ 2 /m 2 ) are summed over in the PDF evolution [16]. Accordingly, the higher order coefficient functions [8] in Eq. (3) take asymptotic forms [17,18] and, together with the PDFs, require matching. This is the standard procedure when changing the description from QCD with n f light flavors and a heavy quark to a theory with n f + 1 light quarks. Thus, a so-called variable flavor number scheme (VFNS) has to describe this transition in the effective number of light flavors and, moreover, a general-mass formalism for a VFNS has to be consistent with QCD factorization, see Refs. [19,20]. However, the VFNS formalism cannot be routinely extrapolated to the region of Q ∼ m, where the power corrections due to the heavyquark mass effects break the factorization and therefore the massive quarks cannot be considered as partonic constituents of the nucleon. In this kinematical region, one strictly applies QCD with n f light flavors and one heavy quark, thus working in the so-called fixed flavor number scheme (FFNS). The VFNS can be used for the whole kinematics of the existing DIS data only if it is matched in order to provide a smooth transition to the FFNS at small Q. Details of this transition cannot be derived within the VFNS framework and are therefore subject to model assumptions (see Ref. [7] for the review of the modern state of art and history of this modeling). On the other hand, in the FFNS the large logarithms ln(µ 2 /m 2 ) are contained in the higher order corrections to the coefficient functions. Thus, the basic motivation for the use of a VFNS weakens, at least for the realistic DIS kinematics [21].
Near threshold, for s ≃ 4m 2 or equivalently η ≪ 1, higher order perturbative corrections are much enhanced. The coefficient functions in Eq. (3) exhibit large double logarithms α s l ln 2l β, with β = 1 − 4m 2 /s being the velocity of the heavy quark and these Sudakov logarithms can be resummed to all orders in perturbation theory. Currently this has been achieved to the next-to-leading logarithmic (NLL) accuracy and the resummed result can be employed to generate approximate results at NNLO in QCD (see e.g. [22]). In the present work, we focus on the impact of these approximate NNLO corrections due to soft gluons and assess their impact on global fits of PDFs in the FFNS. To that end, it is instructive to express the convolution in Eq. (1) as an integration over the partonic variable η. In this way, we obtain with z(η) = x(1 + 4(η + 1)/ξ) where the integration is bounded by η max = ξ/4(1/x − 1) − 1.
In Fig. 1 we plot the shape of the gluon distribution g(z) ≡ f g (z) for representative kinematics of charm quark production at the HERA collider. The typical values of x and Q employed for xg[z(η)] in Fig. 1 are correlated: The minimal (maximal) value of x corresponds to the minimal (maximal) value of Q. For small scales Q the value of xg[z(η)] is suppressed at large η due to its argument, because z(η) rises with η about linearly (see Eq. (4)) and the gluon PDF decreases with rising argument. Therefore, for small scales Q the region of η < ∼ 1 provides the dominant contribution to the charm structure function F c 2 and the parton kinematics relevant for the coefficient functions are effectively constrained to the threshold region (see Ref. [23]). Recall that F c 2 is dominated by boson-gluon fusion. On the other hand at large virtualities Q the rise of z(η) with η is not so fast. As a consequence, the suppression of the large-η region due to the shape of the gluon PDF is weaker, as can be seen in Fig. 1. However, for the case of bottom-quark production the large-η suppression is stronger due to the larger quark mass. Hence, even at large scales Q the bottom structure function F b 2 is still saturated by parton processes (and coefficient functions) close to threshold. Let us now briefly summarize the threshold approximation to the hard scattering coefficient functions in Eq. (3). To that end, we are following the standard procedure for resumming Sudakov logarithms, see e.g. Refs. [24,25]. In differential kinematics for the one-particle inclusive DIS production of a heavy quark (see [8,22]) the (dominant) gluon coefficient function is given by where c (0,0) 2,g is the Born contribution. We have u 1 = s ′ + t 1 − s 4 and s ′ = s + Q 2 and In a physical interpretation s 4 denotes the additional energy carried away by soft gluon emission above the partonic threshold. At higher orders, Eq. (5) contains plus-distributions of the type α l s [ln 2l−1 (s 4 /m 2 )/s 4 ] + that give rise to the Sudakov logarithms upon integration, i.e. the well known double logarithms α s l ln 2l β of an inclusive formulation. At the differential level (oneparticle inclusive kinematics) the threshold resummation for DIS heavy-quark production has been performed to NLL accuracy in Ref. [22]. Subsequently, the resummed result has been used to generate the factors K (i) (cf. Eq. (5)) at fixed-order perturbation theory through NNLO. These factors K (i) contain the large logarithms.
In the present paper we improve the approximate NNLO results of Ref. [22] by performing a matching at one-loop and by including the NLO Coulomb corrections. This provides us the with the first three powers of Sudakov logarithms at all orders and we arrive at the following expressions through NNLO, , Our improved NNLO approximation in Eq. (9) is exact in the region of phase space s ≃ 4m 2 , where perturbative corrections receive the largest weight from the convolution with the gluon PDF (see discussion above and Refs. [22,23]). In Eq. (5) we have restricted ourselves to the case µ 2 = m 2 . While it is straightforward to allow for general choices µ 2 m 2 to logarithmic accuracy in the threshold resummation formalism one can even derive the exact µ r and µ f scale dependence through NNLO [22] with the help of renormalization group methods. Thus, the functions c (2,1) 2,g and c (2,2) 2,g in Eq.
However, for phenomenological applications the sole knowledge of ln(x)-terms is usually insufficient and additional assumption have to be supplied. Thus, a particular functional form for the coefficient function at small-x has been suggested in Ref. [27], where the subtraction term (−4) at ln(1/z) is motivated by the terms of similar size in the smallx limit of other known coefficient and splitting functions. The factor β times the polynomial (1 − ax/z) 20 suppresses large-z effects by a large power and a is given below Eq. (1).
We plot the small-x asymptotics of c (2,0) 2,g in Fig. 2 where we display Eq. (10) for the original result of Ref. [26] and the model coefficient function of Eq. (12). The variants of Eq. (12) without account of the subtraction term (−4) and the suppression factor of β(1 − ax/z) 20 are also given for comparison. The form of Eq. (12) is very sensitive to the particular choice of these terms. Due to the large-z suppression Eq. (12) vanishes at η 20, and the region of η, where this happens is defined by the power of 20 of the polynomial. Similarly, at large η the form of Eq. (12) is entirely defined by the subtraction term. The presence of subleading terms in the small-x expansion with large numerical coefficients is a well documented feature at higher orders in QCD, see e.g. the example of the three-loop gluon splitting function [15]. Moreover, precise phenomenological predictions based on small-x approximations only are extremely difficult to make, because the convolution in Eq. (4) is non-local. Thus, the small-x terms in the coefficient function are weighted by the PDFs in the large-x region (and vice versa), see [15]. For phenomenological applications our approximate NNLO result Eq. (9) for c (2,0) 2,g is added on top of the exact NLO predictions [8] and supplemented with the exact NNLO scale dependent functions c (2,1) 2,g and c (2,2) 2,g of Ref. [22]. This provides the best present estimate for, say, the nucleon structure function F c 2 for DIS electro-production of charm quarks. We investigate the impact of these NNLO (gluon induced) contributions to F c 2 on the nucleon PDFs extracted from global fits and perform a modified version of the fit of Ref. [6]. That fit is based on global data on inclusive charged-lepton DIS off nucleons supplemented by data for dimuon nucleon-nucleon production (i.e. the Drell-Yan process). Within the FFNS we take into account the NNLO corrections to the QCD evolution [14,15] and the massless DIS and Drell-Yan coefficient functions [9][10][11][12][28][29][30]. The value of m was fixed at 1.25GeV, close to the world average and the factorization scale was selected as Q 2 + 4m 2 . For comparision we also provide NNLO results for the factorization scale at µ = Q. The latter choice naturally leads to larger deviations at small virtualities Q due to a much increased numerical value of α s (note that we have identified µ = µ f = µ r ).  [27] for c (2,0) 2,g at large-η added (dashes). The NLO calculation (dashed-dots) and the data of Ref. [4] are given for comparison.
We also perform two variants of this fit taking into account only the LO corrections of Refs. [1,2] and the NLO corrections of Ref. [8]. The predictions for F c 2 based on these three fits are compared in Fig. 3 to the ZEUS data of Ref. [4]. The latter data are not used in the fits. At the smallest values of x and Q in the plot the predictions rise monotonically with increasing orders of perturbative QCD, thus improving agreement with the data. As we discussed above, in this region the value of F c 2 is not sensitive to the coefficient functions at large η and therefore our predictions (cf. Eq. (9)) can be considered as a good approximation to the full NNLO result for F c 2 . At bigger values of x and Q, as a result of a negative contribution from c (2,1) 2,g at large η the NNLO predictions dip below the NLO ones. In this region of η the value of c (2,0) 2,g was set to zero as our choice of matching the threshold approximation to fixed order perturbation theory (see e.g. [31] for related discussions). Checking the curves of Fig. 3 at large values of x and Q one can conclude that this contribution should be positive in order to improve the agreement with the data, the particular numerical impact depending, of course, on the gluon distribution shape as one can conclude from Fig. 1. In this region (large x and Q), the slope of F c 2 appears to be distinctly flatter in Q, particularly at higher x. However, the kinematics for large values of Q is far from threshold and beyond control of our soft gluon approximation.
The small-x contribution to c (2,0) 2,g modeled in Ref. [27] affects the comparison to data at the lowest Q and at small-x only. In order to assess the impact of the small-x term of Eq. (10) quantitatively, we focus on the lowest bin Q 2 = 2 GeV and illustrate its effect in Fig. 4. At the lowest value in x we do observe a slight sensitivity on the small-x term, which in terms of η(z) corresponds to the region of larger η ≫ 1. The effect of the high-η model of Ref. [27] amounts at most to a 30%-fraction of the dominant NNLO contribution coming from the threshold region at the lowest values of x. However, in this region, the addition of the small-x contribution to c (2,0) 2,g overshoots the data while it vanishes quickly at larger x. Recall in our analysis we set the contribution of the threshold logarithms in c (2,0) 2,g to zero for η > 1. As we discussed above, the ansatz of Eq. (10) has an inherent model uncertainty of 100% since it is driven by ad hoc parameters. Therefore it cannot be used in quantitative comparisons.  Figure 5: The same as Fig. 3 for the charm electro-production cross section data obtained by the EMC experiment [32].
At fixed-target energies the charm contribution to the inclusive sample is much smaller than at the HERA collider. This makes the experimental determination of F c 2 more difficult. The only conclusive data on fixed-target charm electro-production were obtained some time ago by the EMC collaboration [32]. These data are compared to the predictions of our fits in Fig. 5. At smallest values of Q the data are in disagreement with the predictions. Moreover the disagreement increases with the order of the perturbative QCD correction. At bigger values of Q the difference between the NLO and the NNLO predictions is marginal due to the difference in the gluon PDFs obtained in these variants of the fit. However the general agreement between data and calculations is far from ideal. Due to the limited collision energy the EMC data are sensitive to the region of η < 1 only, where our NNLO approximation in Eq. (9) should describe the exact coefficient function c (2,0) 2,g very well. Thus, we see no way to improve the agreement with the EMC data by performing a complete calculation of c (2,0) 2,g . Since the EMC data are unique it seems to be useful to have additional experimental input to clarify this disagreement. Despite the fact that we do not use data for the charm electro-production in the fit our results are sensitive to the details of the description of F c 2 . This is because F c 2 amounts to a substantial contribution, up to 30%, to the inclusive DIS structure functions at the HERA collider energies. The biggest variation in the fitted PDFs due to the NNLO corrections to F c 2 is observed for the sea quarks PDF at small-x. The latter becomes smaller in order to compensate the positive contribution from the NNLO term at small-x (see Fig. 6). The gluon distribution at x ∼ 0.02 becomes bigger by about one standard deviation as a compensation of the negative contribution of the NNLO term and the fitted value of α s (M Z ), which is anti-correlated with the gluon distribution at small-x, goes down by about 1σ. Other PDFs are essentially not affected by the corrections.
Let us summarize: We have improved the perturbative QCD predictions for the heavy-quark DIS structure functions. Our NNLO approximation takes along the first three powers of Sudakov logarithms for the boson-gluon channel γg, which, in an inclusive formulation (performing the integration and keeping the leading terms in β only) corresponds to all logarithmically enhanced terms ln k β, k = 2, . . . , 4. Moreover, we have employed the exact expressions for all scale dependent terms through NNLO [22]. Subsequently, we have applied the NLO QCD corrections [8] to the charm structure function F c 2 and our new approximate NNLO result in a global fit to data for charged-lepton DIS and dimuon production in the Drell-Yan process. Especially the use of the threshold-approximated NNLO result is legitimate because the gluon PDF constrains the parton kinematics to values around s ≃ 4m 2 . Modifying the fit of Ref. [6] in this way we have studied the effects for the determination of PDFs. We have found that our approximate NNLO result for F c 2 gives better agreement between the fitted PDFs and the HERA collider data. The agreement with data extends even down to small values of x. The results of our fit are also in good agreement with ZEUS data [4] on charm electroproduction, which was not used in the fits. Comparing with EMC data [32] we did find disagreement, though, and it would be interesting to get new and independent experimental information in order to resolve it.
On the theory side we could in principle extend the NNLO threshold approximation of Eq. (9) further to include (after integration) the linear logarithm in ln β at two loops and the two-loop Coulomb corrections following the procedure of Ref. [31] for heavy-quark hadro-production. However, we leave this to future research.