1 Introduction

Beyond the discovery of the Higgs boson [1, 2], signs of new physics have yet to appear at the LHC and the Standard Model has so far survived all forms of scrutiny. It has therefore become more likely that the Standard Model continues to describe nature accurately up to very high energy scales. At these very high energies heavy particles like electroweak gauge bosons, Higgs bosons and top quarks can start to appear as constituents of jets [3, 4] or otherwise contribute to radiative corrections. Virtual corrections have been shown to become sizeable even at LHC energies in exclusive observables [5,6,7,8,9,10,11,12,13,14,15,16,17,18]. For instance, corrections to transverse momentum at LHC energies can already reach about 10% for exclusive dijet production [14, 15], and about 20% for single vector boson production [16,17,18], and they can be expected to grow even larger at future collider energies [19, 20].

Work such as [21, 22] has instead focussed on the incorporation of electroweak logarithms in the resummation of inclusive observables at high energies. However, many practical observables are not fully exclusive or fully inclusive. In most cases, the only practical solution is to instead include EW effects in a systematic way in parton shower algorithms as part of a general-purpose event generator such as Pythia [23], Sherpa [24] or Herwig [25]. These produce fully differential final states, and can be set up to include specific phenomena exclusive to the EW sector such the effects of symmetry breaking and gauge boson decay. Furthermore, the inclusion of EW effects in parton shower algorithms enables automatic interleaving with the QCD shower, as well as the interfacing with the models of nonperturbative aspects of the event generation process, such as multi-parton interactions and hadronization. As an example, ATLAS has reported on measurements that are sensitive to the collinear enhancements associated with W radiation in jets [26], and compared them to event generators that incorporate these effects. There, it is also pointed out that these types of effects will play a significant role for several measurements at high energy scales, which will become more abundant as the LHC gathers more data.

Electroweak corrections have been incorporated in parton showers in the past. An implementation [27, 28] is available in Pythia event generator [23] which only includes the radiation of electroweak gauge bosons and does not retain any spin information. The radiation of electroweak gauge bosons was similarly included in the Sherpa event generator [24] to study W emissions in jets [29]. Another approach [30] was employed in ALPGEN [31] where fixed-order matrix element calculations are combined with analytic Sudakov factors to achieve results similar to those of an electroweak parton shower. A more recent work [32] has implemented a final-state electroweak shower in the Pythia \(1 \rightarrow 2\) transverse momentum ordered shower formalism that retains spin information and includes all branchings present in the electroweak sector.

In this paper, extend the antenna-based Vincia parton shower [33,34,35] to include all possible final-state collinear EW branchings, including those associated with Yang–Mills triple vector boson and Higgs couplings. Furthermore, collinear vector boson emissions off the initial state are also included. Vincia is a plugin to the Pythia event generator and already allows for QCD evolution with partons of definite helicity states [36, 37]. This feature is especially important in the electroweak theory due to its chiral nature. The electroweak shower described in this paper will thus be responsible for the electroweak component of the shower evolution, and is interleaved with the default Vincia QCD shower.

The formalism descibed in this paper makes use of the spinor-helicity formalism to compute the large number of spin-dependent branching kernels associated with the electroweak shower. Helicity-dependent QCD antenna functions have previously been computed with comparable methods [38, 39]. We start with a brief overview of the spinor-helicity formalism and the conventions used in the calculation of the branching kernels. In Sect. 3, the spinor-helicity formalism is used to compute branching kernels for all branching processes in the electroweak sector. Section 4 discusses the collinear limits of those branching kernels given in terms of Altarelli–Parisi splitting functions [40]. The results in that section are found to be in agreement with [32]. Section 5 details the implementation of an electroweak shower in the Vincia framework and treats a number of peculiarities exclusive to the electroweak sector such as the presence of bosonic interference and the matching to resonance decays. To show the significance of an electroweak shower, its effects are investigated in Sect. 6 for highly energetic particles, but also at LHC energies. We finally conclude in Sect. 7, discussing a number of missing features that are required to bring the modelling of EW radiation further in line with current QCD showers.

2 The spinor-helicity formalism

Due to the chiral nature of the electroweak theory, it is important to calculate electroweak branching kernels for individual spin states. We choose to perform these calculations using the spinor-helicity formalism using conventions similar to those described in [41]. The branching kernels computed here describe the correct soft and collinear factorization properties of branchings in the electroweak sector, while remaining Lorentz-invariant and independent of a particular representation of the Dirac algebra or explicit forms for fermionic spinors. Furthermore, their analytic nature will be useful in dealing with issues that appear due to gauge dependence for longitudinal polarizations of gauge bosons. We first briefly summarize our conventions and techniques.

2.1 Spinors

Helicity spinors for massive fermions may be defined as

(1)

where \(\lambda \) is the fermion helicity and k is a lightlike reference vector that defines the meaning of the helicity of the fermion. Due to its massive nature, helicity is not a Lorentz-invariant quantity and does not coincide with the chirality of the fermion. The spin vector associated with the spinors defined in Eq. (1) is

$$\begin{aligned} s_{\lambda }^{\mu } = \frac{\lambda }{m} \left( p^{\mu } - \frac{m^2}{p{\cdot }k} k^{\mu } \right) . \end{aligned}$$
(2)

We therefore choose the reference vector

$$\begin{aligned} k = (1,-\vec {e}), \end{aligned}$$
(3)

where \(\vec {e}\) is a unit vector pointing in the direction of \(\vec {p}\). With this choice, the massive helicity spinors retain the usual meaning of helicity as the projection of spin along the direction of motion.

2.2 Polarization vectors

The polarization vectors for a massive vector boson with momentum p are defined as

(4)

where k is again a massless reference vector. Here, \(\epsilon _{\pm }^{\mu }(p)\) are the transverse polarizations and \(\epsilon _{0}^{\mu }(p)\) is the purely longitudinal polarization which only exists for massive vector bosons. By again choosing Eq. (3), it is immediately clear that \(\epsilon _{\pm }^{\mu }(p)\) are purely transverse and \(\epsilon _{0}^{\mu }(p)\) is purely longitudinal.

2.3 Spinor products and amplitude evaluation

Having expressed all massive spinors and polarization vectors in terms of massless spinors, amplitudes for particles with definite helicities can now be calculated very efficiently. We first define the spinor product

$$\begin{aligned} S_{\lambda }(k_a, k_b) \equiv {\bar{u}}_{\lambda } (k_a) u_{-\lambda }(k_b) \end{aligned}$$
(5)

