W production at LHC: lepton angular distributions and reference frames for probing hard QCD

Precision tests of the Standard Model in the Strong and Electroweak sectors play a crucial role, among the physics program of LHC experiments. Because of the nature of proton–proton processes, observables based on the measurement of the direction and energy of final state leptons provide the most precise probes of such processes. In the present paper, we concentrate on the angular distribution of leptons from W→ℓν\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$W \rightarrow \ell \nu $$\end{document} decays in the lepton-pair rest-frame. The vector nature of the intermediate state imposes that distributions are to a good precision described by spherical harmonics of at most second order. We argue, that contrary to general belief often expressed in the literature, the full set of angular coefficients can be measured experimentally, despite the presence of escaping detection neutrino in the final state. There is thus no principle difference with respect to the phenomenology of the Z/γ→ℓ+ℓ-\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z/\gamma \rightarrow \ell ^+ \ell ^-$$\end{document} Drell–Yan process. We show also, that with the proper choice of the reference frames, only one coefficient in this polynomial decomposition remains sizable, even in the presence of one or more high pT\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_T$$\end{document} jets. The necessary stochastic choice of the frames relies on probabilities independent from any coupling constants. In this way, electroweak effects (dominated by the V-A\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$V-A$$\end{document} nature of W couplings to fermions), can be better separated from the ones of strong interactions. The separation is convenient for the measurements interpretation.


