Top-quark pair production in association with a $W^\pm$ gauge boson in the POWHEG-BOX

We present a new Monte Carlo event generator for the production of a top-quark pair in association with a $W^\pm$ boson at hadron colliders in the POWHEG-BOX framework. We consider the next-to-leading-order QCD corrections to the $pp\to t\bar{t} W^\pm$ cross section, corresponding to the $\mathcal{O}(\alpha_s^3\alpha)$ and $\mathcal{O}(\alpha_s\alpha^3)$ terms in the perturbative expansion of the parton-level cross section, and model the decays of $W$ and top quarks at leading order retaining spin correlations. The fixed-order QCD calculation is further interfaced with the Pythia8 parton-shower event generator via the POWHEG method as implemented in the POWHEG-BOX. The corresponding code is now part of the public repository of the POWHEG-BOX. We perform a comparison of different event generators for both the case of inclusive production and the case of the two same-sign leptons signature at the Large Hadron Collider operating at a center-of-mass energy of $13$ TeV. We investigate theoretical uncertainties in the modelling of the fiducial volume stemming from missing higher-order corrections, the different parton shower matching schemes, and the modelling of decays. We find that the subleading contribution at $\mathcal{O}(\alpha_s\alpha^3)$ is particularly sensitive to differences in the matching scheme and higher-order parton shower effects. We observe that in particular jet observables can differ quite visibly although these differences play only a subordinate role in the description of physical observables once all contributions are combined.


I. INTRODUCTION
The production of top-quark pairs in association with electroweak gauge bosons (W, Z, γ) can be measured at the Large Hadron Collider (LHC) and future hadron colliders (HL-LHC, FCC-hh, CppC) in a multitude of decay channels and provides new avenues to test the consistency of the Standard Model (SM) in processes that have been beyond the energy reach of existing colliders or statistically limited until recently. At the same time, these processes represent some of the most important backgrounds for Higgs-boson precision measurements and searches of new physics beyond the Standard Model (BSM). In this context, the hadronic production of W ± bosons in association with top-quark pairs is particularly interesting both from a phenomenological and a theoretical point of view.
On top of its intrinsic interest as a SM process, ttW ± provides very distinctively polarized top quarks that can be used to unveil the imprint of new physics interactions. Indeed, in contrast to top-quark pair production, the top quarks originating from the ttW ± production process are highly polarized and give rise to a large tt charge asymmetry [1,2]. A measurement of the tt charge asymmetry in ttW ± production can then be sensitive to the chiral nature of new physics contributing to the process and can become a unique indirect probe of BSM physics.
The ttW ± processes also represent a very important background to the production of a Higgs boson in association with top quarks in the multi-lepton decay channels [3][4][5][6][7], where it limits the accuracy of the direct measurement of the top-quark Yukawa coupling. In addition, ttW ± is the dominant background in searches for the SM production of four top quarks [8,9]. In general, ttW ± is a background to any search of new physics in signatures with same-sign leptons, missing energy, and b jets, common in many BSM models.
Due to its phenomenological relevance, the process has been studied extensively on the theory side, starting from the first calculation of next-to-leading order (NLO) QCD corrections to the production and decay process in Ref. [10]. Further studies at fixed order include the calculation of the leading NLO electroweak (EW) corrections [11] and the assessment of the impact of formally subleading mixed QCD and EW corrections [12]. Beyond fixed order the pure NLO QCD ttW ± calculation has been matched to parton showers using the POWHEG method [13,14] as implemented in the PowHel framework [15] as well as using the MC@NLO method [16,17] in the MG5 aMC@NLO framework [1] 1 Separately, the resummation of soft gluon emission effects have been studied at the next-to-next-to-leading logarithmic (NNLL) accuracy [25][26][27][28][29].
On the experimental side, the associate ttW ± production has been measured by the LHC experiments both as inclusive cross section [30,31] and as a background in searches for ttH and tttt signals [6][7][8][9] in multi-lepton decay channels. Some of these measurements have resulted in larger values with respect to SM predictions. Because of this, the modelling of the ttW ± processes has come under more thorough scrutiny, with the recent inclusion of off-shell and non-resonant effects at fixed-order NLO QCD [32,33] and by studying the effect of spin-correlations and formally subleading EW corrections in the fiducial volume of specific ttW ± signatures [34]. Also, estimates of the impact of higher-order QCD corrections beyond NLO corrections, estimated via multi-jet merging, on inclusive ttW ± samples, has been presented in Ref. [35]. Furthermore, the ATLAS collaboration recently performed a dedicated comparison of the implementation of ttW ± production in existing NLO partonshower Monte Carlo event generators, including both QCD and EW effects and allowing for multi-jet merging [36].
In this paper we continue the investigation of modelling uncertainties of the ttW ± process, by presenting a study based on a new implementation of the pp → ttW ± production in the POWHEG-BOX framework [37]. We consider next-to-leading-order QCD corrections to the pp → ttW ± cross section, corresponding to the O(α 3 s α) and O(α s α 3 ) terms in the perturbative expansion of the parton-level cross section, and model the decays of W and top quarks at leading order (LO) retaining spin correlations.
As part of our study, we perform a detailed comparison between different NLO partonshower Monte Carlo event generators at both the inclusive and the fiducial level in order to address modelling uncertainties. We compare results obtained with our POWHEG-BOX implementation (interfaced to the PYTHIA8 [38,39] parton shower), with Sherpa (using its parton shower based on Catani-Seymour dipoles [40]), and with MG5 aMC@NLO (also interfaced to PYTHIA8) and present a first study of the consistency between these different frameworks.
Our comparison provides a solid basis on which to develop a more robust estimate of the 1 A comparison of ttW ± differential distributions obtained from MG5 aMC@NLO [18], PowHel [15] and Sherpa [19,20]+OpenLoops [21][22][23] including O(α 3 s α) NLO QCD corrections has been presented in Ref. [24] as a validation of the corresponding Monte Carlo tools. Ref. [24] also provided LHC ttW ± cross sections for √ s = 13 and 14 TeV including O(α 3 s α) and O(α 2 s α 2 ) NLO corrections as obtained from both MG5 aMC@NLO and Sherpa+OpenLoops. residual theoretical uncertainty, and suggests which aspects of the theoretical prediction for hadronic ttW ± production still need improvement.
The paper is organized as follows. In section II we review the POWHEG method to the extent of allowing us to establish our notation. In section III we provide details of our implementation of the ttW ± processes in the POWHEG-BOX framework. In section IV we present theoretical predictions for both inclusive ttW ± production and for a two same-sign leptons signature, comparing results from different Monte Carlo event generators. Finally, we give our summary and outlook in section V.

