Searching for heavy Higgs bosons in the $t \bar t Z$ and $t b W$ final states

In the context of two-Higgs doublet models, we explore the possibility of searching for heavy Higgs bosons in the $t \bar t Z$ and $t bW$ final states. We develop realistic analyses strategies and in the case of the $t \bar t Z$ channel provide a detailed evaluation of the new-physics reach at the 14 TeV LHC. We find that already with an integrated luminosity of $300 \, {\rm fb}^{-1}$ searches for the $t \bar t Z$ signature can provide statistically significant constraints at low values of $\tan \beta$ for heavy Higgs masses in the range from around $450 \, {\rm GeV}$ to $1150 \, {\rm GeV}$. Future searches for heavy Higgses in the $tbW$ final state are also expected to be able to probe parts of this parameter space, though the precise constraints turn out to depend sensitively on the assumed systematics on the shape of the $t \bar t$ background.


Introduction
The most important accomplishment of the LHC Run-1 physics programme has been the discovery of a new spin-0 resonance h with a mass of around 125 GeV in 2012 [1,2]. In the last five years the LHC Higgs programme has matured, providing precise measurements of processes such as pp → h → γγ and pp → h → ZZ * → + − + − (see [3][4][5][6] for the latest LHC results at √ s = 13 TeV) with standard model (SM) rates of around 100 fb and 5 fb, respectively.
The finding that the 125 GeV spin-0 resonance has properties close to the one expected for the SM Higgs [7] implies that if additional Higgs bosons exist in nature such states can only be slightly mixed with the h. An extended Higgs sector has thus to be approximately aligned, either via decoupling or via alignment without decoupling. While in the former case the extra spin-0 particles might be too heavy to be accessible at the LHC, in the latter case the additional Higgs bosons can have masses at or not far above the electroweak (EW) scale without being in conflict with any other observation. In the case of alignment without decoupling, direct searches for extra Higgslike particles are hence particularly well-motivated as they can provide complementary information with respect to the LHC programme of precision Higgs measurements.
All the channels mentioned so far have in common that they only have limited sensitivity to additional Higgses with masses above the top threshold, in particular if the H/A → tt branching ratio is sizeable as it happens to be the case in the MSSM at low and moderate tan β. In order to gain sensitivity to new-physics scenarios of the latter kind the channels pp → H/A → tt, pp → ttH/A → tttt and pp → bbH/A → bbtt have been proposed (see [31][32][33][34][35] for instance) and first experimental searches for the tt [36], tttt [37,38] and bbtt [38] final state have been carried out recently. While, at first sight, all three signatures seem to offer good prospects for probing heavy Higgs bosons, it turns out that in practice they all suffer certain limitations. In the case of pp → H/A → tt, interference effects between the signal and the SM tt background [39][40][41][42][43] represent a serious obstacle, while for what concerns the searches for pp → ttH/A → tttt and pp → bbH/A → bbtt the small signal-over-background ratio is in general an issue. As a result, a very good experimental understanding of the systematic uncertainties plaguing the overwhelming tt background is crucial in order for the tt, tttt and bbtt final states to provide statistically significant constraints at low to moderate values of tan β.
In this work two novel search strategies for neutral Higgs particles with masses above the top-quark threshold are devised. The first strategy exploits the ttZ final state and is based on the isolation of the irreducible ttZ process from other SM backgrounds, followed by the discrimination of the signal from SM ttZ production using the distinctive kinematic features of the new-physics signal. As a compromise between purity and statistics, we consider final states where the Z boson decays into charged leptons, and only one of the two top quarks decays semileptonically. The examined final state thus involves three charged leptons, i.e. a pair of same-flavour leptons compatible with the Z → + − decay and one charged lepton from t → bW → b ν, missing transverse energy (E T,miss ) associated to the neutrinos from top decays and four jets, two of which are produced via bottom-quark fragmentation. The invariant masses of the ttZ and tt systems (m ttZ and m tt ) can be experimentally reconstructed and their distributions are peaked at the masses of the heavy Higgs bosons appearing in H/A → A/HZ → ttZ. To separate signal from background a shape fit to the distribution of the variable ∆m ≡ m ttZ − m tt can be used, since this spectrum is smoothly falling with ∆m in the case of the SM ttZ background. An observable that provides additional information while being well measurable is the Z-boson transverse momentum (p T,Z ). The shape of the p T,Z spectrum of the ttZ signal is in fact predicted to be Jacobian with an endpoint that is related to the structure of the massive two-body phase describing the H/A → A/HZ transition. We will show that by using the experimental informations on ∆m and p T,Z , future searches for the ttZ final state should allow to set unique bounds on the part of the 2HDM parameter space that features heavy Higgses with masses above the top threshold and small values of tan β.
Our second search strategy targets the tbW final state. We point out that there is a number of kinematic handles that can be used to separate the tbW signal from the tt background. In the case of the two-lepton final state, one can exploit the invariant masses of the bb and bb systems (m bb and m bb ) since the corresponding distributions have kinematic endpoints, while in the one-lepton final state the Breit-Wigner peaks in the invariant mass spectra of the tb and tbW systems (m tb and m tbW ) can be harnessed. Based on these observations, we sketch the main ingredients of an actual two-lepton analysis. In our exploratory study, the signal-over-background ratio however turns out to be at most a few percent for the considered 2HDM realisations, making it difficult to determine the precise LHC reach of the proposed two-lepton analysis. Similar statements also apply to the one-lepton case. A full exploration of the potential of the tbW final state is therefore left to the ATLAS and CMS collaborations once they have collected data in excess of 300 fb −1 .
The outline of this article is as follows. Section 2 discusses the structure of the relevant Higgs interactions and the resulting decay modes, while the anatomy of the ttZ and tbW signal is studied in Section 3 and 4, respectively. A concise description of our Monte Carlo (MC) generation and detector simulation is given in Section 5. The actual analysis strategies are detailed in Sections 6 and 7. In Section 8 we present our numerical results and examine the new-physics sensitivity of the ttZ signature at upcoming LHC runs. We conclude in Section 9. Supplementary material is provided in Appendices A and B.

