Search for charged Higgs bosons decaying into a top quark and a bottom quark at s\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \sqrt{\mathrm{s}} $$\end{document} = 13 TeV with the ATLAS detector

A search for charged Higgs bosons decaying into a top quark and a bottom quark is presented. The data analysed correspond to 139 fb−1 of proton-proton collisions at s\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \sqrt{s} $$\end{document} = 13 TeV, recorded with the ATLAS detector at the LHC. The production of a heavy charged Higgs boson in association with a top quark and a bottom quark, pp → tbH+ → tbtb, is explored in the H+ mass range from 200 to 2000 GeV using final states with jets and one electron or muon. Events are categorised according to the multiplicity of jets and b-tagged jets, and multivariate analysis techniques are used to discriminate between signal and background events. No significant excess above the background-only hypothesis is observed and exclusion limits are derived for the production cross-section times branching ratio of a charged Higgs boson as a function of its mass; they range from 3.6 pb at 200 GeV to 0.036 pb at 2000 GeV at 95% confidence level. The results are interpreted in the hMSSM and Mh125\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {M}_h^{125} $$\end{document} scenarios.


JHEP06(2021)145
highest branching ratio for charged Higgs bosons above the top-quark mass [31]. Events with one charged lepton ( = e, µ) and jets in the final state are considered, and exclusive regions are defined according to the overall number of jets, and the number of jets tagged as containing a b-hadron. In order to separate signal from SM background, multivariate analysis (MVA) techniques combining jet multiplicities and several kinematic variables are employed in the regions where the signal rate is expected to be largest. Compared to the ATLAS result using 36 fb −1 of Run 2 data [25], improved limits on the pp → tbH + production cross-section times the H + → tb branching ratio are set by means of a simultaneous fit to the MVA classifier outputs in the different analysis regions, which determines both the contribution from the H + → tb signal and the normalisation of the backgrounds. The improvement is small at low H + mass, where the measurement is dominated by systematic uncertainties, but larger than the simple scaling with the square root of the ratio of inte-

Data and simulation samples
The data used in this analysis were recorded with the ATLAS detector at the LHC between 2015 and 2018 from √ s = 13 TeV pp collisions, and correspond to an integrated luminosity of 139 fb −1 . ATLAS [40][41][42] is a multipurpose detector with a forward-backward symmetric cylindrical geometry and a near 4π coverage in solid angle. 2 It consists of an inner tracking detector (ID) surrounded by a thin superconducting solenoid providing a 2 T axial magnetic field, electromagnetic 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. Lead/liquid-argon (LAr) sampling calorimeters provide electromagnetic (EM) energy measurements with high granularity. A steel/scintillator-tile hadron calorimeter covers the central pseudorapidity range 2 ATLAS uses a right-handed coordinate system with its origin at the nominal interaction point (IP) in the centre of the detector and the 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 .

