Anomalous coupling, top-mass and parton-shower effects in W+W− production

We calculate the process pp→W+W−→e+νeμ−ν¯μ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ pp\to {W}^{+}{W}^{-}\to {e}^{+}{\nu}_e{\mu}^{-}{\overline{\nu}}_{\mu } $$\end{document} at NLO QCD, including also effective field theory (EFT) operators mediating the ggW+W− interaction, which first occur at dimension eight. We further combine the NLO and EFT matrix elements produced by GoSam with the Herwig7/Matchbox framework, which offers the possibility to study the impact of a parton shower. We assess the effects of the anomalous couplings by comparing them to top-mass effects as well as uncertainties related to variations of the renormalisation, factorisation and hard shower scales.


Introduction
Among the most important goals for the next phase of LHC data taking are precision tests of the electroweak symmetry breaking sector and the search for signs of new physics. In this respect, the final state of W + W − plays a prominent role. For example, the continuum pp → W + W − → llνν is the dominant background in the measurement of H → W + W − → llνν. The process pp → W + W − (+jets) also can be a major background for new-physics processes involving missing energy. Therefore it is very important to have good theoretical control on the pp → W + W − cross section.
Final states with two massive vector bosons recently have attracted additional interest due to the fact that a slight excess at about 2 TeV in the search for di-boson resonances has been reported by both ATLAS [1, 2] and CMS [3,4], most pronounced in the hadronic decay channel, which however does not seem to persist in Run II. Further, both the ATLAS and CMS measurements for the W + W − total inclusive cross sections using the leptonic decay channels, at 7 TeV [5, 6] and at 8 TeV [7,8], are about 10-20% higher than the NLO predictions obtained from MCFM [9, 10], which include the gg-initiated sub-process [11]. However, the latter discrepancy has been largely reduced by the NNLO predictions which became available recently [12][13][14]. In addition, it has been noticed that resummation of large logarithms arising from the jet veto condition needs to be taken into account carefully [15][16][17][18][19], and that the discrepancy for the fiducial cross section is only at the 1σ level, such that -1 -JHEP05(2016)106 the effective gg W + W − coupling could affect various distributions. To this aim we include the most general effective field theory (EFT) operators mediating the gg W + W − interaction, which first occur at dimension eight, in our automated setup. This allows us to assess the impact of these operators in various effective coupling scenarios. We refrain from the discussion of dimension six operators involving Higgs bosons, modifying e.g. the ggH, ttH or HW + W − couplings, because they have been discussed extensively in the literature, and mostly are severely constrained by LHC measurements already, see e.g. [59][60][61][62][63][64] and references therein. The operators mediating an effective gg W + W − coupling can be studied in isolation from the dimension 6 operators involving Higgs bosons because the two classes of operators do not mix at the order considered here. The diagrams contributing to the process gg → e + ν e µ −ν µ (see figure 1) comprise of diagrams involving two W propagators ("doubly resonant") as well as diagrams involving only one W propagator ("singly resonant", see figure 1b). Note that the latter are important to maintain gauge invariance. Non-resonant diagrams, i.e. diagrams containing no W propagator, do not contribute to the e + ν e µ −ν µ final state, but they would contribute to the e + ν e e −ν e final state, which we are not considering here. We include massive top-and bottom-quark loops, all other quarks (and leptons) in the hard process are assumed to be massless. The photon-exchange graphs vanish due to Furry's theorem. The Z-exchange diagrams are proportional to (m 2 u − m 2 d )(p 2 3 − p 2 4 ), when summed over up-and down-type contributions. Therefore these diagrams also vanish for massless quarks. For arbitrary invariant masses of the charged lepton -neutrino pairs, i.e. p 2 3 = p 2 4 , and the third generation, where we assume m t = 0 and m b = 0 in the loops, we find that the contributions from doubly-resonant (figure 1a) and singly-resonant (figure 1b) diagrams with internal Z-propagator cancel each other. The only triangle graphs that contribute are thus the Higgs-exchange diagrams (figure 1c) where the amplitude contains the QQH Yukawa couplings. The box diagrams do not involve these couplings and therefore form a gauge-invariant subset.