for lightlike (reference) vectors \(k_a\) and \(k_b\), which obey

$$\begin{aligned} |S_{\lambda }(k_a, k_b)|^2 = 2 k_a{\cdot }k_b. \end{aligned}$$
(6)

Spinor products are not uniquely defined, but one possible representation is given by

$$\begin{aligned} S_{\lambda }(k_a, k_b)&= (\lambda k_a^2 + i k_a^3) \sqrt{\frac{k_b^0 - k_b^1}{k_a^0 - k_a^1}} \nonumber \\&\quad - (\lambda k_b^2 + i k_b^3)\sqrt{\frac{k_a^0 - k_a^1}{k_b^0 - k_b^1}}. \end{aligned}$$
(7)

Using the spinors and polarization vectors of the previous section, all amplitudes can be expressed in terms of these spinor products. Structures like

(8)

may still appear, where \(p_i, p_j, \ldots \) may be massive. Such structures may be expressed in terms of the spinor products Eq. (5) by defining

$$\begin{aligned} {\hat{p}}_i = p_i - \frac{p_i^2}{2 p_i{\cdot } k_i} k_i, \end{aligned}$$
(9)

which is explicitly massless. Eq. (8) may then be written as

$$\begin{aligned} S_{\lambda }(k_a, p_i, p_j, \ldots ,k_b) = S_{\lambda }(k_a, {\hat{p}}_i) \, S_{-\lambda }({\hat{p}}_i, p_j,\ldots , k_b). \end{aligned}$$
(10)

This procedure is then repeated, the next time making \(p_j\) massless by subtracting \((p_j^2 / 2 {\hat{p}}_i {\cdot }p_j) {\hat{p}}_i\), until only spinor products remain.

3 Electroweak branching amplitudes

We now use the spinor-helicity formalism to compute branching amplitudes for the electroweak sector. We first recount the phase space regions where radiative amplitudes factorize into a non-radiative amplitude and a radiative correction. For the moment, we restrict ourselves to splittings in the final state, where we denote the splitting momenta as \(p_{ij} \rightarrow p_i + p_j\). The parent momentum \(p_{ij}\) may be considered to be off-shell by an amount that vanishes in the singular limits.

We first consider the quasi-collinear limit [42,43,44], which may be defined using the Sudakov decomposition

$$\begin{aligned} p_i&= z p_{ij} + \alpha _i n + q_{\perp } \nonumber \\ p_j&= (1-z) p_{ij} + \alpha _j n - q_{\perp } \end{aligned}$$
(11)

where n is a lightlike reference vector usually taken to be the anti-collinear direction, and \(q_{\perp }\) is the spacelike transverse momentum with respect to \(p_{ij}\) and n. These momenta satisfy \(q_{\perp }{\cdot }p_{ij} = 0\) and \(q_{\perp }{\cdot }n = 0\). The variable z is the collinear momentum fraction, and the parameters \(\alpha _i\) and \(\alpha _j\) are fixed by the conditions \(p_i^2 = m_i^2\) and \(p_j^2 = m_j^2\). The off-shellness \(Q^2\) is then given by

$$\begin{aligned} Q^2&= (p_i + p_j)^2 - m^2_{ij} \nonumber \\&= \frac{p_{\perp }^2}{z(1-z)} + \frac{m_i^2}{z} + \frac{m_j^2}{1-z} - m_{ij}^2, \end{aligned}$$
(12)

where \(p_{\perp }^2 = -q_{\perp }^2\). The singular factorization of the radiative amplitude occurs when this off-shellness becomes small with respect to the energy scale of the process. In the massless case, the collinear region is thus defined by \(p_{\perp } \rightarrow 0\). This limit then generalizes to the quasi-collinear limit in the massive case, defined by

$$\begin{aligned} p_{\perp }, m_i, m_j \rightarrow 0 \text { with fixed ratios } \frac{m_i}{p_{\perp }}, \frac{m_j}{p_{\perp }}. \end{aligned}$$
(13)

In this limit, matrix elements symbolically factorize as

(14)

where \(P(\lambda _{ij}, \lambda _{i}, \lambda _{j} , z)\) is the helicity-dependent Altarelli–Parisi splitting kernel [40, 45,46,47,48].

The second singular limit is the soft limit, where

$$\begin{aligned} E_j \rightarrow m_j \text{ and } E_i \gg E_j \end{aligned}$$
(15)

or equivalent with \(p_i\) and \(p_j\) switched, with the soft particle being a gauge boson. The amplitude exhibits the usual eikonal factorization

(16)

where c is some spin-dependent coupling.

Using the spinor-helicity formalism, we compute \(1 \rightarrow 2\) electroweak branchings kernels that capture these soft and quasi-collinear factorization properties. To that end, we compute a branching amplitude

(17)

for every possible \(1\rightarrow 2\) branching in the electroweak sector. The branching kernel is then given by

$$\begin{aligned} B&_{\lambda _{ij}, \lambda _i, \lambda _j}(p_{ij}, p_i, p_j) = \frac{1}{Q^4} |M^{\{ij\} \rightarrow i j}(\lambda _{ij}, \lambda _i, \lambda _j)|^2. \end{aligned}$$
(18)

A comprehensive list of the branching amplitudes is given in Appendix B, which also includes vector boson emission from initial state fermions. Here, we highlight one of these calculations to discuss the treatment of artifacts that may appear for longitudinal polarizations.

We consider vector boson emission from a fermion. The branching amplitudes are straightforwardly found to be

(19)

where v and a are the vector and axial couplings of the fermion and vector boson at hand as defined in Appendix A. Note that this amplitude is applicable to all types of electroweak vector boson emission off fermion branching processes, including the cases where the fermion changes flavour due to a W-emission, or where the vector boson is massless and the axial coupling drops out in case of a photon. For massive vector bosons, inserting the corresponding polarization vector leads to the appearance of terms proportional with

(20)

The term proportional to \(Q^2\) cancels the propagator, leaving a contribution that is not singular in the quasi-collinear or soft limit. We note that similar situations may appear in equivalent calculations of QCD branchings kernels. Nonsingular gauge-dependent terms may appear that do not contribute to the logarithmic accuracy of a parton shower [49, 50] and may be used as parameters for uncertainty estimation [34]. In the case of emission of electroweak gauge bosons, the situation is however not the same. Due to the \({\mathcal {O}}(E/m)\) scaling of the longitudinal polarization, the leftover contributions lead to unphysical unitarity-violation at large energies.

