Search for doubly charged Higgs boson production in multi-lepton final states with the ATLAS detector using proton–proton collisions at \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\sqrt{s}=13\,\text {TeV}$$\end{document}s=13TeV

A search for doubly charged Higgs bosons with pairs of prompt, isolated, highly energetic leptons with the same electric charge is presented. The search uses a proton–proton collision data sample at a centre-of-mass energy of 13 TeV corresponding to 36.1 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {fb}^{-1}$$\end{document}fb-1 of integrated luminosity recorded in 2015 and 2016 by the ATLAS detector at the LHC. This analysis focuses on the decays \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H^{\pm \pm }\rightarrow e^{\pm }e^{\pm }$$\end{document}H±±→e±e±, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H^{\pm \pm }\rightarrow e^{\pm }\mu ^{\pm }$$\end{document}H±±→e±μ± and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H^{\pm \pm }\rightarrow \mu ^{\pm }\mu ^{\pm }$$\end{document}H±±→μ±μ±, fitting the dilepton mass spectra in several exclusive signal regions. No significant evidence of a signal is observed and corresponding limits on the production cross-section and consequently a lower limit on \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$m(H^{\pm \pm })$$\end{document}m(H±±) are derived at 95% confidence level. With \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\ell ^{\pm }\ell ^{\pm }=e^{\pm }e^{\pm }/\mu ^{\pm }\mu ^{\pm }/e^{\pm }\mu ^{\pm }$$\end{document}ℓ±ℓ±=e±e±/μ±μ±/e±μ±, the observed lower limit on the mass of a doubly charged Higgs boson only coupling to left-handed leptons varies from 770 to 870 GeV (850 GeV expected) for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$B(H^{\pm \pm }\rightarrow \ell ^{\pm }\ell ^{\pm })=100\%$$\end{document}B(H±±→ℓ±ℓ±)=100% and both the expected and observed mass limits are above 450 GeV for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$B(H^{\pm \pm }\rightarrow \ell ^{\pm }\ell ^{\pm })=10\%$$\end{document}B(H±±→ℓ±ℓ±)=10% and any combination of partial branching ratios.


Introduction
Events with two prompt, isolated, highly energetic leptons with the same electric charge (same-charge leptons) are produced very rarely in a proton-proton collision according to the predictions of the standard model (SM), but may occur with higher rate in various theories beyond the standard model (BSM). This analysis focuses on BSM theories that contain a doubly charged Higgs particle H ±± using the observed invariant mass of same-charge lepton pairs. In the absence of evidence for a signal, lower limits on the mass of the H ±± particle are set at the 95% confidence level.
Doubly charged Higgs particles can couple to either lefthanded or right-handed leptons. In LRS models, two cases are distinguished and denoted H ±± . Besides the leptonic decay, the H ±± particle can decay into a pair of W bosons as well. For low values of the Higgs triplet vacuum expectation value v , it decays almost exclusively to leptons while for high values of v the decay is mostly to a pair of W bosons [9,12]. In this analysis, the coupling to W bosons is assumed to be negligible and only pair production via the Drell-Yan process is considered. The Feynman diagram of the production mechanism is presented in Fig. 1.
The analysis targets only decays of the H ±± particle into electrons and muons, denoted by . Other final states X that are not directly selected in this analysis are taken into account by reducing the lepton multiplicity of the final state. These states X would include, for instance, τ leptons or W bosons, as well as particles which escape detection. The total assumed branching ratio of H ±± is therefore B(H ±± → e ± e ± ) + B(H ±± → e ± μ ± ) + B(H ±± → μ ± μ ± ) + B(H ±± → X ) = B(H ±± → ± ± ) + B(H ±± → X ) = 100%. Moreover, the decay width is assumed to be negligible compared to the detector resolution, which is compatible with theoretical predictions. Two-, three-, and four-lepton signal regions are defined to select the majority of such events. These regions are further divided into unique flavour categories (e or μ) to where at least one of the lepton pairs is e ± e ± , e ± μ ± , or μ ± μ ± increase the sensitivity. The partial decay width of H ±± to leptons is given by: with k = 2 if both leptons have the same flavour ( = ) and k = 1 for different flavours. The factor h has an upper bound that depends on the flavour combination [23,24]. In this analysis, only prompt decays of the H ±± bosons (cτ < 10μm) are considered, corresponding to h 1.5 × 10 −6 for m(H ±± ) = 200 GeV. In general, there is no preference for decays into τ leptons, as the coupling is not proportional to the lepton mass like it is for the SM Higgs boson.
Additional motivation to study cases with B(H ±± → ± ± ) < 100% is given by type-II see-saw models with specific neutrino mass hypotheses resulting in a fixed branching ratio combination [13,25,26] which does not necessarily correspond to B(H ±± → ± ± ) = 100%.
The ATLAS Collaboration previously analysed data corresponding to 20.3 fb −1 of integrated luminosity which were recorded in 2012 at a centre-of-mass energy of 8 TeV [27]. This study resulted in the most stringent lower limits on the mass of a potential H ±± L particle. Depending on the flavour of the final-state leptons, the observed limits vary between 465 and 550 GeV assuming B(H ±± L → ± ± ) = 100%. The analysis presented in this paper extends the one described in Ref. [27] and is based on 36.1 fb −1 of integrated luminosity collected in 2015 and 2016 at a centre-of-mass energy of 13 TeV. A similar search has also been performed by the CMS Collaboration [28].

