Collinear and TMD parton densities from fits to precision DIS measurements in the parton branching method

Collinear and transverse momentum dependent (TMD) parton densities are obtained from fits to precision measurements of deep inelastic scattering (DIS) cross sections at HERA. The parton densities are evolved by DGLAP evolution with next-to-leading-order (NLO) splitting functions using the parton branching method, allowing one to determine simultaneously collinear and TMD densities for all flavors over a wide range in $x$, $\mu^2$ and $k_t$, relevant for predictions at the LHC. The DIS cross section is computed from the parton densities using perturbative NLO coefficient functions. Parton densities satisfying angular ordering conditions are presented. Two sets of parton densities are obtained, differing in the renormalization scale choice for the argument in the strong coupling alpha_s. This is taken to be either the evolution scale $\mu$ or the transverse momentum $q_t$. While both choices yield similarly good $\chi^2$ values for the fit to DIS measurements, especially the gluon density turns out to differ between the two sets. The TMD densities are used to predict the transverse momentum spectrum of Z-bosons at the LHC.


Introduction
Parton density functions (PDFs) play an essential role for precise predictions of production processes in hadronic collisions, obtained from the factorization of the cross sections in hardscattering process and PDFs, containing a non-perturbative input with perturbatively calculable evolution.The most advanced determination of parton densities come from the application of DGLAP [1][2][3][4] evolution with next-to-leading order (NLO) [5,6] and next-to-next-toleading order (NNLO) [7,8] splitting functions.The collinear parton densities as a function of the longitudinal momentum fraction x and the evolution scale µ 2 are obtained by several groups, for example ABM [9], CTEQ [10], HERAPDF [11], NNPDF [12] and MSTW [13,14].The different groups use the same DGLAP evolution, with ordering in virtuality and the same choice of the renormalization scale, but they differ in, for example, the treatment of heavy flavors, and the experimental data sets which are used for the determination of the starting distributions.
In Refs.[15,16] a new method, the Parton Branching method (PB), was introduced to treat DGLAP evolution.The method applies at exclusive level, and provides an iterative solution of the evolution equations.It agrees with the usual methods to solve the DGLAP equations for inclusive distributions, but it provides also additional features: in addition to the standard ordering in virtuality, angular ordering can be applied with the necessary change in the argument of α s [17,18].The transverse momentum at every branching vertex can be calculated, leading to a natural determination of the transverse momentum dependent (TMD) parton densities.The PB method uses the unitarity formulation of QCD evolution equations [19] and is close in spirit to the works in [20][21][22][23][24][25].As shown in Refs.[16,26], it can be applied to NLO and NNLO splitting functions.
In this article we present a determination of collinear and TMD parton densities at NLO applying the PB method for the parton evolution.The initial parton distributions are determined from a fit to HERA I+II inclusive DIS cross section measurements [11].An early fit was presented in Ref. [16].Here, we present results obtained with angular ordering, both for collinear (integrated, iTMD) and TMD parton densities, and for different choices of the renormalization scale in α s including a full treatment of experimental and model dependent uncertainties.We show an application of these TMDs to the calculation of the transverse momentum of the Z-boson in Drell-Yan (DY) production at the Large Hadron Collider (LHC).

Parton Branching method and evolution equation
The PB method has been described in detail in Refs.[15,16].Here we limit ourselves to recalling its main elements.

