Search for dark matter produced in association with a Standard Model Higgs boson decaying into b-quarks using the full Run 2 dataset from the ATLAS detector

The production of dark matter in association with Higgs bosons is predicted in several extensions of the Standard Model. An exploration of such scenarios is presented, considering final states with missing transverse momentum and b-tagged jets consistent with a Higgs boson. The analysis uses proton-proton collision data at a centre-of-mass energy of 13 TeV recorded by the ATLAS experiment at the LHC during Run 2, amounting to an integrated luminosity of 139 fb−1. The analysis, when compared with previous searches, benefits from a larger dataset, but also has further improvements providing sensitivity to a wider spectrum of signal scenarios. These improvements include both an optimised event selection and advances in the object identification, such as the use of the likelihood-based significance of the missing transverse momentum and variable-radius track-jets. No significant deviation from Standard Model expectations is observed. Limits are set, at 95% confidence level, in two benchmark models with two Higgs doublets extended by either a heavy vector boson Z′ or a pseudoscalar singlet a and which both provide a dark matter candidate χ. In the case of the two-Higgs-doublet model with an additional vector boson Z′, the observed limits extend up to a Z′ mass of 3 TeV for a mass of 100 GeV for the dark matter candidate. The two-Higgs-doublet model with a dark matter particle mass of 10 GeV and an additional pseudoscalar a is excluded for masses of the a up to 520 GeV and 240 GeV for tan β = 1 and tan β = 10 respectively. Limits on the visible cross-sections are set and range from to 0.05 fb to 3.26 fb, depending on the missing transverse momentum and b-quark jet multiplicity requirements.


Introduction
Various astrophysical observations based on gravitational interactions [1,2] strongly support the existence of dark matter (DM) which interacts through neither the strong nor the electromagnetic force.However, the Standard Model of particle physics (SM) provides no suitable DM candidate particle.There are many complementary search strategies for DM including direct-detection [3][4][5][6][7] and indirect-detection experiments [8] as well as searches at particle colliders.
Since DM particles do not interact electromagnetically or strongly even direct-detection experiments have low efficiencies.Therefore, instead of attempting to detect them directly, any DM particles produced in proton-proton collisions at the Large Hadron Collider [9] (LHC) would be deduced from an imbalance in the transverse momentum measured in that collision event ( miss T ).This means that DM particles can only be detected if they are produced in association with visible particles.When there is only one such particle, this gives rise to event topologies referred to as 'mono-' final states, where  refers to the visible particle.Prominent examples of these topologies are the mono-jet [10-12], mono-/ [13-15] and mono-photon [16][17][18] final states.
Mono- topologies are typically dominated by cases where the visible particle is produced from initial-state radiation (ISR), as the couplings between the SM particle  and the initial-state quarks or gluons are much larger than its couplings to the final-state DM particles.A counterexample is the case where the visible particle is a SM Higgs boson (named the mono-Higgs signature).Given that the coupling of the Higgs boson to light quarks and gluons is highly suppressed, a Higgs boson is more likely to be produced through final-state radiation (FSR) or as part of the same process that produces the DM particles.This means that this topology is only sensitive to models where the Higgs boson couples directly to DM or some other beyond-the-SM (BSM) particle involved in DM production.However, in these cases the DM-SM interaction is probed directly [19], potentially providing more information about the structure of the DM-SM coupling in the event of a discovery.There are also many models where the coupling between BSM particles and the Higgs boson is enhanced, for instance models where DM is connected to electroweak symmetry breaking [20,21] or where DM particles couple to the SM only through the Higgs sector (Higgs portal models) [22].These features make the mono-Higgs signature an important part of the LHC DM search programme.
The two-Higgs-doublet model (2HDM) [23] extends the SM with a second Higgs doublet.This predicts a total of five Higgs bosons after mixing: two charged scalars  ± , two neutral CP-even scalars ℎ and , and one neutral CP-odd scalar .Two simplified benchmark signal models are used: the  -2HDM [24] and 2HDM+ [25,26].In all models considered here the mass of the lighter neutral CP-even scalar ℎ is required to match that of the Higgs boson observed at the LHC and the Yukawa couplings are defined according to the Type II 2HDM.
The  -2HDM has an additional heavy vector boson, denoted  , whose coupling to quarks   and mass   are free parameters.The production mechanism for the mono-Higgs signature is shown in Figure 1.DM is introduced into the model as a new fermion which couples to the CP-odd scalar  with a coupling strength denoted by   .The coupling strength and the DM particle mass   are treated as free parameters.This model is used mainly as a benchmark for high-mass resonances.
The 2HDM+ scenario is the simplest renormalisable and gauge-invariant extension of a simplified pseudoscalar mediator model.It adds a new pseudoscalar singlet which mediates the interactions between the SM and a singlet fermion  identified as the DM candidate.The coupling between the pseudoscalar and DM singlets, denoted   , and the mass of the DM   are free parameters.This singlet mixes with the pseudoscalar  from the two Higgs doublets, with the mixing angle  and the mass of the resulting pseudoscalar  being free parameters of the model.A major advantage of the 2HDM+ scenario over simpler models is that it generates a wider variety of experimental signatures which can provide complementary exclusion regions from different types of experiment.There are two main production mechanisms for the mono-Higgs signature in this model, as shown in Figure 2. In the type-II 2HDM considered in this paper the coupling between down-type quarks and the  boson scales with tan , the ratio of the vacuum expectation values of the two Higgs doublets.This means that for low tan  values (tan  5) the gluon-gluon fusion (F) mechanism shown in Figure 2(a) dominates, whereas for higher tan  values the -associated production ( ) shown in Figure 2(b) is dominant.Signal grids are generated where each of the two production mechanisms are used exclusively.For each grid a tan  value is chosen that ensures that the corresponding production mechanism is dominant: tan  = 1 for the F grid and tan  = 10 for the   grid.
Similar analyses were performed using data taken during the years 2015-2016 by ATLAS [27] and CMS [28].Other mono-Higgs analyses were also performed on the same datasets in final states where the Higgs boson decays into a pair of photons in ATLAS [29] or either a pair of photons, a  +  − pair in CMS [30] or a pair of  or  bosons [28].Beyond the large increase in integrated luminosity (from 36 fb −1 to 139 fb −1 ) a number of analysis improvements extend the sensitivity beyond the previous ATLAS search.
In the previous ATLAS search, events were required to have either one or two -jets, whereas in this search the events are divided into regions with either exactly two or at least three -jets.Introducing the exclusive three -jet category improves the sensitivity to the 2HDM+   production mechanism which was not considered in the previous ATLAS search, while the one -jet category does not provide any significant improvement due to large backgrounds with high uncertainties.
The analysis also benefits from using particle-flow objects for jet reconstruction [31], neural-network based -jet [32] and -lepton [33] identification, and variable-radius track-jets [34] to identify boosted Higgs boson candidates.In order to reduce backgrounds containing fake  miss T the likelihood-based  miss T significance S [35] is used.The event selections were reoptimised resulting in an improved sensitivity, especially to highly boosted signals.

