Precise predictions for the angular coefficients in Z-boson production at the LHC

The angular distributions of lepton pairs in the Drell-Yan process can provide rich information on the underlying QCD production mechanisms. These dynamics can be parameterised in terms of a set of frame dependent angular coefficients, Ai=0,…,7, which depend on the invariant mass, transverse momentum, and rapidity of the lepton pair. Motivated by recent measurements of these coefficients by ATLAS and CMS, and in particular by the apparent violation of the Lam-Tung relation A0 − A2 = 0, we perform a precision study of the angular coefficients at Oαs3\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \mathcal{O}\left({\alpha}_{\mathrm{s}}^3\right) $$\end{document} in perturbative QCD. We make predic-tions relevant for pp collisions at s=8\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \sqrt{s}=8 $$\end{document} TeV, and perform comparisons with the available ATLAS and CMS data as well as providing predictions for a prospective measurement at LHCb. To expose the violation of the Lam-Tung relationship we propose a new observable ΔLT = 1 − A2/A0 that is more sensitive to the dynamics in the region where A0 and A2 are both small. We find that the Oαs3\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \mathcal{O}\left({\alpha}_{\mathrm{s}}^3\right) $$\end{document} corrections have an important impact on the pT,Z distributions for several of the angular coefficients, and are essential to provide an adequate description of the data. The compatibility of the available ATLAS and CMS data is reassessed by performing a partial χ2 test with respect to the central theoretical prediction which shows that χ2/Ndata is significantly reduced by going from Oαs2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \mathcal{O}\left({\alpha}_{\mathrm{s}}^2\right) $$\end{document} to Oαs3\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \mathcal{O}\left({\alpha}_{\mathrm{s}}^3\right) $$\end{document}.


