TMD gluon density in nuclei versus experimental data on heavy flavor production at LHC

Analytical expressions for the Transverse Momentum Dependent (TMD, or unintegrated) gluon and sea quark densities in nuclei are derived at leading order of QCD running coupling. The calculations are performed in the framework of the rescaling model and Kimber-Martin-Ryskin (KMR) prescription, where the Bessel-inspired behavior of parton densities at small Bjorken $x$ values, obtained in the case of flat initial conditions in the double scaling QCD approximation, is applied. The derived expressions are used to evaluate the inclusive heavy flavor production in proton-lead collisions at the LHC. We find a good agreement of our results with latest experimental data collected by the CMS and ALICE Collaborations at $\sqrt s = 5.02$ GeV.


Introduction
It is well known that the usual concept of QCD factorization for proton-proton (pp) interactions assumes that the corresponding inclusive production cross section is calculated as a convolution of short-distance partonic cross sections and parton (quark or gluon) distribution functions in a proton (PDFs).Suggesting that there is no hot medium formed in proton-nucleus (pA) collisions, this concept is often extrapolated to pA interactions by replacing usual PDFs with nuclear PDFs (nPDFs) while keeping hard scattering cross sections the same (see, for example, [1][2][3][4]).However, some additional phenomena, so called cold nuclear matter effects 1 , can affect this picture (see [5][6][7][8] and references therein).Detailed knowledge of nPDFs, in particular, gluon distribution in nuclei, is necessary for theoretical description of pA processes studied at modern (LHC, RHIC) and future colliders (FCC-he, EiC, EicC, NICA).Moreover, it is important in discriminating the initial nuclear effects from the subsequent hot medium effects essential for more complex nucleusnucleus (AA) collisions, where nuclear matter can reach extremely high energy densities and temperatures, transforming into quark-gluon plasma (see also review [9] for more information).In this sence, production of charm and beauty flavors in pA interactions is of particular interest, because these processes directly probe the gluon distributions in the colliding particles.
Experimental investigation of deep inelastic scattering of leptons on nuclei performed by the European Muon Collaboration reveals the appearance of a significant nuclear effect [10], which excludes the naive idea of the nucleus as a system of quasi-free nucleons (see also [11,12] for review).From theoretical point of view, there are two main scenarios to determine the nPDFs.In the first of them, nPDFs at some initial (or starting) scale µ2 0 are extracted from global fits to nuclear data using empirical parametrization of their shape and normalizations (see, for example, [1][2][3][4] and references therein).Then, numerical solution of Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) equations [13] is employed to describe the corresponding QCD evolution (dependence on the scale µ 2 ).An alternative approach is based on some nPDF models (see [14] for more information).
Here we follow the rescaling model [15] based on the assumption that the effective size of gluon and quark confinement in the nucleus is greater than in the free nucleon [16].Within the perturbative QCD, this confinement rescaling predicts that PDFs and nPDFs can be connected by scaling the argument µ 2 (see also review [17]).So that, the rescaling model demonstrates the features inherent in both approaches: there are certain relationships between conventional and nuclear PDFs that arise as a result of shifting the values of the scale µ 2 and, at the same time, both densities obey DGLAP equations.Initially, the rescaling model was proposed for the domain of valence quarks dominance, 0.2 ≤ x ≤ 0.8, where x is the Bjorken variable.Later it was extended to small x [18,19] (see also a short review [20]), where certain shadowing and antishadowing effects were found for gluon and sea quark distributions 2 .
The main goal of our study is to derive analytical expressions (at leading order in the QCD coupling α s ) for the Transverse Momentum Dependent (TMD, or unintegrated) gluon and sea quark densities in nuclei (nTMDs) using the rescaling model [15] and popular Kimber-Martin-Ryskin (KMR) formalism [23] (see also [24]).These quantities encode nonperturbative information on hadron structure, including transverse momentum and polarization degrees of freedom.Currently they are widely used in a number of applications to topical issues in high energy physics phenomenology, especially for multi-scale or non-inclusive collider observables (see, for example, review [25] for more information).Next, the calculated TMD gluon density will be applied to investigate the heavy flavor (charm and beauty) production in proton-lead collisions at √ s NN = 5.02 TeV.Such processes are known to be sensitive to the gluon content of the nucleus and provide us with possibility to reconstruct the full map of the latter [26].Of course, they are great of importance to first test of the derived expressions.We will use the k T -factorization [27], or high-energy factorization [28] approach, which turns to be a convenient alternative to explicit higher-order pQCD calculations.In fact, a large piece of NLO + NNLO + ... corrections important at high energies can be effectively taken into account in the form of the TMD gluon density [25].Our predictions are compared with latest experimental data taken by the CMS [29][30][31] and ALICE [32] Collaborations at the LHC.The consideration below extends and continues the line of our previous studies [21,22,33,34].
The outline of our paper is following.In Section 2 we briefly describe our theoretical framework and basic steps of our calculation.In Section 3 we present the numerical results and discussion.Section 4 sums up our conclusions.

