Electroweak Splitting Functions and High Energy Showering

We derive the electroweak (EW) collinear splitting functions for the Standard Model, including the massive fermions, gauge bosons and the Higgs boson. We first present the splitting functions in the limit of unbroken SU(2)xU(1) and discuss their general features in the collinear and soft-collinear regimes. We then systematically incorporate EW symmetry breaking (EWSB), which leads to the emergence of additional"ultra-collinear"splitting phenomena and naive violations of the Goldstone-boson Equivalence Theorem. We suggest a particularly convenient choice of non-covariant gauge (dubbed"Goldstone Equivalence Gauge") that disentangles the effects of Goldstone bosons and gauge fields in the presence of EWSB, and allows trivial book-keeping of leading power corrections in the VEV. We implement a comprehensive, practical EW showering scheme based on these splitting functions using a Sudakov evolution formalism. Novel features in the implementation include a complete accounting of ultra-collinear effects, matching between shower and decay, kinematic back-reaction corrections in multi-stage showers, and mixed-state evolution of neutral bosons (gamma/Z/h) using density-matrices. We employ the EW showering formalism to study a number of important physical processes at O(1-10 TeV) energies. They include (a) electroweak partons in the initial state as the basis for vector-boson-fusion; (b) the emergence of"weak jets"such as those initiated by transverse gauge bosons, with individual splitting probabilities as large as O(30%); (c) EW showers initiated by top quarks, including Higgs bosons in the final state; (d) the occurrence of O(1) interference effects within EW showers involving the neutral bosons; and (e) EW corrections to new physics processes, as illustrated by production of a heavy vector boson (W') and the subsequent showering of its decay products.


Contents
1 Introduction

Electroweak parton showers
Process-independent parton showers in QED and QCD have long served as invaluable tools for particle physics in high energy collisions and decays. By exploiting formal factorizations between hard/wide-angle physics and soft/collinear physics [1][2][3], the extremely complicated exclusive structure of high energy scattering events can be viewed in a modular fashion. The dominant flows of energy and other quantum numbers are modeled with manageable, low-multiplicity matrix elements. These are subsequently dressed with soft/collinear radiation, and hadronization applied to bare color charges. Detailed implementations have varied significantly in specific approach, but showering programs such as PYTHIA [4], HERWIG [5], and SHERPA [6] are now standard workhorses required for describing realistic collider events. They have also found widespread use in modeling the interactions of high-energy cosmic rays [7], as well as the exclusive products of dark matter annihilation and decay [8,9]. Collinear parton showers become a ubiquitous phenomenon for processes at energies far above the mass scales of the relevant final-state particles, such as the electron mass in QED or the confinement scale in QCD. With the upgraded LHC and proposed future accelerators [10][11][12] and a growing suite of instruments sensitive to indirect signals of multi-TeV dark matter [13][14][15], we are now forced to confront processes at energies far above the next known mass threshold in Nature, the electroweak (EW) scale v ≈ 246 GeV (the electroweak vacuum expectation value, "VEV" in short). Consequently, we are entering a phase in particle physics where it becomes appropriate to consider electroweak parton showers, extending the usual SU (3) QCD × U (1) EM showers into the fully SU (3) QCD × SU (2) L × U (1) Y symmetric framework of the Standard Model (SM). In effect, we will start to see electroweak gauge bosons, Higgs bosons, and top quarks behaving like massless partons [16,17], appearing both as constituents of jets [18] as well as of initial-state beam particles. This is in stark contrast to the conventional perspective in which they are viewed as "heavy" particles that are only produced as part of the hard interaction.
The concept of electroweak bosons as partons has a long history, beginning with the "effective-W approximation" [19][20][21]. This picture of electroweak vector bosons radiating off of initial-state quarks is now strongly supported by the experimental observation of Higgs boson production via vector boson fusion (VBF) at the LHC [22]. As we imagine probing VBF-initiated processes at even higher energies, with both the initial weak bosons and their associated tag jets becoming significantly more collinear to the beams, the idea of weak parton distribution functions (PDFs) within protons becomes progressively more appropriate.
Many calculations have further revealed large negative electroweak virtual corrections to a variety of exclusive high-energy processes, wherein real emission of additional weak bosons is not included. Such large "non-emission" rate penalties indicate the onset of the universal, logarithmically-enhanced Sudakov form-factors characteristic of massless gauge theories [23,24]. For example, exclusive di-jet production receives corrections from virtual W/Z exchange that begin to exceed −10% for transverse momenta exceeding 3 TeV [25,26], and grow to approximately −30% at the 10's of TeV energies expected at future hadron colliders. For processes that include weak bosons at the hard event scale, such as γ/Z/W +jets or vector boson pair production, the corrections can quickly grow to O(1) [27][28][29][30][31][32][33]. A process-independent framework for extracting all such log-enhanced electroweak virtual corrections at fixed leading-order has been developed in [34,35], and next-to-leading logarithmic resummation of the gauge corrections has been achieved using SCET formalism in [36][37][38][39][40].
The total rates of real W/Z emissions and other electroweak parton splittings have a direct correspondence with the "lost" event rates encoded in the negative electroweak virtual corrections, with matching logarithmic enhancements in accordance with the Kinoshita-Lee-Nauenberg theorem. Iterating this observation across all possible nested emissions and loops within a given process builds up the usual parton shower picture, allowing formal resummations of the logarithms that would otherwise still appear in well-defined exclusive rates. Many studies have addressed aspects of electroweak parton showering in the past several years [41][42][43][44][45][46][47][48][49]. Parts of the complete shower are already available in public codes and are being tested at the LHC, with ATLAS recently making a first observation of collinearenhanced W/Z radiation within QCD jets [50]. A detailed listing of electroweak collinear splitting functions and PDF evolution equations, restricted to processes that survive in the unbroken limit, has been worked out in [43]. There, the effects of electroweak symmetry breaking (EWSB) are addressed minimalistically by including a hard phase space cutoff and working in a preferred isospin basis. These results and more recent SCET-based calculations have also been adapted for the problem of TeV-scale dark matter annihilation in [51][52][53][54][55][56][57]. For general-purpose applications, recent versions of PYTHIA incorporate radiation of W and Z bosons off of light fermions [47], including a detailed model of how this component of the shower turns off due to W/Z mass effects. A study using SHERPA [48] instead breaks down these emissions into separate transverse (V T ) and longitudinal (V L ) components, coupling in the latter strictly using Yukawa couplings by appealing to the Goldstone-boson Equivalence Theorem (GET) [21,58]. The problem has been approached in different way within ALPGEN [46,59], by multiplying exclusive hard event rates with the fixed-order Sudakov factors of [34,35] and supplementing with exact fixed-order real emission processes. This approach, which is itself a first step towards electroweak shower matching, works well when the soft/collinear phase space enhancements are modest and the need for added accuracy of higher-multiplicity hard event generation balances the added computational complexity. However, a complete matching prescription will also ultimately involve a dedicated parton shower step, especially when convolved with QCD radiation. The simpler, process-independent parton shower approach will also become particularly useful in new physics applications [60,61].

Our approach
Notably, no existing general-purpose parton showering algorithm that is capable of generating fully exclusive events has addressed the full scope of universal collinear electroweak physics. In particular, a complete treatment must include the high-rate of non-Abelian splittings amongst the weak bosons themselves, as well as showers that involve longitudinal/scalar states and many of the sometimes subtle effects of spontaneous symmetry breaking. The goal of the present paper is to outline such an algorithm, providing a comprehensive framework in which all collinear electroweak showering phenomena can be implemented, and including a systematic treatment of EWSB. Towards this end, we derive and tabulate the complete set of electroweak splitting functions in the broken phase, including the massive fermions, gauge bosons, and the Higgs boson. These generalize and unify both the unbroken-phase evolution equations of [43] and the purely broken-phase effects already observed within the effective-W approximation, namely the generation of longitudinal vector boson beams from massless fermions [19][20][21]. We further investigate some of the physical consequences of these various electroweak showering phenomena.
Relative to QED and QCD showers, the complete electroweak parton shower exhibits many novel features. At the level of the unbroken theory at high energies, the shower becomes chiral and the particle content is extended to include an EW-charged scalar doublet. Most of the degrees of freedom contained in this scalar are to be identified with the longitudinal gauge bosons via the Goldstone-boson Equivalence Theorem. Including Yukawa couplings, the set of core splitting function topologies expands from the usual three to seven. EWSB also already makes a subtle imprint here due to the presence of a preferred isospin basis for asymptotic states, leading to interference and self-averaging effects between different exclusive isospin channels. The latter are intimately related to "Bloch-Nordsieck violation" when occurring in the initial state [41,45,62]. As the shower evolves down through the weak scale, it becomes physically regulated by the appearance of gauge boson, scalar, and fermion masses. Unlike in QCD where the shower regulation occurs nonperturbatively due to confinement, or in QED where a small photon mass is sometimes used as an artificial regulator for soft emissions, the electroweak shower exhibits a perturbative transition with genuinely massive gauge bosons. It is possible to describe this transition rather accurately, but doing so requires a careful accounting of symmetry-violating effects beyond simple kinematic suppressions, and a consistent elimination of gauge artifacts. In particular, Goldstone-boson equivalence ceases to hold at relative transverse momenta of order the weak scale, allowing for an additional burst of many "ultra-collinear" radiation processes that do not exist in the unbroken theory, and are highly suppressed at energy scales k T ≫ v. To cleanly isolate these effects, we introduce a novel gauge dubbed "Goldstone Equivalence Gauge" (GEG). This is a particularly convenient choice of non-covariant gauge, allowing a completely transparent view of Goldstone-boson equivalence within the shower, as well as systematic corrections away from it in the splitting matrix elements, organized in a power series in VEV factors. The naively bad high energy behavior of the longitudinal gauge bosons is deleted, and the Goldstone fields allowed to interpolate physical states, at the cost of re-introducing explicit gauge-Goldstone boson mixing.
Our formalism developed here has deep implications and rich applications at TeV-scale energies and beyond. Some aspects include EW parton distribution functions associated with initial state radiation (ISR), multiple emissions in EW final state radiation (FSR), consistent merging of EW decays with EW showering, a quantum-coherent treatment of the Sudakov evolution of γ/Z/h states, as well as modeling of general ultra-collinear processes including, e.g., t R → ht R and h → hh. We also make some preliminary studies of the impact of EW showering on new physics searches in the context of a heavy W ′ decay. Quite generally, we begin to see the emergence of the many nontrivial phenomena of "weak jets" across a broad range of SM and BSM phenomena.
Before proceeding, we also clarify what is not covered in our current treatment. We make exclusive use here of the collinear approximation, which, in physical gauges such as GEG, explicitly factorizes all soft and collinear divergences particle-by-particle, isolating them to 1 → 2 real emission diagrams and self-energy loops. This furnishes a formally leading-log model of EW showering, capturing all double-log effects from the soft-collinear region of gauge emissions, as well as the single-logs associated to all hard-collinear splittings. The former are identical to the double-logs that would be inferred from the collinear limits of the eikonal approximation, whose particle-by-particle factorization can be seen upon application of Ward identities [34,35,45]. However, there are additional single-log soft divergences within gauge emission interferences and virtual exchanges between different particles, which do not factorize so simply. For non-singlet EW ensembles, these contributions lead to global entanglements of isospin quantum numbers between different particles in the event, which are absent in our shower. These isospin entanglements are somewhat analogous to the global kinematic entanglements that occur due to soft gluon emissions/exchanges at NLL level in QCD. Nonetheless, the dominant effects of isospin rearrangements, in particular the Bloch-Nordsieck violation, arise already at the double-log level, and are modeled by our shower up to residual single-log ambiguities. We will address approaches to the NLL resummation of isospin entanglements in a future work [63].
The rest of the paper is organized as follows. We begin in Section 2 with a generic discussion of splitting and evolution formalism with massive particles. We then outline some of the other nontrivial features such as PDFs for massive particles, interference between different mass eigenstates, showers interpolating onto resonances, and back-reaction effects from multiple emissions. In Section 3, we introduce the splitting kernels for the unbroken electroweak theory, namely SU (2) L × U (1) Y gauge theory with massless fermions in SM representations, a single (massless) scalar doublet, and Yukawa interactions. We then proceed in Section 4 to generalize these results to the broken phase. After a discussion of the violation of the Goldstone-boson Equivalence Theorem, we introduce the Goldstone Equivalence Gauge. We then discuss the EWSB modifications to the unbroken splitting functions and present a complete list of ultra-collinear processes that arise at leading-order in the VEV. Section 5 explores some key consequences of electroweak showering in finalstate and initial-state splitting processes, including a discussion of EW parton distribution functions and multiple EW final state radiation. We emphasize the novel features of the EW shower and illustrate some of the effects in the decay of a heavy vector boson W ′ . We summarize and conclude in Section 6. Appendices give supplementary details of Goldstone Equivalence Gauge, the corresponding Feynman rules and illustrative examples of practical calculations, more details on the density-matrix formalism for coherent Sudakov evolution, and a short description of our virtuality-ordered showering program used for obtaining numerical FSR results.

Showering Preliminaries and Novel Features with EWSB
We first summarize the general formalism for the splitting functions and evolution equations with massive particles that forms the basis for the rest of the presentation. We then lay out some other novel features due to EWSB.

Splitting formalism
Consider a generic hard process nominally containing a particle A in the final state, slightly off-shell and subsequently splitting to B and C, as depicted in Fig. 1 (left figure). In the limit where the daughters B and C are both approximately collinear to the parent particle A, the cross section can be expressed in a factorized form [2] dσ X,BC ≃ dσ X,A × dP A→B+C , where dP is the differential splitting function (or probability distribution) for A → B + C. A given splitting can also act as the "hard" process for later splittings, building up jets. The factorization of collinear splittings applies similarly for initial-state particles, leading to the picture of parton distribution functions (PDFs) for an initial state parton B or C, as in Fig. 1 (right figure), We will discuss this situation in the next subsection. While our main focus here is on the leading-log resummation of these splitting effects in a parton shower/evolution framework, at a leading approximation Eqs. (2.1) and (2.2) can also be taken as-is, with a unique splitting in the event and no virtual/resummation effects, in order to quickly capture the tree-level collinear behavior of high energy processes. In our further analyses, we will refer to such a treatment as a "fixed-order EW shower" or "fixed-order EW FSR (ISR)." Integrating out the azimuthal orientation of the B + C system, the splitting kinematics are parametrized with two variables: a dimensionful scale (usually chosen to be approximately collinear boost-invariant) and a dimensionless energy-sharing variable z. Common choices for the dimensionful variable are the daughter transverse momentum k T relative to the splitting axis, the virtuality Q of the off-shell particle in the process, and variations proportional to the daughters' energy-weighted opening angle θE A . Our descriptions here will mainly use k T , as this makes more obvious the collinear phase space effects in the presence of masses. For our numerical results in Section 5, we switch to virtuality, which allows for a simpler matching onto W/Z/t decays. Mapping between any of these different scale choices is however straightforward. The energy-sharing variable z (z ≡ 1 − z) is commonly taken to be the energy fraction of A taken up by B (C). The splitting kinematics takes the form When considering splittings involving massive or highly off-shell particles, various possible definitions of z exist which exhibit different non-relativistic limits. Besides strict energy fraction, a common choice is the light-cone momentum fraction, Our specific implementation in Section 5 uses the three-momentum fraction which makes phase space suppression in the non-relativistic limit more transparent. However, in the relativistic regime, where the collinear factorization is strictly valid, all of these definitions are equivalent, and we do not presently make a further distinction. 1 In the simplest cases, generalizing the collinear splitting function calculations to account for masses is straightforward. Up to the non-universal and convention-dependent factors that come into play in the non-relativistic/non-collinear limits, the splitting functions can be expressed as (2.5) Here, M (split) is the A → B + C splitting matrix-element, which can be computed from the corresponding amputated 1 → 2 Feynman diagrams with on-shell polarization vectors (modulo gauge ambiguities, which we discuss later). This may or may not be spin-averaged, depending on how much information is to be kept in the shower. Depending upon the kinematics, the mass-dependent factors in the denominator act to either effectively cut off collinear divergences at small k T or, in final-state showers, to possibly transition the 1 There is unavoidably some frame-dependence to this setup, as there is in all parton showers that are defined strictly using collinear approximations. A more complete treatment would exhibit manifest Lorentzinvariance and control of the low-momentum region, at the expense of more complicated book-keeping of the global event's kinematic and isospin structure, by using superpositions of different 2 → 3 dipole splittings. Extending our treatment in this manner is in principle straightforward, but beyond the scope of the present work.
system into a resonance region. In cases where interference between different mass eigenstates can be important, this basic framework must be further generalized. Resonance and interference effects are introduced in Section 2.3. On dimensional grounds, |M (split) | 2 goes like either k 2 T or some combination of the various m 2 's. Conventional splitting functions typically scale like dk 2 T /k 2 T , which is exhibited by all of the gauge and Yukawa splittings of the massless unbroken electroweak theory, as to be shown in Section 3. There can also be mass-dependent splitting matrix elements that lead to m 2 dk 2 T /k 4 T type scaling. These splittings are highly suppressed for k T ∼ > m. However, they are much more strongly power-enhanced at low k T , a behavior which we call ultra-collinear. Upon integration over k T , the total rate for an ultra-collinear splitting comes out proportional to dimensionless combinations of couplings and masses, with the vast majority of the rate concentrated near k T ∼ m. Such processes exist in familiar contexts like QED and QCD with massive fermions, for example the helicity-flipping splittings e L → γe R and g → b LbL . They are usually not treated as distinct collinear physics with their own universal splitting functions, though they are crucial for systematically modeling shower thresholds. We choose to treat them on independent footing, since the threshold behaviors of the electroweak shower are highly nontrivial, including processes that are qualitatively different from the massless limit.
In both the conventional collinear and ultra-collinear cases, the remaining z dependence after integrating over k T can be either dz/z or dz×(regular). The former yields additional soft logarithms (again, formally regulated by the particle masses), and appears only in splittings where B or C is a gauge boson.

Evolution equations
When applied to the initial state, the splitting functions outlined in the previous section lead to both initial state radiation (ISR) as well as the dynamical generation of B and C parton distribution functions from a parent A. Considering a generic parton distribution function f i (z, µ 2 ) with a factorization scale µ in k T -space, the leading-order convolution relation is where µ 0 is an input factorization scale. Differentiating with respect to µ 2 and incorporating as well the evolution of the f A leads to the celebrated Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) equation [64][65][66].
Gauge theories such as QED and QCD predict that at high energies the splitting functions dP/dk 2 T go like 1/k 2 T , and thus that the PDFs evolve like ln(Q 2 /µ 2 ). This is the classic violation of the Bjorken scaling law [67]. In the broken electroweak theory, there are also the qualitatively different ultra-collinear splitting functions, which instead go as m 2 /k 4 T . The PDFs arising from these splittings "live" only at the scale k T ∼ m. Instead of evolving logarithmically, they are cut off by a strong power-law suppression at k T ∼ > m. The corresponding PDFs preserve Bjorken scaling, up to contributions beyond leading order. In particular, longitudinal weak boson PDFs are practically entirely determined at splitting scales of O(m W ), even when used as inputs into processes at energies E ≫ m W . 2 Numerical computation of electroweak PDFs with a proper scale evolution do not exist yet in the literature, though the complete unbroken-theory evolution equations appear in [43], and fixed-order results are straightforward to obtain with the simple convolution in Eq. (2.6). In the resummed treatment, contributions from the region k T ∼ m W can perhaps most simply be incorporated as perturbative "threshold" effects, essentially adding in their integrated fixed-order contributions up to some scale (a few)×m W as δ-functions in k T -space. These would include the finite, mass-suppressed contributions from the turnon of f → W T f splittings, as well as the entire ultra-collinear f → W L f contribution. Equivalently at leading-order, they may instead be folded continuously into the DGLAP evolution using the massive splitting functions defined as in Eq. (2.5). This latter approach may also be simpler when alternative scaling variables are used, such as virtuality.
The other qualitatively new electroweak effects in the PDFs concern the treatment of weak isospin. First, the chiral nature of the EW gauge interactions leads to more rapid evolution toward low-x for left-handed fermions than for right-handed fermions. Furthermore, the isospin non-singlet nature of typical beam particles yields an additional interesting subtlety. In QED and color-averaged QCD evolution, the soft-singular limits of, e.g., q → gq at a given scale become indistinguishable from q → q with no splitting. Indeed, this allows for the balancing of real and virtual IR divergences as z is formally taken to zero at fixed k T , conventionally encoded in the plus-prescription. However, following this prescription for the electroweak evolution of fermion PDFs at k T ≫ m W leads to unregulated divergences in isospin-flipping transitions, such as u L ↔ d L via arbitrarily soft W ± emission. This is a manifestation of the so-called Bloch-Nordsieck violation effect [41,45,62]. Regulation and resummation of this effect requires the introduction of some form of explicit cutoff z ∼ > k T /E in the evolution equations when formulated in (k T , z) space, in order to avoid non-collinear emission regions [43]. 3 The net effect is a gradual, controlled merging of the u L and d L PDFs (or e L and ν L PDFs in the case of electron beams) into a common "q L " ("ℓ L ") PDF. Unlike conventional PDF evolution, implementing the z cutoff in this way necessitates extending the arguments of the PDFs to explicitly include the (CM-frame) beam energy. While this is not a major complication, we do point out that different choices of 2 This observation persists even in the presence of QCD corrections. We can imagine that a quark is first evolved to large kT (and hence large space-like virtuality Q) from multiple gluon emissions, and then splits into an on-shell quark and space-like longitudinal vector boson. The former emerges as an ISR jet and the latter participates in a hard interaction. We would find (e.g., using Goldstone Equivalence Gauge, introduced in Section 4.1) that the collinear-enhanced piece of the scattering amplitude carries a net suppression factor of O(m 2 /Q 2 ), which cannot be compensated by integration over the collinear emission phase space. 3 In QED and QCD, these non-collinear emissions are implicitly and "incorrectly" integrated over in the plus-prescription. However, in the limit E ≫ kT , the numerical impact of doing so is of sub-leading importance.
scaling variables may yield the same non-collinear regulation without requiring the extra energy argument. A particularly simple choice would be the energy-weighted angle θE A . We defer a detailed study of these issues to future work [63]. We caution that this treatment of the initial state using PDFs remains strictly valid only within the leading-log, collinear approximation. Soft W ± virtual exchanges between the isospin non-singlet beams will induce single-log entanglements that do not factorize between the individual beams, and even more complicated entanglements emerge when we also consider isospin-exclusive final states. The proper generalization for the initial state is from running PDFs to running quantum-ensemble parton luminosities defined for pairs of beams. But it is also possible to define a scheme where these beam-entanglement effects are selectively treated at fixed-order, and PDF resummation still suffices [63]. (The entanglement effects actually wash out as the scale is raised and the isospin ensembles become incoherent.) However, these PDFs will still likely reference the global beam setup via the aforementioned non-collinear cutoff.
Even applying the conventional factorization at leading-log, some of the PDFs must also still be treated as matrices [43]. This is particularly relevant for the photon and transverse Z-boson PDFs, which develop sizable off-diagonal contributions. Indeed, the naive concept of independent "photon PDF" and "Z PDF" at k T ≫ m Z is necessarily missing important physics, as γ and Z are not gauge eigenstates. We outline the appropriate treatment in Section 2.3.2 and Appendix C.
The same splitting functions that govern ISR and PDF generation also serve as the evolution kernels for final-state radiation (FSR). This integrates to the well-known Sudakov form factor ∆ A (t) characterizing the possible time-like branchings of parent A at scales below t ∼ log(k T ) or log(Q) where the allowed z range is determined by kinematics. Practically, we perform the evolution starting at a high k T or virtuality scale characterized by the CM-frame energy of the hard partonic process, and running continuously down through the weak scale with the proper mass effects. The Sudakov factor, evaluated in small t steps, functions as a survival probability for A, upon which the usual Markov chain monte carlo is constructed. (See, e.g., [68].) If A does not survive at some step, it is split into a state B + C. This splitting acts as the "hard" process that produced particles B and C, and Sudakov evolution is continued on each of those particles. The "resolution" scale t 0 can be any scale well below m W , at which conventional QED and QCD showers can take over. Of course, the basic framework leaves many details unspecified, and allows for a great deal of freedom in specific implementation. For example, besides the choice of evolution variable, one must also specify a treatment of kinematic reshuffling. We elaborate on some additional aspects of our own implementation of final-state showers below and in Appendix D. We will generally refer this treatment of Sudakov formalism as the "full EW shower" or "full EW FSR", in contrast to the fixed-order splitting calculations in Eqs. (2.1) and (2.2).