General features
The method is based on introducing a soft-gluon resolution scale z M into the QCD evolution equations to separate resolvable and non-resolvable emissions, and treating these via, respectively, the resolvable splitting probabilities P (R) ba (α s , z) and the Sudakov form factors Here a, b are flavor indices, α s is the strong coupling at a scale being a function of µ 2 to be specified in Section 3, z is the longitudinal momentum splitting variable, and z M < 1 is the soft-gluon resolution parameter.For easier reading we use the notation ∆ a (µ 2 ) = ∆ a (z M , µ 2 , µ 2 0 ).The form factors eq. ( 1) have the interpretation of probabilities for nonresolvable branchings between the evolution scales µ 0 and µ.The functions P (R) ba (α s , z) have the structure where the first term on the right hand side contains the pole singularity in the soft-gluon radiation region z → 1 and the second term contains logarithmic terms and analytic terms for z → 1.The coefficients k b and R ba in eq. ( 2) have the perturbation series expansions The explicit expressions for the n = 1 (LO) and n = 2 (NLO) contributions in the expansions in eq. ( 3) are given in [16].The n = 3 (NNLO) contributions can be read from [7,8] and are used for NNLO calculations in the PB method in [26].The integrals appearing in the Sudakov form factors eq. ( 1) are positive at LO, NLO and NNLO, while the functions eq. ( 2) can be negative at NLO and NNLO.The positivity of the integrals in eq. ( 1) is essential for the application of the PB method.
The PB method allows one to take into account simultaneously soft-gluon emission in the region z → 1 and transverse momentum q ⊥ recoils in the parton branchings along the QCD cascade.Its advantage is twofold: on one hand, in collinear distributions additional QCD features can be studied such as the color radiation's angular ordering, determined by soft-gluon interferences, and its effects on factorization and renormalization scales; on the other hand, the method can be applied to obtain transverse momentum dependent (TMD) distributions.
The PB evolution equations for TMD parton densities A a (x, k, µ 2 ) are given by [16] A in terms of the ∆ a form factors, eq. ( 1), and ba functions, eq. ( 2).The scale in α s is a function of q 2 , as discussed in Section 3.These equations can be solved by an iterative Monte Carlo method.In this method every resolvable branching is reconstructed explicitly and the full kinematics at each branching is taken into account.The PB method allows to solve eq. ( 4) in an easy and direct way, with the possibility to include, for example, also heavy quark masses and soft-gluon coherence conditions.
The collinear parton densities f a (x, µ 2 ) are related to the TMD densities by and are described as integrated TMD (iTMD).The evolution equations for iTMD densities analogous to eq. ( 4) can be written as These equations have been shown to be equivalent to DGLAP evolution equations at NLO [15,16,20,21] and NNLO [26] for α s = α s (µ 2 ) and z M → 1.

PB method and determination of initial distribution
The PB method has been implemented in the xFitter package [27] to allow fits to be made to cross section measurements.A full Monte Carlo solution of the evolution equation for every new set of initial parameters would be too time consuming to be efficient.Instead, a method developed already in [28][29][30] is applied: first, a kernel K int ba x , µ 2 0 , µ 2 is determined from the Monte Carlo solution of the evolution equation for any initial parton * of flavor b evolving to a final parton of flavor a; then this kernel is folded with the non-perturbative starting distribution f 0,b (x, µ 2 0 ), The kernel K int ba includes the full parton evolution from µ 2 0 to µ 2 , as in eq. ( 6), with Sudakov form factors and splitting probabilities, and is determined with the PB method.In eq. ( 7) the kernel K int ba depends on x, µ 2 0 and µ 2 for the k t -integrated (iTMD) distributions.To include also the transverse momentum k t , we define a new kernel K ba x , k 2 t,0 , k 2 t , µ 2 0 , µ 2 for the TMD distributions, with 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.
Technically, the results of the kernel evolution are stored in a grid of size 50×50(×50) (for the TMD densities).The grid spacing is logarithmic (µ 0 < µ < 14000 GeV and 0.01 < k t < 14000 GeV), the x range is divided into 5 subregions with logarithmic spacing: subregions of 10 bins are defined with the boundaries 10 −6 , 0.01, 0.1, 0.4, 0.9, 1 which is optimized to ensure appropriate behavior for large x, where the parton densities (and the kernel) are varying rapidly.
In Fig. 1 we show the result of convoluting the starting distribution (here taken to be the benchmark parameterization of Ref. [31]) with the kernel as given in eq. ( 7) for the integrated distribution, and compare this with the prediction from a standard evolution program (QCDNUM) for different values of the evolution scale µ 2 .The kernel is evolved using NLO splitting functions with resolution scale parameter z M , separating resolvable from nonresolvable branchings, set to the value z M = 0.99999.Very good agreement is observed over the whole range.Only the quark distribution shows differences at very large x of the order of a few percent, which come from the finite grid spacing in x when storing the kernel (changing to a uniform logarithmic grid spacing in x leads to significantly larger deviations at large x).In most of the phase space region relevant for high precision physics at HERA and the LHC the differences are at the per mille level.Figure 1: Comparison of the results from the convolution in eq.( 7) with the prediction from QCD-NUM [32] using the same input distributions, for d-quarks (left) and gluons (right) at different values of the evolution scale µ 2 starting from µ 2 0 = 1.9 GeV 2 with α s (µ 2 ).The lower panels show the ratio of the parton density with the one predicted by QCDNUM.The evolution is performed with NLO DGLAP splitting functions and using z M = 0.99999.

