Search for doubly and singly charged Higgs bosons decaying into vector bosons in multi-lepton final states with the ATLAS detector using proton-proton collisions at $\sqrt{s}$ = 13 TeV

A search for charged Higgs bosons decaying into $W^\pm W^\pm$ or $W^\pm Z$ bosons is performed, involving experimental signatures with two leptons of the same charge, or three or four leptons with a variety of charge combinations, missing transverse momentum and jets. A data sample of proton-proton collisions at a centre-of-mass energy of 13 TeV recorded with the ATLAS detector at the Large Hadron Collider between 2015 and 2018 is used. The data correspond to a total integrated luminosity of 139 fb$^{-1}$. The search is guided by a type-II seesaw model that extends the scalar sector of the Standard Model with a scalar triplet, leading to a phenomenology that includes doubly and singly charged Higgs bosons. Two scenarios are explored, corresponding to the pair production of doubly charged $H^{\pm\pm}$ bosons, or the associated production of a doubly charged $H^{\pm\pm}$ boson and a singly charged $H^\pm$ boson. No significant deviations from the Standard Model predictions are observed. $H^{\pm\pm}$ bosons are excluded at 95% confidence level up to 350 GeV and 230 GeV for the pair and associated production modes, respectively.


Introduction
Experimental signatures with two leptons with the same electric charge (same-charge) or multi-lepton final states are extensively exploited in searches for physics beyond the Standard Model (BSM physics) at hadron colliders. In many models, heavy BSM particles may be produced in proton-proton collisions and decay into multiple massive Standard Model (SM) electroweak gauge bosons or top quarks. Subsequent decays of these gauge bosons into final states with leptons can occur with considerable branching ratios. These final states are favourable to the search for new phenomena since the yields predicted within the SM are generally low and the experimental effects are well understood. At the Large Hadron Collider (LHC) [1], signatures with two leptons of the same charge, or three or four leptons with a variety of charge combinations have been used by the ATLAS and CMS experiments to explore the landscape of possible SM extensions and their phenomenology [2][3][4][5][6][7][8][9][10][11][12][13][14][15][16][17][18][19][20][21]. Among these proposed extensions, the addition to the SM sector of one weak gauge triplet of scalar fields with a non-zero vacuum expectation value of the neutral component is a compelling way to account for neutrino masses through the type-II seesaw mechanism [22].
Guided by such a model, this paper presents a search for doubly and singly charged Higgs bosons. In the considered type-II seesaw model, after the electroweak symmetry breaking, there are seven scalar bosons: ±± , ± , 0 (CP odd), 0 (CP even) and ℎ 0 (CP even). The vacuum expectation value is assumed to be 100 MeV, the lightest CP-even Higgs boson is the discovered SM Higgs boson, and the chosen model parameters are compatible with the measured SM Higgs boson properties and electroweak precision data. The theoretical and experimental constraints help in choosing allowed values for the parameters of the model from which the predominant production and decay modes follow. Only the pair production of doubly charged ±± bosons via a virtual photon or a boson, and the associated production of a doubly charged ±± boson and a singly charged ± boson via a ± boson are considered. For the ±± pair production mode the ±± boson mass is at least 100 GeV lower than the ± boson mass, while for the associated production of ±± and ± bosons the mass difference between the doubly and singly charged Higgs bosons is at most 5 GeV. Signal processes were simulated for different ±± boson mass hypotheses between 200 and 600 GeV. The ±± boson is considered to decay into a ± ± pair and the ± boson into a ± pair. The four gauge-boson final states are identified by the presence of charged leptons (electrons or muons), missing transverse momentum and hadronic jets. The analysis uses three final states (channels) classified according to the number of leptons: two same-charge (SC) leptons (2ℓ sc ), three leptons (3ℓ), and four leptons (4ℓ). Electrons and muons from -lepton decays are also considered.
An earlier ATLAS analysis using this signature to target the pair production of doubly charged ±± bosons at √ = 13 TeV with an integrated luminosity of 36 fb −1 placed lower limits of between 200 and 220 GeV on the ±± boson mass [23]. In this paper the search is extended and refined, using a larger data set of collisions at √ = 13 TeV recorded by the ATLAS detector corresponding to an integrated luminosity of 139 fb −1 . In addition, the production of a ±± boson in association with a ± boson is included, extending the initial search and allowing a search for doubly and singly charged Higgs bosons in a previously unexplored region of the considered model's parameter space. Section 2 presents the theoretical framework while Section 3 gives an overview of the experimental set-up. The different Monte Carlo simulations used in the analysis are described in Section 4, followed by a description of the event reconstruction in Section 5. Signal regions are built from a set of requirements on the kinematic properties of the selected events, in order to isolate the signature of BSM processes from SM backgrounds, as presented in Section 6. Section 7 describes the estimation of the backgrounds. The systematic uncertainties are detailed in Section 8. Event yields in data are then compared with the estimated contributions from SM processes, and the results are discussed in Section 9. Conclusions are presented in Section 10.  Figure 1: Feynman diagrams illustrating production and subsequent decay of the ±± and ± bosons: (a) ±± pair production and (b) ±± and ± associated production.