ATLAS detector
The ATLAS detector [36] at the LHC covers nearly the entire solid angle around the collision point. 1 It consists of an inner tracking detector surrounded by a thin superconducting solenoid, electromagnetic and hadronic calorimeters, and a muon spectrometer incorporating three large superconducting toroidal magnets.
The inner-detector system (ID) is immersed in a 2 T axial magnetic field and provides charged-particle tracking in the range || < 2.5.The high-granularity silicon pixel detector covers the vertex region and typically provides four measurements per track, the first hit normally being in the insertable B-layer installed before Run 2 [37,38].It is followed by the silicon microstrip tracker, which usually provides eight measurements per track.These silicon detectors are complemented by the transition radiation tracker (TRT), which enables radially extended track reconstruction up to || = 2.0.The TRT also provides electron identification information based on the fraction of hits (typically 30 in total) above a higher energy-deposit threshold corresponding to transition radiation.
The calorimeter system covers the pseudorapidity range || < 4.9.Within the region || < 3.2, electromagnetic calorimetry is provided by barrel and endcap high-granularity lead/liquid-argon (LAr) calorimeters, with an additional thin LAr presampler covering || < 1.8 to correct for energy loss in material upstream of the calorimeters.Hadronic calorimetry is provided by the steel/scintillator-tile calorimeter, segmented into three barrel structures within || < 1.7, and two copper/LAr hadronic endcap calorimeters.The solid angle coverage is completed with forward copper/LAr and tungsten/LAr calorimeter modules optimised for electromagnetic and hadronic measurements, respectively.
The muon spectrometer (MS) comprises separate trigger and high-precision tracking chambers measuring the deflection of muons in a magnetic field generated by the superconducting air-core toroids.The field integral of the toroids ranges between 2.0 and 6.0 T m across most of the detector.A set of precision chambers covers the region || < 2.7 with three layers of monitored drift tubes, complemented by cathode-strip chambers in the forward region, where the background is highest.The muon trigger system covers the range || < 2.4 with resistive-plate chambers in the barrel, and thin-gap chambers in the endcap regions.
Interesting events are recorded by the first-level trigger system implemented in custom hardware, followed by selections made by algorithms implemented in software in the high-level trigger [39].The first-level 1 ATLAS uses a right-handed coordinate system with its origin at the nominal interaction point (IP) in the centre of the detector and the -axis along the beam pipe.The -axis points from the IP to the centre of the LHC ring, and the -axis points upwards.Cylindrical coordinates (, ) are used in the transverse plane,  being the azimuthal angle around the -axis.
The pseudorapidity is defined in terms of the polar angle  as  = − ln tan(/2).Angular distance is measured in units of Simulated event samples corresponding to the 2HDM+ signal were generated at LO using M -G 5_aMC@NLO v2.6.7 interfaced to P 8.244 with the A14 tune.Samples were generated separately for the F and   production modes, with the former generated in the 4-flavour scheme setting tan  = 1 and the latter in the 5-flavour scheme setting tan  = 10.Following the recommendations of the LHC Dark Matter Working Group [25], the mass of the DM candidate was set to   = 10 GeV, the Yukawa coupling between the DM candidate and the pseudoscalar  was set to   = 1 and the Higgs quartic couplings were set to  3 =  1 =  2 = 3.The chosen value of   ensures that the  →   branching ratio is significant for all values of   used.The pseudoscalar mixing angle was set to sin  = 0.35 and the alignment limit was assumed.The samples were generated with varying   between 250 and 2000 GeV and   between 100 and 600 GeV.
For all signal samples the masses of the heavy Higgs bosons were considered degenerate (  =   =   ± ).The mass of the lightest CP-even Higgs boson was set to match that of the Higgs boson discovered at the LHC, i.e.  ℎ = 125 GeV [49].
Background events from the production of a single weak vector boson ( = , ) in association with jets or of a pair of weak bosons (diboson) were simulated using S 2.2.1 [50], with S 2.2.2 used for -initiated diboson production.For the +jets samples, next-to-leading-order (NLO) matrix elements for up to two jets and LO matrix elements for up to four jets were calculated using the C [51] and O L [52,53] libraries.For the  q-initiated diboson samples, NLO matrix elements for up to one additional jet and LO matrix elements for up to three additional jets were used, while the -initiated processes were generated using LO matrix elements for up to one additional jet.The samples were matched with the S parton shower [54] using the MEPS@NLO prescription [55][56][57][58].
The /+ℎ samples were generated using P B v2.The cross-sections of the  q-initiated processes were calculated at NNLO QCD and NLO EW accuracy, using the P MiNLO procedure [78,79].The cross-sections of the -initiated processes were calculated at NLO+next-to-leading-logarithm (NLL) accuracy in QCD [80][81][82].
The  t samples were generated using M G 5_aMC@NLO v2.3.3 at NLO.The cross-sections were calculated at NLO QCD and EW accuracies as provided by Ref. [76].
All samples were generated using the NNPDF3.0NLOparton distribution function (PDF) set [83] apart from the S /+jets and diboson samples, which were generated using the NNPDF3.0NNLOPDF set along with a dedicated tune developed by the S authors, and the t-channel single-top-quark production samples, which were generated using the NNPDF3.0NLONF4PDF set.The P B samples were interfaced with P 8.230 for the parton shower and hadronisation.Of these, the  t, single-top-quark and  tℎ samples used the NNPDF2.3LOPDF set [83] and the A14 tune [84], while the /+ℎ samples used the CTEQ6L1 PDF set [85] and the AZNLO tune [86].For the P B top-quark samples the decays of -and -hadrons were simulated using E G 1.6.0[87].The  t samples were interfaced with P 8.210 using the same tune and PDF set as the P B  t, single-top-quark and  tℎ samples, and using E G 1.2.0 for the decays of -and -hadrons.
In order to simulate the effect of additional   collisions in the same and neighbouring bunch crossings (pile-up) all samples were overlaid with multiple   collisions simulated with P 8.186 using the NNPDF2.3LOPDF set and the A3 tune [88].The response of the detector was modelled with a detector simulation [89] based on G 4 [90].

