Collider constraints on electroweakinos in the presence of a light gravitino

Using the GAMBIT global ﬁtting framework, we constrain the MSSM with an eV-scale gravitino as the lightest supersymmetric particle, and the six elec-troweakinos (neutralinos and charginos) as the only other light new states. We combine 15 ATLAS and 12 CMS searches at 13 TeV, along with a large collection of ATLAS and CMS measurements of Standard Model signatures. This model, which we refer to as the ˜ G - EWMSSM, exhibits quite varied collider phenomenology due to its many permitted electroweakino production processes and decay modes. Characteristic ˜ G -EWMSSM signal events have two or more Standard Model bosons and missing energy due to the escaping gravitinos. While much of the ˜ G -EWMSSM parameter space is excluded,


Introduction
Although supersymmetry (SUSY) was not invented to address shortcomings of the Standard Model (SM) of particle physics or cosmology, it addresses them in various aspects.Inflation, dark matter, the cosmic matterantimatter asymmetry, neutrino masses, patterns of fermion families, gauge and Yukawa couplings, naturalness, and more, can all be accommodated if supersymmetry is a symmetry of nature that is broken near the TeV scale; see for example Refs.[1][2][3][4][5] for reviews.Consequently, a major goal of the Large Hadron Collider (LHC) is to search for superpartners.So far, the LHC experiments have found no concrete evidence for SUSY and the impact of the null results in simple SUSY scenarios has been well explored (see e.g. the global fits in Refs.[6][7][8][9][10][11][12][13][14][15][16]).For example, in our previous work [17], we investigated the collider constraints on the electroweakino sector of the Minimal Supersymmetric Standard Model (MSSM).Gravitinos, however, are an interesting and often ignored possibility in SUSY collider phenomenology.
Motivated by GMSB, in this work we consider the electroweakino sector and an approximately massless gravitino LSP, with the other superpartners decoupled.The next-to-lightest supersymmetric particle (NLSP) must then be a neutralino or a chargino, though the latter is unusual in the MSSM parameter space [46,47].The electroweakinos, χ0 1,2,3,4 and χ± 1,2 , may decay to a gravitino and an SM particle.Naively, one might expect this to proceed slowly through gravitational interactions.However, as the gravitino acquires goldstino interactions through the super-Higgs mechanism [48,49], the decay may be prompt when m 3/2 1 keV [50].The neutralino decays χ0 → {h, H, A, Z} G and the chargino decays χ± → {H ± , W ± } G could be kinematically allowed depending on the mass spectrum, whereas the neutralino decays χ0 → γ G are guaranteed to be allowed and dominate for the lightest neutralino, χ0 1 , across much of parameter space [51].We thus assume that the electroweakinos may decay promptly through any kinematically open channel to an SM particle and a gravitino.
Direct LHC production of gravitino pairs, or associated production of a gravitino and another superpartner, can only reach detectable rates if m 3/2 ≪ 1 eV [52,53].For scenarios with electroweakinos within LHC reach and an eV scale gravitino, which is the focus of our study, the dominant gravitino production mode is through the prompt decay of the NLSP.This gives rise to distinctive collider signatures, such as two gravitinos that carry away missing energy and two energetic photons.Whilst the NLSP always decays promptly to a gravitino, an eV scale gravitino implies that the heavier electroweakinos decay predominantly to lighter ones [50], unless the mass degeneracy between the electroweakinos is severe (see below).Production of heavier electroweakinos will therefore typically result in multi-step decay chains that terminate with the decay of the NLSP to the gravitino.
The phenomenological impacts of electron-positron collider [51,[54][55][56][57][58][59][60], Tevatron [50,51,53,[61][62][63][64][65][66][67] and LHC [52,[68][69][70][71][72][73][74][75][76] searches on these scenarios have been previously studied.Reference [68], for example, establishes limits on the electroweakino sector using light gravitino pair-production via electroweakino decay in the context of GMSB in the MSSM.This study shows that while LHC searches specifically designed for such scenarios are important, other LHC searches and measurements provide useful complementary constraints.Using the GAMBIT software [77,78], we here go be-yond previous works by performing the first global fit of electroweakinos in the presence of a light gravitino.We include up-to-date results from LHC Run 2, described in Sec.3.1, and for the first time in a global fit we check that our models are allowed by a suite of measurements of SM-like final states using Contur [79,80]; see Sec. 3.2 for further details.Lastly, we include constraints from the Large Electron-Positron collider (LEP); see Sec. 3.3.We do not include Tevatron searches as these constraints are in general superseded by LHC results, and performing event simulations for Tevatron searches in addition to LHC searches would greatly increase the computational expense of our study.
Whilst a gravitino LSP could play the role of dark matter (DM), and there are strong constraints that we do not consider [81][82][83][84], each of these requires some additional assumptions.It was originally thought that to avoid over-closing the Universe it must be that m 3/2 1 keV [85].Although this constraint is weakened when one considers inflation [86,87], non-thermal production of gravitinos and the NLSP decays to gravitinos are both constrained by the measured abundance of DM.There are, furthermore, constraints from cosmic structure [88] and big-bang nucleosynthesis [89,90], however, the latter does not apply to our scenario where the NLSP decays promptly.We choose not to include constraints from the dark matter properties of the gravitino in this work, in order to explore electroweakinos more generally without making any limiting assumptions about cosmology.
A recent motivation for studying the possibility of light electroweakinos in this scenario is the surprising result from the CDF measurement of the W boson mass [91], which gives a value considerably above both the SM prediction and above existing experimental results.See Ref. [92] for a review of the SM value and a summary of the experimental status.Light electroweakinos, in particular light winos and Higgsinos, are known to result in significant positive corrections to the W mass [93][94][95].However, given the current uncertainty about the interpretation of the new result and its compatibility with other recent measurements, e.g.Ref. [96], we will not use this as a constraint on our model.

Model
The model under consideration in this study is a variant of the MSSM where all supersymmetric states except the electroweakinos and a quasi-massless gravitino are decoupled.This model, henceforth G-EWMSSM, differs from the model in our previous study [17] by the addition of the light gravitino.As discussed in the introduction, a very light gravitino can be motivated in certain supersymmetry breaking scenarios, e.g.gauge mediation.
The general neutralino can be any linear combination of the neutral gauginos ( B, W 0 ), and the neutral Higgsinos ( H0 u , H0 d ), where N ij are the mass eigenvectors indicating the weight of each field component in the gauge basis, (ψ 0 ) T = ( B, W 0 , H0 d , H0 u ).The corresponding bilinear terms in the Lagrangian density are where the neutralino mass matrix, M N , is given by and M 1 , M 2 and µ are the gaugino and Higgsino softbreaking bilinear couplings, respectively, which are free parameters in our model.Further, we have s β = sin β and c β = cos β, g and g are the SU (2) L and U (1) Y gauge couplings, and v is the electroweak VEV.Amongst these, only the ratio tan β = v u /v d is not fixed by data and remains an additional free parameter in our model.The general chargino eigenstates correspond to the charged Higgsinos ( H+ u , H− d ), and gauginos ( W + , W − ).The corresponding bilinear terms in the Lagrangian density are where the chargino mass matrix, M C , is given by The gravitino mass m 3/2 depends on the dynamics of the supersymmetry breaking, but for the purpose of our study we fix it to m 3/2 = 1 eV, similar to what is commonly assumed in ATLAS and CMS searches, see for example Ref. [97].In terms of the collider phenomenology, this makes the gravitino effectively massless and ensures prompt decays of the NLSP.We do not set the mass to exactly zero since the limit m 3/2 → 0 corresponds to no supersymmetry breaking.The exact choice for the small gravitino mass has very little impact on the results as long as m 3/2 = 0.The one small exception is for a wino-like chargino around the W mass or lower, where the gravitino mass may dictate whether the chargino decays directly to the gravitino or via the neutralino NLSP.However, scenarios with such light Fig. 1: Branching ratios for the lightest neutralino as a function of µ, with M1 = 900 GeV, M2 = 400 GeV, and tan β = 1 (top) or tan β = 10 (bottom).The wino and Higgsino NLSP regions are shown in red and green, respectively.The pink line (dash dot) shows the combined branching ratio for decays to all states other than on-shell (Z, h, γ) + G.The thin, grey bar marks a parameter region where m χ± charginos are in any case heavily constrained, independent of this decay.Since we do not consider direct production of gravitinos, where the cross section would be low and the signature difficult to disentangle from backgrounds, the LHC phenomenology of this model is dominated by the production and decay of the light electroweakinos.The hierarchy of M 1 , M 2 and µ, and to some extent the value of tan β, determines their gaugino and Higgsino components, production cross sections and branching ratios.
A chargino NLSP will decay promptly to the gravitino and a (possibly off-shell) W boson.However, having a chargino NLSP is only possible in narrow regions of parameter space; see Fig. 1 for an example.
Throughout most of parameter space the lightest neutralino is the NLSP.In general, a neutralino NLSP has three possible decay modes: χ0 1 → {γ, Z, h} G.In the limit, m 3/2 m { χ,Z,h} , the decay widths take the form [82,98]: Here s W , c W , s α and c α are the sines and cosines of the Weinberg angle θ W and the mixing angle α between the CP -even neutral Higgs states, and .
In Fig. 1 we show representative branching ratios for the lightest neutralino, using the full expression for the widths from Refs.[82,98,99], including also decay modes through off-shell bosons in the total width.The plots use values of µ picked to illustrate the generic behaviour in the different wino NLSP (red) and Higgsino NLSP (green) regions (see below for further discussion), and two different values of tan β, which cover the impact of tan β on decays to Z and h.The bino NLSP region (low M 1 values) is much simpler and not illustrated since here dominantly χ0 1 → γ G, with some small branching ratio to Z G.We see that the dominant decay mode of the lightest neutralino depends strongly on the relative ordering of the masses M 1 , M 2 , and µ, and the size of tan β.
To make our presentation more systematic, we now discuss the properties of these three major phenomenological regions in terms of the ordering of the gaugino, M 1 and M 2 , and Higgsino, µ, masses.
Wino NLSP: With |M 2 | < |M 1 |, |µ|, the two lightest electroweakinos, χ± 1 and χ0 1 , are a charged and neutral wino with relatively large LHC production cross sections.The lightest neutralino decays as χ0 1 → {Z, γ} G, see for example the wino NLSP region (red) of Fig. 1 with µ > M 2 .For the lightest chargino, when m χ± 1 m W the small mass difference between the wino-like chargino and neutralino leads to decays directly to the gravitino and an on-shell W , χ± 1 → W ± G.For smaller chargino masses we have instead decays to two fermions (via an off-shell W ), together with the gravitino or lightest neutralino χ0 1 .
1 , χ0 2 and χ± 1 , are dominantly Higgsino and have somewhat smaller production cross sections compared to the wino scenario.Pure Higgsinos do not decay to photons at tree level, so in this case the decays χ0 1 → {Z, h} G are typically dominant, unless the NLSP mass is so small that the available phase space becomes limiting, or even that these decays go off-shell.In this case decays to photons become important again, especially at low masses, along with three-body final states with two opposite-sign SM fermions at intermediate masses.The relationship between the branching ratios to Higgs and Z final states is determined by the sign of µ and the value of tan β.In particular we note that taking µ < 0 and tan β → 1 suppresses the Z G channel, due to cancellation between the N 13 and N 14 terms in Eq. ( 6).This interplay of decays is again illustrated in Fig. 1 in the Higgsino NLSP region (green) with |µ| < M 2 .The heavier neutralino and the chargino typically decay to the lightest neutralino and SM fermions in three-body decays, instead of the gravitino, due to the generically larger mass differences between the lightest electroweakinos in the Higgsino scenario [47].
1 and the direct pair production cross section at the LHC is small.Most of the production is then likely to be from decays of the heavier, winoor Higgsino-dominated electroweakinos, depending on the hierarchy of M 2 and µ.The bino NLSP decays dominantly as χ0 1 → γ G.The overall pattern that can be deduced from the above discussion is that the model predicts events with a pair of bosons picked from {h, Z, W, γ}, along with missing energy from the escaping gravitinos, possibly with one or both bosons being off-shell if the mass of the NLSP is below 125 GeV.Additional bosons may also be produced from the decays of heavier electroweakinos into the NLSP.In addition to the classic signature of di-photons plus missing energy, we see that this model features events with final state SM fermions from the decays of the massive bosons, meaning that many LHC searches are relevant for the model.
Apart from the addition of the light gravitino LSP, our implementation of the G-EWMSSM model in GAM-BIT is identical to our implementation of the EWMSSM model described in detail in Ref. [17].In particular, the Higgs mass, which in this study only matters for event kinematics, is set by hand to 125.09 GeV.

Collider likelihoods
The total likelihood function explored in our global fit consists of likelihoods for LHC searches for new particles, LHC measurements of SM signatures, and LEP crosssection limits for electroweakino production.We describe each of these likelihoods below.

LHC searches
The likelihood contribution from LHC searches is based on passing simulated signal events through our emulations of the 13 TeV ATLAS and CMS searches in Refs. .Reproducing a collider search to sufficient accuracy can be challenging, e.g.due to limited available information about technical details of the analysis, or due to limitations in the tool-chain used for fast event simulation.In some cases we can therefore only incorporate a subset of the signal regions defined by the search.In Appendix B we provide a short description of each search, and point out which signal regions our signal simulation includes.
For all the included LHC searches we have the background uncertainty for each signal region, but in many cases there is no public information on how these uncertainties are correlated.We then take a conservative approach and, for each search, construct a likelihood function that only uses the signal region i with the best expected sensitivity for the given G-EWMSSM parameter point.Our likelihood function for each of these searches is then constructed from a simple product of a Poisson and a Gaussian factor, where s i , b i and n i are, respectively, the expected signal yield, expected background yield and observed yield for the given signal region i.The Gaussian factor with the nuisance parameter γ i is introduced to account for the uncertainty in the total predicted yield, and we therefore set the width σ i by adding in quadrature the uncertainties of s i and b i .For our parameter scans we need a likelihood function that only depends on the predicted signal yield.Thus, for each sampled G-EWMSSM parameter point we profile L 1SR search (s i , γ i ) over the nuisance parameter γ i : where γi is the γ i value that maximises L 1SR search (s i , γ i ) for a given s i .
CMS have for a number of their searches published covariance matrices for the background uncertainties, following the simplified likelihood approach [127,128].For these searches we can generalise Eq. ( 8) to a likelihood that utilises the information in all signal regions.A search with n SR signal regions is then described by the likelihood function Here Σ is the n SR × n SR covariance matrix for the nuisance parameters γ i .We construct Σ by taking the covariance matrix provided by CMS and adding in quadrature our signal yield uncertainties along the diagonal.
To obtain a likelihood that only depends on the set of signal yields s we, for each G-EWMSSM point, profile L search (s, γ) over the set of n SR nuisance parameters, We also note that for the searches in Refs.[101, 106, 109-111] ATLAS have published the information required to fully utilise all signal regions, through the full likelihood framework [129].We will make use of these likelihoods in future GAMBIT studies.
The LHC experiments often present results for multiple categories of final states in a single publication, e.g. the CMS multilepton search for charginos and neutralinos in Ref. [123], which presents results for searches in 2-lepton, 3-lepton and 4-lepton final states.In these cases we follow the same approach as in [17] and treat the results for the different final states as approximately independent searches, meaning that for each final state category we include a separate likelihood contribution of the form given in Eqs.(9) or (11). 1  Similar to the approach in Refs.[17,131,132], we normalise the likelihood function for each LHC search with the corresponding background-only (s = 0) likelihood.The log-likelihood contribution from each search therefore takes the form of a log-likelihood difference (12) Treating the searches as independent, what we consider as the combined log-likelihood from all the LHC searches is where ∆ ln L j is the log-likelihood contribution from search j.A positive value for the ∆ ln L searches (s) indicates that the combined set of G-EWMSSM signal predictions s(θ) for parameter point θ gives an overall better agreement with current LHC search results than the background-only assumption does.This happens when the predicted G-EWMSSM signals can help accommodate data excesses in some searches, without conflicting strongly with the results of the other searches.We will present the result of our global fit as profile likelihood maps in different G-EWMSSM mass planes.For each plane we show the 1σ (68.3%) and 2σ (95.4%) confidence regions, derived using the likelihood ratio L(θ)/L(θ best-fit ), where θ best-fit is the highest-likelihood G-EWMSSM parameter point.Therefore, if the bestfit point can explain some excesses in the search data 1 A new method for identifying non-overlapping combinations of signal regions from large collections of LHC searches was recently presented in Ref. [130].We plan to implement this method in GAMBIT and use it in future studies.
(∆ ln L searches (s) > 0), the G-EWMSSM parameter regions outside the 2σ contour should not be considered "excluded" in the same sense as for an exclusion limit from an LHC search.Rather, these parameter regions simply provide a significantly worse fit to the combined data compared to that of the best-fit point.It is then interesting to also ask a different question: What G-EWMSSM parameter regions are excluded by the combination of LHC searches, when judged relative to the background-only expectation?A simple way to estimate this is to replace ∆ ln L searches (s) in Eq. (13) with This log-likelihood penalises G-EWMSSM parameter points that give a joint prediction in worse agreement with data than the background-only prediction, while all other points are assigned the same log-likelihood of 0. We note that the maximum value ∆ ln L cap searches (s) = 0 can be obtained in two different ways: The first case is when none of the included searches are sensitive to the given G-EWMSSM parameter point, i.e. the limit s → 0. This is typically what happens for high-mass scenarios, due to small production cross-sections.The second case is when a G-EWMSSM scenario fits the results from some LHC searches sufficiently better than the SM does, enough to offset any likelihood penalty from tensions with other LHC analyses.In Sec. 5 we will present results both for the "full likelihood" (L searches ) case and the "capped likelihood" (L cap searches ) case.This is the same approach as was taken in Refs.[17,131,132].

LHC measurements of SM signatures
The complexity of the phenomenology of the model means that the possibility that it may produce events which could contribute to well-measured SM-like final states must also be taken into account.This is the scenario for which Contur [79,80] is designed.Via Rivet [133], Contur has access to an extensive library of measurements from the LHC experiments, mostly corrected for detector effects and thus not requiring explicit detector simulation.Simulated events are passed through Rivet and projected into the fiducial phase space of the measured cross sections.In the release of GAMBIT accompanying this paper, we have interfaced Contur and Rivet to the GAMBIT ColliderBit module.
As binned unfolding of detector effects requires statistically stable bin populations, a χ 2 test has proven indistinguishable from Poisson log-likelihood differences for measurement interpretations.The χ 2 is evaluated and used as the log-likelihood difference between the "signal-injection" hypothesis and the SM null hypothesis, in this application assuming the data to be equal to the SM: with the log-likelihood difference then evaluated as ∆ ln L meas (s) = ln L meas (s) − ln L meas (s = 0).The set of active bins is conservatively selected to avoid acceptance overlaps, as described in Sec.4.1, and y i and ∆y i are the bin values and uncertainties respectively.The experimental uncertainties are taken into account in the χ 2 construction, but are treated as uncorrelated in the version of Contur (2.3.0)used here.
The set of 13 TeV analyses used by Contur in this analysis are those described in Refs. .These cover final states with (multiple) jets, isolated photons and leptons, as well as missing energy.When discussing our results in Sec. 5 we will highlight the analyses with the greatest impact.

Cross-section limits from LEP searches
In addition to the above LHC searches and measurements that are implemented at the event level, we include LEP searches and measurements that were published as upper limits on particular electroweakino production cross-sections.See [17,182] for general details of our treatment of LEP searches.First, there are searches for electroweakinos that we applied in [17] that we re-interpret as searches for gravitinos.Specifically, we consider searches for pair production of charginos that each decay into SM particles and a stable neutralino, χ± → SM + χ.In our gravitino model, the chargino may decay into SM particles and a gravitino, χ± → SM + G.This leads to an identical signature as both the gravitino and a stable neutralino only contribute to missing energy.
Second, we include a multi-photon and missing energy search by L3 at √ s = 207 GeV [183].In our model, neutralinos can be pair produced at LEP and can each decay to a photon and a gravitino, giving a signature of missing energy and two photons.The number of observed events in the search was less than expected from SM backgrounds, leading to strong constraints on the e + e − → χ0  1 → G Gγγ cross section as a function of the gravitino and neutralino masses for masses less than about √ s/2.We apply the 95% limits shown in Fig. 6c of [183] following the treatment described in [182].The impact of this constraint on our G-EWMSSM model is limited, however, as our assumption of decoupled selectrons typically leads to a very small e + e − → χ0 We perform our study of the G-EWMSSM with the GAMBIT 2.4 global fit framework [77,184], utilising the SpecBit, DecayBit, ColliderBit and ScannerBit modules [182,185,186].To compute the chargino and neutralino mass spectrum at one-loop level, SpecBit employs a FlexibleSUSY [187,188] spectrum generator which uses SARAH [189,190] and routines from SOFT-SUSY [191,192].A more detailed discussion of this spectrum computation is given in [17].
For this study we have extended DecayBit with the capability to compute decay widths for a neutralino or chargino decaying to final states with a gravitino.The implementation is based on analytical expressions given in Refs.[82,98,99].To compute neutralino and chargino decays into final states with a lighter neutralino or chargino, DecayBit uses SUSY-HIT 1.5 [193], which includes the packages SDECAY [194] and HDECAY [195].
We simulate LHC events with electroweakino production at √ s = 13 TeV using ColliderBit's parallelised interface to Pythia 8 [196,197] and native fast detector simulator BuckFast [182]. 2 Due to the cost of computing higher-order production cross-sections, we use the cross-sections computed by Pythia 8 at leading-order plus leading-logarithmic (LO+LL) accuracy.As we will see in Sec. 5, the lowest-mass scenarios not disfavoured by current results are scenarios where the lightest electroweakinos are Higgsinos with masses around 200 GeV.For such scenarios the production cross-sections at nextto-leading order with next-to-leading-logarithmic corrections (NLO+NLL) can be up to 30% higher compared the LO+LL cross-sections [198], so this choice is somewhat conservative.
For each parameter point included in our final scan results we generate 16 million LHC events to evaluate the impact of the LHC searches.The main reason that such a high number of events is needed is that for many of the searches we do not have the information needed to allow a proper statistical combination of all the signal regions in the search.As discussed in Sec.3.1, for these searches the conservative approach is, for each sampled parameter point, to identify the signal region with the best expected sensitivity, and only use this signal region when computing the likelihood contribution from the given search.Many searches will for large parts of the G-EWMSSM parameter space have several signal regions with low and near identical expected sensitivities.Thus, the signal region choice, and through it the likelihood value, becomes highly sensitive to Monte Carlo noise. 3s a post-processing step, we generate a further 100 000 events at each sampled parameter point, which are then passed to first Rivet and then Contur using the new ColliderBit interface.This enables evaluation of whether the parameter point in question would have led to significant but unnoticed collective deviations from the SM expectation in existing measurements.Since LHC measurements have much higher acceptances than LHC searches, we here need fewer simulated events to ensure sufficiently small Monte Carlo uncertainties and a stable identification of the most sensitive measurements.Contur tests the full set of measurements for each parameter point, evaluating the expected likelihood ratio for each measurement.As is usual with Contur, to account for statistical correlations between measurements and avoid double-counting of BSM effects, these measurements are divided into non-overlapping "analysis pools" based upon the run period, experiment and final state.Only the most sensitive measurement from each pool is used, and the set of pool-likelihoods is then combined to provide an overall Contur likelihood, which in Col-liderBit is then combined with the likelihoods for the LHC searches and the LEP cross-section limits.The likelihood provided by Contur in this post-processing step had a significant impact on the final results, which will be discussed in detail in Section 5.4.

Scanning strategy
With the gravitino mass fixed at 1 eV, the collider phenomenology of our model is determined by the mass parameters M 1 , M 2 and µ, and the dimensionless tan β parameter.We restrict our attention to scenarios where the electroweakino masses are all 1 TeV.This is due to the substantial computational cost of accurately mapping out the profile likelihood function across wide,

Parameter
Range/value Sampling priors The "hybrid" prior refers to a prior that is flat on |x| < 10 GeV, and logarithmic elsewhere.
many-dimensional parameter regions where the likelihood function is mostly flat -especially when MC event simulation is performed for each scan point.The high detectability of final states with photons and missing energy ensures that current LHC searches can exclude specific scenarios of electroweakino production where the masses of the produced electroweakinos are close to or beyond 1 TeV.These are typically scenarios with production of a dominantly wino chargino-neutralino pair and a large BR( χ0 1 → γ G) [113,124].But as we will see, within the general electroweakino parameter space explored here, there are still large, unconstrained parameter regions with all electroweakinos 1 TeV.
In Tab. 1 we summarise our choices for the scan input parameters.The MSSM parametrisation we use is implemented in the GAMBIT MSSM model hierarchy as MSSM11atQ_mA_mG (Appendix C), which has 11 free parameters.For the six parameters not listed in Tab. 1 we use the following fixed values: the trilinear couplings A d3 = A e3 = A u3 = 0; the gluino mass parameter M 3 = 5 TeV; the pseudo-scalar Higgs mass m A = 5 TeV; and the squared soft sfermion mass parameters m 2 l = m 2 q = (3 TeV) 2 .The parameters are defined at an input scale Q = 3 TeV.The specific values for these fixed parameters are not important, as they simply ensure that all superpartners except the gravitino and the electroweakinos are decoupled from the collider phenomenology.
In order to obtain accurate profile likelihood maps we must ensure that the parameter space is explored in sufficient detail.We therefore combine the parameter samples from multiple scans using different combinations of the priors (metrics) listed in Tab. 1 to scan the parameters.The "hybrid" prior in Tab. 1 combines a logarithmic prior for |x| > 10 GeV with a flat prior for |x| < 10 GeV (x = M 1 , M 2 , µ).As the physics is invariant under a global sign change for M 1 , M 2 and µ, we follow the common approach in the literature of re-stricting M 2 to positive values.All scans are performed with the differential evolution sampler Diver 1.0.4[186], interfaced via ScannerBit.We run Diver in the jDE mode (self-adaptive rand/1/bin evolution), which is based on Ref. [200].The final combined data set consists of around 3.1 × 10 5 parameter samples.

Best-fit scenarios
In Fig. 2 we show our fit result in terms of the profile likelihood function across the (m χ0 2 , m χ0 1 ) and (m χ± 1 , m χ0 1 ) planes.We will present most of our results in one or both of these planes as they are well suited for mapping out the key phenomenological aspects of the high-likelihood scenarios.For reference, in Appendix A we provide profile likelihood maps in terms of the input parameters.
We find that the G-EWMSSM scenarios in best agreement with current LHC searches and measurements are scenarios where the lightest electroweakinos are dominantly Higgsino, i.e. scenarios with |µ| < |M 1 |, M 2 , corresponding to the Higgsino NLSP region (green) in Fig. 1.As the µ parameter largely controls the mass of three Higgsino states, these scenarios have near-degenerate masses for χ0 1 , χ0 2 and χ± 1 , explaining why the best-fit region falls along the diagonals of the (m χ0 2 , m χ0 1 ) and (m χ± 1 , m χ0 1 ) planes.For the best-fit point, marked by a white star in Fig. 2, the three Higgsinos have masses m χ0 1 = 169.9GeV, m χ0 2 = 178.9GeV and m χ± 4 at m χ0 4 = 788.1 GeV.The scenarios allowed at 2σ confidence level (CL) relative to the best-fit point, all predict such a trio of near-degenerate Higgsinos with masses no less than about 140 GeV and no greater than about 500 GeV.
The scenarios within the 2σ region in Fig. 2 are largely scenarios with negative µ parameter, |µ| < M 2 , |M 1 |, and tan β 5, with the highest-likelihood solutions favouring tan β values close to 1.For such scenarios, the dominant and subdominant decay modes for the lightest neutralino are the χ0 1 → h G and χ0 1 → Z G channels, respectively -see e.g. the region around µ ∼ −300 GeV in the branching ratio plots in Fig. 1.
Low branching ratios for decays to γ G final states ensure that the scenarios in the 2σ region escape the otherwise highly constraining photons + E miss T searches.Many of these scenarios also have sizeable branching ratios for χ0

through χ0
2 → Z * χ0 1 , as often assumed in LHC searches for Higgsino production.Similarly, many scenarios in the higher-mass part of the 2σ region (m χ± 1 300 GeV) have large branching ratios for direct decays of χ± 1 to the gravitino, through χ± 1 → W ± G.By tuning the branching ratios BR( χ0 4 the model can partly fit small excesses in the ATLAS and CMS leptons + E miss T searches and the ATLAS b-jets + E miss T search.(The preference for a small signal contribution in b-jet final states in part explains the preference for tan β ∼ 1, since this increases the branching ratio for χ0 1 → h G, see Sec. 2.) In combination, this produces a weak preference for the lower-mass end of the diagonal in Fig. 2, at masses around 170 GeV. 5e found a preference for low-mass electroweakino scenarios also in our EWMSSM fit in [17].The EWMSSM parameter regions favoured in that study allow for electroweakino decay chains that produce multiple on-shell Z, h and W bosons, and terminate in a binodominated χ0 1 that provides the missing energy signal.The favoured low-mass scenarios in the G-EWMSSM predict a similar collider phenomenology, but now with the gravitino rather than a bino-like neutralino terminating the decay chains.However, in the present study the preference for these low-mass scenarios is weaker, as the previously-observed data excesses are less pronounced in the now larger ATLAS and CMS data sets.

Non-excluded scenarios
Assuming that these small data excesses are just background fluctuations rather than a true BSM signal, it is interesting to consider what electroweakino mass combinations the current combined data clearly exclude in the G-EWMSSM.We investigate this in Fig. 3 by showing profile likelihood plots where we use the capped likelihood, L cap searches (Eq.14), as described in Sec.3.1.To understand the structures visible in Fig. 3, we first consider Fig. 4, where we show the Higgsino, wino and bino components of the lightest neutralino for the highest-likelihood point in each bin across the (m χ0 2 , m χ0 1 ) plane.This allows us to identify which of the three NLSP scenarios discussed in Sec. 2 are preferred in different parts of the mass plane.We see clearly that the preferred scenarios along the diagonal are scenarios with a mostly Higgsino NLSP (left panel), as discussed above.Moving away from the diagonal, towards higher m χ0 2 , the best-fitting scenarios are wino NLSP scenarios (middle panel).We note that around m χ0 1 , m χ0 2 ∼ 400 GeV, the current collider data prefers a fairly even wino/Higgsino admixture for the χ0 1 .Finally, at even higher χ0 2 -χ0 1 mass splittings, the best possible fits are obtained for bino NLSP scenarios (right panel).We will in the following use the term profile-likelihood surface to refer to the set of parameter samples that appear in figures like Fig. 4, where for each bin in the given plane we visualise some property of the highestlikelihood parameter sample belonging to that bin.For the interpretation of these figures it is important to remember that apparent discontinuities, such as the boundaries between the yellow and black regions in Fig. 4, typically result from the projection done by the profile likelihood procedure: two neighbouring bins in a mass plane can have their respective highest-likelihood points coming from very different parts of the fourdimensional G-EWMSSM parameter space.So for instance the black region in the right-hand panel of Fig. 4 does not imply that there are no parameter samples that predict the given χ0 1 and χ0 2 masses and a binoa consequence of our scan settings: Since we restrict our study to the parameter space that has all electroweakino masses below 1 TeV, having the lightest neutralino mass close to 1 TeV will correspond to parameter points with |M1| ∼ M2 ∼ |µ| ∼ 1 TeV.dominated χ0 1 , only that there for these mass predictions exist other parameter points that give a better fit to data and for which the χ0 1 is dominantly wino or Higgsino.We can now go back and reconsider Fig. 3. Along the diagonals of the two mass planes, we see the allowed scenarios with Higgsino-dominated χ0 1 , χ0 2 and χ± 1 .This region extends all the way up towards the edge of our scan range, corresponding to masses around 1 TeV.In addition, there are three other non-excluded scenarios visible.

First, in the (m χ0
2 , m χ0 1 ) plane, we find an allowed horizontal region at around m χ0 1 ≈ 450 GeV, with winodominated and mass degenerate χ0 1 and χ± 1 .Second, in the region of m χ0 1 450 GeV and m χ0 2 , m χ± 1 800 GeV, we see solutions with a lonely, light, bino-dominated χ0 1 .Lastly, in the (m χ0 2 , m χ0 1 ) plane around m χ0 1 > 700 GeV and away from the diagonal, we see a region of solutions allowed at 2σ, where again the χ0 1 and χ± 1 are mostly wino, though with non-negligible Higgsino components.
Before we explore these findings further, let us briefly compare them with the capped-likelihood results from our analysis of the EWMSSM [17].In [17] we found that essentially no combinations of χ± 1 and χ0 1 masses could be conclusively ruled out by the combination of LHC search results at the time of that study.The conclusion is markedly different in the present G-EWMSSM study, where only four distinct scenarios for electroweakinos below 1 TeV remain viable.There are several factors contributing to this result: (1) the overall stronger constraining power due to the now larger LHC data sets; (2) the diminishing of the data excesses that in [17] helped improve the fit for low-mass solutions in the EWMSSM; (3) the additional constraining power in the present study, coming from our inclusion of LHC measurements in addition to direct BSM searches; and (4) the distinctive G-EWMSSM collider signatures, in particular the photon signatures, that result in strong constraints on large parts of the G-EWMSSM parameter space.

Impact of different searches
To understand our results in greater detail, we will in the following discuss the contributions from the LHC searches and measurements that most strongly influence the fit result.To aid this discussion we consider Figs. 5, 6 and 7: In Fig. 5 we show the total log-likelihood difference ln L(s) − ln L(s = 0).The various solutions in Figs. 2 and 3 are visible as regions of greater likelihood.In Fig. 6 we consider the profile likelihood surface for the (m χ0 2 , m χ0 1 ) plane and break the total log-likelihood down into contributions from photon searches, lepton searches, other searches and measurements of SM-like final states.Finally, in the six panels of Fig. 7 we show the total electroweakino LHC production cross-section and a selection of relevant branching ratios across the (m χ0 2 , m χ0 1 ) profile likelihood surface.The top-left panel of Fig. 6 shows that for the scenarios with a bino NLSP (see Fig. 4, right), the most constraining LHC analyses are the photons + E miss T searches.This can be understood from the fact that for these scenarios the dominant χ0 1 decay mode is χ0 1 → γ G (Fig. 7, top right), while the heavier wino-or Higgsinodominated electroweakinos, which here dominate the production cross-section, decay via the χ0 1 rather than directly to a G final state (Fig. 7, bottom right).Towards larger masses for the heavier electroweakinos the production cross-section diminishes (Fig. 7, top left) enough to leave an allowed region at m χ0 2 , m χ0 1 ) plane, where the highest-likelihood scenarios are wino NLSP scenarios (see Fig. 4, middle), the most important contributions to the profile likelihood surface come from the leptons + E miss T searches (Fig. 6, top right), and searches for jets + E miss T final states, with or without leptons (Fig. 6, bottom left).This is largely explained by the fact that the dominant decay modes of the now wino-dominated and near mass-degenerate χ0 1 and χ± and χ± 1 → W ± G, respectively (Fig. 7, middle right and bottom left).Thus, χ± 1 χ0 1 production will for these scenarios typically give rise to the same collider signatures as the commonly studied SUSY scenarios where wino-dominated χ± 2 are produced and decay to final states with a stable, light χ0 1 through χ0 2 → Z χ0 and χ± 1 → W ± χ0 1 .However, while χ± 1 is the most important production mode for these G-EWMSSM scenarios, relevant signal contributions can also arise from production of some of the heavier, Higgsino-dominated electroweakinos.Towards low m χ0 1 (m χ0 1 200 GeV), phase space suppression of the χ0 1 → Z G decay makes χ0 1 → γ G the dominant decay mode for χ0 1 (Fig. 7, top right).Here the photons + E miss T searches contribute strongly to the total log-likelihood, as does the measurements of SM signatures, to be discussed in more detail below (Fig. 6, top left and bottom right).At around m χ0 1 ∼ 450 GeV, the reduction in the production crosssection with increasing mass (Fig. 7, top left), combined with a balancing of the χ0 1 → γ G and χ0 1 → Z G branching ratios (Fig. 7, top right and middle right) means that the combined constraining power of the searches is sufficiently weakened so that a horizontal band in the mass plane avoids exclusion at the 2σ level.This is also partly due to the model fitting some weak excesses in leptons + E miss T and photons + E miss T searches (light blue bands in Fig. 6, top left and top right).However, towards even higher m χ0 1 , the ATLAS search for E miss T + boosted bosons [100] gains sensitivity (Fig. 6, bottom left) and the total likelihood therefore drops below the 2σ threshold for m χ0 1 between ∼ 500 GeV and ∼ 700 GeV.

Impact of measurements
The present study is the first to include LHC measurements of SM signatures in a many-parameter BSM global fit.It is therefore interesting to explore what impact these likelihood contributions have on our results.The log-likelihood contribution ∆ ln L meas on the (m χ0 2 , m χ0 1 ) profile-likelihood surface is shown in the bottom-right panel of Fig. 6.The contribution is significant in the regions with wino-or Higgsino-dominated χ0 2 , m χ0 1 ) profile likelihood surface discussed above is by definition made up of the overall least constrained parameter sample within each (m χ0 2 , m χ0 1 ) bin.To get a more complete picture of the constraining power of the SM signature measurements, it is interesting to also look at ∆ ln L meas across the surface of parameter samples that are most strongly constrained by this loglikelihood contribution.This is shown in the top-left panel of Fig. 8.For the G-EWMSSM scenarios where the SM signature measurements have their largest sensitivity, they rule out scenarios that have both m χ0 2 and m χ0 1 below ∼ 500 GeV, and scenarios towards higher m χ0 2 when m χ0 1 150 GeV.The three other panels in Fig. 8 show the individual log-likelihood contributions from the pools of measurements that contribute most strongly to the combined ∆ ln L meas in the upper-left panel: ATLAS measurements of the pp → ZZ → 4l cross-section (top right) [141, 159, 174]; ATLAS measurements of final states with two different flavour leptons and missing energy, with or without jets (bottom left) [140,149,162,165], where the dominant contribution is coming from the pp → W + W − cross-section measurements in [162,165]; and an ATLAS measurement of the pp → Z(→ l + l − )γ + X production cross-section [180] (bottom right).
In Fig. 9 we show the χ0 1 composition for the parameter samples contributing to Fig. 8. From Figs. 8  and 9 we see that the ZZ cross-section measurements most strongly constrain low-mass scenarios where the χ0 1 is dominantly Higgsino or a wino-Higgsino mixture.These G-EWMSSM scenarios combine a large total electroweakino production cross-section, 7 with significant 7 A balanced wino-Higgsino mixture for a low-mass χ0 1 implies that M2 for these points typically is within ∼ 100 GeV of |µ|.This means that at least four of the five heavier electroweakino states will have masses not too much larger than m χ0 1 .
branching ratios for some of the decays χ0 i → Z G and/or χ0 i → Z χ0 j .The measurements of W + W − production cross-sections exclude low-mass scenarios with winodominated χ0 1 .Here the strongest W + W − signal contribution comes from the production of pairs of light, wino-dominated χ± 1 , which decay as χ± 1 → W ± G. Finally, the Z(→ l + l − )γ + X cross-section measurement constrains scenarios with bino-dominated χ0 1 .These scenarios typically have a large BR( χ0 1 → γ G) and a nonnegligible BR( χ0 1 → Z G), such that production of any pair of electroweakinos that decay to χ0 1 's can result in signal contributions to the measured cross-section.
Since the best-fit region predicts light Higgsinos, at masses around 170 GeV, the LHC searches and measurements performed at √ s = 8 TeV can also be relevant.A full investigation of the impact of 8 TeV results is beyond the scope of this study, as it would effectively double the computational cost of our parameter scans.However, to gauge the possible impact, we generate 100 000 events at 8 TeV for each of our 100 highest-likelihood parameter points.We pass the events through Rivet and Contur to compute a log-likelihood contribution from the collection of √ s = 8 TeV measurements in Rivet.The result of this is illustrated in Fig. 10, where we show the change in the total log-likelihood for each point when the contribution from 8 TeV measurements is added.In the left-hand panel we show the points close to the best-fit point at m χ0 1 ∼ 170 GeV.Of our 100 highest-likelihood points, some also belong to the higher-mass region, at m χ0 1 280 GeV, shown in the right-hand panel.For the best-fit points in the low-mass region, including the 8 TeV measurements reduces the total log-likelihood by around 0.2 units.As expected, there is a smaller impact on points in the higher-mass region.

Scenarios with a chargino NLSP
In contrast with the EWMSSM, the G-EWMSSM admits the possibility of a chargino as the lightest electroweakino.Such a scenario was highlighted in Fig. 1 where the gray band signals a sudden drop in branching ratio due to m χ± 1 < m χ0 1 .While rare for MSSM-like electroweakino mass matrices, and featuring small mass differences, our scan identified still-viable parameter regions with m χ± 1 < m χ0 1 , shown in Fig. 11.We find that in these cases, the points with the highest likelihoods have Higgsino-like electroweakinos, with only small splittings for the χ± 1 , χ0 1 and χ0 2 , with masses preferred to be in the region of 400-500 GeV.Here, the decay mode for χ± 1 is always χ± 1 → W G. Hence, the detectable signal for χ± 1 χ± 1 pair production is two on-shell W bosons and some missing energy from the gravitinos.For the χ0 1 , the dominant decay modes are χ0

and χ0
1 → h G due to the dominant Higgsino component.The detectable signal for χ0 1 χ± 1 production would then be on-shell W Z or W h plus missing energy from the gravitinos.Finally, χ0 2 decays to soft SM fermions and the χ0 1 or χ± 1 .Thus, the production of χ0

Conclusions
In this study we have investigated the current viability of the G-EWMSSM, the simplest realisation of a light supersymmetric electroweak sector together with a nearly massless gravitino LSP.We have confronted the G-EWMSSM with a comprehensive selection of the relevant Run 2 searches at the LHC, relevant past searches at LEP, and, we have, for the first time in a global fit, used a broad set of SM measurements at the LHC to constrain the model by building a new interface between GAMBIT and Contur.
Our best-fit region for the model is where |µ| < |M 1 |, M 2 , and is characterised phenomenologically by a trio of relatively light degenerate Higgsinos in the mass range of 140-500 GeV, with a best fit point around 170 GeV.Due to the collective effect of small excesses over multiple ATLAS and CMS searches we find closed 2σ contours in the parameter space, but we emphasise that this is a model-specific best-fit region and does not constitute a measure of goodness-of-fit.
Our main result is that the bulk of the G-EWMSSM parameter space with electroweakino masses below 1 TeV is excluded by collider searches and measurements.The four exceptions, classified according to the nature of the lightest electroweakinos, are: i) degenerate Higgsinos from 140 GeV and up, ii) a region of degenerate winos around 400-500 GeV allowed at the 2σ level, iii) degenerate winos above 700 GeV, and iv) a 'lonely' bino from 62 GeV and up, decoupled from heavier Higgsinos and winos lying above 800 GeV.
For Run 3 of the LHC the degenerate Higgsino region, i), will be challenging to test fully.Drawing from the lessons learnt in this study, the measurement of SM multi-lepton signatures will continue to be important to exclude parameter space at the low-mass end of the region.Potential improvements to searches sensitive to the important χ0 1 → h G decay (see Fig. 7, middle left), will also improve the reach.However, fully excluding this still very viable region will need future e + e − or muon colliders operating at high enough centre-of-mass energies.
On the other hand, the surviving wino band, ii), with masses around 450 GeV seems to be fully excludable with the slightly higher Run 3 centre-of-mass energy and more data, in particular since its survival is already marginal.For the same reason it should also be possible to push the remaining wino region, iii), to somewhat higher masses with higher cross sections and more data.
For the lonely bino region, iv), the search for pair production of light binos decaying to photons is also hampered by low production cross sections.However, we expect some impact here with increasing statistics in Run 3 and beyond to the High-Luminosity LHC, in particular on the parts of parameter space where there is bino production through the decay of heavier electroweakinos, which could realistically be pushed out beyond 1 TeV.
We emphasise the still open interesting possibility of a reverse mass hierarchy of charginos and neutralinos, with m χ± 1 < m χ0 1 , with distinct signal predictions for LHC Run 3 searches.Although the base production cross section is not so high given their Higgsino nature, the preferred region of this scenario should be within reach of Run 3 statistics and the slightly higher centreof-mass energy, when considering all final states W W , W Z and W h.
We make all our generated parameter samples available from Zenodo for further study [201].