Other novel features in EW showering
There are several additional novel features in EW showering beyond those encountered in the standard formalism. We outline a few relevant to our later discussions and also propose concrete schemes for their implementations.

Mass effects
Besides the basic kinematic modifications and the emergence ultra-collinear splitting phenomena, the existence of a mass scale m W,Z ∼ gv and m f ∼ y f v requires some special treatments as we approach kinematic thresholds and the boundaries of turnoff regions. An immediate complication is that final-state weak showering smoothly connects onto the on-shell weak decays of top quarks, W/Z bosons, and (to a much lesser extent) Higgs bosons. The shower describes the highly off-shell behavior of these particles, including resummed logarithmically-enhanced effects. But the effect of the pole is nonetheless visible, encoded in the last term in the denominator of Eq. (2.5). Within the resonance region, the dominant behavior is more correctly captured by the standard Breit-Wigner line-shape governed by the physical width Γ, which involves a very different kind of resummation. However, a few Γ above the peak, both descriptions can be expanded perturbatively and yield numerically similar predictions. 4 It is therefore straightforward to define a wellbehaved matching prescription. This is easiest to formulate within a virtuality-ordered shower: Halt the shower at some matching scale Q match = m+(a few)Γ, and if the state has survived to this point, distribute its final mass according to a Breit-Wigner resonance below Q match . The exact choice of matching scale here is not crucial, as long as it is within the region where the Breit-Wigner and shower predictions are comparable. For other shower ordering variables, such as k T , we can instead run the shower down to its nominal kinematic limit, but not integrating z within the region that would yield Q < Q match . In either case, the parton shower may be restarted on the resonance's decay products.
Another place where mass effects can become important is in multiple emissions. In massless showers, sequential splittings are dominantly very strongly-ordered in scale, and as a consequence a given splitting rate can be computed without regard to the subsequent splittings while still capturing the leading behavior. However, in showers with massive particles, a large fraction of the available phase space for secondary splittings may require nontrivial kinematic rearrangements within the preceding splittings. For example, a W boson might nominally be produced with a kinematic mass m W via emission off of a fermion. If the W subsequently splits into a W and a Z boson at a virtuality Q ≫ m W , there is a chance that the off-shell W now sits near a suppressed region (i.e., dead cone) for emission off of the mother fermion. In order to avoid badly mis-modeling such cases, secondary splittings can be weighted according to the relative rate modification that would be incurred on the previous splitting. This back-reaction factor depends in detail on how kinematic arrangements are done in the shower. Generally, a given (z, Q) or (z, k T ) parametrizing the mother splitting will be mapped onto a new (z * , Q * ) or (z * , k * T ) for producing the off-shell daughter. The required back-reaction factor is the ratio of the new differential splitting function to the original one, multiplied by the Jacobian for the change of variables. For a final-state shower sequence A * → B * C → (DE)C, for the nested splitting we can use a splitting function multiplied by the back-reaction factor: The simplest implementation would compute this factor independently for each daughter branch, assuming an on-shell sister and neglecting possible correlations in the potentially fully off-shell final configuration A * → B * C * . But a more thoroughly correlated weighting scheme could be pursued if deemed numerically relevant. The above prescription also generalizes beyond massive showers, wherein it has a sizable overlap with the effects of standard angular vetoing. We further show below how back-reaction factors can be conveniently applied for a complete treatment of mixed neutral bosons, wherein an "on-shell" kinematic mass is not necessarily determined at their production. The above back-reaction effects can be particularly important for ultra-collinear emissions, as these occur almost exclusively at the boundaries delineated by finite-mass effects. For example, the prototypical ultra-collinear emission is f → W L f ′ with massless fermions [19][20][21]. It proceeds only via a delicate balancing between a suppression factor m 2 W /E 2 in the squared splitting matrix element and a strong 1/k 4 T power enhancement from the fermion propagator that gets cut off at k T ∼ m W , controlled by the form of the denominator in Eq. (2.5). Within a final-state shower, if either the W L or its sister f ′ is set far off-shell by a secondary splitting at some scale Q (possibly a QCD splitting), that cutoff moves out to k T ∼ Q but the original production matrix element stays approximately the same, and the total rate picks up an additional relative power suppression factor of O(m 2 W /Q 2 ). 5 Roughly speaking, ultra-collinear processes can only occur near the "end" of the weak parton shower as it passes through the weak scale, or conversely near the "beginning" of weak PDF evolution. Such a feature is essentially built-into k T -ordered parton evolution. The back-reaction correction ensures that it is also enforced in showers built on other ordering variables, such as virtuality, while still allowing further low-scale showering such as q → gq and W L → γW L .