Heavy Higgs interactions and decays
The addition of the second Higgs doublet in 2HDMs leads to five physical spin-0 states: two neutral CP-even ones (h and H), one neutral CP-odd state (A), and the remaining two carry electric charge of ±1 and are degenerate in mass (H ± ). Following standard practice, we identify the 125 GeV resonance discovered at the LHC with the h field, denote the angle that mixes the neutral CP-even states by α, and define tan β to be the ratio of the Higgs vacuum expectation values (VEVs).
The tree-level couplings of the Higgses h, H to EW gauge bosons satisfy in all 2HDMs with a CP-conserving Higgs potential the relations where V = W, Z and we have used the shorthand notation s β−α ≡ sin (β − α) and c β−α ≡ cos (β − α). Notice that in the so-called alignment limit, i.e. α → β − π/2, the interactions of h with EW gauge bosons resembles those in the SM while the couplings between H and W-boson or Z-boson pairs vanish identically. The consistency of the LHC Higgs measurements with SM predictions requires that any 2HDM Higgs sector is close to the alignment limit, meaning that small values of c β−α are experimentally favoured. The combinations of mixing angles appearing in (2.1) also govern the interactions between two Higgses and one EW gauge boson. Explicitly one has Notice that the first relation leads to a suppression of the A → hZ decay rate in the alignment limit.
In contrast, the decay rate H → AZ (A → HZ) is unsuppressed for c β−α → 0 and can be large if this  . Like g HAZ also g HH ± W ∓ and g AH ± W ∓ are non-vanishing in the alignment limit, and in consequence the decays H → H ± W ∓ and A → H ± W ∓ are phenomenologically relevant if they are open.
In order to tame dangerous tree-level flavour-changing neutral currents the Yukawa interactions in 2HDMs have to satisfy the natural flavour conservation hypothesis [44,45]. Depending on which fermions couple to which Higgs doublet, one can divide the resulting 2HDMs into four different types. While the Higgs couplings to light fermions turn out to be model dependent, the couplings of h, H and A to top quarks take in all four cases the generic form where we have introduced the abbreviation t β ≡ tan β. These expressions imply that in the alignment limit the coupling of h to top quarks becomes SM-like, while the top couplings of H, A are both t β suppressed. The only charged Higgs coupling to fermions relevant to our work is the one to right-handed anti-top and left-handed bottom quarks. This coupling resembles the form of g Att , and in consequence the charged Higgs decays dominantly via H + → tb if this channel is open. The magnitudes of 2HDM couplings that involve more than two Higgses depend on the precise structure of the full scalar potential. For what concerns the coupling g Hhh that describes the selfcoupling between a H and two h, it turns out that it is homogenous in c β−α , and therefore vanishes in the alignment limit in pure 2HDMs (see [46] for example). In fact, in the limit c β−α → 0 and M H ± > M H > v, M h with v 246 GeV the Higgs VEV and assuming that the quartic couplings λ i that appear in the scalar potential are of order 1, the g Hhh coupling behaves approximately as g Hhh ∝ c β−α M 2 H ± /v. It follows that for a sufficiently large mass splitting M H ± − M H > 0, the partial decay width Γ (H → hh) ∝ g 2 Hhh /M H can be numerically relevant in pure 2HDMs. In contrast, in the MSSM the trilinear Hhh coupling scales as g Hhh ∝ M 2 Z /v s 4β in the limit α → β − π/2. The coupling g Hhh is hence non-zero in the alignment (or decoupling) limit of the MSSM, but since Htt M H , the branching ratio of H → hh is always small for Higgs masses M H sufficiently above the top threshold.
The above discussion suggests that close to the alignment limit the decay pattern of the heavy Higgses H and A is rather simple in all 2HDMs. To corroborate this statement we show in Figure 1 the branching ratios of H and A for four type-II 2HDM benchmark models. The different benchmarks thereby cover values of c β−α that range from the pure alignment limit c β−α = 0 to the case of maximally allowed misalignment, which amounts to around c β−α = 0.15 in the type-II 2HDM after LHC Run-1 (see for instance [9]). Our calculation of the branching ratios is based on the formulas and results given in [47][48][49][50][51][52]. From the upper left panel one observes that for M H 600 GeV the decay H → tt almost fully saturates the total width of H, while for M H 600 GeV the decay mode H → AZ becomes important quickly and even dominant for M H 800 GeV. A similar picture arises in the case of the A with A → tt and A → HZ representing the two dominant decay modes for M A 600 GeV. This feature is illustrated in the upper right panel in Figure 1.
Notice that to obtain the latter plots we have fixed M H ± = M H and M H ± = M A , respectively. These choices are well-motivated, because only in these two cases [53][54][55][56][57] can the H or the A have a sizeable mass splitting from the rest of the non-SM Higgses without being in conflict with EW precision measurements. The left (right) panel shown in the lower row of Figure 1  fourth largest branching ratio for heavy CP-even Higgses H. We add that the results shown in the latter panel correspond to the choice λ 3 = 3, where λ 3 is the quartic coupling that multiplies the term |H 1 | 2 |H 2 | 2 in the 2HDM scalar potential, and H 1 and H 2 denote the two Higgs doublets in the Z 2 basis.