The ATLAS search for electroweak production of charginos and neutralinos in final states with two boosted, hadronically-decaying bosons and missing transverse momentum [100]:
This search (ATLAS_2BoostedBosons) targets the pair production of electroweakinos, where each of them is assumed to decay into the LSP and an on-shell W , Z or SM Higgs boson.The mass difference between the produced electroweakinos and the LSP is assumed to be at least 400 GeV.The analysis is optimised on three different scenarios: 1) a baseline MSSM scenario where the produced electroweakinos and the LSPs can be either binos, winos or Higgsinos, 2) a general gauge mediationinspired scenario in which the LSP is a gravitino and the heavier particles are Higgsinos and 3) a scenario with an axino LSP, where the heavier particles are assumed to be Higgsinos.Various simplified models are considered in each case.The analysis is performed  jets we do not include the signal regions that rely on this.

The ATLAS search for gluino and squark production in final states with jets and missing transverse momentum [101]:
This is the flagship ATLAS supersymmetry search for squarks and gluinos (ATLAS_0lep), targeting events with multiple jets and significant missing transverse momentum.Although it is optimised on models of squark and gluino production, similar final states can be produced by electroweakino production with subsequent cascade decay processes that produce hadronically-decaying gauge bosons.We implement the optimised single-bin signal regions that are designed to present the ATLAS results in a model-independent way (2j-1600, 2j-2200, 2j-2800, 4j-1000, 4j-2200, 4j-3400, 5j-1600, 6j-1000, 6j-2200 and 6j-3400).The signal region selections include requirements on the multiplicity and transverse momenta of the jets in each event, the angular separation between the jets and the missing transverse momentum vector, the aplanarity, E miss T / √ H T and m eff .