Parton densities obtained from fits to inclusive HERA DIS measurements
The most recent and most precise measurements of the lepton-proton DIS cross section over a wide range in x and Q 2 were performed at HERA with a combination of the measurements from the H1 and ZEUS collaborations [11].These measurements are the basis for any determination of parton densities.In Ref. [11] a fit to the inclusive DIS measurements was performed using DGLAP at LO, NLO and NNLO, resulting in the HERAPDF2.0parton distributions.These fits were performed with QCDNUM [32] within the xFitter framework [27] using a starting scale µ 0 = 1.9 GeV 2 and the renormalization and factorization scales set to The light quark matrix elements were taken from QCDNUM, the heavy-quark contributions were obtained within the general-mass variable-flavor scheme RTOPT [33][34][35] for neutral current, while for charged current interactions the zero-mass approximation from QCDNUM was used.The following parameterizations are used for the different parton flavors: The quark-number sum rules and the momentum sum rule can be used to constrain the normalization parameters, A u v , A d v , A g , A g .The B parameters are set B Ū = B D for the sea distributions.The strange-quark distribution is parameterized as a d-type sea with an xindependent fraction, f s , xs = f s x D at µ 2 0 with f s = 0.4.A further constraint was applied by setting A total of 1145 data points of neutral-current and charged-current deep-inelastic cross section measurements were used in the range of 3.5 < Q 2 < 50000 GeV 2 and 4 • 10 −5 < x < 0.65.
The same data sets, kinematic ranges and hard-scattering coefficient functions, including the heavy-quark treatment, are used for the fits described here.We use NLO DGLAP splitting functions [5,6] as well as NLO coefficient functions [36] for light quarks.For heavy quarks we apply the general-mass variable-flavor scheme RTOPT [33][34][35] for neutral current, while for charged current interactions the zero-mass approximation is used.
In the next section we determine the free parameters of the initial distributions given by eq. ( 10) via fits to the HERA DIS data in the range of Q 2 > 3.5 GeV 2 using NLO DGLAP splitting functions within the PB method using z M = 0.99999.The PB method allows the explicit calculation of the kinematics at every branching vertex (see Fig. 2 left).Once the physical meaning of the evolution scale is specified in terms of kinematic variables, the transverse momenta of the propagating and emitted partons can be calculated.In Ref. [15] it was pointed out that angular ordering gives transverse momentum distributions which are stable with respect to variations of the resolution parameter z M .In angular ordering, the angles of the emitted partons increase from the hadron side towards the hard scattering, as shown in Fig. 2 right.The transverse momentum q t i can be calculated in terms of the angle Θ i of the emitted parton with respect to the beam directions from q t,i = (1 − z i )E i sin Θ i .Associating the "angle" E i sin Θ i with µ i gives In the following, we use the PB method to determine collinear (iTMD) and transverse momentum dependent (TMD) parton densities using NLO DGLAP splitting functions for two different scenarios: first we only apply the angular ordering condition for the calculation of the transverse momentum and keep the evolution scale µ 2 i as the argument in α s (Set 1); in a second scenario (Set 2), we use (in eqs.(1,4,6)) the transverse momentum |q 2 t,i | as the argument in α s , as suggested in Ref. [17,18].An additional parameter q cut needs to be introduced in α s (max(q 2 cut , |q 2 t,i |)) to avoid the non-perturbative region, since with large z the scale |q 2 t,i | = (1 − z i ) 2 µ 2 i can become very small.We take the default choice for this parameter to be q cut = 1 GeV, and we estimate the model dependence with a variation around the default choice.
In the first case, the integrated parton density, and the initial parameters, will be the same (up to numerical precision) as the ones obtained by HERAPDF2.0,and we use this as a benchmark for the whole method.In the second case, even the integrated parton distributions differ, because of the different scale in α s .In both cases a reasonably good fit is obtained with χ 2 /ndf ∼ 1.2, as for HERAPDF2.0.In Tab. 1 results of the fits are given.The starting scale µ 2 0 is chosen differently for the 2 scenarios: for Set 1 we chose (as in HERAPDF) µ 2 0 = 1.9 GeV 2 while for Set 2 we chose µ 2 0 = 1.4 GeV 2 , which gave the best χ 2 /ndf .In the appendix we show results obtained from a fit when µ 2 0 = 1.9 GeV 2 is chosen instead of µ 2 0 = 1.4 GeV 2 .The distributions agree within their uncertainties.The values of the parameters at the starting scale µ 2 0 are given in Tab.

