Impact of the PDFs on the Z and W lineshapes at LHC

The parton distribution functions (PDFs) of the proton play a role in determining the lineshape of Z and W bosons produced at the LHC. In particular, the mode of the gauge boson virtuality is shifted with respect to the pole due to the dependence of the partonic luminosity on the boson virtuality. The knowledge of this shift contributes to the systematic uncertainty for a direct measurement of the boson mass. A detailed study of the shift and of its systematic uncertainty due to the limited knowledge of the PDFs is obtained using a tree-level model of Z and W boson production in proton-proton collisions at s=13TeV\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\sqrt{s}=13~\hbox {TeV}$$\end{document}. A Monte Carlo simulation is further used to validate the tree-level model and study the dependence of the shift on the transverse momentum of the gauge bosons. The tree-level calculation is found to provide a good description of the shift. The systematic uncertainty on the lineshape due to the PDFs is estimated to be below one MeV in the phase-space relevant for a future high-precision mass measurement of the gauge boson masses at the LHC.


Introduction
The masses of gauge bosons are some of the most relevant observables to test the electroweak theory and they are measured with high precision. The Z mass was extracted at LEP from a fit to the lineshape, i.e. from the cross section of the process e + e − → γ, Z → f f measured with different beam energies spanning center-of-mass energies near √ s ∼ M Z . At LEP, the lineshape was distorted by the initial state radiation of the colliding electrons, which was theoretically very well understood, resulting in a systematic uncertainty of less than 0.1 MeV on the Z mass. The Z mass was ultimately measured with an uncertainty of 2.1 MeV [1] where the largest contribution was the precision of the beam energy calibration. a e-mail: valerio.bertacchi@pi.infn.it (corresponding author) The experimental situation at CERN's Large Hadron Collider (LHC) differs by at least two aspects. First, protonproton collisions at a center-of-mass-energy √ s result in partonic f f collisions with a broad distribution of partonic center-of-mass-energies √ŝ = √ x 1 x 2 s, where x 1 and x 2 are the fractions of the proton momentum carried by the interacting partons, as described in the empirical parton density functions (PDFs). The center of mass of the initial state is unknown and changes on an event by event basis covering all relevant energies. The lineshape of the boson can be measured using the decay product kinematics, like the invariant mass of the μ + μ − pairs in Z decays.
Secondly, at the LHC the gauge bosons are produced at a rate that is several orders of magnitude larger than at LEP. The ATLAS and CMS experiments already collected about of 400 million W and 40 million Z bosons each during the Run 2 of the LHC. By the end of Run 3, 10 8 Z leptonic decays will be available and a factor of 10 more at the end of the High Luminosity LHC program. These numbers must be compared to the 10 million Z bosons produced at LEP in all decay channels. This unprecedented number of W and Z bosons produced at the LHC offers new opportunities for precise measurements, but it also forces consideration of sources of systematic uncertainty which may have been neglected so far, for example those related to the modeling of the virtuality of the gauge bosons [2]. For instance, by analyzing the dileptonic Z decays of Run 3 of the LHC a statistical-only precision of about 10 −5 GeV on M Z might be achievable. This level of precision would demand a control of the dilepton mass lineshape at the sub-MeV level.
A thorough assessment of the systematic uncertainty on the lineshape is important because the Z mass is used to calibrate the muon momentum scale of the detectors [3,4]. In addition this information may be useful for a new, more precise, measurement of the Z mass in case the momentum scale of the detector could be calibrated to a relative uncertainty of 10 −5 using independent experimental information, like the J/ψ mass which is known with a relative uncertainty of 10 −6 [5], and assuming that final state radiation effects can be understood to this level of precision.
At first order, the distribution of the virtuality Q of a gauge boson (V ) originates from the convolution of a relativistic Breit-Wigner with the partonic luminosity function, see e.g. Ref. [6]. The latter is a function of the dimensionless parameter τ = Q 2 /s. The non-trivial dependence of the partonic luminosity on τ implies a distortion of the lineshape compared to a pure Breit-Wigner. Given the narrowness of the electroweak gauge boson width Γ V , this effect can be treated, in first approximation, as a shift Δ V of the mode of the distribution compared to M V . The limited knowledge of the PDFs introduces an uncertainty on Δ V , which contributes directly to the model uncertainty in the extraction of M V from the dilepton mass distribution.
The goal of this paper is to assess the size of this shift and of its uncertainty due to the limited knowledge of the PDFs. This shift can be regarded as a proxy of the systematic uncertainty on the Z mass extracted from the fit to the dilepton mass distribution.
The W boson lineshape is also distorted by the same effect. However, this is mostly of academic interest since the invariant mass of the leptonic final state cannot be measured in W decays due to the presence of the neutrino. The traditional measurement of the W mass at hadron colliders uses non-Lorentz invariant quantities (e.g. the transverse mass or the lepton transverse momentum) whose distributions have a dependence on the PDFs which induces a systematic uncertainty much larger than the effect discussed in this paper [7]. However, the larger than ever number of W decays collected by the LHC opens up possibilities for novel measurements that might be less sensitive to the PDFs [8,9]. Likewise, the increased coverage in rapidity planned by the upgrades of the LHC experiments may have a strong impact in reducing the aforementioned PDF uncertainty [10,11]. For the above reasons, the results of this paper may become relevant also for a future W mass measurement at LHC. This paper is organized as follows. In Sect. 2 a tree-level calculation of the shift Δ V is presented. In Sect. 3, a Monte Carlo (MC) simulation of Drell-Yan production is used to validate the tree-level model and extend the study to the full phase-space. The results are summarized in Sect. 4.