Anatomy of the tt Z signature
The discussion in the last section singles out the ttZ and tbW final states as promising to search for the presence of heavy Higgs particles. Prototypes of Feynman diagrams that lead to the former signal in 2HDMs are shown in the upper row of Figure 2. In the graph on the left-hand side a H is produced in association with a Z boson from a top-quark box, while in the right diagram the H is emitted from a top-quark triangle and then decays via H → AZ → ttZ. Graphs where the role of the neutral Higgses H and A is interchanged also contribute to the ttZ signature in 2HDMs but are not explicitly shown in the figure.
In order to understand the anatomy of the ttZ signal in the 2HDM context, one first has to notice that the upper right Feynman diagram in Figure 2 allows for resonant ttZ production if the two conditions M H > M A + M Z and M A > 2m t are satisfied. Once the channels H → AZ and A → tt are kinematically accessible the triangle graph therefore always dominates over the box contribution displayed on the upper left of the latter figure. In fact, the dominance of the triangle contribution allows one to estimate the signal strength s (pp → ttZ). Since in the narrowwidth approximation (NWA) the signal strength factorises into production and decay and given that BR (A → tt) 100% for the parameters of interest, one obtains in the case of pp → H → AZ → ttZ the following approximate result (3.1) If the ttZ signature arises instead through A → HZ, the role of H and A has to simply be interchanged. The total H, A production cross sections appearing in (3.1) are easy to calculate at leading order (LO). In the exact alignment limit and assuming that t β is not too large, we obtain at √ s = 14 TeV the following expressions These approximations work to better than 20% for Higgs masses in the range of [400, 1000] GeV. They imply that the total production cross section of a A is always larger than that of a H if these particles have the same mass. For M H = M A = 600 GeV the relations (3.2) predict for example an enhancement factor of around 1.6. We emphasise that the formulas given in (3.2) serve mostly an illustrative purpose and have only been used to obtain the approximate signal strengths for pp → ttZ and pp → tbW production as shown in Figures For relative decay widths below the quoted value the NWA should be applicable.
The resonant contributions not only enhance the ttZ signal cross section, but also lead to interesting kinematic features that one can harness to discriminate signal from background. Firstly, since both heavy Higgses tend to be on-shell in the production chain pp → H → AZ → ttZ, the invariant masses m ttZ and m tt of the ttZ and tt systems show characteristic Breit-Wigner peaks at The difference ∆m between m ttZ and m tt can therefore be used to determine the mass splitting of the heavy Higgses. In the considered case, the ∆m distribution of the ttZ signal will for instance be peaked at Second, since the four-momenta of the decay products A and Z that enter H → AZ are fixed by H being preferentially on-shell, also the p T,Z spectrum will have a characteristic shape. In fact, it is straightforward to show that the p T,Z distribution of the resulting ttZ signal is a steeply rising function of p T,Z with a cut-off at that is smeared by the total decay width Γ H of the heavy Higgs H. Needless to say, that the same line of reasoning and formulas similar to (3.3), (3.4) and (3.5) apply when one considers the process Figure 4 shows the ∆m and p T,Z distribution that we obtain from a MadGraph5_aMCNLO [58] simulation of the new-physics contribution to the ttZ final state. The displayed results have been obtained in the context of the type-II 2HDM by employing a modified UFO implementation [59] of the 2HDM discussed in [60]. The chosen parameters are M H = M H ± = 800 GeV, M A = 500 GeV, c β−α = 0 and t β = 1. The red predictions correspond to the pure new-physics signal (S), while the black distributions take into account the interference between the signal process and the background from SM ttZ production (S+I). All relevant box and triangle diagrams have been included in our simulation. The distinctive kinematic features of the signal discussed earlier are clearly visible in the two panels with the ∆m distribution peaked at about 300 GeV and an edge in the p T,Z spectrum  at around 230 GeV. One also observes that, in contrast to the case of tt production [39][40][41][42][43], signalbackground interference leads only to minor distortions of the shapes of the most interesting ttZ distributions. Although the interference effects are observed to be small (roughly of the size of the statistical uncertainties in the shown example), we will include them in Section 8 when determining the sensitivity of the ttZ signature in constraining the parameter space of 2HDMs.