II. REVIEW OF THE POWHEG FRAMEWORK
Matching parton-shower Monte Carlo event generators with fixed-order perturbative calculations to achieve NLO accuracy for inclusive observables, has proven fundamental to describe LHC data. Two major strategies for this matching are commonly employed, the MC@NLO [16,17] and POWHEG [13,14] methods. In this article we employ the latter using the POWHEG-BOX framework [37] to study event simulation associated to ttW ± production at hadron colliders.
Within the POWHEG method the matching of NLO matrix elements to parton showers is achieved by generating the hardest emission first with NLO QCD accuracy, while subsequent emissions are modelled by the parton shower. Starting with the NLO fixed-order cross section where dΦ n denotes the n-particle Lorentz-invariant phase space measure and B(Φ n ), V (Φ n ), and R(Φ n+1 ) are the differential Born, virtual, and real cross sections, one introduces a jet function F (Φ n+1 ) to split the real radiation contribution as: and uses it to define the soft R s (Φ n+1 ) and hard R h (Φ n+1 ) real contributions according to: The jet function F (Φ n+1 ) is a real function which takes values between 0 and 1, and should approach smoothly 1 in the infrared (soft and collinear) limits of the (n + 1)-particle phase space Φ n+1 . The precise functional form of F (Φ n+1 ) is in principle arbitrary, but a judicious choice is in certain cases necessary to avoid large matching corrections, which are formally subleading. Below we discuss standard choices made in the POWHEG-BOX, and in section IV we will study their impact on the process at hand.
To generate the first (hardest) emission while keeping the NLO accuracy for inclusive observables, a one-step parton shower is introduced according to where the real emission phase space is factorized as dΦ n+1 = dΦ n dΦ r in terms of the underlying Born phase space with n final-state particles (dΦ n ) and the phase space of the radiated particle (dΦ r ), the B(Φ n ) function is defined by and we have introduced the modified Sudakov form factor ∆(Φ n , p T ), which is defined as: In Eq. (4) the parton shower infrared cutoff scale is denoted with p min T and p T = p T (Φ r ) is the transverse momentum of the emitted particle.
The modified Sudakov form factor ensures that no double counting of real radiation in R s is produced, while radiation from R h is treated as an independent contribution. This highlights the importance of the choice of F (Φ n+1 ) for the POWHEG method. As already mentioned, after the first (hardest) emission, subsequent splittings can be generated by a standard parton shower without affecting the NLO accuracy for inclusive observables.
In the POWHEG-BOX the function F (Φ n+1 ) is written as the product of two functions [37,41] according to where by default F damp (Φ n+1 ) and F bornzero (Φ n+1 ) are set identical to 1, such that F (Φ n+1 ) = 1 and consequently R s = R and R h = 0. When enabling a non-trivial F damp in the POWHEG-BOX, the corresponding function takes the form where P ij (Φ r ) ⊗ B(Φ n ) is an approximation to R(Φ n+1 ) based on the factorization properties of the real amplitudes in the soft and collinear limits. The dimensionless parameter h bornzero controls how much phase space outside of the singular limits is associated to R s 2 .
In particular, it ensures that if B(Φ n ) vanishes in certain regions of phase space, no large contributions from the first term on the right-hand side of Eq. (4) will be produced. While historically this was the main reason to introduce F bornzero , this function also helps to identify and distinguish different enhancement mechanisms of the real matrix elements, which should not be interpreted as due to QCD splittings. As we will see, this plays a crucial role for the O(α s α 3 ) contributions in our present study.
As a final remark, we would like to highlight that the role played by the POWHEG damping functions is different from the role played by the initial shower scale µ Q in the MC@NLO method. In the POWHEG method the damping functions are used to control the degree of resummation for non-singular contributions, while the initial shower scale is always assigned by the POWHEG framework to be the transverse momentum of the parton splitting independent of whether that splitting is associated to R s , the soft or to The default value of h bornzero is 5 in the POWHEG-BOX.
have a substantial impact on the generated event sample but only a mild impact on the subsequent parton shower evolution. In contrary, the initial shower scale µ Q of the MC@NLO method directly controls the available phase space for subsequent parton shower emissions and therefore can have a strong impact on the shower evolution. The impact of the damping function and the initial shower scale is formally of higher order and does not spoil the NLO accuracy of the predictions. Nonetheless, these higher-order corrections can become sizable.