ATLAS detector
The ATLAS detector [29] at the LHC is a multi-purpose particle detector with a forward-backward symmetric cylindrical geometry and an almost 4π coverage in solid angle. 1 It consists of an inner tracking detector (ID) surrounded by a thin superconducting solenoid providing a 2 T axial magnetic field, electromagnetic (EM) and hadronic calorimeters, and a muon spectrometer. The inner tracking detector covers the pseudorapidity range |η| < 2.5. It is composed of silicon pixel, silicon micro-strip, and transition radiation tracking detectors. A new innermost layer of pixel detectors [30] was installed prior to the start of data taking in 2015. Lead/liquid-argon (LAr) sampling calorimeters provide electromagnetic energy measurements with high granularity. A hadronic (steel/scintillator-tile) calorimeter covers the central pseudorapidity range (|η| < 1.7). The end-cap and forward regions are instrumented with LAr calorimeters for both EM and hadronic energy measurements up to |η| = 4.9. The muon spectrometer surrounds the calorimeters and features three large air-core toroidal superconducting magnets with eight coils each. The field integral of the toroids ranges between 2 to 6 Tm across most of the detector. The muon system includes precision tracking chambers and fast detectors for triggering. A two-level trigger system is used to select events [31] that are interesting for physics analyses. The first-level trigger is implemented as part of the hardware. Subsequently a software-based high-level trigger executes algorithms similar to those used in the offline reconstruction software, reducing the event rate to about 1 kHz.