Signal model
A type-II seesaw model, in which the scalar sector of the SM is extended with a hypercharge = 2 scalar triplet Δ [24][25][26][27][28][29], is studied. The electroweak symmetry breaking (EWSB) is achieved by requiring the neutral components of the SM scalar doublet and scalar triplet Δ to acquire vacuum expectation values, d and t respectively. After the EWSB, the mixing between these fields results in seven scalar bosons: ±± , ± , 0 (CP odd), 0 (CP even) and ℎ 0 (CP even). The ℎ 0 boson is the SM Higgs boson, and a small mixing between the CP-even scalars is naturally obtained in the model. In addition, the triplet-neutrino Yukawa term provides non-zero neutrino masses proportional to the vacuum expectation value of the triplet t . In this type-II seesaw model, constraints from electroweak precision measurements lead to an upper bound on t of approximately 2.5 GeV, which is significantly lower than the electroweak scale and matches the need for small values suggested by the natural association of t with the neutrino masses. Additional, theoretical constraints are the absence of tachyonic modes, requirements for the vacuum structure and stability of the potential, and unitarity requirements. A detailed discussion can be found in Ref. [29]. The theoretical and experimental constraints help to choose allowed values for the parameters of the model from which the predominant production and decay modes follow. Some type-II seesaw models also extend the weak sector by introducing two scalar triplets and right-handed gauge bosons R (left-right symmetric models). The model studied in this analysis corresponds to the left-handed triplet discussed in Refs. [30][31][32][33][34].
Two scenarios with t equal to 100 MeV are explored in this paper. In the first scenario, the mass difference between ±± and ± bosons is > 100 GeV, the ± boson being the heavier one. Large mass splitting could still be allowed despite the | ± − ±± | < 40 GeV constraint deduced from the electroweak precision data [35] if radiative corrections reduce the tree-level contribution coming from the triplet vacuum expectation value. Only the pair production of ±± bosons via the diagram shown in Figure 1(a) ( → * / * → ±± ∓∓ ) is considered. The triplet vacuum expectation value is taken to be t = 100 MeV such that only the doubly charged ±± boson decay into a pair of bosons with the same charge, ±± → ± ± , is possible [28]. The leptonic decays ±± → ℓ ± ℓ ± are suppressed with increasing t [36,37] and not considered in this paper.
In the second scenario, the mass of the ± boson is chosen to be at most 5 GeV different from the ±± boson mass. The low mass difference ensures that the charged Higgs bosons do not decay into each other with a non-negligible branching fraction, and complies with the additional | ± − ±± | < 40 GeV constraint deduced from the electroweak precision data. In addition, the chosen mass difference between the doubly and singly charged Higgs bosons maximises the signal amplitudes. Only the associated production of ±± and ± bosons via the diagram shown in Figure 1(b) ( → ± * → ±± ∓ ) is considered. The production cross sections for ±± pair production with the same mass settings as for associated production of ±± and ± bosons can be large. However, this production mode is exactly the one described in the first scenario and it is not considered in the second scenario. This choice is motivated by the objective of the search, which is to study a characteristic signature with a benchmark parameter point and report its cross section. The ±± boson decays into a pair of bosons with the same charge, with a branching ratio of 100%. Only the bosonic decays of the singly charged bosons ( ± → ± ) are considered and, depending on the ± boson mass, the branching ratio varies between 40% and 60% (see Table 1). The branching fraction for the ±± = 300 GeV mass hypothesis is lower than for the neighbouring points because of the high dependency of this quantity on the mass of the scalar triplet Δ [29]. Depending on the ± boson mass, other decay modes are ± →¯and ± → ± ℎ 0 [28]. Studies at Monte Carlo generator level show that after the selection of at least two same-charge leptons or at least three leptons and no -jet (where the -jets are selected with 70% efficiency), the contribution from the other possible ± decays is negligible. Similar conclusions are also reached after examining the results obtained for the various control regions used in the analysis. The effect on the mass below which a charged Higgs boson is excluded is negligible.
Pair production of ± bosons is also possible, albeit with a much smaller production cross section than for ±± pair production. Therefore, ± pair production is not considered in this paper. The cross section for single ±± production via vector-boson fusion ( → ± * ± * → ±± ) is proportional to t , and hence negligible.
For the ±± pair production mode, the mixing between the CP-even scalars is taken to be 10 −4 and the remaining five couplings in the potential are adjusted to obtain a given ±± boson mass hypothesis while requiring ℎ 0 to have a mass of 125 GeV. Similar settings are also used for the associated production mode. The branching fraction times cross-section calculation for the pair production of ±± bosons and the associated production of ±± and ± bosons is performed for on-shell and bosons, and therefore only the region ±± > 200 GeV is considered in the present analysis.
Extensive searches for leptonic decays have been performed at various experiments [38][39][40][41][42][43], excluding doubly charged ±± bosons with masses up to about 870 GeV. The CMS Collaboration published results for the ±± → ± ± decay mode in the context of single ±± production through vector-boson fusion with detector (ID) surrounded by a thin superconducting solenoid providing a 2 T axial magnetic field, electromagnetic (EM) and hadron calorimeters, and a muon spectrometer (MS). The inner tracking detector covers the pseudorapidity range | | < 2.5. It consists of silicon pixel, silicon microstrip and transition radiation tracking detectors; the innermost layer is 33 mm from the beamline [48][49][50]. Lead/liquid-argon (LAr) sampling calorimeters provide electromagnetic energy measurements with high granularity. A steel/scintillator-tile hadron calorimeter covers the central pseudorapidity range (| | < 1.7). The endcap and forward regions are instrumented with LAr calorimeters for EM and hadronic energy measurements up to | | = 4.9. The muon spectrometer surrounds the calorimeters and is based on three large air-core toroidal superconducting magnets with eight coils each. The field integral of the toroids ranges between 2.0 and 6.0 T·m across most of the detector. The muon spectrometer includes a system of precision tracking chambers and fast detectors for triggering. A two-level trigger system [51] is used to select events. The first-level trigger is implemented in hardware and uses a subset of the detector information to keep the accepted rate below 100 kHz. This is followed by a software-based trigger that reduces the accepted event rate to 1 kHz on average depending on the data-taking conditions.

