Automated Predictions from Polarized Matrix Elements

The anticipated experimental resolution and data cache of the High Luminosity Large Hadron Collider will enable precision investigations of polarization in multiboson processes. This includes, for the first time, vector boson scattering. To facilitate such studies, we report the automation of polarized matrix element computations in the publicly available Monte Carlo tool suite, MadGraph5_aMC@NLO. This enables scattering and decay simulations involving helicity-polarized asymptotic or intermediate states, preserving both spin-correlation and off-shell effects. As demonstrations of the method, we investigate the leading order production and decay of polarized weak gauge bosons in the process $pp \to j j W^+_\lambda W^-_{\lambda'}$, with helicity eigenstates $(\lambda,\lambda')$ defined in various reference frames. We consider the Standard Model at both $\mathcal{O}(\alpha^4)$ and $\mathcal{O}(\alpha^2 \alpha_s^2)$ as well as a benchmark composite Higgs scenario. We report good agreement with polarization studies based on the On-Shell Projection (OSP) technique. Future capabilities are discussed.

An impeding factor to fully utilizing these data, however, are the limited number fully differential, SM and beyond the SM (BSM) predictions for polarization observables. While incredible efforts are underway to develop precise predictions up to next-to-next-to-leading order (NNLO) in QCD and/or next-to-leading order (NLO) in EW, these are largely restricted to only a handful of SM processes [50,51]. Likewise, the direct simulation of polarized parton scattering in hadron collisions using public Monte Carlo (MC) tools is found almost exclusively at leading order (LO) and again restricted to certain scattering topologies [52]. Though the availability of such public tools has led to a number of complimentary investigations on the production of polarized EW bosons at the LHC [52][53][54][55][56][57][58].
In the present work, we report the development of a scheme to model polarized parton scattering in hadron collisions and its implementation into the publicly available * event generator MadGraph5_aMC@NLO (dubbed mg5amc). By "polarized parton scattering" we specifically mean 2 → n-scattering and 1 → n-decay processes that are determined from polarized matrix elements (MEs). That is to say, matrix elements where some or all external states are in a definite helicity eigenstates and where spin-averaging or spin-summing is truncated or not present. (For simplicity, we refer to all short-distance particles, including massive, colorless EW states, as partons throughout this work.) The method enables the LO simulation of tree-level scattering and decay processes involving external states in fixed helicity eigenstates in an arbitrary reference frame. This includes massless QCD partons, heavy quarks, all leptons, the EW gauge sector, and states up to spin 3/2 and 2. When using the narrow width approximation (NWA), spin correlations of decaying polarized resonances are maintained through the decomposition of fermionic and bosonic propagators into their respective transverse, longitudinal, and auxiliary ("scalar") components; the last of which necessarily vanishes in the on-shell limit. Extension to new physics scenarios is achieved when used with Universal FeynRules Object (UFO) [59][60][61] libraries.
Our work continues as follows: In section (2) we summarize our computational setup. In section (3) we describe our formalism for constructing polarized MEs and its implementation into the program mg5amc. (More technical, implementation details and checks are reported in appendix (A).) We then investigate in section (4) the production and decay of polarized W + W − pairs from EW (section (4.2)) and mixed EW-QCD (section (4.3)) processes pp → jjW + W − in the SM as well as a benchmark CH scenario (section (4.1)). There we compare our methods to the so-called on-shell projection (OSP) technique [62][63][64] and report good agreement with past studies [52][53][54][55][56]. In Sec. 5 we summarize our results.
Throughout this study, we focus on EW and QCD processes at LO. We also report that fully differential event simulation up to NLO in QCD with parton shower matching is also possible for processes involving polarized, color-singlet final states. However, we report such investigations in a companion paper [65].

Computational Setup
In this section, we briefly summarize the computational framework used in our study. In particular, we describe the components of our MC tool chain and their relevant tunings needed for reproducibility. Standard Model and CH model input parameters used in our case studies are also listed. While we heavily utilize and partially expand on the MC suite mg5amc, a full characterization of it is outside the scope of this work and is available elsewhere [66,67]. The description of our method for constructing polarized MEs in hadron collisions and its implementation are given in section (3) and appendix (A).

Monte Carlo Framework and Tuning
We simulate parton scattering with polarized and unpolarized MEs in √ s = 13 TeV pp collisions at LO in perturbation theory using the software suite mg5amc [66]. Working in the so-called HELAS basis [68], tree-level matrix elements are evaluated numerically using helicity amplitudes that are capable of handling massive states [69,70], and with QCD color decomposition based on color flow [71]. Decays of unstable, resonant states are handled using the spin-correlated NWA, as implemented in MadSpin [72].

Standard Model Inputs
For SM inputs, we work in the n f = 4 massless quark scheme, approximate the Cabbibo-Kobayashi-Masakawa (CKM) matrix to be diagonal with unit entries, and take (2.1) In mg5amc, this corresponds to importing the internal sm model library. For SM case studies, we use the NNPDF 2.3 LO parton distribution function (PDF) set with α s (µ) = 0.119 (lhaid=246800) [73]. We set our collinear factorization (µ f ) scale to the MadGraph5 default.

Composite Higgs Inputs
Besides the SM, we also investigate VBS in the context of a CH scenario. For that we use the Higgs Characterization Model libraries of Ref. [74], which provides a general parametrization of the Higgs boson's spin and couplings. We limit ourselves to a somewhat generic CH situation, where the couplings of a SM-like Higgs are rescaled by an overall factor. In Ref. [74], this factor is identified as kSM with kSM=1 denoting the SM limit. Throughout this work, we fix the SM (or SM-like) Higgs mass to m H = 125 GeV. All the other parameters are set as described above with the exception of using the NNPDF 3.1 NLO+LUXqed (lhaid=324900) PDF set [75], with PDF running handled using LHAPDF v6.1.6 [76].

Parton and Hadron Scattering with Polarized Matrix Elements
We now describe the scattering formalism underlying our implementation of polarized parton scattering into the event generator mg5amc. We start in section (3.1) with building meaningful definitions of parton † scattering with polarized MEs in unpolarized hadron collisions, noting instances of reference frame-dependence that are not usually present in standard MC computations. In section (3.2), we describe our treatment of decaying polarized resonances. Additional details related to technical implementation and usage are reported in appendix (A). Physics demonstrations are deferred to section (4).