Dataset and simulated event samples
The data used in this analysis were collected at centre-ofmass energy of 13 TeV during 2015 and 2016, and correspond to an integrated luminosity of 3.2 fb −1 in 2015 and 32.9 fb −1 in 2016. The average number of pp interactions per bunch crossing in the dataset is 24. Interactions other than the hardscattering one are referred to as pile-up. The uncertainty on the combined 2015 and 2016 integrated luminosity is 3.2%. Following a methodology similar to the one described in Ref.
[32], this uncertainty is derived from a preliminary calibration of the luminosity scale using x-y beam-separation scans performed in August 2015 and May 2016.
Signal candidate events in the electron channel are required to pass a dielectron trigger with a threshold of 17 GeV on the transverse energy (E T ) of each of the electrons. Candidate events in the muon channel are selected Footnote 1 continued z-axis along the beam pipe. The x-axis points from the IP to the centre of the LHC ring, and the y-axis points upwards. Cylindrical coordinates (r, φ) are used in the transverse plane, φ being the azimuthal angle around the z-axis. The pseudorapidity is defined in terms of the polar angle θ as η = −ln tan(θ/2). Angular distance is measured in units of R ≡ ( η) 2 + ( φ) 2 . Rapidity is defined as where E denotes the energy and p z is the momentum component along the beam direction. using a combination of two single-muon triggers with transverse momentum ( p T ) thresholds of 26 and 50 GeV. The single-muon trigger with the lower p T threshold also requires track-based isolation of the muon according to the isolation criteria described in Ref. [33]. Events containing both electrons and muons (mixed channel) are required to pass either the combined electron-muon trigger or any of the triggers used for the muon channel or the electron channel. The combined trigger has an E T threshold of 17 GeV for the electron and a p T threshold of 14 GeV for the muon. Events with four leptons are selected using a combination of dilepton triggers. In general, single-lepton triggers are more efficient than dilepton triggers. However, single-electron triggers impose stringent electron identification criteria, which interfere with the data-driven background estimation.
An irreducible background originates from SM processes resulting in same-charge leptons, hereafter referred to as prompt background. Prompt background and signal model predictions were obtained from Monte Carlo (MC) simulated event samples which are summarised in Table 1. Prompt background events mainly originate from diboson (W ± W ± / Z Z / W Z) and tt X processes (tt W , tt Z, and tt H). They also provide a source of reducible background due to charge misidentification in channels that contain electrons. 2 As described in Sect. 5, the modelling of charge misidentification in simulation deviates from data and consequently charge reconstruction scale factors are derived in a datadriven way and applied to the simulated events to compensate for the differences. The highest-yield process which enters the analysis through charge misidentification is Drell-Yan (qq → Z /γ * → + − ) followed by tt production. MC samples are in general normalised using theoretical crosssections referenced in Table 1. However, yields of some MC samples are considered as free parameters in the likelihood fit, as described in Sect. 7.
Another source of reducible background arises from events with non-prompt electrons or muons or with other physics objects misidentified as electrons or muons, collectively called 'fakes'. For both, electrons and muons, this contribution originates within jets, from decays of light-flavour or heavy-flavour hadrons into light leptons. For electrons, a significant component of fakes arises from jets which satisfy the electron reconstruction criteria and from photon conversions. MC samples are not used to estimate this background because the simulation of jets and hadronisation has large uncertainties. Instead, a data-driven approach is used to assess this contribution from production of W +jets, tt and 2 The probability of muon charge misidentification is negligible because muon tracks are measured both in the inner detector and in the muon spectrometer which provides a much larger lever arm for the curvature measurement. NLO [45] multi-jet events. The method is validated in specialised validation regions. The SM Drell-Yan process was modelled using Powheg-Box v2 [37][38][39] interfaced to Pythia 8.186 [34] for parton showering. The CT10 set of parton distribution functions (PDF) [40] was used to calculate the hard scattering process. A set of tuned parameters called the AZNLO tune [42] was used in combination with the CTEQ6L1 PDF set [51] to model non-perturbative effects. Photos++ version 3.52 [52] was used for photon emissions from electroweak vertices and charged leptons. The generation of the process was divided into 19 samples with subsequent invariant mass intervals to guarantee a good statistical coverage over the entire mass range.
Higher-order corrections were applied to the Drell-Yan simulated events to scale the mass-dependent cross-section computed at next-to-leading order (NLO) in the strong coupling constant with the CT10 PDF set to next-to-next-toleading order (NNLO) in the strong coupling constant with the CT14NNLO PDF set [41]. The corrections were calculated with VRAP [53] for QCD effects and Mcsanc [54] for electroweak effects. The latter are corrected from leadingorder (LO) to NLO.
A sample of Z → ee events was generated with Sherpa 2.2.1 [50], in addition to the Powheg prediction, to measure the probability of electron charge misidentification, as explained in Sect. 5. The electron p T spectrum is a crucial ingredient for the estimate of this probability and was found to be better described by Sherpa than by Powheg, especially for invariant masses of the electron pair close to the Z boson mass. Sherpa uses Comix [55] and OpenLoops [56] to calculate the matrix elements up to two partons at NLO and up to four partons at LO in the strong coupling constant. The merging with the Sherpa parton shower [57] follows the ME+PS@NLO prescription in [58].
The tt process was generated with the NLO QCD event generator Powheg-Box v2 which was interfaced to Pythia 8.186 for parton showering. The A14 parameter set [36] was used together with the NNPDF2.3 [35] PDF set for tuning the shower. Furthermore, the PDF set used for generation was NNPDF3.0 [43]. Additionally, top-quark spin correlations were preserved through the use of MadSpin [59]. The predicted tt production cross-section is 832 +20 −30 (scale) ±35 (PDF + α S ) pb as calculated with Top++2.0 [44] to NNLO in perturbative QCD, including soft-gluon resummation to next-to-next-to-leading-log order. The top-quark mass was assumed to be 172.5 GeV. The scale uncertainty results from independent variations of the factorisation and renormalisation scales, while the second uncertainty is associated with variations of the PDF set and α S , following the PDF4LHC [60] prescription using the MSTW2008 68% CL NNLO [61], CT10 NNLO [62], and NNPDF2.3 PDF sets.
Single-top-quark events produced in W t final states were generated by Powheg-Box v2 with the CT10 PDF set used in the matrix element calculations. Single-top-quark events in other final states were generated by Powheg-Box v1. This event generator uses the four-flavour scheme for the NLO QCD matrix element calculations together with the fixed four-flavour PDF set CT10f4. The parton shower, hadronisation, and underlying event were simulated with Pythia 6.428 [46] using the CTEQ6L1 PDF set and the corresponding Perugia 2012 tune (P2012) [47]. The top-quark mass was set to 172.5 GeV. The NLO cross-sections used to normalise these MC samples are summarised in Ref. [45].
The tt W , tt Z, and tt H processes were generated at LO with MadGraph v2.2.2 [63] and MadGraph v2.3.2 using the NNPDF2.3 PDF set. Pythia 8.186 was applied for shower modelling configured with the A14 tune [36], as explained in more detail in Ref. [64]. They were normalised using theoretical cross-sections summarised in Ref. [49].
Diboson processes with four charged leptons, three charged leptons and one neutrino, or two charged leptons and two neutrinos were generated with Sherpa 2.2.1, using matrix elements containing all diagrams with four electroweak vertices. They were calculated for up to three partons at LO accuracy and up to one (4 , 2 +2ν) or zero partons (3 +1ν) at NLO QCD using Comix and OpenLoops. The merging with the Sherpa parton shower [57] follows the ME+PS@NLO prescription. The NNPDF3.0NNLO [43] PDF set was used in conjunction with dedicated parton shower tuning by the Sherpa authors.
Diboson processes with one boson decaying hadronically and the other one decaying leptonically were predicted by Sherpa 2.1.1 [50]. They were calculated for up to three additional partons at LO accuracy and up to one (Z Z) or zero (W W , W Z) additional partons at NLO using Comix and OpenLoops matrix element generators. The merging with the Sherpa parton shower [57] follows the ME+PS@NLO prescription. The CT10 PDF set was used in conjunction with a dedicated parton shower tuning. The Sherpa 2.1.1 diboson prediction was scaled by 0.91 to account for differences between the internal electroweak scheme used in this Sherpa version and the G μ scheme which is the common default. Similarly, loop-induced diboson production with both gauge bosons decaying fully leptonically was simulated with Sherpa 2.1.1. The prediction is at LO accuracy while up to one additional jet is merged with the matrix element.
Additional diboson samples for W Z and Z Z production were generated with Powheg-Box v2 to estimate theoretical uncertainties. Pythia 8.186 provided the parton shower. The CT10 PDF set was used for the matrix element calculation while the parton shower was configured with the CTEQL1 PDF set. The non-perturbative effects were modelled using the AZNLO [42] tune. particles. This choice resulted in a good statistical coverage for all possible decay channels. The production of the H ±± was implemented only via the Drell-Yan process. Originally, the cross-section at √ s = 14 TeV was calculated with NLO accuracy by the authors of Ref. [9]. Subsequently, a rescaling to √ s = 13 TeV with the CTEQ6 PDF [65] set was provided. The cross-sections and corresponding K -factors are summarised in Table 2.
Since this analysis exclusively targets the leptonic decays of the H ±± bosons, the vacuum expectation value of the neutral component of the left-handed Higgs triplet (v L ) was set to zero in order to exclude H ±± → W W decays. The decay width of the H ±± particle to leptons depends on the h couplings. These were set to the value h = 0.02 in all Pythia 8.186 samples. This setting corresponds to a decay width that is negligible compared to the detector resolution. The h τ and h τ τ couplings were fixed at zero. There are 23 MC samples with different H ±± particle masses, starting from 200 GeV up to 1300 GeV in steps of 50 GeV. The ATLAS detector is expected to have the best H ±± mass resolution in the electron-electron final states. Resolutions around 30 GeV for masses of 200-500 GeV and 50-100 GeV for higher masses can be achieved with the event selection defined in Sect. 4. Furthermore, the H ±± mass resolution in electron-muon final states varies from 50 to 150 GeV and from 50 to 200 GeV in muon-muon final states.
For all simulated samples except those obtained with Sherpa, the EvtGen v1.2.0 program [66] was used to model bottom and charm hadron decays. The effect of the pile-up was included by overlaying minimum-bias collisions, simu-lated with Pythia 8.186, on each generated signal and background event. The number of overlaid collisions is such that the distribution of the average number of interactions per pp bunch crossing in the simulation matches the pile-up conditions observed in the data. The pile-up simulation is described in more detail in Ref. [67].
The response of the ATLAS detector was simulated using the Geant 4 toolkit [68]. Data and simulated events were reconstructed with the default ATLAS software [69] while simulated events were corrected with calibration factors to better match the performance measured in data.