Mixed-state evolution
Thus far, the shower formalism that we have presented neglects the possibility of interference between different off-shell intermediate particle states contributing to a specific splitting topology. Traditionally in QED and QCD showers, such interference leads to subleading effects associated with the unmeasured spin and color of intermediate particles [69]. However, the full electroweak theory at high energies presents us with cases where different mass and gauge eigenstates can also interfere at O(1) level, most notably the neutral boson admixtures γ/Z T and h/Z L [43]. All other particles in the SM carry (approximately) 5 When the WL is off-shell, we would naively compensate by using an off-shell gauge polarization, yielding However, the appropriate treatment, discussed in more detail in Appendices A and B, uses on-shell polarization factors throughout. Additional non-collinear corrections might still be present, but are more appropriately viewed as contributions to 1 → 3 splittings. New soft logarithms might also arise in these processes, but new collinear logs will not. conserved charge or flavor quantum numbers that can flow out into the asymptotic state, and therefore they do not tend to interfere in this manner. Interferences originating from CKM/PMNS flavor violations should be small and difficult to observe, and we neglect them for simplicity.
Showering involving superpositions of different particle species can be described using density matrix formalism. Let us consider the simpler case of final-state showers for illustration. The initial value of the density matrix is set proportional to the outer product of production amplitudes: , tracing out over other details of the rest of the event. 6 Here, the indices run over the particle species. The probability for an initial mixed quantum state to subsequently split into a specific exclusive final state must be computed by generalizing the splitting functions to Hermitian splitting matrices dP ij . The exclusive splitting rates are then computed by tracing against the normalized density matrix, 7 . (2.10) Representing the propagator matrix as D ij , and the amputated splitting amplitudes as M (split) i , this modifies Eq. (2.5) to the more complete, yet more complicated form Note that large interference effects can persist even in the massless limit with unmixed propagators. A full treatment, including the Sudakov evolution for ρ ij and the explicit form of the propagators for γ/Z T and h/Z L systems, is given in Appendix C.
Handling the kinematics and decays of mixed states requires some additional steps. "On-shell" kinematics cannot be defined a priori, and we cannot collapse onto mass eigenstates or a showered final-state with well-defined mass until the coherent Sudakov evolution has run its course. A simple prescription is to first produce a mixed boson with its minimum possible kinematic mass (zero for γ/Z T , m Z for h/Z L ) in order to fully fill out the phase space. Splittings that occur before reaching the resonance are weighted by a back-reaction factor as per Eq. (2.9). If the state survives un-split down to the heavier resonance's matching threshold, we can decide to project onto a specific mass eigenstate according to the relative probabilities encoded in the surviving density matrix. The backreaction factor may once again be employed here, implemented as a veto probability for the heavier resonance. (The factor will typically come out less than one for a sensibly-defined change of variables.) If the veto is thrown, the splitting that produced the mixed state is undone, and its mother's evolution continued. This prescription especially becomes relevant when evolving near kinematic thresholds or suppressed regions, for example where Z boson emission would be suppressed but photon emission allowed.
For the mixed γ/Z T system, if a photon is projected out, we can restart a pure QED parton shower (γ → ff ) with virtuality constrained below the Z boson's Q match scale at ≈ 100 GeV. Interference effects below the matching scale can also be incorporated by coherently adding both the γ and Z contributions within the Z resonance region. This requires delineating as well a lower virtuality boundary, ideally at a scale O(1) smaller than m Z . Depending on the integrated probability in this region (modulo the back-reaction veto), we would either create an ff state with an appropriately-distributed mass, or again set the state to a photon and continue running a pure QED shower, now constrained below the Z resonance region.
We also comment that a fully consistent treatment here would require minor changes to the standard output formats of hard event generators. The standard practice of immediately collapsing onto mass eigenstates is equivalent to assuming trivial Sudakov evolution, and cannot formally be inverted such that a proper coherent parton shower can be applied. In particular, only one specific linear combination of γ/Z T states participates in the highrate non-Abelian splittings to W ± T W ∓ T . While collapsing onto mass eigenstates is required to obtain well-defined hard event kinematics, a simple remedy here would be to supply for these particles their production density matrices, using some appropriately-mapped massless kinematics.

Splitting Functions in Unbroken
Before working out the complete set of electroweak splitting functions in the broken phase, it is important to first consider a conceptual limit with an unbroken SU (2) L × U (1) Y gauge symmetry with massless gauge bosons and fermions, supplemented by a massless complex scalar doublet field H without a VEV. This last ingredient is the would-be Higgs doublet. This simplified treatment in the unbroken phase is not only useful to develop some intuition, but also captures the leading high-k T collinear splitting behavior of the broken SM electroweak sector. Some aspects of electroweak collinear splitting and evolution at this level have been discussed, e.g., in [43].
Anticipating electroweak symmetry breaking, we adopt the electric charge basis in weak isospin space. The corresponding SU (2) L bosons are W ± and W 0 , and the hypercharge gauge boson we denote as B 0 . Gauge boson helicities are purely transverse (T ), and are averaged. 8 For the scalar doublet, we decompose as where φ ± , φ 0 will later become the electroweak Goldstone bosons and h the Higgs boson. However, at this stage, we will keep the neutral bosons h and φ 0 bundled into the complex scalar field H 0 , as they are produced and showered together coherently. In the absence of the VEV, the doublet carries a perturbatively-conserved "Higgs number," which may also be taken to flow through RH-chiral fermions in the Yukawa interactions. 9 We denote a generic fermion of a given helicity by f s with s = L, R (or equivalently s = ∓). We do not always specify the explicit isospin components of f at this stage, but implicitly work in the usual (u, d)/(ν, e) basis. Isospin-flips (including RH-chiral isospin where appropriate) will be indicated by a prime, e.g. u ′ = d. Effects of flavor mixing are ignored. The U (1) Y and SU (2) L gauge couplings are respectively taken to be g 1 ≈ 0.36 and g 2 ≈ 0.65 (here evaluated near the weak scale, though in general run to a scale of O(k T )). For compactness we often represent a generic gauge coupling by g V . We represent the gauge charge Q of a particle p coupling to gauge boson V by Q V p , and we give the complete list of the gauge charges for the SM fermions and scalars in Table 8 in Appendix B.1.
The splitting functions that involve only fermions and gauge bosons closely follow those of QED and QCD. Fermions with appropriate quantum numbers may emit transverse SU (2) L and U (1) Y gauge bosons with both soft and collinear enhancements, yielding total rates that grow double-logarithmically with energy. At this stage, fermion helicity coincides with the corresponding chirality, and is strictly conserved in these processes. The SU (2) L bosons also couple to one another via their non-Abelian gauge interactions, and similarly undergo double-logarithmic soft and collinear splittings W 0 → W + W − and W ± → W ± W 0 . This is in direct analogy to g → gg in QCD, except that here we do not sum/average over gauge indices. All of the electroweak gauge bosons may also undergo single-log collinear splittings into fermion pairs, similar to g → qq or γ → ff .
The results can be cast into a familiar form. We write the probability function of finding a parton B inside a parton A with an energy or momentum fraction z in terms of the collinear splitting kernels for A → B as P BA (z). Stripping the common g 2 /8π 2 and 1/k 2 T factors, as well as group theory factors that depend on the gauge representations (hyper-charges or SU (2) L quadratic Casimirs and Dynkin indices), we are left with (1 + z 2 )/z, but it is not independent and can be derived from P V f with z ↔z. The factor of 1/2 in P f V , relative to the standard form in QED with the electric charge stripped (or in QCD with the SU (3) Dynkin index stripped), is due to the fact that we treat each chiral fermion individually. Interference between different gauge groups is a subtlety that is absent in the coloraveraged SU (3) QCD × U (1) EM shower, and arises here from the fact that we have fixed a 9 We have expanded the neutral scalar field as H 0 ∝ h − iφ 0 , adopting a phase convention such that h and φ 0 fields create/annihilate their respective one-particle states with trivial phases, and H 0 annihilates the one-particle state |H 0 ∝ |h + i|φ 0 . Treating h and φ 0 as independent showering particles would be analogous to adopting a Majorana basis instead of a Dirac basis for the fermions in QED or QCD. An incoherent parton shower set up in such a basis would not properly model the flow of fermion number and electric charge. Analogously, H 0 and H 0 * particles carry well-defined Higgs number that we choose to explicitly track through the shower. This leads to correlations between spins and electric charges within asymptotic states. preferred gauge basis for asymptotic states instead of summing over gauge indices. Within different exclusive isospin channels in this basis, exchanges of B 0 and W 0 can exhibit O(1) interference, and thus must be described using density matrices, which have briefly been discussed in Section 2.3.2. In a truly massless theory, the physical preparation and identification of states in any preferred weak isospin basis is actually impossible, since arbitrarily soft W ± can be radiated copiously at no energy cost and randomize the isospin. 10 Our preferred basis here only becomes physical once we turn on the electroweak VEV and cut off the IR divergences. But the tendency for states to self-average in isospin space will persist at high energies.
Beyond these, the major change is the introduction of the scalar doublet. 11 First, the scalars may themselves radiate SU (2) L and U (1) Y gauge bosons. The soft-collinear behavior is identical to their fermionic counterparts, but the hard-collinear behavior is different. Second, the electroweak gauge bosons can split into a pair of scalars, again in close analog with splittings to fermion pairs. Third, fermions with appreciable Yukawa couplings to the scalar doublet can emit a scalar and undergo a helicity flip. Finally, the scalars can split into a pair of collinear, opposite-chirality (same-helicity) fermions. The corresponding splitting function kernels are found to be -s fs (z) =z/2, derived from P V H and P Hf , respectively. 12 The splittings W 0 /B 0 → H 0 H 0( * ) can also be conveniently represented by the final-state hφ 0 , in what will ultimately become hZ L in mass/CP basis. Here the final-state bosons are entangled, but the effects of that entanglement are subtle and only become relevant if both bosons undergo secondary splittings and/or hard interactions. In practice, we will simply take the expedient of collapsing the final state to hφ 0 .
The complete set of splitting functions is summarized in Tables 1 through 3. The tables are organized according to the spin of the incoming particles: polarized fermions with helicity s, transverse gauge bosons (V T ), and scalars. Each table is further subdivided according to the spins of outgoing particles, all together corresponding to seven unique core splitting functions. The various table entries associated to a specific set of incoming and outgoing spins provide the remaining coupling and group theory factors. All of the splitting functions have a conventional collinear logarithmic enhancement dk 2 T /k 2 T , and those involving emission of a massless gauge boson have an additional soft logarithmic 10 Absent the quark chiral condensate at O(100 MeV), massless SU (2)L would also technically confine in the IR, so that asymptotic states would anyway be isospin-singlet bound states, making the situation even more analogous to QCD. 11 We neglect all 1 → 3 splittings coming from either the scalar quartic or the scalar-gauge 4-point. These may feature single-logarithmic collinear divergences, but are expected to be rather highly numerically suppressed due to an additional O(1/16π 2 ) phase space factor. 12 Note that transitions involving the scalars must conserve the Higgs number introduced earlier in this section. For example, we may have Table 1: Chiral fermion splitting functions dP/dz dk 2 T in the massless limit, with z (z ≡ 1 − z) labeling the energy fraction of the first (second) produced particle. The fermion helicity is labelled by s. Double-arrows in Feynman diagrams indicate example fermion helicity directions. Prime indicates isospin partner (u ′ s = d s , etc, independent of s). Yukawa couplings are labelled by the participating RH-helicity fermion. The state H 0 * is the "anti-H 0 ", produced when the RH fermion is down-type and in the initial-state, or up-type in the final-state. Processes with B 0 and W 0 implicitly represent the respective diagonal terms in the neutral gauge boson's density matrix, whereas [BW ] 0 indicates either of the off-diagonal terms (see text). Anti-fermion splittings are obtained by CP conjugation. The conventions for the couplings are given in B.1. Table 2: Transverse vector boson splitting functions dP/dz dk 2 T in the massless limit, where allowed by electric charge flow. N f is a color multiplicity factor (N f = 1 for leptons, N f = 3 for quarks). Other conventions as in Table 1. Table 3: Scalar splitting functions dP/dz dk 2 T in the massless limit via gauge couplings and Yukawa couplings. The symbol H in the column headings represents the appropriate state φ + , H 0 for the given splitting, and H ′ represents the SU (2) L isospin partner (e.g., H 0′ = φ + ). Anti-particle splittings are obtained by CP conjugation. Other conventions as in Tables 1 and 2.  Figure 2: Fixed-order differential emission rate for W ± bosons off a massless fermion at E f = 10 TeV: (a) k T distribution at z = 0.2, (b) z distribution at k T = m W /2. The different curves correspond to massless transversely-polarized W ± T (dotted curves), massive transversely-polarized W ± T (solid curves), and massive longitudinally-polarized W ± L (dashed curves). enhancement dz/z. (The latter are the only emissions that preserve the leading particle's helicity in the soft emission limit.) To represent the off-diagonal terms for the neutral gauge bosons (either in production or splitting, where appropriate), we use the symbol [BW ] 0 . Otherwise, processes involving B 0 or W 0 alone implicitly represent the respective diagonal term in the density matrix.