Monte Carlo event simulation
Monte Carlo (MC) event generators were used to simulate the signal and background events produced in the proton-proton collisions. For the ±± pair production and the ±± and ± associated production signal processes, the events at particle level were generated with M G 5_aMC@NLO [53] using leading-order (LO) matrix elements (ME) and NNPDF3.0lo parton distribution functions (PDF) [54]. The events were subsequently showered using P 8 [55] with a set of tuned parameters called the A14 tune [56]. Signal processes were simulated for different mass hypotheses between 200 and 600 GeV, and presented in Table 1. For the pair and associated production, the cross section decreases rapidly with the charged Higgs boson's mass as shown in Table 1. The event samples were normalised using calculations at relativistic limit. Angular distance is measured in units of Δ ≡ √︁ (Δ ) 2 + (Δ ) 2 . The magnitude of the momentum in the plane transverse to the beam axis is denoted by T . Table 2: Summary of the event generators, parton shower models, and PDF sets used for the simulation of the background event samples. The notation is used to refer to an electroweak gauge boson or / * . In the final column, 'default' refers the to default parameter set provided with the event generator.
Process Generator ME accuracy PDF Parton shower and hadronisation Parameter set next-to-leading order (NLO) in perturbative QCD [52]. The NLO -factor increases the expected event yields by a factor of 1.25, independently of the mass of the charged Higgs bosons.
A summary of the MC event generators used to simulate the background is presented in Table 2 and further details are given below. The notation is used to refer to an electroweak gauge boson or / * .  [57] generator depending on the process, including off-shell effects and Higgs boson contributions, where appropriate. The and processes were simulated using matrix elements at NLO accuracy in QCD for up to one additional parton and at LO accuracy for up to three additional parton emissions. The electroweak production of dibosons in association with two jets ( -EW jj) was generated using LO accuracy matrix elements. The matrix element calculations were matched and merged with the S parton shower based on Catani-Seymour dipole factorisation [58,59] using the MEPS@NLO prescription [60][61][62][63]. The virtual QCD corrections were provided by the O L library [64,65]. The production of triboson ( ) events was simulated with the S v2.2.2 generator using factorised gauge-boson decays. All diboson and triboson processes were generated using the NNPDF3.0nnlo PDF set [54], along with the dedicated set of tuned parton-shower parameters developed by the S authors. The cross sections from the event generator were used for the normalisation.
The production of an electroweak gauge boson or virtual photon in association with jets ( +jets) was simulated with the S v2.2.1 generator using NLO matrix elements for up to two partons, and LO matrix elements for up to four partons calculated with the Comix [58] and O L libraries. The event samples were generated using the NNPDF3.0nnlo PDF set, along with the dedicated set of tuned parton-shower parameters developed by the S authors. The +jets event samples were normalised to a NNLO prediction.
Higgs boson production in association with a vector boson was simulated at LO with P 8.186 and E G [66] using the A14 tune, along with the NNPDF2.3lo PDF set. The Monte Carlo prediction was normalised to cross sections calculated at NNLO in QCD with NLO electroweak corrections for¯/ → , and at next-to-leading-logarithm accuracy in QCD with NLO EW corrections for → [67][68][69][70][71][72][73]. The production of¯events was modelled using the P B v2 [74][75][76][77][78] generator at NLO with the NNPDF3.0nlo PDF set [54]. The events were interfaced to P 8.230 using the A14 tune and the NNPDF2.3lo PDF set [54]. The decays of bottom and charm hadrons were performed by E G v1.6.0. The production of¯, and events was modelled using the M G 5_aMC@NLO v2.3.3 generator at NLO with the NNPDF3.0nlo PDF. The events were interfaced to P 8.210 (P 8.212 for ) using the A14 tune and the NNPDF2.3lo PDF set. The decays of bottom and charm hadrons were simulated using the E G v1.2.0 program. The simulated ,¯, and events were normalised to cross-section calculations at NLO accuracy in QCD. In addition, NLO EW corrections were included for¯,¯and .
Backgrounds from top-quark pair production (¯) and production were estimated at NLO accuracy in QCD using the hvq program [74] in P B v2. The event samples were generated using the A14 tune and NNPDF3.0nnlo PDF set. The interference between¯and production is neglected as it has a negligible impact on the analysis. The¯event sample was normalised to the cross-section prediction at NNLO in QCD including the resummation of next-to-next-to-leading logarithmic (NNLL) soft-gluon terms calculated using T ++2.0 [79][80][81][82][83][84][85]. The inclusive cross section was corrected to the theory prediction calculated at NLO in QCD with NNLL soft gluon corrections [86,87].
The production of¯¯,¯,¯and¯events was modelled using the M G 5_aMC@NLO v2.3.3 generator at NLO with the NNPDF3.1nlo PDF [54]. The events were interfaced with P 8.230 using the A14 tune and the NNPDF2.3lo PDF set. The¯¯,¯and¯contributions were normalised to theoretical cross sections at NLO in QCD including electroweak corrections [88]. The¯production was normalised to a cross section calculated to LO in QCD.
Backgrounds such as → + + − − or → ± ± ∓ with 2ℓ sc , 3ℓ or 4ℓ final states, which have the same signature as the considered signal, have very small production cross sections times branching ratio and their contribution in the preselection, control and signal regions was found to be negligible.
The signal and background events were passed through the G 44 [89] simulation of the ATLAS detector [90] and reconstructed using the same algorithms as are used for the data. The effect of multiple proton-proton interactions in the same or nearby bunch crossings (pile-up) is accounted for using inelastic proton-proton interactions generated by P 8 [91], with the A3 tune [92] and the NNPDF2.3lo PDF set [93]. These inelastic proton-proton interactions were added to the signal and background event samples and weighted such that the distribution of the average number of proton-proton interactions in simulation matches that observed in the data.