III. DETAILS OF THE CALCULATION
In this section we present our implementation of the pp → ttW ± process in the POWHEG-BOX including NLO corrections. We start by discussing the perturbative orders in α n s α k that we are considering.  At tree level the ttW ± final state can be generated in hadron collisions via qq → ttW ± subprocesses, as illustrated in Fig. 1. As illustrated in Fig. 3 The possible coupling combinations contributing to the LO and NLO ttW ± cross section.
The links indicate how a given NLO order originates from a corresponding LO order via QCD or EW corrections. The terms crossed out vanish by color structure. In this study the orders corresponding to the shaded bubbles are neglected.
In the case of ttW ± , the identification in terms of QCD and EW corrections can be restored to a very good approximation by considering the hierarchy of the different leading and subleading orders of the NLO cross section. As shown in Refs. [12,43] The enhancement of the O(α s α 3 ) terms originates from NLO QCD real corrections from the qg channel that opens at NLO QCD, and more specifically from the kind represented by 3 Incidentally, notice that the O(α 2 s α 2 ) and O(α 4 ) NLO corrections include photon-initiated contributions of the form qγ → ttW ± q which however amount to a small fraction of these already subleading corrections [11]. the right-hand side diagram in Fig. 4. In these particular contributions, the parton density enhancement of qg versus qq is largely amplified by the combined effect of several factors, from the t-channel kinematic, to the rescattering of W and Z longitudinal components, and the presence of a large top-quark Yukawa coupling.  Fig. 3). As such, the fixed-order NLO QCD calculation can be consistently interfaced with a QCD parton-shower event generator.
In the following we will then focus solely on the contributions at the perturbative orders O(α 2 s α) and O(α 3 s α), which we will denote from now on as 'ttW ± QCD', and the orders O(α 3 ) and O(α s α 3 ), which we will denote as 'ttW ± EW'. Results from our implementation of pp → ttW ± in the POWHEG-BOX have the same level of theoretical accuracy as the ones presented in Ref. [34], and can be directly compared to the ones that can be obtained via analogous NLO QCD parton-shower Monte Carlo event generators such as Sherpa or MG5 aMC@NLO + PYTHIA8. Indeed, we will show a corresponding comparison among these tools in section IV.
Next, we will discuss in section III A the implementation of pp → ttW ± including the aforementioned orders of NLO QCD corrections in the POWHEG-BOX framework. Furthermore, in section III B we will give some further details on the modelling of a fully realistic final state by including decays of unstable particles while keeping spin correlations at LO accuracy.
A. NLO corrections to the production of ttW ± The implementation of a new process in the POWHEG-BOX requires to provide processspecific ingredients such as the LO and NLO virtual and real matrix elements as well as the parametrization of the Born-level phase space, while all process-independent parts, such as the subtraction of infrared singularities, are automated. In our implementation of ttW ± , all tree-level matrix elements including spin-and color-correlated Born matrix elements are taken from the MadGraph 4 [44,45] interface that is provided within the POWHEG-BOX. The finite remainders of the virtual loop corrections interfered with the Born matrix elements are computed by the fairly new one-loop provider for QCD and EW corrections NLOX [46,47] that uses OneLOop [48] for the evaluation of scalar Feynman integrals. The parametrization of the ttW ± phase space has been directly modelled on the POWHEG-BOX implementation of ttH [49].
We performed several checks to validate our implementation. For example, virtual amplitudes have been successfully compared at a few phase-space points against Recola [50,51] and MadGraph5 [52]. At the same time, total inclusive cross sections at fixed order have been cross checked with MG5 aMC@NLO [18], while a full validation at the differential level has been performed by comparing to an independent calculation obtained from Sherpa [19,20] in conjunction with either the version of the Blackhat library of Ref. [53] or with the OpenLoops program [21][22][23].
On the other hand, as explained in section II, in matching the fixed-order calculation with parton-shower event generator within the POWHEG-BOX framework, we have implemented the following choices for the jet function F (Φ n+1 ) of Eq. (7). First, similar to the case of bbW ± [54], the leading-order matrix element at O(α 2 s α) has vanishing Born-level configurations. Therefore, the usage of F bornzero (Φ n+1 ) (see Eq. (9)) is mandatory for ttW ± QCD production. This is also applied to ttW ± EW production, where it plays a critical role given the presence of the strongly enhanced real matrix element discussed before (due to the diagrams of the right of Fig. 4). Finally we have enabled the damping function F damp (Φ n+1 ) of Eq. (8) to further suppress hard radiation using a dynamic value of h damp . More specifically, our default choices for the damping parameters will be: where is evaluated on the underlying Born kinematics.
In order to disentangle the impact of the jet function F (Φ n+1 ) from possible parton-shower corrections, we investigate the differential distributions of POWHEG-BOX events without taking into account the parton-shower evolution. We show a few representative observables in Fig. 5 at the level of Les Houches events (LHE) [55,56]. On the left we show transverse momentum distributions for the leading jet, the top-quark pair, and the W ± boson for inclusive ttW ± QCD production, while on the right we show the transverse momentum and pseudorapidity distributions of the leading jet, as well as the transverse momentum distribution of the W ± boson for inclusive ttW ± EW production. In all cases, we compare POWHEG-BOX predictions with and without the effect of the jet function F (Φ n+1 ) (labeled as 'LHE -damping' and 'LHE -no damping' respectively) to the corresponding fixed-order differential distributions.
For inclusive ttW ± QCD production we observe large shape differences between the fixedorder NLO prediction and the POWHEG-BOX results without the jet function. In hadronic observables such as the transverse momentum of the hardest jet the differences reach up to a factor of 3 at a transverse momentum of around 600 GeV. A similar behavior can be seen in the transverse momentum of the top-quark pair, where the deviations are up to +70% at the end of the plotted range. Even in non-hadronic observables such as the transverse momentum of the W ± boson moderate differences at the level of 10% are visible. However once the damping mechanism is taken into account the predictions recover the tails of the NLO fixed-order distributions, which are reliably described by fixed-order matrix elements.
Though considerable improvement is achieved for the transverse momentum distribution of the hardest jet with the inclusion of the jet function, still a small difference of about 10% is observed for large p T (j 1 ). On the other hand, shape differences for low p T values such as in the case of the transverse momentum of the leading jet are attributed to the resummation of  (7)) on differential distributions. The LHEdamping curves correspond to our default choice of parameters according to Eq. (10) soft and collinear QCD splittings via the Sudakov form factor and are thus expected to be different from the fixed-order result. In the case of inclusive ttW ± EW production the impact of the jet function is more dramatic. Here, the NLO QCD corrections are dominated by real radiative contributions that are enhanced by the t-channel EW scattering of tW → tW shown in Fig. 2 and Fig. 4. Therefore, the resummation of the full real matrix element yields unphysical results as non-factorizing contributions are resummed and thus the jet function F (Φ n+1 ) has to be used to restrict the resummation to singular QCD splittings that indeed factorize. As a consequence of choosing a trivial jet function F (Φ n+1 ) = 1, the shapes of the transverse momentum and the rapidity distributions of the leading jet are described poorly over the entire phase space. Also the transverse momentum of the W ± boson shows large discrepancies in the tail of the distribution. Nonetheless, using our default jet function F (Φ n+1 ) parameters (see Eq. (10)) we observe an excellent agreement between the NLO fixed-order computation and the POWHEG-BOX results. Furthermore, the impressive agreement of the transverse momentum distribution of the leading jet down to low values of p T (j 1 ) just highlights how radiation attributed to the hard contribution R h (see Eq. (4)) dominates the production of that leading jet, while enhancements from resummation effects are nearly negligible.
B. Decay modelling of the ttW ± system In order to be able to study a broader range of exclusive observables, which depend strongly on spin correlations from production to decaying particles, we include the decay of the ttW ± system in our implementation. Typically, parton-shower event generators decay unstable particles during the shower evolution. This approach is the simplest but it cannot preserve spin correlations in the decays, as each particle is decayed independently. In the following we briefly discuss our implementation of the decays in the POWHEG-BOX that preserves spin correlations at least to LO accuracy, which can have a sizable impact on the top-quark decay products since the emission of the W ± boson in the initial-state polarizes the top quarks [1,34].
Our approach follows closely the method of Ref. [57]. This method has been adopted in the POWHEG-BOX already for several processes (see for example Refs. [49,54,58,59]). The basic idea can be summarized as follows. Starting from an on-shell phase space point for the ttW ± momenta one performs a reshuffling of the momenta to allow for off-shell virtualities of the unstable particles. Afterwards, momenta of decay products are generated uniformly in the decay phase space and finally these momenta are unweighted against the fully decayed matrix element by constructing a suitable upper bounding function.
All previous implementations in the POWHEG-BOX have in common that only the decay of a top-quark pair has been taken into account. In order to allow these reshuffles for more than two unstable particles, we introduce a new momentum mapping in a process-independent and Lorentz-invariant way.
The POWHEG-BOX, after the single-step parton shower process (see section II), generates on-shell momentum configurations Φ OS n that either refer to a ttW ± or a ttW ± j final state. We start by generating independent virtualities v 2 for each top quark and W ± boson according where m f and Γ f are the particle mass and decay width. We constrain the generated virtualities to the window |v − m f | < 5 Γ f . These virtualities will be imprinted onto the momenta of the on-shell phase space point Φ OS n by the repeated application of the mapping presented in Appendix A, where we choose to preserve the momentum of the light jet in the case of a real radiation event. This amounts to choosing Q = p t + pt + p W as the total available momentum in the mapping of Appendix A and excluding the jet from the Lorentz boost. In the following we will call the off-shell phase space configuration simply Φ n .
Afterwards, momenta of the decay products are uniformly generated as a sequence of 1 → 2 decays. This allows also to include off-shell W bosons in the top-quark decays. In the last step we apply the hit-and-miss technique (see e.g. Ref. [58]) on the so obtained final state momenta using an upper bounding function , constructed according to Ref. [57], such that where M undec (Φ OS n ) is the LO matrix element for the undecayed process, while M dec (. . . ) is a so-called decay-chain matrix element, which corresponds to a LO matrix element for the fully decayed process where only diagrams with the resonance structure of interest are kept.
The result of this procedure is such that all spin correlations between unstable particles and decay particles are kept with LO precision.
For completeness, we briefly discuss the extension to the ttW ± process of the method to choose particular decay signatures, which has been previously used in POWHEG-BOX implementations. We present details for the ttW + process, as the ttW − process is treated in an analogous way. All top quarks decay into a W boson and a b-quark, and we now have to take into account the decay of three W bosons. Based on the branching ratios Br(W → i ν i ) and Br(W → q iq i ) we can construct a density matrix ρ for the decay probabilities where X i (X i ) represent the possible (charge conjugated) final states, namely The total branching ratio is then given by which can be less than 1 if, for example, only particular decay channels are selected.
To choose a particular decay channel we simply perform a hit-and-miss procedure on the components of the probability density matrix ρ. To allow for hadronic decays we take into account the CKM mixing of the first two generations, which is parametrized by For example, when the decay W + → ud has been chosen we simply generate a random number r ∈ [0, 1] and if r ≤ |V us | 2 = sin 2 θ c we choose the decay W + → us, otherwise we keep the decay W + → ud.
At last we would like to mention that with the procedure outlined above the symmetry factors for identical final state particles are obtained in the correct way. For example, for a calculation employing full off-shell matrix elements it is clear that due to symmetry factors. Since in our methodology we generate decays subsequently we have to ensure that the probability to generate the e + µ − τ + final state is twice as large as the one for the e + µ − e + final state. This is trivially ensured in our procedure, since and where we have used Br(W → i ν i ) = 1/9.
We performed several cross checks on the modelling of the decays. For instance, we checked that the correct branching ratios are obtained from inclusive event samples that take into account all possible decay modes. Also, we compared at the differential level unshowered leading-order events for a particular decay mode with events obtained through the same procedure by MG5 aMC@NLO in conjunction with MadSpin [60] and found perfect agreement.