Splitting Functions in Spontaneously Broken
While the parton shower formalism of the electroweak theory in the symmetric phase has much in common with that of SU (3) QCD × U (1) EM , care needs to be taken when dealing with the broken phase and systematically accounting for the effects of the VEV (v). In a sense, we must extract the "higher-twist" effects of the broken electroweak theory in terms of powers of v/E. Although the regulating role of v in the shower is somewhat analogous to that of Λ QCD , the electroweak theory remains perturbative at v, and the unbroken QED shower continues into the deep infrared regime. The interplay between gauge and Goldstone degrees of freedom within the shower can also seem obscure, both technically and conceptually.
Most immediately, the splitting functions of the unbroken theory, already detailed in Section 3, must be adjusted to account for the physical masses of the gauge bosons, Higgs boson, and top quark. To large extent, these constitute simple modifications, folding in the kinematic effects discussed in Section 2. As a straightforward example, in Fig. 2 we illustrate the fixed-order emission rate for W ± bosons off a massless fermion at E f = 10 TeV. Both the collinear and soft singularities of the massless theory (dotted curves) become regulated with m W ≈ 80 GeV (solid curves), as seen in the transversely-polarized boson k T distribution in Fig. 2(a) and the z distribution in Fig. 2(b). 13 Indeed, giving the gauge bosons a mass is a common trick for regulating QCD and QED calculations. In the electroweak theory, such regulated splitting functions become physically meaningful. Figure 2 also shows a contribution from longitudinal gauge boson radiation off of a massless fermion (dashed curves). This is a good example of an "ultra-collinear" process which emerges after EWSB at leading power in v/E. In this case it has a splitting probability of the form The rate is seen to be significant in the region k T ∼ m W , and it can be larger than the conventional transverse emissions in the ultra-collinear region k T ∼ < m W as seen in Fig. 2(a). We further show in Fig. 2(b) the z distribution at k T = m W /2, where we can see the dominance of the longitudinal polarization (dashed curve) over the transverse polarization (solid curve) for all values of z at weak-scale values of k T . Here we have defined z as three-momentum fraction, employed a strict kinematic cut-off z > k T /E, and multiplied the splitting rate by the W velocity to account for non-relativistic phase space suppression.
Considering emissions from light initial-state fermions, the ultra-collinear origins of these longitudinal weak bosons leads to quite distinctive PDFs [19][20][21]. Due to the existence of an explicit mass scale m W ∼ gv, the resulting PDFs exhibit Bjorken scaling [67]. In other words, they do not run logarithmically and do not exhibit the usual scaling violations of conventional PDFs in massless gauge theories. Consequently, the ISR jets associated with their generation are constrained to the region k T ∼ m W even for arbitrarily-energetic hard processes. This observation has led to the concepts of "forward-jet tagging" [70][71][72] for the W L W L scattering signal and "central-jet vetoing" [73] for separating the f → W T f ′ backgrounds.
Such processes have no analogs in the unbroken theory. A naive application of the Goldstone-boson Equivalence Theorem (GET) [21,58] would have instructed us to identify longitudinal vector bosons with the eaten scalars from the Higgs doublet, and would have predicted zero rate because massless fermions have vanishing Yukawa couplings. More generally, we expect to see a variety of large effects of EWSB at k T ∼ v, beyond simple regulation of the unbroken-theory splitting functions. These will involve not only the broken-phase masses of the SM particles, but also broken-phase interactions such as scalarvector-vector and the scalar cubics.
The more general role of Goldstone boson equivalence and its violations within the parton shower are rather subtle. We expect that the high-k T showering of longitudinal gauge bosons should closely follow the behavior of the scalars in the unbroken theory. But even this simple identification is obscured by longitudinal polarizations that diverge with 13 Note that in the region z ∼ < mW /E, the W s are non-relativistic, and collinear splitting function language ceases to be strictly appropriate or reliable. This region could more rigorously be matched onto universal soft Eikonal factors, e.g. as in [34,35]. But in practice, our treatment here still yields approximately correct rates for splitting angles ∼ < 1 when the splitting is defined in the hard scatter frame.
energy and by the gauge/Goldstone boson propagators with gauge-dependent tensor and pole structure. For processes with multiple emissions, as well as with the introduction of the novel ultra-collinear emissions, complete isolation and removal of non-collinear gauge artifacts can appear rather complicated. We are thus compelled to seek out a more efficient treatment, such that the bad high energy behavior of the longitudinal gauge bosons is alleviated and the key features of EWSB are made more transparent.

Longitudinal gauge bosons and Goldstone Boson Equivalence
The standard form for the polarization vector of an on-shell longitudinal gauge boson W with a four-momentum k µ where we define the light-like four-vector The second term in Eq. (4.2) is of the order m W /E W , which could seemingly be ignored at very high energies in accordance with the GET. However, there are caveats to this picture, and understanding how pseudo-scalars and longitudinal vector bosons behave as both external and intermediate states requires some care.
In the simplest approach, one would keep only the leading contribution, k µ W /m W . When contracted into scattering amplitudes, this piece effectively "scalarizes" the longitudinal vector boson, realizing the GET. This can often be seen at the level of individual Feynman diagrams. For example, in the decay of a heavy Higgs boson with m h ≫ 2m W , the vertex g m W hW µ W µ simply leads to a scalar interaction (m 2 In other cases, such as in couplings to fermion lines, the naively bad high-energy behavior ∝ E W /m W is fully cancelled thanks to Ward identities, up to possible chirality-flip effects that go like m f /E W . This reproduces the Yukawa couplings of the unbroken theory. When longitudinal and Goldstone bosons appear as off-shell intermediate states, it is also possible to show that neither the naively badly-behaved structure k µ k ν /m 2 W (in unitarity gauge) nor spurious gauge/Goldstone poles (in more general gauges) can lead to new collinear behavior at zeroth-order in the VEV. The unbroken shower emerges as expected as long as k T ≫ m W .
The major complication to the GET picture is that the naively sub-leading effects from EWSB can dominate in the relativistic ultra-collinear regime. Even if the k µ W /m W piece of an emitted gauge boson is removed by Ward identities, the O(m W /E W ) remainder of ǫ µ L (W ) can still receive a compensating ultra-collinear power-enhancement in the region k T ∼ m W . There may also be comparable EWSB contributions lurking within off-shell propagators, including as well the propagators of Higgs bosons and massive fermions.
Disentangling all EWSB effects in an ultra-collinear parton splitting can be accomplished by isolating and removing all parts of a 1 → 2 splitting amplitude that go like (Q 2 − m 2 )/m 2 W , where Q 2 and m 2 are respectively the squares of the four-momentum and pole mass of the off-shell particle in the splitting. Once multiplied by the propagators, such contributions are explicitly not collinear-enhanced, and would need to be combined with other non-collinear (and hence non-universal) diagrams from a hard process. Their extraction can generally be accomplished via manipulations between kinematic quantities, polarization vectors, and couplings. However, carrying out this extraction procedure process-by-process can be tedious, especially when multiple gauge bosons and/or nested collinear emissions are involved, and the effects of EWSB are often not immediately obvious. Within the gauge/Goldstone boson sector, we expect that the k µ W /m W piece of the longitudinal polarization vector must generally reproduce the Goldstone scalar couplings, whereas the effects of EWSB are captured by the remainder term in Eq. (4.2). A more convenient approach for tracking EWSB effects would be to keep the Goldstone scalar contributions manifest, and treat the remainder polarization as a separate entity.
We point out that such a division can be enforced by judicious gauge-fixing. We do so here via a novel gauge which we call Goldstone Equivalence Gauge (GEG). GEG is defined by generalizing off-shell the light-like four-vector n µ that appears in Eq. (4.2) and using it to perform the gauge-fixing in momentum-space. Taking W µ to represent any specific real gauge adjoint, with contraction of gauge indices left implicit, we adopt the gauge-fixing term (dropping here and below the "W " subscript on energy/momentum variables) Taking the ξ → 0 limit effectively introduces an infinite mass term for the gauge polarization associated with the collinear light-like directionn µ ≡ (1,k), aligned with the large components of relativistic momentum modes. This reduces the naive number of dynamical gauge degrees of freedom from four to three. The transverse modes (xy or helicity ±1) are as usual, except that they gain a mass term after spontaneous symmetry breaking. The remaining gauge degree of freedom "W n " explicitly mixes into the Goldstone boson, and becomes associated with exactly the remainder polarization in Eq. (4.2). GEG is essentially a hybrid of Coulomb gauge [74] and light-cone gauge [75], incorporating both the rotational-invariance of the former and the collinear boost-invariance of the latter, while isolating spurious gauge poles/discontinuities away from physical regions. 14 This approach can be contrasted with the more commonly-used R ξ gauges, in which individual splitting diagrams often exhibit unphysical gauge artifacts scaling as 1/v, Goldstone fields live purely off-shell, and Goldstone equivalence can become obscured.
Canonically normalizing such that the gauge remainder field W n interpolates a longitudinal boson state with unit amplitude at tree level, its interaction vertices carry the polarization factor The Goldstone field remains an integral part of the description here, but in a manner quite different from that in R ξ gauges. In particular, it interpolates onto the same external particle as the remainder gauge field. This particle, which may alternately be viewed as a "longitudinal gauge boson" or as a "Goldstone boson", takes on a kind of dual identity in interactions. Processes involving creation/annihilation of this particle are computed by coherently summing over Feynman diagrams interpolated by both remainder gauge fields and Goldstone fields. 15 More details and example calculations are presented in Appendices A and B. However, we can summarize here the key features of GEG that are relevant for parton shower physics: • Gauge artifacts proportional to E/m W are deleted from the description of the theory at the outset, and appear neither in external polarizations nor in propagators. Physical longitudinal gauge bosons are no longer interpolated by a gauge boson field W L and its associated O(E/m W ) polarization vector ǫ µ L , and no propagating component of the gauge field serves a proxy for the eaten Goldstone bosons in high-energy interactions via "scalarization." Instead, only a remainder gauge field W n may still interpolate longitudinal gauge bosons. But it does so via the suppressed O(m W /E) polarization vector ǫ µ n in Eq. (4.5).
• The high-energy equivalence between longitudinal gauge bosons and Goldstone bosons becomes trivially manifest at the level of individual Feynman diagrams. This is because the Goldstone fields behave almost identically as in the unbroken theory at high energies (v/E → 0). The equivalence extends off-shell, encountering neither the usual fake gauge nor Goldstone poles. All propagators exhibit the physical pole at m W or m Z with positive residue. This greatly simplifies the interpretation of an "almost on-shell" boson as an intermediate state in a shower.
• Departures from Goldstone boson equivalence become organized in a systematic power expansion in v/E factors. This allows general ultra-collinear splitting processes to be viewed as simple sums of well-behaved 1 → 2 Feynman diagrams. EWSB contributions in splitting matrix elements can come from remainder-longitudinal gauge insertions, fermion mass terms in spinor polarizations, and a small set of standard EWSB three-point vertices.
As a final remark of this section, we would like to point out that the GET has been shown to be valid including radiative corrections [78][79][80]. Given the close relation between the GET and GEG, we suspect that GEG should also be adequate in dealing with radiative corrections.

Modifications to unbroken-phase splitting functions
The unbroken-phase splitting functions governed by the gauge and Yukawa couplings given in Tables 1 to 3 of Sec. 3 are still valid for k T 's and virtualities far above the masses of all of 15 For a different but related approach, see [77]. the participating particles, provided we make the identification between pseudo-scalars and longitudinal gauge bosons in accordance with the GET. Indeed, in Goldstone Equivalence Gauge, this correspondence is completely transparent. The splitting matrix elements can be used largely unchanged as long as all of the particles are also relativistic, with corrections that typically scale as O(g 2 v 2 /E 2 ).
At k T 's and virtualities approaching the physical masses, EWSB causes these splitting functions to either smoothly shut off or to transition into resonance decays. The modifications are captured by the propagator and kinematic effects outlined in Section 2. In particular, the propagator modifications effectively rescale the unbroken-phase splitting functions of Tables 1-3 as Soft (1/z type) singularities also generally become regulated, though in the 1 → 2 collinear splitting function language this regulation is somewhat convention-dependent. For k T 's far above the physical masses, soft singularities are anyway constrained by kinematics: z,z ∼ > k T /E A . For lower k T 's, such that non-relativistic splitting momenta can be approached, the k T suppression also sufficiently regulates any soft-singular behavior. But additional soft phase space factors can also be applied to reduce artificial spikes in the differential splitting rates. Minimalistically, this involves the product of velocities of the outgoing products in final-state showers, and for initial-state showers involves the product of the onshell daughter's velocity and the space-like daughter's "velocity". We have seen a simple example in Fig. 2(b).
For the neutral boson states, the propagator factors become matrices. These may be conveniently diagonalized by rotating from the interaction basis B 0 /W 0 and H 0 /H 0 * to the mass basis γ/Z T and h/Z L . The former requires the usual rotation by θ W in gauge space. The latter is accomplished by a U (2) rotation into the standard CP-eigenstates. The showering must still be performed coherently in order to capture nontrivial effects such as the flow of weak isospin and Higgs number. The full treatment is detailed in Appendix C. One residual complication is that the off-diagonal terms in the splitting function matrices are proportional to products of different propagator factors. E.g., for a γ/Z T state, the appropriate modification factor for dP γZ would use instead We also note that our convention here is to align the phases of external Z L states with those of the eaten scalar φ 0 . Consequently, terms like dP hZ L are pure imaginary. The above modifications do not explicitly address possible running effects in the masses. Indeed, the numerical impact of the mass terms in the shower is anyway highly suppressed except at splitting scales of O(v). Still, some cases, such as kinematics with k T ∼ v but Q ≫ v, might require special care in the inclusion of higher-order radiative corrections. Similar considerations apply to the purely ultra-collinear splitting processes discussed below.  6). The I V f symbol is a shorthand for the "charge" of a fermion in its Yukawa coupling to the eaten Goldstone boson, or equivalently the fermion's axial charge under the vector V . These are normalized to approximately follow the weak isospin couplings, but are defined independently of the fermion's helicity: Other conventions are given in Appendix B.