Event reconstruction
The analysis is performed in collision data recorded by the ATLAS detector between 2015 and 2018. In this period, the LHC delivered colliding beams with a peak instantaneous luminosity up to = 2.1 × 10 34 cm −2 s −1 , achieved in 2018, and an average number of interactions per bunch crossing of 34. After requirements on the stability of the beams, the operational status of all ATLAS detector components, and the quality of the recorded data, the total integrated luminosity of the data set corresponds to 139 fb −1 [94].
Proton-proton interaction vertices are reconstructed from charged-particle tracks with T > 500 MeV [95,96] in the ID. The presence of at least one such vertex with a minimum of two associated tracks is required, and the vertex with the largest sum of 2 T of associated tracks is chosen as the primary vertex. The anti-algorithm [97,98] with radius parameter = 0.4 is used to reconstruct jets up to | | = 4.9. It uses as inputs particle-flow objects, combining tracking and calorimetric information [99]. The jets are calibrated as described in Ref. [100]. Only jets with T > 20 GeV and | | < 2.5 are considered further. Events are removed if they contain jets induced by calorimeter noise or non-collision background, according to criteria similar to those described in Ref. [101]. Jets with T < 60 GeV and | | < 2.4 from pile-up interactions are suppressed using a jet-vertex tagging discriminant [102].
Jets containing -flavored hadrons are identified in the region | | < 2.5 with the DL1r -tagging algorithm based on a recurrent neural network [103,104]. It makes use of the impact parameters of tracks associated with the jet, the position of reconstructed secondary vertices and their compatibility with the decay chains of such hadrons. The -tagging average efficiency is 70%, as measured in¯events.
Electron candidates are reconstructed as tracks in the ID matched to energy clusters in the EM calorimeter, within | | < 2.47 [105]. Only electrons with T > 10 GeV and not in the transition region between the barrel and endcap calorimeters (1.37 < | | < 1.52) are considered. The electron identification is based on a multivariate likelihood-based discriminant that uses the shower shapes in the EM calorimeter and the associated track properties measured in the ID. The electron candidates must satisfy the 'Loose' identification criteria described in Ref. [105]. Signal electrons are required to satisfy the 'Tight' identification [105] for better rejection of non-prompt electrons. The electron identification efficiency varies with increasing T in → events, from 65% at T = 10 GeV to 88% at T = 100 GeV for the Tight operating point, and from 85% at T = 20 GeV to 95% at T = 100 GeV for the Loose operating point.
The longitudinal impact parameter of the electron track, 0 , is required to satisfy | 0 sin | < 0.5 mm, where is the polar angle of the track. The transverse impact parameter divided by its uncertainty, | 0 |/ ( 0 ), is required to be less than five. For all signal electrons there must be no association with a vertex from a reconstructed photon conversion [105] in the detector material. To further reduce the photon conversion background in the 2ℓ sc channel, additional requirements are applied to the signal electrons [4]: i) the candidate must not have a reconstructed displaced vertex with radius > 20 mm whose reconstruction uses the track associated with the electron, ii) the invariant mass of the system formed by the track associated with the electron and the closest track at the primary or a conversion vertex is required to be larger than 100 MeV; this selection is referred to as the photon conversion veto.
For the signal electrons, the identification criteria are complemented by an isolation requirement, which is based on the energy in a cone around the electron candidate calculated using either reconstructed tracks or energy clusters. Electrons with wrongly reconstructed charge (charge-flip) are suppressed using a boosted decision tree discriminant exploiting additional tracks in the vicinity of the electron and track-to-cluster matching variables [105].
Muon candidates are reconstructed in the region | | < 2.5 from MS tracks matching ID tracks. The analysis only considers muons with T > 10 GeV satisfying the 'Medium' quality requirements defined in Ref.
[106]. The muon reconstruction efficiency is approximately 98% in → events. The longitudinal impact parameter of the muon track must satisfy | 0 sin | < 0.5 mm and the transverse impact parameter must satisfy | 0 |/ ( 0 ) < 3. For signal muons the candidate must satisfy calorimeter-and track-based isolation requirements.
Non-prompt electrons and muons from the decays of -and -flavored hadrons are further rejected using a boosted decision tree discriminant based on isolation and secondary vertex information, referred to as the non-prompt-lepton veto [107].
To avoid cases where the detector response to a single physical object is reconstructed as two different final-state objects, e.g. an electron reconstructed as both an electron and a jet, several steps are followed to remove such overlaps, as described in Ref. [13]. The overlap removal procedure is performed using candidate leptons. Table 3: The requirements applied to define the categories of candidate leptons: loose (L), loose and minimally isolated (L * ) and tight (T) leptons. The overlap removal procedure is not applied for the candidate leptons. In the 2ℓ sc channel, the photon conversion veto is required in addition to the loose and tight criteria.
The different lepton selections used in the analysis are summarised in Table 3 The missing transverse momentum, with magnitude miss T , is defined as the negative vector sum of the transverse momenta of all identified objects (muon candidates, electron candidates and jets) and an additional soft term [108,109]. The soft term is constructed from all tracks that are matched to the primary vertex and are not associated with any other objects.

Event selection
Candidate events are selected for read-out using lepton triggers that require the one electron or muon to satisfy identification criteria similar to those used in the offline reconstruction, isolation criteria, and a transverse momentum requirement of T > 26 GeV [110,111]. With increasing T the requirements on identification and isolation become less stringent. All events must contain at least one offline tight lepton with ℓ T > 30 GeV that triggered the event. The event selection proceeds in two steps: the preselection and the signal regions (SRs) selection.
The preselection is defined in Table 4; the three channels (2ℓ sc , 3ℓ and 4ℓ) are defined to be mutually exclusive. Events are selected only if the absolute value of the sum of charges of the leptons is two, one, and two or zero for the 2ℓ sc , 3ℓ and 4ℓ channels, respectively. In the 2ℓ sc channel, the second-highest-T lepton is required to have T > 20 GeV and both leptons are required to be tight. Similarly in the 3ℓ channel, each lepton in the pair of leptons of the same charge is required to have T > 20 GeV and be tight.
Further preselection requirements are based on miss T , the jet multiplicity jets and the number of -jets, -jet . The lowest miss T value is 30 GeV in the 3ℓ and 4ℓ channels, and 70 GeV in the 2ℓ sc channel; this selection helps to reduce the non-prompt lepton, electron charge-flip and backgrounds. In the 2ℓ sc (3ℓ) channel only events with at least three (two) jets are considered. The background from top-quark production is highly reduced by requiring zero -jets in the event. In order to reduce the background from Drell-Yan processes and neutral mesons, the invariant mass of same-flavour opposite-charge lepton pairs is required to Table 4: The preselection criteria for the 2ℓ sc , 3ℓ and 4ℓ analysis channels. The leptons are required to pass the loose (L), loose and minimally isolated (L * ) or tight (T) requirements. The leptons are ordered by decreasing T (ℓ 1 , ℓ 2 , ...) in the 2ℓ sc and 4ℓ channels, while for the 3ℓ channel ℓ 1 and ℓ 2 denote the two same-charge leptons and ℓ 0 denotes the lepton with a charge opposite to the total lepton charge. ℓ denotes the charge of each lepton. SFOC refers to same-flavour opposite-charge lepton pairs. The symbol "-" means no requirement is made. The equal sign (=) is used to emphasise that the selection criterion has to be exactly the given number.
Selection criteria 2ℓ sc 3ℓ 4ℓ At least one offline tight lepton with ℓ T > 30 GeV that triggered the event be greater than 15 GeV for the 3ℓ and 4ℓ channels. The invariant mass of the same-flavour opposite-charge lepton pair must differ from the nominal -boson mass by 10 GeV. For the 2ℓ sc channel, the -boson invariant mass veto is also applied to ± ± events, in order to reduce the contributions originating from electron charge misidentification.
In addition to miss T the following variables are used to define SRs: • The invariant mass of all selected leptons in the event, ℓ , where can be 2, 3 or 4 corresponding to the 2ℓ sc , 3ℓ or 4ℓ channels.
• The invariant mass of all jets in the event, jets . When there are more than four jets in the event, only the leading four jets are used. This variable is only used for the 2ℓ sc channel.
• The distance in -between two same-charge leptons, Δ ℓ ± ℓ ± . It is used for the 2ℓ sc and 3ℓ channels.
In the 4ℓ channel, two such variables can be calculated per event, Δ min ℓ ± ℓ ± and Δ max ℓ ± ℓ ± , denoting the minimum and maximum values, respectively.
• The transverse momentum of the highest-T jet, leading jet T . This variable is used in the 3ℓ channel.
• The transverse momentum of the highest-T lepton, ℓ 1 T . This variable is used in the 4ℓ channel. • The azimuthal distance between the dilepton system and miss T , Δ ℓℓ, miss T . It is only used in the 2ℓ sc channel.
• The smallest distance in -between any lepton and its closest jet, Δ ℓjet . This variable is used in the 3ℓ channel. Table 5: Definition of the signal regions optimised for the study of different ±± and ± mass hypotheses. The selection is applied on top of the preselection defined in Table 4. For the ±± pair production mode, the ±± = 300 GeV signal regions are also used for ±± = 350 GeV. For the ±± and ± associated production mode, the ±± = 200 GeV, 400 GeV and 500 GeV signal regions are also used for ±± = 220 GeV, 450 GeV and 550 GeV, respectively. The variables used are described in Section 6.  • The variable , used for the 2ℓ sc channel to describe the event topology in the transverse plane, and defined using the spread of the angles of the leptons, miss T , and jets as follows:

Charged Higgs
where the R is the root mean square that quantifies the spread, These variables are found to discriminate between the signal and the background. They exploit both the boosted decay topology of the charged Higgs bosons and the high energy of the decay products.
Signal regions are defined for each channel, as summarised in Table 5. The selection criteria defining the signal regions result from a scan of the multidimensional parameter space of the discriminating variables mentioned above. The SRs were designed by optimising the sensitivity for the ±± pair production mode, using the ±± = 200, 300, 400 and 500 GeV mass hypotheses. The same SRs are used to study the ± associated production mode; even though this approach is not optimal, it is still preferred instead of significantly increasing the number of signal regions. Additional mass hypotheses are used in the signal region defined for the lower mass hypotheses, since the signal discrimination power does not vary significantly in this regime. In particular, for the ±± pair production mode the ±± = 300 GeV signal regions are also used for ±± = 350 GeV. For the ±± and ± associated production mode, the ±± = 200 GeV, 400 GeV and 500 GeV signal regions are also used for ±± = 220 GeV, 450 GeV and 550 GeV, respectively. The SRs defined for the 2ℓ sc channel are further divided into , and final states. Events in the 3ℓ SRs are separated into two categories according to whether or not a same-flavour opposite-charge lepton pair exists in the event. This procedure further improves the expected significance by exploring differences in background composition and lepton-flavour composition between signal and background. The number of events observed in data is shown together with the expected signal and estimated background yields in Section 9.

Background estimation
The background sources can be divided into two main categories. One category is populated by the SM events which contain only reconstructed charged (prompt) leptons originating from leptonic decays of and bosons. The second category is formed by non-prompt leptons and charge-flip electrons. The non-prompt-lepton category refers to the leptons that originate from decays of -and -hadrons, or single pions that mimic electron signatures. The electrons from hadron decays into photons which convert into pairs of electrons in the beam pipe or detector material also enter this category. Lepton candidates reconstructed from these different sources share the properties of being less isolated, having larger impact parameters relative to the primary vertex and being less likely to satisfy the lepton identification criteria.
The backgrounds from SM processes with prompt leptons are estimated with MC simulations, except for background from production, for which the normalisation is corrected using data in a dedicated control region. Background from non-prompt leptons and electron charge-flip are estimated using data-based methods. Background from production can contribute in the SRs if the photon converts to an electronpositron pair, and is estimated using MC simulations. Background from production is estimated from simulation if the two bosons have the same electric charge and from data if the two bosons have opposite electric charge, since the latter only contributes through electron charge-flip and non-prompt leptons.

Background from production
The process is a dominant source of background in the 2ℓ sc and 3ℓ SRs. To correct a mismodelling seen in the jet multiplicity distribution [112], a normalisation factor is computed and then applied to the background events containing two or more jets. The normalisation factor is measured in a dedicated control region. It is selected by requiring exactly three tight leptons with T > 20 GeV, at least two jets and no -jet in the event. Finally, there must be at least one pair of same-flavour opposite-charge leptons with a mass compatible with the -boson mass (| oc ℓℓ − | < 10 GeV). The latter criterion ensures that the control region and the 3ℓ SRs are orthogonal.
The normalisation factor is derived in a fit of a first-order polynomial, as a function of jet multiplicity, to ratios of the data event yield (subtracting all non-contributions) to the event yield predicted in MC simulation. The value of the polynomial at jets = 0 is used to scale the predicted yields and found is to be 0.83 ± 0.07, where the uncertainty includes the statistical and systematic components. The different   sources of systematic uncertainty are discussed in Section 8. The jet multiplicity distribution in the control region is shown in Figure 2(a), and for illustration the distribution in the 3ℓ preselection region is also shown in Figure 2(b). The normalisation factor is applied to the MC contribution. The sum of estimated backgrounds agrees with data within the assigned systematic uncertainties.

