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 an important role, among the physics goals of LHC experiments. Because of the nature of proton-proton processes, observables based on the measurement of the direction and energy of leptons provide the most precise signatures. In the present paper, we concentrate on the angular distribution of leptons from W to l nu 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 polynomials 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 in the final state of neutrino escaping detection. There is thus no principle difference with respect to the phenomenology of the Z/gamma to l^+ l^- Drell-Yan process. We show also, that with the proper choice of the coordinate frames, only one coefficient in this polynomial decomposition remains sizable, even in the presence of one or more high p_T 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 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 eg. [6,7,8], a program of measurements in the domain of Electroweak (EW) and Strong (QCD) interactions is on-going. This is the keystone for establishing the Standard Model as a fundamental 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 intermediate Z and W bosons represent the primary group of measurements of the second 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 the discovery of the W boson, its hadronic production both in pp and pp collisions, mass and the width, have been measured to great precision [14]. To complete physical information on the production process, measurements pursued the boson's differential distribution. The measurements rely on outgoing leptons of the W → ℓν decays in the W -boson rest frame. Because EW 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 inprinted in the angular distributions of outgoing leptons.
In principle, the same standard formalism of the Drell-Yan production Z → ℓℓ [15] can be applied in case of W → ℓν production [16,17]. The angular dependence of the differential cross-section can be written again as 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 lepton in the W rest-frame. The p T , Y denote transverse momenta and rapidity of the intermediate 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σ d p 2 T dY d cosθdφ dσ U+L d p 2 T dY (2) [(1 + cos 2 θ) + 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 cross-section (a convention used in several papers of the 80's).
In case of W boson, (θ, φ) define the orientation of the charged lepton from W → ℓν . 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 it 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 leptonneutrino 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 θ) with the following relations between coefficients; It has been estimated that 1% uncertainty on α 2 corresponds to a shift of the measured m W in pp collision, determined by fitting the transverse mass distribution, of approximately 10 MeV. The α 1 measures the forwardbackward leptonic 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 transvers momentum. The measurements confirmed SM expectations, that α 2 decreases with increasing W boson transvers 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 two-fold 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 sea-quarks and gluon contributions to the structure functions can be neglected. Such 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 tests 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 displays therefore 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 in principle particularly interesting as it is connected to the massive character of the gauge boson [25]. The measurements [28,29] by ATLAS and CMS 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 proportional to A 4 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 Eq. (5) and (6), hold in the helicity frame.
2 Very similar arguments can be made also for the case of 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 asymmetric between left and right handedness. Contrary, 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 being not equal to 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 helicity-frame as function of W-boson transverse momenta is not giving a complete picture on the QCD dynamics of the production process. Already in the first papers [16,33] the point was made, that measurement of the complete set of coefficients is not possible, due to limitations related to the reconstruction of lepton neutrino momenta, leading to a two-fold ambiguity in the determination of the sign of cos θ.
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 Section 3, we discuss variants for the frames of the θ, φ angles definition. In Section 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 Section 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 compatred 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 therefore on the choices of PDF structure functions, QCD factorisation and normalisation scale, or EW scheme used. However, numerical results are sensitive to particular choices.

Kinematical selection
The kinematical selection needs to be applied in the experimental analysis. The limited coverage in the phasespace is needed for the efficient triggering, detection and background suppression. It inevitably reshapes angular distributions of the outgoing leptons. The minimal set of selection, in the context of LHC experiments is to require that in the laboratory frame, the transverse momenta of charged lepton p ℓ T > 25 GeV and pseudorapidity |η ℓ | < 2.5. As the typical 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 selection will define fiducial phase-space of the measurement. Similar selection was used e.g. in measurement [39]. In Fig. 1 we show as an example the pseudorapidity 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.

Solving equation for neutrino momenta
For the leptonic decay mode, W -bosons have the disadvantage with respect to the Z-bosons because 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 Eq. (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 PDG 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 [40], 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 PDG W by e.g. transverse mass m W T can be beneficial. No convincing alternative was found. Replacing m PDG 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 randomly, 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% and in the experimental analysis can be considered as a part of other events losses due to kinematical selection cuts, like thresholds on the lepton transverse momenta or pseudorapidity 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 lepton-pair, 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 leptonpair rest frame. The polar axis (z-axis) is defined in the lepton-pair rest-frame such that it is bisecting the angle between the momentum of one of the proton and inverse of the momentum of the second 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,41,42].
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 [43] 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 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. Figure 2 shows correlation plots between cos θ gen , φ gen calculated using generated neutrino momenta, and cos θ, φ calculated using neutrino momenta from formula (7)    Correlations are shown when correct or wrong solution for 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 PDG W used for solving Eq. (7). 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 PDG W for calculating neutrino momenta p ν z , taking the correct solution for p ν z . We compare the two results. The losses due to the nonexistence 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 PDG 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.
To illustrate the effect of folding into fiducial phase-space, 2D distributions of (cos θ, φ) are shown in Figure 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 PDG 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. Given what we ob-θ cos served in Figure 4 only the second options seems feasible. The technique of templated shapes constructed from Monte 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, we histogram 2D distribution in (cos θ, φ) using our events weighted with where A 8 = 1.0 and P 8 = 1 + cos 2 θ. We obtain a completely flat distribution in (cos θ, φ), where θ, φ are calculated using the generated neutrino momentum, see left plot of Figure 6. This completes our technical test and we can continue the construction of templates.
We fold now events weighted with wt ΣAiPi into fiducial phase-space of the measurement: for the neutrino momentum reconstruction we use m W = m PDG 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. Right plot of Figure 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 ΣAiPi to our events, to model the shape of the P i (cos θ, φ) polynomial in the measurement fiducial phase-space. In Figure 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 the remaining ones 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 Fig. 7, 17 and pseudo-data distributions shown in Fig. 5. We perform a multi-parameter log-likelihood fit in each p W T bin; parameters of the fit are the angular coefficients A i (p W T ). Results of the fitting procedure are shown in Figure 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 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. The same procedure has been repeated for W + → ℓ + ν and results are shown in Appendix C.
We have also performed the fit using templates and pseudo-data distributions prepared with only correct or only wrong solutions for the neutrino momenta 3 . In both cases the fit returned nominal values of the A i coefficients, just 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 that simplified version of the method used in [13] to measure A i coefficients in case of Z → ℓℓ decays, can be applied in case of W → ℓν decays and allows for the measurement of a complete set of angular coefficients in this case as well.

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 [44]. A slightly different variant was successfully used in the Photos Monte Carlo program [45] 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  Figure 6: The 2D distribution of cos θ and φ of charged lepton from W − → τ − ν. On left side distribution of the full phase-space, with generated neutrino momentum used, and events weighted wt ΣAiPi . On right, the same distribution is shown, but: m PDG 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 ΣAiPi is calculated with generated neutrino momenta, as it should be.  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 pseudodata and templates are statistically correlated.
used to approximate the directions and energies of incoming partons 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 [46]. It provides a framework for separating out long-distance effects in hadronic collisions. In consequence, 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 [47] 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. [48,49,50]. We should underline limitations of separating interactions into Electroweak and QCD part. Limitations are well known, since more than 15 years now, see e.g. [51].
Let us come back now to Eq. (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 following 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, and of more interesting nature of the multi-boson couplings, if such is 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,47].

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.

Results with LO simulation
We use samples of events generated with the MadGraph5_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) 12 3 samples, 2M events each, generated with thresholds of p j T > 1, 50, 100 GeV respectively. The incoming partons distributed accordingly to PDFs (using CTEQ6L1 PDFs [52] linked through LHAPDF v6 interface) remain precisely 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. 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) emission. 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 [53] and inclusion of Sudakov form factors [54], 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 Section 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 QCD NLO in the Collins-Soper frame. This confirms the robustness of our conclusions that only the Mustraal frame retains Born like A i coefficients for high p W T configurations.

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 restframe was recently confirmed by experimental publications for Z → ℓℓ process [55,31,13] 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 . 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     Figure 11: The A i coefficients calculated in Mustraal frame for pp → τ ± ν + 1 j processes generated with MadGraph. 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. θ cos