Operators parametrizing the gg W + W − coupling
The first operators that mediate a four-boson interaction between the two gluons and the two W bosons occur at dimension eight.
The gluonic field-strength tensor is defined as the field-strength tensor of the W is defined as We can write the SU(2) fields W I in terms of the physical fields: A CP-odd operator can be introduced via the dual field-strength tensors which are defined byG a,µν = 1 2 µνρσ G a ρσ ,W I,µν = 1 2 µνρσ W I ρσ . (2.4) Based on these field-strength tensors we can build the dimension-eight operators contributing to the Lagrangian: where Λ denotes the scale below which the EFT description is valid.
Combining the SM Lagrangian with the one including the effective couplings, we have Figure 2. The three types of contributions to the squared matrix element.
Combining higher-order QCD corrections from L SM with the part containing the higherdimensional operators means that we are performing a simultaneous expansion in α s /2π and in c i /Λ 4 . This requires a careful assessment of the relative importance of the various terms in such an expansion and of the range of validity of the effective theory. We will come back to this issue in section 3.
As in the SM the gg W + W − coupling is loop induced, we will calculate the following contributions to the gg (→ W + W − ) → e + ν e µ −ν µ cross section, depicted schematically in figure 2: Note that the last term above is suppressed by 1/Λ 8 . In the following we list the Feynman rules corresponding to the dimension-eight operators. All momenta are considered to be incoming:

Loop-induced processes in GoSam
The virtual amplitudes for the process pp (→ W + W − ) → e + ν e µ −ν µ , as well as the spin-and colour-correlated tree amplitudes, have been generated with the program GoSam [56,57], which is an automated package to generate one-loop amplitudes. It is based on a Feynman-diagrammatic approach using Qgraf [65] and Form [66,67] for the diagram generation, and Spinney [68] and Form to produce optimised Fortran90 code. For the reduction of the one-loop amplitudes the user can choose between three different reduction libraries (or a combination thereof): Ninja [69][70][71], Golem95 [72][73][74] or Samurai [75,76], where Samurai and Ninja are based on integrand reduction techniques [77][78][79][80], while Golem95 uses tensor reduction [81] and tensorial reconstruction [80]. The default setup uses Ninja in combination with Golem95 as a rescue system for numerically problematic phase-space points. The scalar basis integrals have been evaluated using Golem95 and OneLOop [82]. We use the complex mass scheme [83] throughout our calculation, which implies that a non-zero top quark width is included in the propagators forming the top quark loops. For numerical values of masses and widths we refer to section 3.

JHEP05(2016)106
The calculation of loop-induced processes is straightforward with GoSam, as it is based on a Feynman-diagrammatic approach. However, we have improved the rescue system for loop-induced processes: in the standard case where a tree-level amplitude exists, the accuracy at a given phase-space point is assessed by comparing the value for the coefficient of the single infrared (IR) pole with the exact value (which is obtained from the universal IR behaviour of the subtraction terms). As this coefficient must be zero in the loop-induced case, we have implemented an accuracy check which tests the precision of this zero. In more detail, we have added the following options to the GoSam input card: • PSP_chk_li1: allows to set the desired precision of the pole part (which should be zero) in comparison to the finite part. If the pole part is at least PSP_chk_li1 orders smaller than the finite part, the point is accepted.
• PSP_chk_li2: for loop-induced processes, this option is used instead of PSP_chk_th2. It is the threshold to declare a phase-space point as "bad", based on the precision of the pole in comparison to the finite part. Points with precision less than this threshold are directly reprocessed with the rescue system (if switched on), or declared as unstable.
According to the verbosity level set, such points are written to a file and not used when the code is interfaced to an external Monte Carlo program in accordance with the updated BLHA standard [84].
• Similarly, PSP_chk_li3 is used instead of PSP_chk_th3 as threshold for the rotation test in the loop-induced case.
• PSP_chk_li4 sets the minimum pole precision for the points which already have been reprocessed by the rescue system. If a rescued point gives a pole coefficient which is at least PSP_chk_li4 orders smaller than the finite part, the point is accepted.
For the Standard-Model case, we have validated our results by comparing to MCFM [11].

Extended BSM support in GoSam
We have implemented various new features in GoSam which facilitate the calculation of corrections beyond the Standard Model, as well as the interference between SM loop corrections and BSM effects described within an EFT framework. Among the new features are: • the import of BSM model files in UFO format [85] in combination with the updated BLHA standard [84] for the definition of new couplings has been implemented, • BSM-SM interference terms can be calculated by specifying the orders in the corresponding couplings, • sub-process specific settings in the GoSam input card are possible, • the coefficients of the effective couplings multiplying the higher-dimensional operators can be modified by the user in the GoSam input card.