Electron charge-flip background
In the 2ℓ sc channel, a background contribution is expected from events with opposite-charge lepton pairs where the charge of one of the leptons is misidentified. The charge-flip background is only significant for electrons and is mainly due to interactions of the electron with the ID material.
The misidentification rate is measured using a large data sample of dielectron events originating mainly from → + − decays selected by requiring two tight electrons with an invariant mass between 80 and 100 GeV. The sample contains mostly opposite-charge dielectron pairs, with a small fraction of same-charge dielectron pairs. The fraction of same-charge dielectron events is used to extract the charge misidentification rate as a function of the electron T and , using the method described in Ref. [105]. This rate is found to range between 0.01% and 4%, where higher values are obtained at large rapidities due to the larger amount of material traversed by the electrons. The statistical uncertainty of this estimate varies between 2% and 26% and is taken as a systematic uncertainty in the charge misidentification rate. The background in both the opposite-charge and same-charge samples is estimated in a cubic polynomial function fit of the high (100-120 GeV) and low (60-80 GeV) ℓℓ sidebands. The uncertainty in the background is estimated by varying the sidebands and propagates to a 3% uncertainty in the estimated charge misidentification rate. The final systematic uncertainty combines all the sources mentioned.
The charge-flip background in a given region is estimated by selecting data events with opposite-charge dilepton pairs, but otherwise identical selection, and weighting them by the probability that the charge of the electrons is misidentified. Another source of systematic uncertainty is estimated by comparing, in simulated +jets,¯, and events, the number of same-charge events estimated from opposite-charge events with the prediction; it accounts for differences between the charge misidentification rates in different processes and is found to be approximately 10%. For this test, the misidentification rate was measured using → + − MC simulations, using the same method as for the measurement performed with data.
To estimate the electron charge-flip background in the regions defined with loose electrons the methodology described above is applied. For this estimation, dedicated electrons misidentification rates measured with a → + − sample selected by requiring one tight and one loose electrons are used.

Non-prompt-lepton background
The composition of the non-prompt-lepton background in the SRs varies considerably among the analysis channels. Therefore, the methods to estimate these contributions are different for the 2ℓ sc , 3ℓ and 4ℓ channels. In the 2ℓ sc and 3ℓ channels the non-prompt-lepton background is estimated using a fake-factor method [23], while the simulation is corrected with scale factors measured in data for the 4ℓ channel.
Non-prompt-lepton background estimate for the 2ℓ sc channel The estimate of the background from non-prompt leptons assumes that these contributions can be extrapolated from a control region enriched in non-prompt leptons with a so-called fake factor. The control region is selected using the kinematic requirements of the preselection or the signal regions but alternative lepton selection criteria. The latter means that at least one of the selected leptons is required to satisfy the loose but not the tight lepton requirements.
The fake factors are calculated separately for electrons and muons in control regions with kinematic selections designed to enhance their content in non-prompt leptons. This is achieved by applying the preselection requirements of the 2ℓ sc channel, except that miss T must be lower than 70 GeV. Only events with electron and muon same-charge pairs are then used. The fake factor is defined as ratio of the number of events in the control region with all the selected leptons required to pass the tight signal requirements, to the number of events in the same region but with one of the selected leptons satisfying alternative lepton requirements. The measurement is performed as a function of the lepton T . The fake factor dependency on the lepton was also checked and found to be negligible. The SM processes with prompt leptons and the charge-flip contributions are subtracted in the control region. For the electron and muon fake-factor measurements, the lepton with the second highest T is assumed to be the non-prompt one.
The measured fake factors are 0.03 ± 0.01 for electrons and muons up to T = 40 GeV, and increase to 0.16 ± 0.05 and 0.09 ± 0.02, respectively, for electrons and muons with T > 60 GeV. The uncertainties are statistical only. A systematic uncertainty of 20% (10%) in the electron (muon) fake factor is estimated by studying the variation of the fake factor with miss T . For this measurement two miss T bins are considered, < 70 GeV and > 70 GeV. This uncertainty accounts for the different compositions of the non-prompt leptons in the control region and the SRs. Another source of systematic uncertainty accounts for how often the non-prompt lepton is actually the one with the highest lepton T and not the one with the second Table 6: The selection criteria that define the control regions enriched in non-prompt leptons used to determine the MC scale factors for the 4ℓ channel. The leptons are required to pass the loose and minimally isolated (L * ) requirements. The transverse mass, T , is calculated as the invariant mass of the vector sum of the transverse momentum of the non-prompt-lepton candidate and the missing transverse momentum. The symbol "-" means no requirement is made. The equal sign (=) is used to emphasise that the selection criterion has to be exactly the given number. Non-prompt-lepton background estimate for the 3ℓ channel The same method as the one employed for the 2ℓ sc channel is used in the 3ℓ channel. Here the opposite-charge lepton ℓ 0 passes the loose selection, and it is assumed to be prompt, an assumption that is found to be valid in MC simulation. The control region used to calculate the fake factors uses the 3ℓ preselection requirements, except that exactly one jet is required. Only events with electron and muon same-charge pairs are then used. The muon fake factor is found to be 0.03 ± 0.01 and the electron fake factor is found to be 0.02 ± 0.01, where the uncertainties are statistical only. A systematic uncertainty of 15% is estimated by measuring the lepton fake factor in a control region enriched in events from¯production. This uncertainty accounts for the different compositions of the non-prompt leptons in the control region and the SRs. Uncertainties in the subtraction of the SM processes with prompt leptons are found to be approximately 55% (45%) for the electron (muon) fake factor. Another 20% uncertainty comes from generator-level studies performed with MC simulations to test the assumption that ℓ 0 is a prompt lepton. The fake factors' dependency on the lepton T is also studied, and the deviations from the nominal fake factors are within the statistical uncertainty. When all sources of systematic uncertainty are combined, the total systematic uncertainty of the electron (muon) fake factor is found to be 60% (50%).

Non-prompt-lepton background estimate for the 4ℓ channel
There are too few data events to use the fake-factor method in the 4ℓ channel. Instead, the non-prompt-lepton yields are estimated from simulation corrected using scale factors that are extracted in control regions enriched in non-prompt leptons. The main sources of non-prompt leptons are -and -hadron decays. For electrons a small component from light-flavour hadrons is also present, while for muons this component is negligible. The scale factor is measured separately for non-prompt leptons from -hadron decays and light-flavour hadrons.
Two data samples enriched in non-prompt leptons originating from +jets and¯events are used to study how well the MC simulations describe the non-prompt leptons originating from light-flavour (LF) and heavy-flavour (HF) hadrons, respectively. These control regions are defined in Table 6. Three scale factors, HF , LF and HF , are obtained by solving the system of equations: Data | +jets − Prompt | +jets = HF HF | +jets + LF LF | +jets ,