IV. PHENOMENOLOGICAL RESULTS
In this section we present and discuss numerical results for pp → ttW ± obtained with the new POWHEG-BOX implementation described in this paper and compare them with analogous results obtained from MG5 aMC@NLO and Sherpa. After reviewing the general settings for the input parameters of our study in section IV A, we will consider the case of inclusive ttW ± production in sec. IV B and the case of a specific signature with two same-sign leptons and jets (usually referred to as 2 SS) in section IV C. In the following we will denote by ttW ± the sum of both ttW + and ttW − production.

A. Computational setup
In our study we consider ttW ± production at the LHC with a center-of-mass energy of √ s = 13 TeV. All results presented in this section have been obtained using the NNPDF3.0 [61] (NNPDF30-nlo-as-0118) parton distribution functions as provided by LHAPDF [62]. We have not performed a detailed study of the PDF uncertainty associated with this production mode since it does not directly affect either the comparison between fixed-order and parton-shower results or the comparison between different NLO parton-shower Monte Carlo event generators. Of course such uncertainty should be included in future more comprehensive assessments of the overall theoretical uncertainty on this production mode.
The necessary Standard Model parameters have been chosen to be: in terms of which the electromagnetic coupling is defined as: The central values of the renormalization and factorization scales are set to: with The theoretical uncertainties associated with this choice of scales are estimated via the 7-point envelope that corresponds to µ R and µ F assuming the following sets of values: We notice that in quoting the dependence of the results presented in this section from scale variation, the previously defined scale variation has been applied only to the hard matrix elements, while scale choices in the parton shower have not been altered.
In order to assess the nature and size of theoretical uncertainties present in the modelling of the ttW ± process, we will compare results from three different NLO parton-shower event generators, POWHEG-BOX, MG5 aMC@NLO, and Sherpa, using the setups described next. In the following figures we will label the results obtained using different tools accordingly.
• POWHEG-BOX: We generate results using the POWHEG-BOX implementation presented in this paper by performing the parton-shower matching to PYTHIA8 [39] (v. 8.303) using the POWHEG method as described in sections II and III A, and modelling the decay of the ttW ± final state as discussed in section III B. We employ by default the damping parameters shown in Eq. (10), that is: where the dynamic damping parameter h damp is evaluated on the underlying born kinematics. The impact of different choices of damping parameters is estimated by the 5-point envelope that corresponds to choosing the values: For PYTHIA8 we use the A14 shower tune and turn off all its matrix element corrections (MEC) to decay processes.
For PYTHIA8 we use the A14 shower tune with the standard MG5 aMC@NLO parameter settings, which also do not include MEC to decay processes.
• Sherpa: Finally, we use the Sherpa (v2.2.10) parton-shower event generator to produce a third set of independent results. For ttW ± QCD production we use the MC@NLO matching procedure and Sherpa's parton shower, which is based on Catani-Seymour dipole factorization [40]. The corresponding one-loop matrix elements for these studies are taken from the OpenLoops [21][22][23] program. To produce predictions with Sherpa's public version for ttW ± EW production we employed a truncated shower merging procedure (MEPS/CKKW) setup [63,64], including LO matrix elements for the processes We do not include in any of the previous setups non-perturbative corrections like those from hadronization or multiple parton interactions. Notice that, even though the hard NLO computation is performed in a massless 5 flavor scheme the used parton showers in this study adopt by default a non-zero bottom mass.
Combining the effect of renormalization and factorization scale variation with the variations induced by different choices of matching schemes and parton-shower codes enables us to address the numerical impact of higher-order corrections that are inherent to each approach. We want to emphasize however that the variation of damping parameters and the initial shower µ Q have different effects as explained at the end of section II. Nonetheless, we will refer to the sensitivity of the predictions on these parameters as an uncertainty.
The showered events are finally passed via the HepMC [66,67] interface to an analysis routine written in the Rivet framework [68,69].
In ancillary material that we make available with this document we provide the necessary files to reproduce our results. They include run cards for POWHEG-BOX and Sherpa as well as our Rivet analyses.