JHEP06(2021)145
(|η| < 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. Only runs with stable colliding beams and in which all relevant detector components were functional are used.
A two-level trigger system, with the first level implemented in custom hardware and followed by a software-based second level, is used to reduce the trigger rate to around 1 kHz for offline storage [43]. Events in this analysis were recorded using single-lepton triggers. To maximise the event selection efficiency, multiple triggers were used, either with low transverse momentum (p T ) thresholds and lepton identification and isolation requirements, or with higher p T thresholds but looser identification criteria and no isolation requirements. Slightly different sets of triggers were used for 2015 and 2016-2018 data due to the increase in the average number of pp interactions per bunch crossing (pile-up). The minimum p T required by the triggers was increased to keep both trigger rate and data storage within their limits. For muons, the lowest p T threshold was 20 (26) GeV in 2015 (2016-2018), while for electrons, triggers with a minimum p T threshold of 24 (26) GeV were used [44]. Simulated events are also required to satisfy the trigger criteria.
Signal and background processes were modelled with Monte Carlo (MC) simulation samples. The pp → tbH + process followed by H + → tb decay was modelled with Mad- Graph5_aMC@NLO [31] at next-to-leading order (NLO) in QCD [45] using a four-flavour scheme (4FS) implementation with the NNPDF2.3NLO [46] parton distribution function (PDF). Parton showers (PS) and hadronisation were modelled by Pythia 8.212 [47] with a set of underlying-event (UE) parameters tuned to ATLAS data and named the A14 tune [48]. Dynamic QCD factorisation and renormalisation scales, µ f and µ r , were set to 1 3 i m(i) 2 + p T (i) 2 , where i runs over the final-state particles (H + , t and b) used in the generation. Only the H + decay into tb is considered. For the simulation of the tbH + process the narrow-width approximation was used. This assumption has negligible impact on the analysis for the models considered in this paper, as the experimental resolution is much larger than the H + natural width [49]. Interference with the SM ttbb background is neglected. A total of 18 H + mass hypotheses are used, with 25 GeV mass steps between a H + mass of 200 GeV and 300 GeV, 50 GeV steps between 300 GeV and 400 GeV, 100 GeV steps between 400 GeV and 1000 GeV, and 200 GeV steps from 1000 GeV to 2000 GeV. The step sizes were chosen to match the experimental mass resolution of the H + signal.
The production of tt + jets events was modelled using the Powheg-Box [50][51][52][53] v2 generator in the five-flavour scheme (5FS), which provides matrix elements (ME) at NLO in QCD, with the NNPDF3.0NLO PDF set [54]. The h damp parameter, which controls the transverse momentum of the first additional emission beyond the Born configuration, was set to 1.5m t [55], where m t is the mass of the top quark. Parton showers and hadronisation were modelled by Pythia 8.230 [56] with the A14 tune for the UE. The scales µ f and µ r were set to the default scale m 2 t + p 2 T,t . The sample was normalised to the Top++ 2.0 [57] theoretical cross-section of 832 +46 −51 pb, calculated at next-to-next-to-leading order (NNLO) -3 -

JHEP06(2021)145
in QCD including resummation of next-to-next-to-leading logarithmic (NNLL) soft gluon terms [58][59][60][61]. The generation of the tt + jets events was performed both inclusively of additional jet flavour, and also with dedicated filtered samples, requiring b-or c-hadrons in addition to those arising from the decays of the top quarks. Events generated with no extra b-hadrons were taken from the unfiltered sample and merged with the tt + jets events from the filtered sample, taking the appropriate cross-section and filter efficiencies into account. Single-top t-channel production was modelled using the Powheg-Box v2 generator at NLO in QCD, generated in the 4FS with the NNPDF3.0NLOnf4 PDF set [54]. The scales µ f and µ r were set to m 2 b + p 2 T,b following the recommendation in ref. [62]. Single-top tW and s-channel production was modelled using the Powheg-Box v2 generator at NLO in QCD, generated in the 5FS with the NNPDF3.0NLO PDF set. The scales µ f and µ r were set to the default scale, which is equal to the top-quark mass. For tW associated production, the diagram removal scheme [63] was employed to handle the interference with tt production [55]. All single-top events were showered with Pythia 8.230.
Production of vector bosons with additional jets was simulated with the Sherpa 2.2.1 generator [64]. Matrix elements with NLO accuracy for up to two partons, and with leading-order (LO) accuracy for up to four partons, were calculated with the Comix [65] and OpenLoops [66,67] libraries. The default Sherpa PS algorithm [68] based on Catani-Seymour dipole factorisation and the cluster hadronisation model [69] was used. It employs the dedicated set of tuned parameters developed by the Sherpa authors for this version, based on the NNPDF3.0NNLO PDF set. The NLO ME of a given jet multiplicity were matched to the PS using a colour-exact variant of the MC@NLO algorithm [70]. Different jet multiplicities were then merged into an inclusive sample using an improved CKKW matching procedure [71,72], which is extended to NLO accuracy using the MEPS@NLO prescription [73]. The merging cut was set to 20 GeV.
The production of ttV events, i.e. ttW or ttZ, was modelled using the MadGraph5_ aMC@NLO 2.3.3 generator, which provides ME at NLO in QCD with the NNPDF3.0NLO PDF set. The scales µ f and µ r were set to the default scale 1 performed in the G µ -scheme [75], ensuring an optimal description of pure electroweak interactions at the electroweak scale. All samples were generated using the NNPDF3.0NNLO PDF set, along with the dedicated set of tuned PS parameters developed by the Sherpa authors.
Other minor backgrounds (tHjb, tHW , tZq, tZW and four top quarks) were also simulated and accounted for, even though they contribute less than 1% in any analysis region. All samples and their basic generation parameters are summarised in table 1.
Most of the samples mentioned above were produced using the full ATLAS detector simulation [76] based on Geant4 [77], and the rest were produced using fast simulation [78], where the complete Geant4 simulation of the calorimeter response is replaced by a detailed parameterisation of the shower shapes, as shown in table 1. For the observables used in this analysis, the two simulations were found to give compatible results. Additional pile-up interactions, simulated with Pythia 8.186 using the A3 set of tuned parameters [55], were overlaid onto the simulated hard-scatter event. All simulation samples were reweighted such that the distribution of the number of pile-up interactions matches that of the data. In all samples the top-quark mass was set to 172.5 GeV, and the decays of b-and c-hadrons were performed by EvtGen v1.2.0 [79], except in samples simulated by the Sherpa event generator.

Object reconstruction and event selection
Charged leptons and jets, including those compatible with the hadronisation of b-quarks, are the main reconstructed objects used in this analysis. Electrons are reconstructed from -5 -JHEP06(2021)145 energy clusters in the electromagnetic calorimeter associated with tracks reconstructed in the ID [80], and are required to have |η| < 2.47. Candidates in the calorimeter transition region (1.37 < |η| < 1.52) are excluded. Electrons must satisfy the tight identification criterion described in ref.
[81], based on shower-shape and track-matching variables. Muons are reconstructed from either track segments or full tracks in the MS which are matched to tracks in the ID. Tracks are then re-fit using information from both detector systems. Muons must satisfy the medium identification criterion [82]. Muons are required to have |η| < 2.5. To reduce the contribution of leptons from hadronic decays (non-prompt leptons), both the electrons and muons must satisfy isolation criteria. These criteria include both track and calorimeter information, and have an efficiency of 90% for leptons with a p T greater than 25 GeV, rising to 99% above 60 GeV, as measured in Z → ee and Z → µµ data samples [80,82]. Finally, the lepton tracks must point to the primary vertex of the event, 3 the longitudinal impact parameter must satisfy |z 0 | < 0.5 mm and the transverse impact parameter significance must satisfy |d 0 |/σ d 0 < 5 (3) for electrons (muons).
Jets are reconstructed from three-dimensional topological energy clusters [83] in the calorimeter using the anti-k t jet algorithm [84] with a radius parameter of 0.4. Each topological cluster is calibrated to the electromagnetic scale response prior to jet reconstruction. The reconstructed jets are then calibrated with a series of simulation-based corrections and in situ techniques based on 13 TeV data [85]. After energy calibration, jets are required to have p T > 25 GeV and |η| < 2.5. Quality criteria are imposed to identify jets arising from non-collision sources or detector noise, and any event containing such a jet is removed [86]. Finally, to reduce the effect of pile-up, jets with p T < 120 GeV and |η| < 2.4 are matched to tracks with p T > 0.5 GeV, thus ensuring they originate from the primary vertex. This algorithm is known as the jet vertex tagger (JVT) [87]. To identify jets containing bhadrons, referred to as b-jets in the following, the MV2c10 tagger algorithm [88], which combines impact parameter information with the explicit identification of secondary and tertiary vertices within the jet into a multivariate discriminant, is used. Jets are b-tagged by requiring the discriminant output to be above a threshold, providing a specific b-jet efficiency in simulated tt events. A criterion with an efficiency of 70% is used to determine the b-jet multiplicity in this analysis. For this working point and for the same tt sample, the c-jet and light-flavour-quark or gluon jet (light-jets) rejection factors are 8.9 and 300, respectively [89].
To avoid counting a single detector signal as more than one lepton or jet, an overlap removal procedure is applied. First, the closest jet within ∆R y = (∆y) 2 + (∆φ) 2 = 0.2 of a selected electron is removed. 4 If the nearest jet surviving that selection is within ∆R y = 0.4 of the electron, the electron is discarded. Muons are removed if they are separated from the nearest jet by ∆R y < 0.4, which reduces the background from semileptonic decays of heavy-flavour hadrons. However, if this jet has fewer than three associated tracks, the 3 Events are required to have at least one reconstructed vertex with three or more associated tracks which have p T > 400 MeV. The primary vertex is chosen as the vertex candidate with the largest sum of the squared transverse momenta of associated tracks. 4 The rapidity is defined as y = 1 2 ln E+pz E−pz , where E is the energy and pz is the momentum component along the beam pipe. muon is kept and the jet is removed instead; this avoids an inefficiency for high-energy muons undergoing significant energy loss in the calorimeter.
The missing transverse momentum (of size E miss T ) in the event is computed as the negative vector sum of the p T of all the selected electrons, muons and jets described above, with a correction for soft energy not associated with any of the hard objects. This additional 'soft term' is calculated from ID tracks matched to the primary vertex to make it resilient to pile-up contamination [90]. The missing transverse momentum is not used in the event selection, but included in the multivariate discriminant used in the analysis.
Events are required to have exactly one electron or muon, with p T > 27 GeV, within ∆R < 0.15 of a lepton of the same flavour reconstructed by the trigger algorithm, and at least five jets, at least three of which must be b-tagged. The total event acceptance for the H + signal samples ranges from 2% (at 200 GeV) to 8.5% (at 1000 GeV). Above 1000 GeV, the acceptance decreases due to the boosted topology of the events, which fail the requirement on jet multiplicity. At 2000 GeV, the acceptance is 6%. The selected events are categorised into four separate regions according to the number of reconstructed jets (j) and b-jets (b) in the event, in order to improve the sensitivity of the fit and constrain some of the systematic uncertainties. The analysis regions are 5j3b, 5j≥4b, ≥6j3b and ≥6j≥4b, where XjYb means that X jets are found in the event, and among them Y are b-tagged. In addition, the ≥5j2b region is used to derive data-based corrections, which are implemented to improve the level of agreement between simulation and data.

Background modelling
With the isolation criteria applied both at the trigger and analysis level, as well as the purity-enhancing identification criteria used for electrons and muons (section 3), the background due to non-prompt leptons is expected to be negligible. To confirm this assumption, the ratio (N Data − N total MC )/N total MC was checked and found to not decrease when moving from a loose to a tight isolation selection. Such behaviour shows that the non-prompt-lepton background, which is not present in the simulation, provides a negligible contribution to the data, as expected given that non-prompt leptons are unlikely to be isolated in data. If data and the MC predictions differed due to a mismodelling of the other backgrounds, tighter isolation requirements would remove events in data and MC simulation alike. All backgrounds in this analysis are estimated using the simulation samples described in section 2.
To define the background categories in the likelihood fit (section 7), the tt + jets background is categorised according to the flavour of the jets in the event. Generator-level particle jets are reconstructed from stable particles (mean lifetime τ > 3 × 10 −11 s) using the anti-k t algorithm with a radius parameter R = 0.4, and are required to have p T > 15 GeV and |η| < 2.5. The flavour of a jet is determined by counting the number of bor c-hadrons within ∆R = 0.4 of the jet axis. Jets matched to one or more b-hadrons, of which at least one must have p T above 5 GeV, are labelled as b-jets; c-jets are defined analogously, only considering jets not already defined as b-jets. Events that have at least one b-jet, not including heavy-flavour jets from top-quark or W -boson decays, are labelled -7 -

JHEP06(2021)145
as tt + ≥1b; those with no b-jets but at least one c-jet are labelled as tt + ≥1c. Finally, events not containing any heavy-flavour jets, aside from those from top-quark or W -boson decays, are labelled as tt + light.
After the event selection, tt + jets constitutes the main background. It is observed that the simulation of tt + jets does not properly model high jet multiplicities nor the hardness of additional jet emissions, and data-based corrections are applied to the simulation [91,92]. Given that the additional jets in the tt + jets sample are simulated in the parton shower, the mentioned mismodelling is expected to be independent of whether the additional jets are b-tagged or not. Therefore, data and MC predictions are compared and reweighting factors are derived in a sample with at least five jets and exactly two b-tagged jets, and then applied in the 3b and ≥4b regions. The level of agreement between data and simulation in these regions improves to the point where the remaining differences are well within the model's systematic uncertainty. The reweighting factors are expressed as: , where x is the variable mismodelled by the simulation. In this context, tt + light, tt + ≥1b and tt + ≥1c, as well as W t single-top contributions, are included in the tt sample. For the range of H + masses considered in this analysis and assuming the observed upper limits on the cross-section times branching ratio published in ref.
[25], signal events contribute less than 1% to the ≥5j2b region and are neglected. Weights are calculated from the number of jets distribution (R(nJets)) first and subsequently applied to the H all T distributions 5 to derive the reweighting factors in the 5j2b, 6j2b, 7j2b and ≥8j2b regions (R(H all T )). Thus, events are weighted by the product R(nJets) × R(H all T ) depending on their jet multiplicity and H all T value. Figure 2 shows R(nJets) in the ≥5j2b region and R(H all T ) in the 5j2b, 6j2b, 7j2b and ≥8j2b regions. Among various functions tried, a hyperbola plus a sigmoid functional form was found to be the best fit to the H all T weight distributions. After the reweighting, agreement between simulation and data in the analysis regions improves, as can be seen, for example, in figure 3, which shows the leading jet's p T distribution before the fit, both before and after applying the reweighting. The final tt + ≥1b and tt + ≥1c normalisation factors and their uncertainties, which account for the remaining mismodelling observed after applying the reweighting, are not applied. These normalisations are extracted from the fit to data, as described in section 7.

Analysis strategy
To enhance the separation between signal and background, a neural network algorithm (NN) is used. Its architecture is sequential with two fully connected layers of 64 nodes, and is implemented with the Python deep learning library, Keras [93]. The activation function used is the commonly employed 'rectified linear unit' and the loss function is the 'binary cross-entropy'. Batch normalisation [94] is performed to speed up the learning process, 5 H all T is defined as the scalar sum of the transverse momenta of all jets and the lepton in the event.  dropout [95] is applied at a 10% rate, and the Adam algorithm [96] is used to optimise the parameters.
All signal samples are used in the training against all background samples, which are weighted according to their cross-sections. The training is performed separately in each analysis region, but the separate trainings include all H + mass samples, and use the value of the H + mass as a parameter [97]. For signal events the parameter corresponds to the mass of the H + sample they belong to, while for background events a random value of the H + mass, taken from the distribution of signal masses, is assigned to each event. A total of 15 variables are used in the NN: • kinematic discriminant D defined in the text, • scalar sum of the p T of all jets, • centrality calculated using all jets and leptons, • p T of the leading jet, • p T of fifth leading jet, Data Data Data Data Data Data  Comparison of the predicted leading jet p T and data before the fit in the four analysis regions before (top) and after (bottom) the reweighting was applied. The uncertainty bands include both the statistical and systematic uncertainties. Since the normalisations of the tt + ≥1b and tt + ≥1c backgrounds are allowed to vary in the fit, no cross-section uncertainties associated with these processes are included. The lower panels display the ratio of the data to the total prediction. The hatched bands show the uncertainties before the fit to the data, which are dominated by systematic uncertainties. The χ 2 /ndf and the χ 2 probability are also shown. Statistical uncertainties on MC predictions and data are uncorrelated across bins, while systematic uncertainties on the predictions are correlated.
• invariant mass of the jet triplet with maximum p T , • invariant mass of the untagged jet-pair with minimum ∆R (not in 5j≥4b), • average ∆R between all b-jet pairs in the event, • ∆R between the lepton and the pair of b-jets with smallest ∆R, • second Fox-Wolfram moment calculated using all jets and leptons [98], • number of jets (only in ≥6j3b and ≥6j≥4b regions), and • number of b-jets (only in 5j≥4b and ≥6j≥4b regions).
The variables are chosen to provide the best discrimination against the tt + ≥1b background. Among them, the kinematic discriminant, scalar sum of the p T of all jets, centrality, and leading jet p T provide the most discrimination. The centrality is computed as the scalar sum of the p T of all jets and leptons in the event divided by the sum of their  energies. The kinematic discriminant is a variable reflecting the probability that an event is compatible with the H + → tb hypothesis rather than the tt hypothesis, and is defined as where P H + (x) and P tt(x) are probability density functions for x under the signal hypothesis and background (tt) hypothesis, respectively. The event variable x indicates the set of the missing transverse momentum and the four-momenta of the reconstructed lepton and the jets [25]. Figure 4 shows the predicted NN output distributions in the four analysis regions for selected H + signal samples and the SM background. These distributions are used in a fit to extract the amount of H + signal in data. The separation of the H + signal from the background is most difficult for low H + masses because the two processes have very similar kinematics and topology. The kinematic discriminant has large separating power at low H + masses, whereas at higher masses, where the topologies of the H + signal and the background are no longer alike, other variables, such as the scalar sum of the p T of all jets, provide the largest separation.

Systematic uncertainties
Various sources of experimental and theoretical uncertainties are considered in this analysis. They may affect the overall normalisation of the processes, the shapes of the NN output distributions, or both. All the experimental uncertainties considered, with the exception of that in the luminosity, affect both normalisation and shape in all the simulated samples. Uncertainties related to the modelling of the signal and background affect both normalisation and shape, with the exception of cross-section uncertainties, which only affect the normalisation of the sample considered. Nonetheless, the normalisation uncertainties modify the relative fractions of the different samples, leading to a shape variation in the -11 -

JHEP06(2021)145
final NN output distributions. A single independent nuisance parameter (NP) is assigned to each source of systematic uncertainty in the statistical analysis. Some of the systematic uncertainties, in particular most of the experimental uncertainties, are decomposed into several independent sources. Each individual source has a correlated effect across all analysis regions and signal and background samples.
The uncertainty of the integrated luminosity for the full Run-2 dataset is 1.7% [99], obtained using the LUCID-2 detector [100] for the primary luminosity measurements. A variation in the pile-up reweighting of the simulated events described in section 2 is included to cover the uncertainty in the ratio of the predicted and measured inelastic cross-sections in a given fiducial volume [101].
Uncertainties associated with charged leptons arise from the trigger selection, the lepton reconstruction, identification and isolation criteria, as well as the lepton momentum scale and resolution. The reconstruction, identification and isolation efficiency of electrons and muons, as well as the efficiency of the trigger used to record the events, differ slightly between data and simulation, which is compensated for by dedicated correction factors (CFs). Efficiency CFs are measured using tag-and-probe techniques in Z → + − data and simulated samples [82,102], and are applied to the simulation to correct for the differences. The effect of these CFs, as well as of their uncertainties, are propagated as corrections to the MC event weight. Additional sources of uncertainty originate from the corrections applied to adjust the lepton momentum scale and resolution in the simulation to match those in data. The impact of these uncertainties on σ(pp → tbH + )×B(H + → tb) is smaller than 1%.
Uncertainties associated with jets arise from the efficiency of pile-up rejection by the JVT, from the jet energy scale (JES) and resolution (JER), and from b-tagging. Correction factors are applied to correct for differences between data and MC simulation for JVT efficiencies. These CFs are estimated using Z(→ µ + µ − ) + jets with tag-and-probe techniques similar to those in ref. [87]. The JES and its uncertainty are derived by combining information from test-beam data, collision data and simulation [85]. Additional uncertainties are considered, related to jet flavour, quark/gluon fraction, pile-up corrections, η dependence, high-p T jets, and differences between full and fast simulation. The JER uncertainty is obtained by combining dijet balance measurements in data and simulation [103]. The b-tagging efficiencies in simulated samples are corrected to match efficiencies in data. Correction factors are derived as a function of p T for b-, c-and light-jets separately in dedicated calibration analyses. For b-jet efficiencies, tt events in the di-lepton topology are used, exploiting the very pure sample of b-jets arising from the decays of the top quarks [89]. For c-jet mistag rates, tt events in single-lepton topology are used, exploiting the c-jets from the hadronically decaying W bosons, using techniques similar to those in ref. [104]. For light-jet mistag rates, the so-called negative-tag method similar to that in ref.
[105] is used, but using Z+jets events instead of di-jet events.
All the uncertainties described above on energy scales or resolutions of the reconstructed objects are propagated to the missing transverse momentum. Additional uncertainties in the scale and resolution of the soft term are considered, which account for the disagreement between data and simulation of the p T balance between the hard and the -12 -

JHEP06(2021)145
soft components. A total of three independent sources are added: an offset along the hard component p T axis, and the smearing resolution along and perpendicular to this axis [106, 107]. Since the missing transverse momentum is not used in selection but only in the event reconstruction, the associated uncertainties have an impact smaller than 1%.
The uncertainty in the H + signal due to different scale choices is estimated by varying µ f and µ r up and down by a factor of two. The uncertainties from the modelling of the PDF are evaluated replacing the nominal NNPDF2.3NLO PDF set by a symmetrised Hessian set, PDF4LHC15_nlo_30, following the PDF4LHC recommendations for LHC Run 2 [108].
The modelling of the tt + jets background is one of the largest sources of uncertainty in the analysis, and several different components are considered. The 6% uncertainty for the inclusive tt production cross-section predicted at NNLO+NNLL includes effects from varying µ f and µ r , the PDF, α S , and the top-quark mass [109]. This uncertainty is applied to tt + light only, since the normalisation of tt + ≥1b and tt + ≥1c are allowed to vary freely in the fit. Besides normalisation, the tt + light, tt + ≥1b and tt + ≥1c processes are affected by different types of uncertainties: the uncertainties associated with additional Feynman diagrams for the tt + light are constrained from relatively precise measurements in data [110]; tt + ≥ 1b and tt + ≥ 1c can have similar or different Feynman diagrams depending on the flavour scheme used for the PDF, and the different masses of the b-and the c-quarks contribute to additional differences between these two processes. For these reasons, all uncertainties in the tt + jets background modelling are assigned independent NP for tt + light, tt + ≥1b and tt + ≥1c. Systematic uncertainties in the acceptance and shapes are extracted by comparing the nominal prediction with alternative MC samples or settings. Such comparisons would significantly change the fractions of tt + ≥1b and tt + ≥1c. However, since the normalisation of these sub-processes in the analysis regions is determined in the fit, these alternative predictions are reweighted in such a way that they keep the same fractions of tt + ≥1b and tt + ≥1c as the nominal sample in the phase-space selected by the analysis.
The uncertainty due to initial state radiation (ISR) is estimated by using the Var3cUp (Var3cDown) variant from the A14 tune [48], corresponding to α ISR The nominal PowhegBox+Pythia sample is compared with the PowhegBox+Herwig sample to assess the effect of the PS and hadronisation models, and to the MadGraph5_aMC@NLO sample to assess the effect of the NLO matching technique. The nominal PowhegBox+Pythia 5FS prediction for the dominant tt + ≥ 1b background, in which all the additional partons are produced by the PS, is compared with an alternative PowhegBox+Pythia 4FS sample, in which the bb pair is generated in addition to the tt pair at the ME level. An uncertainty resulting from the comparison of the shapes of the two models is included. Finally, the weights derived in section 4 to improve the agreement of the simulation with data are varied within their statistical uncertainties, in a correlated way among the three tt + jets components. All the sources of systematic uncertainty for the tt + jets modelling are summarised in  A 5% uncertainty is considered for the cross-sections of the three single-top production modes [111][112][113][114][115]. Uncertainties associated with the PS and hadronisation model, and with the NLO matching scheme are evaluated by comparing, for each process, the nominal PowhegBox+Pythia sample with a sample produced using Powheg-Box+Herwig [116] and MadGraph5_aMC@NLO+Pythia, respectively. As mentioned in section 4, the W t single-top mode is included in the reweighting procedure, and thus the same uncertainties used for tt are applied here. The uncertainty associated to the interference between W t and tt production at NLO [63] is assessed by comparing the nominal PowhegBox+Pythia sample produced using the "diagram removal" scheme with an alternative sample produced with the same generator but using the "diagram subtraction" scheme [62,63].
The predicted SM ttH signal cross-section uncertainty is +5.8% −9.2% (QCD scale) ±3.6% (PDF + α S ) [13,[117][118][119][120][121]. Uncertainties of the Higgs boson branching ratios amount to 2.2% for the bb decay mode [13]. For the ISR and FSR, the amount of radiation is varied following the same procedure as for tt. The nominal PowhegBox+Pythia sample is compared with the PowhegBox+Herwig sample to assess the uncertainty due to PS and hadronisation, and to the MadGraph5_aMC@NLO sample for the uncertainty due to the NLO matching.
The uncertainty of the ttV NLO cross-section prediction is 15%, split into PDF and scale uncertainties as for ttH [13,122]. An additional ttV modelling uncertainty, related to the choice of PS and hadronisation model and matching scheme, is assessed by comparing the nominal MadGraph5_aMC@NLO+Pythia samples with alternative samples generated with Sherpa.
An overall 50% normalisation uncertainty is considered for the four-top-quarks background, covering effects from varying µ f , µ r , PDF and α S [45,123]. The small background

JHEP06(2021)145
tZq is assigned a 7.9% and a 0.9% uncertainty accounting for µ f and µ r scales and PDF variations, respectively. Finally, a single 50% uncertainty is used for tZW [45].
An uncertainty of 40% is assumed for the W +jets normalisation, with an additional 30% for W + heavy-flavour jets, taken as uncorrelated between events with two and more than two heavy-flavour jets. These uncertainties are based on variations of the µ f and µ r scales and of the matching parameters in the Sherpa samples. An uncertainty of 35% is applied to the Z+jets normalisation, uncorrelated across jet bins, to account for both the variations of the scales and matching parameters in the Sherpa samples and the uncertainty in the extraction from data of the correction factor for the heavy-flavour component [54,124]. Finally, a 50% normalisation uncertainty in the diboson background is assumed, which includes uncertainties in the inclusive cross-section and additional jet production [125].

Results
A binned maximum-likelihood fit to the data is performed simultaneously on the NN output distributions in the four analysis regions, and each mass hypothesis is tested separately. A profile-likelihood-ratio test and the CL S method [126][127][128] are used to quantify the level of agreement with the background-only hypothesis or background-plus-signal hypothesis and to determine exclusion limits. The parameter of interest is the product of the production cross-section σ(pp → tbH + ) and the branching ratio B(H + → tb). In addition, two unconstrained factors, for which no prior knowledge is assumed, are used to normalise the tt + ≥1b and tt + ≥1c backgrounds. These normalisation factors range from 1.2 to 1.6 (0.2 to 1.8) with a typical uncertainty of 0.2 (0.6) for the tt + ≥1b (tt + ≥1c) background, depending on the H + mass hypothesis used in the fit. The fitted tt + heavy-flavour jets backgrounds in the different fits are consistent within two standard deviations. All systematic uncertainties are implemented as nuisance parameters with log-normal constraint terms. There are about 170 nuisance parameters considered in the fit, the number varying slightly across the range of mass hypotheses. A summary of the systematic uncertainties with similar sources grouped together is given in table 3. Depending on the particular H + mass hypothesis, the total systematic uncertainty is dominated by the uncertainties in the modelling of the tt background, in particular tt + ≥1b and tt + ≥1c, jet flavour-tagging uncertainties, and jet energy scale and resolution. Table 4 shows the event yields after the background-plus-signal fit under the 200 GeV and 800 GeV H + mass hypotheses. The corresponding post-fit distributions of the NN output in each analysis region are shown in figure 5 for the 200 GeV and 800 GeV H + mass hypotheses. After the fit, good agreement between the data and simulation is found in the input variables to the NN. No significant excess above the expected SM background is observed in all regions and mass intervals and upper limits on the cross-section times branching ratio are derived as function of the H + mass.
The 95% confidence level (CL) upper limits on σ(pp → tbH + ) × B(H + → tb) obtained using the CL S method are presented in figure 6. Uncertainties in the predicted H + cross-       neutralinos (M 125 h (χ)) or light staus (M 125 h (τ )). In both cases a charged Higgs boson of sufficiently high mass is allowed to decay into the supersymmetric particles. Finally, the value of tan β in both the M 125 h (alignment) scenario, characterised by one of the two neutral CP-even scalars having couplings like those of the SM Higgs boson, and the M 125 h 1 (CPV) scenario, which includes CP violation in the Higgs sector, is already constrained to be in the range 1-20 by previous searches at the LHC [36]. Uncertainties in the predicted H + cross-sections or branching ratios are not included in the limits. For all scenarios except the hMSSM, Higgs boson masses and mixing (and effective Yukawa couplings) have been calculated with the code FeynHiggs [129][130][131][132][133][134][135]. Whereas in the hMSSM the branching ratios are computed solely with HDECAY [136,137], all other scenarios combine the most precise results of FeynHiggs, HDECAY and PROPHECY4f [138,139].
In the context of these scenarios, tan β values below 1 are observed to be excluded at 95% CL for H + masses between 200 and

Conclusion
A search for charged Higgs bosons is presented using a data sample corresponding to an integrated luminosity of 139 fb −1 from proton-proton collisions at √ s = 13 TeV, recorded with the ATLAS detector at the LHC. The search for pp → tbH + is performed in the H + mass range 200-2000 GeV. A neural network that combines jet multiplicities and several kinematic variables is built in the regions where the signal rate is expected to be largest. The neural network output depends on the H + mass, and a fit to the data is performed simultaneously on the neural network output distributions in the analysis regions, separately for each signal mass hypothesis.
No significant excess above the expected Standard Model background is found and observed (expected) upper limits at 95% confidence level are set on the σ(pp → tbH + ) production cross-section times the branching ratio B(H + → tb), which range from σ × B = 3.6 (2.6) pb at m H + = 200 GeV to σ × B = 0.036 (0.019) pb at m H + = 2 TeV. Compared to the previous ATLAS search for tbH + production followed by H + → tb decays with 36 fb −1 , the observed σ × B limits improved by 5% to 70%, depending on the H + mass, apart from the lowest one. In the low H + mass region the measurement is dominated by systematic uncertainties, while in the high H + mass region the use of tighter lepton triggers and refined b-tagging techniques, along with the H + mass-independent training of the neural network, leads to an improvement beyond the simple scaling with the square root of the ratio of integrated luminosities.
In the context of the hMSSM and several M 125 h scenarios, some values of tan β, in the range 0.5-2.1, are excluded for H + masses between 200 and 1200 GeV. For H + masses between ∼200 and ∼750 GeV, values of tan β > 34 are also excluded. Compared to previous results of the same search channel, this analysis excludes a broader region of large tan β. Additionally, an extended region of low tan β and low and high H + masses is also excluded.  [105] ATLAS collaboration, Calibration of light-flavour b-jet mistagging rates using ATLAS proton-proton collision data at √ s = 13 TeV, ATLAS-CONF-2018-006 (2018). [