Collinear Parton Densities (iTMD)
The fits to HERA measurements are performed using χ 2 minimization, as in the case of the HERAPDF fits, implemented in xFitter [27].The definition of χ 2 includes systematic shifts, a treatment of correlated and uncorrelated systematic uncertainties.In total 162 systematic uncertainties plus procedural uncertainties from the combination of H1 and ZEUS are treated as correlated uncertainties.
Central The experimental uncertainties of the resulting parton densities are determined with the Hessian method [37] (as implemented in xFitter ) with ∆χ 2 = 1.The model dependence of the PDF fits is obtained by varying charm and bottom masses and the starting scale of the evolution µ 2 0 .For Set 2 also the parameter q cut is varied.The central values and the range of variation is given in Tab. 3.
In Fig. 3 the Ū -type quark and gluon densities are shown as functions of x for different values of the evolution scale µ 2 = Q 2 including the experimental uncertainties (red band) and the uncertainties coming from the model dependence (yellow band).For Set 2 the uncertainty of the parameter q cut is shown as the green band.The results of Set 1 are identical to the ones obtained in HERAPDF 2.0.Although the fits (Set 1 and Set 2) to HERA I+II data are of similar quality, the resulting parton distributions, especially for the gluon, are significantly different.With increasing evolution scale, however, they become more and more similar.
In Fig. 4 the total uncertainties (experimental and model) of the parton densities are shown.The uncertainties of Set 2 for the gluon distribution at large x become large.We have investigated a possible bias coming from the chosen form of the parameterization by including additional terms for the gluon density: The obtained χ 2 of the fit does change by at most 1 unit, the resulting gluon distribution does not change visibly.Details of the bias study are given in the appendix.
In Fig. 5 we show predictions for the inclusive DIS cross section and the inclusive charm cross section obtained from the two different parton distributions, and compare them with the measurements from HERA [11,38].While the inclusive DIS cross section is well described, the prediction using Set 2 differs from inclusive charm measurement at low Q  [11], lower row: inclusive charm production [38].The dashed lines include the systematic shifts in the theory prediction.small x.For values x > 0.001 all predictions agree reasonably well with the data.It has been checked explicitly that including the charm measurements in the fits does not significantly change the fit result (the charm data have too large an uncertainty compared to the precise inclusive measurements).In Fig. 5 the predictions including the systematic shifts are also shown, visually showing that the quality of the two different fits is similar.

Transverse Momentum Dependent Parton Densities (TMD)
Within the PB method both collinear and TMD densities can be determined, as the transverse momentum is calculated at every step of the branching process.TMD parton densities can be obtained via the PB method once the relationship between kinematical variables and evolution scale µ is specified, and the transverse momentum at each individual branching is calculated with eq.(11).The parameters for the starting distributions are obtained for the collinear parton densities by a fit to inclusive DIS cross section measurements, as described previously.The TMD parton densities are then obtained from a convolution of the TMD kernel with the starting distribution as given in eq. ( 8).The starting distribution is taken from the collinear iTMD described in Sec.3.1.In Fig. 6 we show the TMD parton densities for ū-quarks and gluons as a function of the transverse momentum k t = k 2 for different values of the evolution scale µ = 10, 100, 1000 GeV and different values of x for Set 1 and Set 2. One can clearly see that both sets give identical results for larger k t , while they are different for small k t , a consequence of the different scale choices for the argument of α s .
In Fig. 7 the parton densities for all flavors are shown as a function of k t at x = 0.01 and for different values of the evolution scale µ = 10, 100, 1000 GeV.The large scales are relevant for phenomenology at the LHC, and it is interesting to observe that the transverse momenta extend to very large values, up to the values of the factorization scales (for µ = 1 TeV the transverse momenta extend to k t ∼ 1 TeV).However, the large k t values are suppressed compared to smaller ones.The different quark flavors show a different behavior at small k t , coming essentially from the no-branching probability times the starting distribution (first term in eq. ( 4)), while they are very similar at larger k t , a result of perturbative splittings (second term in eq. ( 4)).In Fig. 8 the gluon and ū densities as a function of the transverse momentum are shown for µ = 100 GeV and x = 0.01 together with the uncertainty bands obtained from the fits.The panels show the uncertainties coming from the experimental sources as well as the total uncertainty coming from experimental and model sources separately.Although only collinear splitting functions are used, and the fit was obtained with collinear parton densities, a k t dependence of the uncertainties is obtained.At small k t essentially the first term in eq. ( 4) contributes without any resolvable branching and the uncertainty comes from the starting distribution at x, while at large k t several branching may have occurred and therefore the uncertainty comes from the starting distribution at x/z x.The experimental uncertainties are small over the whole range, while the model dependent uncertainties dominate.
The parametrization of the intrinsic transverse momentum distribution is another uncertainty.With the fit to inclusive DIS data, this distribution cannot be further constrained.In Fig. 9 we show the TMD distribution for gluon and ū for Set 1 and Set 2 at µ = 10(100) GeV and x = 0.01 when q 0 in B(k 2 t,0 , µ 2 0 ) is varied from q 0 = 0.25 GeV to q 0 = 1 GeV.We do not include the variation of q 0 as a systematic uncertainty, since it is not constrained by the fit (in future we plan to use also Z-boson transverse momentum spectra, which would constrain q 0 ).
The resulting TMD parton densities, PB-NLO-2018-Set1 and PB-NLO-2018-Set2, including uncertainties (as well as with variation of q 0 ) are available in TMDLIB [39].The TMD- PLOTTER [40,41] interface allows easy and fast comparison to other TMDs, once they are made publicly accessible and available in TMDlib.