The ATLAS search for top squarks in the jets plus missing transverse momentum final state [102]:
This search (ATLAS_0lep_stop) seeks to uncover evidence of stop production in final states with four or more jets plus missing transverse momentum.Five sets of signal region are defined in the analysis, targeting different stop simplified models, with a range of different included sparticles and sparticle mass differences.The six SRA and SRB regions employ top-mass reconstruction to increase sensitivity to models in which the stop produces a top quark, which makes them less relevant for the scenario considered in this paper.The five SRC regions use recursive jigsaw variables to target regions with a small t1 − χ0 1 mass difference, the details of which are highly-dependent on the treatment of initial state radiation in the Monte Carlo generator used to model LHC events.We do not include these SRC regions due to known deficiencies of the Pythia initial state radiation model in this region.The two SRD regions are optimised for direct top squark production where both top squarks decay via t → b χ± 1 .At least five jets are required, two of which must be b-tagged, and further requirements are placed on the jet transverse momenta and the scalar sum of the transverse momenta of the two jets with the highest b-tag weights.Finally, the SRE signal region is designed for models with highly boosted top quarks.Requirements on the jet mass of reclustered fat jets are used, alongside requirements on the main is optimised on simplified models of stop production with decays that produce one lepton (through a real or off-shell leptonically-decaying W boson), and also on a dark matter model with a spin-0 mediator produced in association with two top quarks.All signal region are required to have exactly one signal lepton, and 2, 3 or 4 jets.Five regions labelled tN are optimised for the decay pattern t → t χ0 1 , using selections on variables such as the am T2 variable, the transverse mass m T formed from the lepton and missing transverse momentum, the H miss T,sig 8 and the mass of a reconstructed hadronic top quark.Note that we do not include three signal regions that use a boosted decision tree in the definition of the signal region, since this is very difficult to reproduce outside of the ATLAS collaboration.Two additional signal regions, bWN and bffN, are dedicated to the three-body ( t → bW χ± 1 ) and four-body ( t → bf f χ± 1 ) decay searches.Six signal regions target various t → b χ± 1 scenarios: three are optimised on a simplified model that assumes m χ± 1 = 2m χ0 1 (labels bC2x_diag, bC2x_med, bCbv), and three are designed to search for the case of a Higgsino LSP, in which the χ± 1 , χ0 2 and χ0 1 are close in mass (labels bCsoft_diag, bCsoft_med, bCsoft_high).In the latter case, the signature is characterised by low-momentum leptons or jets from highly off-shell W or Z bosons, and the analysis benefits from a dedicated soft lepton reconstruction.Finally, three extra signal regions (DM_low_loose, DM_low, DM_high) are optimised on the dark matter mediator model, with the analysis using similar variables to the regions targeting the decay t → t χ0 1 .