The model
This section provides a short description of the calculation steps for the TMD gluon and sea quark densities in a proton and nuclei and a brief review of the k T -factorization formulas for heavy flavor production.

Conventional (collinear) PDFs in a proton
As it was argued [35], the HERA small-x data can be well interpreted in terms of so-called double asymptotic scaling approximation, which is related to the asymptotic behaviour of DGLAP evolution (see also [36][37][38]).In this approximation, flat initial conditions for conventional gluon and sea quark densities in a proton f a (x, µ 2 ) at some scale Q 2 0 could be used: , where a = q or g [36,37,39].At leading order (LO) of perturbative QCD, small-x expressions for f a (x, µ 2 ) read3 : where and I ν (σ) and Ĩν (σ) are combinations of modified Bessel functions (at s ≥ 0, i.e. µ 2 ≥ Q 2 0 ) and usual Bessel functions (at s < 0, i.e. µ 2 < Q 2 0 ): Here N f = 4 is the number of active (massless) quark flavors and β 0 = 11 − 2N f /3 is the first coefficient of the QCD β-function in the MS-scheme.The parameters A a and Q 2 0 have been extracted [39] from a fit to HERA data on the proton structure function

Rescaling model and nuclear PDFs
In the rescaling model [15], the structure function F A 2 (x, Q 2 ) and, consequently, the valence part f A V (x, µ 2 ) of the quark density in a nucleus A are modified at intermediate and large x values, 0.2 ≤ x ≤ 0.8, as follows where the new scale µ 2 A,V is related to µ 2 by [18] So that, kernel modification of the main variable s defined in (2) depends on the µ 2independent parameter δ A V having small values [18].Next, since rise of sea quark and gluon densities increases with increasing values of µ 2 , the small-x PDF asymptotics (1) were applied [18] to the small x region of the EMC effect.In fact, in the case of nuclei, the evolution scale is less than µ 2 , that can directly reproduce the shadowing effect observed in global fits.So, one can assume that where f ± a (x, µ 2 ) are given by (1).Thus, there are two free parameters µ 2 A,± which should be determined from the analysis of experimental data for the EMC effect at low x.Corresponding values of s A ± turned out to be where δ A ± can be presented as [18] −δ and In particular, for 208 Pb we have δ Pb + = −0.34 and δ Pb − = −0.78.

Kimber-Martin-Ryskin approach
The KMR approach is a formalism to construct the TMD gluon and quark distributions from conventional (collinear) PDFs.The key assumption here is that the transverse momentum dependence of the parton densities enters only at the last of QCD evolution (namely, DGLAP).The KMR procedure is believed to take into account effectively the main part of next-to-leading logarithmic (NLL) terms α n s ln n−1 µ 2 /Λ 2 QCD compared to the leading logarithmic approximation (LLA), where terms proportional to α 2 s ln n µ 2 /Λ 2 QCD are taken into account.
In the integral formulation of KMR approach, the TMD gluon and quark densities at the leading order4 of α s can be written as [23] where D a (x, µ 2 ) = f a (x, µ 2 )/x are the conventional PDFs in a proton, P aa ′ (z) are the usual unregulated leading order DGLAP splitting functions and a, a ′ = q or g.The Sudakov form factors T a (µ 2 , k 2 T ) enable one to include logarithmic loop corrections and have the following form The cut-off parameter ∆(k 2 T ) = |k T |/(µ + |k T |) imply the angular-ordering constraint5 specifically to the last evolution step to regulate soft gluon singularities.Following [39], everywhere below we use the phenomenological infrared modification of QCD coupling which effectively increases its argument at small scales, namely, , where m ρ is the ρ meson mass ('freezing' treatment), see [40] and discussions therein.