Anatomy of the t bW signature
Two example diagrams that gives rise to a tbW signal through the exchange of a charged Higgs boson are displayed in the lower row of Figure 2. In the left graph a H + and a W − are radiated off a box diagram with internal top and bottom quarks, while in the diagram on the right-hand side a H is emitted from a top-quark triangle which then decays via H → H + W − . In both cases the charged Higgs boson decays to a tb pair. Notice that diagrams with H − or A exchange also lead to a tbW signal. These contributions while not explicitly shown in the lower row of Figure 2 are all included in our analysis.
The tbW final state can be resonantly produced via In such a case triangle diagrams provide the leading contribution to the tbW signal strength. In Figure 5   600 GeV and M A,H 400 GeV. Since for M H = M A the total production cross section σ (pp → A) is bigger than σ (pp → H), one again notices a small asymmetry between the signal strengths pp → A → H ± W ∓ → tbW and pp → H → H ± W ∓ → tbW with the former being always slightly larger than the latter.
Like the case of the ttZ signal, also the kinematic distributions of the tbW signature have distinctive features that can be exploited to tame SM backgrounds. Figure 6 shows an assortment of invariant mass distributions that can serve this purpose. The displayed results have been obtained in the type-II 2HDM using MadGraph5_aMCNLO. The choice of parameters is M H = 800 GeV, M A = M H ± = 400 GeV, c β−α = 0 and t β = 1. The left panel in the upper row of the figure depicts the invariant mass m bb of the bb system in pp → tbW → bbW + W − . One sees that the m bb distribution has a sharp edge at around 320 GeV, which corresponds to the kinematic endpoint [61,62] Similarly, also the invariant mass m bb of the bb final-state configuration that appears in the tbW channel from the sequential decay H ± → tb → bbW → bb ν has a kinematic endpoint. It is located at [61,62]   The associated edge in the m bb spectrum arises at about 390 GeV, a feature that is evident in the upper right panel of Figure 6. Notice that bb final states also arise from the leptonic decay of the W bosons involved in H/A → H ± W ∓ . The corresponding invariant mass m bb has a very soft endpoint at m bb M H/A − M W and no edge because the lepton is not emitted directly from the backbone of the whole decay chain. Since experimentally one can separate the two cases (see Section 7), an example of a m bb signal distribution has not been depicted in Figure 6.
In the lower left panel the invariant mass m tb of the tb system is depicted. As expected, this distribution shows a Breit-Wigner peak at The invariant mass m tbW of the tbW final state is displayed in the lower right panel of Figure 6. The two mass peaks at Realise that not only the process gg → H ± W ∓ → tbW (example diagrams are shown in the lower row of Figure 2) but also graphs corresponding to gg → H/A → tt → tbW contribute to the tbW signature in 2HDMs. To separate the charged Higgs contributions to the tbW channel from the neutral Higgs contributions associated to tt production, we employ the so-called diagram removal (DR) procedure [63]. In this scheme the tbW final state is defined by removing from the tbW scattering amplitude all doubly resonant diagrams, i.e. graphs in which the intermediate top quarks can be on-shell. Singly resonant contributions are on the other hand kept. The DR procedure is also applied to the SM amplitudes, and as a result only tb-fusion (but no top-fragmentation) diagrams contribute at LO in QCD to the tbW final state.
Based on the DR definition of the tbW final state, we have studied the impact of signalbackground interference. The red curves in Figure 6 correspond to the pure new-physics signal (S), while the black distributions in the four panels take into account the interference between the signal process and the background from tbW production within the SM (S+I). Comparing the two sets of histograms, we observe that the kinematic features in the m bb , m bb , m tb and m tbW distributions are always less pronounced for the S+I predictions than the S results. Since the size of the signalbackground interference typically exceeds the statistical uncertainties expected in future LHC runs, a rigorous assessment of the prospects of the tbW final state to search for heavy Higgses should be based on MC simulations that include interference effects between the new-physics signal and the SM background.