Tree level study
A simplified model of Drell-Yan production is first considered based on a minimal set of tree-level diagrams. This approximation amounts to considering just one Feynman diagram per quark-antiquark pair, as illustrated in Fig. 1a. Besides accounting already for the bulk of the total cross section (about 80% for a 20 GeV threshold on the transverse momentum of the extra parton at √ s = 13 TeV), these dia- grams are also expected to be the most sensitive to the PDFdependent shift under study. Indeed, they are the only 2 → 1 diagrams contributing to the amplitude, whereas higher-order diagrams are at least 2 → 2, see e.g. Fig. 1b. As such, they include additional invariants besides Q. The existence of these extra scales is expected to dilute the sensitivity of the lineshape to the details of the PDFs. This assumption will be validated by a MC analysis of pp → V + X production discussed in Sect. 3. In the following, the NNPDF3.0 [12] PDF set will be used to evaluate the PDFs relevant for W and Z production in proton-proton collisions at √ s = 13 TeV. Within the tree-level approximation, the doubledifferential cross section for pp → V (→ ) + X , as a function of the quark momentum fractions x 1,2 , is given by where N C is the number of QCD colours, M V and Γ V are the mass and width of the resonance, BR V →ab are the relevant branching fractions, and the sum at the right-hand side runs over the different combinations of quark flavours contributing to the process under study. For simplicity, the scaledependence is omitted from the quark PDF f i (x). In Eq. 1 a relativistic Breit-Wigner function with fixed width has been assumed, which is also the functional form used for the MC simulation discussed in Sect. 3. The opportunity of using a running width scheme [5] has been studied as well. However, since the interesting feature in Eq. 1 concerns the core of the distribution, where Q ∼ M V , the results would not change because the two schemes differs only when |Q − M V | Γ V . In order to express the double-differential distribution of Eq. 1 as a function of Q and y the canonical transformation is performed. By combining Eqs. 1 and 2, the singledifferential distribution dσ V /d Q, conditional on the rapidity y, is obtained: where x 1,2 = √ τ e ±y and the constants C V i j include terms that depend on y but not on Q. The fact that Γ V /M V 1 and that f i are smooth functions in the relevant range of Bjorken x values (10 −3 x 10 −1 ) can be exploited to perform a Taylor expansion of the right-hand side of Eq. 3 aroud Q = M V : where the flavour-dependent terms have been factored out of BR V →q i q j in the form of the square of the V matrix elements. The latter should be interpreted as the usual CKM matrix for the case of W production, and as T 3 i and Q i are the weak isospin and electric charge of quark i, respectively. In Eq. 4, the following auxiliary functions of y have been introduced: where f ( f ) are the first (second) order derivatives 1 of the PDFs with respect to x. The constants H V and K V defined in Eq. 4 represent the appropriate average of the auxiliary functions of Eq. 5 over the flavour space. The validity of the Taylor expansion of Eq. 4 has been assessed by comparing the lineshape from Eqs. 3 and 4 at different values of y: the relative difference between the two is found to be below 0.5% for Q in a range of ±2 GeV around M V .
In Eq. 4, the contribution of the PDFs to the lineshape is fully encoded in the constants H V and K V . The mode Q 0 of the lineshape can be readily calculated from Eq. 4: where the approximation Γ 2 represents the displacement of the mode Q 0 from M V . Part of it is simply due to the Jacobian factor from the transformation of Eq. 2, and does not depend on the PDFs. The right-hand side of Eq. 6 depends on K V only at higher order in Γ V /M V because it enters as the coefficient of a quadratic correction to the Breit-Wigner functions, which is symmetric around M V . Both H V and K V are functions of y and M V , albeit the dependence on the latter is negligible in the range of experimental uncertainty on M Z (∼ 2 MeV) and M W (∼ 12 MeV) compared to the PDF uncertainties on Δ V . The shift Δ V determined from Eq. 7 is plotted in Fig. 2 as a function of the boson rapidity y for Z and W ± production. The error bars correspond to the RMS of the distribution obtained by sampling 100 replicas in the NNPDF30_nlo_nf_5_pdfas set from the LHAPDF library [13]. A comparison with NNPDF31_nlo_pdfas, which includes Drell-Yan measurements from the 8 TeV run of the LHC [14], has been performed as well, showing consistent results both in the central value and in the uncertainty. Numerical values are reported in Table 1 for three representative values of y. A negative shift with typical size |Δ V | ∼ 13 MeV is observed in the central region |y| 3, steeply increasing at larger rapidity values. This behaviour can be understood qualitatively in terms of the valence quark density xu V and xd V , which are typical benchmarks in PDF fits [12]. Indeed, for any derivable and positive-definite func- 1 All the derivatives are evaluated as f = 1 , with a step of Δx = 10 −5 . This step has been varied in Δx ∈ [10 −6 , 10 −4 ] and the result is found to be stable in this range.
holds that The left-hand side of Eq. 8 is of the same form of the terms that appear in the definition of H V (see the second line of Eq. 5). The valence quark densities feature a local maximum at x ∼ 10 −1 , which corresponds to |y| ∼ 3 at Q ∼ 90 GeV. By identifying f in Eq. 8 with u V or d V one can easily see that the terms (x f ) / f vanish around |y| ∼ 3, thus giving the smallest shift, whereas they steeply decrease at higher rapidity values since f → 0 and (x f ) becomes negative. The relative PDF uncertainty on Δ V is found to be in the 5% ballpark, ranging from 0.3 MeV at |y| ∼ 0 to 1 MeV at |y| ∼ 3.5. A further inspection of Fig. 2 shows that in the case of the Z boson, the shift happens to lie between the shifts for the W boson of opposite charge, as an effect of the different partonic favours probed by the gauge bosons.
Finally, the approximation which leads to Eq. 7 has been numerically validated. The resulting K V values range between 5 and 15, depending on y, with a relative uncertainty below 3%.