JHEP05(2016)106
For the process considered here we use a model file in UFO format [85] which we create by extending the SM Lagrangian of FeynRules [86] with the EFT operators outlined in section 2.2. The SM parameters, which are provided per default with the UFO model file, are then overwritten in accordance to the updated BLHA standard [84] by the ones listed in section 3. In order to be able to compute the pure SM NLO QCD contributions as well as the additional interference between the SM one-loop contribution and the effective coupling with only one model file, we restrict the one-loop contributions to allow only for SM couplings, with the help of GoSam's diagram filter facilities, while the tree-level contributions are allowed to include the effective coupling.

Interface to Herwig7 and computational setup
Herwig7 features the full simulation of particle collision events up to the particle level, i.e. perturbative as well as non-perturbative physics. The perturbative part provides the simulation of hard processes at NLO QCD (including several built-in LO and NLO matrix elements, LH event file input as well as the fully automated assembly of NLO QCD calculations for almost all Standard-Model processes, utilizing various interfaces to several external matrix-element providers), shower Monte Carlo algorithms, as well as the corresponding LO and NLO matching procedures (dedicated matrix element corrected shower plug-ins and built-in matched cross sections, as well as a fully automated matching machinery). For details we refer to the Herwig7 release notes [34].
The capabilities of Herwig7 for integration, unweighting and sub-process parallelization, as well as the steering at the level of input files, are significantly improved. By virtue of the Matchbox framework new integrator modules were introduced, which provide for an efficient, automated multi-channel phase-space sampling: the one which we employ, for the process considered here, is based on the standard sampling algorithm contained in the ExSample library [87].
Based on the Matchbox framework, Herwig7 facilitates the automated setup of all sub-processes and ingredients necessary for a full NLO QCD calculation in the subtraction formalism: an implementation of the dipole subtraction method based on the approach by Catani and Seymour [88], interfaces to various external matrix-element providers, or plug-ins to various in-house calculations for the hard sub-processes or to built-in colour and helicity sub-amplitudes.
In the case where the necessary one-loop and tree-level parts are obtained from external matrix-element providers there exist several possibilities: either at the level of colourordered sub-amplitudes or at the level of squared matrix elements through the updated BLHA standard interface [84]. A more detailed description of the interface between Her-wig++ and GoSam has been given in [89], for the example of Z+jet production. The interface between Herwig7/Matchbox [34,58] and GoSam-2.0 [57] is fully automated for one-loop QCD corrections, and extensions thereof to handle loop-induced processes and additional contributions from EFT operators are employed for the process considered here.
Fully automated matching algorithms are available, inspired by MC@NLO [40] and Powheg [90] (referred to as subtractive and multiplicative matching respectively), for the systematic and consistent combination of NLO QCD calculations with both shower variants -7 -

JHEP05(2016)106
in Herwig7 (facilitating two coherent shower algorithms -an angular-ordered parton shower [91] as well as a dipole shower [92], including the simulation of decays with full spin correlations). For the process studied here, we eventually combine our fixed-order NLO result (supplemented with loop-induced and EFT contributions) with the Herwig7 angular-ordered parton shower [91] through the subtractive (i.e. MC@NLO-like) matching algorithm based on [40,58].