The ATLAS search for top squarks in final states with two opposite-charge leptons and missing transverse momentum [104]:
This search (ATLAS_2lep_stop) is optimised on similar models of direct stop production to the 0 lepton and 1 lepton searches.Events are required to have exactly two light leptons (electrons or muons) of opposite charge, with an invariant mass outside of the Z boson mass window in the case of same flavour leptons.A series of discriminating variables are constructed from the missing transverse momentum and p T 8 H miss T,sig is defined as , where H miss T is the negative vectorial sum of the momenta of the signal jets and the lepton, M = 100 GeV is an offset parameter, and the denominator is computed from the per-event jet energy uncertainties.
values of the leading leptons and jets, with other useful variables including a variant of m T2 and the super-razor variables first defined in Ref. [203].Various signal regions are optimised for 2-body, 3-body and 4-body stop decays.For the case of 2-body decays, the ATLAS analysis also defines a set of seven inclusive signal regions (labelled SR2bInc) intended to provide less model-specific sensitivity.Our implementation of the search uses this set of inclusive signal regions.

The ATLAS search for top squarks in events with a Higgs or Z boson [105]:
This search (ATLAS_2OSlep_Z) is optimised on various simplified models of top squark production in which a top squark decays to produce a Higgs or Z boson.Top squark decays involving Z bosons are targeted using a 3-lepton selection, with at least one same-flavour-opposite-sign pair (SFOS) whose invariant mass is consistent with Z boson mass.Further selections are placed on the transverse momenta of the three leading leptons, the jet multiplicity, the b-jet multiplicity, the transverse momenta of the leading jet and b-jet, the missing transverse energy, a variant of m T2 and the transverse momentum of the SFOS pair.Events containing Higgs bosons are targeted using a 1-lepton event selection, with further selections placed on the jet and b-jet multiplicity, the transverse mass formed from the lepton and the missing transverse momentum, and the missing transverse energy significance.In addition, a Higgs tagger built from a neural network is used to identify Higgs boson candidates, and events must contain at least one of them.Due to the difficulty of reproducing this Higgs tagging with sufficient accuracy, our implementation of this search covers only the 3-lepton final states (labels SRZ1A, SRZ1B, SRZ2A, SRZ2B).