Such problematic terms originate from the scalar component of the longitudinal polarization of the electroweak gauge bosons, which originate from their corresponding Goldstone bosons. In the above example, this Goldstone boson couples to the fermion through its Yukawa, which in this case appears as its kinematic off-shell mass. Due to the analytic nature of the spinor-helicity formalism, the offending terms are straightforwardly identified and corrected. Equivalently, the unitarity-violating term of Eq. (20) cannot survive if all other Feynman diagrams of a complete branching process are included, while all other terms are associated with a propagator factor \(1/Q^2\) which cannot be cancelled elsewhere. In Eq. (19) the unitarity-violating terms may thus be removed by the replacement

(21)

The replacements required for all other branching amplitudes are described in Appendix B.

In leading-colour QCD showers, soft gluon interference may be approximately incorporated by angular ordering [51], or by \(2 \rightarrow 3\) dipole-like functions [52,53,54] or antennae functions [34, 35, 55]. In the EW sector, no leading-colour approximation exists and many eikonal terms may contribute to the soft limit of an EW gauge boson emission. An algorithm that includes the full multipole structure of QED radiation in the Vincia shower was described in [56, 57]. The inclusion of the all soft interference effects is accomplished by defining a \(n \rightarrow n+1\) branching kernel and sectorizing the phase space into multiple regions to regulate the intricate soft and collinear singular structure. An extension of this algorithm requires the computation of an equivalent spin-dependent \(n \rightarrow n+1\) branching kernel for soft-enhanced EW gauge boson radiation, as well as an implementation that involves a similar sectorization of phase space. This is beyond the scope of this paper and is left for future work.

4 Collinear limits of the branchings amplitudes

In this section we discuss the behaviour of the branching amplitudes in the quasi-collinear limit described by Eq. (13). In this limit, the reference vectors simplify to

(22)

The branching amplitudes can be expressed in terms of the energy sharing variable z by replacing

$$\begin{aligned} p_i \rightarrow z \, p_{ij} \text{ and } p_j \rightarrow (1-z) \, p_{ij}. \end{aligned}$$
(23)

The only two remaining spinor products in the branching amplitudes are related by

$$\begin{aligned} S_{-\lambda }(k, p_j, p_i, k) = -S_{-\lambda }(k, p_i, p_j, k). \end{aligned}$$
(24)

Up to a phase factor, they are

(25)

where

$$\begin{aligned} {\tilde{Q}}^2 = Q^2 + m_{ij}^2 - \frac{m_j^2}{1-z} - \frac{m_i^2}{z}. \end{aligned}$$
(26)

Tables 1, 2, 3, 4, 5, 6 and 7 contain the collinear limits of all electroweak branching amplitudes. These limits are related to the Altarelli–Parisi splitting kernels by

$$\begin{aligned} P(\lambda _{ij}, \lambda _{i}, \lambda _{j} , z) = |M(\lambda _{ij}, \lambda _{i}, \lambda _{j})|^2 \end{aligned}$$
(27)

where M is the branching amplitude. Note that we used the replacement of Eq. (25) for the sake of notation, leading to a missing phase factor that is irrelevant for the computation of the splitting kernels. The splitting functions found here agree with the results of [32].

Table 1 Branching amplitudes for vector boson emission off a fermion. For the antifermion, the interchange \((v-\lambda a) \leftrightarrow (v+\lambda a)\) is applied
Table 2 Branching amplitudes for vector boson splitting to fermions
Table 3 Branching amplitudes for vector boson emission off a vector boson
Table 4 Branching amplitudes for Higgs emission off (anti)fermions
Table 5 Branching amplitudes for Higgs splitting to fermions
Table 6 Branching amplitudes for Higgs emission off a vector boson
Table 7 Branching amplitudes for Higgs splitting to vector bosons

The quasi-collinear limits of the branching amplitudes listed in Tables 1, 2, 3, 4, 5, 6 and 7 show a rich landscape of splitting modes that are the result of several physical effects. For instance, in the case of vector boson emission off a fermion, listed in Table 1, we find that the splitting functions for a transverse emission without a spin flip correspond with the familliar spin-summed form