B. Inclusive NLO+PS observables
In this section we study the on-shell inclusive production of ttW ± and the impact of the parton shower evolution on inclusive observables. To this end we focus on differential distributions that can also be computed at fixed-order and compare the NLO fixed-order results to the corresponding results obtained by parton-shower matching in different frameworks.
For this analysis we keep the top quarks and W ± bosons stable. No cuts are applied on the top quarks and W ± bosons. Jets are formed using the anti-k T jet algorithm [70] with a resolution parameter of R = 0.4 as implemented in FastJet [71,72]. We require a minimal transverse momentum of p T > 25 GeV for all jets. Note that we do not distinguish between light and b-flavored jets in obtaining the results presented in this section.
In the following we will discuss first the ttW ± QCD and ttW ± EW contributions to the ttW ± inclusive cross section separately, as they show different overall features. Afterwards, we will study the impact of the EW contribution on the QCD+EW combined prediction for a few representative observables.
ttW ± QCD contribution We start our discussion by focusing on the dominant ttW ± QCD contribution. Fig. 6 shows the transverse momentum of the top-quark pair (p T (tt)) as well as of the W ± boson   the fixed-order result. Only results obtained from MG5 aMC@NLO show a small shape difference with respect to the fixed-order NLO prediction. However, within the estimated theoretical uncertainties of 10 − 15%, predictions from all generators considered agree well with each other. Again, the matching uncertainties are small compared to the scale uncertainties, as it is expected for inclusive NLO observables, and are below 5%.
Next we turn to less inclusive observables such as the transverse momentum (p T (j 1 )) and the pseudorapidity (η(j 1 )) of the hardest jet, as shown in Fig. 7. These observables are only accurate to LO since they do not receive one-loop corrections and can thus be more affected by the parton shower evolution. For example, in the case of the transverse momentum of the hardest jet, modifications to the spectrum are expected for small transverse momenta due to the Sudakov resummation in the parton shower, while the high-energy tail should be well described by using fixed-order matrix elements. Indeed all predictions including partonshower effects differ from the fixed-order curve by up to 21 − 30% for small transverse momenta, where POWHEG-BOX predictions differ the most. On the other hand, in the high-  and Sherpa give a slightly softer spectrum than the fixed-order curve by roughly 10% at the end of the plotted spectrum. In addition to the shape differences we also notice a severe reduction of the theoretical uncertainties in the beginning of the spectrum where the partonshower is expected to dominate. While the fixed order has nearly constant uncertainties of the order of 30 − 35% over the whole range, the parton-shower based predictions show only a variation of below 10% at the beginning of the spectrum while the uncertainty grows up to 35% at high p T once the real matrix elements dominate the spectrum again. Over the whole plotted range the scale uncertainties dominate over the matching related ones. However, the latter are larger for MG5 aMC@NLO as compared to the POWHEG-BOX and amount to roughly 15%.
The plot on the right-hand side of Fig. 7 illustrates the pseudorapidity of the hardest jet.
Here, we observe that including parton-shower effects generates overall positive corrections of the order of 13 − 20%. These corrections can be simply attributed to the fact that after the parton shower evolution the number of events with at least one jet is higher than in the corresponding fixed-order NLO computation. To be precise, MG5 aMC@NLO predicts 11%, Sherpa 13%, and the POWHEG-BOX 15% more events with at least one hard jet. The small differences between MG5 aMC@NLO and Sherpa can be attributed to different shower dynamics, while the difference between the POWHEG-BOX and MG5 aMC@NLO are related to the parton shower matching scheme (as shower differences are minor between these two). Also notable is the reduction of the theoretical uncertainties by a factor of two, from ±30% at fixed-order to ±15% in the central rapidity bins once the predictions are matched to a parton shower. As we can see from the bottom panel the matching uncertainties are negligible over the whole spectrum.
ttW ± EW contribution Moving on to ttW ± EW contributions to the inclusive ttW ± signature, we illustrate in the left-hand side plot of Fig. 8 the transverse momentum of the top quark (p T (t)). As this observable is already accurate to NLO we find the expected good agreement between fixed-order prediction and the POWHEG-BOX as well as MG5 aMC@NLO. On the other hand, the Sherpa prediction captures the shape of the distribution well but underestimates the normalization by 13% as can be seen from the corresponding bottom panel. The reason for the discrepancy is that the Sherpa result is based on merging of tree-level matrix elements for ttW ± and ttW ± j and therefore misses higher-order loop corrections, as well as parts of the dominant real radiative corrections below the merging scale. The estimated theoretical uncertainties from scale variations is around ±20 − 25% for fixed-order, POWHEG-BOX, and  the normalization by about −12%. Nevertheless, within the scale uncertainties of each prediction, which amount to 16% at the beginning of the spectrum and 21% at the end of the plotted range, all distributions agree well with each other. Moreover, we find that the observable is very stable with respect to matching related parameters, as these uncertainties are almost negligible.
In Fig. 9 we further consider the transverse momentum (l.h.s.) and the pseudorapidity (r.h.s.) of the hardest jet. Contrary to the ttW ± QCD predictions we observe large differences between the various predictions. For the transverse momentum distribution of the leading jet we find that the POWHEG-BOX predicts a slightly softer spectrum than the fixed-order result with +12% corrections in the first bin, while for the remaining spectrum the curve is nearly a constant −5% below the fixed-order one. It is the parton shower evolution that softens the spectrum slightly and thus these corrections are of higher or- prediction in the very forward region, while it is closer to the POWHEG-BOX for the central rapidity bins. The significant shape differences also lie mostly outside of the estimated uncertainty bands. Let us remind the reader that, as illustrated in Fig. 5, the POWHEG-BOX reproduces the fixed-order distribution well when no shower evolution is taken into account.
Therefore, the significant corrections in the central rapidity region can be attributed to formally higher-order corrections generated by the parton shower. This can be also deduced from the dependence of the MG5 aMC@NLO results on the initial shower scale. We observe a large shape differences of around 20% in the central rapidity bins by varying the shower scale, which is comparable to the scale uncertainties.  Therefore, in order to have a visible effect of the EW contribution on differential distributions one has to focus on phase space regions where the EW production mode is enhanced with respect to the QCD one. In the following we show two representative observables that illustrate the impact of the EW contribution on combined differential distributions using our POWHEG-BOX implementation.
In Fig. 10 we show again the pseudorapidity of the hard- est jet, for the three considered Monte Carlo generators on the left-hand side while on the right-hand side the same distribution is shown as predicted by the ttW ± QCD contribution only and the combined ttW ± QCD+EW contribution. The overall agreement between the different predictions is good with shape differences of only 10% in the very forward region.
As can be deduced from the right plot of Fig. 10 the EW contribution becomes sizable in the forward region by modifying the shape of the distribution by nearly 80% and thus the observed shape differences between the various predictions are related to the modelling discrepancies in the EW contribution. Even though the theoretical uncertainties are dominated by missing higher-order corrections and amount to roughly ±15% the inclusion of the EW contribution represents a systematic shift of at least +5%. As a second example   we show in Fig. 11 the inclusive cross section as a function of the number of jets (N jets ), both light and b jets. For the first three bins we find excellent agreement between the three generators. Afterwards the predictions diverge and the POWHEG-BOX predicts the smallest while Sherpa the largest cross section for high jet multiplicities with differences as large as 60% for the jet bin with 6 or more jets. However, while the first three bins are dominated by uncertainties stemming from missing higher-order corrections the remaining bins are mostly affected by parton-shower effects whose uncertainty can be large. Indeed, while in the case of POWHEG-BOX the matching uncertainty does not exceed 5%, in the MG5 aMC@NLO case the spectrum is dominated by the dependence on the initial shower scale µ Q that can easily account for the aforementioned shape differences. Finally, we want to notice that even if the impact of the inclusive ttW ± EW contribution starts at +9% and increases with the number of jets to about a +20% correction on top of the QCD prediction for six and more jets, in this region the uncertainty of the QCD contribution overshadows this correction.
To summarize our findings of this section, we can say that ttW ± QCD inclusive production is rather robust with respect to matching uncertainties and different parton-shower algorithms. On the other hand, the ttW ± EW contribution is very sensitive to different matching procedures as sizable higher-order corrections can be generated by the parton shower even for inclusive observables. This however, might not be surprising all together in this particular case, since the NLO corrections in ttW ± EW production are highly dominated by real radiation matrix elements and therefore the description of this process is essentially only at LO accuracy. This can also be seen from the fact that the Sherpa prediction, which in the EW case is based only on tree-level matrix elements, still recovers the main features of many observables. In addition, the EW production channel is mediated by t-channel exchanges of color singlets and could be thus also sensitive to radiation patterns of parton shower implementations [73]. At last, we have to note that even though we observe large differences in the modelling of the ttW ± EW contribution these differences are much less visible once the QCD and EW contributions are combined.