Monte Carlo simulation study
The tree-level calculation of Sect. 2 has been validated by using a MC simulation of pp → V + X production. Besides corroborating the tree level model, the MC analysis will also allow us to extend the study to the full phase-space, which includes the contribution of other diagrams. Given the similarity between neutral-and charged-current Drell-Yan production, and the observation that W boson production, split by charge, can serve as a good proxy also for the Z boson (see Fig. 2), and given the much larger MC sample available, the analysis has been restricted hereafter to the special case V = W ± About 8×10 7 events in the final state W ± → μ ± ν μ have been generated using the MG5_aMC@NLO [15] program interfaced with Pythia8 [16]. The dilepton mass is reconstructed using the muon momentum before QED final state radiation. The MC simulation is NLO accurate for observables inclusive in additional QCD radiation and it assumes M MC W = 80.419 GeV and Γ MC W = 2.047 GeV. As already anticipated in Sect. 2, the tree-level prediction is expected to be reproduced in the limit q T → 0, where q T is the transverse momentum of the boson. Indeed, in this regime the relative contribution of the tree-level 2 → 1 diagrams, which produce the boson at rest in the transverse plane, is enhanced compared to higher-order 2 → 2 diagrams. In contrast, a reduction of the shift in the large q T region is expected, where gluon-initiated diagrams dominate, thus reducing the sensitivity of the boson virtuality on the partonic luminosity.