Analytical expressions for TMDs and nTMDs
Using conventional PDFs given by (1) as an input for KMR procedure, one can derive the analytical expressions for TMD gluon and quark densities in a proton.After some algebra we have [33] where and Figure 1: The TMD gluon densities in a proton (solid curves) and nuclei (dashed curves) calculated as a function of k 2 T for different values of x and µ 2 .Note that results for µ 2 = 100 GeV 2 are multiplied by a factor of 100.The ratios f A g (x, k 2 T , µ 2 )/f g (x, k 2 T , µ 2 ) are presented also. Here T ) corresponds to the angular ordering constraint and 'frozen' treatment of QCD coupling is implied.
Following [33,34], the analytic expression for TMD parton densities (12) has to be modified at large x in the form: that is in agreement with similar modifications of conventional PDFs (see [18,19] and references therein).The value of β a (0) can be estimated from the quark counting rules or extracted from the data.Below we set β g (0) = 3.03, which was derived [33] from the best description of LHC data on inclusive b-jet production in pp collisions 6 .Performing calculations in a similar way, for TMD gluon and sea quark densities in nuclei we have where s → s A ± (k 2 T ) and k 2 A,± could be easily derived from The TMD gluon densities in a proton and nuclei calculated according to ( 12) and ( 16) with the 'frozen' scenario for strong coupling and appropriate treatment of β g (0) are shown in Fig. 1 as a function of gluon transverse momentum k 2 T for different values of x and hard scale µ 2 .

Cross section of heavy flavor production
In the k T -factorization approach, the heavy flavor production in pp or pp collisions is dominated by the direct leading-order off-shell gluon-gluon fusion subprocess where Q = c or b and four-momenta of corresponding particles are given in the parentheses.The contribution from the quark induced subprocesses is of almost no importance due to comparatively low quark densities.Corresponding cross section is calculated as a convolution of the off-shell (dependent on non-zero virtualities of the incoming gluons) partonic cross section and TMD gluon distribution in a proton [27,28].As it was already mentioned above, here we extend the consideration into pA collisions by employing the master factorization formula: where the initial off-shell gluons have fractions x 1 and x 2 of the parent proton and nucleus longitudinal momenta and azimuthal angles ϕ 1 and ϕ 2 .As usual, the off-shell partonic cross section reads The analytic expression for the off-shell gluon-gluon fusion amplitude | Ā2 (g * g * → Q Q)| is known for quite a long time (see, for instance, [27,28,42]).Below we use it with the derived formulas ( 12) and ( 16) for TMD gluon densities in a proton and nuclei, respectively.In all other respect our calculation is generally identical to that performed previously [43].