JHEP11(2017)003
Distributions for the angular coefficients can be extracted experimentally through fits to the measured final-state lepton kinematics, which can then be compared to the corresponding predictions obtained in pQCD. The measurement of these angular coefficients is therefore interesting in its own right and much effort has been devoted to their precise determination. Moreover, such a measurement also plays an important role in the determination of the W-boson mass M W at hadron colliders. Indeed, a precise extraction of M W requires control of the Monte Carlo samples used to describe the kinematic distribution of leptons resulting from W-boson decay. The approach to generating these samples (and/or reweighting them) can in part be validated by using the Z-boson production process as a case study, where the predicted values for all relevant angular coefficients can be directly compared to data.
On the theoretical side, the angular coefficients have been computed in pQCD up to O(α s ) [10][11][12][13] and O(α 2 s ) [14][15][16] for non-vanishing transverse momenta p T,Z of the Z boson. For the inclusive Drell-Yan process, the O(α 2 s ) corrections are available in the parton-level generators DYNNLO [22] and FEWZ [23], which retain the full kinematical information of the final state and allow for a direct comparison to data in the fiducial region. These fixed-order predictions have been further matched to parton showers at NNLO in ref. [24], where a comparison to the angular coefficients has also been performed. Using the results obtained with DYNNLO and FEWZ, a detailed comparison to all available hadron collider and fixed target data has been carried out in ref. [25]. Studies of the Lam-Tung relation in the context of the intrinsic transverse momentum of the parton have also recently been considered in ref. [26].
Experimentally, a number of the angular coefficients were determined in fixed target experiments by the NA10 [27], E615 [28], and FNAL E866/NuSea [29,30] collaborations using a variety of beams (pions, protons) and targets (tungsten, deuterium). It is worth noting that the kinematical range probed in these fixed-target experiments was restricted to small invariant masses of the lepton pairs, typically from a few GeV up to ∼ 15 GeV. In this regime, photon-exchange in the Drell-Yan process is by far the dominant contribution and only the parity-even angular coefficients could be determined.
At high-energy colliders such as the Tevatron and the LHC, on the other hand, leptonpair invariant masses around the Z-boson mass are considered, which are dominated by Z-boson exchange and also allow for the study of the parity-odd angular coefficients. The measurement of angular coefficients at hadron colliders were performed by the CDF [31] collaboration in pp collisions at a centre-of-mass (CoM) energy of √ s = 1.96 TeV, and more recently by the CMS [32] and ATLAS [33] collaborations in pp collisions at √ s = 8 TeV. Each of these analyses were performed in an invariant-mass window around the Z-boson resonance and in the Collins-Soper reference frame [10]. Most notably, both ATLAS and CMS observe for the first time clear evidence for the violation of the Lam-Tung relation in Z-boson production. 1 The new results from ATLAS and CMS are therefore particularly JHEP11(2017)003 interesting as they demonstrate the violation of the Lam-Tung relation at energies never probed before.
However, compared to the fixed-order O(α 2 s ) prediction for lepton-pair production using the fixed order parton level code DYNNLO [22], a "significant deviation" is reported by the ATLAS collaboration [33] for the difference (A 0 − A 2 ) in the region with p T,Z > 20 GeV. Although less significant, a similar trend is also observed in the CMS data [32] where the O(α 2 s ) prediction for the Drell-Yan pair production reaction is obtained using the parton-level generator FEWZ [23]. Both experiments observe that the data exceeds the corresponding theory prediction for this observable. A tension is also observed in the p T,Z spectrum for the angular coefficient A 2 , where the data tends to undershoot the theory prediction. It is worth noting that although both FEWZ and DYNNLO yield predictions which are accurate at next-to-next-to-leading order (NNLO) for the inclusive Z-boson production cross section, in analogy to the case for the p T,Z or φ * η distributions studied in refs. [34,35], these codes provide NLO-accurate predictions for the angular coefficients A i and LO-accurate predictions for the difference (A 0 − A 2 ) (since A 0 = A 2 up to O(α s ) by virtue of the Lam-Tung relation).
The purpose of this work is to reassess the compatibility of the LHC data to theory by providing predictions for the phenomenologically most important angular coefficients in high-mass lepton pair production at O(α 3 s ), while focussing on the kinematic region with p T,Z > 10 GeV, where many of the angular coefficients start to acquire non-vanishing values. This accuracy is achieved through the calculation of the Z+jet process at O(α 3 s ) [36] at finite p T,Z without requiring any resolved jets in the final state.
The layout of the paper is as follows. In section 2, the theoretical formalism for decomposing Z-boson production in terms of angular coefficients and spherical harmonics is discussed. We further propose a new observable ∆ LT , which is particularly suited to study the violation of the Lam-Tung relation in the low-p T,Z regime. Numerical predictions for these angular coefficients are provided in section 3, along with a detailed comparison to the available LHC data. In addition, predictions for the LHCb experiment are provided, for which no measurement is available at present. A summary of our findings and concluding remarks are presented in section 4.