Introduction
The main purpose of the LHC experiments [1,2] is to search for effects of New Physics. This program continues after the breakthrough discovery of the Higgs boson [3,4] and measurement of its main properties [5]. In parallel to searches of New Physics, see e.g. [6][7][8], a program of precise measurements in the domain of Electroweak (EW) and Strong (QCD) a e-mail: Z.Was@cern.ch interactions is on-going. This is the keystone for establishing the Standard Model as a consistent theory. It is focused around two main directions: searches (setting upper limits) for anomalous couplings and precision measurements of the Standard Model parameters. Precision measurements of the production and decay of Z and W bosons represent the primary group of measurements of the latter domain, see e.g. [9][10][11][12][13]. The study of the differential cross-sections of W production and decays is essential for understanding open questions related to the electroweak physics, like the origin of the electroweak symmetry breaking or the source of the CP violation.
Since discovery, the W boson mass and width, have been measured to a great precision in pp and pp collisions [14]. To complete physical information on the production process, measurements pursued the boson's decay differential distributions too. The measurements rely on outgoing leptons of the W → ν decays in the W -boson rest frame. Because electroweak interaction of the decay vertex are known with much better precision than the QCD interaction of the production, the measurement predominantly tests dynamics of QCD imprinted in the angular distributions of outgoing leptons.
In principle, the same standard formalism of the Drell-Yan production Z → [15] can be applied in the case of W → ν production [16,17]. The angular dependence of the differential cross-section can be written again as dσ dp 2 T dY dΩ * = Σ 9 α=1 g α (θ, φ) where the g α (θ, φ) represent harmonic polynomials of the second order, multiplied by normalisation constants and by dσ α which denote helicity cross-sections, corresponding to nine helicity configurations of W matrix elements. The angle θ and φ in dΩ * = d cos θ dφ are the polar and azimuthal decay angles of the charged lepton in the W rest-frame. The p T , Y denote transverse momenta and rapidity of the inter-mediate W boson in the laboratory frame. The z-axis of the W rest frame can be chosen along the W momentum of the laboratory frame (the helicity frame), or constructed from the directions of the two beams (the Collins-Soper frame [18]). We rewrite Eq. (1) explicitly, defining polynomials and corresponding coefficients dσ dp 2 T dY d cos θ dφ + 1/2 A 0 (1 − 3 cos 2 θ) + A 1 sin(2θ) cos φ + 1/2 A 2 sin 2 θ cos(2φ) + A 3 sin θ cos φ + A 4 cos θ + A 5 sin 2 θ sin(2φ) + A 6 sin(2θ) sin φ + A 7 sin θ sin φ] where dσ U +L denotes the unpolarised differential crosssection (a convention used in several papers of the 1980s). In case of W boson, (θ, φ) define the orientation of the charged lepton from W → ν decay. The coefficients A i ( p T , Y ) are related to ratios of corresponding cross-sections for intermediate state helicity configurations. The full set of A i coefficients has been explicitly calculated for pp → W (→ ν) + 1 j at QCD NLO in [16,17]. The first term at Born level (no jets): (1 + cos 2 θ) results from spin 1 of the intermediate boson. The dynamics of the production process is hidden in the angular coefficients A i ( p T , Y ). This allows to treat the problem in a model independent manner. In particular, as we will see, all the hadronic physics is described implicitly by the angular coefficients and it decouples from the well understood leptonic and intermediate boson physics. Let us stress, that the actual choice of the orientation of coordinate frames represents an important topic; we will return to this later.
The understanding of how QCD corrections affect lepton angular distributions is important in the measurement of the W mass (m W ), independently of whether leptonic transverse momentum or transverse mass (m W T ) of the W are used. In fact, the first measurements of the angular coefficients explored this relation in the opposite way. Assuming the mass of the W boson measured by LEP, from the fit to transverse mass distribution of the lepton-neutrino system m W T , information on the angular orientation of the outgoing leptons was extracted.
The cross-section has been parametrised [19] using only the polar-angle (i.e. integrating over azimuthal angle) as dσ d cos θ ∼ (1 + α 1 cos θ + α 2 cos 2 θ) (3) with the following relations between coefficients; It has been estimated that a 1% uncertainty on α 2 corresponds to a shift of the measured m W in pp collision, determined by fitting the transverse mass distribution, of approx-imately 10 MeV. The α 1 measures the forward-backward decay asymmetry. The measurements of α 2 at 1.8 TeV pp collisions have been conducted by D0 and CDF experiments and published in [20,21]. It was based on the data collected in 1994-1995 by Fermilab's Tevatron Run Ia. The fit to m W T was performed in several ranges of the W boson transverse momentum. The measurements confirmed standard model (SM) expectations, that α 2 decreases with increasing W boson transverse momentum, which corresponds to increase of the longitudinal component of the W boson polarisation. The ratio of longitudinally to transversely polarised W bosons in the Collins-Soper W rest frame increases with the W transverse momentum at a rate of approximately 15% per 10 GeV.
With more data collected during Fermilab Tevatron Run Ib, the measurement of the W angular coefficients was performed using a different technique; through direct measurement of the azimuthal angle of the charged lepton in the Collins-Soper rest-frame of the W boson [22]. The strategy of this novel measurement was documented in a separate paper [23]. Because of the ambiguity on determining the sign of cos θ (due to neutrino momenta escaping detection) which was not resolved, only the measurement of the coefficients A 2 and A 3 was performed and angular coefficients were measured as function of the transverse momentum of the W boson. The measurement was performed specifically for the W − bosons; angular coefficients of the W + were obtained by CP transformation of Eq. (2).
The pure V − A interactions of W ± without QCD effects, lead for pp collisions to α 2 = 1.0 and α 1 = 2.0, thus to pure transversely polarised W boson. This assumes that the W boson is produced with no transverse momenta, and seaquarks and gluon contributions to the structure functions can be neglected. Such a simple parton-model could guide intuition for the pp collisions at Tevatron, but had to be revisited for the pp collisions at LHC.
The dominance of quark-gluon initial states, along with the V-A nature of the coupling of the W boson to fermions, implies that at the LHC W bosons with high transverse momenta are expected to exhibit a different polarisation as the production mechanism is different at low p W T and high p W T [24,25]. W bosons produced with low p W T , and therefore moving generally along beam axis, exhibit a left-handed polarisation [26]. This is because the W -boson couples, in the dominant production diagram, to the left-handed component of valence quarks, and to the right-handed one of the sea anti-quarks. At high p W T , the situation becomes more complex due to contributions of higher-order processes. Of special interest, to quantify the validity of the QCD predictions, becomes the behavior of polarisation fractions as function of p W T . It was recently pointed out in [27], that events with high p W T can test the absorptive part of the scattering amplitudes and hence offer a non-trivial test of perturbative QCD at one and higher-loop levels. In all p W T ranges, the production at LHC therefore displays new characteristics: asymmetries in charge and momentum for W bosons and their decay leptons.
The LHC experiments pursued measurement techniques different than Tevatron. With 7 TeV data of pp collisions, the helicity frame and not the Collins-Soper frame was used. The interest was not to measure A i coefficients directly but rather the helicity fractions, f 0 , f L , f R . The helicity state of the W boson becomes a mixture of the left and right handed states, whose proportions are respectively described with fractions f L and f R . The f 0 denotes the fraction of longitudinally polarised W bosons, which is possible at higher transverse momenta, due to a more complicated production mechanism. This state is particularly interesting as it is connected to the mass of the gauge boson [25]. The measurements [28,29] by CMS and ATLAS experiments established that W bosons produced in pp collisions with large transverse momenta are predominantly left-handed, as expected in the Standard Model.
In the standard notation of the helicity fractions, the following relations with A i 's of Eq. (2) are valid The difference between left-and right-handed fraction is pro- Note, that even if Eq. (2) is valid for any definition of the W-boson rest frame, the A i ( p T , Y ) are frame dependent. The relations Eqs. (5) and (6), hold in the helicity frame.
Very similar arguments can be made also for the case of the Z production. However, the different characteristic of couplings have to be considered: the coupling of Z-boson to quarks does not involve the chirality projector 1 2 (1 − γ 5 ), but is asymmetric between left and right handedness. The analysing power of Z leptonic decays is severely affected by the coupling to right-handed leptons, being similar to the coupling to left handed leptons. As a consequence the angular coefficients f L , f R , f 0 can no longer be interpreted directly as polarisation fractions of the Z boson. The respective matrix transformation, involving left-and right couplings of Z boson to fermions, relates them to the Z -boson polarisation fractions [30].
For the case of Z → channel the measurement of the complete set of A i 's coefficients in the Collins-Soper frame was recently performed at 8 TeV pp collisions by the CMS Collaboration [31] and the ATLAS Collaboration [13]. The precision of the measurement by the ATLAS Collaboration allowed to clearly show that the violation of the Lam-Tung sum rule [32] i.e. A 0 = A 2 , is much stronger than predicted by NLO calculations. It has shown also an evidence of A 5 , A 6 , A 7 deviation from zero.
As of today, the situation with the measurement of A i coefficients for W → ν production in hadronic collisions is far from satisfactory. Measuring only some coefficients like α 2 in the Collins-Soper frame or f L , f R , f 0 in the helicityframe as function of W-boson transverse momenta does not give a complete picture on the QCD dynamics of the production process. In the early papers [16,33], the point was made, that measurement of the complete set of coefficients may not be possible, due to limitations related to the reconstruction of lepton neutrino momentum: in particular a two-fold ambiguity in the determination of the sign of cos θ defined in the next section.
In the present paper we argue, that following the strategy outlined in [13], one can design a measurement which allows to measure the complete set of coefficients also in the case of W → ν in pp collision. Then, we move to the discussion of the reference frames used for W → ν decay and demonstrate that the Mustraal [34] frame introduced and detailed for LHC in [35] will be interesting in the case of W → ν production as well.
Our paper is organized as follows. Section 2 is devoted to the presentation of the strategy which allows to measure complete set of the A i 's coefficients in case of W → ν process. We follow this strategy and show a proof of concept for such measurement. In Sect. 3, we discuss variants for the frames of the θ, φ angles definition. In Sect. 4 we collect numerical results for the A i 's coefficients in the case of pp → ν + 1 j generated with QCD LO MadGraph5_aMC@NLO Monte Carlo generator [36] and QCD NLO Powheg+MiNLO Monte Carlo generator [37,38]. We elaborate on possible choices of the coordinate frame orientation. We recall arguments for introducing the Mustraal frame [35], (where the orientation of axes is optimized thanks to matrix element and next-to-leading logarithm calculations) and compare the Collins-Soper and Mustraal frames. We demonstrate that, similarly to the Z → case discussed in [35], with the help of probabilistic choice of reference frames for each event, the results of formula (3.4) from [34] are reproduced and indeed only one non-zero coefficient in the decomposition of the angular distribution is needed. Finally, in Sect. 5 we conclude the paper.
To avoid proliferation of the figures, we generally present those for W − → − ν only, while the corresponding ones for W + → + ν are deferred to Appendixes A-C.