MC generation and detector simulation
In our study we consider throughout pp collisions at √ s = 14 TeV. We generate the signal samples using a modified version of the Pseudoscalar_2HDM UFO together with MadGraph5_aMC@NLO and NNPDF23_lo_as_0130 parton distribution functions [64]. Compared to the UFO presented in [60] our new implementation is able to calculate the interference between the loop-induced ttZ and tbW signals and the corresponding tree-level SM backgrounds. The obtained parton-level events are then decayed and showered with PYTHIA 8.2 [65] which allows us to study the fully interfered signals and backgrounds at the detector level.
Our ttZ analysis will address the three-lepton final state, with two opposite-sign same-flavour leptons from the Z-boson decay and one lepton from the semileptonic decay of one of the two top quarks. For the description of the SM backgrounds to this final state, SM processes involving at least three leptons coming from the decay of EW gauge bosons are simulated. Most of the backgrounds are generated at LO with MadGraph5_aMC@NLO. The dominant irreducible background is ttZ which is generated with up to an additional jet. The tbW background is instead generated with up to two additional jets. The dominant diboson background, i.e. WZ, is simulated with up to three additional jets. The minor backgrounds considered are tZ and tWZ both of which are obtained using MadGraph5_aMC@NLO. In each case the decay of the top quarks and the EW gauge bosons is performed by MadGraph5_aMC@NLO. The reducible ZZ background is generated at next-to-leading order (NLO) with POWHEG BOX [66].
A potential significant background in the case of the ttZ channel arises from processes where two leptons are produced in the decays of EW gauge bosons, and a third lepton is either the result of a misidentification in the detector or the decay of a B meson. The latter is experimentally strongly suppressed by requiring the leptons to be isolated. The estimate of these backgrounds requires a profound understanding of the detector performance, and indeed the ATLAS and CMS use data-driven techniques rather than MC simulations to determine them. A recently published search for ttZ by the ATLAS experiment [67] shows that the requirement of having at least four jets of which two are identified as coming from the fragmentation of bottom quarks reduces the background from misidentified leptons to a level well below the other backgrounds.
In the case of the tbW analysis, the background evaluation is performed through the generation of SM processes involving at least two leptons coming from the decays of EW gauge bosons. The backgrounds from tt [68], tW [69], WW, WZ and ZZ production [70,71] were all generated at NLO with POWHEG BOX. The Z + jets sample is generated at LO with MadGraph5_aMC@NLO, considering up to four jets for the matrix element calculation. The latter MC code is also used to simulate the ttV backgrounds with V = W, Z at LO with a multiplicity of up to two jets, and the tZ and tWZ backgrounds at LO. As for the ttZ analysis, we do not consider final states where one or both of the leptons are either fake electrons from jet misidentification or real non-isolated leptons from the decay of heavy flavours.
All partonic events are showered with PYTHIA 8.2 and the SM backgrounds are normalised to their NLO cross section calculated either with MadGraph5_aMC@NLO or with POWHEG BOX where relevant. The simulated analyses are performed on experimentally identified electrons, muons, photons, jets and E T,miss which are constructed from the stable particles in the generator output. Jets are constructed by clustering the true momenta of all the particles interacting in the calorimeters, with the exception of muons. An anti-k t algorithm [72] with a parameter R = 0.4 is used, as implemented in FastJet [73]. Jets originating from the hadronisation of bottom quarks (b-jets) are experimentally tagged in the detector (b-tagged). The variable p T,miss with magnitude E T,miss is defined at truth level, i.e. before applying detector effects, as the vector sum of the transverse momenta of all the invisible particles (neutrinos in our case). The effect of the detector on the kinematic quantities used in the analysis is simulated by applying a Gaussian smearing to the momenta of the different reconstructed objects and reconstruction and tagging efficiency factors. The parametrisation of the smearing and of the reconstruction and tagging efficiencies is tuned to mimic the performance of the ATLAS detector [74,75] and is applied as a function of the momentum and the pseudorapidity of the physical objects. The discrimination of the signal from the background is significantly affected by the experimental smearing assumed for E T,miss . To simulate this effect, the transverse momenta of unsmeared electrons, muons and jets are subtracted from the truth E T,miss and replaced by the corresponding smeared quantities. The residual truth imbalance is then smeared as a function of the scalar sum of the transverse momenta of the particles not assigned to electrons or jets. The same techniques have also been employed in [76,77]. 6 Analysis strategy for the tt Z signature In the case of the ttZ channel the generated signal and background events are preselected by requiring exactly three charged leptons (electrons or muons) with a pseudorapidity of |η | < 2.5. The leading lepton must have p T, > 25 GeV, while the other two are required to satisfy p T, > 20 GeV. At least one pair of leptons of opposite charge and same flavour must be present, and the invariant mass m of this pair must meet the requirement |m − M Z | < 15 GeV. In case the event includes more than one such lepton pair, the pair with the invariant mass closest to the nominal value of M Z is selected as the Z-boson candidate. All events furthermore need to contain four jets with p T, j > 20 GeV and |η j | < 2.5, of which two must be tagged as bottom-quark jets (b-tagged).
Notice that the large jet multiplicity and the requirement of having two b-tagged jets leads, on the one hand, to a strong reduction of the WZ and tWZ SM backgrounds, and on the other selects all objects needed for a full reconstruction of the event. The leptonically decaying W boson is reconstructed from the charged lepton not assigned to Z → + − and the amount of E T,miss by solving the W-boson mass constraint. If the solution for the longitudinal momentum of the neutrino is imaginary, the real part of p z,ν is taken. The reconstructed jets are assigned to the different W and top decays, by choosing the assignment which gives the best compatibility with the decay of two top quarks in terms of reconstructed masses. If in the process two solutions are found for the reconstructed leptonic W decay, the one giving the best mass compatibility is selected.
In Figure 7 the distributions of ∆m and p T,Z for the SM backgrounds and two type-II 2HDM benchmark models after applying the selections described above are displayed. Our benchmarks correspond to M H = M H ± = 600 GeV, M A = 400 GeV (M H = M H ± = 1000 GeV, M A = 500 GeV), c β−α = 0 and t β = 1 and are indicated by the dashed (dotted) black lines. In the case of the ∆m spectrum, one observes from the left panel that the sum of the SM backgrounds is a steeply falling distribution, while both new-physics ttZ signals exhibit a Breit-Wigner peak. In fact, as expected from (3.4) the peaks are located at around 200 GeV and 500 GeV. Our results for the p T,Z distributions are presented in the right panel of the latter figure. In agreement with (3.5) the two 2HDM benchmark models lead to spectra that show distinctive Jacobian peaks with edges at roughly 150 GeV and 370 GeV. The SM backgrounds are in contrast again smoothly falling and featureless. Notice that a measurement of the p T,Z distribution in ttZ production does, unlike a measurement of the difference ∆m of invariant masses, not require the full reconstruction of the final state. As a result, p T,Z is less subject to experimental uncertainties than ∆m. In order to stress the experimental robustness of the proposed ttZ signature, we will in our sensitivity study consider both the ∆m and the p T,Z distribution as final discriminants.