Phenomenological studies
For our phenomenological studies of the process pp → e + ν e µ −ν µ + X we choose a centerof-mass energy of 13 TeV. We also study the behaviour of the BSM effects when going from 8 TeV to 13 TeV.
We use the MMHT2014nlo68cl_nf4 [93] parton distribution functions (PDFs), with 4 massless quark flavours 1 in the initial state, and we set α s (M Z = 91.1876 GeV) = 0.12 in accordance with the PDF set we use.
Our default scale choice for the 13 TeV results is a dynamic scale, The mass of the top quark has been set to M t = 174.2 GeV, the Higgs mass to M H = 125.7 GeV. We further use a non-zero top width of Γ t = 1.4 GeV, and a Higgs width of For the electroweak input parameters we follow the SLHA [94] scheme 2 for a set of SM low-scale input parameters to fix the electroweak sector: the electroweak input parameters we choose are M Z = 91.1876 GeV, α = α(M Z ) = 1/128.91 and G F = 1.16637 · 10 −5 GeV −2 , from which M W and sin 2 (θ w ) are subsequently derived. Furthermore we choose Γ Z = 2.4952 GeV and Γ W = 2.085 GeV.
The events are analysed using an in-house analysis for Rivet-2.4 [95] interfaced by Herwig7. The W bosons are directly reconstructed from their leptonic decay products (not via the built-in WFinder function of Rivet). Our cuts on the analysis level are as follows. To mimic W -identification cuts, we employ a cut on the invariant mass of each same flavour lepton-neutrino pair of 60 GeV ≤ m lν l ≤ 100 GeV. We have verified that the exact boundaries of this mass window do not have a major impact on our results. Therefore, these cuts should match within the experimental uncertainties to cuts where a reconstructed transverse mass m lν l T is used for the W -identification. We further require the net transverse momentum of the lepton-neutrino pair to be larger than 10 GeV, and a minimum p T of 25 GeV for each identified lepton and for the missing transverse energy of the event. The identified leptons are required to be in the rapidity range −3 ≤ y l ≤ 3.
For numerical stability we also employ cuts at the generator level, which are of course less restrictive than the cuts employed at the analysis level.
In the following three subsections we will first concentrate on the fixed-order results, and then discuss the impact of a parton shower in section 3.4.
1 For the b quarks circulating in the loops we use m b = 4.2 GeV. We have found that the effect of finite b-quark masses in the loops is below 0.1%.
2 This scheme is frequently used per default by many FeynRules/UFO models.