Angular coefficients in W → ν production
The production of vector bosons at LHC displays new characteristics compared to the production at Tevatron due to proton-proton nature of the collision: asymmetries in charge and momentum for vector bosons and their decay leptons. Large left-handed polarisation is expected in the transverse plane. Contrary to the case of pp collisions, the angular coefficients in pp collisions of the W + and W − are not related by CP transformation, due to absence of such symmetry in the proton structure functions. Only quarks can be valence, while both quarks and anti-quarks may be non-valence.
For the numerical results presented in this section we use a sample of 4M events pp → τ ν + 1 j generated at QCD LO with MadGraph5_aMC@NLO Monte Carlo [36], with minimal cuts on the generation level, i.e. p j T > 1 GeV, and default initialisation of other parameters. 1 The purpose of presented results is not so much to give theoretical predictions on A i 's but to illustrate the proof of concept for the proposed measurement strategy. Therefore we will not elaborate on the choices of PDF structure functions, QCD factorisation and normalisation scale, or EW scheme used. However, numerical results are sensitive to particular choices. The experimental precision of reconstruction of hadronic recoil necessary for the neutrino momentum reconstruction can be impressive as it is the case of W mass [39], but it is less accurate for other measurements where such precision was not essential.

Kinematic selection
Kinematic selections need to be applied in the experimental analysis. Limited coverage in the phase-space is due to the need for the efficient triggering, detection and background suppression. It inevitably reshapes angular distributions of the outgoing leptons. The minimal set of selections, in the context of LHC experiments, is to require (in the laboratory frame), the transverse momenta and pseudo-rapidity of charged lepton p T > 25 GeV and |η | < 2.5 respectively. Typically selection to suppress background from the multi-jet events, we require neutrino transverse momenta p ν T > 25 GeV and the transverse mass of the charge-lepton and neutrino system m T > 40 GeV. This set of selections will define fiducial phase-space of the measurement. Similar selection was used e.g. in measurement [40]. In Fig. 1 we show as an example the pseudo-rapidity distribution of the charged lepton from W ± → ± ν decay, in the full phase-space and in the fiducial phase-space as defined above. Clearly, the distributions are different between W + → + ν and W − → − ν processes. 1 In principle any other lepton flavour could have been used for presentation of numerical results. Our choice to generate τ ν final states is motivated by the planned extensions of the work.  For the leptonic decay mode, W -bosons have the disadvantage (with respect to the Z -bosons) that the decay kinematics cannot be completely reconstructed due to the unobservability of the outgoing neutrino. On the other hand, we can profit from a simplification: the electroweak interaction does not depend on the virtuality of the intermediate state. The transverse components of the neutrino's momentum p ν x , p ν y can be approximated from missing transverse momentum balancing the event. The longitudinal component p ν z can only be calculated up to a twofold ambiguity when solving the quadratic equation on the invariant mass of the lepton-neutrino system m W , assuming its value is known.
Let us recall the corresponding simple formulas: where Equation (7) has two solutions. Moreover solutions exist only if Δ = (b 2 − 4a · c) is positive. It requires also, that the mass of the W boson, m W , is fixed, usually m P DG W is taken (no smearing due to its width). The solution for the neutrino momentum allows to calculate its energy, completing the kinematics of massless neutrino Some studies of the past [41], investigated if a better option can be designed than taking one of the two p ν z solutions randomly, with equal probabilities. In particular in case that solutions do not exist, if replacing the m P DG W by e.g. transverse mass m W T can be beneficial. No convincing alternative was found. Replacing m P DG W with the transverse mass was creating spikes in shapes of angular distributions that are difficult to control. Similar effect (i.e. spiky distortions of the angular distributions), was caused by favoring some solutions of the neutrino momenta e.g. by selecting the one in the most populated regions of the multi-dimensional phase space, or taking the bigger of the two, etc.
In the analysis which will be outlined below, we propose to: -Use nominal PDG value for m W to solve the equation for the neutrino momenta p ν z .
-Choose with equal probabilities, one of the two solutions for the neutrino momenta p ν z . This solution will be called a random solution.
We estimated that the loss of events due to Δ < 0 is on the level of 10%. In the experimental analysis this loss can be considered as a part of other events losses due to kinematical selection cuts, like thresholds on the lepton transverse momenta or pseudo-rapidity bounds due to limited detector acceptances.