Theoretical preliminaries
We consider the inclusive production of lepton pairs through the decay of an intermediate gauge boson, p(p 1 ) + p(p 2 ) → V (q) + X → (k 1 ) +¯ (k 2 ) + X as depicted in figure 1. The cross section for this process can be written as the contraction of a lepton tensor (L µν ) describing the final-state decay with a hadronic tensor (H µν ) that describes the production sub-process, namely L µν H µν . The lepton tensor in this context takes the role of an analyser, providing a probe of the structure of H µν . Note that the definition of the hadronic tensor includes the convolution with the PDFs as well as the integral over any degrees of freedom associated with the hadronic recoil "+X". As a result, H µν only depends on the four-momenta p 1 , p 2 , and q. Based on Lorentz-and gauge-invariance, the -3 -
It is interesting to note that lepton-pair production satisfies an analogous relation to the Callan-Gross relation in deep-inelastic scattering (DIS) known as the Lam-Tung relation [11][12][13], This relation, formulated in a covariant manner, is frame independent and characteristic of the spin-1 2 nature of the quark. It has been further shown [13] that eq. (2.2) is not affected by O(α s ) QCD corrections, 3 which follows as a direct consequence of the vector-coupling of the spin-1 gluon to quarks [37]. However, this relation has been shown to be violated at O α 2 s [16]. As such, the Lam-Tung relation offers a unique opportunity to study the pQCD predictions of the underlying dynamics encoded in H µν in more detail than through rate measurements alone.
To further elucidate the Lam-Tung relation, let us consider the kinematics of this process in the lepton-pair rest frame where the final-state lepton momenta can be expressed in terms of the angles θ and φ: where so far the orientation of the coordinate axes remains unspecified. The only nonvanishing entries of the hadronic tensor (2.1) are the space-space components H ij (H µ0 = H 0ν = 0) , where dσ unpol. is the unpolarised cross section. We note that the first term inside the parenthesis equal to (1 + cos 2 θ) is not accompanied by a separate angular coefficient, as its normalisation is described by dσ unpol. that has been extracted as a global pre-factor in eq. (2.5). The unpolarised cross section is given by the trace of the hadronic tensor and for Z exchange explicitly reads with α denoting the fine-structure constant and g ± the chiral Z-boson couplings to charged leptons. 4 These couplings are defined according to where θ w is the weak mixing angle. The form factors H 1,...,9 in eq. (2.1), or equivalently the nine non-vanishing components H ij of the hadronic tensor, are directly related to the eight angular coefficients A 0,...,7 and the unpolarised cross section. Explicitly, the A i are given by where the proportionality factors c ± arise from the fact that dσ unpol. has been removed as a prefactor in the definition (2.5) and are given by The unpolarised cross section (2.6) completes eq. (2.8) as the ninth linearly independent combination of the H ij .
Let us now choose a specific reference frame by defining the direction of the axes in the lepton-pair rest frame. To this end, we consider the Collins-Soper frame [10] shown  The definition of the Collins-Soper [10] angles in the di-lepton rest frame. in figure 2: the z-axis is chosen as the external bisector of the incoming beam directions, e CS z ∼ ±( p 1 − p 2 ), where the positive z-direction is aligned with the z-direction of the lepton pair in the laboratory frame. The x-axis lies in the hadron plane orthogonal to the z-axis and points in the direction ofê CS x ∼ −( p 1 + p 2 ). Lastly, the y-axis is chosen to complete a right-handed Cartesian coordinate system and is orthogonal to the hadronic event plane. The four-momenta of the incoming hadrons in this reference frame are given by 5 (2.10) Returning to the Lam-Tung relation (2.2), one can derive the corresponding relation in terms of the angular coefficients A i in the Collins-Soper frame where we have used for the non-vanishing diagonal components of the hadronic tensor. We observe that the Lam-Tung relation is equivalent to A 0 − A 2 = 0. Note that the result of eq. (2.11) is not JHEP11(2017)003 frame independent but only holds if both the z-and x-axis in the lepton-pair rest frame lie in the hadronic event plane. This condition enters in the step where the form factors H i are expressed in terms of the diagonal H ij components using eq. (2.12) and can be understood by inspecting the covariant formulation of eq. (2.2) in the lepton-pair rest frame: the only form factors that contribute to the trace of the hadronic tensor on the r.h.s. of eq. (2.2) are H 1,...,4 . The tensor structures multiplying H 2,3,4 only involve momenta lying inside the hadronic plane and it is solely the tensorg µν multiplying the form factor H 1 which has a non-vanishing component orthogonal to it. The Lam-Tung relation therefore distinguishes the direction perpendicular to the hadronic plane and can be interpreted as a statement about the current-current correlation of the hadronic tensor in this direction. 6 Making use of the completeness of the spherical harmonics, the angular coefficients appearing in the decomposition provided in eq. (2.5) can be extracted through the projectors where . . . denotes taking the (normalised) weighted average over the angular variables θ, φ and is defined as . (2.14) The dominant angular coefficients are A 0,...,4 , while A 5,6,7 vanish at O(α s ) and only receive small O α 2 s corrections from the absorptive parts of the one-loop amplitudes in Z + jet production. We therefore will not discuss the coefficients A 5,6,7 in the following. In the case of pure γ * exchange, the relevant coefficients are the parity-conserving coefficients A 0,1,2 . A 3 and A 4 , on the other hand, are odd under parity and proportional to the product of vector-and axial-vector-couplings of the gauge boson to the fermions. As such, they are sensitive to the relative rate of incoming down-and up-type quark fluxes as well as the weak mixing angle s w . All the coefficients A i vanish in the limit p T,Z → 0 with the exception of A 4 , which is finite in this limit and directly related to the forward-backward asymmetry.
One of the goals of this work is to assess the compatibility of the observed extent of the Lam-Tung violation with that expected in predictions based on pQCD. This can be done by directly studying the p T,Z distribution for the difference of the angular coefficients A 0 and A 2 . Here, we propose a new observable which has the benefit that the strong suppression of the individual angular coefficients in the low-p T,Z region is absent. In addition, the dependence on the unpolarised cross section