Ultra-collinear broken-phase splitting functions
The remaining task is to compute all of the ultra-collinear splitting functions, proportional to the EWSB scale like in Eq. (4.1). Generalizing the standard massless-fermion f → W L f ′ calculation [19][20][21], we include the splittings involving arbitrary particles in the SM. The electroweak VEV (v), to which all of these splitting functions are proportionate, has been explicitly extracted, as well as universal numerical factors, the kinematic factork 4 T as in Eq. We present these "purely broken" splitting functions in Tables 4−6, using similar logic as in Section 3, though now working exclusively in mass basis for the neutral bosons. Unlike conventional collinear splittings, ultra-collinear splittings do not lead to collinear logarithms. Instead, integrating the emissions at a fixed value of z yields a rate that asymptotes to a fixed value as the input energy increases. However, they are also unlike ordinary finite perturbative corrections, in that they are highly collinear-beamed, and subject to maximally large Sudakov effects from the conventional parton showering that can occur at higher emission scales.
Ultra-collinear emissions of longitudinal gauge bosons, when formed by replacing a transverse boson in any conventional gauge emission by a longitudinal boson, retain softsingular behavior ∼ 1/z. (Within GEG, the 1/z factors within the splitting matrix elements become regulated to 2E W /(E W + k W ).) Fully integrating over emission phase space, these still lead to single-logarithmic divergences at high energy. This result might seem at odds  Table 4 and in Appendix B. with smoothly taking the unbroken limit. For f → W L f ′ , as we dial v to zero at fixed fermion energy, the emission rate for longitudinal bosons grows unbounded. However, the spectrum of those bosons has a median energy fraction z ∼ m W /E f , and also tends to zero. Moreover, in theories where the fermion has a gauge-invariant mass, such as QED, the nominal ultra-collinear region k T ∼ < m W becomes subsumed by the usual emission dead cone at k T ∼ < m f .
Many of the other (soft-regular) splitting functions are close analogs of the unbroken splittings, but with "wrong" helicities. For example, there are processes where a fermion emits a transverse gauge boson but undergoes a helicity flip, and also where a fermion emits a Higgs boson without flipping its helicity. There are also new processes such as h → hh where such an identification is not possible. Schematically, all of these processes can be viewed as arising from 1 → 3 splittings in the unbroken theory, where one of the final-state particles is a Higgs boson set to its VEV.
To make Tables 4−6 more compact, and to make closer contact with practical applications, we have made one additional simplification by neglecting neutral boson interference effects for outgoing particles. E.g., for an ultra-collinear process such as t s → (h/Z L )t s (helicity non-flipping scalar emission), we treat the outgoing Higgs and longitudinal Z states incoherently. For final-state radiation, such a treatment is easily justified, since, as discussed in Section 2.3.1, the particles produced out of an ultra-collinear splitting have suppressed secondary showering. And for PDF evolution starting from an initial-state composed exclusively of light matter, there are simply no available ultra-collinear processes where such interference effects can occur (e.g., there is GET-violating q s → Z L q s , but not q s → hq s ). At higher scales, where heavier particles begin to populate the PDFs, further ultra-collinear splittings are again suppressed. Note, however, that we retain interference effects for incoming neutral bosons, which can remain important for final-state splittings like γ/Z T → W ± L W ∓ T . We also re-emphasize that interference effects for outgoing particles should still be retained for the conventional splitting functions, even in the broken phase. This is particularly important for the generation of the mixed γ/Z T PDF.

Shower Implementation and Related New Phenomena
We are now in a position to implement the splitting formalism and to present some initial physics results. Our studies here involving PDFs have been generated using simple numerical integration techniques. Our studies involving final-state radiation, which provide much more exclusive event information, have been generated using a dedicated virtuality-ordered weak showering code. Some technical aspects of this code can be found in Appendix D. We do not presently study the more technically-involved exclusive structure of weak ISR radiation. More detailed investigations of specific physics applications will appear in future work [63].
We first show some representative integrated splitting rates for an illustrative set of electroweak splitting processes in Table 7, at incoming energies of 1 and 10 TeV, as well as the leading-log asymptotic behavior. We have mainly focused on examples from Sections 3 and 4 that exhibit single-or double-logarithmic scaling with energy. Unless otherwise noted, the rates are summed/averaged over spins and particle species. (For instance, q = u L , u R , d L , d R , and f denotes all twelve fermion types of either spin.) The symbols in the parentheses denote the conventional collinear-enhanced (CL), infrared-enhanced (IR) and ultra-collinear (UC) behaviors, respectively. Radiation of a V T boson exhibits the usual CL+IR double-log behavior. Notably, the largest splitting rates occur for V T → V T V T , due to the large adjoint gauge charge. Splittings of this type occur with roughly 35% probability at 10 TeV, a factor that is enormous for an "EW correction" and which clearly indicates the need for shower resummation. We also see the analogous UC+IR process V T → V L V T , which only grows single-logarithmically, but which still represents a sizable fraction of the total splitting rate (even more so if we focus on low-k T regions, similar to Fig. 2). Similarly, the other ultra-collinear channels are smaller but not negligible.
We next present our numerical results for various exclusive splitting phenomena, paying special attention to the novelties that arise in the EW shower.  Table 7: Representative electroweak splitting behaviors and integrated fixed-order splitting probabilities for an illustrative set of processes at two parent energies E = 1, 10 TeV. The symbols in the parentheses denote the collinear (CL), infrared (IR), and ultra-collinear (UC) behaviors, respectively.

Weak boson PDFs
We first revisit the classic calculation of weak boson PDFs within proton beams [19,20]. The basic physical picture has been dramatically confirmed with the observation of the Higgs boson signal via vector boson fusion at the LHC [22]. It is anticipated that at energies in the multi-TeV regime, the total production cross section for a vector boson fusion process V 1 V 2 → X can be evaluated by convoluting the partonic production cross sections over the gauge boson PDFs, originated from the quark parton splittings q → W ± q ′ , q → γ/Zq. 16 A useful intermediate object in this calculation is the parton-parton luminosity, consisting of the convolutions of the PDFs from each proton. We write the cross section in terms of the parton luminosity of gauge boson collisions as and can approximate this luminosity at fixed-order using the concept of weak boson PDFs of individual quarks within the proton: 16 It should be noted that a formal factorization proof for electroweak processes in hadronic collisions is thus far lacking. For instance, it is not presently demonstrated whether contributions from gauge boson exchanges between the two incoming partons are factorizable. Nonetheless, we expect that the factorized PDF approach should furnish a reliable and useful calculation tool at very high energies at leading order, as indicated by simple scaling arguments [81,82].
Here, τ = s/S is the ratio of the partonic and hadronic energies squared, and τ low and τ high the kinematic boundaries (e.g., defining a bin in a histogram). We assume τ low ≫ 4m 2 W /S. The objects f V ∈q are evaluated at fixed-order as where the upper boundary of the k T integration is of order the partonic CM energy. For example [19,20], where the PDFs have been integrated up to k 2 T = s/4, assumed to be much larger than m W .
We emphasize that in deriving these illustrative fixed-order weak boson PDFs, we have not resummed the logarithmic enhancement, which remains explicit in Eq. (5.4) for the transverse bosons. There are also corresponding double-and single-log EW enhancements in the virtual corrections for the sourcing quarks, arising from integrating over both z and k T , which we have not accounted for. While these are of formally higher-order concern in determining the weak boson PDFs, they would also be required for an all-orders resummation of the leading-order effects. (We comment on other novel EW effects on the quark PDFs at the end of this subsection.) A related issue is that there are factorization scales implicit in the definition of the sourcing quark PDFs. Since the weak coupling and log(E/m W ) factors are together still below O(1) size at planned future machines, the choice of factorization scale might also seem to be of strictly higher-order concern. However, the interleaving of the much faster QCD evolution complicates the situation somewhat, especially at a large value of the energy fraction z. We have already noted above that the longitudinal W/Z PDFs would not continue to be sourced above m W , as their ultra-collinear generation is constrained to the region k T ∼ m W . It is therefore important to fix a factorization scale of O(m W ) for the quark PDFs from which the fixed-order W L PDFs are derived, even for processes where √ s ≫ m W [83]. However, the transverse W/Z PDFs are sourced continuously at all scales.
Higher-order calculations and/or full solution of the mixed QCD/EW DGLAP equations would be required to more fully resolve the issue of scale choices for the transverse bosons.
Here we simply fix the scale for the sourcing quark PDFs to be the geometric mean of √ s and m W (e.g., O(1 TeV) in a 10 TeV process). 17