Collins-Soper rest-frame
For the Drell-Yan productions of the lepton-pair in hadronic collisions, the well known and broadly used Collins-Soper reference frame [18] is defined as a rest-frame of the leptonpair, with the polar and azimuthal angles constructed using proton directions in that frame. Since the intermediate resonance, the W -or Z -boson are produced with non-zero transverse momentum, the directions of initial protons are not collinear in the lepton-pair rest frame. The polar axis (zaxis) is defined in the lepton-pair rest-frame such that it is bisecting the angle between the momentum of one of the protons and the reverse of the momentum for the other one. The sign of the z-axis is defined by the sign of the lepton-pair momentum with respect to z-axis of the laboratory frame. To complete the coordinate system the y-axis is defined as the normal vector to the plane spanned by the two incoming proton momenta in the W rest frame and the x-axis is chosen to set a right-handed Cartesian coordinate system with the other two axes. Polar and azimuthal angles are calculated with respect to the outgoing lepton and are labeled θ and φ respectively. In the case of zero transverse momentum of the lepton-pair, the direction of the y-axis is arbitrary. Note, that there is an ambiguity in the definition of the φ angle in the Collins-Soper frame. The orientation of the x-axis here follows convention of [16,42,43].
For the Z → + − production, the formula for cos θ can be expressed directly in terms of the momenta of the outgoing leptons in the laboratory frame [44] cos with where E i and p z,i are respectively the energy and longitudinal momentum of the lepton (i = 1) and anti-lepton (i = 2) and p z ( + − ) denotes the longitudinal momentum of the lepton system, m( + − ) its invariant mass. The φ angle is calculated as an angle of the lepton in the plane of the x and y axes in the Collins-Soper frame. Only the four-momenta of outgoing leptons and incoming proton directions are used. That is why the frame is very convenient for experimental purposes.
In case of W ± → ± ν production we follow (in principle) the same definition of the frame. We use the convention that the θ and φ angles define the orientation of the charged lepton, i.e. anti-lepton of W + production and lepton in case of W − production. We calculate cos θ for the chosen solution of neutrino momenta with formula (9) and φ from the event kinematics as well.
However, because neutrino momentum remains unobserved, one of the two solutions to Eq. (7) has to be chosen and m P DG W must be used instead of ± ν pair mass. Figure 2 shows correlation plots between cos θ gen , φ gen calculated using generated neutrino momenta, and cos θ , φ calculated using neutrino momenta from formula (7)    The correlation plots of cos θ gen and φ gen calculated using generated neutrino momenta and cos θ and φ calculated using m W = m P DG W for solving Eq. (7). The plots for correct (two top) or wrong (two bottom) solution for neutrino momenta are shown. Correlation plots are prepared for the full phase-space of W − → τ − ν process neutrino momenta 2 are selected. The cos θ − cos θ gen and φ − φ gen , can be anti-correlated in case of correct and wrong solutions, the effect is much stronger for wrong solutions and the cos θ variable. We observe also inevitable migrations between bins due to the approximation m W = m P DG W used for solving Eq. (7). The following two figures will demonstrate how elements of reconstruction based on Eq. (7) affect the obtained angular distributions. Figure 3 shows cos θ and φ distributions of the charged lepton from W → ν decays in the Collins-Soper rest frame. We use the generated W boson mass m gen W of a given event or the fixed PDG value m P DG W for calculating neutrino momenta p ν z , taking the correct solution for p ν z . We compare the two results. The losses due to the non-existence of a solution of Eq. (7) are concentrated around cos θ = 0 but are uniformly distributed over the full φ range. Figure 4 demonstrates the variation of cos θ and φ distributions for charged lepton of W − → − ν decays when m W = m P DG W is used for solving Eq. (7) and the selection of the fiducial regions applied. In each case, distributions are shown for correct, wrong and random solution for p ν z . Selection of the fiducial region enhances modulation in the φ distribution. Corresponding distributions for W + → + ν decay are shown in Appendix A. Now we are ready to illustrate the effect of folding into fiducial phase-space, 2D distributions of (cos θ, φ) are shown in Fig. 5: (i) for events in the full phase-space, when generated neutrino momenta are used, (ii) in the fiducial phase-space when m W = m P DG W is used for Eq. (7) and random solution of neutrino momenta is taken. Clearly, original shapes of distributions are significantly distorted, but still, as we will see later, basic information on the angular correlation of the outgoing charged lepton and the beam direction is preserved. In particular, it is non trivial that the information is preserved despite approximate knowledge of the neutrino momentum. Moreover, the information is carried by both, correct and wrong, solutions for neutrino momenta. These observations are essential for the analysis presented in our paper.