Object definitions Primary Vertex
Primary vertices are constructed using at least two ID tracks with  T > 500 MeV [91].The primary vertex with the largest sum of squared track transverse momenta (  2 T ) is selected as the hard-scatter vertex, henceforth only referred to as the primary vertex.

Jets
Jets are reconstructed using the anti-  algorithm [92,93].The analysis considers three types of jets to better match the different event topologies.Small-radius (small-) jets are constructed from particle-flow objects formed from ID tracks and calorimeter energy clusters [31] using a radius parameter of  = 0.4.This radius parameter is designed to capture jets initiated by a gluon, light quark or -quark.Small- jets are classified as central (|| < 2.5) or forward (2.5 < || < 4.5).Central small- jets are required to have  T > 20 GeV and forward small- jets  T > 30 GeV.In order to remove the impact of jets predominantly formed from particles from pile-up vertices, central small- jets are required to pass the 'Tight' jet vertex tagger (JVT) [94] working point (WP). 2n topologies where the Higgs boson decay ℎ →  b cannot be resolved into two small- jets, large-radius (large-) jets with a radius parameter of  = 1.0 are used, constructed from calorimeter energy clusters calibrated using the local hadronic cell weighting (LCW) scheme [95].This radius parameter is chosen so that a single large- jet should capture all jets produced in the decay of a boosted heavy object, such as a Higgs boson.To reduce the impact of pile-up, these jets are then 'trimmed', removing any  = 0.2 subjets which have less than 5% of the original jet energy [96].In order to identify subjets originating from -hadrons within the large- jets, jets are also constructed from ID tracks, using a variant of the anti-  algorithm with a radius parameter that shrinks as the  T of the proto-jet increases [34].These are referred to as variable-radius (variable-) track-jets and are matched to the large- jets by ghost association [97].The radius parameter is set to  = 30 GeV/ T , with minimum and maximum values of 0.02 and 0.4, respectively.The reduced radius at high  T allows the algorithm to reconstruct separate jets from closely spaced -hadrons, such as in highly boosted ℎ →  b decays.
Both the small- and large- jet energies are calibrated using a sequence of simulation-derived corrections.Small- jets additionally have an area-based energy subtraction applied to reduce the impact of pile-up, as well as a series of additional data-derived corrections [98].
Central small- jets and variable- track-jets containing -hadrons are identified using the DL1 tagger [32].This multivariate algorithm uses the impact parameters of ID tracks as well as information about secondary vertices and reconstructed flight paths of -and -hadrons within the jet.For both classes of jet a WP is chosen which tags jets containing -hadrons with 77% efficiency in  t events.The decays of these -hadrons can produce muons which are vetoed when building particle-flow objects and therefore not included in the energies of either the small-or large- jets.In order to correct for this, the four-momenta of non-isolated muons falling inside these jet cones can be added into the jet, improving the resolution of their four-momenta.For small- (large-) jets, this is done for the muon (two muons) closest to the jet axis.This correction is only used when calculating the mass of the Higgs boson candidate  ℎ and was shown in Ref. [99] to improve the resolution of this measurement.Correcting  ℎ in this way improves the ability of the fit to separate the signal from the major backgrounds.

Leptons
Leptons are divided into 'baseline' and 'signal' categories.Events with baseline leptons are vetoed in the analysis search regions and signal leptons are used to define control regions to constrain background components.
Electrons are reconstructed from a track which is coincident with a cluster built from energy deposits in the calorimeter [100].They are then identified using a multivariate likelihood technique, using several features including the shape of the measured shower, the track quality and the distribution of energy within the calorimeter [101].For this analysis, the 'LooseAndBLayer' WP [101] is used for both the signal and baseline electrons.Isolation selections are also applied to distinguish between electrons produced in the initial collision or decays of / bosons or -leptons (prompt) and those produced in decays of other objects [101].Requirements are placed on the energy of calorimeter clusters and the  T of tracks measured in isolation cones around the electron.For signal electrons with  T < 200 GeV and all baseline electrons the total energy of clusters within Δ = 0.2 of the electron, excluding the electron cluster, must be less than 20% of the  T of the electron.The total  T of tracks matched to the primary vertex that lie within a cone whose size is set to the smaller of Δ = 10 GeV/ T and 0.2, excluding the electron track, must be less than 15% of the  T of the electron.For signal electrons with  T > 200 GeV the total energy of clusters within Δ = 0.2 of the electron is required to be less than the smaller of 0.015 ×  T and 3.5 GeV.For these electrons, no track-based isolation selection is applied.Both the baseline and signal electrons are required to have || < 2.47.Baseline electrons are required to have  T > 7 GeV and signal electrons  T > 27 GeV.In order to ensure that they are compatible with the primary vertex, the track from which the electron is reconstructed is required to have ( 0 ) < 5 and | 0 sin | < 0.5 mm, where ( 0 ) is the significance of the transverse impact parameter,  0 is the longitudinal impact parameter, and  is the polar angle of the track.
Muons are reconstructed by matching track segments formed in the MS to a track from the ID [102].Identification is performed through selections on the qualities of the tracks used in the reconstruction, as well as their compatibility, for example, in the measurements of  T in the MS and ID.Similarly to electrons, isolation selections are also applied.Baseline muons are required to pass the 'Loose' identification WP and signal muons are required to pass the 'Medium' identification WP [102].For baseline muons, the total energy of clusters within Δ = 0.2 of the muon is required to be less than 30% of the  T of the muon.For signal muons, the total  T of tracks within Δ = 0.2 of the muon's primary track is required to be less than 1.25 GeV.Both the baseline and signal muons are required to satisfy || < 2.5, ( 0 ) < 3 and | 0 sin | < 0.5 mm.Baseline muons are required to have  T > 7 GeV and signal muons  T > 25 GeV.
Hadronically decaying -lepton reconstruction is seeded from  = 0.4 anti-  jets built using the LCWcalibrated clusters [103].As hadronic -lepton decays yield either one or three charged pions the jets are required to have either one or three tracks within Δ = 0.2 of the jet axis.A recurrent neural network (RNN) classifier is used to identify the -leptons [33].The inputs to the RNN are built from the clusters and tracks associated with the -lepton.All -leptons are required to pass the 'VeryLoose' WP [33] and have || < 2.5 and  T > 20 GeV.As there is no dedicated -lepton control region, no signal -lepton selection is defined.