Analysis strategy for the t bW signature
In the case of the pp → H ± W ∓ → tbW → bbWW signal the dominant QCD backgrounds are tt production and tW production in association with a b-jet. By vetoing events where the observed W bosons and b-jets are kinematically compatible with the decay of two top quarks the overwhelming tt background can be reduced by approximately two orders of magnitude, making it comparable to the tW background in size. After this selection the signal is however still two orders of magnitude smaller than the background. Notice that this is in contrast to the ttZ channel where the signal and the background are of the same size after background suppression. To improve the signal-over-background ratio in the case of the tbW signal, one needs to exploit the decay kinematics of the heavy Higgses by identifying the decay products of the top quark in the signal events. The invariant mass of the top quark with the additional b-jet will be peaked at M H ± , while the invariant mass of the two b-jets and the two W bosons equals M H or M A depending on which mass is larger. Experimentally the signal can be looked for in events with two, one or zero isolated charged leptons resulting from W → ν. In the following we will sketch a possible analysis procedure for the two-lepton final state. Given the small signal-to-background ratio for the irreducible backgrounds, we however expect that our conclusions will be valid for the one-lepton final state as well.
For the two-lepton case, the reconstruction of mass peaks is not possible due to the presence of two neutrinos in each event. However, given the presence of multi-step sequential decays leading to undetected neutrinos, the invariant mass distributions of the visible decay products are bounded from above [78]. In the case of top decays, the invariant mass m b of the resulting b-quark and lepton must be lower than m max b = m 2 t − M 2 W 153 GeV. Thus exactly two opposite-sign leptons ( 1 and 2 ) with p T, 1 > 30 GeV and p T, 2 > 25 GeV and exactly two b-tagged jets (b 1 and b 2 ) with p T,b > 30 GeV are required in the event, and events are selected in which none of the two possible pairings among b-jets and leptons is compatible with the decay of two top quarks. A convenient way of rejecting events compatible with two top decays consists in introducing the observable where the minimisation runs over all pairs { j a , j b } of distinct jets inside a predefined set of test jets. Based on the number of b-tagged jets in the event, the set of test jets is defined as follows. If the event includes one or two b-tagged jets, an additional test jet is considered, chosen as the non-btagged jet with the highest b-tagging weight and p T, j > 25 GeV. If three b-tagged jets are found, they are all taken as test jets. The requirement m t b > 180 GeV suppresses the tt background by approximately two orders of magnitude. A dangerous background is also also due to the production of a Z boson in association with b-jets. Vetoing same-flavour lepton pairs compatible with a Zboson decay and requiring some E T,miss reduces this background to roughly the level of the signal, albeit with large uncertainties. To avoid this possible issue we completely remove the latter and all other backgrounds including a real Z boson by requiring that the two selected leptons have different flavours. After these selections the remaining background consists of approximately one half of tt and one half of tW events.
A further separation of signal from background can be achieved by exploiting the fact that in the case of the signal the invariant mass m bb of the two b-tagged jets as well as the invariant mass m bb of the two b-tagged jets with the lepton from a top decay are bounded from above. See (4.1) and (4.2). In order to illustrate this point we show in the left panel of Figure 8 the distribution of m bb for two signal samples with M A = M H ± = 400 GeV (dashed black line) and M A = M H ± = 600 GeV (dotted black line), respectively. The remaining 2HDM parameters are set to M H = 800 GeV, c β−α = 0 and t β = 1 and the background has been scaled down by a factor of 100 for better visibility. Upper cuts on m bb matching the kinematic endpoint of the signal for different values of M H ± will improve the signal-over-background ratio, bringing it to a level of at most 3% over the parameter space relevant for this analysis. The variable m bb is less effective as it has a less sharp edge, and suffers from an ambiguity in the choice of the lepton.
The final experimental handle is the fact that the invariant mass of the tbW system will peak at the mass of the H/A bosons for the bulk of the signal (cf. Figure 6) whereas for the background m tbW has a broad distribution centred at around 400 GeV. However, if the charged Higgs mass is not known, the observable m tbW cannot be reconstructed because of the two undetected neutrinos. A large literature on the reconstruction of the mass of new particles in events with two invisible particles in the final state is available (see e.g. [79] for a review). In our exploratory study, we employ the variable m eff as an estimator of m tbW , which is defined as [78,80] m eff ≡ a= 1 , 2 ,b 1 ,b 2 p T,a + E T,miss . In the right panel of Figure 8 we show the m eff distribution for M H = 600 GeV (dashed black line) and M H = 800 GeV (dotted black line), fixing the other parameters to M A = M H ± = 400 GeV, c β−α = 0 and t β = 1. For better visibility the background has been scaled down by a factor 50 after applying the cut m bb < 280 GeV. The significance of the signal can be extracted from a shape fit to the m eff spectrum for signal and background. Given the difference in shape between signal and background, and the large number of kinematic handles, it should be possible to extract a significant signal for a signal sample corresponding to 3 ab −1 of integrated luminosity if the shape of the background can be experimentally controlled to a level below 2%. In these conditions a reliable evaluation of the coverage in parameter space can only be performed by the experimental collaborations. It is however worth noting that, due to the shape of the m bb background distribution, the maximum sensitivity of the tbW analysis is expected to arise for M A = M H ± 400 GeV, making the tbW coverage complementary to that of the ttZ search.