The ATLAS search for charginos and sleptons in final states with two leptons and missing transverse momentum [106]:
This search (ATLAS_ 2OSlep_chargino) is optimised on simplified models of slepton and chargino production, targeting chargino pair production with decays to lightest neutralinos and W bosons, chargino cascade decays through sleptons to lightest neutralinos, and the direct production of slepton pairs.Events are required to have exactly two opposite-charge light leptons with an invariant mass greater than 100 GeV.Selected events must also have no b-tagged jets, and large values of E miss T and E miss T significance.Further discrimination comes from the use of the m T2 variable.Four sets of signal regions are defined (labels SR-SF-0J, SR-SF-1J, SR-DF-0J, SR-DF-1J) based on whether the leptons have the same or a different flavour, and whether the events have 0 or 1 non-b-tagged jets.From this, the ATLAS analysis defines a total of 16 inclusive signal regions, intended for more model-independent sensitivity, and a set of 36 signal regions with fine-grained binning in m T2 , to maximise sensitivity to the simplified model studied by ATLAS.In our study we use the inclusive signal regions.
The ATLAS search for bottom and top squarks in final states with two b-tagged jets and missing transverse momentum [107]: This search (ATLAS_2b) is optimised on various simplified models of stop and sbottom production, targeting final states with 2 b-tagged jets, large missing transverse momentum and either zero leptons or one lepton.A long list of discriminating variables is used, including the minimum ∆Φ between any of the leading jets and the missing transverse momentum vector, H T (defined as the scalar sum of the p T values of a subset of the jets in the event), m eff , ratios of the missing transverse energy with m eff and √ H T , the contranverse mass, am T2 and others.Three zero lepton signal regions and three one lepton signal regions are defined.Due to challenges in reproducing the cuts based on am T2 our study only uses the zero lepton signal regions (labels 0L_SRA350, 0L_SRA450, 0L_SRA550, 0L_SRB, 0L_SRC).