JHEP11(2017)003
appearing in the denominator of eq. (2.14) cancels in the ratio between the two coefficients. Consequently, this observable may help to expose the extent of the Lam-Tung violation in this region. In section 3, we shall compare our predictions to the available ATLAS and CMS data for the p T,Z distributions of both (A 0 − A 2 ) and ∆ LT . In the latter case, the data will be re-expressed in terms of ∆ LT .

Numerical predictions
In this section, we provide a comparison of the predictions for a set of angular coefficients to the available ATLAS and CMS data in pp collisions at √ s = 8 TeV. While the LHCb collaboration has not yet performed a measurement of the angular coefficients, previous measurements of the Z-boson p T spectrum [7][8][9] and forward-backward asymmetry [21] indicate that there is potential for such a measurement in the forward region. We therefore also provide predictions in the LHCb fiducial region at √ s = 8 TeV. In all cases (ATLAS, CMS, LHCb), the angular coefficients are defined in the Collins-Soper reference frame [10]. In addition to the angular coefficients, we also provide absolute predictions for the unpolarised p T,Z distributions. 7 Special attention is also given to the difference (A 0 − A 2 ), where the quality of theoretical description with respect to the observed distributions is quantified by means of a χ 2 test. Furthermore, we also present a comparison to data for the new observable ∆ LT .
The measurements of the angular coefficients are performed differentially in p T,Z and for various rapidity intervals, where in all cases an invariant-mass window for the leptonpair final state is imposed around the Z-boson resonance. For non-vanishing values of p T,Z , the LO prediction for this distribution can be obtained from the O(α s ) tree-level Z + jet process, where the transverse momentum of the Z boson is balanced with that of a single final-state QCD parton. The NLO QCD and EW corrections to this process have been computed in refs. [38,39], and more recently the NNLO QCD corrections to this process have been completed [36,40]. In this work, we employ the calculation of ref. [36] based on the antenna subtraction formalism [41][42][43][44][45][46][47][48][49] to provide NNLO-accurate QCD predictions for the p T,Z distributions of the angular coefficients. This process is implemented in the flexible parton-level Monte Carlo generator NNLOjet.
The predictions are provided in the G µ -scheme, where we take the following set of numerical inputs: M os Z = 91.1876 GeV, Γ os Z = 2.4952 GeV, M os W = 80.385 GeV, Γ os W = 2.085 GeV, and G µ = 1.16638 · 10 −5 GeV −2 . In the extraction of the corresponding numerical values for α and s 2 w , we additionally include the dominant one-and two-loop universal corrections to the ρ-parameter [50] which relate M W −M Z interdependence present beyond tree-level. Including these contributions leads to the effective values of α eff. = 0.007779, s 2 w,eff. = 0.2293. As a baseline PDF set, we use the central member of PDF4LHC15_nnlo_30 [51][52][53][54][55][56], and extract α s from the grid provided with the PDF set -corresponding to α s (Q = M os Z ) = 0.118.