Templated shapes and extracting A i 's coefficients
The standard experimental technique to extract parameters of complicated shapes is to perform the multi-dimensional fit to distributions of experimental data using either analytical functions or templated shapes [45]. Given what we observed in Fig. 4 only the second options seems feasible. The technique of templated shapes constructed from Monte 2 As correct we denote solution which is closer to the generated p ν z value, as wrong the other one. In practice, the two-fold ambiguity for solution of Eq. (9) is present. We will return to this point later, in Sect. 2.3. It is important for the observables we advocate in the paper. θ cos Carlo events, elaborated in [13] for the A i 's measurement in Z → case, is followed here and shortly described below.
We use the Monte Carlo sample of W ± → ± ν events and extract angular coefficients of Eq. (2) using moments methods [17]. The first moment of a polynomial P i (cos θ, φ), integrated over a specific range of p T , Y is defined as follows: Owing to the orthogonality of the spherical polynomials of Eq. (2), the weighted average of the angular distributions with respect to any specific polynomial, Eq. (10), isolates its corresponding coefficient, averaged over some phase-space region. As a consequence of Eq. (2) we obtain: We extract coefficients A i using generated neutrino momenta to calculate cos θ and φ. As a technical test, 2-dimensional histogram of (cos θ, φ) distribution obtained from our events weighted with where A 8 = 1.0 and P 8 = 1 + cos 2 θ is used. By construction, thanks to Eqs. (11) and (10), weighted with (12) sample, feature unchanged Y , p T distribution, but matrix element dependence of angular distribution of leptons in lepton pair rest-frame is completely removed.
If averages for (10) are taken for sub-samples in appropriately narrow bins of Y and p T this feature holds precisely for configurations of up single high p T , thus degrading predictions of the Monte Carlo simulation results, to at worst NLO (NLL) level. We have found that for numerical results binning in p T alone is sufficient. Indeed flat distribution in (cos θ, φ), where θ, φ are calculated using the generated neutrino momentum, see top plot of Fig. 6, is obtained. This completes our technical test and we can continue the construc-  Fig. 6 The 2D distribution of cos θ and φ of charged lepton from W − → τ − ν. On top, distribution of the full phase-space, with generated neutrino momentum used, and events weighted wt Σ Ai Pi . Bottom, the same distribution is shown, but: m P DG W is used for solving Eq. (7), randomly one of the solutions for p ν z is taken and fiducial selection is applied. The weight wt Σ Ai Pi is calculated with generated neutrino momenta, as it should be tion of templates. Further refinements are straightforward but require substantially more CPU time as binning in more than one or even more than two variables is then necessary.
We fold now events weighted with wt Σ Ai Pi into fiducial phase-space of the measurement: for the neutrino momentum reconstruction we use m W = m P DG W and take one of the solutions at random, then we recalculate θ, φ angles and finally we apply the kinematical selection of the fiducial phase-space. Bottom plot of Fig. 6 shows how the initially flat distribution is distorted by this folding procedure.
We can now model any desired analytical polynomial shape of the generated full phase-space folded into fiducial phase-space of experimental measurement. It is enough to apply wt i = P i · wt Σ Ai Pi to our events, to model the shape of the P i (cos θ, φ) polynomial in the measurement fiducial phase-space. In Fig. 7, we show 2D distributions modeling polynomials P 0 (cos θ, φ) and P 4 (cos θ, φ) in the full and fiducial phase-space as an example. Distributions for P 1 (cos θ, φ), P 2 (cos θ, φ) and P 3 (cos θ, φ) are shown in the Appendix C.
We can now proceed with the fit of a linear combination of templates to distributions of the fiducial phase-space pseudo- data. We bin in p W T both templates shown in Figs. 7, 21 and pseudo-data distributions shown in Fig. 5. We perform a multi-parameter log-likelihood fit 3 in each p W T bin; the only parameters of the fit are the 8 angular coefficients A i ( p W T ) (of a given p W T bin). Results of the fitting procedure are shown in Fig. 8. The black points represent fitted values of A i 's with their fit error, black open circles are the generated values of the A i 's (which we extracted with moments method 3 The Root framework [46] was used for this purpose. All eight A i ( p W T ) coefficients were fitted simultaneously. Statistical errors were calculated. Only for selected (non-zero), coefficients statistical errors (pools) are presented in plots, such as Fig. 8. For the remaining ones they were found to be of the same order. Non diagonal elements of correlation matrix [47] were all smaller than 0.7. described above). Bottom panels show difference between fitted and true values divided by their errors (so called pulls distributions). Pulls are small because of the samples correlations. We confirm closure of the method, i.e. extracted coefficients are equal to their nominal value for analysed events sample. However, estimation of the statistical errors for all A 0 to A 7 of the method, and in particular how it would compare with the case when ν momentum would be measured directly and full phase space was available for the measurement is premature. Corresponding discussion should include discussion of experimental systematic error as well.
The same procedure has been repeated for W + → + ν and results are shown in Appendix C.  In the bottom panels shown are pulls (difference between gen-erated and fit value, divided by the statistical error of the fit). Pulls are smaller than one could expect. This is because events of pseudo-data and templates are statistically correlated We have also performed the fit using templates and pseudo-data distributions prepared with only correct or only wrong solutions for the neutrino momenta. 4 In both cases the fit returned nominal values of the A i coefficients, confirming that both solutions of neutrino momenta carry the same information on the angular correlations.
In this proof of concept for the proposed measurement strategy we have not discussed possible experimental effects like resolution of the missing transverse energy which will be used to reconstruct neutrino transverse momenta or e.g. background subtraction which can limit precision of the measurement.
-We have shown, that despite the neutrino escaping detection, one can define, under some assumptions, the equivalent to the rest frame of lepton-neutrino system and preserve sensitivity for the complete set of angular coefficients of the decomposition. -Then, we have shown with a simplified version of the method used in [13] to measure A i coefficients in case of Z → decays, applied in case of W → ν decays allows for the measurement of a complete set of angular coefficients in this case as well.