Scattering Formalism with Polarized Matrix Elements
Preliminaries: Scattering with Unpolarized Matrix Elements In unpolarized proton collisions, the scattering observableÕ built from the n-body final state B at momentum transfers ( Q 2 ) much larger than the nonperturbative QCD scale (Λ NP ) is governed by the Collinear Factorization Theorem [77][78][79][80][81], For protons m = 1, 2, with 4-momenta P m = ( √ s/2)(1, 0, 0, ±1), the above stipulates that inclusive, hadron-level observables (dσ/dÕ) that are functions of external momenta, i.e.,Õ = g(p 1 , . . . , p n ), can be expressed as the product of probabilities (convolution) for (a) finding partons i and j ∈ {q, q, g, γ}, with q ∈ {u, d, c, s}, in proton m, which is described by f ; (b) the renormalization group (RG) scale evolution of i and j from a proton to the hard scattering process, described by ∆; and (c) the exclusive, parton-level hard scattering process ij → B, governed by dσ/dÕ. Here, τ = Q 2 /s is the hard threshold at which ij → B proceeds, and for τ < τ 0 = min{Q 2 }/s, the production of B is kinematically forbidden.
More specifically, f k/p (ξ m , µ f ) is the collinear PDF, which for momentum fraction 0 < ξ m < 1, represents the likelihood of parton k in proton m possessing a momentum p k = (ξ m E m , 0, 0, ±ξ m E m ). Using the DGLAP evolution equations [82][83][84], f can be RG-evolved to the collinear cutoff / factorization scale µ f . This accounts for (resums) an arbitrary number of initial-state emissions that are produced in association with k and carry a relative transverse momentum p T < µ f . Factors of (1 ↔ 2) and (1 + δ ij ) account for identical beam and identical initial parton symmetrization.
The Sudakov factor ∆ ij (z, µ f , µ r , µ s ) accounts for (resums) soft and/or collinear emissions of massless partons carrying a momentum fraction z = Q 2 /ξ 1 ξ 2 s, away from the (ij) system prior to the hard ij → B scattering process. Through various RG evolutions between µ f , the UV renormalization scale µ r , and the Sudakov cutoff/ factorization scale µ s , ∆ ensures that eq. (3.1) remains RG scale-independent [85]. In our notation, ∆ additionally records i → k parton depletion and k → i parton buildup for hard scattering partons i and j. In general-purpose, MC event generators, ∆ ij (z) ≈ δ ij δ(1 − z) + O(α s ) can be identified as the parton shower and principally models collinear parton emissions, though developments to expand this domain are ongoing [86][87][88].
When built from the n-body final state B, parton-level scattering observables dσ/dÕ are derivable from the fully differentiated scattering rate dσ/dP S n , where dP S n is the separately Lorentz-invariant, n-body phase space measure given by 3) can be expressed in terms of perturbative matrix elements by the usual expression: Here s k = 1/2 and N k c are the helicity and SU(3) c color symmetrization factors for massless parton k = i, j. For massive spin-1 states, s k = 1 and the 2Q 2 flux factor is scaled by the kinematic Källen function. After summing over all external helicity and color polarizations (dof), |M| 2 is the (squared) Lorentz-invariant matrix describing ij → B scattering. The total parton-level ij → B cross section (σ) is recoverable upon integration over dP S n σ = dσ = dP S n dσ dP S n .
(3.6) While eq. (3.1) is formally proved for only a handful of processes [80], we make the strong but standard assumption that the relation, with appropriate modifications, broadly holds for other processes, including heavy quark and multijet production. For conciseness, we omit insertion of fragmentation functions (J) into eq. (3.1) for exclusive such final states.

Scattering Helicity-Polarized Partons at the Parton Level
In building the expression for unpolarized parton scattering in eq. (3.5), one takes the crucial step of averaging over discrete spacetime and internal quantum numbers for initial-state (IS) partons but only sum discrete degrees of freedom (dof) for final-state (FS) partons. This leads to the familiar IS dof-averaged and FS dof-summed, squared matrix element In dropping all summations over all external helicity eigenstates and fixing the helicities of all external partons in the ij → B process, which we denote generically as withλ representing the set of n helicity eigenstates, one can define the totally helicitypolarized, IS color-averaged and FS color-summed squared matrix element as Generically, we use the terms "polarized matrix elements" and "helicity-polarized matrix elements" interchangeably to mean M(i λ j λ → Bλ). The label "totally helicity-polarized" qualifies that all external partons are in a fixed helicity state, as oppose to instances where only a subset of external partons are in a fixed helicity state. Such configurations corresponds to "partially helicity-polarized" matrix elements and can be constructed analogously. For example: for an unpolarized i, j, and a totally polarized B, the dof-averaged and colorsummed, squared matrix element is Unambiguously, eqs. (3.7)-(3.10) are related by reintroducing helicity averaging / summing: Other configurations, such as with totally or partially polarized IS partons with unpolarized FS partons, can be also constructed so long as helicity averaging factors are consistently accounted. Subsequently, these permutations need not be discussed further. Given a definition for squared matrix elements describing parton scattering with fixed, external helicity polarizations as in eq. (3.9), one can construct scattering observables as done for unpolarized parton scattering. To do this, we promote the fully differentiated scattering rate dσ/dP S n for unpolarized parton scattering in eq. (3.5) by using instead the totally polarized squared matrix elements in eq. (3.9). Explicitly, the fully differentiated scattering rate for the totally polarized partonic process i λ + j λ → Bλ is (3.14) Likewise, for unpolarized IS partons but a totally polarized FS B, the fully differentiated scattering rate is given by It is clear that the relation among eq. (3.14), eq. (3. 16), and the unpolarized case proceeds identically to that established in eq. (3.12). Finally, upon phase space integration one obtains parton-level, total cross sections and differential observables as defined in eq. (3.3).
Unlike eq. (3.7), the polarized expressions of eq. (3.9) and eq. (3.10) are not guaranteed to be Lorentz-invariant quantities as they technically possesses uncontracted Lorentz indices. To be precise: in standard construction of helicity amplitudes within the SM, e.g. Ref. [89], spin-1/2 spinors u m (p, λ f ), v m (p, λf ), spin-1 polarization vectors ε ρ (p, λ V ), and their conjugations, each carry two indices ‡ : one to denote a component within their Lorentz group representation, e.g., m = 1, . . . , 4, and ρ = 0, . . . , 3, and a second to denote their helicity polarization, e.g., λ f , λf = ±1 and λ V = 0, ±1. (Such statements hold also for tensor fields h ρσ , etc., but need not to be discussed further as the conclusions are the same.) The Lagrangian-based formulation of quantum field theory, and hence Feynman rules, leads to scattering amplitudes that are manifestly reference frame-independent for only the first type of index when all such indices are contracted. That is to say, when all u, v spinors are acted upon byū,v, and all ε µ , ∂ µ , σ µν , . . . are acted upon by ε µ , ∂ µ , σ µν , . . . , in some appropriate permutation. Lorentz invariance is only achieved when all indices of the first type are contracted and all indices of the second type are summed. Since helicity polarizations are reference frame-dependent, one must stipulate a reference frame when using eq. (3.9) or its variations. While conceptually simple, for MC event generators this introduces a technical restriction on exploiting Lorentz invariance that is often used in computing matrix elements for unpolarized parton scattering.
In mg5amc, this technicality is managed by exploiting the separately Lorentz-invariant nature of the phase space volume measure given in eq. (3.4). To summarize: A point in phase space is first generated for computing a polarized ME in the same manner as for an unpolarized ME. External momenta are then Lorentz boosted to a definite reference stipulated by the user or to a default option; see appendix (A.3). Helicity amplitudes are then evaluated numerically in this frame. Upon completion of phase space integration, weighted or unweighted events are written to file in standard Les Houches format [90]. In the present implementation, phase space cuts on momenta are applied in the partonic c.m. frame, with the exception of rapidity cuts, which are applied in the lab frame § . In principle, it is also possible to apply phase space cuts in an reconstructable reference frame in mg5amc using the dummy_fct.f capabilities.