JHEP11(2017)003
As discussed in section 2, the theoretical predictions for the coefficients A i can be obtained by computing the normalised expectation values of the spherical harmonics according to eq. (2.14). To assess the theoretical uncertainty in the extraction of the coefficients through this method, various solutions are possible. In this paper, for comparison with LHC data, we choose to perform an independent variation of factorisation (µ F ) and renormalization (µ R ) scales in both the numerator and denominator of this expression. The scales µ num. R , µ num. F , µ den. R , and µ den.
F are each independently varied by factors of 1 2 and 2 about the transverse energy of the lepton pair, with the constraint that all pairs of these uncorrelated scales satisfy 1 In total this corresponds to 31 possible combinations and the associated uncertainty is obtained as the envelope around the central scale µ num.
An alternative approach is to correlate the scale uncertainties between numerator and denominator. However, this treatment can lead to an underestimation of the uncertainty due to missing higher-order effects. For example, at LO the renormalisation scale dependence is fully encapsulated in the strong coupling α s (µ) which entirely cancels if the scales between numerator and denominator are correlated.
To further demonstrate this point, we show the impact of these two different approaches in figure 3 where the p T,Z distribution for A 2 is evaluated at NLO and NNLO. The distributions are obtained with an invariant-mass cut of 80 < m < 100 GeV on the lepton-pair final state and inclusively with respect to rapidity of the lepton pair. It is clearly seen that at NLO, these two prescriptions result in substantial differences in the scale uncertainty bands, with the correlated approach yielding considerably smaller uncertainty bands for p T,Z above 20 GeV. At NNLO however, the scale uncertainty estimate of the A i coefficients obtained with either choice gives similar results in the low-pt region (p T,Z < 80 GeV) where the NNLO effects are largest. Throughout this work, all distributions are obtained with the uncorrelated prescription discussed above.
It is worth commenting that the angular coefficients are evaluated differentially in p T,Z and in multiple kinematic regions, corresponding to the various experimental setups. In addition, each of these coefficients are computed through the projectors in eq. (2.13), which are highly oscillating functions with respect to the leptonic kinematics, and consequently their stable numerical evaluation is rather challenging. This is particularly true for the difference (A 0 − A 2 ) for which large non-local cancellations occur at the level of the angular coefficient as opposed to the integrand.

Comparison to ATLAS data
The ATLAS measurements have been performed with an invariant-mass cut of 80 < m < 100 GeV on the lepton-pair final state, and distributions for the angular coefficients in p T,Z were extracted for a range of different accessible rapidity regions. The data has also been presented integrated in y Z , obtained after performing an extrapolation to the full phasespace region. The comparison of the theoretical predictions is performed with respect  to this data, referred to as y Z inclusive, for which the measurements are most precise. The predictions of the various p T,Z distributions are provided for the kinematic range of p T,Z ∈ [11.4, 600] GeV, and compared to the available data in this region. Furthermore, in the region of p T,Z < 85.4 GeV we provide our predictions with a coarser choice of p T,Z bins with respect to the data, which are obtained by pair-wise combining neighbouring bins. Before continuing, it is also important to highlight that we perform the comparison to the ATLAS data which is obtained prior to the regularisation procedure outlined as part of the experimental analysis -more detail can be found in appendix C of [33]. Our motivations for doing so are as follows. Firstly, the regularisation procedure introduces large bin-to-bin correlations for the distributions of A i coefficients meaning that a visual comparison to the regularised data can be misleading as large correlations are hidden from view. Secondly, we wish to quantify the agreement between theory and data by performing a χ 2 test, which requires knowledge of the bin-to-bin correlations between the different A i coefficients in p T,Z (this is particularly important if these correlations are large, which is the case for the regularised data). However, to our knowledge, a well-defined covariance matrix for the regularised version of the A i coefficients is not available.
The p T,Z distributions for the angular coefficients A 0 (upper left), A 1 (lower left), A 2 (upper right), and the unpolarised cross section (lower right) are shown in figure 4. The ATLAS data is represented by black points, and is compared to theoretical predictions at LO (blue), NLO (green), and NNLO (red). In the lower panel of each plot, the same distributions are shown normalised with respect to the central NLO prediction.
The NNLO corrections are observed to have an important impact on each distribution, substantially reducing the scale uncertainties in all cases. With respect to the central value at NLO, the NNLO corrections to A 0 are negative and typically below 5% in magnitude.  In the case of the A 2 distribution, the corrections are also negative and most sizeable in the region of p T,Z ∈ [10, 50] GeV. The description of the observed A 2 distribution is visibly improved at NNLO, while the NLO predictions systematically overestimate the data. It should be noted that the y-axis ranges for both A 0 and A 2 distributions are fixed to the same values to allow a straightforward visual comparison of the relative impact of the NNLO corrections in each case (this is also the reason why they are placed in neighbouring positions within the figures). In the case of the A 1 distribution, the corrections are positive at low p T and change sign to become negative at large p T , resulting in a modified shape of the distribution. The size of the corrections across the whole p T,Z range vary between +10% at low values of p T,Z and −5% in the last p T,Z bin shown. In figure 5, the same comparison is performed for the parity-violating angular coefficients A 3 (left) and A 4 (right). The NNLO corrections reduce the scale uncertainty of the prediction, while having little impact on the central value. With respect to the experimental precision, the NNLO corrections to these distributions (which are well described by the central NLO prediction) are phenomenologically unimportant for a comparison to data. As discussed in section 2, these particular coefficients are sensitive to the product of vector-and axial-vector-couplings of the Z-boson to the initial-state quarks. The corresponding predictions for these distributions are therefore sensitive to a combination of the input value of s 2 w as well as to the relative contribution of up-and down-type quark initiated processes. To the accuracy of the experimental distributions for these coefficients, our choice of input parameter scheme (including universal corrections via the ρ-parameter) provides a consistent description of the data. However, if the precision of future measurements of these coefficients improves, it would be important to revisit this comparison while including possibly also the effect of electroweak corrections and to assess the impact of PDF uncertainties on these distributions. We note that while a measurement of these coefficients is sensitive to the weak mixing angle, a more direct extraction of this parameter is possible through the measurement of the forward-backward asymmetry in lepton-pair production. Indeed, such a measurement has already been performed by the ATLAS collaboration [20].