C. Two same-sign leptons signature
In this section we focus on the experimental signature of two same-sign leptons in association with additional jets, usually denoted as "2 SS". The final state is selected by requiring exactly two same-sign leptons 5 with p T ( ) > 15 GeV and |η( )| < 2.5. Jets are formed using the anti-k T jet algorithm with a separation parameter of R = 0.4. We further require that jets (light as well as b jets) fulfill p T (j) > 25 GeV and |η(j)| < 2.5. Finally, we require to have at least two light jets as well as two tagged b jets.
Results shown in this section correspond to the sum of ttW ± QCD and ttW ± EW production modes which contribute to the 2 SS signature defined above. Emphasis will be placed on identifying the major sources of theoretical uncertainty from scale variation and parton-shower matching following the procedure discussed in section IV A.
For the range of parameters described in section IV A and the selection cuts above, we obtain fully consistent fiducial cross sections in all the frameworks considered in our study, 5 We exclude τ leptons here, since these typically form different signatures at hadron collider detectors.
Thus we focus on: e ± e ± , e ± µ ± , µ ± µ ± . We start by considering the fiducial cross section as a function of the number of light   (N l−jets ) or bottom (N b−jets ) jets, as shown in Fig. 12. In the case of light jets, shown on the left side of Fig. 12, we observe that MG5 aMC@NLO and Sherpa predictions align really well in all jet bins. Comparing to the POWHEG-BOX we note that predictions for at least 2, 3, and 4 light jets agree well among all generators, while starting from at least 5 light jets onward the POWHEG-BOX generates a softer spectrum and the deviation grows from 10% up to 47% for events with at least 8 light jets, although the shape differences are within the estimated theoretical uncertainties. Scale uncertainties are slightly asymmetric and start, after symmetrization, at 12% for at least two jets and increase to 25% for the POWHEG-BOX Turning to other exclusive observables we compare the modelling of fiducial differential distributions for the considered event generators starting with the invariant mass spectrum of the two hardest jets M j 1 j 2 depicted on the left of Fig. 13. It is remarkable to see the pronounced Breit-Wigner shape of the W boson, which suggests how the hardest jets in an event are predominantly generated by the hadronic decay of one of the W bosons. More information can be derived from a closer inspection of the p T distribution of these jets, as we will discuss in the following. With respect to the POWHEG-BOX curve, the MG5 aMC@NLO prediction for the M j 1 j 2 distribution shows a slightly different shape with deviations of the order of 5%. MG5 aMC@NLO is slightly more off-shell as can be seen by the depletion of events in the resonance region by about 5%. Sherpa on the other hand is considerable more on-shell with 32% more events on the W resonance. The scale uncertainties are nearly constant for MG5 aMC@NLO and of the order of 7% in the plotted range. For the POWHEG-BOX prediction scale uncertainties are constant and of the order of 8% below M j 1 j 2 ≈ M W , and increase mildly to 10% at the end of the spectrum. The matching uncertainties are for the MC@NLO or POWHEG scheme both at most 4%. The small impact of the shower scale variation can be understood from the fact that PYTHIA8 preserves the momentum of resonant decaying particles. Therefore, the available phase space for further radiation is naturally limited by a shower scale µ Q M W , which is considerable smaller than µ Q = H T /2 and thus the radiation pattern of a hadronic decaying W can to a large extent be independent of the variation of the initial shower scale we have chosen.
On the right-hand side of Fig. 13 we also show the overall transverse momentum distribution of the two hardest light jets p T (j 1 j 2 ). All three predictions agree remarkably well over the whole plotted range with only minor shape differences, with deviations of up to 5%. The Additionally, on the right-hand side of Fig. 15 we show the ∆R bb separation of the two hardest b jets. As can be seen from the shoulder in the distribution, the two leading b jets tend to be generated most of the time in a back-to-back configuration, which is typical for top-quark decays, whereas jets originating from a collinear g → bb splitting in the Scale uncertainties are the smallest around the peak at ∆R bb ≈ 3 where they amount to ±10%, while towards the beginning and the end of the spectrum they increase to about 18%. Matching uncertainties are for the POWHEG-BOX as well as MG5 aMC@NLO below 5% over the whole spectrum.