Validation
The data-based methods employed to estimate the backgrounds are validated by comparing the event yields in data with the combined predictions for these backgrounds, added to MC predictions for SM processes with prompt signal leptons. Distributions of selected variables in the 2ℓ sc , 3ℓ and 4ℓ channels are shown in Figures 3, 4 and 5, respectively, after the preselection requirements from Table 4 are applied. Good agreement is observed in both normalisation and shape, demonstrating that the background contributions are well estimated. The expected distributions of both signal models are shown for ±± = 300 GeV to illustrate the discrimination power of the selected variables.
The signal contamination at the preselection level was studied for ±± = 300 GeV and ±± = 200 GeV mass hypotheses, corresponding to the ±± pair production mode and the ±± and ± associated production mode, respectively. These mass hypotheses were selected because they are close to the sensitive mass range of the current analysis. A maximum value of 8.5% (2%) ±± ∓∓ ( ±± ∓ ) signal contamination was found in all the individual 2ℓ sc , 3ℓ and 4ℓ channels at preselection level.

Systematic uncertainties
Uncertainties in the signal and background yields arise from experimental uncertainties and from the theoretical accuracy of the prediction of the SM background yields. The experimental uncertainties arise from the luminosity determination, modelling of pile-up interactions, the reconstruction and identification of electrons, muons and jets and from the uncertainties associated with the data-based methods that are used to estimate the non-prompt lepton and charge-flip electron backgrounds.    impact of these uncertainties on the background event yields estimated from MC simulation is found to be lower than 2%.
The uncertainties related to event reconstruction include the lepton [106,115] and jet [100]   from SM processes with prompt leptons and the signal yields in the SRs from these systematic uncertainties is found to be lower than 10%. The dominant contribution comes from the jet energy scale component.
The uncertainties in the efficiencies of the electron [116] and muon [106] reconstruction, identification and trigger are also included. Their impact on the estimated yields in the SRs for the signal and the background from SM processes with prompt leptons is lower than 5%. The uncertainties in the -jet identification are found to be negligible.
The uncertainties in the background yield from electron charge-flip are in the range of 20-30% for the 2ℓ sc and channels. The uncertainties in the non-prompt-lepton contributions in the SRs are between 25% and 90% and include the uncertainties described in Section 7.3, and statistical uncertainties in the control regions where the contributions are estimated.
For the background from , and multiboson processes, and , and¯processes (¯, ,¯), the systematic uncertainties due to higher-order QCD corrections are evaluated by varying the renormalisation and factorisation scales independently by factors of two and one-half, and removing combinations where the variations differ by a factor of four. The uncertainties due to the PDF and the s value used in the PDF determination are evaluated using the PDF4LHC prescription [117]. In the SRs defined for ±± = 300 GeV, the theory uncertainty in the background yields from and processes varies between 15% and 40%. For the background from it is approximately 10%, for the background from it is approximately 35% and for the and¯processes it is approximately 14%.
The total uncertainty in the estimated yields is 9% and includes the statistical uncertainty as well as two sources of theoretical uncertainty. Systematic uncertainties due to higher-order QCD corrections are evaluated using the same prescription as for the other diboson processes, and are found to be 3%. The second source is evaluated by comparing the results obtained with the linear fit model (Section 7.1) in the control region, and in a region defined with the control region requirements except that three jets must be in the event. This uncertainty is found to be 8.4%. To validate the assigned uncertainty, several checks are performed. The linear fit function is changed to a quadratic one to check the quality of the fit model. The fit parameterisation is studied in the -enriched region and in a region defined with the control region requirements except that three jets must be in the event. The differences between the obtained results are found to be covered by the total uncertainty. The choice of parton shower model is studied in the control region by comparing event samples simulated with S v2.2.1 and M G 5_aMC@NLO, and the differences are found to be covered by the total uncertainty. For the background in regions with lower jet multiplicities, the uncertainty is estimated with the same methodology as for the other diboson processes described above.
An uncertainty of 50% is assigned to the other backgrounds (¯,¯¯,¯and ), This large value is assigned to cover uncertainties from missing higher-order corrections and the PDF sets. Since these processes produce a larger number of jets at the first order of the perturbative expansion, they are less sensitive to parton shower modelling uncertainties.
The relative uncertainty in the background yields obtained in a fit of background to the observed data is shown in Figure 6 for all SRs. The statistical uncertainties originate from the limited number of preselected and opposite-charge data events used in the fake-factor method and the charge-flip electron background estimate, respectively, as well as the effect of the limited number of simulated events for SM processes with prompt leptons. The total uncertainty is computed using all sources of uncertainty, and they are treated as uncorrelated. The uncertainties range from 10% to 30% and are dominated by the statistical uncertainties in the non-prompt-lepton estimate and the theory uncertainties. An exception is the 2ℓ sc SR defined for ±± = 300 GeV, where the uncertainties from most sources are of similar size. The uncertainties associated with the charge-flip background are small in all 2ℓ sc SRs. In the 4ℓ channel, the statistical uncertainty of the non-prompt leptons comes from the limited number of events in the +jets,ā nd MC simulation samples. The theoretical uncertainties in the predicted signal yields arise from the parton shower model, missing higher-order corrections, and parton distribution functions. The systematic uncertainty due to the parton shower model is evaluated by comparing event samples simulated by P 8 using the A14 tune with samples from H using the H7--underlying-event tune [118], and is found to be less than 5%. The uncertainties due to the PDFs are found to be less than 5%. The uncertainty due to missing higher-order corrections is less than 10% [52]. When those uncertainties are combined in quadrature, an overall uncertainty in the signal yields of approximately 10% is obtained for the signal normalisation.