Comments on systematic errors
Once we have collected numerical results, let us address the issues of systematic errors. The only difference of principle for our measurements, between W and Z decays, is the neutrino momentum reconstruction, which is by far less precise than the measurement of the second lepton of Z .
Let us start with the consequences of the use of m P DG W (instead of m W the invariant mass of the lepton-neutrino system) discussed in Sect. 2.2. Substantial complications arise: the two-fold ambiguity for the solution appears and for some other events solution for the neutrino momentum does not exist at all. On the other hand we have found the numerical impact of this replacement for the measurement of A i coefficients to be rather small. This is because m P DG W − m W is typically of the order of W width, thus an order of magnitude smaller than momenta of its decay products. Also the use of wrong instead of correct solution for neutrino momenta is not introducing bias in the result for A i coefficients. Detector effects on the reconstruction are much larger and must be taken into account for the experimental measurement.
We do not attempt to estimate expected experimental error, discussion of Ref. [48], taught us how unrealistic it can be, 4 With experimental data one can use random solution for p ν z only. The other two possibilities correct and wrong are for tests and to better understand two-fold ambiguity of neutrino reconstruction, which could be expected to be the reason for some of the A i coefficients to remain non-observable. despite the best efforts, once it is confronted with experimental analysis [49]. It can be simply misleading.
In Footnote 3, we have pointed that the multi-parameter fit provides numerical answer for all angular coefficients A 0 −A 7 simultaneously. It is justified to ask if the 8 parameter fit will provide constraints for all A 0 −A 7 . We have verified that it is the case; all non-diagonal elements of correlation matrix were smaller than 0.5. The evaluation of 8 dimensional correlation matrix for the fit parameters was delegated to specialized algorithms present in the Root framework [46]. This is an essential cross check, specialized algorithms for such purposes were developed already long time ago [47]. This encouraging result needs to be confirmed once all experimental effects are introduced. If this holds, the fit will not feature unacceptably strong correlations.

Angular coefficients and reference frames
In this paper, we concentrate on the numerical analysis of tree level parton-parton collisions into a lepton pair and accompanying jets, convoluted with parton distributions, but without parton showers. Even though such approach is limited, it provides input for general discussions. Such configurations constitute parts of the higher order corrections, or can be seen as the lowest order terms but for observables of tagged high p T jets.
For the choices of the reference frames to be discussed here, let us point out that in the limit of zero transverse momenta, all coefficients except A 4 vanish. The ( p T , Y ) dependence of the A i coefficients differs with the choice of the reference frame.
So far, we have introduced and discussed angular coefficients in the Collins-Soper frame only. Let us present now the variant of the reference frame definition we are also going to use, i.e. the Mustraal frame.

The Mustraal reference frame
The Mustraal reference frame is also defined as a rest frame of the lepton pair. It has been proposed and used for the first time in the Mustraal Monte Carlo program [34] for the parametrization of the phase space for muon pair production at LEP. The resulting optimal frame was minimising higher order corrections from initial state radiation to the e + e − → Z /γ * → ff and was used very successfully for the algorithms implementing genuine weak effects in the LEP era Monte Carlo program KORALZ [50]. A slightly different variant was successfully used in the Photos Monte Carlo program [51] for simulating QED radiation in decays of particles and resonances. The parametrization was useful not only for compact representation of single photon emissions but for multi-emission configurations as well.
Recently in [35], the implementation of the Mustraal frame has been extended to the case of pp collisions and studied for configurations with one or two partons in the final state accompanying Drell-Yan production of the lepton pairs. The details of the implementation of this phase space parametrization have been discussed in context of Z → events and the complete algorithm how to calculate cos θ and φ angles was given. There is no need to repeat it here.
Let us point out that unlike the case of the Collins-Soper frame, the Mustraal frame requires not only information on 4-momenta of outgoing leptons but also on outgoing jets (partons). The information on jets (partons), is used to approximate the directions and energies of incoming partons. This is used for the calculation of weights (probabilities) with which each event contributes to one of two possible Born-like configurations. Each configuration requires different cos θ , φ definition. This does not have to be very precise but can introduce additional experimental systematics, and requires attention. No dependence on coupling constants or PDF's is introduced in this way.

QCD and EW structure of angular correlations
The measurement of the angular distribution of leptons from the decay of a gauge boson V → where V = W, Z or γ * , produced in hadronic collisions via a Drell-Yan-type process h 1 + h 2 → V + X provides a detailed test of the production mechanism, revealing its QCD and EW structure.
The predictive power of QCD is based on the factorisation theorem [52]. It provides a framework for separating out long-distance effects in hadronic collisions. Consequently, it allows for a systematic prescriptions and provides tools to calculate the short-distance dynamics perturbatively, at the same time allowing for the identification of the leading nonperturbative long-distance effects which can be extracted from experimental measurements or from numerical calculations of Lattice QCD.
The question of the input from the Electroweak sector of the Standard Model is important, especially for distributions of leptons originating from the intermediate Z /γ * state. We have addressed numerical consequences of this point recently in [53] in the context of τ lepton polarization in Drell-Yan processes at the LHC. A wealth of publications was devoted during last years to this issue, see e.g. [54][55][56]. We should underline limitations of separating interactions into Electroweak and QCD part. Limitations have been well known for more than 15 years now, see e.g. [57].
Let us come back now to Eqs. (1) and (2) and discuss the structure of cross-section decomposition into harmonic polynomials multiplied by angular coefficients. The A i coefficients represent ratios of the helicity cross-sections and following the conventions and notations of [16,17], the fol-lowing coefficients constructed from couplings, appear in A i 's: .
where v i , a i , i = q, denote vector and axial couplings of intermediate boson to quarks and leptons.
In case of W boson the EW sector at leading order is simply a (V − A) coupling only. At higher order and higher p W T the more complicated structure (potentially also more interesting) of the multi-boson couplings, if present, may be revealed. In case of the Z /γ * , the sensitivity to the EW sector is much richer from the physics point of view, in particular for A 3 and A 4 coefficients, and we have discussed it recently in [35,53].