Scattering Helicity-Polarized Partons at the Hadron Level
To finally define a version of polarized parton scattering in unpolarized hadron collisions that can be implemented in MC event generators, we argue that the Factorization Theorem of eq. (3.1) can be extended as desired. While a full, field-theoretic derivation is beyond ‡ The two do not have a one-to-one correspondence. For example: for a scalar field φ, the field operator ∂µφ is in a vector representation but possesses only a single (trivial) helicity state.
§ Note that most of the observables defined in run_card.dat are invariant under boosts along the zdirection and are thus the same in the lab frame or the partonic c.m. frame.
this work, principle tenets are already established ¶ in Ref. [80] and references therein. We start by noting that PDFs describing unpolarized partons k out of unpolarized hadrons P can be defined to all orders in α s as a transition amplitude given by [78,80]: Here,Ô kP (x) denotes the composite field operator that extracts parton k with momentum p k = ξP from P. The integral is a Fourier integral that takes the amplitude forÔ kP (x) into momentum space. As eq. (3.17) is defined at the momentum transfer scale Λ NP , a first principle determination of Ô kP , and hence f , is not possible with perturbative methods. That said, in real scattering experiments, massless, initial-state partons are very nearly on their mass shell, indicating that the underlying dynamics of eq. (3.17) occur on a different time scale, τ ∼ 1/Λ NP , and not on the time scale of the hard process, τ ∼ 1/Q. Hence, the dynamics of IS partons i and j are effectively decoupled from the hard scattering process ij → B. Thus, f are factorizable, i.e., can be written as eq. (3.1), up to corrections of the order O(Λ t NP /Q t+2 ), for t > 0. Since f are factorizable, they can be RG-evolved [85] to a cutoff / factorization scale µ f Λ NP , using perturbative methods (DGLAP evolution), and subsequently entered into real scattering computations.
It follows then that IS partons can be approximated as asymptotic states in definite helicity eigenstates. For massless partons, this becomes a matter of splitting the operator O kP (x) for unpolarized k into two orthogonal pieces using chiral projection operators: Here,Ô k λ P is the operator that extracts k with helicity λ from (unpolarized) P. Such partitioning is possible since chiral and helicity eigenstates are identical for massless particles. Consistently, one can decompose the PDF in eq. (3.17) into left-handed (LH) and right-handed (RH) helicity components: The PDF of eq. (3.20) describes the density of a hadron that effectively contains twice as many helicity-polarized parton species as an "unpolarized" PDF. As DGLAP evolution is derivable from perturbative methods, it is possible to decompose them into helicity components as done, for example, in Refs. [91][92][93][94][95]. Alternatively, one can pragmatically bypass a numerical extraction of f k λ /p (ξ, µ f ) by noting that massless SU(3) c ⊗U(1) QED is a parity-invariant theory. In such theories and for unpolarized hadrons, PDFs of partons with opposite helicities are equal, i.e., f k L /p (ξ, µ f ) = f k R /p (ξ, µ f ) [80]. In other words: while it is possible for, say, an u L quark to split into a gluon that splits into an u R quark, such helicity ¶ More specifically, established for only a few inclusive processes at leading power approximations. We assume consistently that the theorem holds for other processes in which perturbative QCD is valid.
depletion and buildup wash out. One can then introduce a normalization factor N = 1/2, and extract PDFs for polarized partons from PDFs for unpolarized partons using: Imposing parity invariance means that the generation of new polarized PDF sets are not needed in real MC simulations. One only needs a wrapper routine to implement eq. (3.21).
Using identical arguments for splitting the DGLAP evolution equations into orthogonal helicity components, the perturbative (in the coupling sense) component of the Sudakov factor ∆ can also be split into permutations of i's and j's helicities λ and λ : A generic implementation of eq. (3.22) is less clearcut than for PDFs. The difference stems from the fact that IS partons, before Sudakov evolution, propagate along the beam axis and therefore possess an azimuthal rotation symmetry. This means that IS partons from polarized and unpolarized parton PDFs transmit the entirety of their polarization information along the beam axis and is captured entirely by matrix elements. However, Sudakov evolution, particularly as implemented via parton showers, injects relative transverse momentum into external partons through IS and FS radiation. In general, this "kick" breaks preexisting rotational symmetry and induces azimuthal spin correlation. Proposals for how to enforce azimuthal spin correlation in MC simulations appear throughout the literature [96][97][98], and their implementation are under active investigation [99][100][101].
Taken all together, a consistent description of parton scattering via polarized ME in unpolarized hadron collisions emerges. Combining the totally polarized and fully differentiated, parton-level scattering rate for i λ + j λ → Bλ in eq. (3.14), with the polarized PDFs of eq. (3.20) and the polarized Sudakov factor in eq. (3.22), the fully differentiated, hadron-level scattering rate for the production of Bλ from partons i λ and j λ is Accounting for all parton species, including those in different helicity states, the production of Bλ from spin-averaged IS partons, in terms of IS states in definite helicity states is Between the second and third lines, we split the single summation over helicity-polarized parton species into a double sum over unpolarized parton species and parton helicities. We then exploit that massless, IS parton species cannot contribute to a scattering requiring the opposite helicity, e.g., f u R fū R ⊗ ∆ u RūR ⊗σ u LūR = 0. Such helicity inversion is proportional to parton masses, and hence vanishing. (We reiterate that in this notation, factorizable k λ → i λ parton buildup and i λ → k λ parton depletion are handled internally by polarized Sudakov evolution.) This allows us to rewrite the helicity-dependent PDFs and Sudakov factor in terms of their helicity-independent counterparts, and demonstrates that the normalization factor N for polarized parton densities in eq. (3.21) can be identified as the spin-averaging symmetry factor in unpolarized parton scattering. Moreover, one sees that the Factorization Theorem of eq. (3.1) is recovered after a summation over FS helicity polarizationsλ, and therefore shows consistency with the above construction. For leading order processes, we report the implementation of eq. (3.26), for helicity polarizations defined in a reconstructable reference frame, at FO into the event generator MadGraph5_aMC@NLO. Importantly, the FO stipulation implies that the Sudakov factor is expanded to zeroth order, i.e., ∆ ij (z) ≈ δ ij δ(1 − z), leading to the simpler relationship: As helicity information is recorded in MC event files at LO, parton-level polarizations can then be passed to a parton shower as desired. We report also the implementation of eq. (3.27) for polarized, colorless, external states, with or without additional, unpolarized QCD partons and heavy quarks, e.g., pp → Z λ + nj, pp → W λ Z λ , or e + R e − L → Z λ + tt, at NLO in QCD. Details are reported in the companion paper Ref. [65].

Decays of Helicity-Polarized Resonances
The polarization features introduced into mg5amc extend also to unstable resonances. In the default usage of mg5amc, the production and decay syntax trigger the so-called spincorrelated NWA [102]. Whereas the usual (spin-uncorrelated) NWA factorizes matrix elements, for example, for qq, gg → tt → tbW − into the product of two decoupled amplitudes, mg5amc instead first generates the helicity amplitude for the 2 → 2 scattering process qq, gg → tt, but replaces the outgoing v(p t , λ t ) spinor with a fermionic Breit-Wigner (BW) propagator for the internal t and the appropriately contracted 1 → 2 decay current. Likewise, for e + e − → W + W − with W + → e + ν e and W − → e − ν e , mg5amc replaces the outgoing polarization vectors ε * µ (p W + , λ W + ) and ε * ν (p W − , λ W − ), which describe W + and W − respectively in the 2 → 2 scattering amplitude for e + e − → W + W − , are each replaced by a bosonic BW propagator and a contracted 1 → 2 decay current.
To propagate the polarization of an unstable resonance to its decay products, we consider modifying this procedure by inserting a "spin-truncated" propagator in lieu of a normal BW propagator. For fermion F and antifermion F with fixed helicity λ, new propagators are defined by denominators with a BW pole structure but a numerator given by the outer product of spinors at helicity λ. Explicitly, the replacement is The origin of this structure stems from the condition that the full propagator is the coherent sum of the spin-truncated propagator over all helicity states λ. That is, For massive gauge bosons, we introduce a similar spin-truncated propagator given by For gauge bosons, the relation of the spin-truncated propagator to the full propagator is different due to gauge theory redundancies, i.e., using 4-component vectors to describe quantities possessing only two or three degrees of freedom. For massive gauge bosons, the full propagator is recovered from Π λ µν by summing over both transverse polarizations, the longitudinal polarization at a given virtuality, and an auxiliary (or scalar) polarization that rapidly vanishes in the on-shell limit [89,103]. (For massless gauge bosons, there is a cancellation between the longitudinal and auxiliary components.) In the unitary gauge for massive spin-1 states, the full and spin-truncated propagators are related explicitly by We report that these propagator decompositions has been implemented in MadGraph5 and MadSpin. While the transverse and longitudinal polarization vectors are defined according to the HELAS convention [68], we set as our auxiliary (or scalar) polarization vector (3.36)