Overlap removal
In order to avoid the same detector signals being interpreted as different objects, an overlap removal procedure is applied as follows.If any object is rejected at one step it is not considered in later steps.First, if any two electrons share a track the electron with the lower  T is removed.Next, any -leptons within Δ = 0.2 of an electron or muon are removed.Then, any electrons which share a track with a muon are removed.If any small- jet is within Δ = 0.2 of an electron it is removed, and then any electron within a cone of  T -dependent size around a small- jet is removed.If any small- jet with fewer than three tracks has an associated muon or is within Δ = 0.2 of one it is removed, and then any muons within a cone of  T -dependent size around a small- jet are removed.Next, any small- jets within Δ = 0.2 of a -lepton are removed.Finally, any large- jets within Δ = 1.0 of an electron are removed.
Track-jets do not participate in the overlap removal as they are only used for -tagging.

Missing transverse momentum
The missing transverse momentum (with magnitude  miss T [104]) is defined as the negative vector sum of the transverse momenta of all the observable objects in the event, plus a soft term including ID tracks matched to the primary vertex but not to any of the other objects.The  miss T reconstruction uses the baseline electrons and muons as well as all small- jets, and employs a separate overlap removal procedure which takes into account detector signals from each object included [104].In control regions a modified definition of  miss T is used in which electrons and muons are treated as invisible,  miss T, lep.invis., to imitate the kinematics of the  →  ν background process.
Object mismeasurements, especially of jets, are the main source of fake  miss T .Therefore, the  miss T significance (S) [35] is defined to assess the likelihood that the  miss T is really due to invisible particles or is more likely to come from mismeasurements.It is calculated using the expected resolutions of all objects which enter the  miss T calculation and the correlations between them.

Event selection
The basic target final-state topology is a Higgs boson decaying into two -quarks produced with a significant imbalance in the measured transverse momentum.Events are divided into non-overlapping regions designed either to be enriched in the signal process (signal regions) or in a significant background process (control regions).Control regions differ from signal regions primarily through requiring the presence of one or two lepton(s), whereas signal regions veto events containing baseline leptons.
As the angle between the two -jets produced in the Higgs boson decay is inversely proportional to the  T of the Higgs boson, in cases where the Higgs boson is significantly boosted it can become difficult to reconstruct the two -quarks as separate jets.This motivates splitting the analysis into 'resolved' regions in which the decay products of the Higgs boson are reconstructed as two separate jets and 'merged' regions in which the entire Higgs boson decay is reconstructed as a single jet.
In -associated production within the 2HDM+ benchmark model, the Higgs boson and DM particles are produced with an extra pair of -quarks from gluon splitting.Therefore, to enhance sensitivity to these models, all regions are further split into those requiring exactly two -jets and those requiring ≥3 -jets (referred to as 2 -tag and ≥3 -tag, respectively).

Common selections
Events which do not have a reconstructed primary vertex are rejected.Events are also rejected if found to contain any jets with properties consistent with beam-induced backgrounds, cosmic-ray showers or noisy calorimeter cells [105].
Events are required to have  miss T > 150 GeV and are vetoed if they contain a baseline -lepton.In order to further reduce the background from -lepton decays, events are also vetoed if they have any small- jets with Δ(jet,  miss T ) < 22.5°where the track multiplicity in the jet is between 1 and 4.This selection is referred to as the 'extended -lepton veto'.A further source of background is  miss T arising from either leptonic heavy-flavour decays in a jet or a jet which is severely mismeasured.In these cases, the  miss T tends to be aligned with the jet; therefore, events where any of the up to three leading small- jets have Δ(jet,  miss T ) < 20°are rejected.For control regions, these requirements use  miss T, lep.invis.rather than  miss T .Only loose selections are placed on the mass of the Higgs boson candidate ( ℎ , defined in the following sections) because it is used as the discriminating variable for the final fit.The range is 50 GeV <  ℎ < 280 GeV for the resolved regions and 50 GeV <  ℎ < 270 GeV for the merged regions, with the lower limit chosen to be the lowest calibrated large- jet mass and the upper limit chosen to be significantly larger than the Higgs boson mass, with the precise value being determined by the  ℎ binning used in the fit, which depends on the available sample size.

Signal regions
The signal region selections for both the merged and resolved regions are summarised in Table 1.All signal region events are required to have passed the primary  miss T trigger [44].In order to reduce the contribution of SM processes producing  miss T through the decay  → ℓ, events are rejected if they contain any baseline electron or muon.