Comparison to CMS data
A similar measurement of the angular coefficients has also been presented by the CMS collaboration [32]. In this case, the angular coefficients have been measured differentially within rapidity bins of |y Z | ∈ [0.0, 1.0] and |y Z | ∈ [1.0, 2.1], and with an invariant-mass window of 80 < m < 100 GeV on the lepton-pair final state. In the following, we perform For the parity-violating angular coefficients A 3 and A 4 , we see that the NNLO corrections do not alter the shapes of these distributions for the CMS kinematical setup (as was the case for the y Z -inclusive distributions for ATLAS shown in figure 5). Figure 8 shows that these distributions are well approximated by the central NLO predictions. The NNLO corrections to these distributions, as compared to the accuracy of the data, are phenomenologically unimportant.

Predictions for LHCb
The LHCb collaboration has not yet performed a measurement of the angular coefficients in Z-boson production. Such an analysis would however be of interest to provide a probe of the Z-boson production mechanisms which may be enhanced at forward rapidities, and is also an important stepping-stone towards performing an extraction of M W within the LHCb acceptance [57]. We therefore provide predictions for the p T,Z distributions for the LHCb fiducial region of y Z ∈ [2.0, 4.5], placing an invariant mass selection of 80 < m < 100 GeV on the lepton-pair final state. No other cuts are placed on the lepton-pair final state as it is assumed that this will be corrected for in the experimental analysis. The predictions for A 0 (upper left), A 1 (lower left), A 2 (upper right), as well as the unpolarised cross section (lower right) are provided in figure 9. Each p T,Z distribution is provided in the region p T,Z ∈ [10.5, 270] GeV, guided by the choice of binning taken in the recent LHCb measurement of Z-boson production at 13 TeV [9].
The predicted shapes of the distributions within the LHCb acceptance are similar to what is observed at more central rapidities. In addition, the NNLO corrections for each of the angular coefficients are also observed to be of similar size to those at more central rapidities. Although not shown here, the NNLO corrections to A 3 and A 4 were also computed for this kinematic setup and found to be negligibly small. We can therefore conclude that the NNLO corrections to the p T,Z spectrum for A 0 , A 1 , and A 2 should be included when performing a comparison to data, while the central NLO prediction -16 -

JHEP11(2017)003
for the corresponding A 3 and A 4 distributions are likely sufficient. If an experimental determination of the A 3 and A 4 distributions is achievable at LHCb with smaller relative uncertainties as compared to ATLAS and CMS, it would be important to include the effects of electroweak corrections and to also assess the impact of PDF uncertainties on these distributions.

