Z boson production in proton-lead collisions at the LHC accounting for transverse momenta of initial partons

We perform a calculation of inclusive $Z$ boson production in proton-lead collisions at the LHC taking into account the transverse momenta of the initial partons. We use the framework of $k_T$-factorization combining transverse momentum dependent parton distributions (TMDs) with off-shell matrix elements. In order to do it we need to construct appropriate TMDs for lead nuclei which is done using the parton branching method. Our computations are compared with data from CMS taken at $\sqrt{s}=5.02$ TeV. The results are in good agreement with the measurements especially the transverse momentum distribution of the $Z$ boson.


Introduction
The production of Z bosons in hadron-hadron collisions is described in lowest order (LO) calculations as the annihilation of a pair qq → Z. In collinear factorization, the initial quarks do not carry any intrinsic transverse momentum, and therefore the p T of the Z boson vanishes. When higher order corrections in perturbative quantum chromodynamics (pQCD) are included, the Z boson receives a significant transverse momentum corresponding to the additional emission. However, in collinear factorization the transverse momentum spectrum of the Z boson at Λ QCD < p T < O(10)GeV cannot be well described by a fixed order calculation, and a resummation to all orders of soft gluon radiation is needed. Several calculations for this soft gluon resummation exist for pp and pp collisions [1][2][3][4].
In a different approach of k T -factorization, originally developed for small x physics [5,6], the parton densities depend in addition on the partons' transverse momenta. Such transverse momenta come from the intrinsic motion of the partons inside the hadrons but also from the perturbative evolution of the partons from a small scale to the hard scale of the process. The hard process is in general calculated with off-shell initial partons. In the past years, significant progress has been made by the calculations of hard processes not only for initial gluons but also for initial quarks [7]. The transverse momentum dependent parton densities (TMDs) for protons were recently obtained from precision fits to deep-inelastic cross section measurements within the partonbranching (PB) approach [8][9][10].
In this paper we are in particular interested in exploring the transverse momentum structure of the partonic content of lead nucleus at relatively large values of its longitudinal momentum. 1 To achieve this we extend the PB approach to the case of heavy nuclei, in particular to lead nucleus, and apply the newly constructed nuclear TMDs (nTMDs) together with off-shell matrix elements to calculations of Z boson production in pPb collisions at the LHC. 2 The interest is motivated by experiments at CERN where proton-lead and leadlead collisions are studied. The precise knowledge of the partonic structure of the lead nucleus and the factorization used will allow to increase the precision of the theoretical description of the initial state of proton-lead and lead-lead collisions. The interest in knowing the precise transverse momentum structure of the nucleus lies in the following: when the quark gluon plasma is produced in lead-lead collision the propagating jet broadens and gets kicks from the plasma constituents. Therefore, its transverse spectrum changes. Moreover the initial state effects (initial state shower) lead to decorrelations. The precise knowledge of the transverse momentum distribution of partons in lead nucleus will therefore be beneficial for a more precise determination of final-state effects due to jet-plasma interactions. In order to demonstrate the usefulness of the newly obtained nTMDs for lead nucleus, we calculate the cross section for the rapidity and p T spectrum of Drell-Yan pairs with an intermediate Z/γ * boson state. Furthermore, such a final state, being a colorless particle, gives the opportunity for particularly interesting investigations complementary to results obtained in studies of jet final states in [12,13].