Numerical results for Collins-Soper and Mustraal frames
Let us now present numerical results for the angular coefficients A i and compare predictions in the Collins-Soper and Mustraal frame for W − production. Most of results for W + are delegated to Appendix C. We will present numerical results for coefficients A 0 to A 4 only, the other ones are consistent with zero over all kinematical range, for both Collins-Soper and Mustraal frames with similar statistical errors. One can conclude that for both choices of frames all A i coefficients are constrained. Thus, in principle interpretation of results in language of different helicity amplitudes of the scattering process is possible.

Results with LO simulation
We use samples of events generated with the MadGraph-5_aMC@NLO Monte Carlo [36] for Drell-Yan production of W + 1 j with W → τ ν and 13 TeV pp collisions. Lowest order spin amplitudes are used in this program for the parton level process. To better populate higher p W T bins we merged (adjusting properly for relative normalization) 3 samples, 2 × 10 6 events each, generated with thresholds of p j T > 1, 50, 100 GeV respectively. The incoming partons distributed accordingly to PDFs (using CTEQ6L1 PDFs [58] linked through LHAPDF v6 interface) remain precisely  Fig. 9 The A 0 and A 4 coefficients calculated in Collins-Soper (black) and in Mustraal (red) frames for pp → τ − ν + 1 j process generated with MadGraph. Further plots, but for A 1 , A 2 , A 3 and A 0 − A 2 are delegated to Fig. 16 of Appendix B collinear to the beams. At this level, jet (j) denotes outgoing parton of unspecified flavour. Figure 9, collects results for angular coefficients A i of the processes with W − → τ − ν in the final state. We show sets of five angular coefficients A 0 − A 4 only; the remaining ones A 5 − A 7 are close to zero over the full p W T range for both definition of frames; Collins-Soper and Mustraal. For the A 3 and A 4 coefficients we show their absolute value, as in the convention we have adopted the signs depend on the sign of the charge of the W boson (so the charge of lepton). In case of W + production, both A 3 and A 4 are negative, see Appendix C.
In both frames, at low p W T the only non-zero coefficient is A 4 , and is of the same value. Similarly as we observed in case of Drell-Yan Z → process [35], the only significantly non-zero coefficient in the Mustraal frame at higher p W T remains A 4 , while A 0 , A 2 and A 3 are rising steeply for higher p W T in the Collins-Soper frame.  Fig. 10 The A 0 and A 4 coefficients calculated in Collins-Soper frame for pp → τ ± ν + 1 j processes generated with MadGraph. Further plots, but for A 1 , A 2 , A 3 and A 0 − A 2 are delegated to Fig. 17 of Appendix B Figures 10 and 11 show the comparison of predicted coefficients for W + and W − , respectively in Collins-Soper and Mustraal frames. The noticeable difference of the A 4 coefficients at low p W T directly reflects different compositions of the structure functions in pp collisions to produce W + and W − which enters the average over couplings shown in Eq. (13). This difference is present for both the Collins-Soper and the Mustraal frame. For the Collins-Soper frame we observe also sizable A 3 coefficient above p W T = 100 GeV. As stated already in [16], the theoretical uncertainties due to the choice of the factorisation and renormalization scales are very small for A i representing the cross section ratios. Also most of the uncertainties from the choice of structure functions and factorisation scheme cancel in the ratios.

Results with NLO simulation
So far, we have discussed results for samples of fixed order tree level matrix elements and of single parton (jet) emis-  Fig. 11 The A 0 and A 4 coefficients calculated in Mustraal frame for pp → τ ± ν + 1 j processes generated with MadGraph. Further plots, but for A 1 , A 2 , A 3 and A 0 − A 2 are delegated to Fig. 18 of Appendix B sion. In general, configurations with a variable number of jets and effects of loop corrections and parton shower of initial state should be used to complete our studies. We have performed this task partially only, with the help of 10M weighted W + + j and W − + j events, with W ± → τ ± ν generated with Powheg+MiNLO Monte Carlo, again for pp collisions at 13 TeV and the effective EW scheme. The PowhegBox v2 generator [37,38], augmented with MiNLO method for choices of scales [59] and inclusion of Sudakov form factors [60], by construction achieves NLO accuracy for distributions involving finite non-zero transverse momenta of the lepton system. Two jet configurations are thus present.
In Fig. 12 results for A i 's coefficients for W − → τ − ν are shown, extracted using moments method [17] described in Sect. 2.4. Comparisons of results using Mustraal and the Collins-Soper frames feature again the pattern, observed already for QCD LO W + 1 j events generated with MadGraph. As predicted for Z → and higher p Z T , the Lam-Tung relation A 0 = A 2 [32] is again violated at