Event reconstruction and selection
Events are required to have at least one reconstructed primary vertex with at least two associated tracks with p T > 400 MeV. Among all the vertices in the event the one with the highest sum of squared transverse momenta of the associated tracks is chosen as the primary vertex.

Event reconstruction
This analysis classifies leptons in two exclusive categories called tight and loose, defined specifically for each lepton flavour as described below. Leptons selected in the tight category feature a predominant component of prompt leptons, while loose leptons are mostly fakes, which are used for the fake-background estimation. All tracks associated with lepton candidates must have a longitudinal impact parameter with respect to the primary vertex of less than 0.5 mm.
Electron candidates are reconstructed using information from the EM calorimeter and ID by matching an isolated calorimeter energy deposit to an ID track. They are required to have |η| < 2.47, p T > 30 GeV, and to pass at least the LHLoose identification level based on a multivariate likelihood discriminant [70,71]. The likelihood discriminant is based on track and calorimeter cluster information. Electron candidates within the transition region between the barrel and endcap electromagnetic calorimeters (1.37 < |η| < 1.52) are vetoed due to limitations in their reconstruction quality. The track associated with the electron candidate must have an impact parameter evaluated at the point of closest approach between the track and the beam axis in the transverse plane (d 0 ) that satisfies is the uncertainty on d 0 . In addition to this, electron candidates are classified as tight if they satisfy the LHMedium working point of the likelihood discriminant and the isolation criteria described in Ref. [70]. This is based on calorimeter cluster and track isolation, which vary to obtain a fixed efficiency for selecting prompt electrons of 99% across p T and η. Electrons are classified as loose if they fail to satisfy either of the identification or the isolation criteria.
Muon candidates are selected by combining information from the muon spectrometer and the ID. They satisfy the medium quality criteria described in Ref. [33] and are required to have p T > 30 GeV, |η| < 2.5 and |d 0 |/σ (d 0 ) < 10. Muon candidates are classified as tight if their impact parameter satisfies |d 0 |/σ (d 0 ) < 3.0 and they satisfy the most stringent isolation working point of the cutbased track isolation [70]. Muons are loose if they fail the isolation requirement.
Jets or particles originating from the hadronisation of partons are reconstructed by clustering energy deposits in the calorimeter calibrated at the EM scale. The anti-k t algorithm [72] is used with a radius parameter of 0.4, which is implemented with the FastJet [73] package. The majority of pile-up jets are rejected using the jet-vertex-tagger [74], which is a combination of track-based variables providing discrimination against pile-up jets. For all jets the expected average transverse energy contribution from pile-up is subtracted using an area-based p T density subtraction method and a residual correction derived from the MC simulation, both detailed in Refs. [75,76]. In this analysis, events containing jets identified as originating from b-quarks are vetoed. They are identified with a multivariate discriminant [76] that has a b-jet efficiency of 77% in simulated tt events and a rejection factor of ≈ 40 (≈ 20) for jets originating from gluons and light quarks (c-quarks).
After electron and muon identification, jet calibration, and pile-up jet removal, overlaps between reconstructed particles or jets are resolved. First, electrons are removed if they share a track with a muon. Secondly, ambiguities between electrons and jets are resolved. If a jet is closer than ( y) 2 + ( φ) 2 = 0.2 the jet is rejected. If 0.2 < ( y) 2 + ( φ) 2 < 0.4 the electron is removed. Finally, if a muon and a jet are closer than ( y) 2 + ( φ) 2 = 0.4, and the jet features less than 3 tracks, the jet is removed. Otherwise the muon is discarded.