A Distributions for W + → ℓν
In this Appendix, we collect plots corresponding to the ones shown in Section 2 but for W + → ℓ + ν. Figure 13 shows cos θ and φ distributions of charged lepton in the Collins-Soper rest frame. We use generated W boson mass in a given event m W = m gen W and compare to the case when fixed PDG value m W = m PDG W for calculating neutrino momenta p ν z , taking correct solution for p ν z is used. The losses due to nonexisting solution of Eq. (7) are concentrated at cos θ = 0 but are uniformly distributed over the whole φ range. Figure 14 shows variation of cos θ and φ distributions for charged lepton from W + decays when m W = m gen W or m W = m PDG 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 Figure 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 PDG W is used for solving Eq. (7) with random solution of neutrino momenta taken).

B Further plots for templated shapes and extracting A i 's coefficients.
Figure 16 (right) shows (cos θ, φ) distribution of W + → ℓ + ν events, where θ, φ are calculated using generated neutrino momentum and events are weighted with wt ΣAiPi . Right plot of Figure 16 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 Figure 6. This is due to different distributions e.g. the pseudorapidity of charged lepton and thus different acceptance of the fiducial selection. Figure 17 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 PDG W , take randomly one of the solutions to recalculate θ, φ angles from Eq. (7) and to apply kinematical selection of the fiducial phase-space. Figure 18 collects results of the multi-likelihood fit of W + → ℓ + ν, displayed are A i coefficients as function of p W T .

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 19 and 20 show predictions for A i 's coefficients for W + → τ + ν and events generated with MadGraph Monte Carlo for these sub-samples. Figure 21 shows predictions for A i 's coefficients of W + → ℓ + ν and for processes generated with Powheg+MiNLO.  Figure 16: As Fig. 6 but for W + → τ + ν events. The 2D distribution of charged lepton cos θ and φ. On left side distribution of the full phase-space, with generated neutrino momentum used, and events weighted wt ΣAiPi . On right, the same distribution is shown, but: m PDG 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 ΣAiPi 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.  Fig. 12 but for W + → ℓ + ν events. The A i coefficients calculated in Collins-Soper (black) and in Mustraal (red) frames for pp → τ + ν + 1 j process generated with Powheg+MiNLO.