Numerical results and discussion
We are now in a position to present our numerical results.First we describe our input and the kinematic conditions.After we fixed the TMD gluon densities in a proton and nuclei, the cross section (19) depends on the renormalization and factorization scales, µ R and µ F .We take them to be equal to transverse mass of the leading produced heavy quark, , where the ξ parameter is altered from 1/2 to 2 around its default value of 1 to estimate theoretical uncertainties of our calculations.The quark masses are taken as m c = 1.4 GeV and m b = 4.75 GeV.The calculations were made with the one-loop formula for QCD coupling α s with Λ QCD = 143 MeV, Q 2 0 = 0.43 GeV 2 and A g = 0.77, A q = 0.99 [39,44] for 'frozen' treatment of α s .In case of heavy jet production, we fully associate the produced heavy quark with the final jet.
We start from charm and beauty jet production at proton-lead collisions at √ s NN = 5.02 TeV.The CMS data on c-jet transverse momentum spectra [31] refer to the pseudorapidity region (in nucleon-nucleon center-of-mass7 frame) of |η CM | < 1.5.The data on b-jet production [30] were measured at four different pseudorapidity differences: −2.The uncertainty band corresponds to estimates made as described in the text.Experimental data are from [30,31].
these dataset correspond to jet transverse momentum 55 < p T < 400 GeV.Our predictions are shown in Fig. 2 in comparison with the CMS data [30,31].The shaded bands (which, in fact, are quite narrow) represent our theoretical uncertainties estimated as discussed above.One can see that overall description of the data is reasonable good in all the pseudorapidity regions.The shape and absolute normalization of both measured charm and beauty jets are reproduced well, except, may be, last p T bin (p T ∼ 300 GeV), where effects of parton showers and/or hadronization can play a role8 .
To test lower transverse momenta we turn to data on D and B meson production in proton-lead collisions at the same energy.So, ALICE collaboration provided us with D 0 production data [32] measured in the kinematical region 0 < p T < 24 and −0.96 < y CM < 0.04.The D 0 rapidity spectra were measured at three different regions of transverse momenta: 2 < p T < 5, 5 < p T < 8 and 8 < p T < 16 GeV, where rapidity region was extended then to −1.265 < y CM < 0.335.The CMS data on B + meson production [29] were taken at 10 < p T < 60 GeV and |y lab | < 2.4.To convert c or b quarks to D 0 or B + mesons we employ the standard Peterson fragmentation function with corresponding shape parameters ϵ c = 0.06 and ϵ b = 0.006, which are often used in the NLO pQCD calculations.Following to [45], we set branching ratios B(c → D 0 ) = 0.559 and B(b → B + ) = 0.408.The results of performed calculations are shown in Fig. 3.One can see that our predictions for all observables are consistent with the data within the theoretical and experimental uncertainties.There is only disagreement with the ALICE data on D 0 transverse momentum distribution at very low p T ≤ m(D 0 ).However, we would like to point out that in order to obtain the reliable predictions at low p T the large logarithmic terms ∼ α n s ln n m Q /p T have to be resummed.It can be done using, for example, special soft gluon resummation technique [46,47] or Collins-Soper-Sterman framework [48][49][50] (see also [25] for more information).This issue is not considered here.We note also that problematic low p T region, p T ≤ m(D 0 ), is excluded from D 0 rapidity measurements, see Fig. 3 (lower panels).Another important observable is the nuclear modification factor R Q pA , describing the nuclear medium influence on the production dynamics: where A = 208 is the number of nucleons in the Pb nucleus.Our results for R Q pA (with Q = c or b) are summarized in Table 1.All the estimated values are close to unity and our predictions agree with the experimental results within uncertainties.
Thus, our calculations show that the derived expressions for TMD gluon densities in a proton and nuclei provide a reasonably good description of LHC data for heavy flavor production in proton-lead collisions at the LHC, which is extremely sensitive to the gluon content of colliding particles.It is important for future investigations of proton-nucleus and nucleus-nucleus interactions in the TMD-based framework.Note that the consideration could be improved by taking into account Gross-Llewellyn-Smith and Gottfried sum rules for valence and nonsinglet quark distributions and momentum conservation for singlet quark and gluon densities similar to that as it was already done [44].We plan to perform such derivation in a forthcoming study.Also we plan to incorporate the calculated nTMDs into the Monte-Carlo event generator pegasus [51] that significantly extends its feasibilities.

Conclusions
We have derived analytical expressions for the TMD gluon and sea quark densities in nuclei at leading order of QCD running coupling.The calculations are performed in the framework of the rescaling model and Kimber-Martin-Ryskin prescription, where the Bessel-inspired behavior of parton densities at small Bjorken x values, obtained in the case of flat initial conditions in the double scaling QCD approximation, is applied.Then we apply the obtained expressions to inclusive heavy flavor production in proton-lead collisions at the LHC, which is known to be sensitive to the gluon content of colliding particles.The calculations were performed in the framework of the k T -factorization QCD approach and based on the dominating off-shell gluon-gluon fusion g * g * → Q Q subprocess.The analysis covers the transverse momentum and rapidity spectra of charm and beauty jets as well as D 0 and B + mesons.We find a good agreement of our results with latest experimental data collected by the CMS and ALICE Collaborations at √ s = 5.02 GeV with the theoretical and experimental uncertainties.It is important for future studies of proton-nucleus and nucleus-nucleus interactions within the TMD-based approaches.

Figure 2 :
Figure 2: Differential cross sections of c (left) and b (right) jets production as functions of transverse momenta of the leading jet measured at √ s NN = 5.02 TeV.The uncertainty band corresponds to estimates made as described in the text.Experimental data are from [30, 31].

Figure 3 :
Figure 3: Differential cross sections of D 0 (left) and B + (right) jets production as functions of transverse momenta of the meson (upper panels) and their rapidity (lower panels) measured at √ s NN = 5.02 TeV.The uncertainty band corresponds to estimates made as described in the text.Experimental data are from [29, 32].

Table 1 :
Calculated nuclear modification factors R Q pA in comparison with the LHC data.The estimated theoretical uncertainties come from the scale variation as it is described in the text.