Assessing the violation of the Lam-Tung relation
As highlighted in section 2, the Lam-Tung relation is expected to be violated, i.e. (A 0 − A 2 ) = 0, starting at order O(α 2 s ) in the framework of pQCD. Measurements of this violation therefore provide an important test of the Z-boson production dynamics. As discussed in section 2, the extent of the breaking of the Lam-Tung relation can be assessed by measuring the p T,Z distribution for the difference of the angular coefficients A 0 and A 2 , or equivalently through the normalised observable ∆ LT . The latter has the benefit of better exposing the violation of the Lam-Tung relation in the lower p T,Z range, where the angular coefficients A 0 and A 2 are individually relatively small. In the following, we discuss the corrections to (A 0 − A 2 ) and quantify the consistency of the data and the predictions by performing a χ 2 test. We then present results for ∆ LT and perform a comparison to data, where the data points for the latter are obtained by re-expressing ∆ LT in terms of the measured angular coefficients.
Before comparing to data, it is important to comment on the expected accuracy of our theoretical predictions for these two observables. As for the individual angular coefficients, our theoretical predictions for (A 0 − A 2 ) are obtained from the computation of the production process for Z + jet at O(α 3 s ). While the O(α 3 s ) contributions comprise genuine NNLO corrections to the individual A i coefficients as demonstrated throughout this section, the prediction degrades to an NLO-accurate description for the difference (A 0 − A 2 ) and ∆ LT . 8 For consistency with the rest of the paper, we will continue to refer to the corrections of order O(α 3 s ) as "NNLO corrections" and similarly label the figures in this section as NNLO predictions. Figure 10 shows the p T,Z distribution for (A 0 − A 2 ), where the ATLAS data is represented by black points, and is compared to LO (blue), NLO (green), and NNLO (red) theoretical predictions. The NNLO corrections are observed to be large and positive, amounting to +40% at moderate p T,Z values, and provide an improved description of the ATLAS data. It is worth noting that while a reduction of the absolute scale uncertainties is already observed at NNLO with respect to NLO, the relative uncertainty is in fact reduced by almost a factor of two across the shown p T,Z range. This is a reflection of the fact that the computation of the Z+jet-production process at O(α 3 s ) at finite p T,Z used to predict the difference (A 0 − A 2 ) is only NLO accurate and therefore yields corrections and remaining scale uncertainties which are typical for NLO effects. For most of the p T,Z range, and for both experimental setups, it is found that the NLO and NNLO predictions for (A 0 − A 2 ) are consistent within uncertainties. In figure 10, we have also chosen to include the regularised ATLAS data (indicated by the grey fill). As discussed towards the start of this JHEP11(2017)003 section, large bin-to-bin correlations are introduced in the regularisation procedure of the ATLAS data (see for example figure 24 of ref. [33]). We believe that this demonstrates how a visual comparison of the theory prediction with respect to the regularised data (in this case at least) can lead one to overestimate the disagreement between theory and data.
As an alternative to a visual comparison, the quality of the theoretical description of the data can be quantified by performing a χ 2 test according to where O i exp and O i th. are respectively the central value of the experimental and theoretical predictions for data point i, and σ −1 ij is the inverse covariance matrix. In this case, we consider the unregularised ATLAS data for the angular coefficients A 0 and A 2 (and their correlations, as provided through the covariance matrix) within the range of p T,Z ∈ [11.4, 600] GeV, corresponding to a total of 38 data points. To perform this comparison, the central theory predictions for these angular coefficients are evaluated with the same binning choice as the data. This test indeed demonstrates that the NLO predictions give a poor description of the data in the considered p T,Z range, a point that was also highlighted in the experimental analysis [33]. This tension is largely reduced with the inclusion of the NNLO corrections, and from closer inspection of figure 4, can be mainly attributed to the large negative corrections to the A 2 distribution. The corresponding p T,Z distributions for the CMS measurement are shown in figure 11 for the rapidity bins |y Z | ∈ [0.0, 1.0] (left) and |y Z | ∈ [1.0, 2.1] (right). The NNLO corrections to (A 0 − A 2 ) exhibit a similar behaviour for the CMS kinematic selections, and again improve the description of data. This agreement can also be quantified by performing a χ 2 test, where in this case the test is performed directly on the (A 0 − A 2 ) distribution as no covariance matrix for these A i coefficients is publicly available. In total 14 data points are considered, corresponding to seven p T,Z bins for each rapidity selection. The results, assuming uncorrelated bins, are Similar to the findings for the ATLAS data, the description of the CMS data is substantially improved at NNLO. As discussed previously, it is also informative to express the data in terms of the new obserable ∆ LT as defined in eq. (2.15). This comparison is performed in figures 12 and 13 for the ATLAS and CMS measurements, respectively, where the data has been re-expressed in terms of this quantity. 9 It is found that the extent of the Lam-Tung violation observed 9 We omit the lower panels with the K-factors in these figures, as they are almost identical to the case of (A0 − A2) shown in figures 10 and 11 due to the small corrections to the A0 coefficient.  in data is consistently described by the NNLO predictions. While there is some tendency for the data to prefer a stronger Lam-Tung violation for p T,Z > 40 GeV, more precise data is required to confirm this behaviour.