Polarized Vector Boson Scattering in the SM and beyond
Exploring EW VBS is a key step to understanding the SM, and in particular the underlying mechanism of EW symmetry breaking (EWSB). More specifically, VBS is sensitive to whether EWSB is described by more than just the SM Higgs sector due to inevitable disturbances of strong cancellations in amplitudes involving longitudinally polarized weak bosons [1][2][3][4]. As the first observations of VBS were at last achieved by the ATLAS and CMS collaborations during Run II of the LHC program [43-47], their use as a direct probe of new physics is now possible. In a more general setup, VBS is sensitive to peculiar new physics that can be described by dimension-6 and dimension-8 effective operators [104][105][106], assuming the usual decoupling limit [27]. Such new physics would manifest as the anomalous production of EW states in specific helicity configurations, and hence motivates one to investigate means to experimentally disentangle EW boson polarizations.
In this section we investigate VBS production of polarized weak bosons at the √ s = 13 TeV LHC, within the polarized mg5amc framework. We consider EW and mixed EW-QCD production of W + W − boson pairs with helicities (λ, λ ) and two partons at LO, We start in section (4.1) with discussing high energy VBS in the context of Composite Higgs (CH) models. This class of models manifest as an enhancement of scattering rates involving longitudinal weak bosons. Subsequently, we illustrate the mg5amc polarization framework by focusing on the reference frame-dependence of observables, e.g., polarization fractions, built from polarized states and with phase space cuts on particle kinematics. In section (4.2) we extend the study to observables built from the decay of W + λ → µ + ν µ and W − λ → e −ν e , when eq. (4.1) proceeds at O(α 4 ). To do this, we use the MadSpin framework in conjunction with helicity-polarized samples generated from mg5amc. We give special attention to angular observables that are sensitive to the polarization of the parent particle W ± . The same process was studied by the Phantom MC collaboration [52][53][54]56], using the on-shell projection (OSP) technique. The agreement we find not only serves as a check of the two methodologies but also as a basis for future studies. In section (4.3) we repeat this exercise but for when eq. (4.1) proceeds at O(α 2 α 2 s ). Throughout this section, we summarize the new or relevant syntax for mg5amc and MadSpin needed for our study.

Vector Boson Scattering in Composite Higgs Models
In this section we investigate CH models [17][18][19][20][21][22][23][24], in high energy VBS using polarized parton event generation. Promising, modern incarnations of these scenarios predict that the Higgs coupling to weak gauge bosons are rescaled by a common (dimensionless) factor a, and can be described by the effective interaction Lagrangian [25,26] L ⊃ m 2 The presence of a away from unity disrupts fine cancellations in SM amplitudes describing longitudinal weak boson scattering, and leads to amplitudes growing with the invariant mass of the which can potentially be observed at the LHC. Direct measurements of Higgs couplings constrain a at the 95% CL to be a 0.9 [107]. Indirect EW precision data also require a 0.98 [108], but can be relaxed if additional assumptions are satisfied [109].
To quantify the impact of a CH scenario on VBS, we make use of the NLO Higgs Characterization UFO model [74] described in section (2) and focus on the LO EW process for helicity states λ, λ = 0, ±1. Throughout this analysis we do not make the so-called Vector Boson Fusion Approximation, which considers only genuine W W/ZZ → W W scattering diagrams, and is known to neglect significant interference effects [110][111][112]. We instead include all interfering diagrams at O(α 4 ), including other VBS topologies, like γγ → W W , and non-VBS contributions. The mg5amc syntax to model the production of unpolarized and polarized W + W − pairs in unpolarized pp collisions is, respectively, import model HC_UFO-CH generate p p > j j w+ w-QCD=0 QED<=4 generate p p > j j w+{X} w-{Y} QCD=0 QED<=4 Here, one should replace X and Y by all permutations of 0 (longitudinal helicity) and T (transverse helicities) . For event generation, we consider three benchmark scenarios: a = We note that the syntax T coherently sums over both LH and RH helicity states; see appendix (A.1) for details. For further details regarding the usage of mg5amc, we refer readers to Refs. [66,113]. 0.8, a = 0.9, and the SM limit of a = 1. In the UFO model, a is identified as the parameter kSM and can be set in param_card.dat or at runtime with the command set kSM A To define helicity polarizations for the W + λ W − λ pairs, we consider two reference frames: (i) The rest frame defined by the two initial-state partons in the 2 → 4 process, which we label as the partonic c.m. (p-CM) frame. (ii) The rest frame of the W + W − system, which we label as the W W c.m. (W W -CM) frame. Both frames can be specified using the new me_frame selector tag in mg5amc's run_card.dat input file (see appendix (A.3) for details). By momentum conservation, the p-CM frame can be built from either summing the two IS partons' momenta or the four FS partons' momenta. This corresponds to the syntax The W W -CM frame is built most directly from the W + W − itself, and corresponds to [5,6]

= me_frame
To obtain total and differential cross sections * * , we impose the following generator-level phase space cuts, which as described in appendix (A.3), are applied in the p-CM frame: The cuts serve several purposes: First, they regulate collinear and soft singularities from interfering VBS and non-VBS diagrams. Second, they correspond to typical, analysis-level selection cuts that enhance the VBS topology over interfering EW diagrams. Third, they enhance the appearance of new physics. This follows from nonzero a coefficients leading to an enhancement in the scattering amplitude that grows with increasing M (W W ).
In table (1) we show the generator-level cross sections [fb] for unpolarized W + W − production and each W + λ W − λ helicity polarization configuration (λ, λ ), defined in the p-CM, for the CH benchmark scenarios (a = 0.8 and a = 0.9) and the SM limit (a = 1.0). We also report the ratio between the CH and SM rates as well as the polarization fraction f λλ , defined as the ratio of the (λ, λ ) helicity configuration to the unpolarized rate: As expected, nonzero a largely impacts the longitudinal (λ, λ ) = (0, 0) state, which displays roughly a 30% (13%) increase in cross section for a = 0.8 (a = 0.9) over the SM prediction. However, in the absence of more stringent selection cuts, the changes in f λλ indicate that a percent-level determination of polarization fractions would be needed to observe such disturbances. In table (2) we show the equivalent results for polarizations defined in the W W -CM. Only a slight difference is noticed and thus need not be discussed further. Due to rounding errors, the sum of f λλ obfuscates that the sum of the polarization configurations reproduces the unpolarized rate. Satisfying closure requirements in 2 → 4 parton scattering in pp collisions represents a highly nontrivial check of our method.
Turning to differential information, we show in figure (1) the invariant mass distribution of the (W W )-system, dσ/dM (W W ), according to the various helicity configurations for the SM. In the lower panel, we show the differential ratio with respect to the SM, i.e., for the CH scenarios a = 0.8 (dashed lines) and a = 0.9 (solid lines). In figure (1(a)), helicity polarizations are defined in the p-CM frame and in the W W -CM frame in figure (1(b)). We observe explicitly in the lower panel the growing behavior of the CH cross section relative to the SM prediction with increasing M (W W ). We note that tighter selection cuts, such as M (W W ) > 625 − 825 GeV can further enhance the ratio R[M (W W )], though perhaps at a high cross section cost. An alternative possibility is to extract the polarizations via observables built from the W + W − decay products, which we discuss in the next section.