The ATLAS search for Higgsinos in final states with at least three b-tagged jets [108]:
This search (ATLAS_3b) targets Higgsino production and decay in gauge-mediated supersymmetry scenarios, in which each Higgsino is assumed to decay to a Higgs boson and a gravitino.Two complementary analyses, targeting high-and low-mass signals, are performed.For the high-mass analysis, events with at least three b-tagged jets are selected, and jet pairs are assigned to two Higgs candidates.For the low-mass analysis, events with four b-jets are analysed by grouping the jets into Higgs candidates.Selections are placed on a number of kinematic variables including m eff , E miss T , the mass of the Higgs boson candidates, angular variables and the minimum transverse mass formed with the missing transverse momentum vector and any of the leading four jets.Due to some overlaps between the signal regions for the low-mass and high-mass analyses, our analysis only uses the low-mass signal regions.From this analysis we have implemented all the 46 signal regions optimised for exclusion.
The ATLAS search for chargino-neutralino pair production in final states with three leptons and missing transverse momentum [109]: This search (ATLAS_3lep) is optimised on two scenarios of electroweakino production.In the first, a χ± 1 and χ0 2 are produced (both wino-dominated), with subsequent decay to a bino-dominated χ0 1 .In the second, the χ± 1 , χ0 2 and χ0 1 are pure Higgsino states, and are therefore typically more mass degenerate (although an arbitrary mass hierarchy is assigned in order to define a parameter plane in which to optimise the analysis).The analysis has three dedicated selections to cover different mass regimes and assumptions, including an on-shell W Z selection, an off-shell W Z selection and a W h selection.All consider final states with exactly three leptons, possible ISR jets and E miss T .Events with at least one SFOS pair are divided into three bins of the SFOS pair invariant mass, m ll , covering the regions below, on and above the Z mass.Each m ll bin is further divided into E miss T and m T bins, where the transverse mass m T is defined using the lepton that is not in the SFOS pair (and which can therefore be assumed to arise from a W boson decay).Events are further separated by their jet multiplicity, and by two different variants of H T , defined as the scalar sum of the transverse momenta of jets or leptons depending on the definition.Signal regions for events with a different-flavour-opposite-sign (DFOS) lepton pair are defined separately, using selections on the jet multiplicity, E miss T significance, transverse momentum for the third-leading lepton, and the ∆R between the DFOS leptons and the same-flavour-same-sign lepton that is nearest in φ.Our implementation includes 39 of the 41 signal regions targeting on-shell W Z or W h production, leaving out the two regions SR-Wh-DFOS-1 and SR-Wh-DFOS-2 for which some cuts rely on object resolution variables that are not available in our fast event simulation framework.
The ATLAS search for gluino, electroweakino or slepton production in final states with four or more leptons [110]: This search (ATLAS_4lep) is optimised on different R-parity violating and R-parity conserving SUSY scenarios that can produce leptonrich final states.The simplified model used for the R-parity conserving scenarios is a model with Higgsino production and a gravitino LSP, thus highly relevant for our study.Events are required to have four or more leptons (electrons, muons or hadronically-decaying taus).The signal regions for four-lepton events are classified by whether the events have four light leptons and zero taus, three light leptons and one tau, or two light leptons and two taus.Further selections are placed on the number of b-tagged jets, the presence or absence of a Z boson, E miss T and m eff .A further signal region is also defined, requiring at least five light leptons, subject to no kinematic requirements.Our implementation includes all the zero-tau signal regions, i.e. the five-lepton region (label SR5L) and the seven regions with four light leptons (the regions with the SR0 label).In particular, this includes the the two signal regions SR0-ZZ-loose-bveto and SR0-ZZ-tight-bveto designed to target the Higgsino plus gravitino scenario.
The ATLAS search for squarks and gluinos in final states with same-sign leptons and jets [111]: This search (ATLAS_MultiLep_strong) is optimised on simplified models of gluino and squark production, covering the case of both R-parity conversation and R-parity violation.Events are required to have two same-sign leptons and may contain additional leptons.In the case of the R-parity conserving search, large missing transverse momentum is required.Five signal regions are defined using the number of leptons and their relative electric charges, the number of jets and the number of b-tagged jets.Key kinematic variables used include the effective mass m eff , E miss T and its ratio to m eff and the invariant mass of same-sign electron pairs (which reduces contamination from Z → e + e − decays where the charge of one electron is mismeasured).We include all five signal regions in our analysis.
The ATLAS search for gauge-mediated supersymmetry in final states with photons, jets and missing transverse momentum [112]: This search (ATLAS_PhotonGGM_1photon) is optimised on a simplified model in which pair-produced gluinos decay to neutralinos, which in turn decay to a gravitino, at least one photon and jets.Three signal regions are defined which target the cases of large, medium and small mass differences between the gluino and neutralino, and all of them veto leptons in the selected events.Further selections are placed on the transverse momentum of the leading photon, the jet multiplicity, the angular separations of the jet and photon momenta with the missing transverse energy vector, E miss T , H T and a variable called R 4  T , defined as the ratio of the scalar sum of the p T for the four leading jets, and the scalar sum of the p T for all signal jets in the event.Our implementation includes all three signal regions.

The ATLAS search for photonic signatures from gauge-mediated supersymmetry models [113]:
This search (ATLAS_PhotonGGM_2photon) is optimised on models of both strong and electroweak sparticle production, and targets final states with either a single photon and multiple jets, or two photons, plus significant missing transverse momentum in both cases.Discriminating variables include m eff , E miss T , a variant of H T , the angular separation between photons and the missing transverse momentum vector and R 4 T (the scalar sum of the transverse momenta of the four leading jets divided by the scalar sum of the transverse momenta of all jets in the event).Our implementation contains all signal regions from the paper, but in our analysis we only use the two-photon signal regions (labels SRaa_SL, SRaa_SH, SRaa_WL, SRaa_WH), since the one-photon signal regions largely overlap with (and are superseded by) the signal regions in ATLAS_PhotonGGM_1photon above.

The ATLAS search for exotic decays of the Higgs boson to at least one photon and missing transverse momentum [114]:
This search (ATLAS_ Z_photon) targets exotic decays of the Higgs boson to, for example, a gravitino and a lightest neutralino, with the neutralino subsequently decaying to a gravitino and a photon, This generates a final state with a single photon plus missing transverse energy, and one can reduce SM backgrounds by looking for events with a Higgs boson produced in association with a Z boson.One can also generate final states with two photons if the Higgs boson decays to a pair of neutralinos that subsequently decay.Events are selected in the analysis if they have at least one photon, moderate E miss T and two opposite-sign electrons or muons with an invariant mass within 10 GeV of the Z boson mass (and no additional leptons).A γE miss T system is defined by performing the vector sum of the photon momentum or momenta and the missing momentum vector in the transverse plane.The search relies on two discriminating variables that quantify the angular separation of the two-lepton and γE miss T systems and the p T asymmetry of the two systems, and the analysis has only one signal region.
The CMS search for gluino and squark production in final states with multiple jets and missing transverse momentum [115]: This is the flagship CMS search for gluino and squark production (CMS_0lep), and is the CMS equivalent of the ATLAS search presented in Ref. [101].Search regions are defined in a four-dimensional space of variables given by the total number of jets, the number of b-tagged jets, the scalar sum of the jet p T values (H T ), and the magnitude of the vector p T sum of the jets (H T miss ).In total, there are 174 exclusive signal region bins.For our analysis, we implement 12 aggregate search bins which are presented in an appendix of the paper, and which are constructed from the original search bins after taking correlations into account. .In our analysis we use both signal regions, and account for the correlated background uncertainties using the correlation coefficient provided by the CMS collaboration.

The
The CMS search for stop production in final states with one lepton [117]: This search (CMS_ 1lep_stop) is optimised on simplified models of stop production, and targets events with a single isolated electron or muon, jets and large missing transverse momentum.Two sets of signal regions are defined; one for a large range of t1 − χ0 1 mass splittings, and one for compressed spectra.For the first set of 27 signal regions, events are selected based on the number of jets, the E miss T , the invariant mass of the lepton and the closest b-tagged jet and a special "topness" variable.For the second set of 4 signal regions, at least five jets are required, and with the highest p T jet must not be b-tagged since it is expected to arise from initial state radiation.Further selections are placed on the angular separations of the missing transverse momentum and the lepton/jets, plus the p T of the lepton.Our analysis makes use of a smaller set of six aggregated signal regions that are provided in an appendix of the paper.
The CMS search for stop production and dark matter in final states with two opposite-charge leptons [118]: This search (CMS_2lep_stop) is optimised on models of stop production, and on dark matter models with a scalar or pseudo-scalar mediator in which the mediator is produced in association with a pair of top quarks.Events are selected if they have exactly two leptons with opposite charge and, in the case of a same-flavour lepton pair, the invariant mass of the lepton pair must not be close to the Z mass.Events must also have at least two jets, at least one b-tagged jet and a moderate amount of missing transverse momentum.Signal regions are defined in bins of three variables: two variants of m T2 , and E miss T , giving 13 signal regions, which are further split into 26 regions by separating events with same-flavour and different-flavour lepton pair.Our analysis includes all 26 signal regions, plus makes use of the covariance matrix provided by the CMS collaboration to account for the correlated background uncertainties.
The CMS search for charginos and neutralinos in final states with two low-momentum opposite-charge leptons [119]: This search (CMS_ 2lep_soft) is optimised on simplified models of chargino and neutralino production where the mass difference between the mass-degenerate χ0 2 and χ± 1 and the χ0 1 is small, such that decays proceed via off-shell W and Z bosons.A separate series of signal regions targets stop production.Selected events in both cases must contain two opposite-charge leptons (of either the same or different flavour) with a low transverse momentum, moderate E miss T and at least one jet.No b-tagged jets must be present and further selections are applied to variables such as the invariant mass and transverse momentum of the dilepton pair, E miss T /H T , H T and the transverse masses form from the leptons and the missing transverse momentum.For the electroweakino search, signal regions are defined in bins of E miss T and the dilepton invariant mass.For the stop search, signal regions are defined in bins of E miss T and the transverse momentum of the leptons.We implement all of the signal regions in our analysis, and treat the correlated background uncertainties using the covariance matrices provided by CMS.