Fit to the MC sample
The analytical study of Sect. 2 shows that the shift in the MC sample has to be extracted from a statistical analysis of the dilepton mass distribution dσ MC W /d Q. A crucial part of this task is to choose the correct functional form for dσ MC W /d Q, capable of modelling the lineshape without introducing a bias in the estimator of Δ V . Motivated by the tree-level study, an ansatz function of the same form of Eq. 4 has been chosen: The choice α = 1 defines the baseline function, which will be referred to as the modified Breit-Wigner. Indeed, this functional form, which is identical to Eq. 4 from the tree-level study, will be explicitly validated by checking that the estimator of M W and Γ W is consistent with the input values of the MC simulation M MC W and Γ MC W . As a further validation of this choice, two alternative instances of the parametric family of functions in Eq. 9 has been considered. The first is obtained by the choice α = H = K = 0, which reduces to a Breit-Wigner. This function is formally incorrect to model the dilepton mass distribution since it does not account for the Jacobian factor proportional to Q. However, it is a useful benchmark since it is symmetric around Q = M W so that the mass estimator also matches the mode. The second alternative function is obtained by choosing α = 1 and H = K = 0. This function, which will be referred to as Breit-Wigner with Jacobian, would be correct in the absence of the PDF distortion to the Q distribution. It peaks at Q ≈ M W + Γ 2 8M W , which is always larger than M W .
Three statistical analyses of the simulated events have been performed. The first analysis is inclusive in the phasespace of the boson and allows us to benchmark the different fit functions with the largest possible statistical precision. The second analysis is differential in the boson rapidity and is expected to reproduce, at least qualitatively, the ydependence from the tree-level model, as shown in Fig. 2. However, the comparison can only be approximate, since the latter predicts the transverse momentum q T to be identically zero, whereas the MC simulation generates a physical spectrum of transverse momenta. The third analysis is performed in bins of q T , and is inclusive in y. It allows us to both validate the tree-level calculation of Sect. 2 by extrapolating to q T → 0, and to study the dilution effect at larger values of q T . The fit parameters of Eq. 9 are determined by minimizing a χ 2 constructed using the event counts in each bin of the histogram in the range [79, 82] GeV and the value of the fit function at the center of the bin. Figure 3 shows the result of the three functional fits for the inclusive analysis for the W + sample. The W − sample shows very similar results for both the inclusive and the differential analyses, which are thus omitted for brevity. The quality of the fit improves dramatically when using the modified Breit-Wigner, with a reduced χ 2 of about 1.0 compared to 1.6 and 3.9 for the alternative functions. The best-fit value of M W when using the baseline function is consistent with the MC input within 1σ (1 ± 1 MeV), whereas the alternative functions depart from it by −9.9 ± 0.4 MeV and −21.5 ± 0.4 MeV, respectively. Likewise, the best-fit value of Γ W is consistent with the MC input values within 1.5σ (8 ± 5 MeV) for the baseline function, while it departs from it by 18 ± 1 MeV and 16 ± 1 MeV for the alternative functions. The residual discrepancy on Γ W when using the modified Breit-Wigner is ascribed to higher order terms in the power expansion, not included in Eq. 9. In order to validate this assumption, a toy MC has been used to verify that such a discrepancy is indeed consistent with the neglected terms that mostly enter through the tails of the distribution, whereas M, H and K are seen to be robust.
The best-fit value of M W from the differential analysis in the W boson rapidity are reported in Fig. 4a for the W + sam- The same applies to the Breit-Wigner with Jacobian function. For the latter, the discrepancy is even more pronounced. Indeed, the Jacobian factor contributes via a positive bias to the peak position. By neglecting the PDF term, which pulls in the opposite direction, the estimator of M W is thus shifted to even lower values compared to M W . The modified Breit-Wigner function correctly reproduces the input value M MC W in all bins of |y|, including the high |y| regimes, where the alternative functions perform rather poorly.
Finally, the results of the analysis differential in the W boson transverse momentum are shown in Fig. 4b. The modified Breit-Wigner is seen to correctly reproduce the input mass value for all bins of q T , whereas the two alternative functions disagree, especially at low transverse momenta. For q T in excess of about 40 GeV the statistical error is too large to discriminate among the models.