Polarized W Bosons in EW Production of jjW + W −
In weak boson decays to charged leptons, it is well-known that the polarization of a parent boson is imprinted on the kinematics of its decay products. This follows from stable fermions being effectively massless compared to the EW scale. This is especially true of angular observables, which also feature particular sensitivity to the (V − A) structure of bosonic couplings to matter. These observables therefore serve as a test of the SM's chiral structure and, for example, a probe of the coupling structure of new physics.
Here we investigate the production of W + W − λ pairs, via the pure EW process  at LO. The process is defined with an unpolarized W + boson and a polarized W − λ boson with helicity λ = 0, T . We propagate the polarization of the W − λ to its decay products using MadSpin as described in section (3.2). As a high-level check, we also propagate the W − λ polarization using the OSP method [62][63][64]. In the context of VBS, the OSP technique has been used in Refs. [52][53][54][55][56] and is implemented in mg5amc under an unsupported, standalone development branch for the purpose of this work. In short, the method amounts to setting the momenta of the W ± bosons, k ± , to their mass-shell values (k 2 ± = M 2 W ) in the numerators of matrix elements for the full 2 → 6 EW process qq → qq e − ν e µ + ν µ . The virtuality k 2 ± in the denominator of propagators is allowed to float. Non-resonant diagrams are neglected. In practice, the k 2 ± are restricted to the neighborhood of k 2 ± = M 2 W using phase space cuts, thereby approximating the spin-correlated, NWA employed by MadSpin.
For the process in eq. (4.7), we define the polar angle θ as the angle between the W − flight direction in the p-CM frame and the e − flight direction in the W − rest frame, i.e., Here, p W is the 3-momentum of the W − in the p-CM frame and p e is the 3-momentum of the e − in the W − rest frame. Similarly, an azimuthal angle φ can be defined as the opening angle between the W − boson's production plane (defined by the W − and beam direction) and its decay plane. Analytically, this is given by where the two unit vectorsv 1 andv 2 are defined to be, Here P i is the 3-momentum of any of the IS partons in the W − rest frame. For definiteness, we choose the i that makes the smallest opening angle with W − 's flight direction in the p-CM frame. If one defines the z-direction along p W , we can identifyv 1 =ŷ as the unit vector in the y-direction,v 2 =x, and we reproduce the coordinate system of Ref. [114].
In this convention, the matrix element M λ , with λ defined in the p-CM frame, describing eq. (4.7) depends on θ and φ through the W − λ → e −ν e decay, and specifically through the angular dependence of the e − spinor. Hence, the θ and φ dependence of M λ scales as Now, at beam-symmetric experiments such as the LHC, since the momenta of quarks are typically larger than of antiquarks, and since the EW gauge couplings to fermions are chiral, the emission rates of RH, LH, or longitudinal W, Z bosons off initial-state parton lines differ. The relative emission rates are further skewed by non-Abelian gauge and Higgs couplings. Ultimately, this leads to asymmetric production of polarized W − λ for not just single boson production [115] but also in multiboson processes [65]. Allowing for an asymmetric production of W − λ , and in the notation of Ref. [116], the differential cross section for inclusive W − production in terms of the angles θ and φ can be written as [117]: Here, f λ can be interpreted as the likelihood of producing W − λ with helicity λ in the inclusive process. The g λλ result from interference among the different W − λ helicity polarizations but vanish upon integration over φ ∈ [−π, π]. The precise form of eq. (4.13) is somewhat arbitrary as trigonometric identities and redefinitions relate coefficients here to other parameterizations; see, e.g., Ref. [118]. After integrating, one obtains the familiar expression: In the following we use eq. (4.13) and eq. (4.14) as guiding relationships to explore the polarization of W − λ → e − ν e decays in eq. (4.7). We also consider the ability to extract the coefficients f λ using the MC methods developed and reported in section (3). Throughout the following we use the OSP technique as a benchmark to quantify our results.

Production and Decay of Polarized W Bosons in mg5amc
To simulate eq. (4.7) in the SM using mg5amc+MadSpin, we first use the syntax reported in section (4.1) to generate the subprocess qq → qq W − λ W + . The mg5amc commands are , of the pure EW process pp → jjW + W − λ , with W + → µ + ν µ and W − λ → e −ν e , in the SM for unpolarized W + and W − helicity polarization λ, defined in the p-CM frame, assuming generator-and analysis-level cuts of eq. (4.15) and eq. (4.23), and using the MadSpin and OSP decay schemes for W ± . import model sm generate p p > w+ w-{X} j j QED<=4 QCD=0 One should replace X in the generate command with 0 or T for λ = 0 or λ = T . We consider the polarizations defined in the p-CM frame, and thus in run_card.dat set me_frame = [1,2]. For comparison, we consider also the unpolarized process, which is simulated by generate p p > w+ w-j j To regulate infrared poles in the matrix element and enhance the VBS topology over interfering topologies, we require events to fulfill the generator-level cuts: We relax cuts relative to eq. (4.4) since we are not strictly interested in isolating the pure VBS topology. After event generation, unpolarized and polarized W + W − λ pairs are decayed using MadSpin. As described in appendix (A.2), the MadSpin syntax is the same for unpolarized and polarized W λ . Therefore in the madspin_card.dat file, one only needs: decay w+ > mu+ vm decay w-> e-veT he process with identical final states and kinematic cuts is also simulated using our implementation of the OSP method. In table (3), we report generator-level cross sections [fb] and polarization fraction (f λ ) [%] for the full 2 → 6 process using both the MadSpin and OSP decay schemes. We report good agreement in generator-level normalizations † † . † † As implemented, we also report comparable MC efficiency between the MadSpin and OSP techniques.

Leptonic Observables from Polarized W Boson Decays
We now turn to kinematic observables built from final-state charged leptons in the EW process of eq. (4.7). Throughout this section we present in upper panels of plots overlapping distributions for unpolarized jjW + W − production (black; dash-double dot), trans-  versely polarized W − λ production (green), longitudinally polarized W − λ production (blue), and polarization-summed W − λ production (red). W + W − decays are treated using the Mad-Spin (solid) and OSP (bar) methods. For unpolarized production, we only use MadSpin. To quantify potential disagreement between the two decay techniques and unless specified, for each observable we also report in lower panels of plots the OSP-to-MadSpin ratio of the polarized and polarization-summed curves. In summary, we find good agreement with the shape and normalization between the MadSpin and OSP samples. Differences are consistent with MC statistics and therefore demonstrate strong checks of both the methods.
We start with figure (2(a)), which shows the polar distribution cos θ as defined in eq. (4.8). As anticipated from eq. (4.14), we observe that the longitudinal component of W − λ exhibits a polar dependence that behaves as dσ(W λ=0 )/d cos θ ∼ sin 2 θ, while the transverse modes are given as a coherent sum of left and right contributions. We see a preference for dσ(W λ=T )/d cos θ > 0, indicating that f L > f R , and consistent with above arguments that the production of W − λ=−1 is preferred over W − λ=+1 at the LHC. In figure (2(b)) we show the azimuthal distribution φ as defined in eq. (4.9). Notably, the longitudinal mode exhibits a flat behavior and the transverse modes oscillate. This follows from eq. (4.13), which shows that the λ = 0 polarization is only sensitive to φ through λ = T interference terms; these are neglected in polarized computations. Consistently, for λ = T , we observe the cos 2φ behavior that originates from the λ = ±1 interference, which is modeled since λ = ±1 modes are summed coherently. In comparison to the unpolarized sample, which sums all polarizations coherently, and the semi-incoherent sum of λ = 0 and T polarizations, we find that the difference is small for both the MadSpin and OSP schemes. This suggests that the g L0 , g R0 interference terms are small, and that the interference is dominated by g RL . The difference between the unpolarized curve and the MadSpin sum quantifies the interference between λ = 0 and T modes, and is small.
To extract the polarization fractions f λ from the distribution in figure (2(a)), we use the Legendre expansion technique as used by Ref. [54] for VBS, which is related to the moment method used by Refs. [114,118] for pp → W ± λ + nj. We start by noting that the polar distribution of eq. (4.14) can be written in term of first Legendre polynomials, with After explicit integration, the polarization fractions in terms of Legendre coefficients are: We can extract the values of α l for l = 0, . . . , 2, from our simulated predictions (or from data) by performing a sum over each histogram bin and by approximating the α l integral: In the first line, we make the change of variable x = cos θ for clarity. In the second, we partition the integral into a large number of disjoint integrals over continuous ranges (bins), such that the bin widths satisfy |x k+1 − x k | 1. We then factor the normalized histogram weight g(x) using the Mean Value Theorem, and thereby approximate g(x) as a constant g(x * ) for x * ∈ [x k , x k+1 ]. This allows us to evaluate the integrals exactly.
Choosing x * = x k , i.e., the bin starting boundary, we obtain the following fractions The uncertainty we report is statistical, but other theory uncertainties, e.g., scale uncertainties, can be propagated in a straightforward manner. The reconstructed distribution in cos θ is shown in figure (3(a)), together with the simulated expression. For nearly the entire domain of cos θ, we report a good reproduction of the MC simulation using Legendre polynomials. A large disagreement is observed at the boundaries, near cos θ = ±1, for the λ = 0 distribution. This can be attributed simply to the fact that the distribution dσ(W − λ=0 ) itself vanishes smoothly at the endpoints, resulting in an ill-defined ratio.