Event selection
In this search, events are classified in independent categories, called analysis regions, which serve different purposes. The so-called control regions are used to constrain free background parameters in the statistical analysis detailed in Sect. 7. The background model is validated against data in validation regions. Both the control and validation regions are designed to reject signal events. A dedicated selection targeting signal events is utilised to define the signal regions. The selection criteria utilised for each region are summarised in Table 3. The main variable that defines the type of the region is the invariant mass of same-charge lepton pairs. Invariant masses are required to be above 200 GeV in signal regions and below 200 GeV in most control and validation regions.
The lepton multiplicity in the event is used to define the analysis regions. Events with two or three leptons are required to contain exactly one same-charge lepton pair, while four-lepton events are required to feature two samecharge pairs where the sum of all lepton charges has to be zero. An exception is the opposite-charge control region (OCCR) where exactly two electrons with opposite charge are required. In all regions, events with at least one b-tagged jet are vetoed, in order to suppress background events arising from top-quark decays. In regions with more than two leptons, events are rejected if any opposite-charge sameflavour lepton pair is within 10 GeV of the Z boson mass (81.2 GeV < m( + − ) < 101.2 GeV). This requirement is applied to reject diboson events featuring a Z boson in the final state, and is inverted in diboson control regions, where at least one Z boson is present. Furthermore, the Z boson veto is not applied in four-lepton control and validation regions to increase the available number of simulated diboson events.
The invariant mass of the same-charge lepton pair is used in the final fit of the analysis for the two-and three-lepton regions. In the OCCR, the invariant mass of the oppositecharge lepton pair is used. A lower bound of 60 GeV on the invariant mass is imposed in all regions to discard low-mass events which would potentially bias the background estimation of the analysis while maximising the available number of events.
In the electron and mixed channels the lower bound is increased to 90 GeV in the three-lepton regions and to 130 GeV in the two-lepton regions. The motivation for increasing the lower mass bound in regions containing electrons is the data-driven charge misidentification background correction, where the Z → ee peak is used to measure the charge misidentification rates (described in Sect. 5). Differences between data and MC simulation in the dielectron Table 3 Summary of all regions used in the analysis. The table is  Electron channel

Fig. 2
Dielectron mass distributions for opposite-charge (black) and same-charge (red) pairs for data (filled circles) and MC simulation (continuous line). The latter includes a correction for charge misidentification. The hatched band indicates the statistical error and the luminosity uncertainty summed in quadrature applied to MC simulated events same-charge Z → ee peak (see Fig. 2) were minimised by construction following the methodology described in Sect. 5, and the Z → ee peak was therefore not used in the fit. In the two-lepton regions, this bound is set to 130 GeV to completely remove the Z peak region. In the three-lepton regions, where this effect is not as strong, the bound is relaxed to 90 GeV to reduce the statistical uncertainty of the sample. As the charge misidentification background is not present in the muon channel, there is no need to increase the lower mass bound there.
In the mixed channel, events are further divided into two categories, where the same-charge pair features differentflavour leptons or not, indicated by e ± μ ± ∓ and e ± e ± μ ∓ or μ ± μ ± e ∓ , respectively.
In order to maximise the sensitivity in two-lepton and three-lepton signal regions (SR1P2L and SR1P3L), additional requirements are imposed on same-charge lepton pairs, regardless of the flavour. These exploit both the boosted decay topology of the H ±± resonance and the high energy of the decay products. The same-charge lepton separation is required to be R( ± , ± ) < 3.5. Their combined transverse momentum has to be p T ( ± ± ) > 100 GeV. 3 Finally, the scalar sum of the leptons' transverse momenta is required to be above 300 GeV in the signal regions. In SR1P2L and SR1P3L, the signal selection efficiency combined with the detector acceptance varies greatly with the assumed branching ratio into light leptons. It is the highest for B(H ±± → ± ± ) ≈ 60% where about 40% of signal events are selected either in SR1P2L or SR1P3L. For B(H ±± → ± ± ) = 100%, about 25% of signal events are selected in either of the regions.
In the four-lepton signal region (SR2P4L), the fit variable is the average invariant mass of the two same-charge lepton pairsM ≡ (m ++ + m −− )/2. A selection on the variable M/M ≡ |m ++ − m −− |/M is applied to reject background where the two same-charge pairs have inconsistent invariant masses. The M/M requirement is optimised for different flavour combinations which generally feature different mass resolutions. This selection corresponds to M values which are required to be below 15-50 GeV for M = 200 GeV, 30-160 GeV forM = 500 GeV, and 50-500 GeV forM = 1000 GeV. In the 2P4L signal region, the fraction of signal events that are selected is approximately 50% for the B(H ±± → ± ± ) = 100% case and lower for branching ratios into light leptons below 100%.
The same-charge validation region (SCVR) is used to validate the data-driven fake-background estimation and the charge misidentification effect in the electron channel. The three-lepton validation region (3LVR) is used to validate the SM diboson background and fake events with three reconstructed leptons with different proportions across channels. The four-lepton validation region (4LVR) is used to validate the diboson modelling in the four-lepton region. Furthermore, the diboson control region (DBCR) is used to constrain the diboson background yield in each channel while the opposite-charge control region is used to constrain the Drell-Yan contribution in the electron channel only. The four-lepton control region (4LCR) is used to constrain the yield of the diboson background in four-lepton regions.