$$\begin{aligned} P_{f \rightarrow f'V} \propto \frac{1 + z^2}{1-z}, \end{aligned}$$
(28)

in correspondence with the QCD equivalent of gluon emission off a quark. The presence of fermion and vector boson masses induces a shift of \(1/Q^2 \rightarrow {\tilde{Q}}^2 / Q^4\) in the propagator structure. The fermionic mass corrections in \({\tilde{Q}}^2\) also appear for gluon emission and reproduce the mass contributions to the eikonal factor in the soft limit. For W and Z emission, a vector boson mass correction is also present. In case of a transverse emission with a spin flip, the splitting functions are suppressed with a relative factor \(m^2/Q^2\), where m may the the vector boson mass or either of the fermion masses. Again, similar mass corrections appear in QCD [58], and cause such modes to be suppressed at off-shellness \(Q^2\) much larger than the electroweak scale.

The splitting function for longitudinal vector boson emission with a spin flip behaves like the scalar splitting function

$$\begin{aligned} P_{f \rightarrow f\varphi } \propto (1-z). \end{aligned}$$
(29)

In the unbroken Standard Model, this splitting mode indeed corresponds with the emission of a Goldstone boson, and the proportionality constant is given by the Yukawa coupling to the Goldstone. The spin flip mode is in this case not mass suppressed, reflecting the behaviour of scalar emissions in the unbroken phase. On the other hand, the mode without spin flip is mass suppressed, and terms that scale as \(1/m_j\) appear, originating from the scalar piece of the longitudinal polarization, as well as terms that scale as \(m_j\), corresponding with the vector piece.

5 The electroweak shower implementation

We use the branching kernels computed in the previous section to implement an electroweak shower in the Vincia framework, which was set out in [33,34,35,36,37, 59,60,61]. Here we first provide a brief summary before continuing with a description of some details specific to the electroweak shower.

Parton showers are constructed as a Markov chain of emissions that are distributed according to an approximation to the radiative matrix element and an associated Sudakov factor [62]. In most modern showers, these branchings are kinematically modelled as \(2 \rightarrow 3\) processes while the dynamics vary according to the shower model. If both branching partons are in the final state, momenta are labelled as \(I, K \rightarrow i,j,k\), and the Vincia shower approximation to the matrix element may be written as

$$\begin{aligned} |M_{n+1}|^2 d\varPhi _{n+1} = |M_n|^2 d\varPhi _{n} \times a(s_{ij}, s_{jk}) \, d\varPhi _{\mathrm{ant}}, \end{aligned}$$
(30)

where

$$\begin{aligned} d\varPhi _{\mathrm{ant}}^{\mathrm{FF}} = \frac{1}{16 \pi ^2} m_{IK}^2 \lambda ^{-\frac{1}{2}}(m_{IK}^2, m_I^2, m_K^2) ds_{ij} ds_{jk} \frac{d\varphi }{2\pi }, \end{aligned}$$
(31)

where \(s_{ij} = 2 p_i{\cdot }p_j\) and \(\lambda \) is the Källén function. Equation (31) represents an exact factorization of the radiative phase space. An associated kinematic map is defined between the pre-branching and post-branching momenta that conserves total momentum and is soft- and collinear-safe [59]. The branching kernels are so-called antenna functions \(a(s_{ij}, s_{jk})\) that capture the leading collinear and soft singularities associated with QCD emissions. The equivalent expressions to Eq. (30) for radiation from the initial state, as well as the definition of the kinematic maps, may be found in [60]. In the electroweak sector, we instead use the branching kernels computed in Sect. 3. These kernels contain the correct singular behaviour in the soft-collinear and collinear limits, and are simultaneously usable for the treatment of bosonic interference effects and recoiler selection to be discussed below.

Vincia supports QCD evolution of partons with definite helicity [36, 37], making for a natural framework for the inclusion of an electroweak shower. Interference effects between intermediate spin configurations have previously been incorporated in parton showers, such as in the Herwig [25, 63] parton shower [58] and in Deductor [64, 65]. The Vincia parton shower currently makes no attempt to incorporate such interferences, and correspondingly the same is true for the electroweak sector.

5.1 Ordering and resonance showers

The electroweak shower includes a number of branchings that would normally be associated with the decay of resonances by Pythia [66], for which Vincia is a plugin. In particular, the Standard Model particles that have such decay-like branchings are Z, \(W^{\pm }\), Higgs and top quark. With the inclusion of an electroweak shower, the decay modes of the resonances are now also all present as shower branchings. The shower enables highly-energetic resonances to branch and disappear at scales much higher than their width, where they should indeed be treated as any other non-resonant particle. At scales close to the resonance width, the Breit–Wigner character of the resonance decay should instead dominate the distribution.

Here, we set up a method of matching the parton shower distribution to a Breit–Wigner smoothly. However, a full-fledged incorporation of resonance decays at high scales in the Pythia event generator framework is a much more involved issue. In the absence of EW shower branchings, Pythia always associates the scale of a resonance decay with the width of the resonance. As the widths of all SM resonances are close to \(\varLambda _{\mathrm{QCD}}\), the approach of Pythia is to factorize the QCD shower of the hard scattering from showers in resonance decay systems. When resonances are instead able to decay at much larger scales, this approximation is no longer valid and a more general solution is required.

Most generally, resonances should be allowed to decay at any scale during the shower evolution of the hard scattering. Once a decay occurs, a factorized shower in the resonance system may be performed using Vincia’s resonance shower described in [61]. Afterwards, the showered decay products may be joined with the rest of the hard system and the evolution can be continued. As these modifications require substantal restructuring of both the Pythia and Vincia framework, they will be treated separately in [67].

Here, we do however proceed to set out the matching of the EW shower to the Breit–Wigner shape through a simple sampling procedure. To that end, we must first choose a suitable ordering scale. An issue specific to resonance branchings is that there are regions of phase space where the off-shellness \(Q^2\) is negative. As more negative values of off-shellness should correspond with shorter-lived resonances, we are led to define an ordering scale

$$\begin{aligned} |Q^2| = |s_{ij} + m_i^2 + m_j^2 - m_{ij}^2|, \end{aligned}$$
(32)

and equivalent for initial state branchings. For branchings that are not of the resonant decay type, \(Q^2\) is strictly positive and the ordering scale corresponds to off-shellness ordering, which indeed regulates the soft-collinear and collinear singularities.

Note that, even without the presence of EW masses, the choice of ordering scale is by no means unique and many shower algorithm employ different functional forms. The inclusion of EW scale masses leads to further ambiguities, since the virtuality or transverse momentum of the branching is not the only scale considered to be small in the quasi-collinear limit. In [27] the effects of shifting the ordering scale to account for the presence of gauge boson masses were found to lead to essentially identical results, motivating Eq. (32) as it leads to simpler matching to the Breit–Wigner shape.

For most types of non-resonant electroweak branchings the phase space is naturally cut off due to the masses of the post-branching momenta. Beyond resonance branchings, photon emission is the only remaining branching that is not cut off naturally. They are instead cut off at the same scale as the QCD shower approximately equal to \(\varLambda _{\mathrm{QCD}}\), and QED radiation at lower scales is included as is described in [57].

We set up the Breit–Wigner distribution by making the replacement

$$\begin{aligned} \frac{1}{Q^4} \rightarrow \frac{1}{Q^4 + m^2 \varGamma ^2}, \end{aligned}$$
(33)

in the relevant branching kernels computed in Sect. 3, where \(\varGamma \) is the width of the resonance. These kernels are then normalized to represent a probability distribution from which off-shellness scales, kinematics and post-branching spin states may be sampled. Because the ordering scale of the shower is given by Eq. (32), it is straightforward to define a matching scale \(|Q^2_{\mathrm{Match}}|\) where the shower is matched to the Breit–Wigner. However, while it may be possible to pick the matching scale to ensure the distribution is continuous, it will in general not be smooth. Furthermore, the shape of the shower distribution at scales close to the resonance width depends on the starting scale and the antenna mass, so such a choice would in any case not always be continuous. We thus ensure the distribution is smooth by sampling the value of the matching scale for every resonance from the distribution

$$\begin{aligned} {\mathcal {P}}(|Q^2_{\mathrm{Match}}|) \propto \frac{m^2 \varGamma ^2 |Q_{\mathrm{Match}}^2|}{(Q_{\mathrm{Match}}^4 + m^2 \varGamma ^2)^2}. \end{aligned}$$
(34)

The shower contribution to the off-shellness spectrum is then effectively multiplied by a factor

$$\begin{aligned} \int _0^{|Q^2|} \! \! \! \! \! \! \! \! d|Q^2_{\mathrm{Match}}| \, {\mathcal {P}}(|Q^2_{\mathrm{Match}}|) = \frac{Q^4}{Q^4 + m^2 \varGamma ^2}, \end{aligned}$$
(35)

which ensures a smooth suppression of the shower kernel. In this treatment, the total distribution is dominated by the shower at hight scales, while the Breit–Wigner dominates at low scales.

The implementation within the parton shower framework is relatively straightforward. When a shower of a resonance is initiated, a matching scale \(|Q^2_{\mathrm{Match}}|\) is sampled from Eq. (34), serving as a local cutoff scale. If during showering the branching scale drops below the matching scale, a new off-shellness scale is instead drawn from the Breit–Wigner distribution. We emphasize that this solution serves as an approximate means of matching the shower to a Breit–Wigner, and a more sophisticated method that closely matches that of Pythia will be developed in [67].

5.2 Recoiler selection

While in the QCD portion of the Vincia shower the colour structure dictates the pairings of branching partons I and K, no such guidance exists in the electroweak sector. Furthermore, the branching kernels only describe the soft-collinear and collinear singularities associated with the branching of particle I. The choice of recoiler K does not contribute to the accuracy of the shower in these limits. It may however be chosen probabilistically to minimize the physical consequences of the recoil on previously generated branchings. For the branching of particle i, the probability to select a spectator j from a pool of N available ones is

$$\begin{aligned} {\mathcal {P}}_{ij} = \frac{\left| M^{x \rightarrow i j} \right| ^2}{\sum _{j'=1}^N \left| M^{x \rightarrow i j'} \right| ^2}. \end{aligned}$$
(36)

That is, all candidate spectators j have a probability to be assigned as a recoiler for particle i if the pair ij may have been produced by an electroweak branching. All of the contributions in the denominator of Eq. (36) thus correspond with possible shower histories that contribute to the current state. The selection of a spectator is then more likely if the shower history where the current brancher and that spectator were created by a previous branching.

We clarify the choice of Eq. (36) through an example, which may be expressed in terms of diagrams as

(37)

The spectator for the splitting of the vector boson is chosen to be either of the other external legs based on the probabilities that the vector boson was emitted by either of those legs. When the vector boson splits, it is brought off its mass shell by transferring some momentum of the spectator to the vector boson. Because the vector boson momentum is modified, the emission kernel that it was produced with is no longer correct. In the strong-ordering phase space region where \(Q^2_{\mathrm{emit}} \gg Q^2_{\mathrm{split}}\) this type of mismodelling is absent. However, the shower covers all of phase space, including regions where the scales \(Q^2_{\mathrm{emit}} > Q^2_{{\mathrm{split}}}\) are of comparable size. Recoil effects of previous brancings may be especially relevant for branchings that involve masses of the order of the electroweak scale.

The Vincia \(2 \rightarrow 3\) kinematic map conserves the invariant mass of the original two-particle system. This means that the probability of Eq. (36) ensures the propagator structure of the emitter and vector boson pair that was most important in the emission process is most often conserved. Note that a similar procedure has in the past been used to select a recoiler for gluon splitting [55, 59], which also features only collinear singularities. In [32] this effect is referred to as ‘kinematic back-reactions‘ and is accounted for as a multiplicative factor of the branching kernels.

5.3 Bosonic interference

A unique type of interference effect appears in the electroweak sector due to the existence of multiple neutrally charged bosons [68]. It is possible to treat such interference effects comprehensively by evolving density matrices that contain mixed states of neutral bosons [32], but such procedures quickly become computationally prohibitive. We instead opt for a simpler approach that incorporates the most imporant physical effects while preventing the shower to become a bottleneck in the Monte Carlo event generation chain [69].

Interference between neutral bosons occurs when the electroweak shower produces a neutral boson, and it subsequently disappears by splitting. The heavy neutral bosons will always split due to the matching to a Breit–Wigner distribution described in Sect. 5.1, but photons may survive the showering procedure. As such, interference effects are corrected for by applying an event weight

$$\begin{aligned} w_{\mathrm{BI}} = \sum _x \frac{|M^{x \rightarrow x' b_1} M^{b_1\rightarrow i j} + M^{x \rightarrow x' b_2} M^{b_2\rightarrow i j}|^2}{|M^{x \rightarrow x' b_1} M^{b_1\rightarrow i j}|^2 + |M^{x \rightarrow x' b_2} M^{b_2\rightarrow i j}|^2} \end{aligned}$$
(38)

after the splitting, where we have dropped the helicity indices for readability. Equation (38) corrects the branching kernels for the interference between two bosons \(b_1\) and \(b_2\). These bosons may be either a photon and a transversely polarized Z boson, or a Higgs boson and a longitudinally polarized Z boson. Interference between spin states is not taken into account in correspondence with the rest of the shower algorithm. The weight of Eq. (38) then sums over all possible emitters x of the bosons \(b_1\) and \(b_2\), accounting for all possible shower histories. This includes the different spin states of the same particle, but note that these contributions are summed over incoherently. The weight has the property \(0 \le w_{\mathrm{BI}} \le 2\), and since the overal rate of electroweak boson emissions is moderate even at high energies, there is little danger of wildly fluctuating weights leading to inefficiencies.

5.4 Overestimate determination

The implementation of electroweak radiation in the shower formalism through the Sudakov veto algorithm [70,71,72] requires finding overestimates of the associated branching kernels. Due to the large number of types of branchings in the electroweak sector, it is desirable to automate this procedure.

For final-state particles with an electroweak charge, a recoiler is selected through the procedure outlined in Sect. 5.2. The branching phase space is then given by Eq. (31). All final-state electroweak branching kernels are overestimated by a parameterized function

$$\begin{aligned} {\mathcal {O}}^{\mathrm{FF}}&= c_1^{{\mathrm{FF}}} \frac{1}{|Q^2|} + c_2^{\mathrm{FF}} \frac{1}{|Q^2|} \frac{E_{IK} (E_{IK} + |\vec {p}_{IK}|)}{s_{ij} + s_{ik} + m_i^2} \nonumber \\&\quad + c_3^{{\mathrm{FF}}} \frac{1}{|Q^2|} \frac{E_{IK} (E_{IK} + |\vec {p}_{IK}|)}{s_{ij} + s_{jk} + m_j^2} + c_4^{{\mathrm{FF}}} \frac{m_I^2}{Q^4}, \end{aligned}$$
(39)

where \(|Q^2|\) is given by Eq. (32). The term multiplying \(c_1^{{\mathrm{FF}}}\) reflects the \(1/|Q^2|\) behaviour of most branching kernels, while the second and third terms incorporate the soft behaviour associated with vector boson emission, containing ratios \(E_{IK}/E_i \sim 1/z\) and \(E_{IK}/E_j \sim 1/(1-z)\) in terms of shower variables. The term multiplying \(c_4^{{\mathrm{FF}}}\) represents the mass corrections that may be present for massive branchers. The contribution of post-branching masses are typically negative, and therefore do not improve the overestimate much.

Initial-state branchings are only allowed to recoil against other initial states. In this case, the antenna phase space is

$$\begin{aligned} d\varPhi _{{\mathrm{ant}}}^{{\mathrm{II}}} = \frac{1}{16\pi ^2} \frac{x_A^2}{x_a^2} \frac{x_B^2}{x_b^2} \frac{1}{s_{AB}} ds_{aj} ds_{bj} \frac{d\varphi }{2\pi }, \end{aligned}$$
(40)

where A branches to a and j, and B is the recoiler. The electroweak shower currently only implements vector boson emission from fermions in the initial state, which are treated as massless by Vincia. The ordering scale is crossed into the initial state to give

$$\begin{aligned} Q^2 = s_{aj} - m_j^2. \end{aligned}$$
(41)

An absolute value qualification is not required here since resonance type branchings do not occur in the initial state. A sufficient overestimate is

$$\begin{aligned} {\mathcal {O}}^{{\mathrm{II}}}&= c_1^{{\mathrm{II}}} \frac{1}{Q^2}\frac{s_{ab}}{s_{AB}} \nonumber \\&\quad + c_2^{{\mathrm{II}}} \frac{1}{Q^2} \frac{x_A^2 s_{ab}^2}{x_A s_{bj}(s_{ab} - s_{bj}) + x_B s_{aj}(s_{ab} - s_{aj})}. \end{aligned}$$
(42)

The factor \(s_{ab}/s_{AB}\) accounts for the additional factor of 1/z that shows up in the Altarelli–Parisi splitting kernels when crossed to the initial state. The second term represents the soft enhancement \(1/(1-z)\) that may appear for vector boson emissions.

The parameters \(c_1^{{\mathrm{FF}}}\) through \(c_4^{{\mathrm{FF}}}\), \(c_1^{{\mathrm{II}}}\) and \(c_2^{{\mathrm{II}}}\) are automatically determined for all possible branchings in the electroweak shower. To do that, brancher-recoiler pairs are generated from antennae with randomly chosen invariant masses. Branchings are then generated with a distribution \(1/|Q^2|\) for the final state or \(s_{ab}/s_{AB} \, 1/Q^2\) for the initial state to roughly model the branching kernel behaviour. For every event i, the value of the branching kernel \(B_i\) as well as the terms \(A_{ij}\) multiplying the parameters \(c_j\) are stored. The problem of finding suitable values for the overestimate parameters can then be formulated as

$$\begin{aligned}&\text{ Minimize } \sum _{i=1}^n (A{\mathbf {c}})_i - B_i \nonumber \\&\text{ subject } \text{ to } (A{\mathbf {c}})_i \ge B_i \text{ and } {\mathbf {c}} \ge 0. \end{aligned}$$
(43)

The minimization condition minimizes the average difference between the branching kernel and its overestimate. The constraints ensure the overestimate is larger than the branching kernel for all samples and the parameters are positive definite. Equation (43) is an instance of a linear programming problem, for which many libraries are available. We make use of the Python [73] package PuLP [74].

Fig. 1
figure 1

Branching spectra of a 1 TeV \(\tau _{-}\) to \(\tau \gamma \) (left) and \(\tau Z\) and \(\nu _{\tau } W^{-}\) (right)

Fig. 2
figure 2

Branching spectra of a 1 TeV top quark to \(b W^+\) (left), tZ (right) and \(t \gamma /h\) (bottom)

Fig. 3
figure 3

Branching spectra of a 1 TeV \(W^+_+\) to fermions (left), \(W^+Z\) (right) and \(W^+ \gamma /h\) (bottom)

Fig. 4
figure 4

Branching spectra of a 1 TeV Higgs boson to fermions (left) and VV/hh (right)

Fig. 5
figure 5

Differential rates of the shower histories \(e_- \rightarrow e_- \gamma / Z_T \rightarrow e_- X\) (top) and \(e_+ \rightarrow e_+ \gamma / Z_T \rightarrow e_+ X\) (bottom) without (left) and with (right) the bosonic interference correction. The showers are initiated from a 10 TeV electron with a neutral recoiler

Fig. 6
figure 6

Electroweak shower approximation of electroweak virtual corrections to exclusive dijet production and exclusive \(W^+\) + jet production at center-of-mass energy \(\sqrt{s} = 14\) TeV as a function of the transverse momentum of the hard scattering process \(p_{\perp }^{{\mathrm{Hard}}}\). The solid line corresponds to the parton shower prediction without the QCD shower, while the dashed line shows the effect of interleaving with QCD radiation

Fig. 7
figure 7

Average number of weak boson emissions in exclusive dijet production and exclusive \(W^+\) + jet production at center-of-mass energy \(\sqrt{s} = 14\) TeV as a function of the transverse momentum of the hard scattering process \(p_{\perp }^{{\mathrm{Hard}}}\)

5.5 Overview of the shower algorithm

We conclude this section with a short description of the complete shower algorithm. Branching kernels are constructed using the formalism described in Sect. 3 for all possible electroweak branchings and all helicity configurations. Overestimates are found using the optimization algorithm of Sect. 5.4, where the post-branching helicities are summed over. This leaves a total of 277 types of final-state branchings, of which 74 are resonance decays, and 90 types of initial-state branchings.

As the shower initializes, a recoiler is selected for all final-state particles that have electroweak charge, making use of the selection probability described in Sect. 5.2. Initial-state branchers are always paired with the other initial-state particle, recoiling against the entire event.

While the shower runs, electroweak branchings compete against the QCD branchings generated by Vincia. The overestimates are used to generate trial branchings which are accepted or rejected through the usual Sudakov veto algorithm. For resonance decay branchings, the procedure outlined in Sect. 5.1 is used to match the shower to a Breit–Wigner distribution. We make use of the same kinematic maps as Vincia, and first-order running of the electroweak coupling constant is incorporated as part of the veto procedure. Wherever applicable, the bosonic interference weight Eq. (38) is included. After accepting a branching, a helicity state is selected with probability

$$\begin{aligned} {\mathcal {P}}_{\lambda _I, \lambda _i, \lambda _j} = \frac{B_{\lambda _I, \lambda _i, \lambda _j}(p_I, p_i, p_j)}{\sum _{\lambda _i,\lambda _j}B_{\lambda _I, \lambda _i, \lambda _j}(p_I, p_i, p_j)}, \end{aligned}$$
(44)

and equivalent for resonance decay branchings. The QCD shower and the electroweak shower run interleaved until the QCD cutoff scale is reached, after which only QED radiation is simulated.

6 Results

The inclusion of electroweak radiation in parton showers opens up a rich field of showering phenomena, in particular at energies well above the electroweak scale. Rather than considering specific sets of observables for specific processes, in this section we first consider radiation spectra of several particles at energies compatible with future colliders, showing overal branching rates of the electroweak shower, the relative importance of the large number of emission modes and the effects of bosonic interference. We conclude by showing results for hadron collision processes at LHC energies and comparing with the Pythia electroweak shower [27].

6.1 Branching spectra

In this section, we show electroweak branching spectra of several highly energetic particles as a function of the invariant mass, which closely corresponds with the ordering scale. Figures 12, 3 and 4 show this for the first branching of a left-handed \(\tau \) and top, a transverse \(W^+\) boson and a Higgs. Note that multiple emission rates are related to these single branching spectra through the usual generalized Poisson statistics associated with shower algorithms due to the multiplicative property of the Sudakov factor.Footnote 1 All particles are produced at an energy of 1 TeV together with a recoiler that is uncharged under electroweak interactions. For photon emissions, a cutoff around \(\varLambda _{{\mathrm{QCD}}}\) is imposed. All other branchings are automatically regulated by the particle masses or by matching to the Breit–Wigner distribution.

Figure 1 shows the branching spectrum of a negative-helicity \(\tau \). The two dominant photon production channels are those where the \(\tau \) helicity is conserved. The mass-suppressed spin-flip mode only contributes at very small invariant masses, as is to be expected from the branching kernel behaviour of \(m_{\tau }^2/Q^4\). The other spin-flip mode is highly suppressed in the collinear limit as is indicated in Table 1. For the emission of other vector bosons, the spin-flip contributions do not become sufficiently enhanced to show up before the kinematic limit is reached. The longitudinal vector boson emission channels have a characteristic form which looks very similar for the \(W^{-}_{0}\) and the \(Z_0\) channels, and which becomes comparable to the transverse channels at scales close to the kinematic limit.

Figure 2 shows the branching spectrum of a negative-helicity top. The left graph displays the resonance branchings as generated by the sampled matching procedure outlined in Sect. 5.1. The right and bottom graphs show all other branchings that are not of the resonance decay type. Spin-flip modes now show up for \(t \rightarrow b W^{+}\), \(t \rightarrow t Z\) and \(t \rightarrow t \gamma \) due to the large top mass, and they show the expected \(m_t^2/Q^4\) scaling with the emission scale. The ‘natural‘ mode of spin-flip Higgs emission is relatively flat compared with the fermion mass scaling mode of Higgs emission without spin flip.

Figure 3 shows the branching spectrum of a transverse \(W^{+}\). Resonance peaks only appear for decays to negative-helicity states due to their small masses. The branchings \(W^+_+ \rightarrow t {\bar{b}}\) with a spin-flipped top do occur on the other hand. The \(W^+_+ \rightarrow W^+ Z\) and \(W^+_- \rightarrow W^+ \gamma \) channels are dominated by the all-positive helicity configuration because of its \(1/z(1-z)\) scaling in the collinear limit as can be seen in Table 3. The modes to opposite transverse helicities are almost identical for the \(W^+ Z\) channels due to symmetry in the collinear limit and almost identical mass, but they are widely different for low scales in the \(W^+ \gamma \) channels. This is caused by the \(z^3/(1-z)\) and \((1-z)^3/z\) scaling of the collinear limits, where the photon can attain a very small collinear momentum fraction while that of the \(W^+\) is constrained by its mass. The single-longitudinal channels in \(W^+ Z\) are also almost identical for very similar reasons. The \(W^+_0 Z_0\) is a mode that is related to the Goldstone bosonic part of the \(W^+\) and Z, and it can be seen to be very similar to the \(W_0^+ h\) channel. On the other hand, the \(W^+_+ h\) mode differs significantly from the \(W^+_+ Z_0\) channel because it is dominated by the vectorial part of the longitudinal polarization.

Figure 4 shows the branching spectrum of a Higgs. The only significant resonance decay channels are \(b_{\pm } {\bar{b}}_{\pm }\) and \(\tau _{\pm } {\bar{\tau }}_{\pm }\) as may be expected due to the coupling to the fermion mass and the Higgs spin zero nature. On the other hand, the mass-suppressed \(t_{\pm } {\bar{t}}_{\mp }\) channel is comparable with the natural \(t_{\pm } {\bar{t}}_{\pm }\) channel. All channels to \(W^+ W^-\) and ZZ are almost identical since their branching kernels only differ in the gauge boson mass and a factor of \(1/c_w\) in the coupling. Also included is the \(h \rightarrow hh\) cubic Higgs coupling which is proportional to the Higgs mass \(m_h\), or equivalently the Higgs self-coupling \(\lambda \). This is the only branching where it makes an appearance, and it can be seen to provide a significant contribution to the total branching rate.

6.2 Bosonic interference

We now consider the effect of the application of the bosonic interference factor described in Sect. 5.3. Figure 5 shows rates for the shower histories \(e_- \rightarrow e_- \gamma / Z_T \rightarrow e_- X\) and \(e_+ \rightarrow e_+ \gamma / Z_T \rightarrow e_+ X\) using a similar setup as in the previous subsection, but starting from a 10 TeV source electron. Multiple interesting features appear when the bosonic interference weight Eq. (38) is included. The most striking difference occurs for the \(W^+ W^-\) channel, where the bosonic interference causes an increase in case of the \(e_-\), but a major decrease in case of the \(e_+\). This may be understood by considering the structure of the interfering branching amplitudes. Factoring out coupling constants and other kinematic components, the interference is proportional to

$$\begin{aligned} \frac{1}{M_{WW}^2} + \frac{c_w}{s_w} \frac{1}{4 s_w c_w} (1 - 4 s_w^2 - \lambda _e) \frac{1}{M_{WW}^2 - m_z^2 + i m_z \varGamma _z}, \end{aligned}$$
(45)

where the factor \(c_w/s_w\) comes from the ZWW-coupling and \(\lambda _e\) is the electron helicity. The second term in brackets interferes destructively with the photon contribution for sufficiently large values of \(M_{WW}\), and the remaining terms in the Z contribution cancel for \(\lambda = 1\).

The effects of the bosonic interference factor on charged fermion rates close to the Z peak may be understood through a similar argument. The rates close to the Z peak are significantly affected by the simplified and preliminary method of matching to resonance decays as described in Sect. 5.1, and will be improved upon in [67].

6.3 Electroweak corrections to proton collision processes

We finally consider the parton shower predictions of electroweak corrections to some common proton collision processes at LHC energies and compare with the Pythia electroweak shower [27]. Since the weak vector bosons produced by the electroweak shower at high energies are massive and thus observable, they may provide a rich environment for phenomenological studies including kinematic effects on the hard scattering, jet substructure due to vector boson decay inside jet cones and external high-energy jet and lepton production.

With the goal of examining the general significance of electroweak Sudakov effects in common LHC processes, we generate dijet and \(W^{+}\) plus jet events at \(\sqrt{s} = 14\) TeV using the default tune of Pythia 8.2 [23] and the NNPDF2.3 sets [75]. Figure 6 shows the approximate electroweak virtual corrections as predicted by the Pythia electroweak shower, which only incorporates vector boson emission from fermions, and the Vincia electroweak shower as a function of the transverse momentum of the hard scattering. The virtual corrections may be estimated by counting the events that contain at least one weak vector boson emission. The probability for the shower to produce no additional weak bosons is given by the Sudakov factor

$$\begin{aligned} \varDelta _{{\mathrm{EW}}} = 1 - {\mathcal {O}}(\alpha ), \end{aligned}$$
(46)

and thus the \({\mathcal {O}}(\alpha )\) corrections are given by the probability for at least one weak boson emission. Virtual corrections to these processes were calculated in for example [14, 15] for exclusive dijet production and in [17, 18] for vector boson production.

For dijet production, the results of the showers are very similar. In the case of \(W^+\) plus jet production the substantial difference between the showers is caused by the absence of the Yang–Mills vector boson coupling in the Pythia shower. Furthermore, we find that the contribution to the weak boson emission rates of the initial-state quarks is significantly smaller than that of the final state. Because at large x and high scales the PDFs are predominantly quark-like and the hard scattering is dominated by \(q {\bar{q}}' \rightarrow W^+ g\), and the phase space for initial state radiation is small at large x, virtual corrections decrease at large values of transverse momentum.

Also shown in Fig. 6 are the results of interleaving the electroweak shower with the QCD showers of Pythia and Vincia. In the strongly ordered limit shower branchings are unaffected by subsequent branchings, but subleading effects due to the kinematics and the creation of weakly charged quarks still lead to minor differences.

Figure 7 shows the average number of weak boson emissions of the showers. In \(W^+\) plus jet production the first purely vector boson branching is always \(W^+ \rightarrow W^+ Z\) explaining the large increase in the Z boson emission rates. Similarly, Pythia’s \(W^{+}\) rate is small since the final-state quark is always down-type. The increase in the Vincia shower is thus caused entirely by secondary emissions from prior weak vector boson emissions.

7 Conclusion and discussion

The effects of weak corrections in parton shower algorithms are known to become significant already at LHC energies, in particular with the upcoming luminosity upgrade, and will be even more relevant at future colliders. One of the major challenges of the construction of such a shower is the calculation of the relevant branching kernels, which in this paper was done using the spinor-helicity formalism. Compared with QCD, the electroweak theory involves many theoretical subtleties that have to be handled carefully. One major issue is the chiral nature of the electroweak theory, which forces the shower to be helicity-dependent and leads to a large number of possible types of branchings. In particular, the scalar components of longitudinal polarizations lead to unphysical, unitarity-violating contributions that have to be treated carefully. The collinear limits of the computed branching kernels are found to be in agreement with the results of [32]. The electroweak shower also includes many branchings that would usually be considered to be decays of resonances, in which case the distribution follows a Breit–Wigner peak. A strategy to match the parton shower to a resonance decay was proposed, but this may likely be improved upon by a better understanding of the interplay between the virtual corrections contained in the Sudakov factor and the decay width. A more sophisticated treatment of this matching is beyond the scope of this paper, and will be the topic of future study [67]. Further electroweak effects added to the shower include a recoiler selection procedure that compensated for recoiler effects of previous branchings and treatment of bosonic interference effects. Results were shown that quantify the general size of electroweak shower corrections at future collider energies and at LHC energies.

Several features that are currently lacking from the electroweak shower were already pointed out. They include topics such as soft and spin interference effects, although such issues have similarly not been fully solved in the QCD sector of commonly used shower codes. However, further issues particular to the electroweak sector still remain. One is the inclusion of the CKM quark-mixing matrix [76,77,78], which would lead to an even larger number of possible electroweak branchings. The impact of the CKM matrix is however not expected to be large since the off-diagonal terms of the third-generation row and column are close to zero. The other off-diagonal terms are not as small, but they mix quarks that are treated as massless in Vincia anyway.

One other peculiar property of the electroweak theory is the appearance of Bloch–Nordsieck violation [10, 79]. The parton shower formalism is fundamentally based on the principle of unitarity and the cancellation of infrared divergences between real and virtual corrections. Since the electroweak vector bosons are massive, divergences associated with their emission are mass-regulated. The flavour-changing nature of W-boson emission from the initial state spoils the exact cancellation of the infrared divergences and some mass-regulated logarithms may be left-over.

There is no straightforward method to incorporate these violations in the shower formalism, since they explicitly break unitarity. We note that, while Bloch–Nordsieck violations are not particularly significant at the LHC [10, 15], they will be important at future collider energies and a comprehensive treatment will be necessary.

Finally, hard processes initiated by vector bosons have been considered for a long time [80, 81]. PDF sets with QED corrections have been available for some time [82,83,84,85], and recent progress was made towards PDFs with complete electroweak corrections [86,87,88]. The current shower implementation only allows for the emission of vector bosons from the initial state. The calculation of the other required initial-state branching kernels is in principle as straightforward as the calculation of those available already, but an implementation in the Pythia framework is likely not simple.