Gluon-induced contributions and effects of higher-dimensional operators
We start the discussion of the results with the gg-initiated processes, i.e. contributions which are either loop induced in the SM, or require dimension-eight operators to contribute at the tree level. In the following we will distinguish three different contributions, see eq. (2.7). One is the pure Standard-Model loop-induced contribution, where the matrix element is the square of the gg-initiated one-loop amplitude. In the plots this contribution is denoted as gg_SM. The second contribution is the interference term between the SM one-loop amplitude and the EFT tree-level amplitude, which means it is a linear term in the higher-dimensional operators. This contribution is labeled as gg_Interf. Finally the third contribution stems from the square of the EFT tree-level amplitude and is therefore quadratic in the higherdimensional operators. This contribution is labeled as gg_Eff2. This distinction allows us to separate the term linear in the higher-dimensional operators from the quadratic one and study their behaviour independently.
Consequently gg_SM+Interf and gg_Eff2+Interf denote the combination of gg_SM or gg_Eff2 with gg_Interf respectively, while gg_All finally denotes the combination of all contributions to the gg initial state (cf. figure 2 or eq. (2.7)).
For the numerical results we have set unless stated otherwise. Let us first focus on the invariant-mass distribution of the W -boson pair, shown in figure 3a. The invariant mass is defined via the momenta of the decay products: (3. 2) The most striking feature is the fact that the interference term of the SM loop-induced gg → W + W − amplitude with the amplitude induced by the dimension-eight operators is negative. As this cannot be displayed in a logarithmic plot, we display it with the sign switched, denoted by gg_NegInterf. Certainly, the sign of (some of) the dimension eight operators could also be negative, which could entail a positive sign of the interference term. However, for the sake of discussing a concrete example, we will stick to the case of positive signs in fron of the dimension eight opersators in the following.
As expected, the term linear in the dimension-eight operators dominates over the pure EFT contribution (gg_Eff2) which is quadratic in these operators. This is the typical behaviour as the quadratic terms receive an additional suppression of a factor c i /Λ 4 . The contribution of the quadratic term increases with the center-of-mass energy (here √ŝ = m W W ) and will eventually dominate over the linear term. The point where this happens depends on the setup and in particular on the chosen value for the anomalous coupling constants. In our example this happens at about 500 GeV (where the yellow and purple curves cross). As the linear term is negative, this is related to the point where the sum of the two higher-dimensional contributions (gg_Eff2+Interf) becomes positive, which happens a bit earlier at about 400 GeV (where the green and the dashed red curves cross). While the SM contribution drops rapidly as m W W is increasing, the dimension-eight contributions -9 -JHEP05(2016)106  increase and start to dominate at around m W W ∼ 500 GeV. This is the expected behaviour as the contribution from a dimension-eight operator can increase maximally with E 4 /Λ 4 , so an increase ∼ s 2 is the "worst case scenario". We will derive in section 3.1.1 that the actual scaling is milder.
The scale-variation uncertainties (shown as shaded bands in the plots) are relatively large, due to the fact that the results for the gg-initiated subprocess are leading-order predictions. The NLO corrections to this subprocess, calculated in ref. [14], lead to a reduction of the scale uncertainties by about a factor of two, and to a K-factor of about 1.5 for the central scale choice. However, this channel increases the total NNLO corrections only by about 2% [14].
Considering the fact that the gg-initiated channel including the EFT contributions only constitutes a relatively small contribution to the total pp → W W cross section in the region where the EFT approach does not yet conflict with unitarity, we believe that a calculation of the NLO corrections to the dimension eight contribution would not lead to a valuable gain in precision.
The fact that the term linear in the dimension-eight operators is negative (for our choice of positive c i values) leaves us with a phenomenologically interesting constellation. In the region where the linear term is dominant we get a decrease in the cross section compared to the SM prediction, and with increasing invariant mass we observe a (partial) cancellation between the linear and the quadratic term. At the point where linear and quadratic terms are of the same magnitude, we recover exactly the Standard-Model contribution. This means that putting experimental constraints on these types of couplings would be more difficult, and signs of new physics would be masked: in the low energy region the effects of the dimension-eight operators are anyway suppressed by a factor of 1/Λ 4 , and for larger -10 -JHEP05(2016)106 energies we find (partial) cancellation between the linear and the quadratic term. However, we again emphasize that the sign of the dimension-eight operators is not necessarily fixed to be positive. A negative value would lead to a constructive interference between the linear and the quadratic term rather than to a cancellation. The bounds on negative values will therefore be much more restrictive than on positive values.
The region where the quadratic term becomes equally important and eventually dominates over the linear term has to be interpreted with care. For weakly coupled EFT scenarios, the two terms being equally important means that the suppression of the quadratic terms by the additional factor c i /Λ 4 is compensated by a factor of s 2 . In other words s ∼ Λ 2 , which means that we are probing the scale of New Physics and which is the point where the EFT approach becomes invalid, as it is based on the assumption that it is a low-energy effective theory and that the scale of New Physics is much higher than the scale we are probing. It is this assumption that allows us to be confident that operators of lower dimensions are more important compared to higher-dimensional operators. If higher-dimensional operators were not sufficiently suppressed, it would not be justified to neglect dimensionten operators, whose linear terms are actually less suppressed than the quadratic term of a dimension-eight operator. And even worse, in the case of s ∼ Λ 2 all higher-dimensional operators could contribute equally and there is no physically meaningful expansion anymore. On the other hand, in strongly interacting EFT scenarios with c i > 1, the EFT expansion could still be meaningful even if the squared terms dominate over the linear terms, in case we are in a situation where c 2 i (s/Λ 2 ) 2 > c i s/Λ 2 , while still s 2 /Λ 4 < 1. More details, in particular with regards to a possible violation of unitarity, will be discussed in section 3.1.1.
In figure 3b we display the observable ∆R W W = (y 1 − y 2 ) 2 + (φ 1 − φ 2 ) 2 , the separation in rapidity and azimuthal angle between the two W bosons. Here we see that the effects of the higher-dimensional operators lead to an enhancement of the distribution over the whole range. Note that the region below ∆R W W = π is not populated because we only show the fixed-order results in this subsection.
In figures 4a and 4b we show the same observables calculated at a center-of-mass energy of √ s = 8 TeV. We observe that the on-set of the BSM effects in the m W W distribution is around 550 GeV. However, the relative size of the BSM contributions with respect to the SM contribution is much larger at √ s = 13 TeV, as can be seen by comparing figures 3b and 4b. We have verified that the qualitative behaviour in comparing 13 to 8 TeV stays the same if we also use a fixed scale (M W ) for the √ s = 13 TeV case. In figure 5 we compare the two scale choices µ r = µ F = M W (fixed) and µ r = µ F = m W W = (p e + + p µ − + p νe + pν µ ) 2 (dynamic), including the scale uncertainty band, obtained as usual by varying by a factor of two up and down from those central scale choices. The fixed scale M W , being relatively low, leads to a larger value of α s and therefore an increase in the cross section. Since the bands do not overlap, this also means that the scale variations by factors of two up and down are not sufficient to capture the uncertainties correctly.
Another interesting distribution is the relative azimuthal angle between the two charged leptons, ∆φ e + µ − , shown in figure 6a. The contributions from the higher-dimensional operators lead to more highly boosted W bosons, and therefore the associated leptons are more -11 -JHEP05(2016)106  likely to be "back-to-back". A similar behaviour can be seen in the ∆R distribution of the leptons (see figure 6b). We also show the various contributions to the transverse momentum of the positron from the W + decay in figure 7a and the invariant mass of the charged leptons in figure 7b. The BSM effects lead to a clear enhancement in the p e ⊥ spectrum, which becomes quite substantial already for p e ⊥ values as low as ∼ 60 -100 GeV. In the m e + µ − distribution, the effect of the higher-dimensional operators is also clearly visible for energies larger than about 150 GeV to 190 GeV already (taking scale-variation uncertainties into account). However, as this concerns only the gg-initiated contribution, the effect will be washed out once all subprocesses contributing to the pp initial state, plus higher-order corrections, are taken into account, as will be discussed in section 3.3. In order to investigate differences between the three dimension-eight operators, we will now consider them one at a time, always setting the coupling constant of the two others to zero respectively. Let us note that the operators are orthogonal to each other, as can be -13 -  seen from symmetry considerations in eq. (2.8). Therefore, interference terms between the operators will vanish in any case.