Background composition and estimation
Prompt SM backgrounds in all regions are estimated using the simulated samples listed in Sect. 3. Prompt light leptons are defined as leptons originating from Z , W , and H boson decays or leptons from τ decays if the τ has a prompt source (e.g. Z → τ τ ). MC events containing at least one nonprompt or fake selected tight or loose lepton are discarded to avoid an overlap with the data-driven fake-background estimation. Prompt electrons in the remaining simulated events are corrected to account for different charge misidentification probabilities in data and simulation.
Electron charge misidentification is caused predominantly by bremsstrahlung. The emitted photon can either convert to an electron-positron pair, which happens in most of the cases, or traverse the inner detector without creating any track. In the first case, the cluster corresponding to the initial electron can be matched to the wrong-charge track, or most of the energy is transferred from one track to the other because of the photon. In case of photon emission without subsequent pair production, the electron track has usually very few hits only in the silicon pixel layers, and thus a short lever arm on its curvature. Because the electron charge is derived from the track curvature, it could be incorrectly determined while the electron energy is likely appropriate as the emitted photon deposits all of its energy in the EM calorimeter as well. For a similar reason high-energy electrons are more often affected by charge misidentification, as their tracks are approximately straight and therefore challenging for the curvature measurement. The modelling of charge misidentification in simulation deviates from data due to the complex processes involved, which particularly rely on a very precise description of the detector material. A correction is obtained by comparing the charge misidentification probability measured in data to the one in simulation. The charge misidentification probability is extracted by performing a likelihood fit on a dedicated Z → ee data sample (see Fig. 2 denoting the expected number of same-charge pairs in bin (i, j), where i and j indicate the kinematic configuration of the two electrons in the pair, given the charge misidentification probabilities P i and P j . N i j SC is the measured number of same-charge pairs. The formula for the negative log likelihood used in the likelihood fit is given in Eq. 1: The charge misidentification probability is parameterised as a function of electron p T and η, P( p T , η) = σ ( p T )× f (η).
The binned values, σ ( p T ) and f (η), are free parameters in the likelihood fit. To ensure the proper normalisation of P( p T , η), the area of the distribution describing f (η) was The fake-lepton background is estimated with a datadriven approach, the so-called 'fake factor' method, as described in Ref. [27]. The b-jet veto significantly reduces fake leptons from heavy-flavour decays. Most of the fake leptons still passing the analysis selection originate from inflight decays of mesons inside jets, jets misreconstructed as electrons, and conversions of initial-and final-state radiation photons. The fake factor method provides an estimation of events with fake leptons in analysis regions by extrapolating the yields from the so-called 'side-band regions'. For each analysis region a corresponding side-band region is defined. It requires exactly the same selection and lepton multiplicity except that at least one lepton must fail to satisfy the tight identification criteria. The ratio of tight to loose leptons is measured in dedicated 'fake-enriched regions'. It is determined as a function of lepton flavour, p T , and η, and referred to as the 'fake factor' (F( p T , η, flavour)). It describes the probability for a fake lepton to be identified as a tight lepton. The definitions of the fake-enriched regions for the electron and muon channels are reported in Table 4. In the measurement of the fake factor, a requirement on the unbalanced momentum in the transverse plane of the event, Table 4 Selection criteria defining the fake-enriched regions used to measure the ratio of the numbers of tight and loose leptons, the so-called fake factor, for the electron and muon channels E miss T , is imposed to reject W + jets events and to further enrich the regions with fake leptons. The fake factor method relies on the assumption that no prompt leptons appear in the fake-enriched samples. This assumption is not fully correct with the imposed selection. Therefore, the number of residual prompt leptons in the fake-enriched regions is estimated using simulation and subtracted from the numbers of tight and loose leptons used to measure the fake factors.
The number of events in the analysis regions containing at least one fake lepton, N fake , is estimated from the side-bands. Data are weighted with fake factors according to the loose lepton multiplicity of the region: with N data SB denoting the number of data events in the sideband, N L ,i is the loose lepton multiplicity in the i-th event of the side-band region and l indicates the loose lepton. The contamination of prompt leptons in the side-band region is subtracted using simulated events, denoted by N MC SB . Dedicated two-lepton and three-lepton validation regions, defined in Table 3, are used to verify the data-driven fakelepton estimation in regions as similar to the signal regions as possible. They are designed to contain only a negligible number of signal events. Orthogonality between signal and validation regions is ensured by requiring the invariant mass of the same-charge lepton pair m( ± ± ) to be less than 200 GeV in the validation regions. Furthermore, diboson modelling and the electron charge misidentification backgrounds are tested. Each background estimation is validated in the corresponding regions, defined to be enriched in the given contribution.

Systematic uncertainties
Several sources of systematic uncertainty are accounted for in the analysis. These correspond to experimental and theoretical sources affecting both background and signal predictions. All considered sources of systematic uncertainty affect the total event yield, and all except the uncertainties on the Distribution of dilepton mass for data and SM background predictions in three-lepton validation regions: a the three-electron validation region, b the three-muon validation region, c the 3LVR with an electron-muon same-charge pair (e ± μ ± ∓ ), and d the 3LVR with a same-flavour same-charge pair (e ± e ± μ ∓ or μ ± μ ± e ∓ ). The hatched bands include all systematic uncertainties post-fit, with the correlations between various sources taken into account luminosity and cross section also affect the distributions of the variables used in the fit (Sect. 7). The cross-sections used to normalise the simulated samples are varied to account for the scale and PDF uncertainties in the cross-section calculation. The variation is 6% for diboson production [77], 13% for tt W production, 12% for tt Z production, and 8% for tt H production [49]. The theoretical uncertainty in the Drell-Yan background is estimated by PDF eigenvector variations of the nominal PDF set, variations of PDF scale, α S , electroweak corrections, and photon-induced corrections. The effect of the PDF choice is considered by comparing the nominal PDF set to several others, namely CT10NNLO [62], MMHT14 [78], NNPDF3.0 [43], ABM12 [79], HERAPDF2.0 [80,81], and JR14 [82]. An envelope is constructed by taking into account the largest deviations from the nominal choice. The predominant prompt background, arising from diboson production, is assigned an additional theoretical uncertainty by comparing the nominal Sherpa 2.2.1 prediction with the Powheg prediction. This uncertainty varies from 5 to 10%. Furthermore, the theoretical uncertainty in the NLO cross-section for pp → H ++ H −− is reported to be about 15% [9]. It includes the renormalization and factorization scale dependence and the uncertainty in the parton densities. Lastly, the theoretical uncertainty in the simulated pp → H ++ H −− events is assessed by varying the A14 parameter set in Pythia 8.186 and choosing alternative PDFs CTEQ6L1 and CT09MC1 [83]. The impact on the signal acceptance is found to be negligible.  Fig. 6 Relative uncertainties in the total background yield estimation after the fit. 'Stat. Unc.' corresponds to reducible and irreducible background statistical uncertainties. 'Yield fit' corresponds to the uncertainty arising from fitting the yield of diboson and Drell-Yan backgrounds. 'Lumi' corresponds to the uncertainty in the luminosity. 'Theory' indicates the theoretical uncertainty in the physics model used for simulation (e.g. cross-sections). 'Exp.' indicates the uncertainty in the simulation of electron and muon efficiencies (e.g. trigger, identification). 'Fakes' is the uncertainty associated with the model of the fake background. Individual uncertainties can be correlated, and do not necessarily add in quadrature to the total background uncertainty, which is indicated by 'Total Unc.' A significant contribution arises from the statistical uncertainty in the MC samples and data sideband regions. Analysis regions have a very restrictive selection and only a small fraction of the initially generated MC events remains after applying all requirements. The statistical uncertainty varies from 5 to 40% depending on the signal region. Experimental systematic uncertainties due to different reconstruction, identification, isolation, and trigger efficiencies of leptons in data compared to simulation are estimated by varying the corresponding scale-factors. They are at most 3% and less significant than the other systematic uncertainties and MC statistical uncertainties. The same is true for lepton energy or momentum calibration.
The experimental uncertainty related to the charge misidentification probability of electrons arises from the statistical uncertainty of both the data and the simulated sample of Z /γ * → ee events used to measure this probability. The uncertainty ranges between 10 and 20% as a function of the electron p T and η. Possible systematic effects were investigated by altering the selection requirements imposed on the invariant mass used to select Z /γ * → ee events analysed to measure the misidentification probability. The effects estimated with this method are found to be negligible compared to the statistical uncertainty.
The experimental systematic uncertainty in the datadriven estimate of the fake-lepton background is evaluated by varying the nominal fake factor to account for different effects. The E miss T requirement is altered to consider variations in the W + jets composition. The flavour composition of the fakes is investigated by requiring an additional recoil-ing jet in the electron channel and changing the definition of the recoiling jet in the muon channel. Furthermore, the transverse impact parameter criterion for tight muons (defined in Sect. 4.1) is varied by one standard deviation. Finally, in the fake-enriched regions, the normalisation of the subtracted simulated samples, to remove the prompt lepton component, is altered within its uncertainties. This accounts for uncertainties related to the luminosity, the cross-section, and the corrections applied to simulation-based predictions. The statistical uncertainty in the fake factors is added in quadrature to the total systematic error. The uncertainty ranges between 10% and 20% across all p T and η bins.
The total relative systematic uncertainty after the fit (Sect. 7), and its breakdown into components, is presented in Fig. 6. All experimental systematic uncertainties discussed here affect the signal samples as well as the background.

Statistical analysis and results
The statistical analysis package HistFitter [84] was used to implement a maximum-likelihood fit of the dilepton invariant mass distribution in all control and signal regions, and thē M distribution in four-lepton regions to obtain the numbers of signal and background events. The likelihood is the product of a Poisson probability density function describing the observed number of events and Gaussian distributions to constrain the nuisance parameters associated with the systematic uncertainties. The widths of the Gaussian distributions correspond to the magnitudes of these uncertainties, whereas between various sources taken into account. The notation ± ± ∓ indicates that the same-charge leptons have different flavours and ± ± ∓ indicates that same-charge leptons have the same flavour, while the opposite-charge lepton has a different flavour Table 5 The number of predicted background events in control regions after the fit, compared to the data. Uncertainties correspond to the total uncertainties in the predicted event yields, and are smaller for the total than the sum of the components in quadrature due to correla-tions between these components. Due to rounding the totals can differ from the sums of components. Background processes with a negligible yield are marked with the en dash (-)  Table 6 The number of predicted background events in two-lepton and four-lepton validation regions (top) and three-lepton validation regions (bottom) after the fit, compared to the data. Uncertainties correspond to the total uncertainties in the predicted event yields, and are smaller for the total than the sum of the components in quadrature due to correlations between these components. Due to rounding the totals can differ from the sums of components. Background processes with a negligible yield are marked with the en dash (-)   Table 7 The number of predicted background events in two-lepton and four-lepton signal regions (top) and three-lepton signal regions (bottom) after the fit, compared to the data. Uncertainties correspond to the total uncertainties in the predicted event yields, and are smaller for the total than the sum of the components in quadrature due to correlations between these components. Due to rounding the totals can differ from the sums of components. Background processes with a negligible yield are marked with the en dash (-)   Poisson distributions are used for MC simulation statistical uncertainties. Furthermore, additional free parameters are introduced for the Drell-Yan and the diboson background contributions, to fit their yields in the analysis regions. Fitting the yields of the largest backgrounds reduces the systematic uncertainty in the predicted yield from SM sources. The fitted normalisations are compatible with their SM predictions within the uncertainties. The diboson yield is described by four free parameters, each corresponding to a different diboson region: electron channel, muon channel, mixed channel, and the four-lepton channel. After the fit, the compatibility between the data and the expected background was assessed. For various branching ratio assumptions, 95% CL upper limits were set on the pp → H ++ H −− cross-section using the CL s method [85].

Fit results
The observed and expected yields in all control, validation, and signal regions used in the analysis are presented in Fig. 7 and summarised in Tables 5, 6, 7. No significant excess is observed in any of the signal regions. Correlations between various sources of uncertainty are evaluated and used to estimate the total uncertainty in the SM background prediction. Two-and four-lepton signal regions are presented in Fig. 8 and three-lepton signal regions are presented in Fig. 9. In the four-lepton signal region only one data event is observed.
It is an e + μ + e − μ − event with invariant masses of 228 and 207 GeV for the same-charge lepton pairs. The likelihood fit to the two-, three-, and four-lepton control and signal regions was designed to fully exploit the pair Observed 95% CL limit Expected 95% CL limit Observed 95% CL limit Expected 95% CL limit   production of the H ±± boson with its boosted topology and lepton multiplicity. For B(H ±± → ± ± ) = 100% the production cross-section is excluded down to 0.1 fb, corresponding to 3-4 signal events, which is the theoretical limit of a 95% CL exclusion. Some representative cross-section upper limits as a function of the H ±± boson mass are presented in Fig. 10, for different combinations of the branching ratios for decay into light-lepton pairs. The final result of the fit is a lower limit on the twodimensional grid of the H ±± boson mass for any combination of light lepton branching ratios that sum to a certain value. The fit was performed for values of B(H ±± → ± ± ) from 1% to 5% in 1% intervals, and from 10% to 100% in 10% intervals. Expected limits for B(H ±± → ± ± ) = 100% are presented in Fig. 11 for H ±± L and in Fig. 12 for H ±± R . Results of the fit are presented in Figs. 13 and 14 for H ±± L and H ±± R , respectively. Here, three specific decay scenarios to only e ± e ± , μ ± μ ± , and e ± μ ± , are considered and the minimum limit for each value of B(H ±± → ± ± ) is given. The minimum limit is obtained by taking, for each value of B(H ±± → ± ± ), the least stringent limit for any combination of branching ratios that sum to B(H ±± → ± ± ). The lower mass limits for these four cases are similar, which indicates that the analysis is almost equally sensitive to each decay channel.

Conclusion
The ATLAS detector at the Large Hadron Collider was used to search for doubly charged Higgs bosons in the same-charge dilepton invariant mass spectrum at high values, using e ± e ± , e ± μ ± and μ ± μ ± final states as well as final states with three or four leptons (electrons and/or muons). The search was performed with 36.1 fb −1 of data from proton proton collisions at √ s = 13 TeV, recorded during the 2015 and 2016 data-taking periods. No significant excess above the Standard Model prediction was found. As a result of the search, lower limits are set on the mass of doubly-charged Higgs bosons. These vary between 770 and 870 GeV for the H ±± L mass and for B(H ±± → ± ± ) = 100% and above 450 GeV for B(H ±± → ± ± ) ≥ 10% for any combination of partial branching ratios. The observed lower limits on the H ±± R mass vary from 660 to 760 GeV for B(H ±± → ± ± ) = 100% and are above 320 GeV for B(H ±± → ± ± ) ≥ 10%.