Leptonic observables
We now focus on leptonic observables, such as the transverse momentum p T ( ) and the invariant mass M spectrum of the two same-sign lepton pair that are shown in Fig. 16.
These observables are only indirectly affected by QCD corrections, because the leptons will only recoil against further emissions in the parton shower evolution. Nonetheless, these observables are crucial to illustrate the dynamical correlations between the decay of the  and increase up to 25% at the end of the spectrum. In addition, the spectrum is not very sensitive to the parton shower evolution as the corresponding uncertainties are below 5% for both MG5 aMC@NLO and the POWHEG-BOX.  Additionally, we show in Fig. 17 the azimuthal angle between the two same-sign leptons ∆φ which is sensitive to spin correlations. Let us note at this point that the observed spin correlations in our case have a different origin than those typically studied in topquark pair production [75][76][77][78][79][80]. In top-pair production the dilepton pair originating from the top-quark decays are spin correlated. In our case however, since we are considering signatures with same-sign leptons, only one lepton emerges from a top-quark decay, while the other one is produced by the decay of the prompt W boson. At the same time, it is exactly the presence of the prompt W boson that fully polarizes the top quarks [1], and as a consequence ∆Φ depends on the spin correlation between them. To highlight the importance of spin correlations, in the right-hand side plot of Fig. 17 we show POWHEG-BOX predictions that either include or neglect spin correlations in the decay modelling for the 2 SS final state. We observe that spin correlations lead to a shift in the spectrum, where events are shifted from the upper end of the spectrum at ∆Φ ≈ π to the opposite end.
The induced corrections to the distribution can reach up to 10% for ∆Φ /π 0.15. Similar effects are found in other leptonic observables (either dimensionless or dimensionful).