JHEP05(2016)106
In figures 8 and 9 the effects of the individual operators are shown for two observables, the angle between the decay planes of the W bosons, cos (Ψ eν,µν ), and their invariant mass, m W W . For these comparisons we have always set one of the c i coefficients to the value 0.3 and the other two to zero respectively. Looking at the decay planes of the W bosons in figure 8, we find that the first two operators, O 1 and O 2 , show the same angular dependence, whereas the angular dependence of the third operator O 3 , which contains the dual fieldstrength tensorW I,µν , is different, which is seen to be more prominent for √ s = 13 TeV. To distinguish between the first and the second operator, the invariant-mass distribution is also a suitable observable, as can be seen in figure 9a. Here the first operator leads to a stronger decrease of the distribution in the region around m W W ∼ 500 GeV. Therefore, the combination of these two observables in principle allows for a distinction between the three operators. However, it should be noted that this can only be a qualitative discussion, as the impact of the dimension-eight operators strongly depends on the size (and on the sign) of the coupling constants c i . The overall size of the BSM effects for our default choice of the anomalous couplings is discussed in section 3.3, where we combine all sub-processes contributing to the pp initial state.

Unitarity bounds
In the context of higher-dimensional operators it is also important to talk about unitarity. As the effects of these operators grow with increasing center-of-mass energy, they will -14 - eventually violate unitarity. For the case of stable W bosons, i.e. for 2 → 2 scattering, a unitarity bound on the total cross section can be derived along the lines of ref. [96]. In more detail, we can start from eq. (48) of ref. [96] (but use total angular momentum J = 0 for the gg-initiated case), where the bound for an inelastic 2 → 2 scattering amplitude T in , summed over the final state helicities λ 3 , λ 4 , is given by To obtain the bound for the total cross section, we include the flux factor 1/(2ŝ), average over initial state colours and helicities, and sum over the colour and helicity configurations contributing to the cross section. This leads to where we have used λ 1 ,λ 2 ∈{+,−} = 4. The sum over the colour states contributing to the amplitude is given by δ ab δ ab = N 2 c − 1 (the trace of the identity matrix in the adjoint representation). Therefore, with N c = 3, we have (3.6) -15 -