The CMS search for supersymmetry in final states with two opposite-sign same-flavour leptons and missing transverse energy [120]:
This search (CMS_2OSlep) is optimised on various simplified models of gluino, squark, slepton and electroweakino production, and targets three potential signatures: 1) an excess of events with a lepton pair, whose invariant mass is consistent with the Z boson mass, 2) a kinematic edge in the invariant mass distribution of the lepton pair and 3) non-resonant production of two leptons.A set of strong production signal regions is defined using selections on the jet and b-jet multiplicities, H T , a variant of m T2 and E miss T , with the signal regions binned in the latter of these variables.On-Z electroweak production signal regions are defined using selections on the jet and b-jet multiplicities, the dijet invariant mass (sometimes defined using the b-jets), two different variants of m T2 and E miss T .In addition, a set of boosted signal regions is defined by requiring that there is a large radius jet with p T > 200 GeV, consistent with a hadronically-decaying gauge boson.For the dilepton edge search, a first approach is based on a fit to the dilepton invariant mass using events that pass selections on m T2 and E miss T .A second approach uses counts in various bins of E miss T , after applying other selections on m T2 and the jet and b-jet multiplicities.An additional requirement is placed on a novel variable that characterises how "t t-like" the events are.Finally, a set of slepton search regions are defined using selections on the jet and b-jet multiplicities, m T2 , E miss T and the ratio of the transverse momenta of the sub-leading lepton with the leading jet.Our implementation includes all signal regions except the edge fit regions, due to difficulty of implementing these outside of the CMS collaboration.We use the covariance matrices provided by CMS to account for correlated background uncertainties.

The CMS search for charginos and stops in final states with two opposite-charge leptons [121]:
This search (CMS_2OSlep_chargino_stop) is optimised on various simplified models of chargino and stop production and decay.Events are selected if they contain two opposite-charge electrons or muons, plus missing transverse momentum.For events with a same-flavour lepton pair, the invariant mass of the dilepton pair must not be close to the Z mass.For the chargino search, signal regions are defined in bins of the E miss T , number of b-tagged jets, number of jets, same-flavour or different-flavour status of the leptons and m T2 .For the stop search, an extra requirement is added on the number of "ISR jets", defined as jets with p T >150 GeV and no b-tag.We implement all of the chargino and stop regions, and make use of the covariance matrices provided by the CMS collaboration.
The CMS search for beyond-Standard Model physics in final states with jets and two samesign or at least three charged leptons [122]: This search (CMS_2SSlep_stop) is optimised on various simplified models of gluino production and decay, including both R-parity conserving and violating processes.Six exclusive categories of events are defined using preliminary selections on the lepton multiplicity and charge, plus E miss T .For events with two leptons, the leptons must have the same sign.For events with at least three charged leptons, separate categories are defined for the cases where there either is or is not an opposite-sign same-flavour pair with an invariant mass consistent with the Z boson mass.The main search is performed using a very large number of binned regions, and we instead implement the set of 17 inclusive signal regions (labels ISR1-ISR17) that are designed for easier interpretation of the results.These are defined in bins of the minimum transverse mass formed from either of the leptons and the missing transverse momentum, E miss T , the jet and b-tagged jet multiplicities and H T .
The CMS search for electroweak production of charginos and neutralinos in final states with three or four leptons, up to two hadronicallydecaying τ leptons or two same-sign light leptons [123]: This search (CMS_MultiLep) is optimised on a broad range of simplified models of chargino and neutralino production, including χ2 0 χ1 ± production with decays to lightest neutralinos via intermediate sleptons or with Higgs, W or Z bosons, and χ1 0 χ1 0 production with decays to gravitinos and Z or Higgs bosons.A number of event categories are defined in the analysis, based on the lepton multiplicity, whether there is a same-sign lepton pair or same-flavour-opposite-sign lepton pair, and the hadronically-decaying τ lepton multiplicity.The analysis is performed using a large number of binned signal regions for each category.
The key variables used are various variants of m T2 , the transverse momentum of the dilepton system in 2 lepton events, E miss T , the dilepton invariant mass of same-flavour-opposite-sign lepton pairs, the transverse mass formed from the trilepton system, dilepton system or a single lepton and the missing transverse energy vector, H T , the ∆R separation between leptons, and the invariant mass formed from a light lepton and a hadronically-decaying τ lepton.Our analysis includes the signal regions labelled SS01-SS20 (two same-sign leptons), A01-A64 (three leptons, one SFOS pair), B01-B03 (three leptons, no SFOS pair), G01-G05 (four leptons, two SFOS pairs) and H01-H03 (four leptons, one or zero SFOS pairs).
The CMS search for gauge-mediated supersymmetry in events with at least one photon and missing transverse momentum [124]: This search (CMS_photon) is optimised on various simplified models of gauge-mediated supersymmetry breaking, in which neutralinos and charginos decay to produce gravitinos and photons, as well as Z, W and Higgs bosons.Events are required to contain at least one high-energy photon plus large missing transverse momentum.The signal regions feature common selections on the transverse mass formed from the missing transverse momentum and the photon, and the E miss T .After these selections, four signal regions are defined as bins in S γ T , the scalar sum of E miss T and the p T of all photons in the event.We implement all four signal regions in our analysis.
The CMS search for gauge-mediated supersymmetry in events with two photons and missing transverse momentum [125]: This search (CMS_2photon) is optimised on simplified models of squark and gluino production, with cascade decays terminating in lightest neutralinos that always decay to a photon and a gravitino.Events are required to have two photons and large missing transverse momentum, and to satisfy various selections on their electromagnetic activity.Six signal regions are defined by selecting different ranges of E miss T , and all of these are implemented in our analysis.
The CMS search for gauge-mediated supersymmetry in events with a photon, an electron or muon and large missing transverse momentum [126]: This search (CMS_1photon_1lepton) is optimised on simplified models of gauge-mediated supersymmetry breaking that include gluino and quark production, plus direct production of neutralinos and charginos.Events are required to have at least one photon and at least one electron or muon.In the case of more than one light lepton, the lepton with the highest transverse momentum is used in the analysis.Selections are placed on the transverse mass formed from the lepton plus missing transverse momentum, and the E miss T , and multiple signal regions are defined in bins of E miss T , the transverse momentum of the photon and H T (defined as the scalar sum of the jet p T values in the event).The CMS search defines 18 signal regions per lepton channel (i.e.electron or muon), and we implement all regions in our analysis.

Appendix C Code extensions
In this Appendix we describe the extensions to the GAM-BIT framework introduced for this study, and how to use them.See Refs.[77,184] for a general introduction to the GAMBIT framework and Ref. [182] for an introduction to the ColliderBit module.