Resolved regions
The resolved regions are defined by selecting events with  miss T < 500 GeV.Events in the resolved regions are required to have at least two -tagged small- jets, with the two with the highest  T forming the Higgs boson candidate.The combined  T of this two-jet system ( T ℎ ) is required to be greater than 100 GeV and its mass is corrected for nearby muons as described in Section 4 to form  ℎ .
The dominant background in the resolved region is  t production where one top quark decays leptonically, but the lepton is either not reconstructed or not correctly identified.In these cases, all the  miss T in the event (beyond that from mismeasurement) originates from the decay of one of the two  bosons, and therefore the transverse mass of the  miss T and the corresponding -jet should be approximately bounded from above by the top-quark mass.Here the transverse mass is defined as where  ,min In order to suppress contributions from multĳet backgrounds, the object-based  miss T significance is required to satisfy S > 12.After this selection, data-driven estimates of the remaining multĳet contribution to the signal regions were found to be substantially smaller than the expected statistical uncertainty of the data so the impact of multĳet processes is not included in the background estimation.The studied signal models typically have fewer reconstructed jets than the dominant backgrounds.Therefore, events with exactly two -tagged jets (2 -tag) are required to have at most four small- jets, where only central jets are counted.
For events with at least three -tagged jets (≥3 -tag), this requirement is relaxed to at most five small- jets to ensure a sufficient sample size in the corresponding control regions.

Merged regions
The merged regions are defined by selecting events with  miss T > 500 GeV.At least one large- jet is required, and the two leading variable- track-jets associated with the leading large- jet are required to be -tagged.This large- jet is defined to be the Higgs boson candidate and its mass is corrected for nearby muons as described in Section 4 to form  ℎ .Events are separated into those which have no additional -tagged variable- track-jets (2 -tag selection) and those which have at least one such -tagged variable- track-jet not associated with the Higgs boson candidate (≥3 -tag selection).
The merged 2 -tag selection is split into two  miss T bins, 500 GeV <  miss T < 750 GeV and  miss T > 750 GeV, while no further splitting is done for the ≥3 -tag selection.

Background modelling and control regions
The dominant backgrounds in the signal regions consist of  t and / bosons produced in association with heavy-flavour jets.The /+jets backgrounds are subdivided according to the true flavour of the jets that constitute the Higgs boson candidate. 3If the flavour of either (or both) those two jets is a -quark, the event is considered to be /+HF background, where HF stands for Heavy Flavour.In the 2 -tag resolved regions the dominant backgrounds are  t and +HF, with the latter becoming more important as the  miss T increases.The 2 -tag merged regions are dominated by +HF.Both the resolved and merged ≥3 -tag regions are dominated by  t, where the extra -jet typically is a mis-tagged jet originating from a hadronic  boson decay.At higher  miss T values the +HF background becomes important again.
The  t, +HF and +HF contributions are modelled using simulation with their normalisations corrected from data by using background-enriched control regions besides the signal regions.Smaller backgrounds are taken directly from simulation.These include the production of a  or  boson in association with light jets or at most one jet containing a -hadron, single top-quark production (dominated by production in association with a  boson) and diboson processes.Small contributions also arise from  t processes in association with vector bosons or a Higgs boson.Another background contribution stems from vector-boson production in association with a Higgs boson ( ℎ), which mimics the signal due to the presence of a Higgs boson peak in association with jets.In the case where the vector boson is a  boson decaying to neutrinos this is an irreducible background.Similarly, the diboson decay   →  is nearly irreducible due to small difference in  and ℎ mass peaks (compared to the Higgs candidate mass resolution).The multĳet background is negligible in all regions after the requirements on the object-based  miss T significance  and on Δ(jet,  miss T ), and thus not further considered.The total background estimates in all regions, including uncertainties, are determined in a simultaneous fit to all regions, which is described in Section 6.
Top-quark pair production and +HF processes contribute in the signal regions if leptons in the decays are either not identified or outside the kinematic acceptance.The main contribution arises from decays involving hadronically decaying -leptons.As the shape of the event variables is the same for all lepton flavours, the kinematic phase space of the signal region can be closely approximated by control regions requiring an isolated signal muon (1-muon control regions).In this case, to better approximate the signal regions which veto the presence of any leptons,  miss T, lep.invis. is used as proxy for  miss T .Also, any other variable using  miss T in its calculation, e.g. miss T significance and  b,min/max T , is constructed using  miss T, lep.invis. .This ensures that the  miss T -related quantities in the control regions correspond to those in the signal regions.Otherwise, the 1-muon control regions are defined by the same criteria as the signal regions.
As the momentum of the  boson does not depend on its decay mode, the main background in the signal regions,  →  in association with heavy-flavour jets, can be closely modelled by  → ℓ + ℓ − events.This means that the normalisation of  → +HF contribution can be corrected by measuring  → ℓ + ℓ − +HF events.To select these events, control regions requiring exactly two baseline electrons or muons with opposite charge are defined (2-lepton control regions).These events are collected using primary triggers selecting an isolated electron or muon [106,107].Further, one of the electrons or muons is required to be a signal electron or muon with  T > 27 GeV or  T > 25 GeV, respectively.The invariant mass of the leptons is required to be consistent with the mass of the  boson within 10 GeV.While keeping all other criteria of the signal regions, an additional criterion of  < 5 is imposed to suppress a remaining contribution of  t processes.To be similar to the signal regions,  miss T, lep.invis. is used as proxy for  miss T and in the calculation of any other variable using  miss T .