JHEP05(2016)106
The derivation of the unitarity bound on the total cross section in eq. (3.6) is based on a 2 → 2 scattering process. To obtain an estimate for the unitarity bound including the decay of the W bosons one could integrate the 2 → 2 process numerically and rescale the result with the branching ratios of the two W bosons decaying into leptons. This procedure, however, does not take into account the effect of the cuts on the leptons, and therefore we refrain from showing a unitarity bound in the plots for the distributions.
Unitarity arguments can also be employed to calculate an upper bound on the absolute value of the anomalous coupling constants. To do so we use the ansatz to require unitarity of the amplitude for a given set of helicities and project the amplitude onto partial waves [97]. Looking at the scattering a + b → c + d with the corresponding helicities λ a , . . . , λ d , the partial wave decomposition reads with λ = λ a − λ b and µ = λ c − λ d , and where θφλ c λ d |T (E)|00λ a λ b denotes the transition matrix element. Its unitarity must hold for each partial wave independently, i.e.
|T J | ≤ 1 . (3.8) Therefore we project the full amplitude onto single partial waves, where the strongest constraints typically come from the lowest order partial waves. In the case where λ = µ = 0, the d J functions reduce to the Legendre polynomials, i.e. d J 00 (θ) = P J (cos θ). Usually it is assumed that the strongest constraints stem from longitudinally polarized W bosons, as in the limit of large momentum k the longitudinal polarization vector behaves like lim (3.9) Projecting onto the 0th partial wave we find c 1 Λ 4 , (3.10) For the third operator, O 3 , the contribution vanishes for longitudinally polarized W bosons. It is also interesting to note that the contributions from the first two operators increase more mildly with energy than naively expected. For dimensional reasons the denominator in eq. (3.10) could be ∼ s 2 , which in turn would mean that the amplitude itself could be ∼ s 2 . However we find only a behaviour which grows like ∼ s.
The fact that the third operator vanishes for longitudinal polarizations, and that we do not find the strongest possible increase with the center-of-mass energy, suggests to also consider the situation where the W bosons are transversely polarized. Projecting these amplitudes onto the 0th partial wave we find (3.11) This means that the strongest constraints for energies above the weak scale come from transversely polarized W bosons and, in order to maintain unitarity, one can roughly assume c i Λ 4 π s 2 . The range in m e + µ − has been limited to 300 GeV here for better visibility of the SM/BSM transition region. We have verified that for very large values of m e + µ − , the yellow and blue curves merge again, as expected.

Impact of heavy-quark loop contributions
We have taken both bottom-and top-quark masses into account for the quark loops mediating the SM gg → W + W − interaction. The effect of the bottom-quark mass is negligible (of the order of the Monte Carlo integration error), while top-quark mass effects have a considerable impact on the partonic cross section in the gg-initiated channel. Figure 10 shows the effect of massive top-quark loops on the invariant-mass distribution of the charged leptons and on the R-separation between the charged leptons. We observe that top-mass effects decrease the m e + µ − distribution by more than 30% below values of m e + µ − ∼ 250-300 GeV. The SM contribution with M t set to zero (gg_SM_massless) is of the same size as the SM+BSM result with masses taken into account (gg_All) at m e + µ − ∼ 200 GeV, which means that neglecting the top-quark mass in the SM calculation could potentially "fake" BSM effects.
It should be noted here that in the case of massless top quarks also the Yukawa coupling between the Higgs boson and the top quark vanishes.
This result is contrasted to a calculation where the top-quark loops have been dropped altogether, 3 shown in figure 11. This has a considerable impact on the m e + µ − distribution beyond about 150 GeV, however, the effect is much less pronounced than in the case where top-quark loops are taken into account but the top-quark mass is neglected (cf. figure 10).
3 It should be noted here that omitting the top quarks also eliminates almost all contributions involving bottom quarks. Only the diagrams where a Higgs boson couples to a b-quark pair remain, which are numerically negligible. Therefore, omitting the top quark loops basically means excluding the third quark generation. Even though the mass effects are below the 10% level once the full pp initial state including the NLO corrections is taken into account (see figure 12), this study demonstrates that massive top-quark loops should be taken into account to describe measurements of highly boosted W bosons.

Combination of gluon-and quark-initiated channels
In this section we compare the gg-initiated contribution to the full process pp (→ W + W − ) → e + ν e µ −ν µ at NLO QCD, considering first the results at the fixed-order level. Shower effects will be discussed in section 3.4. The important questions here are how visible the effects of the anomalous couplings are if the quark-initiated Standard-Model contributions are added, and to what extent the effects of the higher-dimensional operators are hidden in the theoretical uncertainties.
In figure 13 we show the invariant-mass distribution of the W -boson pair including all SM as well as EFT contributions. The effects of scale variations are plotted as well, where the scale uncertainty bands have been obtained by varying by a factor of two up and down from the dynamic scale choice µ R = µ F = m W W . We show the SM NLO contribution with and without the EFT contributions, and in comparison to that the effects of the higher-dimensional operators in the gg-initiated contributions alone. This allows to directly assess the impact of the anomalous couplings. The loop-induced, gg-initiated Standard-Model contribution leads to an O(10%) increase over the quark-or quark-gluon-initiated NLO result. We therefore observe that in the full result, combining all channels, a visible deviation from the SM prediction begins to show at larger m W W values, of about 700 GeV, while in the gg-initiated contribution, shown in figure 3a, the deviation already starts to be visible at about 500 GeV to 600 GeV (taking scale-variation uncertainties into account).