Impact of Selection Cuts on Polarized Distributions
As noted in the double differential distribution in eq. (4.13), interference terms between different W − λ polarizations possess a dependence on the azimuthal angle φ. Hence, observables integrated over φ are insensitive to the interference between the different polarization terms. For such observables, the incoherent sum of transverse and longitudinal contributions agrees with the unpolarized prediction. However, it is known that realistic experimental conditions are not totally inclusive with respect to φ due to kinematical selection cuts, which are motivated by detector acceptance or analysis criteria. To study the impact of realistic selection cuts on our modeling of polarized matrix elements as well as the residual size of possible polarization interference, we consider the following selection cuts applied to the decay products of W − λ :  The selection cuts here are applied at the analysis-level and are in addition to the generatorlevel cuts of eq. (4.15). We report the resulting cross section for both the MadSpin and OSP samples in table (3). Observed differences between the two are consistent with MC statistical uncertainty and rounding errors. We find that while there is a 20−25% reduction in cross section, the W − λ polarization fractions remain essentially the same. We do, however, see the emergence of a sub-percent discrepancy between the incoherent sum of helicitypolarized cross sections and the full, unpolarized cross section when the selection cuts of eq. (4.23) are applied to decay products.
In figure (2(c)) and figure (2(d)), we show the same polar and azimuthal observables described above and shown in figure (2(a)) and figure (2(b)) but after applying the selection cuts of eq. (4.23). In comparing to the distributions without cuts, a large impact can be noted, namely a total depletion of events in the neighborhood of cos θ = −1. This results in an increased forward-backward asymmetry and stems from the softer nature of "backwards" flying e − originating from W − L decays, which are selected out by the p l T cut in eq. (4.23) [114,118]. As the concavity of the λ = T curve, and hence the helicity polarization-summed curve as well, now changes multiple times as a function of cos θ, the functional form of eq. (4.14) does not serve as a good description of figure (2(c)). Therefore, to recover polarization fits as reported in eq. (4.22), modern unfolding techniques are necessary; see, for example, Refs. [119][120][121][122][123] and references therein. Such techniques "correct" reconstructed distributions/observables for real, detector-level and analysis-level acceptance efficiencies, enabling more direct comparisons to truth-level, MC predictions [119,120]. This, however, comes at the cost of introducing systematic uncertainties stemming from imperfect modeling of underlying physics and detector response. The availability of a polarized MC event generator, which is a main result of this work, can significantly help to reduce such uncertainties. For example: the ability to generate specific helicity samples provides a means to directly model detector responses to kinematic regions that are strongly suppressed in the SM or exhibit ultra low detector acceptance efficiencies, say from a forward-backward asymmetry, and subsequently help ameliorate singularities that may otherwise appear in unfolding response matrices.
For the azimuthal distribution, we observe similar shapes to the generator-level cut curves, albeit with larger maxima and minima differences. Lastly, to reiterate, we observe good agreement between the OSP and MadSpin samples.
In figure (4)(a) and figure (4)(b) we show the p e − T distributions before and after the cuts on the decay products (eq. (4.23)) respectively. We observe a small difference between the incoherent sum of polarizations with respect to the unpolarized simulation, which we attribute to the interference between longitudinal and transverse polarizations in some restricted region of phase space. In figure (4)(c) and figure (4)(d) we show the invariant mass of the di-lepton system m(e − , µ + ) before and after the cuts in eq. (4.23) respectively. Unlike the p e − T , no difference between unpolarized and polarized samples can be observed. Turning to more reconstructed objects, we show in figure (5)(a,b) the m jj distribution and in (c,d) the M (W + W − ) distribution, assuming only (a,c) generator-level cuts and (b,d) with the cuts of eq. (4.23). We find that both before and after eq. (4.23) the observables are insensitive to interference and that the incoherent polarization sum describes the unpolarized distributions well. We attribute this insensitivity to the fact that interference effects appear first at the W − λ decay level, though the angle φ as defined in eq. (4.9). By working at the W − λ level, we are inclusive with respect to φ, leading to a washout of interference effects. By identical arguments, an insensitivity to interference can be found in figure (6) In all distributions and cross sections we find good agreement between the MadSpin and OSP method. We also find the interference effect between transverse and longitudinal polarization channels to be small, and the incoherent sum of polarization describes the distributions we consider to a good degree. The largest difference, although still small, is observed in the p T (e − ) distributions. The difference remains negligible even after applying selection cuts defined in eq. (4.23). Of course, this observation somewhat follows the fact that this process is dominated by transverse modes, it is hard to access the effect of the longitudinal bosons, and interference with the transverse modes are even less accessible.