Extraction of Δ W
Since the fit reproduces well the true values for M W and Γ W , we fix the values of M and Γ to the MC input values in Eq. 6 and repeat the fit with A, H, K as the only free parameters.
For the inclusive sample (Fig. 3) the value of the shift is found to be: The first uncertainty is statistical-only while the second is the estimation of the systematic uncertainty from the PDFs. The latter is estimated from the RMS of the shifts determined using the first 100 replicas, as described in Sect. 2. In these fits the parameters M and Γ have been left free, since the uncertainty on the PDFs would be otherwise over-constrained by the imposed knowledge on the mass and the width of the resonance.
The fitted values of Δ W for the differential analyses are shown in Fig. 5 in bins of |y| and q T , separately for W + and W − . The variation of Δ W with the boson rapidity is shown in Fig. 5a. It agrees well with the tree-level expectation of a flat shift in the central rapidity region followed by a rapid decrease at larger rapidity values. However, the shift in the central region is found to be smaller by a factor of about two, like for the inclusive results. Such a difference has been interpreted as the result of the dilution from higher-order diagrams. Indeed, in the limit q T → 0, the measured shift gets closer to the tree-level result as shown by Fig. 5b, while it vanishes for q T in excess of about 40 GeV. A simple linear extrapolation to q T → 0 yields limiting values of (q T → 0 extrapolation) Δ W + = −10.1 ± 0.5 (stat.) ± 0.2 (PDF) MeV, Δ W − = −10.0 ± 0.6 (stat.) ± 0.2 (PDF) MeV. (11) Although reasonably close to the tree-level expectation, this result still disagrees with it by roughly 30%. This residual difference is interpreted as a pure next-to-leading-order correction to the leading-order prediction, stemming from collinear gluon emission and from gluon-initiated diagrams which contribute to the small-q T regime. The relative PDF uncertainty is found to agree reasonably well with the expectation from the tree-level model averaged over the W boson rapidity.
As a cross-check of this result, the fit has been repeated varying the range symmetrically by ±10%. The results for Δ W are stable, with a maximum discrepancy of 5%, which is within the uncertainty of the parameter. The fit has been also repeated after changing the renormalization and factorization scales in the matrix elements of the MC simulation by factors of 0.5 and 2, respectively. The results are again found to be stable within the PDF uncertainty.
Due to the agreement between the MC studies above and tree-level model of Sect. 2, the results for Δ W are expected to be valid also for the Z boson case with Δ W − < Δ Z < Δ W + and a similar PDF-related uncertainty.

Conclusions
In this paper, the impact of the PDFs on the lineshape of gauge bosons at the LHC has been investigated. Given the narrow width of the electroweak gauge bosons, the PDF impact can be treated, to a first approximation, as a shift Δ V of the mode of the dilepton mass spectrum from the boson mass M V . The origin of such shift can be traced back to the dependence of the partonic luminosity on the virtuality Q of the gauge boson. This effect is particular important for a possible future precise measurement of the Z boson mass since it directly affects the extraction of M V from the kinematics of the dilepton final state. It has been first studied analytically by using a tree-level model of Drell-Yan production and then validated by a statistical analysis of a MC simulated sample. The tree-level calculation agrees reasonably (a) (b) Fig. 5 The shift Δ W ± in bins of the W boson rapidity y (left) and transverse momentum q T (right). For the latter, a linear fit in the range [0, 40] MeV is performed to extrapolate the result to q T → 0.
The shaded boxes correspond to the PDF systematic uncertainty, as described in the text. On the right side of each plot, the equivalent scale for the H W parameter is reported well with the MC study in the phase-space where the two are expected to be comparable. The results of this study prove that the PDF uncertainty on Δ V is below one MeV all over the phase space relevant for future mass measurements at the LHC.