Numerical results
Based on the search strategy outlined in Section 6, we now study the sensitivity of future LHC runs to the ttZ signature. To evaluate the upper limit on the ratio of the signal yield to that pre-dicted in the 2HDM framework, a profiled likelihood ratio test statistic applied to the shapes of the ∆m and p T,Z distributions is used. The CLs method [81] is employed to derive exclusion limits at 95% confidence level (CL). The statistical analysis has been performed by employing the RooStat toolkit [82]. A systematic uncertainty on the absolute normalisation of the SM background (signal) of 15% (5%) is assumed. This choice of uncertainties is in accordance with the uncertainties obtained by ATLAS and CMS for existing searches in similar final states. For the signal, the main uncertainty is generated by the impact on the selection efficiency of uncertainties on the measurement of quantities such as e.g. the energy scale and resolution for jets and E miss T . In the case of the background there is in addition an important contribution to the total uncertainty that is associated with the procedure used to obtain the background estimate, which is typically achieved through a mixture of MC and data-driven techniques. Since we perform a shape analysis, the obtained fit results have reduced sensitivity to the absolute normalisation uncertainties, and are essentially determined by the uncertainties on the prediction of the shape of the distribution of the fitted variable for the SM background. The magnitude of these uncertainties is difficult to forecast, as they include different factors, such as the shape distortion from uncertainties on energy and efficiency determinations, or theoretical uncertainties associated to the simulation of the background. A variety of techniques are used by the experiments to control shape uncertainties, including the usage of appropriate control regions and the profiling of experimental uncertainties. In the case of the ttZ final state, shape uncertainties of a few percent seem to be an achievable goal, and we will determine the LHC reach, assuming a representative value of 5% for the latter uncertainty.
The results given in the following are for integrated luminosities of 300 fb −1 and 3 ab −1 , corresponding to the LHC Run-3 phase and the high-luminosity option of the LHC (HL-LHC), respectively. As the LHC experimental community is still working on the detailed assessment of the impact of the high pileup on the detector performance in the HL-LHC phase, we assume for simplicity the same detector performance for the two benchmark luminosities. The analysis based on the ∆m variable, relying on an accurate measurement of E miss T and the momenta of jets will likely be affected by pileup. In contrast, we expect only a minor impact on the p T,Z variable, built from two high-p T leptons. Under the assumptions presented above, we find that a shape analysis using ∆m leads to only marginally better results than a fit to p T,Z . We are therefore convinced that the conclusions of this study are valid also in the presence of a much higher pileup than the one experienced in the ongoing LHC run. In this section, we will only show the results of our ∆m shape fit. A comparison of the performance of the ∆m and p T,Z fits is provided in Appendix A.
The results of our sensitivity study are displayed in Figures 9 and 10. In the two panels of the first figure, we show the 95% CL exclusion limits in the M H -M A plane that derive from a shape fit to the ∆m observable introduced in (3.4). The red (yellow) contours illustrate the constraints that follow from 300 fb  Although our ttZ analysis is performed keeping effects due to off-shell H/A production and decay, and due to the interference with the SM background (see Section 3), it ignores possible modifications of the H/A line shape [85][86][87]. The latter effects have been studied in [88,89], where it was found that for a heavy Higgs boson different treatments of its propagator can lead to notable changes in the inclusive production cross sections compared to the case of a Breit-Wigner with a fixed width, as used in our work. In consequence, the exclusion limits in the upper left and lower right corner of the right plot in Figure 9 carry some (hard to quantify) model dependence related to the precise treatment of the H/A propagators.
In Figure 10 we show furthermore the 95% CL exclusion contours in the M A -t β plane for the type-II 2HDM scenarios with M H = M H ± = M A + 200 GeV and λ 3 = 6. In both panels the results of our ∆m shape fit are given for four different values of c β−α . Notice that for the chosen parameters there are no issues with perturbativity and vacuum stability, and that the A is sufficiently narrow for the NWA to hold. From the results shown on the left-hand side one observes that with 300 fb −1 of √ s = 14 TeV data, all values of t β 2.5 can be excluded for A masses close to the top threshold in the exact alignment limit, i.e. c β−α = 0. If the Higgs sector is not perfectly aligned, the branching ratio H → AZ is reduced see (2.2) , and as a result the bounds in the M A -t β plane become weaker. For instance, for the choice c β−α = 0.15, we find that the reach is decreased by roughly a factor of 1.5 compared to the case of full alignment. One also sees that M A values up to around 450 GeV can be excluded with 300 fb −1 for c β−α ≤ 0.15 and t β = O(1). With 3 ab −1 of integrated luminosity it turns out that the obtained limits can be pushed to values of M A and t β that are higher by approximately 30% than the corresponding 300 fb −1 bounds. This statement is illustrated by the coloured contours that are displayed in the right panel of Figure 10. The discovery reach corresponding to the latter figure are provided in Appendix B.
The constraints on the type-II 2HDM parameter space presented in this section should be compared to the bounds that have been derived in [32][33][34]. These analyses have considered the pp → H/A → tt, pp → ttH/A → tttt and pp → bbH/A → bbtt channels, and there seems to be a consensus that future searches for pp → ttH/A → tttt should provide the best sensitivity to neutral Higgses with masses M H,A > 2m t at low values of t β . While a one-to-one comparison with the exclusions obtained in [32][33][34] is not possible, we note that the limits derived in our work appear to be more stringent than the bounds reported in the latter articles. In this context it is important to realise that the reach of the pp → ttH/A → tttt searches has found to be strongly dependent on the systematic uncertainty of the normalisation of the tt background. The ttZ analysis strategy proposed by us does in contrast not rely on knowing the absolute size of the relevant backgrounds to the level of a few percent, since the search gains its discriminating power from shape differences. We therefore expect future ttZ searches to lead to the most robust coverage of the 2HDM parameter space with M H,A > 350 GeV, |M H − M A | > M Z and t β = O(1).