Nuclear TMDs
In the PB approach, the parton density is evolved with the DGLAP evolution equation from a small scale (where the initial parton density is parametrized) to the scale of the hard process using an iterative procedure. In this way, every single splitting process during the evolution is calculated, and kinematic constraints in each parton splitting step are treated. Once a physical meaning is given to the evolution scale, the transverse momentum of the partons involved in each splitting can be calculated, and a transverse momentum dependent (TMD) parton density can be obtained.
We apply DGLAP splitting functions at next-to-leading order (NLO) and we use angular ordering for the parton evolution, which relates the evolution scale µ i to the transverse momentum of the emitted parton q 2 t,i = (1−z i ) 2 µ 2 i , where z i is the fractional momentum in the splitting process (details of the procedure are described in Ref. [10]).
The nuclear TMD, A Pb a (x, k 2 t , µ 2 ), is obtained by a convolution of the starting distribution A Pb 0,b (x , k 2 t,0 , µ 2 0 ) with the evolution kernel K ba x , k 2 t,0 , k 2 t , µ 2 0 , µ 2 as described 2 Whenever we refer to Z boson production we mean a production of a lepton pair via both Z and in Ref. [10] (with k 2 t = k 2 ): The evolution of the kernel starts at x 0 = 1 at µ 2 0 . In general, the starting distribution A 0 can have flavor and x dependent k t,0 distributions, for simplicity we use here a factorized form: where the intrinsic k t,0 distribution is given by a Gauss distribution with σ 2 = q 2 0 /2 for all flavors and all x with a constant value q 0 = 0.5 GeV. The q 0 value is assumed to be the same as in the proton TMD fit [10].
The starting distribution f Pb 0,b (x, µ 2 0 ) for lead nucleus is taken to be one of the available collinear nuclear PDFs (nPDFs), e.g. the nCTEQ15 [14]. We always produce two sets of nTMDs which differ by the initial scale µ 2 0 and the argument of the strong coupling, α s . We refer to these sets as Set1 and Set2, where Set1 features µ 2 0 = 1.9 GeV 2 and the scale of α s coinciding with the DGLAP evolution variable in [10]. On the other hand Set2 has µ 2 0 = 1.4 GeV 2 and the α s scale is chosen to be the transverse momentum |q 2 t,i | of the splitting process, as suggested in Refs. [15][16][17] and also applied in [10]. This choice is motivated by accounting for coherence effects effectively introducing higher order corrections. For the purpose of practical calculations in Sec. 3 we used three different collinear nPDFs to produce the corresponding nTMDs. We adopted the two most commonly used nPDFs nCTEQ15 [14] and EPPS16 [18]. Additionally we did comparisons with distributions obtained in Ref. [19] which allow for more reliable description of the low-x region. However, the nTMDs based on these nPDFs (PB-gluon D c ncteq1568CL Pb) turned out to be rather similar to the results obtained based on the nCTEQ15 nPDFs which is why we show them only later when discussing results for Z boson production.
We now present the obtained nTMDs. In Fig. 1 we show the k t integrated distributions of the up quark at the scale µ = 10 GeV and µ = 100 GeV as a function of x. First in Figs. 1a and 1b we plot the distributions obtained form the nCTEQ15 starting distribution evolved with the PB method using different scales in α S according to the Set1 and Set2 prescriptions (PB-nCTEQ15FullNuc 208 82), and compare them with the original nCTEQ15 collinear PDFs. One can observe clear differences between the Set1 and Set2 distributions, but as expected after integration over the k t the Set1 distributions reproduce the collinear PDFs. In the computations that will follow in Sec. 3 we mostly concentrate on the Set2 distributions which are preferred by the phenomenological applications [10]. In Figs. 1c and 1d we   also compare the PB-nCTEQ15FullNuc 208 82 Set2 distributions with another Set2 TMDs obtained by using the EPPS16 lead nPDFs as starting distributions [18] (PB-EPPS16nlo CT14nl Pb208) and with the PB-NLO ptoPb208 proton TMDs 3 obtained   in Ref. [10], which we will also employ in our calculations.
In Fig. 2 the transverse momentum dependent parton distributions for different quark species (u, d,ū,d) are shown at x = 0.01 and µ = 100 GeV as a function of the transverse momentum. We show here the PB-nCTEQ15FullNuc 208 82 TMDs for Set1 and Set2 and compare them with the corresponding combination of proton TMDs PB-NLO ptoPb208. In Fig. 3 we compare different nuclear TMDs with PB-NLO ptoPb208 distributions restricting only to the Set2 case at x = 0.01 and µ = 100 GeV. We can see that both nuclear TMDs (PB-nCTEQ15FullNuc 208 82 and PB-EPPS16nlo CT14nl Pb208) are quite similar and differ from the PB-NLO ptoPb208 distributions. We can see that for small values of transverse momentum (below 1 GeV)