Conclusions and outlook
Using our calculation of the Z + jet process at NNLO [36], we have computed the p T,Z distributions for the angular coefficients in Z-boson production to O(α 3 s ). We have focussed on the phenomenologically most relevant angular coefficients A i=0,...,4 for pp collisions at √ s = 8 TeV and have compared them with available LHC data. With the theory uncertainties estimated by the uncorrelated variation of the factorisation and renormalisation scales (as described at the beginning of section 3), we find that these coefficients display a good perturbative convergence. In particular, a reduction of scale uncertainties is observed at each successive order and the residual scale uncertainties at NNLO are typically at the level of 5%. The NNLO corrections are observed to have an important impact on the predicted shapes of the distributions for A 0 , A 1 , and A 2 . Of particular note is that the corrections to the A 2 distribution are both large and negative (up to −20%) in the direction of data for both ATLAS and CMS measurements. It is found that the impact of -20 -

JHEP11(2017)003
the NNLO corrections to A 3 and A 4 distributions is small, and that these distributions are well described by the central NLO prediction. Besides comparing our predictions to the available LHC data for these coefficients, we have also provided predictions for the p T,Z distributions for A i=0,1,2 within the LHCb fiducial region, which would allow to probe the Z-boson production mechanism at forward rapidities. We find that the corrections to these distributions exhibit a similar behaviour to that at central rapidities both in size and shape. Particular emphasis has been placed also on testing the consistency between the Lam-Tung violation observed in CMS and ATLAS data with respect to the theory predictions. To this end, we have studied the p T,Z distributions for the observable (A 0 − A 2 ) directly, where the quality of the data-theory comparison has been assessed through a χ 2 test. Here, the inclusion of NNLO corrections leads to a significant improvement in the χ 2 /N dat. values: in the case of ATLAS, the χ 2 /N dat. value reduces from 4.89 at NLO to 1.80 at NNLO; whereas for CMS, a reduction from 1.75 at NLO to 1.01 at NNLO is observed. With respect to the NNLO predictions, no significant deviation is observed for the ATLAS data, and the CMS data is found to be fully consistent.
We further introduced a new observable ∆ LT defined in eq. (2.15), which is designed to better expose the violation of the Lam-Tung relation in the lower p T,Z regime. Expressed through this quantity, it becomes clear that the extent of Lam-Tung violation observed within the range of p T,Z ∈ [10, 40] GeV, where this effect is the strongest, is consistent with the NNLO predictions. There however still remains some tendency for the data to systematically exceed the corresponding predictions at larger p T,Z values. More precise data is required to clarify this situation.
Throughout this work, we have shown how the NNLO QCD predictions obtained via the calculation of the Z + jet process at O(α 3 s ) are essential to provide an adequate description of the p T,Z distributions of several angular coefficients present in Z-boson production. It is therefore likely that a similar statement will also apply to the case of W-boson production. At present, a precise extraction of M W at the LHC relies on an accurate modelling of the corresponding angular coefficients in W-boson production [58] based on the fixed-order O(α 2 s ) prediction. Our studies indicate that this level of theoretical accuracy is inadequate. The O(α 3 s ) corrections to the decay lepton distributions in vector-boson production computed here are providing an important step towards improving the theoretical description of reference quantities necessary for the precise measurement of M W .