Conclusions
In this article, we have proposed to use the ttZ and In fact, the involved couplings g HAZ , g Htt , g Att , g HH ± W ∓ , g AH ± W ∓ and g H ± tb are all non-vanishing for c β−α 1 see (2.2) and (2.3) which corresponds to the so-called alignment limit that is experimentally favoured by the agreement of the LHC Higgs measurements with SM predictions. As a result, appreciable ttZ and tbW rates associated to H/A production turn out to be a rather generic prediction in 2HDMs that feature a SM-like 125 GeV scalar and non-SM Higges that are heavier than about 350 GeV with some of their masses split by around 100 GeV or more.
By analysing the anatomy of the ttZ and tbW signatures in 2HDMs, we have demonstrated that many of the resulting final-state distributions show peaks and/or edges that are characteristic for the on-shell production of a resonance followed by its sequential decay into visible and invisible particles. These kinematic features can be used to disentangle the new-physics signal from the SM background. In the case of the ttZ final state, we found that the difference ∆m ≡ m ttZ − m tt between the masses of the ttZ and tt systems and the transverse momentum p T,Z of the Z boson are powerful discriminants, while in the case of the tbW final state the invariant masses m bb , m bb , m tb and m tbW can be exploited for a signal-background separation. Our MC simulations have furthermore shown that the discussed observables can all be reconstructed and well measured under realistic experimental conditions through either a dedicated three-lepton (∆m and p T,Z ) or twolepton (m bb and m bb ) analysis strategy (see Sections 6 and 7).
Applying our three-lepton analysis strategy to simulated 14 TeV LHC data, we have then presented a comprehensive sensitivity study of the ttZ signature in the 2HDM framework. We have derived various 95% CL exclusion limits on the parameter space of the type-II 2HDM that follow from a shape fit to the ∆m (see Section 8) Figure 9). For the scenarios M H = M H ± = M A + 200 GeV, we have also derived the 95% CL exclusion limits in the M A -t β plane for four different values of c β−α (see Figure 10). For the choice c β−α = 0.1, we found for instance that it should be possible to exclude t β values up to almost 2 for M A = 350 GeV, assuming 300 fb −1 of data. The HL-LHC is expected to improve the quoted LHC Run-3 limits noticeably. The 5σ discovery reach in the M H -M A and M A -t β plane can be found in Appendix B. The constraints obtained in our work are complementary to and in many cases stronger than the exclusions that future LHC searches for the processes pp → H/A → tt, pp → ttH/A → tttt and pp → bbH/A → bbtt are expected to be able to provide (cf. [32][33][34]) on 2HDMs with neutral Higgses with M H,A > 2m t and small values of t β .
In the case of the tbW final state, we have found that for both the two-lepton and one-lepton analysis the signal-over-background ratios do not exceed the level of a few percent. As a result, a reliable evaluation of the coverage of the 2HDM parameter space would require to make strong assumptions about the systematic uncertainties that plague the normalisation and shape of the tt background at future LHC searches. Since we feel that it would be premature to make these assumptions, we hope that the ATLAS and CMS collaborations will explore the tbW signature further. This is a worthwhile exercise, because we expect that at the HL-LHC with 3 ab −1 of integrated luminosity, this channel should also allow to probe parts of the 2HDM parameter space that feature heavy non-SM Higgses and t β values of the order of a few. In fact, the maximum sensitivity of our tbW analysis arises for M A = M H ± 400 GeV, making the tbW coverage complementary to that of the proposed ttZ search.