82/208f
p + (208 − 82)/208f n , with the neutron distribution obtained assuming isospin symmetry. there is a supression of the nTMDs compared to the proton TMDs. On the other hand, when the k t rises above 1 GeV we start to observe an enhancement of nTMDs which presists up to the high k t values (∼ 80 GeV). A similar behaviour is observed for other x values and other flavours. To even further quantify the nuclear modifications of TMDs, in Fig. 4 we provide additional distributions for different µ and k t values as a function of x. At lower k t values we can observe a behaviour similar to the collinear PDFs which exhibit a rather well defined shawoding, anti-shadowing and EMC regions.
As a last part of this section, in Fig. 5, we present the uncertainties of the obtained nTMDs. We plot them for u andū quarks at x = 0.01 and two scales µ = 10 GeV and µ = 100 GeV. We can see that the low k t region features the largest uncertainties: around 10% for up and 15% forū distributions. The size of the error band becomes smaller with the evolution in µ but this reduction is still limited between the two presented scales. The distributions for other quarks exhibit similar uncertainties, on the other, hand the uncertainties for the gluon are substantially larger reaching up to 50% at small k t values (covering also the differences between proton and nuclear distributions).

Results
We present now our predictions for the inclusive Z boson production in pPb collisions at the LHC pPb → (Z/γ * ) → ¯ at √ s = 5.02 TeV and compare them with CMS data [20]. The intermediate vector boson is decaying into a pair of electrons or muons and these two channels are combined and we compare with this combined data. The measurement is done in the fiducial region defined by: p T > 20 GeV, |η lab | < 2.4 and 60 < m < 120 GeV. In our calculations we use leading order (LO) off-shell matrix elements as calculated by the KaTie Monte Carlo generator [7] and the TMDs (and PDFs) discussed in Sec. 2. The factorization and renormalization scales are set to be equal to the Z boson mass, µ = m Z . The considered data is restricted to the Z boson (center of mass) rapidity ranging between y * ∈ (−2.8, 2.0) which at LO corresponds to x values in lead x Pb ∼ (3 * 10 −3 , 0.3) (where negative rapidities correspond to large x Pb and positive y * correspond to small x Pb ). This shows that the probed x values are not very small and we do not need to worry about non-linear effects. In Fig. 6 we first present a comparison of predictions obtained using different PB TMDs for lead and proton [8]. In all cases we use Set2 distributions. We can see a very good description of the data provided by the k T -factorization framework. This is true for all the nuclear TMDs (PB-nCTEQ15FullNuc 208 82, PB-EPPS16nlo CT14nl Pb208, PB-gluon D c ncteq1568CL Pb). We observe only minor 10 100 CMS pPb √ s = 5.02 TeV 60 < m l + l − < 120 GeV PB-EPPS16nlo CT14nlo Pb208 set2 PB-nCTEQ15FullNuc 208 82 set2 PB-gluon D c ncteq1568CL Pb set2 PB-NLO ptoPb208 set2 CMS Data (b) Figure 6: Comparison of predictions for: (a) Z boson (center of mass) rapidity y * , and (b) Z boson transverse momentum, p T , distributions obtained within k T -factorization using different PB nTMDs/TMDs with the CMS data [20].
differences between predictions obtained with different nTMDs which suggest that the details of the nuclear corrections are here of secondary importance compared to the framework itself. On the other hand, we can see that neglecting the nuclear correction entirely (using lead composed out of free-proton TMDs -PB-NLO ptoPb208) undershoots the measured cross sections as a function of the rapidity and the transverse momentum. One should highlight that these are absolute distributions meaning that the LO framework we are using predicts not only the shape of the distributions but also their normalization (the uncertainty of the prediction is discussed at the end of this section). In order to understand the obtained results in more detail we perform an additional comparison of our predictions with LO collinear results as well as with the hybrid approach [21][22][23]. 4 The results are very similar for all the nuclear TMDs and therefore we concentrate only on the ones obtained with the PB-nCTEQ15FullNuc 208 82 distributions. In Fig. 7 we first compare the k T -factorization results obtained with the PB-nCTEQ15FullNuc 208 82 TMDs of Set1 and Set2 with the corresponding collinear predictions (of course, in the collinear case at LO the p T of the Z boson vanishes). We clearly see that the rapidity distribution obtained in the collinear -and k T -factorization with Set1 TMDs are very similar but give much worse description of the measurement compared to the calculation performed in k T -factorization calculation with Set2. The similarity of the collinear and TMD Set1 results is not surprising as the Set1 TMDs 4 Strictly speaking the hybrid approach is valid in the situation of more exclusive observables where one of the partons involved in the scattering has low-x and the other not, e.g. for jet production in the forward region. Here it will serve us rather as a tool for better understanding the k T -factorization results. reduce to the collinear distributions after the integration over the transverse momentum, and for the rapidity distributions such an integration is effectively carried out (with additional corrections coming from the off-shell matrix elements). Interestingly, the results obtained with Set2 are very different and describe the data much better than the Set1. The same holds for the p T distribution where results obtained with Set2 is spot on the data whereas the Set1 results are lower.
To further analyze the situation in Fig. 8 we do an additional comparison using only Set2 nTMDs but in addition to the k T -factorization and collinear results we include also results from the hybrid approach. In the hybrid approach the parton distributions of one beam are given by regular collinear PDFs whereas the second beam is described by the TMDs. Since we are in the central region (|η lab | < 2.4), a priori there is no clear choice on which beam should be described by the TMD and which by the collinear PDF. We present computations for both choices. Rather surprisingly we can see from the rapidity distribution in Fig. 8a that the hybrid calculation with an off-shell parton from the lead beam using the nTMDs gives results very close to the calculation with both partons being off shell (as in the k T -factorization calculation). On the other hand the hybrid calculation with an off-shell parton from the proton beam gives results very close to the collinear calculation. In the p T distribution, the hybrid result with an off-shell parton from the lead beam is also closer to the data, however, at low p T both hybrid results fail to describe the Sudakov suppression and fail to describe the data (which are very well described when both partons are treated off-shell in the calculation in k T -factorization).
In order to estimate missing higher orders it is conventional to perform a scale variation by a factor two up and down. Higher order contributions come from the radiation of additional partons, and therefore will depend on the transverse momentum of the Z boson: at p T = 0 no real parton can be emitted and the scale variation should lead to minimal changes. However, at larger p T the scale variation should give significant effects. We therefore introduce a scale µ 2 = m 2 Z + p 2 T , where we vary only the transverse momentum in order to estimate the effects from missing higher orders. A similar choice of the scale and scale variation was applied in Ref. [24]. In Fig. 9 we show a comparison of calculations using two different choices for the factorization scales µ 2 = m 2 Z (that we used in earlier plots) and µ 2 = m 2 Z + p 2 T ; the predictions are very similar. In Fig. 10 we show the result of the scale variation with µ 2 = m 2 Z + p 2 T , where the dynamical part (p T ) is varied by a factor of two up and down. The scale dependence is visible in the rapidity distribution and at larger transverse momenta (as expected). The size of the variation is similar to what is expected in NLO calculations.
Finally, in Fig 11 we show the effect of the PB-nCTEQ15FullNuc 208 82 Set2 TMD uncertainties on the rapidity and transverse momentum distributions of the Z boson. One can see that the TMD uncertainties are larger than the scale variation. They are also larger than the spread of predictions using different nTMDs (shown in Fig. 6) but not wide enough to cover also the predictions with proton TMDs (also in Fig. 6). This again highlights that the use of nuclear corrections is important but with the current experimental errors their details are secondary. dσ dp

Predictions for 8 TeV data
The data for Z boson production in proton-lead collisions at 8.16 TeV are currently being analyzed by LHC collaborations. As an addition to the results presented above we also provide prediction for such measurements. In Fig. 12 we show both rapidity and transverse momentum distributions of the produced Z boson (lepton pair). We present results for two ranges of the invariant mass of the produced lepton pair: around the Z boson mass peak 60 < m ll < 120 GeV, and for low mass DY pair 15 < m ll < 60 GeV.

Conclusion
We have constructed the first sets of nuclear TMDs for gluons as well as all the quarks. It was done using the parton branching method employing different choices of α S arguments. The obtained nTMDs are provided to all interested colleagues on the TMDlib [25] website https://tmdlib.hepforge.org. We have used the constructed nTMDs in the k T -factorization framework to calculate predictions for Z boson production at pPb collisions at the LHC. The obtained results show very good description of the CMS data [20] for both rapidity and p T distributions. One should highlight here that predictions are not only for the shape of the distributions but also for their normalization.