Figures 3(a) and 3(b)
show the predicted fixed-order luminosities for a variety of possible colliding partons, including quarks as well as polarized W ± bosons and photons, at the 14 TeV LHC and a 100 TeV pp collider. At low scales, the "EW" PDFs are of course wholly dominated by photons. However, at scales above m W , the W ± PDFs are of comparable 17 This calculation uses only QCD evolution for the quark PDFs. The additional impact of electroweak evolution effects on the sourcing of the electroweak PDFs should indeed be small. Note also that mixed processes, such as VT VL → X would generally need a different factorization scale for each sourcing quark PDF. size. This can be seen here by comparing the qγ and qW ± T parton luminosities, as well as the W + T γ and W + T W − T luminosities. Note that in this comparison, we have also derived the photon PDF at fixed-order, sourced from quark PDFs. Attempts at fitting the photon PDFs with LHC data have recently been made [84]. Some recent discussions regarding the factorization scale uncertainties can be found in Ref. [85]. More importantly, a complete description will ultimately require including as well the Z T and mixed γ/Z T PDFs [63].
The PDFs and corresponding parton luminosities for longitudinal gauge bosons can be seen to be significantly smaller than those of transverse bosons. Of course, these nonetheless remain uniquely important for probing the nature of the electroweak sector beyond the Standard Model [21,58,73,[86][87][88]. In Fig. 3(c), we show the ratios of the partonic luminosities at the 100 TeV collider and the LHC dL 100 (s)/dL 14 (s). The increase with energy is largest for W L W L , with an enhancement factor about two orders of magnitude for √ s = 1-4 TeV.
As discussed in Sec. 2.2, some additional novel electroweak effects in the PDFs involve the different gauge interactions of left-handed and right-handed chiral fermions, and the isospin non-singlet nature of typical beam particles. The former leads to more rapid evolution to low-x for left-handed fermions than for right-handed fermions. The latter leads to Bloch-Nordsieck violation [41,45,62]. In PDF language, this appears as a self-correcting instability wherein the two LH isospin components of the beam flip between one another at a progressively increasing double-logarithmic rate, via soft/collinear W ± emissions. Both effects contribute to spontaneous beam polarization. In particular, in unpolarized proton beams the u L and d L PDFs will gradually split off from the u R and d R PDFs, and begin to asymptotically merge together into a common "q L " PDF at high energies. We investigate these phenomena in future work [63].

Final states with multiple gauge bosons
The collinear showering approximation allows us to estimate the leading contributions for multiple EW gauge boson production at high energies. A major component is splittings amongst the gauge bosons themselves via their non-Abelian interactions, in analogy with g → gg splittings in QCD. These have so far received little dedicated study in the electroweak theory within a parton shower framework. For some earlier studies of the fixedorder Sudakov effects in high-p T gauge boson production, see for example [28,29,31]. 18 As a simple illustration of the onset of shower-dominated behavior, we show in Fig. 4(a) a 2D kinematic distribution in fixed-order W ± Z + q/g production at a 100 TeV proton collider, generated with MadGraph5 [89]. A single kinematic cut p T (q/g) > 3 TeV is applied. The horizontal axis is the ∆R separation between the W and Z, and the vertical axis is the relative transverse momentum carried by the W : 2p T (W )/H T with H T defined as the scalar sum of all object p T s. Several features are immediately apparent. Most of the rate is concentrated along a curved band at low ∆R(W, Z), indicating W (q/g) production with a secondary collinear W → ZW splitting, and with enhancements at high (low) relative p T for W (Z) events. A second clear concentration of events occurs at ∆R(W, Z) ≃ π and nearmaximal relative H T indicating W q production with a secondary q → Zq splitting. A third, more subtle concentration is visible at ∆R(W, Z) ≃ π and low relative H T , representing Zq production with a secondary q → W q ′ splitting.
We can show how portions of this distribution arise within an available showering framework by generating V j events within PYTHIA8, and applying its native weak parton shower [47]. This shower currently includes only q → V q splittings, and does not model the 18 As a simple cross-check of our shower framework, we can make a comparison to the pT -dependent EW radiative corrections in W j production, as computed to NLO and approximate NNLO [31]. Since our shower is defined only for FSR, we study W q production and square the inferred Sudakov factor for the final-state quark. This approximately includes the Sudakov contribution of the initial-state quark. We select events without W/Z emissions, but allow final-state photons. At pT = 1 TeV, the EW correction to (NLO,NNLO) order is computed to be −(27, 24)%, whereas our resummed shower Sudakov also predicts −24%. At pT = 2 TeV, the EW correction to (NLO,NNLO) order is computed to be −(42, 34)%, whereas our resummed shower Sudakov predicts −33%.  Figure 4: Event population for exclusive W Z + j production in the plane of 2p T (W )/H T versus ∆R(W, Z) with p T (j) ≥ 3 TeV at a 100 TeV proton collider. (a) 2 → 3 fixedorder W Zj production generated with MadGraph; (b) 2 → 2 dressed with the PYTHIA weak shower, which includes only q → V q splittings; (c) 2 → 2 W j and Zj production dressed with fixed-order FSR splitting functions; (d) 2 → 2 dressed with the full EW FSR shower, including all collinear final-state Sudakov effects. QCD showering is not incorporated. An integrated luminosity of 10 ab −1 is used for illustration.
V → V V splittings responsible for the dominant rate near ∆R(W, Z) ≃ 0. The resulting incomplete distribution is shown in Fig. 4(b).
As a step toward gaining a more complete picture, we show in Fig. 4(c) the same distribution with hard V j events supplied by PYTHIA8 but dressed with our own EW FSR treatment (Appendix D), for the moment using fixed-order splitting functions and without Sudakov evolution effects. Now including V → V V as well as V → V q, the agreement  Figure 5: Normalized rates versus the number of multiple final-state W/Z emissions with a 10 TeV initial state particle, (a) d L -initiated showers for q → V q and V → V V splittings with full EW FSR (solid histogram), q → V q splitting only (long-dashed), and q → V q without back-reaction correction (short-dashed). Output from PYTHIA q → V q weak shower is also included for comparison (dotted histogram). (b) W T -initiated showers for fully constrained FSR (solid histogram), compared with various stages of approximations as labeled.
becomes quite good in all of the collinear-enhanced regions where we expect splitting functions to furnish a reliable description. 19 Besides the simpler generation of high-multiplicity final-states in collinear regions, the advantage of the parton shower is the ability to automatically fold in Sudakov corrections, going beyond fixed-order predictions. We show the result of running the full parton shower evolution Fig. 4(d), including as well important contributions such as V → ff. Exclusive W ± Z(q/g) events are selected as including exactly one each of "on-shell" W and Z, defined as lying within 10Γ of their pole mass, and we allow for multiple photon emissions. While the distribution looks similar to that at fixed-order, the overall rates in the collinear regions are reduced by several tens of percent due to the Sudakov corrections.
While formally any secondary parton splittings involve rate penalties of O(α W ), they become progressively more log-enhanced at high energies. This is again in close analogy to QCD. However, unlike in QCD, individual weak splittings in arbitrarily soft/collinear limits are in principle both observable and subject to perturbative modeling. Figure 5 shows the predicted number of W/Z generated from showering off a highly energetic particle with E = 10 TeV. In this calculation, we keep the weak bosons stable and include only the 19 Physics parameters here and in the MadGraph simulation are evaluated at a fixed scale of mZ for simplicity of comparison, using MadGraph's defaults. The PDF set is CTEQ6L1, evaluated at a factorization scale of 3 TeV. The PYTHIA simulation does not track fermion chirality throughout the hard event, and directly collapses γ/Z states into mass basis instead of providing a gauge-space wave function. We have explicitly corrected for both of these effects in this comparison and below.
splittings f → V f and V → V V . QCD showering is also turned off. We construct "weak jets" by clustering particles with the anti-k T algorithm [90] with R = π/2, and count the contained W/Z bosons. In Fig. 5(a), we show the results for a left-handed chiral fermion (d L ). Roughly speaking, we see that the emission of each additional gauge boson comes with an O(10%) suppression factor, which can be compared to the naive (not log-enhanced) O(1%) suppression typical of adding gauge bosons to lower-energy processes. The solid histogram shows the total rate and the long-dashed histogram indicates the rate with non-Abelian gauge splittings turned off. The difference indicates the large contribution from the gauge boson self-interaction beyond the first emission. As a cross-check, we include as well the prediction from the PYTHIA8 weak shower [47], as shown by the dotted histogram. Our own shower by default includes a back-reaction correction, discussed in Section 2.3.1, which approximates the expected suppression of multiple emissions due to dead cone-like effects for off-shell particles. To make a more direct comparison, we have also switched this off, and plotted the result as the short-dashed histogram. The two showers, both modeling unrestricted q → V q emissions, are then seen to be in close agreement.
In Fig. 5(b), we show the predicted number of W/Z contained in "weak jets" generated from showering off of a highly energetic transversely-polarized W ± boson with E W = 10 TeV. As already indicated in Table 7, the overall emission rates are much higher, close to 40% for the first emission (including both photons and Z bosons). Here we have again considered the effect of turning on/off back-reaction corrections. In addition, from experience with QCD showers, it is known that coherence effects in emission amplitudes lead to effective color-screening and approximate angular-ordering of nested emissions in non-Abelian splittings. To test this, we have also turned on/off a strict angular-ordering veto in our shower simulation. The results, visible in Fig. 5(b), are that both the backreaction correction and the angular ordering can have an O(1) effect at high multiplicities, but that the two effects come with sizable overlap. Splittings with large opening angles tend to exhibit large back-reaction effects, and vice-versa. This observation provides some evidence that modeling of the high-multiplicity region might be made to quickly converge, though more study is required.
It should be noted that at higher energy scales, the production of multiple gauge bosons could be the characteristic signature in many scenarios for physics beyond the SM [91,92].

EW Showers initiated by top quarks
Top quarks are instrumental in searches for new physics related to the EWSB sector, and for exotica such as resonances with large couplings to the third generation, as well as thirdgeneration squarks [93]. High-energy tops can be produced copiously at the LHC and at future accelerators, and multi-TeV top quarks offer a particularly rich laboratory to study the effects of weak showering.
We start by considering splittings that follow the same structure as the top quark's weak decay, t → W + b. Figure 6(a) shows the resulting W b mass spectrum from applying this splitting process to 10 TeV top quarks of left-handed or right-handed helicities. One immediate feature is the transition between shower and decay: the Breit-Wigner peak centered at m t continuously matches onto a high-virtuality shower dominated either by W T emission from left-handed top quarks, or W L emission from right-handed top quarks. 20 The former are simple manifestations of SU (2) L gauge showers with a larger rate (upper curve), whereas the latter are a due to the Goldstone-equivalent Yukawa showers with a smaller rate (middle curve). Ultra-collinear emissions are necessary for properly modeling the shower/decay transition, as shown in more detail in Appendix B (see Fig. 12). We also show the unpolarized top decay with a fixed-width Breit-Wigner without shower (lower curve in Fig. 6(a)). The events are understandably much more constrained to the region M (W b) ≃ m t . It is very important to appreciate the difference, for example since one must properly model the properties of off-shell top quarks in searching for new physics [93][94][95][96][97][98] associated with the top quark as well as the Higgs sector.
Top quarks may also radiate Higgs bosons and, analogously, longitudinal Z bosons. Both of these Yukawa-showering processes occur with similar rates off of left-handed and right-handed tops, and grow single-logarithmically with energy. In Fig. 6(b), we present a 10 TeV right-handed top quark splitting via the EW shower. The rates for t R → ht L and to Z L t L are governed by the Yukawa coupling and essentially the same, due to the GET. The channel t R → Z T t R , shown for reference, is via the gauge coupling of nearly pure B 0 , which is rather small. The other two channels t R → ht R , Z L t R are helicity-conserving scalar emissions and are of the ultra-collinear nature. The integrated splitting rates for all the above channels are of similar size: 20 To improve the matching, we have distributed the "decay" events according to a Breit-Wigner distribution weighted by Γt(Q)/Γt(mt). This constitutes approximately a 30% effect at the given matching scale of 187 GeV. the rates for the ultra-collinear processes are concentrated toward smaller virtualities (and correspondingly smaller k T s). Though the total splitting rate represented in Fig. 6(b) is only a few percent, the fact that top quarks are produced through strong interactions can lead to significant numbers of showered events at a hadron collider. On the other hand, the splitting rates to a Higgs boson are in sharp contrast to the much smaller rate for an on-shell top quark decay to a Higgs boson in the Standard Model [99], of the order 10 −9 .
In considering determination of the top-quark Yukawa coupling in the processes tth/ttZ at high energies [100], the qualitative features shown here should be informative.

EW Showers initiated by neutral bosons
The neutral bosons γ, Z T , h, and Z L contain rich physics at high energies, but their showering requires special treatment due to the presence of sizable interference effects.

γ/Z T coherence
For the γ/Z T system, these interference effects have two aspects: the mass basis is misaligned with the gauge interaction basis, and even when viewed within the B 0 /W 0 interaction basis, the existence of a preferred physical isospin basis for asymptotic states leads to observable coherence between B 0 and W 0 exchanges. A rigorous final-state shower must address both of these aspects simultaneously by using Sudakov evolution based on density matrices, as outlined in Section 2.3.2. More specific details can be found in Appendix C.
As a simple example of the basis alignment issue, consider high energy showering of neutral bosons γ/Z → W + W − . A naive treatment would shower the photon and Z including the triple-vector processes γ → W + W − and Z → W + W − . 21 However, depending on the gauge charges of the initial sources, the interference between these two mass-basis splitting channels can be O(1). In particular, for an energetic γ/Z emitted from a righthanded chiral electron line, the SU (2) L content of the produced neutral gauge bosons is practically zero, suggesting a near absence of collinear W + W − splittings in the final state. We explicitly compute these splittings assuming either an e − L or e − R source, which radiate off 2.5 TeV γ/Z bosons (e.g., via neutral boson pair-production at a 5 TeV e − e + collider). The results are displayed in Fig. 7. Our full EW FSR treatment is labeled as "coherent shower," contrasting with the hypothetical incoherent contributions from individual γ or Z. For the γ/Z produced by left-handed electrons in Fig. 7(a), the W 0 fraction is prominent from the constructive interference between γ/Z, leading to a total splitting rate of roughly 15% (black solid curve) and noticeable Sudakov distortions relative to a simple fixedorder splitting calculation (dashed curve). Fig. 7(b) shows the result for a right-handed electron source, exhibiting the almost complete destructive interference between the γ and Z channels, due to the fact that the produced boson is nearly pure B 0 when viewed in gauge basis. The small residual rate at high virtualities is actually dominated by the unbroken-phase vector-to-scalar splitting In our GEG approach, this is simply computed as a distinct process, rather than due to a delicate cancellation.
Perhaps more subtle are the interference effects between different exclusive isospin channels. Naively, we might expect to be able to treat SU (2) L × U (1) Y in a manner analogous to SU (3) QCD × U (1) EM , wherein the showers of the two gauge groups are simply run independently of one another. However, weak isospin quantum numbers are directly correlated with electric charge, and are therefore usually experimentally distinguishable. (Consider, e.g., the response of a detector to e L versus ν L .) Therefore, weak isospin cannot be summed/averaged like QCD color. As a consequence, observable rate asymmetries arise due to interference between the SU (2) L and U (1) Y gauge boson exchanges. Although a well-known effect, it has never been implemented in a parton shower framework. Again, we illustrate this by the splittings of 2.5 TeV γ/Z neutral bosons, here produced off of a left-handed chiral electron line. This boson may subsequently split into a ℓ − ℓ + or νν pair. The splitting rates with/without interference effects are shown in Fig. 8. 22 Besides the full coherent EW evolution (solid curves), two hypothetical incoherent treatments are shown using γ-Z mass basis (dashed curves) and B 0 -W 0 gauge basis (dotted curves). It is instructive to see that Z → ν LνR contribution alone gives the correct result as seen in Fig. 8(a); B 0 → ℓ − Rl + L alone also gives the correct result at high masses as seen in Fig. 8(c), although it misses substantial destructive interference near m Z due to the unequal γ and Z masses; and ℓ − Ll + R would need coherent treatment in the whole kinematical regime as seen in Fig. 8(b). The same issues of course arise in hadron colliders, though the numerical impact is often smaller because of the healthy admixtures of u/d flavors and LH/RH chiralities, as well as the charge-rearranging effects of hadronization. Nonetheless, we strongly advocate for a consistent treatment based on matrix-valued splitting functions and Sudakovs.

Higgs splitting and h/Z L coherence
Analogous interference effects also occur between the Higgs boson and longitudinal Z boson. In the high-energy gauge theory, these appear as different components of the same complex scalar, and particular linear combinations carry a partially-conserved "Higgs number" that flows through the shower. As a simple illustration, consider high energy production of W + T → (h/Z L )W + L . The coherently mixed h/Z L carries Higgs number of −1, and corresponds to the "anti-Higgs" state H 0 * . This state preferentially splits into W + T W − L (or, equivalently, W + T φ − ), as shown in the top curve of Fig. 9(a), labeled by the W helicities and charges as T + L − . The charge conjugate state W + L W − T (labeled L + T − ) carries the opposite Higgs number and thus is highly suppressed. It arises only at low virtuality, mainly due to the Higgs-Z mass difference. An incoherently-showered admixture of h and Z L would instead distribute probability equally between these two different polarization channels, as shown in the figure with the middle curve. (A similar charge-polarization correlation also occurs in splittings to top quark pairs.) The contributions from the other sub-leading ultra-collinear polarization channels are shown by curves labeled L + L − and T + T − . Though not obvious from the virtuality distributions, we note that coherence effects also significantly influence these channels. In particular, the ultra-collinear splitting H 0 * → W + L W − L inherits the soft divergence from the regular gauge splitting H 0 * → W + T W − L , but only in the limit as the W + L becomes soft. Similarly for the CP-conjugate process. The individual h and Z L incoherent showers, on the other hand, exhibit parts of the soft-singular behaviors of each of their H 0 and H 0 * components. See Table 6.
As a final novel example of neutral boson showering, we consider the purely ultracollinear splitting h → hh. This proceeds through the Higgs cubic interaction that arises after EWSB, and it is the unique 1 → 2 splitting process in the SM that is strictly proportional to Higgs boson self-interaction λ h . Isolating the h component of a general energetic h/Z L state, the total splitting rate comes out to about 0.14% for E ≫ m h . We illustrate in Fig. 9(b) the kinematic distribution ∆R(h, h), for an example initial Higgs energy of 1 TeV. The distribution peaks at roughly 2m h /E, which in this example is close to 0.25. Gener-ally, the majority of the phase space for high-energy production hhX for any X becomes dominated by such collinear configuration. While this ultra-collinear splitting process lacks any log-enhancements, integrating the splitting phase space yields a total rate relative to hX that scales like λ h /16π 2 , whereas the non-collinear regions contribute a relative rate of order λ 2 h /16π 2 × v 2 /E 2 . Therefore the "collinear enhancement" here is E 2 /λ h v 2 ∼ E 2 /m 2 h , rather than a conventional logarithm. Though the splitting rate is still quite small, for a 100 TeV pp collider with 10's of ab −1 integrated luminosity, we expect thousands of such events arising from the (also novel) high-energy production process qV L → q (′) (h/Z L ) at p T ∼ 1 TeV. In future precision Higgs physics [101], accurate description of such Higgs splittings could serve an interesting role.

EW showers by a new heavy state: W ′ example
The possibility of multiple weak boson emissions in the same event, and indeed even from the same parent particle, leads us inevitably to start considering final-states in terms of "weak jets" rather than in terms of individual, well-separated EW-charged particles (possibly dressed with QCD and EM radiation). Besides altering the energy spectra of the particles emerging from a hard interaction, EW emissions can significantly alter the multiplicity and flavor structure of an event. In particular, this new feature could have major consequences for how a new physics signal would be detected and reconstructed.
While it is beyond the scope of this current paper to present detailed examples for physics beyond the SM in high energy collisions [102], we study a simple case for illustration. We consider the decay of a narrow heavy W ′+ resonance into ν L ℓ + R , with a left-handed coupling and M W ′ ≫ m W . Nominally, the resonance is reconstructed from the charged lepton and the missing transverse momentum using the transverse mass variable M T (ℓ, E T ), which gives a Jacobian peak at M W ′ . When multiple EW emissions are taken into account, various new flavor channels open up, as well as additional kinematic handles that can facilitate more accurate resonance reconstruction. For example, in [60], it was pointed out that collinear weak emissions ν → Zν can effectively reveal the neutrino's direction-of-flight when the Z decays visibly. For illustration here, we simply divide up the showered signal by inclusive lepton multiplicity, focusing on channels up to three charged leptons. Quarks and τ -leptons may be present in the secondary W/Z showering/decays, but are ignored here for simplicity. Within each lepton multiplicity channel, we approximately reconstruct the resonance using the "cluster transverse mass" variable M T cl , defined as [103] The result of this analysis is displayed in Fig 10(a), taking M W ′ = 20 TeV. Solid curves are those from the nominal EW shower for 1ℓ+X, 2ℓ+X and 3ℓ+X, where X represents the rest of the particles in the event (mainly neutrinos and quarks). The dotted line shows the result of the naive two-body decay calculation, without the parton shower. To focus on the weakscale contributions, we have terminated the EW shower at a lower virtuality of 50 GeV. The showering reduces the total visible rate within 10% of the nominal peak by about 10% due to the radiation. In this window, the relative contributions from 1-lepton, 2-lepton, and 3lepton are respectively 0.81, 0.13 and 0.06. Although higher lepton multiplicities are rarer,  their M T cl distributions are also more sharply-peaked. It is also instructive to compare these predictions to those of a simple fixed-order splitting calculation, which captures the leading-log corrections but does not resum them. We find that this calculation predicts 9% more 1-lepton events than the full EW shower in the near-peak region.
Like e L and ν L , left-handed top and bottom quarks live together in a weak isospin doublet, and can also convert into one another through soft/collinear W ± emissions. Similar to the Bloch-Nordsieck violation effect discussed above for PDFs, the distinction between t L -and b L -jets therefore becomes somewhat blurred at high energy [62]. This effect, which is double-log enhanced at fixed order, is automatically resummed in the parton shower. Consider again, as a simplified example, a narrow 20 TeV W ′+ resonance, this time decaying to t LbR of 10 TeV each in energy. The final flavor content of two heavy quarks should gradually average out. We show in Fig. 10 the mass spectrum of the two-quark system resulting from the decay plus EW parton shower, individually in tb, bb, tt, and bt channels. (For this purpose, the threshold between the "shower" and "decay" of a top quark is set to m t + 10Γ t .) Respectively, these are dominated by unshowered events, events with a single t → W + b splitting, events with a singleb → W −t splitting, and events with one of each such splitting. The relative rates of the four channels are about 0.77, 0.09, 0.12, and 0.015. Within 10% of the W ′ mass peak, the nominal tb signal would be reduced by almost 30% from purely electroweak effects. Of course, this observation invites "weak jet" reconstructions that add back in the emitted gauge and scalar particles, though inferring the resonance's charge becomes somewhat more complicated.
Finally, we can consider the interplay of EW and QCD radiation, which is shown in Fig. 10(c) for the mass spectra of the quarks when t → gt and b → gb emissions are also turned on. Again the shower is terminated at 50 GeV virtuality to focus on effects at and above the EW scale. The full Standard Model showering leads to dramatic distortions in both mass and flavor distributions. Now the W ′ mass could be more accurately reconstructed by adding back in both-the EW and QCD radiation, which practically may overlap heavily since emitted weak bosons dominantly decay hadronically.

Summary and Conclusions
At very high energies, far above the electroweak scale, the full gauge and Yukawa structure of the Standard Model emerges, leading to an extremely rich set of parton showering phenomena. As this full SM parton shower evolves down in scale, it ultimately passes back through the electroweak scale. There it encounters additional showering phenomena that arise uniquely from EWSB, and then finally transitions back into the SU (3) QCD × U (1) EM gauge showers familiar from the past several decades of theoretical and experimental work.
With an eye towards experiments in the next decade and beyond, in this paper we have attempted lay out the above picture of electroweak showering in a more comprehensive manner. We have systematically presented the electroweak collinear splitting functions in the SM in the SU (2) L × U (1) Y symmetric phase as well as in the broken phase after electroweak symmetry breaking. We discussed their general features in the collinear and soft-collinear regimes and identified the general class of EWSB contributions that are uniquely "ultra-collinear," namely localized at k T ∼ v with appreciable rates, but otherwise absent in conventional showering regimes. Effects of the ultra-collinear part of the shower include counter-intuitive "violations" of the Goldstone-boson Equivalence Theorem. We have also identified a convenient way to isolate EWSB effects within the shower, especially by disentangling contributions from gauge bosons and Goldstone bosons at high energies, using a novel gauge choice which we call Goldstone Equivalence Gauge (GEG). We further implemented the full EW shower in a numerical monte carlo, and showed a number of new results regarding its subtleties and practical impact in SM processes and beyond.
Our main observations and results are as follows: • The splitting functions of the unbroken SU (2) L × U (1) Y theory, presented in Sec. 3, typically act as the leading contributions to showering processes at energies far above the EW scale.
• At splitting scales k T ∼ gv and yv, the unbroken splitting functions become regulated and the new ultra-collinear splitting functions arising from EWSB appear, as presented in Sec. 4. The latter is the analogue of "higher-twist" terms in terms of the formal power counting. While they do not contribute to the leading logarithmic evolution, numerically they can be larger than the unbroken contributions at low k T , and in some cases can also account for a sizable fraction of the integrated splitting rates.
• Goldstone-boson equivalence ceases to hold in the ultra-collinear regime, allowing, e.g., for emission of relativistic longitudinal bosons from massless fermions. This effect is generalized here to all splitting functions in the SM, often involving nontrivial interplays of EWSB effects in gauge, Yukawa, and scalar couplings.
• We introduced the Goldstone Equivalence Gauge (as detailed in Appendix A) that practically as well as conceptually disentangles the effects from the Goldstone bosons and the gauge fields. Utilization of this gauge choice makes the GET transparent and organizes its leading violations in a straightforward diagrammatic expansion (see Appendix B). The concept of a "nearly on-shell" gauge/Goldstone boson as an intermediate state in the shower also becomes unambiguous.
• We implemented a practical EW showering scheme based on the calculated collinear and ultra-collinear splitting kernels in a Sudakov formalism. As discussed in Sec. 2, some additional novel features in the implementation include matching between showering and resonance decay, kinematic back-reaction corrections for multiple emissions of massive particles, and a density matrix treatment for the mixed-state evolution of neutral bosons (γ/Z/h). Our treatment of EW showering is fully self-contained, and far beyond the currently existing monte carlo simulation packages.
• We applied the EW showering formalism to a number of important physical processes at high energies. They include: electroweak partons in PDFs as the basis for vector-bosonfusion; EW FSR as a leading source of multiple gauge boson production, with splitting probabilities at the level of 10s of percent; EW showers initiated by top quarks, including Higgs bosons in the final state; and showers initiated by neutral bosons γ/Z/h, for which care must be taken to obtain meaningful results. The emergence of "weak jets" from high-energy new physics processes was illustrated using a heavy W ′ as an example.
In summary, we have derived the collinear splitting functions for the Standard Model electroweak sector, including the massive fermions, gauge bosons, and the Higgs boson, and implemented a collinear showering scheme in the Sudakov formalism for all SM particles at high energies. We have highlighted many novel features and the necessity to include them for physics explorations in and beyond the SM at high energies, including any physics at future colliders, as well as other processes in high energy regimes much above the electroweak scale.
While our paper has explored collinear EW showering at a new level of detail compared to earlier works, it leaves open several interesting issues that we intend to address in future publications [63]. One such issue is a more comprehensive picture of PDF evolution, folding together QCD and EW effects into a unified set of DGLAP equations that incorporate both quantum coherence effects and ultra-collinear effects, and allowing for a complete QCD+EW ISR showering scheme. Implications for the exclusive structure of multi-TeV VBF events would be particularly interesting to study. We also intend to address issues related to soft wide-angle EW exchanges, which lead to quantum entanglements between the isospins of the beams and the final state at NLL. These entanglements represent a formally subleading aspect of the notorious Bloch-Nordsieck violation, which naively implies double-and single-logarithmic divergences in inclusive cross sections sourced by isospinexclusive initial states. The collinear formalism developed here would allow for simple LL resummation of the soft-collinear, double-logarithmic contributions. (See, e.g., Section 5.5 for simple examples in the final-state shower.) Capturing and resumming the remaining single-log, quantum-coherent contributions, as well as motivating factorization of the initial state at NLL, requires a more advanced formalism that uses the language of quantum ensembles.
Representing a generic gauge adjoint component of a vector field by W µ , we decompose the gauge degrees of freedom as the components of W n (Wn) aligned (anti-aligned) with n µ and the two ±1 helicity (or "xy") transverse modes, collectively W T : withn µ ≡ (1, +k sign(k 0 )). Since W µ is a real vector field here, we have chosen the above definition such that n µ (k) * = n µ (−k). Introducing the gauge-fixing Lagrangian in momentum space as the large light-like component of the on-shell longitudinal polarization, Wn field, ceases to propagate because of its infinite "mass" 1/ξ. This is the key feature for GEG by design.
We are left with three physical degrees of freedom that can propagate. It is interesting to note that GEG respects the rotational symmetry under SO(3) by construction. The surviving polarization states are also invariant (up to a possible rescaling) under boosts collinear to k. Incorporating EWSB, neither the gauge boson mass nor the would-be-Goldstone field φ are folded into the gauge-fixing procedure. The normalization of W n and its associated polarization vector ǫ µ n ∝ n µ can be chosen such that W n will interpolate external particles with unit amplitude: This polarization vector is what remains of the standard longitudinal polarization ǫ µ L (k) upon subtraction of the Goldstone-equivalence term (scalarization term) k µ /m W . Preserving Hermiticity of the W n field also necessitates introduction of a factor of i into the polarization vector, such that (iǫ µ n (k)) * = iǫ µ n (−k). This will also conveniently synchronize the phase of states created by the W n field and the φ field. 24 Accounting for the gauge-Goldstone mixing term, the quadratic Lagrangian can then be expressed as Inverting yields the propagators 24 When working in a complex gauge basis, as for W ± , these polarization phase factors become simply ±1.
In all cases, care must be taken to rigorously define the orientation of momentum flows when computing amplitudes, since ǫ µ n (−k) = −ǫ µ n (k), and the sign is often needed to determine the relative phase between gauge-interpolated and Goldstone-interpolated diagrams. Figure 11: Schematic tree-level collinear factorization for an arbitrary process with a splitting Goldstone/longitudinal in the final state.
These propagators are naively fully Lorentz-invariant, though choosing a polarization basis in the first place has anyway tied us to a specific frame. They share a unique, common pole at k 2 = m 2 W with residue +1. The mixed W n and φ fields interpolate the same particle: the "longitudinal gauge boson" or "Goldstone boson," depending on perspective. 25 Note that the apparent spurious pole at k 2 = 0 in the mixed propagator is purely an artifact of our momentum-dependent field normalization, and does not lead to light-like gauge poles in complete Feynman diagrams. 26 Goldstone boson equivalence in the high-energy limit now emerges trivially, diagramby-diagram. For a process where |k 2 | ≫ m 2 W for all internal gauge/Goldstone lines and E ≫ m W for all external bosons, the mixed propagators and ǫ n factors scale away, leaving over only the Goldstone contributions. In addition, since there are no terms that go like k/m W or E/m W , power-counting of corrections ∝ m W becomes straightforward at the level of individual Feynman diagrams. Upon introduction of complete fermion and scalar sectors, we may generalize to counting VEV factors associated with arbitrary masses and interactions introduced by spontaneous symmetry breaking. Some simple examples for splitting calculations are given in Appendix B.
We can also see how this gauge choice facilitates a factorized picture of longitudinal gauge/Goldstone boson production and splitting in the parton shower, beyond the simple Goldstone-equivalent picture at zeroth-order in the VEV. Fig. 11 illustrates how this works schematically in a final-state shower. A generic hard process produces an off-shell gauge/Goldstone boson of virtuality k 2 with m 2 W ≪ k 2 ≪ E 2 , and this boson subsequently splits. There are four contributing classes of diagrams, corresponding to the four possible propagator exchanges between the production and splitting processes. We would like to approximate this as an on-shell production amplitude multiplied by a universal splitting amplitude. The decomposition is trivial for the leading pure Goldstone exchange diagram, but the other, subleading diagrams involve interplays between the propagators and the offshell polarization vectors ǫ µ n ∝ ( √ k 2 /E)n µ . For the mixed diagrams, the propagator factor m W / √ k 2 can be combined with the polarization factor √ k 2 /E to yield an approximate on-shell polarization proportional to m W /E. Assuming that there is no large back-reaction in the hard production matrix element (at least to O(m 2 W )), contracting with the rescaled off-shell polarization approximately reproduces the on-shell hard process. For the mixed diagram where the gauge field contracts with the splitting process, this decomposition would simply instruct us to compute the splitting amplitude with an effective on-shell ǫ n . The pure gauge exchange does not immediately fit this pattern, but it can be separated into two pieces: The former piece has the correct structure to provide m W / √ k 2 factors to each gauge polarization. The latter piece cancels the √ k 2 's from each polarization vector, but leaves over no poles or mass factors. It therefore produces a non-collinear interaction that goes as 1/E 2 instead of 1/(k 2 − m 2 W ), and can be grouped together with the neglected non-collinear diagrams. We can view all of the remaining collinear contributions as a simple product of on-shell gauge+Goldsone production and gauge+Goldstone splitting matrix elements, connected by the standard scalar propagator i/(k 2 − m 2 W ). Analogous results were obtained for the factorization of logarithmic virtual corrections to external gauge/Goldsone bosons in [74] by working directly in Coulomb gauge, and in [34,35] by invoking the Goldstone Boson Equivalence Theorem in Feynman-'t Hooft gauge. Our own approach directly exhibits the applicability of the Equivalence Theorem in the corresponding real emission processes at tree-level, and extends them beyond the strict Goldstone limit to O(m W /E).

B.1 Lagrangian, couplings, and charge conventions
In Goldstone Equivalence Gauge, each physical longitudinal gauge boson state is interpolated by two fields: V n and φ V , where V = W ± , Z. Unlike, e.g., in R ξ gauges, the relative phases of V n -mediated and φ V -mediated processes must be explicitly kept track of. Here, we first present the Lagrangian of the SM in GEG to set the conventions. Before electroweak symmetry breaking (EWSB), the Lagrangian with the gauge fixing is written as L Ghost =c a n µ D ab µ c b .
The flavor indices are suppressed since we do not consider the effects of flavor mixing. The covariant derivative D µ and SU (2) L field strength component W a µν are defined as The gauge-fixing vector n µ of Eq. (A.1) would here be treated as a differential operator of schematic form (1, −∂ t ∇/ ∂ 2 t ∇ · ∇). This becomes a well-defined operation in momentum-space. We take the formal ξ → 0 limit in what follows.
After EWSB H 0 = v/ √ 2, and particles acquire masses. The neutral gauge fields W µ 3 and B µ mixing to form mass eigenstates Z µ and A µ . Gauge and fermion masses go as with g 1 ≈ 0.36 and g 2 ≈ 0.65 at the weak-scale, y t ≈ 1, and v ≈ 246 GeV. The Higgs field self-coupling is normalized such that such that λ h ≃ 0.52 for m h ≃ 125 GeV.
As for the gauge-fermion interactions in a general basis, we denote them using g V as the gauge coupling constant for a vector boson V = B 0 , W 0 , W ± , γ, Z, where the chirality projection operators are P R/L = 1 2 (1 ± γ 5 ). They are all built up from the underlying U (1) Y and SU (2) L gauge couplings. Specifically, As usual, the weak mixing angle is defined as We denote the gauge charge Q of a particle p (chiral fermion or scalar) under a given gauge boson V by Q V p . 27 We list the full set of charges in Table 8. We now turn to the quadratic Lagrangian terms involving gauge fields and Goldstone fields. The quadratic terms of Z and φ Z Lagrangian are (B.9) 27 For V = W ± , two different components of a left-handed doublet participate, but they can be assigned a common charge of 1/ √ 2, with either flavor plugged in.
This in turn determines the phase factor of the polarization vector Z n . (Though of course our convention choices ultimately have no effect on physical rates.) For W ± µ /φ ± , the unmixed kinetic and mass terms are analogous, and the quadratic mixing term is given by

B.2 External polarizations and propagators
We decompose all fermions and gauge bosons into helicity basis within the hard process CM frame, including off-shell particles. We emphasize that in computing leading-order 1 → 2 splitting functions, all particle polarization states should be set on-shell, since the off-shell corrections are strictly non-collinear. An on-shell polarization can be associated with an off-shell momentum, for example, by adjusting the three-momentum at fixed energy.
The fermion external spinors are as usual, though to facilitate extraction of O(v) effects we Taylor expand in m f /E = (y f / √ 2)(v/E). Explicitly, for fermions moving approximately along the z-axis, possibly offset toward the x-axis by a small angle θ, (B.11) Propagators are also as usual, but given our approximate decomposition into on-shell spin states, they fall into a factorizable form. For a generic off-shell k µ , we can build an effective on-shellk µ by keeping k 0 ≡ E fixed but changing We may then rewrite the propagator as exploiting the fact that the leading correction away from a factorized numerator is set up to cancel the propagator's denominator. We ignore possible coherence effects between different spin channels. Transverse gauge bosons are also assigned their standard polarization vectors with the complex-conjugate ǫ µ * ± used for outgoing bosons. However, the longitudinal gauge/Goldstone sector is treated somewhat unconventionally. Longitudinal gauge bosons can be created by Goldstone/pseudo-scalar boson fields. We set our phase conventions so that these creation and annihilation amplitudes are unity, maintaining continuity with the unbroken theory. However, longitudinal bosons may also still be created by gauge fields, in association with the "remainder" field component V n expanded out in Eq. (A.2). Synchronizing these component fields such that they also create/annihilate external bosons with unit amplitude, their associated polarization vectors then carry nontrivial phases: incoming Z : iǫ µ n ; outgoing Z : (iǫ µ n ) * = −iǫ µ n ; incoming/outgoing : The light-like gauge-fixing vector n µ is defined in Eq. (A.1). The corresponding propagators are given in Eq. (A.6). Photons are subjected to the same gauge conditions, but this has little practical bearing on their showering behavior. As usual, they have purely transverse external polarization states, and only their transverse modes contribute to collinearenhanced physics. As discussed in Section 2.3.2, the transverse photon and Z propagators should be treated coherently within a parton shower. The h and φ 0 /Z n propagators should also be treated coherently. We will see an example of this, including the Z n component, in Appendix B.4.

B.3 Feynman rules for three-point couplings
Feynman rules in GEG are largely similar to those of standard gauges. We list below many of the relevant three-point vertex rules. For brevity, we omit four-point interactions, which do not play a role in 1 → 2 splittings at this order. Wherever explicitly referenced, we reckon all four-momenta as flowing into the vertex. We use the small arrows next to a particle line to indicate the flow of the momenta as well as the electric charge, where relevant. When no arrows labelled for the charged particles, charge conservation is implied at each vertex for the particles involved.
The symbol ⊗ denotes the mass (or v) insertion from the EWSB.

B.4 Example calculations with GEG
Calculations in high energy processes involving longitudinal vector bosons can be complicated in dealing with gauge artifacts, often exhibiting artificial "bad high energy behavior" containing factors of E/v. Here we show some explicit examples to demonstrate how to calculate ultra-collinear splitting amplitudes in GEG, where all such amplitudes are automatically free of such artifacts and are simply proportional to the VEV. We focus in detail on the specific massive fermion splitting t s → W + L b s , where the fermion helicity s = L, R is preserved. This calculation is also trivially adapted to cases where one or both fermion is a massless flavor, such as the usual u L → W + L d L , and is straightforward to extend to Z L boson emission with appropriate replacements of couplings and remainder polarization phases. We also outline below the diagrammatic construction of a few other processes for illustration.
We first reemphasize that the longitudinal gauge boson W + L in GEG should be interpolated by both the Goldstone field φ + and the remainder gauge field W + n , leading us to break up the splitting amplitude as Applying the three-point Feynman rules in Sec. B.3, and taking the exact collinear limit (θ, k T → 0) to extract the leading behavior, we have for the LH process 29 The full LH splitting amplitude is then Plugging this in Eq. (2.5), we have the splitting function . (B.20) As for the RH transition t R → W + L b R , there is no analogous amplitude for W n at O(v) due to the absence of RH charged-currents, so the amplitude is dominated by the Yukawa contribution, and the splitting function is  Table 4. When combined with conventional collinear top quark splittings, the ultra-collinear splittings become important for modeling the approach to the top resonance peak. This includes as well the process t R → W + T b L . We show these individual shower contributions and their continuity with a simple Breit-Wigner model of top decay (weighted by Γ t (M (W b))/Γ t (m t )) in Fig. 12. Here we have taken 10 TeV top quarks of either helicity, zooming into near the top quark pole, and set a decay/shower matching threshold of 187 GeV. All polarizations are measured in "lab frame" (as opposed to the top's rest frame). QCD and other electroweak showering effects are not incorporated.
We have seen above how GEG allows us to organize the amplitude's dependence on EWSB by explicitly decomposing it into individual mass-insertion terms, or equivalently VEV-insertion terms. External-state fermion mass insertions are found by Taylorexpanding the fermion Dirac spinors, and external-state gauge boson mass insertions are found via the remainder polarization ǫ n . For more general processes, there may also be three-point interactions that function as VEV-insertions, such as interactions between the scalars or the hV µ V µ vertices (listed in Sec. B.3). Generally, we may rather straightforwardly construct any ultra-collinear amplitude at O(v) by adding together diagrams with exactly one mass-insertion or EWSB interaction. Besides helping to organize a calculation, this approach serves as a convenient tool for visualizing where different EWSB contributions arise in a given amplitude.  Table 4; • Fig. 13b: W ± T → W ± L Z T , representative calculation for Table 5; • Fig. 13c: Z L → W + L W − L , representative calculation for Table 6; • Fig. 13d: h → W + L W − L , representative calculation for Table 6.

C Coherent Showering
Showering involving superpositions of different particle species can be described using density matrix formalism. The initial value of the density matrix is proportional to the outer product of production amplitudes Figure 13: Representative ultra-collinear splittings with multiple contributing diagrams. The effects of the VEV are indicated schematically by the symbol ⊗.
tracing out over other details of the rest of the event. Here, the indices run over the species. We nominally assign the state its smallest possible kinematic mass (zero for γ/Z, m Z for h/Z L ), and subsequently reweight/veto the splitting probability and adjust the global kinematics as necessary (see Section 2.3.1). This prescription specifically becomes relevant when evolving near kinematic thresholds.
The probability for an initial mixed quantum state to subsequently split into a specific exclusive final state, e.g. γ/Z → e − L e + R or ν LνR , must be computed by generalizing the splitting functions to Hermitian splitting matrices dP ij . The exclusive splitting rates are then computed by tracing against the normalized density matrix: If a boson is not split, the Sudakov evolution of ρ proceeds analogous to mixed-state radioactive decay: dρ ij = − 1 2 channels (ρ ik dP kj + dP ik ρ kj ). (C.2) As usual, this just represents the wave-function running, now applied to multi-component states. The splitting matrices for an initial mixed quantum state are computed from outer products of splitting amplitudes, convolved with the mixed propagators. Representing the propagator matrix as D ij , and the amputated splitting amplitudes as M In more complete generality, a mixed state can split into another mixed state, leading to an enlarged set of indices for the splitting matrices. However, in most cases, the finalstate density matrices are fully determined by the initial-state density matrices, such that in practice a single pair of indices suffices.
While the formalism is basis-independent, we default to some standard bases in our EW shower approach. Within the unbroken phase (Section 3), we present neutral gauge and scalar splitting functions in the interaction basis (B 0 , W 0 ), (H 0 , H 0 * ). In the broken phase (Section 4), we present them in the mass basis (γ, Z), (h, Z L ). The corresponding propagator matrices in the unbroken-phase basis, including the effects of EWSB, are 30 for the gauge bosons (θ W is the weak mixing angle), and for the neutral scalars. In the mass basis, the matrices are diagonal and have entries corresponding to the usual poles: Similar considerations apply in the application and generation of PDFs [43]. The γ/Z and (in principle) h/Z L PDFs should each properly be treated as 2 × 2 matrices, and hard process cross sections sourced by these PDFs computed by tracing against the hard matrix elements. The PDF evolution equations involve matrix-valued splitting functions. In the high-k T /high-virtuality limit, these follow straightforwardly from the splitting functions presented in the Section 3. However, unless working well above the TeV-scale, mass effects can still be important. The above propagator modifications must then be applied at the (spacelike) virtual leg emerging from a splitting.

D Final-State Shower Simulation
In order to facilitate studies of final-state weak showering at the level of exclusive rates, we have programmed a variation of the PYTHIA6 [68] timelike virtuality-ordered parton shower. Basic collinear QCD is included by default, extended to the massive showering formalism outlined in Section 2, and including purely ultra-collinear processes. In addition, the full set of weak showering processes described in this paper has been added, with a number of novel features compared to standard showering programs, outlined in the main text. In particular, see Section 2.3. Here we describe a few additional technicalities of the implementation.
Splitting functions in the virtuality-ordered shower are simple to relate to those in the k T -ordered shower, which we have used by default for most of the presentation. Using the relativistic/collinear approximation for a splitting A → B + C, we get Working in log Q, we can build the translation dP dz d log Q 2 ≃ 1 zz (D. 2) The function in parentheses goes either as k 2 T or as v 2 . For given Q, z, and daughter masses, the former is simple to derive either by inverting the approximate Eq. (D.1) or by using exact kinematics. For the energy-sharing variable z, we use CM-frame threemomentum fraction | k B |/(| k B | + | k C |). To approximately model the phase space effects in the nonrelativistic limit, we further weight the splitting probabilities by a velocity factor | k B || k C |/E B E C . We also suppress splittings at angles larger than θ ≈ π/2, where the collinear shower would be highly untrustworthy.
As in PYTHIA6, the input to the shower is a "hard" partonic configuration with some characteristic virtuality scale, assumed here to be large compared to the weak scale. Evolution is based on a simple recoiler method, whereby particles are showered in pairs. (At the current level, no dipole coherence effects or color/isospin flows are incorporated, nor are they strictly necessary at leading-log level, but they would be possible to include in more advanced approaches.) Each particle in a pair undergoes a trial QCD/EW Sudakov evolution, defined in the hard event's rest frame, and ignoring the possible evolution of its sister. In general, each particle may undergo a 1 → 2 splitting and acquire an offshell mass. Kinematics are then adjusted within the pair's rest frame, by boosting each showered system along the pair's axis to preserve momentum and energy. If the summed masses from the trial evolutions exceeds the original pair's mass, the more off-shell splitting is vetoed, and that particle's evolution restarted. The procedure is easily recursed to build up completely showered events, with the two daughters from a given splitting serving as paired sisters in subsequent evolution.
Kinematic back-reaction effects are also incorporated, as discussed in Section 2.3.1 and parametrized in Eq. 2.9. The kinematic re-arrangments required by setting a daughter offshell through its secondary showering can have a sizable effect on the mother's splitting rate. We introduce this back-reaction factor as an additional weight multiplying the daughter's splitting probability. In our virtuality-ordered implementation, the virtuality of the mother (invariant mass of the daughter pair) remains unchanged, so Q * = Q. The Jacobian for the transformation is then simply |dz * /dz|, and its explicit form is tied to our kinematic prescription above. Within the mother splitting A → B + C, assume that particle B with momentum-fraction z is the one to be set off-shell: B → B * . Within the A rest-frame, the direction of B (C) is held at a fixed angle Θ (π − Θ) relative to A's boost axis from the CM-frame. The angle Θ has a one-to-one mapping to both the old z and the new z * , and is a useful intermediate variable. The symbols A, etc, here are shorthand for various quantities built out of A's velocity β A , and daughter kinematics in its rest-frame: the A-frame three-momentum magnitude of either of the daughters P , and their individual A-frame energies and kinematic masses E B , E C , m B , m C . We have (D.6) Analogous formulas hold with Y * and z * , defining the coefficients A, etc, using the A-frame kinematic quantities redefined with B set off-shell. (Prescriptions yielding simpler analytic formulas than ours almost certainly exist.) The differential splitting function of the mother must also be re-evaluated using the off-shell daughter kinematics. This is much simpler, as there the main effect is just the change in z. Explicit EWSB mass factors for the daughters, which appear in the numerators of the ultra-collinear splitting functions, are not adjusted from their on-shell values. Angular-ordering may also be invoked. If the showering pair was itself produced from a splitting, the event-frame angles of each daughter splitting and mother splitting can be compared, and the former splitting(s) vetoed if it has a larger angle. This veto may be applied selectively depending on the nature of the splitting and its parent splitting.
In our approach, parton shower evolution is automatically matched onto decay for W ± , Z, Higgs, and top. This matching is particularly simple in the virtuality-ordered shower. Particles that survive down to their decay/shower matching scale are assigned masses drawn from a Breit-Wigner distribution and final-state flavors assigned according to known branching fractions. In practice, we also weight the Breit-Wigner distribution accounting for the different available decay phase space at different off-shell virtualities. Similar to a shower splitting, the decays are then further weighted with back-reaction factors, if the decaying particle was itself produced in a splitting. The back-reaction factor here is applied as a simple probabilistic veto.
Finally, we re-emphasize that the neutral bosons γ/Z T and h/Z L are produced and evolved as general quantum mixed states. They are assigned initial kinematic masses of zero and m Z , respectively, and given nontrivial 2 × 2 density matrices that evolve via matrix-valued Sudakov factors. There is one major practical difference in implementing these Sudakovs relative to simple number-valued Sudakovs. In the latter case, a given particle's wavefunction decreases in magnitude as its evolution proceeds, but the surviving probability is an automatic outcome of the differential splitting rates integrated via monte carlo. In practice, these splitting rates are integrated over z with the expedient of overestimator functions, and vetoed-down to the true rates. In the matrix-valued case, however, the wavefunction can also rotate, and capturing this effect using over-estimator functions and a veto algorithm does not appear to be as straightforward. Instead, we use explicit formulas for the z-integrated splitting matrices at each virtuality step. These formulas are necessarily approximate, but we have verified that they yield results similar to what would be obtained by costly brute-force numerical integration.