Application to Z-boson production at the LHC
The transverse momentum spectrum of Z bosons in Drell-Yan (DY) production at small values of transverse momentum q T cannot be described by fixed-order perturbative calculations, and resummation of soft gluon emissions to all orders in α s is needed.See e.g.[42] for a recent discussion.The DY q T spectrum can be described by the CSS method [43][44][45][46] using TMD factorization at small q T [47,48], or by parton showers within Monte Carlo event generators [25].The ATLAS and CMS experiments at the LHC have measured the q T spectrum of the Z-boson [49][50][51].
The TMD distributions obtained from HERA DIS measurements can be used to predict the DY q T spectrum of the Z-boson at LHC energies.Since we are interested in the low-q T region, we use the LO expression for Z production matrix elements.† The transverse momentum of the initial state partons is calculated according to the TMDs and added to the event record in such a way that the mass of the produced DY pair is conserved, while the longitudinal momenta are changed accordingly.This procedure is common in standard parton shower approaches [53,54] and is implemented in the CASCADE package [55,56] (version newer than 2.4.X) where events in HEPMC [57] format are produced, for further processing with Rivet [58].The importance of the proper inclusion of transverse momentum effects from parton showers has been pointed out in Ref. [59,60].With the TMD distributions described here, these effects can be included already at the level of the cross section calculation.
In Fig. 10 (left) we show the predictions for the transverse momentum spectrum of the Z-boson obtained with the two TMD distributions, compared with the measurements of AT-LAS [51].The uncertainties coming from experimental and model sources are shown for both Set 1 and Set 2 with the colored bands (Fig. 10 left); the experimental and full uncertainties are shown for Set 2 in Fig. 10 (right).The difference between the full and experimental uncertainties from the fit is very small.
In general the shape of the spectrum is described by both TMD fits.The TMD Set 2, applying the transverse momentum as the renormalisation scale (instead of the evolution scale µ), provides a significantly better description of the transverse momentum spectrum of the Z-boson, coming from the different k t spectrum of the TMD already visible in Fig. 6.One should note that no adjustment of any parameter is made, and that the TMDs are entirely constrained by the fits to inclusive DIS data.The description of the transverse mo-  mentum spectrum of the Z-boson obtained with the PB-TMD Set 2 is of similar quality as the NLO+NNLL prediction of Ref. [61], however, one should note, that the approach of PB-TMDs is more general and can be applied directly to other processes as well without further modification.

