Heavy-flavor hadro-production with heavy-quark masses renormalized in the ${\overline{\rm MS}}$, MSR and on-shell schemes

We present predictions for heavy-quark production at the Large Hadron Collider making use of the \msbar and MSR renormalization schemes for the heavy-quark mass as alternatives to the widely used on-shell renormalization scheme. We compute single and double differential distributions including QCD corrections at next-to-leading order and investigate the renormalization and factorization scale dependence as well as the perturbative convergence in these mass renormalization schemes. The implementation is based on publicly available programs, ${\texttt{MCFM}}$ and ${\texttt{xFitter}}$, extending their capabilities. Our results are applied to extract the top-quark mass using measurements of the total and differential $t\bar{t}$ production cross-sections and to investigate constraints on parton distribution functions, especially on the gluon distribution at low $x$ values, from available LHC data on heavy-flavor hadro-production.

Beyond fixed-order perturbation theory, the resummation of logarithms in the ratio p T /m, relevant when the transverse momentum of the heavy-quark p T is much larger than its mass m, up to the next-to-leading-logarithmic accuracy, performed through the fragmentation function approach [12], has been combined with NLO predictions [13][14][15]. The resummation of recoil logarithms related to soft gluon radiation from initial state partons, as well as the one of threshold logarithms and high-energy logarithms, have also been worked out (up to various degrees of accuracy) and presented in a number of papers (see e.g. [16][17][18][19][20][21]). The transition from quarks to hadrons is described either by fragmentation functions (FFs) [22][23][24], or by matching NLO matrix elements to parton shower and hadronization approaches [25][26][27][28].
One of the inputs of all aforementioned calculations are the values of the heavy-quark masses. The SM by itself does not predict the values of the quark masses, which are thus fundamental parameters and subject to renormalization. The ultraviolet divergences appearing in the heavy-quark self-energies, to be eliminated by renormalization, require to fix a specific renormalization scheme for relating the bare masses to the renormalized ones. The most common choice is the mass in the on-shell scheme, defining the pole mass by the relation that the inverse heavy-quark propagator with momentum p vanishes on-shell, i.e. for p 2 = m 2 , and it is known at four loops in QCD [29,30]. Another option, also known at four loops [31], is the MS prescription. In complete analogy to the renormalization of the strong coupling constant α S , only divergent terms are absorbed such that the quark propagator becomes finite after wave-function renormalization. Finally, the MSR scheme [32] defines a scale-dependent short-distance mass, that interpolates smoothly between a valid mass definition at low scales much below the mass and the MS mass for scales larger than the mass, using infrared renormalization group flow. Thus, predictions for physical observables in perturbation theory carry scheme dependence through the choice of a particular mass renormalization scheme. For a given observable, the behavior of the truncated expansion in perturbative QCD, including the apparent convergence and the residual scale dependence, can vary significantly between different schemes employed.
In this paper we describe the phenomenological effects of the use of the MS and MSR schemes for the renormalization of the heavy-quark masses in charm, bottom and top production at hadron colliders. We investigate the perturbative convergence in these schemes, by providing comparisons between physical quantities calculated at various levels of accuracy, and we discuss applications concerning the extraction of mass values and parton distribution functions (PDFs) from collider data. The implementation of the MS and MSR schemes for renormalizing the heavy quark masses, as an alternative to the on-shell scheme, is described in Sec. 2. The obtained theoretical predictions for differential distributions including NLO QCD corrections and mass renormalization in these three schemes, are presented in Sec. 3, together with considerations on the convergence of the perturbative expansion in the strong coupling constant. The results are applied in Sec. 4 to investigate the impact of LHC data on possible extractions of PDFs from collider data and on determinations of the heavy-quark mass values in different mass renormalization schemes. In Sec. 4.3, we also use predictions for total cross-sections of charm hadro-production up to NNLO accuracy in QCD. Finally, our conclusions are summarized in Sec. 5.