Summary
The interest in the decomposition of results for measurement of final states in Drell-Yan processes at the LHC, into coefficients of second order spherical polynomials (for angular distributions of leptons in the lepton-pair rest-frame), was recently confirmed by experimental publications for Z → process [13,31,61] and W → ν processes [28,29]. Inspired by those measurements, we have investigated the possibility to apply a strategy similar to [13] for W → ν production and decay. Thus, to contest the statements which are often made in the literature, that the neutrino escaping detection makes measurement of the complete set of angular coefficients not possible. We have shown a proof of concept for the proposed strategy, for the measurement of the complete set of A i . To establish precision level of 1% (required for the W mass measurement at 10 MeV precision level), require experimental details to be included, it is thus out of scope of the present paper. As an example, Monte Carlo events of simulated W → ν +1 j process were analysed and the complete set of A i coefficients was extracted from a fit to pseudo-data distributions in the fiducial phase-space, in agreement with the prediction for this sample calculated with the moments method. The results were cross-checked to hold when QCD NLO effects were taken into account.
In the second part of this paper, we have discussed the optimal reference frame for such measurement. Two frames: Collins-Soper and Mustraal [35] were studied. We have presented predictions for the angular coefficients in those frames as function of W -boson transverse momenta. We have shown that as expected, in case of the Mustraal frame, only one coefficient remains significantly non-zero and constant, almost up to p W T = 100 GeV, where it starts decreasing. Similarly as we argued in [35], this may help to facilitate the interpretation of experimental results into quantities sensitive to strong interaction effects. The longitudinal W boson polarisation seems to appear predominantly as kinematic consequence of the choice of reference frames even in configurations of high p T jets.  Figure 14 shows variation of cos θ and φ distributions for charged lepton from W + decays when m W = m gen W or m W = m P DG W are used for solving Eq. (7). In the second case selection of the fiducial region is applied. In each case distributions are shown for correct, wrong and random solution for p ν z . One can notice, comparing with Fig. 4, that the shapes of the cos θ distributions in case of W + and W − are not mirrored when random or wrong solution for the neutrino momenta are used. Figure 15 shows 2D distributions in (cos θ, φ) for events in: full phase-space with generated neutrino momenta used and in the fiducial phase-space (also m W = m P DG W is used for solving Eq. (7) with random solution of neutrino momenta taken).

Appendix B: Further plots for templated shapes and extracting A i 's coefficients
Let us start the Appendix B with plots for A 1 , A 2 , A 3 and A 0 − A 2 (Figs. 16,17,18,19) supplementing the plots for A 0 and A 4 presented in Sect. 4.  Fig. 16 The A 1 , A 2 , A 3 and A 0 − A 2 coefficients calculated in Collins-Soper (black) and in Mustraal (red) frames for pp → τ − ν + 1 j process generated with MadGraph. Plots for A 0 and A 4 are given in Fig. 9 (GeV)  Fig. 17 The A 1 , A 2 , A 3 and A 0 − A 2 coefficients calculated in Collins-Soper frame for pp → τ ± ν +1 j processes generated with MadGraph. Plots for A 0 and A 4 are given in Fig. 10 (GeV)  Fig. 18 The A 1 , A 2 , A 3 and A 0 − A 2 coefficients calculated in Mustraal frame for pp → τ ± ν + 1 j processes generated with MadGraph. Plots for A 0 and A 4 are given in Fig. 11 Fig. 19 The A 1 , A 2 , A 3 and A 0 − A 2 coefficients calculated in Collins-Soper (black) and in Mustraal (red) frames for pp → τ − ν + 1 j process generated with Powheg+MiNLO. Plots for A 0 and A 4 are given in Fig. 12 Figure 20 (bottom) shows (cos θ, φ) distribution of W + → + ν events, where θ, φ are calculated using generated neutrino momentum and events are weighted with wt Σ Ai Pi . Bottom plot of Fig. 20 show how the initially flat distributions are distorted by this folding procedure. Note that the shape is not the same as in case of W − → − ν shown in Fig. 6. This is due to different distributions e.g. the pseudo-rapidity of charged lepton and thus different acceptance of the fiducial selection. Figure 21 collects examples of 2D distributions for polynomials P 1 (cos θ, φ), P 2 (cos θ, φ) and, P 3 (cos θ, φ) in the full and fiducial phase-space. Now for each event we reconstruct neutrino momenta using m W = m P DG W , take randomly one of the solutions to recalculate θ, φ angles from Eq. (7) and to apply kinematical selection of the fiducial phase-space.  Fig. 6 but for W + → τ + ν events. The 2D distribution of charged lepton cos θ and φ. On top distribution of the full phasespace, with generated neutrino momentum used, and events weighted wt Σ Ai Pi . Bottom, the same distribution is shown, but: m P DG W is used for solving Eq. (7), randomly one of the solutions for p ν z is taken and fiducial selection is applied. The weight wt Σ Ai Pi is still calculated with generated neutrino momenta In the bottom panels shown are pulls (difference between generated and fit value, divided by the statistical error of the fit). Pulls are smaller than one could expect. This is because events of pseudo-data and templates are statistically correlated

Appendix C: Additional plots on A i s coefficients
In the generated sample information on incoming and outgoing partons flavours is stored. We will use this information for tests, to define sub-samples of qq and q(q)G parton level initial states. Figures 23 and 24 show predictions for A i 's coefficients for W + → τ + ν and events generated with MadGraph Monte Carlo for these sub-samples. Figure 25 shows predictions for A i 's coefficients of W + → + ν and for processes generated with Powheg+MiNLO.