Conclusion
The parton branching method has been used to determine a first complete set of collinear and TMD parton densities from fits to precision DIS data over a large range in x and Q 2 as measured at HERA.The parton densities are obtained with NLO DGLAP splitting functions and 2-loop α s with α s (M Z ) = 0.118.The renormalisation scale in the evolution has been chosen to be the evolution scale µ i (Set 1) or the transverse momentum q t i (Set 2).Two different collinear and TMD sets are obtained for these different choices, both giving a similar χ 2 /ndf = 1.2.The obtained parton densities are valid over a wide range in x and scale µ, up to the multi-TeV scale, relevant for LHC physics.Experimental uncertainties of the fit are obtained using the Hessian method with ∆χ 2 = 1 and model dependent uncertainties are determined.
The obtained TMDs are applied to calculate the transverse momentum spectrum of the Z-boson in DY production at LHC energies.Good agreement with the measurement is observed if angular ordering is applied.The uncertainties of the prediction come only from the TMD uncertainties determined in the fit to HERA measurements.
For the first time, precision DIS measurements have been used to obtain both collinear and TMD parton densities, including uncertainties, over a wide range in x and µ values, which are relevant for LHC and future collider phenomenology as well as for low-energy and small-k t physics.A potential bias of the form of the parameterization was checked by extending the original parameterization xg(x) = A g x B g (1 − x) C g − A g x B g (1 − x) C g with additional parameters: In Fig. 12 we show the gluon distribution after fitting C g and including the additional factors D g and E g one after the other.The starting scale is µ 2 0 = 1.4 GeV 2 (as for the original fit Set 2).The obtained χ 2 is larger by 1 unit after including additional terms, the shape of the distribution does not change significantly.The uncertainty band of Set 2 corresponds to the uncertainties coming form the experimental sources, no model or parameterization uncertainty is included.The parton distributions agree within the uncertainties shown, excluding a significant bias from the chosen form of the parametrization.
The mass of the charm quark is set m c = 1.47 GeV, and m b = 4.5 GeV is used for the bottom quark mass.The strong coupling is set to α s (M 2 z ) = 0.118.The parameterized PDFs are the gluon distribution, xg, the valence-quark distributions, xu v , xd v , and the u-type and d-type anti-quark distributions, x Ū , x D. The relations x Ū = xū and x D = x d + xs are assumed at the starting scale µ 0 .

Figure 2 :
Figure 2: Left: Branching process b → a + c.Right: Schematic view of a parton branching process.

Figure 3 :
Figure 3: Parton densities for different values of the scale µ 2 = Q 2 .The different choices for the renormalization scale in α s are shown.The red band shows the experimental uncertainty, the yellow band the model dependence.The green band shows the uncertainty coming from the variation of the parameter q cut in Set 2.

Figure 4 :
Figure 4: Total uncertainties (experimental and model uncertainties) for the two different sets at different values of the evolution scale µ 2 .

Figure 5 :
Figure5: Measurement of the reduced cross section obtained at HERA compared to predictions using Set 1 and Set 2. Upper row: inclusive DIS cross section[11], lower row: inclusive charm production[38].The dashed lines include the systematic shifts in the theory prediction.

Figure 6 :
Figure 6: Transverse Momentum Dependent parton densities (PB-NLO-2018-Set1 and PB-NLO-2018-Set 2) as a function of k t for different scales µ.Upper row shows the densities for ū, lower row the densities for gluons for two different values of x.

Figure 7 :
Figure 7: Transverse Momentum Dependent parton densities (PB-NLO-2018-Set1 upper row and PB-NLO-2018-Set2 lower row) as a function of k t for different scales µ at x = 0.01 for all flavors.

Figure 8 :
Figure 8: Transverse Momentum Dependent parton densities for ū and gluon from Set 1 and Set 2 as a function of k t for µ = 100 GeV at x = 0.01.In the lower panels we show the relative uncertainties coming from experimental uncertainties as well as the total of experimental and model uncertainties.

Figure 9 :
Figure 9: Transverse Momentum Dependent parton densities (ū and gluon) from Set 1 and Set 2 as a function of k t for µ = 10(100) GeV at x = 0.01, when the width of the intrinsic transverse momentum distribution is varied by a factor of two.

Figure 10 :
Figure10: Transverse momentum q T spectrum of Z-bosons obtained from the two TMDs, compared with measurements from[51].Left: comparison of predictions using Set 1 and Set 2 including the full (experimental and model) uncertainties.Right: prediction using Set 2, with experimental and full uncertainties separated (the difference is very small).

Figure 11 : 2 0 = 1 . 4
Figure 11: Comparison of gluon density of Set 2 type obtained at µ 2 0 = 1.4 GeV 2 and µ 2 0 = 1.9 GeV 2 at a scale of Q 2 = 3 GeV 2 .The ratio of the gluon densities is shown with respect to the default Set 2. The uncertainties for the new fit include only those from experimental sources, the uncertainties for Set 2 are the same as in Fig 4.

Figure 12 :
Figure 12: Comparison of gluon densities after fit when additional terms in the gluon parametrization are included.The uncertainty band of Set 2 corresponds to the uncertainties coming form the experimental sources.

Table 1 :
Values of χ 2 for the different fits at NLO.

Table 2 :
Parameter values of the initial distributions at NLO.The parameter C = 25 was fixed, as in HERAPDF2.0.The parameters correspond to a starting scale µ 2 0 = 1.9(1.4)GeV 2 for Set 1 (Set 2).