Results
The observed data event yields and the corresponding estimates for the backgrounds in the SRs are shown in Figure 7. More details of the event yields in the signal region defined for ±± = 300 GeV are given in Table 7. No significant excess over the expected yields is observed in any of the SRs. Table 7 includes the acceptance for the signal from pair production of ±± bosons, PP , and from associated production of ±± and ± bosons, AP . It is defined as the number of selected reconstructed events divided by the number of events at the event generation stage and represents the signal reduction due to phase-space acceptance, branching ratio and detector efficiency. The results in all the SRs are shown in the Appendix, in Table 8. The miss T distribution is shown in Figure 8 for the SRs of the ±± = 300 GeV signal mass hypothesis, where the selection requirement on miss T has been removed. The other criteria summarised in Table 5 used to define the SR are applied. Good agreement between data and the expected event yields is observed for low values of miss T , demonstrating that the background contributions are also well estimated close to the SRs, where the requirements on the kinematic variables are tighter than at preselection.
The statistical interpretation is based on a likelihood ratio test [119] using the CL s method [120]. The signal strength, a free parameter in the fit, modifies the cross section of the signal hypothesis under investigation.  The expected background and the observed data event yields in the signal region defined for the ±± = 300 GeV mass hypothesis. The signal yield is for the corresponding mass point and is normalised to the luminosity of 139 fb −1 . The displayed numbers include all sources of statistical and systematic uncertainties. The overall signal acceptances PP and AP and the observed upper limit on extra contributions to each signal region at 95% confidence level, 95 , are also presented. Selections with , three or four leptons are not affected by the electron charge-flip background, so these contributions are denoted by "-". A separate likelihood function is constructed for every signal hypothesis as the product of the Poisson probability distributions of the six channels of the corresponding signal region. Gaussian distributions are used to constrain the nuisance parameters associated with the systematic uncertainties. The widths of the Gaussian distributions correspond to the magnitudes of these uncertainties. Statistical uncertainties in the backgrounds are estimated using Poisson distributions.
The expected and observed upper limits on the ±± pair production and the ±± and ± associated production cross sections times branching fraction at 95% confidence level (CL) for the five mass hypotheses are shown in Figures 9(a) and 9(b), respectively. They are obtained from the combination of 2ℓ sc , 3ℓ and Data  Figure 9: Observed and expected upper limits on the (a) ±± pair production and (b) ±± and ± associated production cross section times branching fraction at 95% CL obtained from the combination of 2ℓ sc , 3ℓ and 4ℓ channels. The region above the observed limit is excluded by the measurement. The bands represent the expected exclusion curves within one and two standard deviations. The theoretical prediction [29] including the NLO QCD corrections [52] is also shown.
4ℓ SRs. Assuming a linear dependence of the cross-section limit between neighbouring mass points, the observed 95% CL lower limit on the mass of the ±± boson is 350 GeV for the pair production mode and 230 GeV for the associated production mode. To confirm the validity of the linear extrapolations close to the exclusion threshold, upper limits on the ±± pair production cross section times branching fraction are also computed for the complementary mass hypothesis of 350 GeV, as shown in Figure 9(a), and are found to match the extrapolated values well. A similar test was also done for the ±± and ± associated production mode, and the upper limits on the cross section times branching fraction are computed for the 220 GeV, 450 GeV and 550 GeV complementary mass hypotheses. The results are shown in Figure 9(b), and again good matching is found.
A tighter limit on the ±± boson mass is obtained for the ±± pair production mode mainly due to the cross sections times branching fraction for ±± pair production being higher than those for associated production of ±± and ± . Other important reasons are the different branching ratios of ±± ∓ (≈ 16%) and ±± ∓∓ (≈ 26%) decays into 2ℓ sc , 3ℓ and 4ℓ leptons, and the signal acceptance, which is higher for the pair production mode than for the associated production mode (see Tables 7 and 8). The definition of common signal regions also plays a role. As discussed in Section 6, the SRs are not optimal for associated production of ±± and ± . For similar reasons, the upper limits on the charged Higgs boson production cross sections times branching fraction are weaker for associated production of ±± and ± than for ±± pair production.
For associated production of ±± and ± , the limit on the ±± boson mass implies a constraint on the ∓ boson mass which is at most 5 GeV different from the ±± mass. Upper limits on the ±± pair ( ±± and ± associated) production cross section times branching fraction range from 15 fb to 5 fb (40 fb to 10 fb) when the assumed ±± mass varies from 200 GeV to 600 GeV. Observed upper limits at 95% CL on the number of BSM events for each signal region are derived using the CL s prescription and the results are shown in Table 7 for ±± = 300 GeV and for all SRs in the Appendix, in Table 8.

Conclusion
A search for pair production of ±± bosons, and for associated production of ±± and ± bosons, in the context of a type-II seesaw model, is presented. For the ±± pair production mode the ±± boson mass is at least 100 GeV lower than the ± boson mass, while for the associated production of ±± and ± bosons the mass difference between the doubly and singly charged Higgs bosons is at most 5 GeV. Only the bosonic decays ±± → ± ± and ± → ± are studied. Dedicated signal regions in the 2ℓ sc , 3ℓ and 4ℓ channels are defined as a function of different mass hypotheses within the model to look for evidence of ±± and ± bosons in the √ = 13 TeV proton-proton collision data sample collected between 2015 and 2018 with the ATLAS detector at the LHC. The data sample corresponds to an integrated luminosity of 139 fb −1 . The increase in the total integrated luminosity allowed significant improvements in the lepton fake factors measurement and the estimation of the associated uncertainties. The data are found to be in good agreement with the estimated background for all channels investigated. Combining those channels, the model considered is excluded at 95% confidence level for ±± boson masses below 350 GeV for the pair production of ±± bosons, and below 230 GeV for the associated production of ±± and ± bosons. Upper limits on the pair (associated) production cross section times branching fraction range from 15 fb to 5 fb (40 fb to 10 fb) when the assumed ±± mass varies from 200 GeV to 600 GeV. The results obtained for the ±± pair production scenario raise the exclusion limits beyond those from a previous, similar search by ATLAS using a smaller 13 TeV data set by approximately 130 GeV. In addition, the associated production of ±± and ± bosons is also explored, extending the previous search.
The expected background and the observed data event yields in the signal regions are shown in Table 8. The signal yield is for the corresponding mass point and is normalised to the luminosity of 139 fb −1 , and the displayed numbers include all sources of statistical and systematic uncertainties. The overall signal acceptances PP and AP and the observed upper limit on extra contributions to each signal region at 95% confidence level, 95 , are also presented. The data and background yields obtained for ±± = 200 GeV, 300 GeV and 400 GeV are also valid for ±± = 220 GeV, 350 GeV and 450 GeV, respectively. The data and background yields obtained for ±± = 500 GeV signal regions are also valid for ±± = 550 GeV and ±± = 600 GeV. Selections with , three or four leptons are not affected by the electron charge-flip background, so these contributions are denoted by "-".  The ATLAS Collaboration