Implementation of heavy-quark mass renormalization schemes
In this work light quarks are assumed to be massless. For the heavy-quark masses, on the other hand, different renormalization schemes can be adopted and we briefly recall the relevant relations for the above mentioned cases of the MS, MSR and on-shell schemes. Other choices for the mass renormalization are possible. For physical observables inherently connected to the heavy-quark production threshold, for instance, the potential subtracted mass [33] was suggested as a possibility to produce an improved perturbative convergence at energies slightly above the quark-pair production threshold and the 1S mass has been presented [34] as a way of stabilizing the position of the peak of the vector-current-induced total cross-section for tt production in e + e − collisions, as a function of the center-of-mass energy √ s, for √ s ∼ 2m. In boosted kinematics, limited to the case of top-quarks, the topquark jet-mass [35], was introduced in the framework of Soft Collinear Effective Theory. Further mass renormalization schemes are described, e.g., in Refs. [36,37] and references therein.
The on-shell mass coincides with the pole in the propagator of the renormalized quark field and is known up to four loops in QCD [29,30]. Thus, it is the same at all scales and infrared finite to all orders in perturbation theory. This definition of the pole mass m pole , although being gauge invariant, has its short-comings [38,39]. It does not extend beyond perturbation theory, i.e. to full QCD, since it is based on the (unphysical) concept of colored quarks as asymptotic states. Therefore, m pole must acquire non-perturbative corrections, which leads to an intrinsic uncertainty in its definition of the order O(Λ QCD ) related to the renormalon ambiguity [40]. The latter manifests itself as a linear sensitivity to infrared momenta in Feynman diagrams, leading to poorly convergent perturbative series for the observables expressed in terms of m pole .
On the other hand, short-distance mass definitions such as the MS or the MSR schemes are renormalon-free. In general, such short-distance masses m sd are related to the pole mass through the relation m pole = m sd (R, µ R ) + δm pole−sd (R, µ R ) , (2.1) where the term δm pole−sd removes the renormalon and the dependence of the short-mass definition on long-distance aspects of QCD. Here, µ R denotes renormalization scale, connected with the ultraviolet divergences, whereas the scale R is associated with the infrared renormalization group equation (RGE) [32]. In many short-distance mass renormalization schemes, R coincides with the renormalized mass itself, as, for instance in the MS scheme, where R = m(µ R ). However, the possibility to consider other choices of R through the associated RGE allows to improve the stability of the conversion between short-distance mass schemes characterized by different values of R.
In the MS scheme, the renormalized mass of the heavy quark evolves with the RGE in the renormalization scale µ R , governed by the mass anomalous dimensions γ(α S (µ R )), where the perturbative expansion of γ(α S (µ R )) ≡ ∞ i=0 γ i (α S (µ R )/π) i+1 is known at four loops [31]. Precise determinations of the MS masses for charm-and bottom-quarks are summarized by the Particle Data Group (PDG) [22]. For the MS mass of the top-quark, see, e.g., Refs. [22,[41][42][43][44]. The conversion to the on-shell scheme proceeds in the standard manner where the first numerical coefficients c i read [45][46][47], .
Here, ζ denotes the Riemann zeta-function, L ≡ ln(µ 2 R /m(µ R ) 2 ) and the function ∆ accounts for all quarks with masses m i smaller than the heavy-quark one. As the light quarks are taken massless here, i.e., m i = 0, the ∆ term vanishes. The strong coupling α S is evaluated at the scale µ R and renormalized in the MS scheme with the number of active flavors set to n f = n lf + 1 at and above the heavy-quark threshold scale, which is assumed to be equal to its running mass. The number of light flavors is n lf = 3, 4, 5 for charm, bottom and top production, respectively.
For the particular choice m(m) ≡ m(µ R = m(µ R )), i.e. the MS mass renormalized at the specific scale µ R = m(µ R ), the logarithmic terms L cancel and eq. (2.3) evaluates numerically (up to terms O(α 4 S )) as [48] m pole = m(m) 1 + 1.333 α S π + (13.44 − 1.041 n lf ) α S π 2 + 190.595 − 27.0 n lf + 0.653 n 2 The infrared renormalon ambiguity in the conversion in eqs. (2.3), (2.6) manifests itself in practice as factorially growing terms in the perturbative expansion, that spoil convergence. The sizable coefficients in eq. (2.6) indicate the poor convergence of m pole for the case of charm and bottom, when α S at low scales is large. For top-quarks, the convergence is better due to the smaller value of α S at larger scales. Including the four-loop QCD results [31], the residual uncertainty in m pole for top-quarks, including renormalon contributions, is estimated of the order of a few hundred MeV [49], i.e., of the order of O(Λ QCD ). All available relations for scheme changes from m(µ R ) to m pole and vice versa have been summarized in the programs CRunDec [50] and RunDec [51]. While α S in eqs. (2.3), (2.6) is renormalized in the MS scheme, the matrix elements, as well as the PDFs and the associated α S evolution used in the fixed-order massive calculations presented in this paper are all defined with a fixed number of light flavors n lf = 3  Table 1. for charm and bottom production 1 and n lf = 5 for top production, even at scales well above the heavy-quark mass value. Using the standard decoupling relations, it is possible to express the PDFs, α S and the partonic cross-sections in the MS scheme with the same fixed number of flavors, once the matching scale is fixed. In this way, also eqs. (2.3), (2.6) can be re-expressed in terms of α S with the heavy degrees of freedom decoupled. If the decoupling is performed at a scale equal to the MS mass m(m), the coefficients c 1 and c 2 in eqs. (2.4), (2.5) remain identical due to the fact that the leading order coefficient in the decoupling relation for α vanishes. In practice, although the perturbative expansion in eqs. (2.3), (2.6) is known up to four loops [31], we truncate it in this work to two loops (order α 2 S ) for computing the NNLO cross-sections and to one loop (order α S ) for computing NLO cross-sections, respectively, unless stated otherwise. In addition, the evolution of α S as a function of µ R and the corresponding α S values entering in eqs. (2.3), (2.6) and other parts of the fixed-order computation are evaluated retaining three loops for producing NNLO cross-sections and two loops for producing the NLO ones, respectively, unless stated otherwise.
The MSR mass is a specific realization of the short-distance mass introduced in eq. (2.1). It is obtained, e.g., by considering the difference between m pole and m(m), see eq. (2.3), and substituting m(µ R ) with R in the terms proportional to α S to determine the difference between m pole and m MSR (R) as The use of n lf = 3 even for bottom production is justified for bottom-quark production at very low pT (see e.g. the available measurements of B-meson production by the LHCb collaboration [52][53][54] . All values are given in GeV. In the conversion formulas between the expression of masses in different renormalization schemes, we use the coupling constant α S of the effective theory including 5 active flavors in case of top and 3 active flavors in case of charm and bottom, obtained through decoupling from the theory including one additional quark, supposed to be massive. We fix α s (M Z ) n f =5 = 0.118 (α s (M Z ) n f =3 = 0.106) and we evolve α S at four loops in all cases.
where the numerical coefficients a i are given in Ref. [32]. The evolution of the MSR mass with the R scale follows the RGE 8) which is linear in the scale R and where γ MSR (α S (R)) ≡ ∞ i=0 γ MSR i (α S (R)/π) i+1 denotes the R scheme anomalous dimension. In practice, the MSR mass interpolates between the pole and the MS mass m(m). This occurs through the dependence on the scale R, because by construction m MSR (R) → m pole for R → 0 and m MSR (R) → m(m) for R → m(m).
In the following we use what has been called practical MSR mass in Ref. [55], in contrast to the natural MSR mass. For our purposes the numerical differences between those definitions are mostly negligible.
The evolution of the MS heavy-quark masses with renormalization scale is shown in the left panel of Fig. 1. It is calculated at one loop using the CRunDec program [50,51] (n lf = 3 for charm and bottom, n lf = 5 for top). The evolution of the MSR heavy-quark masses with the R scale at one loop is shown in the right panel of Fig. 1. It is obtained by solving the RGE in eq. (2.8): where expressions for the first few coefficients of the anomalous dimension γ M SR can be found in Ref. [32,55]. Here, eq. (2.9) is expanded up to the lowest non-vanishing order of α s . As is visible in Fig. 1 in the perturbative expansion of the anomalous dimension γ MSR is positive, cf. eq. (2.8).
In Table 1

Predictions for differential cross-sections
Predictions for cross-sections of heavy-quark production with different mass renormalization schemes can be obtained from those in the widely used on-shell scheme by substituting eqs. (2.3) and (2.7) in the cross-sections and performing a subsequent perturbative expansion in α s , see Refs. [57,58]. In particular, the cross-sections converted to the MS or MSR mass schemes are determined to NLO accuracy as follows Here σ 0 is the Born contribution to the cross-section (proportional to O(α 2 s )), and the differences m(µ R ) − m pole and m MSR (R) − m pole are calculated up to the lowest nonvanishing order in α s , such that all terms of order O(α 4 s ) are dropped in eq. (3.1), as appropriate for a NLO calculation at order O(α 3 s ). Formulae for scheme conversions up to NNLO haven been given in Refs. [57,58].
We are considering theory predictions for stable heavy quarks (in case of bottom and charm augmented by FFs to describe the final state B-and D-hadrons, as for the applications in Sec. 4). The additional impact of parton showers and the dependence of the quark Cross-section pole mass scheme MS mass scheme MSR mass scheme dσ/dp T M, X M, X X dσ/dy M, X M, X X dσ/dM M, X M d 2 σ/dp T dy M, X X X Table 2. Summary of the capabilities of the MCFM (M) and xFitter (X) frameworks to compute differential cross-sections for heavy-quark hadro-production in different mass schemes.
mass parameter on their cutoff [59] as well as the study of renormalon effects in observables with cuts leading to corrections of order Λ QCD in the extracted quark masses [60] are subject of ongoing theory research.
We have computed double-differential cross-sections as functions of the transverse momentum p T and rapidity y of the heavy quark Q, and single-differential cross-sections as a function of the invariant mass M QQ of the heavy-quark pair in the on-shell, MS and MSR mass renormalization schemes at NLO using both frameworks, MCFM [61,62] with modifications [58,63] and xFitter [64]. In both cases, the original NLO calculations are done in the pole mass scheme [3,4,65]. The modified MCFM program [58] is capable of converting the NLO calculations using a pole mass into those with the heavy-quark mass renormalized in the MS mass scheme, in case of single-differential distributions in p T of the heavy quark, y of the heavy quark and invariant mass M QQ of the heavy-quark pair. On the other hand, the developed xFitter framework implements the calculation of one-particle inclusive cross-sections (i.e., with the other particle integrated over), i.e., it is capable to compute the double-differential cross-sections as a function of p T and y of the heavy quark, but not as a function of the invariant mass M QQ of the heavy-quark pair. It converts the pole-mass NLO cross-sections into the MS and MSR mass schemes for fully differential distributions. The derivative of the Born contribution appearing in eq. (3.1) is calculated semi-analytically in MCFM (see [58]), whereas it is computed numerically in xFitter, which allows for cross-checks of both methods. The differential cross-sections in different schemes which can be computed by MCFM and xFitter are summarized in Table 2. For all crosssections calculated with both programs, MCFM and xFitter, (i.e. all those in the pole mass scheme and the p T and y differential distributions in the MS mass scheme), agreement within one percent accuracy is found 3 . xFitter is also interfaced to other codes, like e.g. aMC@NLO, and thus it can be used for computing cross-sections with heavy-quark masses renormalized in the MS scheme (instead of the pole scheme implemented in the standalone standard version of these codes) for a wider range of processes (e.g. tt+j hadro-production, see Sec. 4.1 for the application of this interface to a phenomenological study).             In the calculations of the differential distributions presented in the following, we use the PDG values [22] for the MS charm-and bottom-quark masses, m c (m c ) = 1.28 GeV and m b (m b ) = 4.18 GeV, and the MS top-quark mass value m t (m t ) = 163 GeV [43]. Alternatively, we use the MS charm-, bottom-and top-quark mass central values extracted from the ABMP16 NLO simultaneous fit of PDFs, α S (M Z ) and MS heavy-quark masses [66]: m c (m c ) = 1.18 GeV, m b (m b ) = 3.88 GeV and m t (m t ) = 162.1 GeV. Although these values are smaller than those quoted by the PDG, they allow for fully self-consistent computa-  tions when used in association with the ABMP16 α S values and PDFs 4 . The pole and MSR mass values are obtained from the previous ones, using a procedure analogous to that adopted for building Table 1, except that the α S (M Z ) values are now those extracted in the Table 1, and α S is evolved at two loops as in the fit. Specifically, for the MSR masses m   quarks R = 1 GeV, in order to avoid using the too small value of m MSR c at R = 3 GeV (see Fig. 1, right panel). For the pole masses, the values from the MS mass conversion at one loop are chosen. The factorization and renormalization scales µ R and µ F are set to µ 0 = 4m 2 Q + p 2 T and the proton is described by the PDF set ABMP16 at NLO. To estimate the theoretical uncertainties, the pair of factorization and renormalization scales, (µ R , µ F ), are varied by a factor of two up and down around the nominal value µ 0 , both simultaneously and independently, and excluding the combinations (0.5,2)µ 0 and (2,0.5)µ 0 , following the conventional seven-point scale variation. All calculations are provided for pp collisions at the LHC at a center-of-mass energy of √ s = 7 TeV.  In Fig. 2 the NLO differential cross-sections are shown together with their scale uncertainties as a function of p T in different intervals of the rapidity y of the charm-quark, and with the charm-quark mass renormalized in the pole, MS and MSR mass schemes. These cross-sections are computed using xFitter. The changes of the cross-sections are in the range of a few percent to ∼ 40%, when using the MS or MSR mass scheme instead of the pole mass scheme, and they are more evident in the bulk of the phase space. However, this is a small effect compared to the size of scale uncertainties at NLO. The latter amount to a factor of ∼ 2 in the bulk of the phase space, decreasing slightly towards large p T values. It turns out that the scale uncertainties are very similar in all mass schemes for variations around the chosen nominal scale µ R = µ F = 4m 2 c + p 2 T . Modifying the value of the charm-quark MS mass, which is set to the PDG value in Fig. 2, to the value extracted in the ABMP16 NLO fit, produce results qualitatively similar, shown in Fig. 3. Differences between predictions in different mass renormalization schemes in Fig. 3 are smaller than in Fig. 2, due to the fact that the ABMP16 MS masses are smaller than the PDG ones.
In Fig. 4 the same comparison of NLO differential cross-sections in the various mass renormalization schemes is presented for bottom-quark production. In this case, the impact of converting the pole mass calculations into the MS or MSR schemes vary from a few percents to 25%, which is still small compared to the NLO scale uncertainties of the order of 50%. With the choice for the nominal scale, µ R = µ F = 4m 2 b + p 2 T , the scale uncertainties are similar in the pole and MSR mass schemes, whereas they are more asymmetric and slightly smaller at low p T in the MS mass scheme. Again, modifying the value of the bottom-quark MS mass, which is set to the PDG value in Fig. 4, to the value extracted in the ABMP16 NLO fit, leads to results qualitatively similar, shown in Fig. 5, with slightly smaller differences (up to ∼ 20% among predictions in different mass renormalization schemes.
Finally, Fig. 6 displays the same comparison for top-quark production. In this case, the impact of converting the pole mass calculations into the MS mass scheme is about 20% at low p T , which is no longer small compared to the NLO scale uncertainties. It decreases towards higher p T values. When converting the cross-sections from the pole to the MSR mass scheme, the impact is below 10% and is within the NLO scale uncertainties for variations around the nominal scale µ R = µ F = 4m 2 t + p 2 T . The scale uncertainties in the MS mass scheme are slightly smaller than in the pole mass scheme, as was already reported previously [58], while the scale uncertainties in the MSR and pole mass schemes are very similar. Again, modifying the value of the top-quark MS mass, which is set to the PDG value in Fig. 6, to the value extracted in the ABMP16 fit, leads to predictions qualitatively similar, shown in Fig. 7.
In general, the differences between predictions in different mass renormalization schemes slightly increase with the rapidity, as can be seen in all Figs. 2-7.
In Figs Fig. 9, whereas for the bottom the MSR and MS mass uncertainty bands, shown in Fig. 11, are enlarged with respect to those in Fig. 10 due to the larger uncertainty accompanying the bottom-mass extraction in the ABMP16 fit (± 0.13 GeV) as compared to the PDG case (± 0.03 GeV).
In Fig. 12 the single-differential cross-sections as a function of the invariant mass M QQ of the heavy-quark pair in the pole and MS mass scheme are shown, as calculated using MCFM (no implementation of the MSR scheme is available for this distribution). The impact of changing from the pole to the MS mass scheme is largest at the lowest values of M QQ close to the production threshold. At a technical level, this is due to the derivative term in eq. (3.1) becoming dominant in this kinematic region. However, this implies that the term δm pole−sd = m pole − m sd in eq.(2.1) for the conversion of m pole to a short distance mass grows parametrically as δm pole−sd ∼ m sd α S , hence is no longer small either. This situation is realized for the MS mass definition and it persists even when changing the MS mass value, as follows from the comparison of Fig. 12 with Fig. 13, where different m(m) values are employed. This excludes the MS scheme from being a suitable mass renormalization scheme for observables very close to threshold, cf. Ref. [58] for a detailed analysis for the top-quark pair invariant mass distribution. Alternative mass renormalization schemes for observables dominated by the production threshold have been mentioned in Sec. 2.
For comparison to current experimental data on pair-invariant mass M QQ distributions in hadro-production, however, this aspect is of minor relevance. For instance, for top production at the LHC [44,67] the size of the lowest M tt bin is large, extending to O(50) GeV above threshold, so that sensitivity to threshold dynamics is significantly reduced and the MS mass scheme is still applicable in analyses of those data.
In Fig. 14 we show the impact of using different PDF sets (together with their α S (M Z ) value) in the rapidity distributions for charm-, bottom-and top-quarks (see also Ref. [63]). We fix the heavy-quark MS masses to the PDG values. Slight changes in the normalization of the distributions can be ascribed to the fact that different PDF fits are accompanied by slightly different values of α S (M z ). On the other hand, larger changes in normalization and in shapes are related to the different behaviour of different PDFs as a function of x and µ F . In particular, in case of charm production, the shape of the rapidity distribution obtained with the central set of the MMHT14 PDF fit [68] for pp collisions at √ s = 7 TeV is much wider with respect to that obtained with the central PDF sets from other widely used fits. This is particularly evident when using the MS heavy-quark mass, instead of the pole mass, in the computation, due to the lower value of the first one with respect to the second one, and is related to the peculiar and very flexible MMHT14 PDF parameterization and the particular behaviour of the gluon distribution at small x. At the scales relevant for the calculation, the MMHT14 NLO central gluon distribution steeply rises for smaller x and displays large uncertainties, in absence of data capable of constraining it for x < 10 −4 in the fit, see also Ref. [69]. On the other hand, in case of top and bottom production, the differences among predictions making use of different PDF sets are smaller than for the charm case, because, for fixed rapidity values, these processes probe larger (x, Q 2 ) values, where more data have been used to constrain the various PDFs. In particular, the predictions obtained by different PDF sets, turn out to be within the scale uncertainty band computed using the ABMP16 NLO PDF nominal set, at least for rapidities away from the far-forward region.
Additionally, in this paper we explore the possibility of using a dynamical scale in the heavy-quark MS mass renormalization, as an alternative to the static value m Q (m Q ) and its variations used in the previous distributions and in Ref [58]. There, the p T distribution of the top-quark at NLO was computed for static central scales µ R = µ F = µ m = m t (m t ), varying them simultaneously by a factor (1/2, 2) around their central value and finding that the scale uncertainty band was reduced with respect to the case when µ R and µ F are varied and µ m is fixed to m t (m t ). In general, we expect that dynamical scales, catching the different kinematics of different events, provide a more accurate description of differential distributions. Thus, in the following we consider the case when the renormalization and mass renormalization scales are chosen dynamically and coincide, i.e.
We fix the central factorization scale to the same value µ F = µ 0 . For this configuration, we compute scale uncertainties, obtained by vary- Another example of dynamical mass renormalization scale choice was shown in Ref. [72], where the MS mass renormalization scheme was studied as an alternative to the pole mass scheme for producing predictions for top-quark related distributions at NNLO. There the tt-pair invariant mass distribution was studied at NNLO, using the MS mass at a scale µ m ∼ M tt /2, setting µ R = µ F = M tt /2 and making a 15-point scale variation of factors (1/2, 2) around the (µ R , µ F , µ m ) central value. Predictions were compared to the case when µ m = m t (m t ), seeing small differences. On the other hand, p T and rapidity distributions were computed using static µ m values. To the best of our knowledge, our paper is the first work where the use of a dynamical µ m scale in computing p T distributions for heavy-quark hadro-production is investigated.
We checked that our NLO predictions are consistent with those reported in Ref. [72], when using their configuration. In Fig. 16 we present the p T distribution of the antitopquark for tt production in pp collisions at √ s = 13 TeV, using as input the NNPDF3.1 NLO PDF set with its α S (M z ) default value and α S evolution, i.e. one of the configurations already considered in Ref. [72], and multiple choices for the (µ R , µ F , µ m ) scales. For fixed µ m = m t (m t ) = 163.7 GeV, we observe that central predictions using µ 0 R = µ 0 F = p 2 T + 4m 2 t (m t ) have larger (µ R , µ F ) uncertainty bands (especially in the peak region) and have smaller absolute values than those using central scales µ R = µ F = p 2 T + m 2 t (m t ) or µ R = µ F = m t (m t ), with differences between central values at the peak amounting to ∼ 10%, as shown in the lower panel. The latter two scale choices can be considered better scale choices (i.e. scale choices leading to a faster perturbative convergence) for tt production than the first one, as also proven by the fact that NNLO corrections, reported in Ref. [72] for the (µ R , µ F ) case in comparison to the NLO ones, are quite small. On the other hand, the central predictions we obtained using µ R = µ F = p 2 T + 4m 2 t (µ R ) are in much better agreement with the previous ones than those with (µ 0 R , µ 0 F ), as shown in the upper panel, and have smaller scale uncertainty bands (not reported in the plot), which shows that the use of m t (µ R ) instead of m t (m t ) in the dynamical scale definition improves the perturbative convergence of the calculation, corresponding to smaller NNLO corrections. The predictions with (µ R , µ F ) are larger than those with (µ 0 R , µ 0 F ) because m t (µ R ) < m t (m t ) and µ R < µ 0 R . The differences at the peak of the p T distribution amount to ∆µ R ∼ −14.5 GeV and ∆m ∼ −7.4 GeV. A similar behaviour emerges when comparing the lower left and right panel of Fig. 15, for which analogous considerations and conclusions apply. On the other hand, if one uses a scale (µ R , µ F ) = p 2 T + m 2 t (µ R ), one finds central predictions only slightly larger than for the case (µ R , µ F ), as also shown in the upper panel of Fig. 16, considering that both (m t (µ R ) − m t (m t )) and (µ R − µ R ) ∼ -0.9 GeV at the peak of the p T distribution.
In summary, the heavy-quark p T -distributions in Figs. 15 and 16 with dynamical renormalization and factorization scales of the type (µ R , µ F ) = p 2 T + κ m 2 Q (µ R ) for some number κ = 1 . . . 4 and the quark masses in the MS scheme m Q (µ R ) evaluated at this dynamical scale directly incorporate the running effects of the mass parameter. If compared to sufficiently precise experimental data, this offers new and complementary ways to test the running, e.g., of the top-quark mass, cf. [44].      Table 3. The values for α S (M Z ) and the top-quark mass in different mass schemes obtained in Ref. [67] and in this work by fitting the CMS data on tt production and the HERA DIS data [73] to theoretical predictions. The fit, model (mod), parametrisation (par) and scale variation (µ) uncertainties are reported. Also the values of χ 2 are reported, as well as the partial χ 2 values per number of degrees of freedom (dof) for the tt data (χ 2 tt ) for 23 tt cross-section bins in the fit. The scale H is defined in the text.

Phenomenological applications
The use of the theory results can be illustrated with a number of applications in phenomenology, determining the strong coupling constant α S (M Z ) and values of the top-quark mass in the different renormalization schemes as well as constraints on PDFs by using available LHC data.

Extraction of m t (m t ) and m MSR
t + α S (M Z ) from differential tt cross-sections at NLO The top-quark mass can be extracted using measurements of the total or differential tt production cross-sections. As an example, we use the recent CMS measurement [67] of normalized triple-differential tt cross-sections as a function of invariant mass and rapidity of the tt pair, and the number of additional jets. These observables provide decent sensitivity to the values of m t (m t ) and m MSR t in a simultaneous fit with α s (M Z ) and the PDFs, i.e. the complete set of input theoretical parameters of fixed-order calculations for stable top-quark pair production. We compare the results with the ones obtained in the CMS analysis [67]. In particular, the distributions of the tt invariant mass and the additional jet multiplicity are sensitive to the top-quark mass through threshold and cone effects [75].  [22,[41][42][43]. The world average labelled as 'PDG2018, appr. NNLO' is based on a single determination of the D0 collaboration [74].
The QCD analysis setup follows the original CMS analysis [67], and the main settings are summarized in the following paragraph 6 . The QCD analysis is done using the xFitter framework [64]. Theoretical predictions for the tt data are obtained at NLO in the pole mass scheme using the MadGraph5 aMC@NLO program [76], interfaced with aMCfast [77] and ApplGrid [78] to store the calculated cross-sections bin-by-bin in the format which is suitable for PDF fits with xFitter. The dependence of the theoretical predictions on the top-quark mass is taken into account by generating several sets of predictions with different values of this parameter and smoothly interpolating them in the fit. The HERA combined inclusive DIS data [73] are included in the fit to provide constraints on the valence and sea quark distributions and to probe the gluon distribution and α s through scaling violations, while the CMS tt data provide direct constraints on the gluon PDF and α s , as well as on the top-quark mass as discussed in Ref. [67].
In our analysis, we convert the NLO calculations for the tt production cross-sections from the pole mass scheme into the MS or MSR mass schemes according to eq. (3.1). Due to the fact that the calculated cross-sections are stored in ApplGrid tables as bin integrated cross-sections, it is not possible to use a dynamic scale µ R = µ F = H = ( i m T,i )/2, defined as one half of the sum of transverse masses m T,i = m 2 i + p 2 T,i of the final-state partons i, since H is not constant within the bin 7 . Instead, we use a static scale µ R = µ F = m pole t , and we perform the extraction of the pole mass with this scale choice. 6 The detailed description of the fitting procedure can be found in Ref. [67], see Section 10 in particular. 7 At the best of our understanding, the publicly available version of MadGraph5 aMC@NLO is not yet capable of computing integrals over bins using a running mass, but only using a pole mass.
As the analysis of triple-differential tt cross-sections requires NLO predictions not only for inclusive tt production (N jet ≥ 0), but even for inclusive tt+1 jet production (N jet ≥ 1), MadGraph5 aMC@NLO is the only public code, among those providing such calculations, that is already interfaced to ApplGrid. In general, also other frameworks implementing NLO QCD corrections could be adopted, even beyond the fixed-order studies considered here, but they are not yet interfaced with ApplGrid.
The fit results obtained using different mass schemes are given in Table 3. The values of χ 2 characterize the fit quality. These values are very similar in all fit variants and illustrate a general good description of the tt data. To estimate uncertainties, we follow the procedure from Ref. [67] and determine fit, model, parametrization and scale variation uncertainties. As in the CMS analysis, the scales are varied coherently in all bins of the measured cross-sections. As shown in Table 3, in the pole mass scheme, switching from the dynamic scale H to the static scale m pole t modifies the extracted pole mass by about 0.6 GeV, a value still smaller than the fit uncertainties amounting to 0.7 GeV, but enlarges the scale uncertainties substantially. Therefore, the larger scale uncertainties obtained in this analysis using the MS or MSR mass schemes, as compared to Ref. [67], are explained by the usage of the static scale in the calculations. Switching from the pole mass m pole t to the MS mass m t (m t ) or the MSR mass m MSR t (3 GeV) does not affect the scale uncertainties significantly. On the other hand, if in the future one would know the value of the MS masses very precisely (from some other measurement), one could use them to get accurate predictions for differential cross-sections with smaller heavy-quark mass uncertainties, while the pole mass would be affected by O(Λ QCD ) uncertainties.
In the light of these observations, it will be worth to implement the transition to the other mass schemes directly in the MadGraph5 aMC@NLO program 8 and in further Monte Carlo integrators/event generators, such that predictions for differential tt cross-sections in association with jets can be obtained in the format which is suitable for PDF fits in different mass schemes and with dynamical scales. The advantages of the latter for the running masses have been illustrated in the previous Sec. 3.
Furthermore, we do not observe any noticeably larger theoretical uncertainty when using the MS running mass instead of the pole mass, as was reported in Refs. [42,43]. Switching from the pole mass scheme to the MSR mass scheme with R = 3 GeV changes the extracted mass value by 0.3 GeV only, which is well within the current experimental and theoretical uncertainties. On the other hand, the value of α s (M Z ) extracted from the fit does not change significantly when using different schemes, as also shown in Table 3. The obtained values are compatible with α s (M Z ) = 0.1191 ± 0.0011 extracted in the ABMP16 fit at NLO [66] within two standard deviations 9 .
The extracted value of m t (m t ) is compared with several other determinations in Fig. 17. In the ABMP16 analysis, the running top-quark mass was determined from measurements 8 At the moment the program MadGraph5 aMC@NLO does not compute directly cross-section integrals using running masses and we have developed the xFitter interface to it to convert the predictions in the pole mass scheme to the MS and MSR schemes. 9 The PDG value of αs(MZ ) = 0.1179 ± 0.0010 is based on comparisons to theory at NNLO and on lattice data.
of total top-quark pair and single-top production cross-sections in a global QCD fit at NNLO [41]. In Ref. [43] ATLAS extracted a m t (m t ) value at NLO from their measurement of tt + 1 jet production cross-sections, while Ref. [42] has obtained m t (m t ) at NLO using the ATLAS measurement of tt + 1 jet production [79] on the basis of LHC Run-1 data. Currently, the world average value of m t (m t ) by the PDG [22] is based on a single determination of this parameter by the D0 collaboration at approximate NNLO [74]. When comparing to the other determinations of m t (m t ) displayed in Fig. 17, it is worth to note that only the results of this work and of the ABMP16 analysis are obtained in a simultaneous fit of m t (m t ), α s (M Z ) and PDFs, preserving correlations among these quantities, while the other determinations were done using a value of α s (M Z ) and a PDF set fixed a-priori.
In line of principle, the applied methodology can be extended to the extraction of the m c (m c ) and m b (m b ) mass values from measurements of charm and bottom production in association with jets at colliders. However, this is a great challenge from the experimental point of view, because measuring jets at low p T , where the sensitivity to the charm-and bottom-quark mass would be particularly large, is hard.

NLO PDF fits with differential charm hadro-production cross-sections
The application of differential distributions for charm hadro-production with the MS mass definition allows for an update of PDF fits which use heavy-flavor measurements from the LHC, to constrain the gluon distribution at low values of the longitudinal momentum fraction x [80][81][82]. In particular, constraints at the lowest x values explored nowadays (x 10 −6 ) can be obtained by considering the charm hadro-production process at high rapidities (|y| 4.5) at the LHC, whereas the bottom hadro-production process at similar rapidities at the LHC is sensitive to slightly larger x values (x 10 −5 ), with a region of sensitivity that partially overlaps with the one of charm data. Because of the large scale dependence of the NLO calculations for charm hadro-production, it is customary to include in such fits only ratios of cross-sections, which are constructed using measurements at different values of rapidity and/or transverse momentum, or at different center-of-mass energies.
As an example, in the PROSA analysis [80] charm and bottom hadro-production crosssections [52,83] as a function of rapidity were used in ratios to the respective cross-section in the rapidity interval 3 < y < 3.5 for each p T bin, together with the inclusive DIS data [84] and the heavy-flavor production DIS data [85,86] from HERA. These ratios feature a reduced scale dependence, but, at the same time, they have reduced sensitivity to the heavy-quark mass. We repeat this PROSA analysis using the MS heavy-quark masses in the calculations of both the DIS structure functions [87] and the charm and bottom hadro-production cross-sections, instead of pole masses, while all other settings are as in Ref. [80]. As a result, we observe only a small impact on the χ 2 value and the fitted PDFs, with a new central PDF that is well within the previously found PDF uncertainties. These small differences are driven mainly by the change in the predictions for the HERA data, because the LHCb data used in the format of normalised cross-sections do not provide any notable sensitivity neither to the heavy-quark mass scheme, nor to the value of the  [85,86,88,89], where the HERA data alone were analyzed to determine the heavy-quark MS masses. The MS masses are also in better agreement with the world average values [22], than the pole masses of eqs. (4.3), (4.4), indicating that the latter carry a significant intrinsic theoretical uncertainty. Therefore in our most recent PDF analysis [90] we have solely adopted heavy-quark running masses.

NNLO PDF fits with total charm hadro-production cross-section
The NLO predictions for differential charm hadro-production at the LHC have very large scale uncertainties (> 100% in some phase space regions), as illustrated in Sec. 3. The lack of theory predictions for differential cross-sections on charm and bottom hadro-production at NNLO prevents including the corresponding existing data in the state-of-the-art PDF fits, which nowadays are mostly provided at NNLO accuracy. In this context measurements of the total charm hadro-production cross-section would be beneficial, because they can be confronted in the PDF fits with the already available inclusive NNLO predictions [6][7][8][9] which have significantly reduced scale uncertainties. However, no such measurements have been performed to date. On the other hand, the ALICE [91,95], ATLAS [92] and LHCb [83,93,94] experiments have provided measurements of charm production in different kinematic regions which cover more than one half of the phase space. One can reliably determine the total crosssection by extrapolating these measurements to the full phase space. The extrapolation procedure is analogous to that adopted for extracting reduced cross-sections for charm production in ep collisions at HERA [89] from measurements in a fiducial phase space. These reduced cross-sections are then routinely used in global PDF fits. In the following, we perform such extrapolations and provide the inferred values of the total cc production crosssection at different center-of-mass energies and their ratios, together with experimental and theoretical uncertainties arising from the extrapolation procedure. We then compare the results to theoretical predictions at NNLO in QCD which are obtained using different PDF sets, and demonstrate how these data can help to reduce PDF uncertainties.
The existing most precise LHC measurements of open charm production are summarized in Table 4. The ALICE measurements at √ s = 5 and 7 TeV cover the central region |y| < 0.5, the LHCb measurements at 5, 7 and 13 TeV provide coverage of the forward region 2 < y < 4.5, and the ATLAS measurement at 7 TeV essentially bridges the gap by providing data at |η| < 2.1. However, while both ALICE and LHCb provide measurements nearly in the full p T range starting from 0 GeV, ATLAS reports the cross-sections only for p T > 3.5 GeV, thus leaving the bulk of the corresponding p T kinematic range unmeasured. Furthermore, it turns out that the most precise data of ALICE and LHCb among all open D-meson data are those for D 0 production, while this final state was not measured by ATLAS.
Given these arguments, we extrapolate ALICE and LHCb measurements of D 0 production at 5 and 7 TeV to the full phase space. In order to maintain the least dependence on the theoretical predictions, both ALICE and LHCb measurements are extrapolated to nearby regions of y, namely to 0 < |y| < 1.5 and |y| > 1.5, respectively: where Here σ ALICE and σ LHCb denote the ALICE and LHCb data on fiducial cross-sections, respectively, and σ NLO in different rapidity ranges are the theoretical predictions. The factor 2 in the second term takes into account that the LHCb data are provided for only 2 < y < 4.5 and need to be extrapolated to 2 < |y| < 4.5. We exploit the symmetry around y = 0 and assume that the cross-sections at 2 < y < 4.5 and −4.5 < y < −2 are equal, as reasonably expected in case of pp collisions. Also the measurements are extrapolated into the full range of p T (not shown in eqs. (4.5), (4.6) for brevity), which implies only a 1% correction for the LHCb data at 7 TeV provided for 0 < p T < 8 GeV, and even smaller corrections for the ALICE data sets. This procedure is used to obtain the total cross-section for D 0 production at collision energies √ s = 5 and 7 TeV, while at 13 TeV we extrapolate solely the LHCb measurement since no other data are available at this energy 10 .
10 Preliminary predictions on D 0 production in pp collisions at √ s = 13 TeV were reported by the ALICE collaboration in a conference proceeding [96] in 2018, but they have neither been confirmed yet nor further refined in a regular article. In addition, the data are presented in plots, but no numerical tables are provided in Ref. [96].
We calculate the total charm production cross-section from the D 0 production crosssection dividing the latter by the fragmentation fraction from Ref. [97]: The factor 2 in eq. (4.7) accounts for the fact that both c andc fragment into charmed hadrons. We assume f (c The uncertainty on f (c → D 0 ), which amounts to 1%, is neglected. We also compute ratios of cross-sections at different center-of-mass energies R 7/5 = σ 7 TeV /σ 5 TeV and R 13/7 = σ 13 TeV /σ 7 TeV , which benefit from a partial cancellation of theoretical uncertainties [98]. The theoretical predictions σ NLO in eqs. (4.5), (4.6) are computed using the MS masses as described in the previous sections. The hard-scattering cross-sections for heavy-quark hadro-production are supplemented with phenomenological non-perturbative FFs to describe the c → D 0 transition. The factorization and renormalization scales are chosen to be µ R = µ F = µ 0 = 4m 2 c (m c ) + p 2 T and varied by a factor of two up and down (both simultaneously and independently for µ R and µ F ) to estimate the scale uncertainties with the conventional seven-point scale variation, leaving out the combinations (µ R , µ F ) = (0.5, 2)µ 0 and (2, 0.5)µ 0 . The MS charm-quark mass is set to m c (m c ) = 1.275 ± 0.030 GeV [22].
The proton is described by the PROSA PDF set [80], which is expected to have a reliable gluon distribution at low x thanks to the heavy-quark data used for its determination. Furthermore, to estimate the PDF uncertainties, the extrapolation is performed using the ABMP16 [66], CT14 [99], MMHT2014 [68], JR14 [100], NNPDF3.1 [71] and HERAPDF FF3A [73] NLO PDF sets. Then, the envelope covering the PROSA PDF uncertainties and the difference obtained using any of the additional PDF sets is constructed. This conservative procedure is essential, because the theoretical calculations for the highest y values involve the gluon PDF at the lowest x values (up to 4 × 10 −8 ), which are not directly covered by data in any of the PDF fits (not even in the PROSA fits which include the charm data up to y = 4.5 as measured by LHCb, for which PDFs at x < 10 −6 and their uncertainties are extrapolated from the results obtained up to x ∼ 10 −6 , using built-in procedures in the LHAPDF library [101]).
All theoretical uncertainties are assumed to be fully correlated for cross-sections in different kinematic regions and at different center-of-mass energies. The robustness of the extrapolation procedure is checked by varying the boundary y between the kinematic regions into which the ALICE and LHCb measurements are extrapolated by ±0.5 (at the same time, this variation tests consistency of the ALICE and LHCb data). As a further check of the method, we have computed predictions for the ALICE and LHCb data using NLO matrix elements matched, according to the Powheg method [102,103], to parton shower and hadronization implemented in PYTHIA8 [104], and found these predictions to be consistent with our NLO + FF predictions within theoretical uncertainties.
The results of the extrapolation are reported in Table 5. The scale, mass, PDF and fragmentation uncertainties are added in quadrature to obtain the total theoretical uncer-  Table 5. Extrapolated total charm production cross-sections and their ratios at different centerof-mass energies together with uncertainties from parametric variations of the scales at NLO, the mass m c (m c ) ± 0.03 GeV, α K ± 1.7, PDFs and the rapidity y ALICE,LHCb ± 0.5. The correlation factor between R 7/5 and R 13/7 is −0.61. α S uncertainties are negligible compared to the PDF ones, computed using as a baseline the CT14 PDF set of eigenvectors at NLO. tainty assigned to the extrapolated results. The experimental uncertainties of the input data are propagated to the extrapolated cross-sections and reported separately. The experimental uncertainties of the input data sets are assumed to be fully uncorrelated 11 . The experimental and theoretical extrapolation uncertainties are approximately of the same size. The total uncertainties are obtained by adding the experimental and theoretical uncertainties in quadrature, and amount to ≈ 10%. Our value for the total charm cross-section at 7 TeV is in agreement with the extrapolated cross-sections reported in Refs. [91,92,105] within uncertainties. While the central values for the extrapolation factors in eqs. (4.5), (4.6) were obtained at NLO, their uncertainties are calculated such that they should cover potential deviations from the unknown true QCD result. Therefore the resulting total cross-sections, with these uncertainties included, can be compared to calculations in any QCD scheme to any order. Furthermore, for determining these extrapolation factors, only the shape of the predictions for the p T and y differential cross-sections is relevant, while a large part of theory uncertainties related to normalization cancels.
The extrapolated cross-sections and their ratios are compared to NNLO predictions obtained using the NNLO PDF sets ABMP16 [41], CT18 [70], MMHT2014 [68], JR14 [100], NNPDF3.1 [71] and HERAPDF [73]. The cross-sections are computed using the Hathor program [106] interfaced in xFitter [64]. The factorization and renormalization scales are chosen to be µ R = µ F = 2m c (m c ) and varied by a factor of two up and down (both simultaneously and independently for µ R and µ F for the 7-point scale variation) to estimate the uncertainties. The MS charm-quark mass is set to m c (m c ) = 1.275 GeV [22].
In Fig. 18 we show the extrapolated cross-sections and their ratios compared to NNLO predictions. For the NNLO predictions, the theoretical uncertainty arising from scale  Figure 18. Comparison of the extrapolated total charm production cross-sections and their ratios with the NNLO theoretical predictions using different PDF sets. Uncertainties from scale variations at NNLO (µ) and PDFs are shown separately. variations and the PDF uncertainty are shown separately. All theoretical predictions agree with the data within uncertainties, but noticeably the MMHT2014, HERAPDF2.0 (and CT14, not plotted in the figure) PDF sets have uncertainties which are larger than both scale and data extrapolation uncertainties for some of the observables. In particular, the MMHT2014 and HERAPDF2.0 predictions for the cross-sections at √ s = 13 TeV are consistent with negative values within uncertainties (see also Ref. [69]). Predictions based on the new CT18 PDFs (and unlike those using the previous PDF set CT14) do not show anymore a large positive uncertainty which greatly exceeds the extrapolated cross-section. These PDF sets could benefit from including in their fits data on charm production crosssections or on their ratios.
Remarkably, also the scale uncertainties appear to be different when using different PDF sets. Among the considered observables, the most conclusive one is R 7/5 for which both data extrapolation and theoretical scale uncertainties are moderate (≈ 10%). Our extrapolated value for this observable can be used in future NNLO PDF fits to constrain the gluon PDF at low x. The other ratio R 13/7 has a larger extrapolation uncertainty suffering from a lack of experimental measurements of charm production in the central rapidity region at 13 TeV. We are confident that this lack will be solved by the data which will appear in forthcoming experimental studies at the LHC (see footnote 10).
As a demonstration that the provided observables can indeed constrain the PDFs, we employ a profiling technique [107] based on minimizing the χ 2 function built from data and theoretical predictions, taking into account both data and theoretical uncertainties arising from PDF variations. The analysis is performed using the xFitter program [64]. We consider the MMHT2014 PDF set at NNLO and the ratios R 7/5 , which exhibits the least scale uncertainties, and R 13/7 . The correlation of R 7/5 and R 13/7 due to the common input of 7 TeV data sets is taken into account. The PDF uncertainties are included in the  Fig. 19, but for the CT18 PDF at Q 2 = 10 GeV 2 (up left) and Q 2 = 100 GeV 2 (up right) and for the CT14 PDF at Q 2 = 10 GeV 2 (down left) and Q 2 = 100 GeV 2 (down right), all at NNLO. χ 2 functional through nuisance parameters, and the values of these nuisance parameters at the minimum are interpreted as optimised (or profiled) PDFs, while their uncertainties determined using the tolerance criterion of ∆χ 2 = 1 correspond to the new PDF uncertainties.
The original and the profiled MMHT2014 gluon PDF are shown in Fig. 19 at the scales Q 2 = 10 and 100 GeV 2 . The profiled distribution exhibits greatly reduced uncertainties at low x, and in this region the distribution is shifted towards larger values of the gluon density. In case of the MMHT2014 set, the original gluon PDF is negative at low x values, while the profiled one remains positive down to at least x ∼ 5·10 −6 , thanks to the constraint realized by adding the ratios of charm data in the PDF fit. We emphasize that the strong impact at low x is obtained as well when working with other PDF sets. As an example, in Fig. 20 the CT14 and CT18 gluon distributions are shown before and after profiling. For these sets the gluon PDF is always positive in the entire x range for all eigenvectors by construction. In case of CT14, adding the aforementioned data strongly reduces PDF uncertainties at low x, whereas the effect is milder for CT18, but still sizable at low Q 2 .

Conclusions
The hadro-production of heavy-quarks is an important class of processes at LHC. Not only for top, but also for bottom and charm, a wealth of very precise high-statistics data has been collected by the LHC collaborations, differential in the relevant kinematic variables, such as the transverse momentum p T , the rapidity y or the pair-invariant mass M QQ of the heavy-quarks (or of the respective heavy hadrons). In comparison to theory predictions in perturbative QCD, these data can be directly used for the extraction of heavy-quark masses, which are typically correlated with the value of the strong coupling constant α S (M Z ). The data also have an impact on fits of fundamental non-perturbative QCD parameters such as PDFs, where they give unique kinematic constraints.
In order to provide meaningful determinations of heavy-quark masses, the value of α S (M Z ) and PDFs, QCD predictions with good accuracy are needed. Currently theory predictions are available at NNLO for top-quark production, also for differential distributions, but not for bottom and charm. In the latter case, the predictions at NLO accuracy are generally not sufficiently precise enough considering the large theoretical uncertainties, which stem predominantly from scale variations. In view of the much smaller experimental uncertainties reached in modern analyses improvements in the theoretical descriptions are clearly required.
One such aspect, which has been studied in this paper, is the choice of a suitable renormalization scheme for the heavy quark masses. We have investigated different heavy-quark mass renormalization schemes with emphasis on the MS and MSR masses as representative short-distance mass definitions. The choice of a particular mass scheme as well as the values for the scales µ R and µ F have an impact on the rate of apparent convergence of the perturbative expansion of the cross-sections. We have investigated a range of dynamical scale choice for the cross-sections and, in case of the MS mass, also for the mass parameter m Q (µ R ). In particular, we have found that dynamical renormalization and factorization scales of the type (µ R , µ F ) p 2 T + κ m 2 Q (µ R ) for heavy-quark p T -distributions with the running of mass m Q (µ R ) included, can lead to significantly reduced residual scale uncertainties. At NLO accuracy in QCD the latter are, however, in general still large for all mass schemes, but theory predictions using MS or MSR masses carry smaller parametric uncertainties in the mass values, being theoretically well-defined and free of renormalon ambiguities.
We have demonstrated these features in extractions of the top-quark MS and MSR masses at NLO from differential distributions measured by CMS, finding good consistency with other determinations. Using differential charm hadro-production cross-sections we have also been able to improve available constraints on PDFs and, using the MS mass scheme, to decrease extrapolation uncertainties when determining total cross-section from open charm data measured in fiducial regions of phase space by the LHC collaborations. In the latter case, ratios of cross-sections are particularly useful observables to cancel residual theoretical uncertainties. In order to carry out theses studies, we have developed software frameworks using the MCFM and xFitter programs to determine differential distributions at NLO in QCD efficiently.
Avenues for theoretical improvements include the obviously needed QCD predictions for bottom and charm hadro-production at NNLO accuracy, possibly combined with the resummation of large logarithms in specific kinematics, but also an improved description of charm-and bottom-quark fragmentation to mesons, an issue which has been side-stepped in the present analysis.
The extended xFitter program, implementing the MSR and MS mass renormalization schemes, as an alternative to the on-shell scheme in heavy-quark hadro-production, is publicly available on the web, and further extensions of the MCFM and Hathor programs used to perform calculations in this paper are available from the authors upon request.