Polarized W Bosons in Mixed EW-QCD Production of jjW + W −
As a final case study, we consider the LO production of the mixed EW-QCD process (4.24) at O(α 2 s α 2 EW ). Aside from its own interesting features, the process is a primary background for the pure EW process jjW + W − λ at O(α 4 EW ). Subsequently, in this section, we discuss the similarities and differences in distributions between EW and mixed EW-QCD production of jjW + W − λ for the same observables discussed in section (4.2). We use both the MadSpin and OSP methods to treat the decays of W bosons. The mg5amc syntax for eq. (4.24) is where the polarization (X) of W − λ is set to T or 0. The MadSpin syntax and all phase space cuts are the same as those reported in section (4.2).
In analogy to table (3), we report in table (4) the cross sections [fb] and helicity fractions (f λ ) [%] for the full process 2 → 6, assuming generator-level cuts of eq. (4.15) and analysislevel cuts of eq. (4.23), using both the MadSpin and OSP decay schemes. Compared to the pure EW process, which shows σ(W − λ=T ) : σ(W − λ=0 ) ratio of about 4 : 1, the mixed EW-QCD process here exhibits a bigger difference of about 6 : 1. This difference can be attributed to the fact that most W ± λ in the mixed EW-QCD case are emitted off massless fermion legs, which only contribute to the W ± λ=T process. The pure EW process permits W ± λ=0 production through diboson and VBS type of scattering topologies. We start our investigation of differential observables sensitive to W − λ polarization with figure (7), where we show the polar (θ) and azimuthal (φ) angular distributions as defined in eqs. (4.8)-(4.9). As in figure (2), panels (a,b) include only generator-level cuts while (c,d) include analysis-level cuts listed in eq. (4.23). In comparison to the pure EW process, we observe a smaller fraction of W ± λ=0 events, consistent with results reported in table (4). In contrast to the EW process shown in figure (2), we observe a milder g T 0 interference pattern in the φ distribution with both sets of phase space cuts. The smaller interference pattern can be attributed to the smaller λ = 0 contribution. Crucially, we note that the shape of distributions are not substantially affected by the production mechanism; only the normalizations are strongly affected, reflecting the processes' different coupling orders. We extract the helicity fractions f λ from figure (7(a)) using the Legendre expansion technique described in eqs. (4.16)-(4.21). For the EW-QCD process, we report f L = 0.5248 ± 0.3%, f R = 0.3307 ± 0.4%, f 0 = 0.1445 ± 2%. (4.25) In comparison to the pure EW process, we observe strong similarities for the production of W − λ possessing LH polarizations, with f L ≈ 50%. For the RH and longitudi- nal polarizations, we see an increase (decrease) of the λ = +1 (λ = 0) modes of about δf λ ∼ +5% (−5%). In figure (3(b)) we show the polar distribution as reconstructed from eq. (4.25) as well as using MadSpin. As in the pure EW process, we find good agreement between the Legendre expansion and MadSpin, except at the boundaries. There the distribution dσ(W − λ=0 ) vanishes and our ratios quantifying disagreements become ill-defined. We report similar shape behaviors between the EW and EW-QCD processes in figure (8), where we show for the EW-QCD process the (a,b) p T (e − ) and (c,d) m( ) distributions, assuming (a,c) only generator-level cuts and (b,d) analysis-level cuts. The analogous distributions for the EW process are shown in figure (4). We find a much smaller g T 0 in- terference pattern for the p T (e − ) curves here, due in part to the smaller W − λ=0 component. Turning to the (a,b) dijet invariant mass m(jj) and (c,d) diboson invariant mass m(W + W − ) in figure (9), we observe large differences with respect to the pure EW process in figure (5). With and without analysis-level cuts, we see that both the dijet and diboson spectra here strongly peak toward smaller invariant masses; the shape is driven mostly by the λ = T modes. The dσ ∼ 1/m k (jj) behavior is typical of s-channel g * → qq splittings and suggests that the mixed EW-QCD process is not driven by valence-valence scattering. This is opposed to the EW process which shows a plateau in the dijet spectrum and a softer peaking of the diboson mass, which are consistent with VBS-like topologies. In both sets of distributions we find that the impact of the g 0T interference is negligible.
Finally, in figure (10) we show in (a,b) the p T (W − ) distributions and in (c,d) the η(W − ), assuming only (a,c) generator-level cuts and (b,d) with analysis-level cuts. For the unpolarized and the incoherent summation curves, we observe little differences between the mixed EW-QCD process here and the pure EW process in figure (6). By individual polarizations, however, we observe that the λ = 0 and λ = T polarizations in the EW-QCD process possess slightly broader peaks than their pure EW counter parts. This feature is hidden because the EW process possesses a relatively larger λ = 0 fraction than the mixed process (see eq. (4.22) and eq. (4.25)), and that the narrower peaks of the λ = 0 and λ = T polarizations in the EW process are more widely separated than in the EW-QCD process. This in turns broadens the polarization-summed curve in the EW process. As a result of this preference for a higher p T , the η(W − ) distributions for both λ = 0 and λ = T in the mixed process are more central than their pure EW counterpart. This is particularly striking when comparing the two λ = 0 curves. In the EW-QCD case, the broad but central single-bump shape is indicative of a moderate recoil against the dijet system, and consistent with process not being driven by valence-valence scattering. In the EW case, the forward, double-bump shape is indicative of forward W − production via VBS.

Conclusions
The SM of particle physics remains the best description to date of how nature functions at small distances and high momentum-transfer scales. This is especially the case in light of a SM-like Higgs boson and the multitude of data collected during Runs I and II of the LHC [37-47]. However, the unambiguous evidence for dark matter and nonzero neutrino masses, as well as theoretical demands to understand the origin of flavor and the stability of the Higgs's mass, require extending the SM. Among the viable solutions are scenarios that predict the production of fermions and EW gauge bosons in high-energy scattering processes that are polarized in a distinctly different manner than that predicted by the SM. Consequently, searches for the anomalous polarization of SM particles represent an important and well-motivated component of the LHC's program.
To facilitate such studies, we report the development of a method for using polarized matrix elements, where some or all external states are in a definite helicity eigenstates and where spin-averaging or spin-summing is truncated or not present, in the publicly available event generator MadGraph5_aMC@NLO. For an arbitrary reference frame, partonic scattering and decay rates of asymptotic states with fixed helicity polarizations can be computed at LO, with little impact on runtime, and supports particle spins up to 3/2 and 2. The helicity polarizations of resonances are transmitted to their decay products via modifications to their propagators. Furthermore, our framework can be used beyond the scope of the LHC from low energy physics to astrophysics.
The scattering formalism underlying our work and main implementation details are given in section (3); technical and usage details are reported in appendix (A). As case studies, we investigated the production and decay of polarized W + λ W − λ pairs in the process pp → jjW + λ W − λ , with helicity eigenstates (λ, λ ) defined in various reference frames. We considereded a benchmark Composite Higgs scenario (section (4.1)) as well as SM production at O(α 4 ) (section (4.2)) and O(α 2 α 2 s ) (section (4.3)). We focused on the helicity polarization decomposition of processes according to their reference frame as well as investigated the impact of typical generator-level and analysis-level selection cuts. In all case studies, we found that accounting for interference between LH and RH W λ bosons is much more important than interference between transverse and longitudinal polarizations. Investigations into the production and decay of polarized EW bosons beyond tree-level are reported in a companion paper [65].

A Polarized Matrix Elements in MadGraph5_aMC@NLO
In this appendix we expand on the formalism reported in Sec. 3 and describe the main features introduced into the event generator mg5amc that allows the modeling of parton scattering with polarized matrix elements. That is, matrix elements where some or all external states are in a definite helicity eigenstates and where spin-averaging or spin-summing is truncated or not present. In appendix (A.1) we describe the new syntax that triggers the creation of scattering amplitudes with a truncated polarization summation. Decays of polarized resonances are described in appendix (A.2). Leading order event generation within a reconstructable reference frame is described in appendix (A.3), while in appendix (A.4) the possibility of event re-weighting of polarized samples is discussed.