Statistical analysis
A binned profile likelihood fit [108,109] is used to obtain background estimates and check the compatibility of the data with the background-only hypothesis as well as to extract upper limits at 95% confidence level (CL) on the signal cross-section.The likelihood function is constructed from a product of Poisson probability functions based on the expected signal and background yields in every region considered in the fit, as shown in Table 2.It contains the parameter of interest, , which multiplies the signal cross-section, as well as floating normalisation factors controlling the background normalisations.Systematic uncertainties are included in the likelihood by nuisance parameters (NP)  which are parameterised by Gaussian or log-normal priors.
Four normalisation factors are used for the backgrounds, scaling the  t and /+HF processes.Of those, two normalisation factors are used for  boson production in association with two -tagged jets or at least three -tagged jets.In the 2 -tag and ≥3 -tag selections,  t and +HF are normalised by a single parameter each, as the production mechanism for  t is the same in the two categories and the contribution of +HF is minor in the ≥3 -tag selection.
The likelihood includes all control and signal regions, as detailed in Table 2, where these regions are binned to increase the sensitivity and improve the determination of the normalisation factors for  t and /+HF.The signal regions are binned in the invariant mass of the Higgs boson candidate as defined in Section 5.The chosen binning in the signal regions is optimised to obtain the best expected sensitivity to the signal models addressed in this paper, while also keeping statistical uncertainties in each bin low.
The 1-muon control regions are split in positive or negative muon charge in the case of the 2 -tag regions.This improves the separation between  t and +HF processes, as +HF processes exhibit a charge asymmetry in   collisions.Only inclusive event yields are used in the other regions, which are dominated by  t, and for the 2-lepton control regions.
The nominal fit results are obtained by maximising the likelihood function with respect to all parameters.Two different fit configurations are used: in the background-only profile likelihood fit, the parameters are determined in a fit to data assuming the presence of no signal.The second fit configuration instead allows for the presence of a specific signal, and is referred to as the model-dependent fit.
The test statistic   is constructed using the profile likelihood [110]:   = −2 ln(L (, θ )/L ( μ, θ)), where μ and θ are the parameters that maximise the likelihood, and θ are the nuisance parameter values that maximise the likelihood for a given .This test statistic is used to measure the compatibility of the background-only model with the observed data and to derive exclusion intervals using the CL prescription [111].

Systematic uncertainties
Signal and background expectations are subject to statistical, detector-related and theoretical uncertainties, which are all included in the likelihood as nuisance parameters.Detector-related and theoretical uncertainties may affect the overall normalisation and/or shape of the simulated background and signal event distributions.
Detector-related uncertainties are dominated by contributions from the jet reconstruction.Uncertainties in the jet energy scale (JES) for small- jets [98] arise from the calibration of the scale of the jet and are derived as function of the jet  T , and also .Further contributions emerge from the jet flavour composition and the pile-up conditions.The 'category reduction' scheme as described in Ref. [98] with 29 nuisance parameters is used.Uncertainties in the jet energy resolution (JER) depend on the jet  T and  and arise both from the method used to derive the jet resolution and from the difference between simulation and data [98], and they are included with eight nuisance parameters.Similarly, uncertainties in the jet energy resolution for large- jets arise from the calibration, the flavour composition and the topology dependence [112].Further uncertainties are considered for the large- jet mass scale [112] and resolution [113].
Uncertainties due to the -tagging efficiency for heavy-flavour jets, including -flavour jets, are derived from  t data [114,115] and are represented by four nuisance parameters.Uncertainties are also considered for mistakenly -tagging a light-flavour jet, with nine nuisance parameters.These are estimated using a method similar to that in Ref. [116].
Uncertainties in the modelling of  miss T are evaluated by considering the uncertainties affecting the jets included in the calculation and the uncertainties in soft term's scale and resolution [117].The pile-up in simulation is matched to the conditions in data by a reweighting factor.An uncertainty of 4% is assigned to this reweighting factor.The uncertainty in the combined 2015-2018 integrated luminosity is 1.7% [41], obtained using the LUCID-2 detector [42] for the primary luminosity measurements.
Scale factors, including their uncertainties, are calculated specifically for this analysis to correct the efficiency of  miss T triggers in simulation to that in data.The uncertainties in the scale factors are at most 1%-2% for low  miss T values.
In the regions requiring the presence of leptons, uncertainties in the lepton identification and lepton energy/momentum scale and resolution are included.These are derived using simulated and measured events with  → ℓ + ℓ − , / → ℓ + ℓ − and  → ℓ decays [100,102].
Modelling uncertainties impact the shape of the  ℎ distribution, the relative acceptance between different  miss T and -tag multiplicity bins and between signal and control regions as well as the overall normalisation of the samples that are not freely floating in the fit.
The theoretical uncertainties are dominated by modelling uncertainties in the  t and +HF backgrounds.For the  t and  processes, the impact of the choice of parton shower and hadronisation model is evaluated by comparing the sample from the nominal generator set-up with a sample interfaced to H 7.04 [118,119].To assess the uncertainty in the matching of NLO matrix elements to the parton shower, the P B sample is compared with a sample of events generated with M G 5_aMC@NLO v2.6.2.For the  process, the nominal sample is compared with an alternative sample generated using the diagram subtraction scheme [77,84] instead of the diagram removal scheme to estimate the uncertainty arising from the interference with  t production.
For the +jet processes, uncertainties arising from the modelling of the parton shower and the matching scheme are evaluated by comparing the nominal samples with samples generated with M -G 5_aMC@NLO v2.2.2.For the diboson processes, the uncertainties associated with the modelling of the parton shower, the hadronisation and the underlying event are derived using alternative samples generated with P B [60][61][62] and interfaced to P 8.186 [120] or H ++.
For all MC samples, the uncertainties due to missing higher orders are estimated by a variation of the renormalisation and factorisation scales by a factor of two, while the PDF and  s uncertainties are calculated using the PDF4LHC prescription [121].
Table 3 gives the impact of the different sources of systematic uncertainties for selected signal models as evaluated in different model-dependent fits.The signal models with lower masses illustrate the impact of the systematic uncertainties in the resolved regions, while the models with larger mediator masses are more impacted by the merged regions.The theoretical uncertainties in the modelling of the  t background, the experimental uncertainties in the calibration of jets and the limited MC sample size show the largest impact.