JHEP05(2016)106
While for values of m W W around 700 GeV the size of the BSM effects is comparable to the size of the scale uncertainties, shown in figure 13, for values of m W W of about 800 -900 GeV, the deviations from the Standard Model due to the higher-dimensional operators start to become clearly visible. On the other hand, the region beyond 1 TeV already probes energies where the EFT approach starts to become invalid.

Parton-shower effects
A realistic description of observables measured in hadronic collisions includes effects from a parton shower. We therefore combine our fixed-order NLO results (supplemented with loop-induced and EFT contributions) with the Herwig7 angular-ordered parton shower [91] through the subtractive (i.e. MC@NLO-like) matching algorithm based on [40,58]. 4 Uncertainties in the shower are mainly quantified by varying the hard shower scale µ Q which provides a reliable estimate of missing logarithmic orders as well as the impact of largeangle, hard and thus unreliably modelled emissions. It also serves as a check to verify the improvements expected from NLO plus parton-shower matching, as studied e.g. in [42,55].
While, as expected, the invariant mass of the reconstructed W + W − system is not affected by the parton shower, a number of other infrared sensitive observables receive large contributions, both through the NLO real radiation and further subsequent partonshower emissions. Typical infrared-sensitive distributions in this case are the R-separation of the two W bosons (where in the zero-jet limit ∆R is composed solely of a difference in rapidity, while ∆φ = π), as well as the transverse momentum of the reconstructed W + W − system, shown in figure 14. Both observables show the expected behaviour with respect to additional radiation; in the region ∆R < π, both the NLO real emission as well as shower emissions off the gg-induced channel contribute. The first contribution includes a small shower uncertainty, as this kinematic range has been improved by the NLO matching. Once the BSM contribution to the gg-channel becomes dominant, pure shower emissions off this sub-process become more important and hence yield a larger uncertainty. Ultimately, NLO QCD corrections, or at least a leading-order multi-jet merging are desirable in this case. Similar features are present in the transverse momentum of the reconstructed W + W − pair. Azimuthal and R-separations of the charged leptons are sensitive to the BSM contribution and very stable against QCD activity, as shown in figure 15.
We finally discuss a few observables which are relevant to the experimental reconstruction of the W + W − final state, particularly lepton-jet separations and the distribution of the transverse momentum of the positron, displayed in figure 16. While shower uncertainties at the level of 10% are observed, the lepton-jet separation is rather stable against QCD activity, and BSM contributions only affect the normalization in the small-∆R region; the experimentally required lepton-jet isolation is thus not introducing any bias. Larger impact is observed on the transverse momenta of the charged leptons (e.g. as shown for the positron in figure 16), which, however, turn out to be rather stable with respect to parton-shower scale variations.

Conclusions and outlook
The production of electroweak gauge-boson pairs is amongst the most important signatures at the LHC. These final states are important Higgs-boson decay channels, and they allow us to study the electroweak sector, with the aim to reveal the mechanism of electroweak symmetry breaking. We have studied the production of a pair of W bosons at NLO QCD, in the light of additional anomalous couplings. In particular, we have also included the gg-initiated (loop-induced) process gg (→ W + W − ) → e + ν e µ −ν µ , which is formally a contribution to the NNLO result, but is enhanced due to the large gluon luminosity at the LHC. In addition to the Standard-Model gg-initiated contribution, we have included gg-initiated contributions stemming from dimension-eight operators which induce a treelevel coupling between gluons and electroweak gauge bosons. This possibility has not been discussed in the literature so far. Their presence leads to an interference between the Standard-Model gg-induced one-loop amplitude and a tree-level amplitude mediated by dimension-eight operators.
We have discussed their effects on a variety of important observables. We have found that the presence of dimension-eight operators can lead to substantial effects in the highenergy tail of the distributions, which can be used by the LHC experiments to constrain the parameter space for the associated effective couplings. Heavy resonances decaying to a W -boson pair would lead to a more distinct signature, unless the energy to produce them directly is not sufficient, in which case a growing slope in the tail of the m W W distribution may also be observed.
Furthermore, we have investigated the importance of heavy (SM) quarks in the loopinduced process, leading to corrections of up to 10%, depending on cuts and center-ofmass energy.
Finally, by combining our fixed-order results with the Herwig7 angular-ordered parton shower, we have studied the effects of a parton shower, including variations of the hard shower scale, on the leptonic observables and on observables related to the reconstructed W + W − system.