A.1 Syntax for Polarized Matrix Elements
In order to fix the helicity polarization of particles in mg5amc, we introduce new syntax commands at the process-definition and event-generation levels. When specifying a scattering or decay process using the usual [66] mg5amc commands, any particle followed immediately (without spacing) by {X} will be polarized in the helicity eigenstate "X". We stress that the notion of helicity polarization is not Lorentz invariant for massive particles. Consequently, using the polarization syntax requires that a reference frame be specified at the time of matrix element evaluation. For massive spin 1/2, 1, 3/2, and 2 particles, we list in table (5) the allowed syntax, the corresponding helicity states in the HELAS basis, and whether the polarization can be transmitted through the propagators of massive particles (see also section (3.2) and appendix (A.2) for details).  Table 5. For a given particle spin, the allowed mg5amc polarization syntax, its helicity state in the HELAS basis, and whether the polarization is transmitted through propagators of massive particles.
At LO, the bracket polarization syntax can be used for any initial-state (IS) or finalstate (FS) particle in any scattering process. Examples of such usage are: generate p p > t t~{R} generate e+{L} e-> w+{0} w-{T} generate z z{R} > w+ w-{0} which respectively describe the Born-level processes: The helicity label 0 denotes a longitudinally polarized massive vector boson; L and R represent LH and RH helicity polarizations for spin 1/2 and 1 particles; and T models transverse polarizations of spin 1 particles as a coherent sum of L and R helicities. Throughout this following, omitting a helicity label expresses an unpolarized particle. The {X} polarization syntax can also be used with multi-particle definitions. For example: to model the diboson process pp → W ± T W ∓ 0 , the following commands are possible: define ww = w+ wgenerate p p > ww{T} ww{0} To avoid polarization definition conflicts, multi-particle definitions consisting of polarized states, e.g., define wwX = w+{T} w-{0}, is not allowed. In standard computations using mg5amc, once a process has been defined, e.g., generate p p > t t∼, the MadGraph sub-program [113,124] will build all helicity amplitudes from ALOHA [70] and HELAS [68] routines, for all contributing sub-channels, e.g., gg, qq → tt, and for all external helicity permutations, e.g., t L t L , t L t R , t R t L , and t R t R . Next, amplitudes are evaluated numerically, squared, and summed. For initial states and identical final states, dof. averaging and symmetry multiplicity factors are then incorporated. When using the polarization features on IS/FS particles, this procedure is changed in two ways:  Table 6. For the gravitino scattering process g rv (p 1 , λ) + g rv (p 2 ) → τ + (p 3 ) + τ − (p 4 ), the squared scattering amplitude |M| 2 as a function of gravitino helicity λ, defined in the partonic c.m. frame at the phase space point provided in the text, as well as the unpolarized, spin-averaged squared matrix element |M| 2 , and the sum of the four squared amplitudes.
• Instead of summing over all helicity polarizations of all external particles, mg5amc only sums over the polarizations allowed in the process definition.
• Averaging symmetry factors over initial state polarizations are modified according to the new number of initial states.
Special attention is needed for processes with identical particles. For initial state particles, the mg5amc convention is that the order of the particles during process declaration matters. The ordering condition is particularly suited for asymmetric beam experiments, where the first state is associated to one beam and the second state to another beam. As an example, consider the production of an unpolarized τ + τ − pair from massless, spin 3/2 gravitino scattering [125] with fixed external momenta: This can be simulated using the syntax import model GldGrv_UFO generate grv{X} grv > ta+ taoutput standalone Polar_grv_grv_tau_tau; launch -f In the above, g rv (p 2 ) is unpolarized and eq. (A.2) is not equivalent to g rv (p 1 )g rv (p 2 , λ) → τ + τ − . In order to recover the unpolarized process g rv g rv → τ + τ − , one has to incoherently sum over the four helicity configurations, λ = ±1, ±3. And as described above, spinaveraging over possible initial states is modified to only account for relevant dof. As a check of modeling asymmetrically polarized, initial-state particles that are identical, we report in table (6)  the squared scattering amplitude |M| 2 of eq. (A.2) for polarizations λ = ±1, ±3, as well as the unpolarized, spin-averaged squared matrix element |M| 2 , and the sum of the four squared amplitudes. One sees precisely that the difference between the spin-averaged result and the summed result is the symmetry factor (2s ggv + 1) = 4, for s ggv = 3/2. For identical, final-state particles, the mg5amc convention demands that each identical particle has a specified polarization. For example, the process implies a sum over all interfering diagrams with one transversely polarized Z and one longitudinally polarized Z. In this sense, pp → Z 0 Z T is equivalent to pp → Z T Z 0 . To recover the full, unpolarized process, pp → ZZ, ones must sum the three helicity configurations: Z 0 Z T , Z 0 Z 0 , and Z T Z T . For identical, final-state particles, a mixed syntax where some identical particles are polarized and others are not, e.g., pp → Z 0 Z, is not supported.

A.2 Decays of Polarized Resonances with MadGraph5 and MadSpin
As described in section (3.2), the helicity polarization features introduced into mg5amc extend to unstable resonances. After specifying a hard scattering or decay process for a massive, polarized final state at LO, one can steer the decay of a resonance to the desired final state in the usual manner [102]. For example: the syntax to model the production and decay of t L or W + 0 W − T pairs at LO is generate p p > t t~{L}, t~> b~wgenerate e+{L} e-> w+{0} w-{T}, w+ > e+ ve, w-> e-veT he {X} syntax changes the standard mg5amc decay protocol by inserting the spin-truncated propagators defined in eqs. (3.30)-(3.31), eq. (3.34) and eq. (3.36), instead of a normal BW propagator. Special care has been taken for the case where the transverse momentum of a spin 1 boson is vanishing in order to consistently adhere to the limit employed by HELAS.
The inclusion of polarized propagators is possible through the extension of the ALOHA package [70] to support non-Lorentz invariant quantities and the auxiliary polarization λ = A, defined in eq. (3.36). The polarization can be called explicitly using the syntax generate p p > z{T} z{A}, z > e+ e-, which describes resonant diboson production qq → Z T Z A , Z λ → e + e − . In principle, the λ = A polarization vector is needed to recover unpolarized events from polarized event samples, particularly in the off-shell region. However, its kinematical structure leads to a highly suppressed or vanishing contribution in practical applications.
In the presence of identically polarized, final state particles, the handling of symmetry factors and optimization of phase space integration requires care. As such, two polarization modes has been implemented: (i) If the user specifies exactly one decay for each polarized particle, like in the following: generate p p > z{X} z{Y}, z > e+ e-, z > mu+ mugenerate p p > z{X} z{Y}, z > l+ l-, z > j j  for the process e + e − → tt at √ s = 1000 GeV, assuming unpolarized particles, totally and partially polarized beams in the partonic c.m. frame using the polbeam1 steering commands (polbeam2=0), totally polarized IS particles in the partonic c.m. frame using the polarization {X} syntax, and the anticipated cross section as derived from the polbeam1 results.
Here we report a statistical error are at 2. 10 −4 pb.

decay w+ > e+ ve
The {X} command is also supported by MadSpin itself. This allows one to force some particles in a decay chain into a fixed helicity polarization that is defined in the same frame as the original, production-level events. Such steering can be called using the commands: This describes the decay of a top quark t into an unpolarized b quark and a transversely polarized W + T boson, which in turns decays to electron-flavored leptons. ‡ ‡

A.3 Event Generation with Polarized Partons
As stressed throughout this text, scattering particles with fixed helicity polarizations requires one to fix a reference frame in order to meaningfully define individual polarization vectors and spinors. For LO processes, this is possible at the event-generation level using the new "matrix element frame" parameter me_frame in the run_card.dat steering file. The parameter is displayed by default only if at least one massive particle is polarized but can technically be used for any processes.

A.4 Event Re-weighting for Arbitrary Reference Frames
A key feature of the Re-Weighting module [126] in mg5amc is the ability to take an event sample defined by one process definition and, within reason, generate a new event sample defined by a second process definition through matrix element re-weighting. It is therefore also possible to use new the polarization syntax in conjunction with the Re-Weighting module, allowing one to study the impact of polarization via re-weighting methods.
In order to have meaningful helicity polarizations one needs a definite reference frame as in previous considerations. By default, the Re-Weighting module will use the frame defined in the run_card.dat file but also allows a user to define an alternative frame. However, Attribute Description status Returns −1 (1) [2] for an initial (final) [intermediate] state particle mother1 Returns the first progenitor particle object.
mother2 Returns the second progenitor particle object if both progenitors have status= −1, otherwise returns mother1. color1 First color index for the leading color associated to the particle color2 Second color index for the leading color associated to the particle  Table 9. List of common Les Houches event file attributes available to the Re-Weighting module; see Ref. [90] for further details. * denotes that the quantity is defined in the lab frame by default.
since the module can interface with a generic Les Houches event file [90,127,128], we have designed a specific syntax for building new frames. The user must simply provide a pythonbased lambda function that selects the particles to include in the Lorentz-boost definition. The fundamental ideas and effects are the same as simulating polarized particle scattering following appendix (A.1); only the procedure for defining a reference frame differs. Particles whose momentum are to be included in the frame definition can be identified through any of the preexisting Les Houches event file attributes [90]. A list of common attributes that can be used is given in table (9). Some examples (and their impact) include: change boost True # use to lab-frame change boost lambda p: p.status==-1 # go to partonic-center-of mass frame change boost lambda p: p.pid in [24,-24] # go to the ww rest-frame