Results
The post-fit background yields are determined in a background-only profile likelihood fit to data in all regions.Figure 3 shows the yields in the 1-muon and 2-lepton control regions.The post-fit normalisation factors for  t and for +HF are found to be 0.93±0.08 and 0.95±0.14,respectively.For  boson production in association with two (at least three) heavy-flavour jets the normalisation factors are determined to 1.41 ± 0.09 (1.85 ± 0.24).An upward scaling of the +HF background relative to the simulation was also   Table 3: Relative importance of the different sources of uncertainty for different  -2HDMs, with the masses of the  boson and the  boson given in the second row, expressed as fractional impact on the signal strength parameter.The fractional impact is calculated by considering the square of the uncertainty in the signal strength parameter arising from a given group of uncertainties (as listed in the left column of the table), divided by the square of the total uncertainty in the signal strength parameter.Due to correlations, the sum of the different impacts of systematic uncertainties might not add up to the total impact of all systematic uncertainties.

Source of uncertainty
Fractional observed in other studies [122], and was attributed to an underestimation of the  →  b rate in S .A larger scaling is observed in the region with ≥3 -tagged jets, dominated by processes with more  →  b splittings.The uncertainty on the +HF normalisation factor increases in the ≥3 -tag region due to lower statistical precision and the smaller contribution of the +HF background.
The distributions of the Higgs boson candidate mass  ℎ after the background-only fit are shown in Figures 4  and 5.The signal to background ratio is higher in the ≥3 -tag signal regions, because the 2HDM+ signal model is shown with tan  = 10, where -associated production dominates.Tables 4 and 5 present the background estimates in comparison with the observed data.No significant deviation from SM expectations is observed, with the largest deficit corresponding to a local significance of 2.3, and the largest excess amounting to 1.6.Figure 6 summarises the total yields in the signal regions as a function of  miss T .The background prediction from simulation is scaled upwards in the fit for lower  miss T values, while the simulation agrees better with the data for large  miss T values.The results are interpreted as exclusion limits at 95% CL in the  -2HDM and the 2HDM+ scenarios in Figures 7 and 8. Considering the  -2HDM case,  masses up to 3 TeV are excluded for  masses of 300 GeV at 95% CL.The exclusion boundaries for the 2HDM+ scenario extend up to   = 520 GeV for   = 1.25 TeV for F production and tan  = 1.This is an improvement of about 200 GeV in   on previous results [123], which reinterpreted the earlier ℎ(→  b) +  miss T analysis using 36.1 fb −1 [27].
The higher exclusion limit at high   , low   , is due to an increase of the cross-section of the  → ℎ process, without resonant  production.It should be noted that with the exact parameter choices adopted in this analysis, the ℎ coupling becomes larger than 4 for   1750 GeV.Moreover, as discussed in Refs.[25,26] values of   1250 GeV (for tan  = 1) or   2150 GeV (for tan  = 10) would not be consistent with the requirement of having a bounded-from-below scalar potential, given the parameter choices discussed in this paper.These constraints can be relaxed substantially if the quartic couplings assume a value closer to the perturbativity limit and also in more general 2HDMs containing additional couplings as discussed in Refs.[124,125].Therefore, the above should not be considered as a strong requirement for the validity of the model predictions.At high   the width of the additional Higgs bosons grows substantially and the theoretical predictions are subject to additional theoretical uncertainties associated with the treatment of the width.Exclusion limits are therefore not shown in the region of very large widths (  > 2200 GeV).
In the case of   production and tan  = 10, the exclusion limits extend up to   = 240 GeV for   = 900 GeV.The inclusion of the ≥3 -tag tag region helps to increase the sensitivity relative to the 2 -tag tag region by about 30-70%.The difference between observed and expected limits arises from data deficits in the ≥3 -tag region, especially the deficit around the Higgs boson peak in the  miss T ∈ [350, 500) GeV region, as shown in Figure 5.The 2HDM+ scenario with   production and tan  = 10 is considered for the signatures discussed in this paper for the first time.
Figure 9 displays the upper limits on the visible cross-section, defined as where (A × ) with the acceptance A and the reconstruction efficiency  quantifies the probability for a certain event to be reconstructed within a window around the Higgs boson mass in a given signal region.The visible cross-section is obtained from the number of signal events in each signal region extracted from a fit to the ( b) distribution as described below, divided by the integrated luminosity.
In contrast to the model-specific exclusion limits in Figures 7 and 8, the upper limits on the visible cross-section are calculated without dependence on a signal model.They only assume that a resonance 500 1000 1500 2000 2500 3000 3500 [GeV]  Figure 8: Exclusion limits for the 2HDM+ signal with (a) tan  = 1 and F production and with (b) tan  = 10 and   production.The solid black line shows the observed limit at 95% CL, the dashed black line the expected limit.The green band gives the ±1 uncertainties of the expected limit, the yellow band the ±2 uncertainties.The hashed area in plot (a) indicates the region where the width of at least one of the Higgs bosons ,  or  ± or of the pseudoscalar  is above 20% of its mass.Extending the exclusion into this region would require additional assumptions within the 2HDM+ model, and as such limits are not shown beyond 2.2 TeV. was produced with a mass close to 125 GeV and decays into a pair of  quarks in association with  miss T .For this purpose, the binning in  ℎ in the signal regions is modified such that all bins under the Higgs boson peak in a range from 90 GeV to 150 GeV are merged into one bin.This bin includes the Higgs boson peak, but excludes important parts of the  boson peak.A simultaneous fit to all control and signal regions is performed, using the modified binning.However, in order not to assume a specific signal model, a signal may only be present in the one bin under the Higgs boson peak.Although all signal regions are fitted simultaneously, the signal contributions in the different signal regions are independent of each other, which is ensured by using different signal strength parameters.