C.1 GAMBIT models with a light gravitino
The first modification is the addition of a new family of models MSSMXatY_mG, mirroring the existing family of MSSM models MSSMXatY in GAMBIT, supplemented with a new parameter mG which codifies the mass of a light gravitino.As in the existing MSSMXatY models, X refers to the number of parameters and Y to the scale at which the parameters are defined (which itself could be a parameter in MSSMXatQ models) (e.g.MSSM30atMGUT has 30 parameters defined at the GUT scale).Models with alternate parametrisation are labelled with increasing alphabetical letters after the number of parameters (e.g.MSSM10batQ) or in the specific case of reparametrisation with m A and µ instead of m Hu and m H d , models are labelled with the suffix _mA (e.g.MSSM30atMGUT_mA.The specific model used in the study is the MSSM11atQ_mA_mG, defined as MSSM11atQ_mA_mG: M1, M2, mu, TanBeta, Ad_3, Ae_3, Au_3, M3, Qin, mA, mG, ml2, mq2 An MSSM parametrisation with 11 parameters plus a gravitino mass mG, of which we vary only four in this study: TanBeta, M1, M2 and mu.We fix the other parameters to mG = 1 eV, Ad_3 = Ae_3 = Au_3 = 0, M3 = mA = 5 TeV and ml2 = mq2 = (3 TeV) 2 , to decouple all superpartners other than the electroweakinos and the gravtino from the low energy phenomenology.The input parameters are defined at a scale Qin = 3 TeV.We choose to use a model parametrised with mA and mu instead of mHu2 and mHd2, as µ controls the Higgsino masses.

C.2 Additions to DecayBit
The GAMBIT module DecayBit takes care of the computation of the decay widths of the various BSM particles [185].As detailed on Section 2, the main channel for the production of gravitinos is through the decays of the light electroweakinos.Therefore, for the purpose of this study we have implemented new capabilities and module functions for the computation of the decays of neutralinos and charginos to gravitinos, which can be seen in Table 3.The capabilities neutralino_i_decay_rates_gravitino, where i runs through the neutralino eigenstates, i.e. i= 1, . . ., 4, compute the decay of each of the neutralinos to a gravitino and a γ, h or Z, following eqs.( 5)- (7).The capabilities chargino_j_decay_rates_gravitino, where j runs through the chargino eigenstates, i.e. j= 1, 2, compute the decay of each of the charginos to a gravitino and a (possibly off-shell) W -boson.
In order to combine the newly added decays to gravitinos with the pre-existing decay channels of neutralinos and charginos, we implemented a new set of module functions that provide the capabilities neutralino_i_decay_rates and chargino_j_decay_rates, respectively.These module functions, also seen in Table 3 are called neutralino_i_decays_all and chargino_j_decays_all, and simply combine the decay tables computed by SUSY-HIT, via the capabilities neutralino_i_decay_rates_SH and chargino_j_decay_rates_SH, with the computation for the decays to gravitinos, via the capabilities introduced above, neutralino_i_decay_rates_gravitino and chargino_j_decay_rates_gravitino.

LHC_measurements histotags_perPool
Table 4: New ColliderBit capabilities for the ColliderBit added in this study, including a new LEP search, restructuring of the event generation capabilities and new capabilities for LHC measurements, using Rivet and Contur as backends.Note that the dependency type is only provided when there is ambiguity and that every Contur related function has a partner function with the prefix Multi_, which allows running different sets of Contur options in the same run.

C.3 Additions to ColliderBit
ColliderBit [182] is the GAMBIT module that handles anything related to collider physics, 9 such as LHC searches, LEP limits and Higgs measurements.It is naturally the most relevant module used in a collider-focused study such as this one.We extended ColliderBit significantly for this study, implementing new LHC analyses (Appendix B), adding a new LEP limit, and introducing the ability to compute the likelihood for BSM models to agree with LHC SM-like measurements, using the new interfaces to Rivet and Contur.
For the most part in this study we reuse the preexisting LEP upper limits on electroweakino production, which we re-interpret for the gravitino model.However, we added a new module function to compute the limit from an L3 multi-photon and missing energy search, which is exclusive of models such as this one.For this purpose we added a new module function, L3_Gravitino_LLike, providing a capability of the same name, to compute the 95% CL upper limit on the cross section.The new capabilitiy and module function can be found in Table 4.
The machinery for computing the likelihood for LHC measurements follows the structure of the existing Collid-erBit, where a Monte Carlo event generator (e.g.Pythia 8) is used to simulate hard scattering events at the LHC, each of those events is passed through a native detector simulation, then a collection of analyses, and a likelihood is computed from the resulting yields.In previous incarnations of ColliderBit, the events generated by Pythia 8 were immediately converted into a HEPUtils::Event type, which is needed for detector simulation and analysis of the events.However, to compute the predicted yields of a given simulated event for LHC SM measurements in Rivet, the event must be provided as a HepMC event (HepMC3::GenEvent type).This is done by splitting the existing module function generateEventPythia, which provides the capability HardScatteringEvent, into three parts, one that returns a Pythia8::Event, carrying the original name, one that converts the event to a HEPUtils::Event, called generateEventPythia_HEPUtils, and one that converts the event into a HepMC3::GenEvent, called generateEventPythia_HepMC. With this, each generated event can be converted to whichever format it is needed in for native ColliderBit LHC search analyses, or for SM measurements using Rivet.This new structure of module functions can be seen in Table 4.
Lastly, some of the most relevant additions to Collid-erBit are module functions to compute the likelihood for LHC SM measurements, using Rivet and Contur.The new module function Rivet_measurements provides 9 With the exception of flavour, which is handled by FlavBit [204].
an eponymous capability and analyses a given HepMC event using Rivet's measurements.Its output, a YODA10 analysis object written as a stream, is then used by the module function Contur_LHC_measurements_from_stream, which provides the capability LHC_measurements and that runs the YODA analysis object through Contur to compute the likelihood.An additional module function, Contur_LHC_measurements_from_file provides the same capability and performs the same computation but starting from a YODA analysis object read from a file.Both of these functions return a Contur_object structure which contains the output from Contur, including the total likelihood value as well as the likelihood value per pool.These respective values are extracted by the module functions Contur_LHC_measurements_LogLike and Contur_LHC_measurements_LogLike_perPool, which provide the capabilities LHC_measurements_LogLike and LHC_measurements_LogLike_perPool, respectively.
For debugging purposes only, there exists an additional module function Contur_LHC_measurements_histotags_perPool, providing the capability LHC_measurements_histotags_perPool that extract the tag of the dominant histogram(s) in each pool.Descriptions of these capabilities and module functions can be found in Table 4. Finally, each of the module functions that use Contur has a partner function that allows running with multiple sets of Contur options simultaneously, and these are called and Multi_Contur_LHC_measurements_LogLike_all.The single variation of the LogLike function is required to provide a single one of the multiple likelihoods to ColliderBit for combination with LHC search likelihoods, whereas the all version -which comes under its own dedicated LHC_measurements_LogLike_Multi capability -can be used to collect and save all the likelihoods for examination after the run.

C.4 New backend interfaces: Rivet and Contur
One the major breakthroughs of this study is the simultaneous combination of LHC searches and measurements in the analysis.This is possible due to the newly developed interfaces to Rivet 3.1.5[205], which contains the extensive library of measurements, and Contur 2.1.1 [79], which takes care of the rigorous statistical combination of the results from individual measurements.
Rivet is written in C++, so its interface to GAMBIT is auto-generated using the package BOSS [77], to allow dynamic loading of Rivet classes.Although some modifications of BOSS were needed to interface Rivet, these will not be documented here, but in a future publication.This interface provides GAMBIT with access to the classes inside of Rivet and some of its global functions.In particular, the only class of interest for the interface is the AnalysisHandler class, which handles the analysis of the passed HepMC events.The usage of this class in ColliderBit was described above.Note that at the time of writing, Rivet 3.1.5and, in particular, the class AnalysisHandler, is not threadsafe.This class must therefore be used within a critical section.In addition to the aforementioned class, GAMBIT also uses the global function addAnalysisLibPath.This function is just called once at the beginning of the scan by the backend initialisation function to inform Rivet of the location of the analysis library.
The interface to the Python package Contur opens the possibility to compute the likelihood for LHC measurements from a YODA analysis object generated by Rivet.Two backend convenience functions were implemented for this purpose: Contur_LogLike_from_file, which reads the analysis object from a YODA file, and Contur_LogLike_from_stream, which reads it from a standard stringstream.Both of these functions provide the same capability Contur_Measurements, but with different signatures, as the former takes just a string with the YODA file path, and the latter a shared pointer to the stream.In addition, both functions take as an argument a list of options for Contur (see [79] for a list of useful options).The value returned by both functions is an object of the class Contur_output, which is a simple class designed to manage the Python dictionary produced by Contur.Lastly, the backend convenience function Contur_get_analyses_from_beam, which provides the capability Contur_GetAnalyses, is used to inform Rivet of the analyses known to Contur.

C.5 Adaptations to Rivet and Contur
Besides the implementation of the interfaces to Rivet and Contur on the side of GAMBIT, minor modifications of each of these packages were necessary to adapt them to the GAMBIT workflow.
Typically, Contur is run on YODA files that have been generated in a separate Rivet run.However, in a high-performance computing environment, the cost of writing and reading from a YODA file at each parameter point would be prohibitive.Therefore, Rivet and Contur were adapted so that the YODA file could be passed between them in memory via stringstream.
In Rivet, this means the AnalysisHandler class received a new overload to the write method to output to stringstream.This is also available in the Python interface, and became available as of Rivet 3.1.4.
The changes to Contur were more significant.The main run function now takes a dictionary of arguments, and if the "YODASTREAM" term contains a Python StringIO -which can be converted to and from a C++ stringstream -then Contur will run on that stream.When run like this, the Contur main will return a dictionary containing various statistics.Which outputs appear in the dictionary is controlled by the "YODASTREAM_API_OUTPUT_OPTIONS" argument, and the options here are summarised by Table 5.These features first appeared in Contur 2.1.1.

4 . 1
Software framework and event generation

∈
(150,200) GeV[108]; ii) a signal region requiring 3 leptons, no opposite-sign, same-flavour lepton pairs and E miss T > 50 GeV [110]; and iii) a signal region requiring ≥ 5 leptons[123].Due to the many different final state combinations of leptons and b-jets that can arise in the decays of 2-4 on-shell and off-shell h, Z and W ± bosons, the best-fit parameter point simultaneously predicts small signal contributions in all of these three signal regions.

1 450
GeV and m χ02 , m χ± 1 800 GeV.In the middle sector of the (m χ0

1 , m χ0 1 )
planes.Here the model obtains positive contributions to ∆ ln L searches from small excesses in leptons + E miss T searches and the ATLAS b-jets + E miss T search (Fig. 6, top right and bottom left).Some examples of the balancing of different branching ratios that these scenarios exhibit, discussed in Sec.5.1, can be seen in the middle left, middle right and bottom left panels of Fig. 7.

1 with m χ0 1 200
GeV, where BR( χ0 1 → γ G) is large.In particular, the SM signature measurements contribute to excluding low-mass scenarios where the constraints from leptons + E miss T searches would otherwise have been largely balanced by positive log-likelihood contributions from the photons + E miss T searches (Fig.6, top panels, m χ0 1 100 GeV).The (m χ0
in two fully-hadronic final states: the qqqq final state arising from W /Z bosons each decaying to light-flavour quarks/antiquarks, and the bbqq final state which arises from a Z or Higgs boson decaying to b b and a W or Z boson decaying to light-flavour quarks/antiquarks.The analysis uses events with at least two large-R jets, and counts the b-multiplicity of each of these jets using a b-tagged track jet procedure.Boosted boson tagging algorithms are then defined to identify various SM boson decays in the two leading large-R jets: W qq (Z qq )-tagging targets W (Z) → qq, whilst Z bb (h bb )-tagging targets Z(h) → bb.V qq -tagging is used to denote the logical OR of W qq -and Z qq -tagging.Signal regions are then defined using the multiciplities of the different boson tags n(W qq ), n(Z qq ), n(V qq ), n(Z bb ) and n(h bb ).Additional background rejection is provided by selections such as a veto on b-jets that do not originate from the boosted boson candidates, lower bounds on the effective mass m eff (defined as the scalar sum of the p T of the two leading large-R jets and E miss T), lower bounds on E miss T , cuts on an event shape variable, and a lower bound on the stransverse mass m T2 constructed from the two leading large-R jets.Our implementation of this search includes the signal regions 4Q-WW, 4Q-WZ, 4Q-ZZ and 4Q-VV.Due to difficulties with reproducing the b-tagging for small radius track discriminating variables H T , E miss T and E miss T / √ H T .The ATLAS search for top squarks in final states with one lepton, jets plus missing transverse momentum [103]: This search (ATLAS_1lep_stop) CMS search for chargino and neutralino production in the WH final state [116]: This search (CMS_1lep_bb) is optimised on simplified models of chargino and neutralino production, with subsequent decays to a lightest neutralino and either a W or Higgs boson.Events are required to have an electron or muon, two b-tagged jets with an invariant mass close to the Higgs boson mass and significant missing transverse momentum.Discriminating variables include the transverse mass of the lepton-neutrino system, the contranverse mass, and the E miss T .The paper defines two signal regions, distinguished by the selection on E miss T

Table 1 :
Ranges and scanning priors for the input parameters.

Table 2 :
The different ATLAS and CMS searches we simulate for our LHC likelihood, with associated short-hand labels.

Table 3 :
Capabilities and module functions added to DecayBit for this study.The Contur_output object.
Runs Contur on YODA file.