Impact of ttW ± EW contributions
We conclude this section by analyzing the impact of the EW production mode of the ttW ± final state on the 2 SS signature under consideration using our POWHEG-BOX implementation.
In this context, we would like to stress that, while we highlight in the following the few special cases where the ttW ± EW contribution becomes significant and cause visible shape modifications, for most observables the inclusion of the ttW ± EW process amounts to a rather flat +10% correction at the differential level. Since the EW production mode enhances the production of light jets in the forward region, as can be seen from the pseudorapidity distribution in Fig. 10, we expect the most severe impact in light-jet observables.   To illustrate this point, we show in the left-hand side plot of Fig. 18 the invariant mass distribution of the leading two light jets M j 1 j 2 as predicted by the ttW ± QCD contribution only and the total ttW ± QCD+EW contribution. The bottom panel shows the ratio with respect to only the ttW ± QCD contribution. Below M j 1 j 2 100 GeV we see that the ttW ± EW contribution is small, about a +8% correction, which can be understood from the fact that this region should be dominated by jets originating from the hadronically decaying W boson instead of jets emitted in the production process. However, above 100 GeV the ttW ± EW contribution starts to grow until it reaches a quite significant +25% correction at the end of the plotted range. Also in the ∆R j 1 j 2 separation between the leading light jets, shown on the right-hand side of Fig. 18, we can see a strong impact of the ttW ± EW contribution.
Below ∆R j 1 j 2 ≈ π the additional contribution is small between 8 − 12% and starts growing rapidly beyond that point. In fact, the ttW ± EW contribution generates corrections of the order of 50% at ∆R j 1 j 2 ≈ 5 which highlights the fact that this contribution preferably populates the forward regions.
To further highlight the impact of the EW contribution we study jets in the forward region. To this end we look for events that have passed the selection cuts specified in section IV A and have additional jets in the forward region defined by: where we do not distinguish between light or b jets.  FIG. 19. The inclusive cross section for the 2 SS fiducial region as a function of the number of hard jets in the central phase space volume (l.h.s.) and in the forward region (r.h.s.). The predictions based only on ttW ± QCD production are given in red while the total ttW ± QCD+EW contributions are shown in blue. The bottom panel shows the percentage change in the shape of the distribution.
In Fig. 19 we show the inclusive cross section as a function of the number of jets in the forward region (r.h.s.) and contrast it with the corresponding plot in the central region (l.h.s.). In the central phase space volume the ttW ± EW contribution only has a very mild impact, and, for events with at least two light jets, this contribution amounts to 11%, while for events with at least eight light jets it grows to 19%. On the contrary, if we look at the inclusive jet multiplicities in the forward region as shown on the right plot in Fig. 19 then we see that for the first bin a 11% corrections is visible that then quickly increases up to 66% for events that have at least four additional forward jets. Last we study the impact on the transverse momentum spectrum of the leading light jet, p T (j 1 ), in the central and forward regions as depicted in Fig. 20. As before, in the central phase space volume we find a rather constant correction of 9−11% over the whole range of the plotted spectrum. The maximal corrections of 11% are obtained around p T ≈ 200 GeV. On the other hand, the hardest forward jet receives large corrections from ttW ± EW production.
Essentially starting from 17% corrections at the beginning of the distribution the ttW ± EW contribution gives rise to 100% corrections for p T 200 GeV.

V. CONCLUSIONS
In this paper we have presented a new NLO parton-shower Monte Carlo event generator for the hadronic production of a top-quark pair in association with a W ± boson taking into account the dominant NLO corrections at O(α 3 s α) and O(α s α 3 ). Decays of unstable particles are included at leading order retaining spin-correlation effects. The POWHEG event generator Wtt dec is publicly available as part of the POWHEG-BOX repository under http://powhegbox.mib.infn.it Motivated by the current tension [6] between the state-of-the-art SM predictions for the ttW ± cross section and the corresponding measurement derived when ttW ± is extracted from a combined signal and background fit in ttH analyses, we performed a detailed generator comparison involving the POWHEG-BOX, MG5 aMC@NLO, and Sherpa. A comparison at the level of on-shell ttW ± production has revealed good agreement for the O(α 3 s α) production mode, while the O(α s α 3 ) contribution is very sensitive to the details of a given generator setup.
Furthermore, a comparison has been made at the fully decayed stage for a fiducial phase space volume corresponding to the 2 SS signature with two same-sign leptons and both light and b jets. We provide this as a proof of concept to estimate in a more robust way the residual theoretical uncertainty on the pp → ttW ± cross section from both fixed-order and parton-shower components, at the inclusive and fiducial level. We found good agreement between all three generators at the differential level. Theoretical uncertainties have been addressed via means of independent variations of the renormalization and factorization scales as well as matching related parameters specific of each generator, such as damping factors (POWHEG-BOX) or the initial shower scale (MG5 aMC@NLO). The investigation of the latter dependence allowed us to explain peculiar shape differences between the generators. In addition, we also investigated the impact of LO accurate spin-correlated decays on the differential distributions and found that leptonic observables are particularly sensitive to spin-correlation effects. Finally, we also quantified the impact of the O(α s α 3 ) contributions at the differential level. For most observables these contributions amount to a flat +10% correction, while for a few observables sensitive to forward jets they can become more sizable.
It is interesting to notice that even though the on-shell modelling of ttW ± EW production largely depends on the matching and parton-shower settings of each generator, these effects are much less visible for full predictions once QCD and EW contributions are combined.
Furthermore, giving the fact that the two same-sign lepton signature is dominated by jets emerging from a hadronic W decay that in all generators is modelled only at LO, for this particular signature it will be very important to improve on the theoretical description of these jets in the future. With respect to this, the modelling of the two same-sign lepton signature could be improved by a fixed-order full off-shell computation similar to Ref. [81], as the corresponding matching to parton showers would be computationally challenging. The narrow-width approximation presents an alternative to include one-loop QCD corrections to the top-quark and W boson decays and could be included in an event generator, as has been already shown in Ref. [82]. On the other hand, other fiducial signatures, like the ones involving three leptons and no hadronic W decay, may require the inclusion of higher-order corrections (like O(α 4 s α) NNLO corrections) to the production process to reach a better control of the corresponding theoretical systematics. In all cases the impact of radiative top-quark decays should be studied carefully. Finally, we want to emphasize the fact that even though for the current center-of-mass energy of the LHC of √ s = 13 TeV the modelling of the subleading O(α s α 3 ) contribution only plays a minor role at the fiducial level it will be crucial to improve on its theoretical accuracy for higher center-of-mass energies, as its radiative contribution will increase [12].
It is clear from our discussion that providing a robust theoretical prediction for hadronic ttW ± production cannot be framed as a unique recipe and care must be taken to analyze the specific characteristics of different observables measured in experiments. Having at our disposal several well tested tools that allow to implement state-of-the-art theoretical calculations in the modelling of collider events is clearly valuable and offers us the possibility of studying the problem in its complexity and identify where improvement is most needed.
which defines y as a measure of the virtuality The upper boundary on the virtuality is given by the value y max for which λ vanishes and is given by In order to preserve momentum conservation the recoiling momenta have to be boosted with the Lorentz transformation