Conclusion
A search for dark matter production in association with a Higgs boson decaying into  b is presented, considering final states with either two or at least three -tagged jets and without leptons.The search uses proton-proton collision data recorded by the ATLAS experiment at the LHC during the data-taking periods in 2015-2018, corresponding to an integrated luminosity of 139 fb −1 .The analysis targets the 2HDM extended by dark matter particles and either a heavy vector boson  ( -2HDM) or a pseudoscalar singlet  (2HDM+) as mediator particle.
Sensitivity to these models is obtained by considering different -jet multiplicity regions, and the cases of boosted or non-boosted Higgs bosons.Further improvements relative to previous searches targeting this signature with a partial dataset from Run 2 of the LHC are obtained through better background rejection of  t processes, better identification of -hadrons from the decay of a boosted Higgs boson, and improved suppression of multĳet processes by using an object-based  miss T significance.
No significant deviation from Standard Model expectations was found.Exclusion limits are set on the  -2HDM and extend up to  -masses of 3 TeV at 95% CL.This is an improvement of about 700 GeV relative to the analysis in Ref. [27] for heavy  mediator masses.In the case of the 2DHM+, the mass of the pseudoscalar  is excluded up to 520 GeV for tan  = 1 and gluon-gluon fusion production, and up to 300 GeV for tan  = 10 in -associated production.Upper limits on the visible cross-section range from 0.05 to 3.26 fb, depending on the signal region.[4]

Figure 1 :
Figure 1: Feynman diagram for the production of the mono-Higgs signature in the  -2HDM.

Figure 2 :
Figure 2: Feynman diagrams for the main production mechanisms of the mono-Higgs signature in the 2HDM+ scenario: (a) gluon-gluon fusion and (b) -associated production.

T
and  ,max T are respectively defined as the  T of the -jet which is closest to (min) or furthest from (max) the  miss T in .Events are required to satisfy  ,min T > 170 GeV and  ,max T > 200 GeV.

The 2
-tag and ≥3 -tag selections are split into three  miss T bins: 150 GeV <  miss T < 200 GeV, 200 GeV <  miss T < 350 GeV and 350 GeV <  miss T < 500 GeV.In the highest  miss T bin the requirement on  T ℎ is tightened to being greater than 300 GeV.This leads to six resolved signal regions: one corresponding to each combination of  miss T bin and number of -tagged jets.The  miss T triggers become fully efficient at an offline  miss T value close to 200 GeV, but the analysis also uses events in the range 150 GeV <  miss T < 200 GeV.In order to correct for Monte Carlo (MC) mismodelling of the  miss T trigger response, the trigger efficiency must be measured in both data and simulation and scale factors calculated to correct the simulation.Given that the  miss T triggers use calorimeter information only, muons are treated as almost invisible particles, meaning that  miss T trigger efficiencies can be measured using events selected by single-muon triggers [106].The scale factors are calculated as a function of  miss T, lep, invis in a region whose selection matches the 150 GeV <  miss T < 200 GeV and 2 -tag region except that all  miss T selections are dropped, exactly one -tagged jet is required and exactly one signal muon is required.Events containing electrons are still vetoed.The scale factors have values in the range 0.95-1.0.

Figure 3 :
Figure 3: Yields in the resolved and merged (a) 1-muon control regions and (b) 2-lepton control regions.The top panel compares the fitted background yields with data, while the bottom panel indicates the ratio of the observed data to the predicted Standard Model backgrounds.The different control region bins included in the fit are indicated on the -axis by first giving the range in  miss T and then the sign of the muon charge (where applicable).

Figure 4 :
Figure 4: Distributions of the Higgs boson candidate mass in the 2 -tag signal regions for different  miss T ranges.The top panel comparies the fitted background yields with data, while the bottom panel indicates the ratio of the observed data to the predicted Standard Model backgrounds.The background yields are obtained in a background-only fit to data.An example signal model from the  -2HDM with parameters (  ,   ) = (1400 GeV, 1000 GeV) and a cross-section of 1.89 fb is displayed for comparison.The signal is scaled by the factors indicated in the legend for better visibility.

Figure 5 :
Figure 5: Distributions of the Higgs boson candidate mass in the ≥3 -tag signal regions.The top panel compares the fitted background yields with data, while the bottom panel indicates the ratio of the observed data to the predicted Standard Model backgrounds.The background yields are obtained in a background-only fit to data.An example signal model from the 2HDM+ with   production and with parameters (  ,   ) = (1000 GeV, 150 GeV), tan  = 10, and a cross-section of 62.7 fb is displayed for comparison.The signal is scaled by the factors indicated in the legend for better visibility.

Figure 6 :
Figure 6:  miss T distributions after requiring (a) two -tagged jets or (b) at least three -tagged jets.The background yields are obtained in a background-only fit to data.Each bin corresponds to one signal region, where the bins at lower  miss T use small- jets and the ones with larger  miss T select events with at least one large- jet.The rightmost bin includes all events with  miss T above the range shown on the plots.The bottom panel compares the data with background estimates as obtained from the fit.Example signal models are overlaid in the top panel.In plot (a) two  -2HDM points are shown, one with parameters (  ,   ) = (1400 GeV, 1000 GeV) and a cross-section of 1.89 fb and the other with parameters (  ,   ) = (2800 GeV, 300 GeV) and a cross-section of 1.14 fb.Both signals are scaled up by a factor of 5 for better visibility.In plot (b) a 2HDM+ is shown with parameters (  ,   ) = (1000 GeV, 150 GeV), tan  = 10 and a cross-section of 62.7 fb.

Figure 7 :
Figure7: Exclusion limits in the  -2HDM model.The solid black line shows the observed limit at 95% CL, the dashed black line the expected limit.The green band gives the ±1 uncertainties of the expected limit, the yellow band the ±2 uncertainties.

Table 1 :
Summary of selections used to define the signal regions used in the analysis.The kinematic variables are defined in the text.-tagged small- jets At least 2 -tagged associated variable- track-jets  T ℎ > 100 GeV if  miss

Table 2 :
Event categories entering the combined fit of the model to the data.The discriminant  ℎ denotes the mass of the light Higgs boson candidate and corresponds either to the dĳet mass    in the regions selecting two small- jets or to the large- jet mass   for the regions requiring a large- jet."Yields" refers to the number of events in a given region.

Table 4 :
Background yields in comparison with data in the 2 -tag signal regions for different  miss