HL-LHC sensitivity to higgsinos from natural SUSY with gravitino LSP

In the realm of natural supersymmetric models, higgsinos are typically the lightest electroweakinos. In gauge-mediated supersymmetry breaking models, the lightest higgsino-dominated particles decay into a $Z$-boson or a Higgs boson ($h$), along with an ultra-light gravitino ($\tilde{G}$) serving as the lightest supersymmetric particle (LSP). This scenario suggests a significant non-resonant $hh$ production. Basing on the recent global fitting results of the $\tilde{G}$-EWMSSM (MSSM with light electroweakinos and an eV-scale gravitino as the LSP) performed by the \textsf{GAMBIT} collaboration, which support a higgsino dominated electroweakino as light as 140 GeV, we develop two simplified models to evaluate their detection potential at the high-luminosity LHC (HL-LHC) with $\sqrt{s} = 14~{\rm TeV}$ and an integrated luminosity of $3000~{\rm fb}^{-1}$. The first model examines the processes where heavier higgsino-dominated states decay into soft $W/Z$ bosons, while the second focuses on direct decays of all three higgsino-dominated electroweakinos into $W/Z/h$ plus a $\tilde{G}$. Our study, incorporating both models and their distinct decay channels, utilizes detailed Monte Carlo simulations for signals and standard model backgrounds. We find that the HL-LHC can probe higgsinos up to 575 GeV, potentially discovering or excluding the natural SUSY scenario in the context of a gravitino LSP. Further, we reinterpret this discovery potential using the GAMBIT global fit samples, and find that the entire parameter space of $|\mu| \leq 500~{\rm GeV}$ with an electroweak fine-tuning measure ($\Delta_{\rm EW}$) under 100 in $\tilde{G}$-EWMSSM is accessible at the HL-LHC.

1 Introduction Supersymmetry (SUSY) is a unique non-trivial extension of relativistic invariance and provides a unified description for fermions and bosons.It is also an elegant framework for new physics (NP) beyond Standard Model (SM) in particle physics [1,2].Specifically, the supersymmetric Higgs mass parameter µ should be at the electroweak scale to ensure that the Higgs fields H u and H d get the observed non-zero vacuum expectation value (vev) v ≈ 246 GeV after electroweak symmetry breaking [3].Since this µ parameter is the mass of higgsinos, it can be measured at the Large Hadron Collider (LHC) and dark matter (DM) detection experiments [4][5][6][7][8].Theoretically, if µ is much greater than electroweak scale, some other parameters need to be fine-tuned and the tuning extent is reasonably quantified by a dimensionless measure [4,[9][10][11][12]] where m Z is the Z-boson mass and p i are the SUSY parameters at electroweak scale.
Of our special interest is the SUSY model with light electroweakinos and not so heavy top squarks (stops), which is called natural SUSY [4,9,10,.Natural SUSY is well motivated both theoretically and experimentally.Firstly, natural SUSY is supported by the Higgs boson properties, which offer compelling perspectives for testing the effects of SUSY.Ref. [37] reviewed the SUSY corrections to the couplings and decays of the SM-like Higgs boson and their dependence on the 19-dimensional parameter space in the phenomenological minimal supersymmetric standard model (MSSM), which showed that the SM-like Higgs state discovered at the LHC can still be described as a SUSY-like Higgs state.Secondly, the electroweakinos and the stops can still be light [38].The searches at the current LHC experiment have pushed the detection limits on the colored sparticles to about 2 TeV, mostly independent on their decay modes and the mass differences with other sparticles.The bounds on third generation squark masses are somewhat weaker, with the detection limit on the stop mass is around 1.2 TeV.About the electroweakinos, the current LHC limits are rather weak.The searching limits on wino-like states are about 1 TeV for a massless lightest supersymmetric particle (LSP), while for higgsino-like states the limits are only about 550 GeV.However, the current stringent LHC constraints on the SUSY particles strongly depend on the simplified models used in the result interpretation [39].A detailed survey [40] of the combined impacts of DM and collider constraints on the electroweakino sector in the MSSM (EWMSSM), mostly using the 36.1 fb −1 data, shows that the µ parameter around 200 GeV is still allowed 1 , mainly due to the compressed mass spectrum.In the framework of natural SUSY, the LHC bounds on the higgsinos are even much weaker since they are quite compressed with the higgsino-like LSP and their productions merely give missing energy [43][44][45].
In the theoretical view, the µ-parameter is indeed predicted at electroweak scale in some SUSY modles like supergravity (SUGRA).Since in SUGRA the gravitino mass m 3/2 is also expected to be above 100 GeV (the gravitino is expected to acquire mass through the super-Higgs mechanism) with ⟨F ⟩ being the SUSY breaking scale and M P = (8πG N ) −1 = 2.4 × 10 18 GeV being the reduced Planck mass, the higgsinos may be lighter than the gravitino and the LSP may be higgsino-like.Because in SUGRA the gravitino couples to every sparticle with a universal coupling in form and its interaction is of gravitational strength, the gravitino has no collider phenomenology.In contrast, in gauge-mediated SUSY breaking (GMSB) models, the gravitino is expected to be much lighter than other sparticles and hence serve as the LSP.
In the presence of an approximately massless gravitino LSP, the GAMBIT collaboration reinterpreted the collider data in the EWMSSM (the MSSM with light electroweakinos), referred as the G-EWMSSM [88].In such a scenario, it is found that the µ-parameter can still be about 100 GeV, with a value of 170 GeV for the best point, indicating a small fine-tuning measure ∆ EW2 .
In the G-EWMSSM, the electroweakinos are produced in pair at the LHC and subsequently decay into final state of an electroweak boson Z/W/h/γ plus a gravitino G or plus a lighter neutralino ( χ0 1 mostly).For these decay modes, the most difficult channel to detect is the final state of di-Higgs plus missing energy.The ATLAS report [90][91][92][93][94] shows that there exist a gap with higgsino mass m χ0 1,2 = m χ± 1 around the top quark mass m t , which is unreachable at LHC Run-II (see Fig. 21 in [38]).Although the di-Higgs channel is the most difficult to detect, it is one of the most important processes to be studied at the LHC since it allows for a measurement of the Higgs trilinear self-coupling λ hhh and for the search of some NP scenario [95][96][97].For such a di-Higgs production, the decay channels are bbW W [98,99], four b-jets [100,101], bbγγ [102,103], bbτ τ [104] and four W -bosons [105].Each channel has its pros and cons.The four b-jets channel and the bbγγ channel are more frequently used than other channels.The former suffers from the challenge of separating it from other processes produced in the LHC collisions.Although photons are easy to pick out amongst other processes, the bbγγ channel is suppressed by its very small branching ratio.According to the current ATLAS and CMS results, the bbτ τ is a rather promising channel, which has a higher production rate than the bbγγ channel and lower background than four b-jets channel.As for the benchmark NP model like G-EWMSSM, the bbτ τ channel can also be used to probe the process such as pp → χ χ → hh G G, which is expected to fill the unreachable gap in the Fig. 6 of Ref. [94].
In the spirit of leaving no stone unturned, in this work we focus on the G-EWMSSM scenario and study the HL-LHC observability.The content is organized as follows.In Sec. 2, we briefly describe the G-EWMSSM scenario, i.e., the MSSM electroweakino sector with a light gravitino.In Sec. 3, we study the HL-LHC sensitivity via scrutinizing the signature of higgsino pairs in five complementary decay channels, including the final states of ZZ, Zh, hh, W Z and W h plus missing energy E miss T .The prospects of HL-LHC searching limits for the simplified models are summarized in Sec. 4. Then we reinterpret the result in the framework of G-EWMSSM with the help of the global-fit data released by GAMBIT collaboration, and discuss the impact on the fine-tuning measure ∆ EW .In Sec. 5, we make a summary.

Theoretical preliminaries
The gravitino, with a mass generated in Eq. (1.2), couples universally to all matter and all force fields.Its relevant interactions can be described by the goldstino field G [106-108] where (ψ, ϕ) is the scalar and fermionic component of the chiral supermultiplets, λ a is the gaugino associated with the gauge field with field strength F a µν .In the MSSM the electroweak gauginos (including bino B and winos W ) mix with the higgsinos ( H0 u , H+ u , H− d , H0 d ).The neutral states form four neutralinos χ0 i , and the charged states form two pairs of charginos χ± i .In the basis of the neutralino and chargino mass terms are given by where the symmetric neutralino mass matrix M N takes a form with M 1 and M 2 being the soft breaking masses of bino and wino fields respectively, g ′ and g being SU (2) L and U (1) Y gauge couplings, and v u and v d being the vevs of Higgs fields Similarly, the chargino mass matrix M C is given by (2.5)

The light higgsino scenario with a light gravitino LSP
For the naturalness of SUSY, the higgsino mass parameter µ can not be too heavy, which can be quantified by ∆ EW defined in Eq. (1.1).Specifically, at the one-loop level the Z-boson mass is given by [109] where m 2 Hu and m 2 H d are the soft SUSY breaking masses of the Higgs fields H u and H d , Σ u and Σ d represent the radiative corrections.Therefore, from Eq (1.1), one can get (2.7) The fine-tuning ∆ EW µ of less than one percent level indicates µ in the range of 100−500 GeV.Here it should be noted that the fine-tuning measure can be changed by the SUSY breaking boundary condition and the radiative correction.For example, in the constrained MSSM with non-universal Higgs and gaugino masses, ∆ EW µ can be reduced by a factor 2 [110].Moreover, the parameter space with the square root of the stop masses √ m t1 m t2 greater than 1.2 TeV and/or the gluino mass heavier than 1.8 TeV can also contribute a 5% finetuning (∆ EW ≳ 20) [13].In this work, we focus on the light higgsinos in the simplified model defined in Sec.2.3, and neglect the colored sparticles [111] by fixing their masses at 3 TeV.
In the light higgsino scenario, the lightest three electroweakinos ( χ0 1,2 and χ± 1 ) are higgsino like, which are nearly degenerated, with masses around the value of µ.In the following, the MSSM-like heavy Higgs states H, A, H ± are also assumed much heavy, and the SM-like h acts as the 125 GeV Higgs state.
In addition to the decay modes in the MSSM, the interaction in Eq. (2.1) leads to all electroweakinos decay into the final states containing a gravitino.In the mass eigenstates, the decay widths of the neutralinos into neutral gauge bosons are [112][113][114] where depict the components of bino, wino and higgsinos in the neutralino χ0 i .Here N ij is the neutralino mixing matrix element with (i, j) being the (mass, gauge) eigenstate labels, θ W is the Weinberg angle, and α is the mixing angle between the CP -even neutral Higgs states.Similarly, the decay width of charginos into gravitino are given by where with U and V being the chargino rotation matrices defined by U XV = m diag χ± .In terms of a typical decay width of electroweakino into gravitino then the decay widths in Eqs.(2.8), (2.9), (2.10) and (2.12) can be written as (2.15) For a heavier electroweakino χ0 i its decay width to LSP G is enhanced by the quintic factor of the electroweakino mass; while for a heavier LSP G, the decay width is suppressed by m 2 3/2 .

The electroweak MSSM with a light gravitino LSP
The scenario of G-EWMSSM has been encoded in the package GAMBIT, which provides a global fit for the parameter space [88,115].Here we show the result in light of current experiments, focusing on the survived samples3 .In the GAMBIT work, the input parameters of the G-EWMSSM electroweak sector are usually chosen as  We plot the components and the masses of the two lightest neutralinos χ0 1,2 in Fig. 2, and plot the decay widths of χ0 2 ( χ± 1 ) → χ0 1 + X versus χ0 2 ( χ± 1 ) → G + X in Fig. 3, in which the color is coded with the mass splitting.From Fig. 2 one can find that when the χ0 1 mass is smaller than 400 GeV, the current experiments favor the higgsino NLSP scenario.There is also a small part of samples predicting χ0 1 to be bino-dominated or bino-higgsino admixture.Pure higgsinos do not decay to photons, which is verified by the fit result.The decay width of χ0 2 (or χ± 1 ) to χ0 1 is suppressed by the mass splitting between higgsinodominated states.For the most compressed case, in which M 1 and M 2 are much heavier than µ, the total decay width Γ χ0 2 → χ0 1 + X can reach 10 −20 GeV, corresponding to a very small mass splitting around 10 −3 GeV.For a mass splitting at order 1 GeV, the decay width of χ0 2 into all final states including χ0 1 is at the same order of Γ 0 defined in Eq. (2.14).

Two simplified models of light higgsino scenario
We consider two simplified SUSY models inspired by gauge-mediated SUSY breaking and by the GAMBIT result.Both of them contain an almost mass-degenerate higgsino system ( χ± 1 , χ0 1 , χ0 2 ) and an nearly massless LSP G, while all other SUSY particles are much heavy and out of the reach of the LHC/HL-LHC.The higgsino system offers four production processes at the LHC: χ+ 1/2 and χ0 1 χ0 2 .Fig. 4 shows the cross sections of chargino-neutralino pair productions for √ s = 14 TeV at NLO plus next-leading-log (NLO-NLL) contributions using the package Resummino [116][117][118] with the parton distribution  1 .The data is taken from Ref. [88,115].
functions (PDFs) given by CTEQ6.6 [119].In the calculations we assume the higgsino-like charginos and neutralinos are degenerate, denoted as m H .

Simplified model A
In the first simplified model, the gravitino mass m 3/2 is set to 0.1 GeV, and the χ± 1 and χ0 2 masses are set to 1 GeV above χ0  Cross section

[pb]
pp, 14 TeV, NLO+NLL The cross sections for higgsino pair productions, assuming χ0 1 , χ0 2 and χ± 1 with a degenerate mass m H .The four production modes χ+ 2 are shown respectively by the orange, red, cyan and green curves, while the black curve represents the sum of these modes.as two free parameters in our analysis for this simplified model.Therefore, there are three searching channels complementary to each other: (1) HH channel.This channel is targeting the Higgs pair from the χ0 1 decay.As discussed in Sec. 1, we focus on the final states of hh decaying into b bτ τ , and a typical process is shown in Fig. 5 (a).
(2) ZH channel.As shown in Fig. 5 (b), this channel focuses on the final state containing two b-jets plus one opposite-sign same-flavor (OSSF) lepton pair.
(3) ZZ channel.This channel is targeting the final state of four-leptons which can form two OSSF pairs, as shown in Fig. 5 (c).
All production processes of the higgsino system (black line in Fig. 4) are considered in this model.

Simplified model B
The second simplified model is similar to the first one, except the gravitino mass m 3/2 is set to 1 eV.Then both χ± 1 and χ0 2 dominantly decay directly into G plus a W ± or Z/h boson.In this simplified model, only chargino-neutralino pair productions ( χ± 1 and χ± 2 ) are considered in our analysis.So this simplified model can be probed via: (5) WZ channel.This channel is well probed in the previous LHC searches [41,42,120], and a similar search strategy focusing on 3ℓ + E miss T final states is used in this work, which is shown in Fig. 5 (e).
Since in the light higgsino scenario we have |µ| ≪ M 1 , M 2 , the bino and wino components of χ0 1 and χ0 2 can be neglected,

which results in a maximal mixing between H0
u and H0 d in χ0 1 and χ0 2 .Therefore, given the complicated model dependence on the parameters like tan β and α, the branching ratio of the lightest neutralino χ0 1 is treated as a free parameter, with BR χ0 (2.17) These relations provide a good approximation according to Eq. (2.9) and Eq.(2.10) in the decoupling limit m A ≫ m Z [121].As a result, the total branching ratio of the WH channel or WZ channel is equal to 50%.The corresponding cross section is shown by dark blue curve in Fig. 4.

Detecting higgsino pair productions at HL-LHC
Now we present our search strategies, especially the distinctive new signature from the di-Higgs events in the bbτ τ final state.As discussed in the preceding section, enhancing this channel can significantly enhance the sensitivity compared to the leading four b-jets signal for the di-Higgs events.Several relevant studies have explored the electroweakino productions with a light gravitino scenario [122][123][124][125][126][127].In this work, we conduct a comprehensive investigation for the decays of the higgsino pairs, which are subsequently combined with all other channels to provide a final assessment of the LHC sensitivity to the µ parameter.We perform Monte-Carlo simulations to estimate both the SM backgrounds and the higgsino pair signals, as well as to examine the efficiency for various electroweakino productions.In our analysis, all events are generated using Madgraph5 aMC@NLO event generator [128,129] and we use Pythia8 [130,131] for parton shower and hadronization.Nextto-leading-order (NLO) cross sections are used for background and signal normalization, calculated using Madgraph5 aMC@NLO [128] and the package Resummino [117], respectively.The events are processed through the fast ATLAS detector simulation using the Rivet-HEP package [132] for parameterized simulation and object reconstruction with the ATLAS RUN2 smearing and efficiency setting [133,134].
We reconstruct jets using the anti-k T clustering algorithm [135] with a distance parameter of 0.4, as implemented in the FASTJET package [136].The jets are required to have p T > 30 GeV in the region |η| < 2.47.The b-jets containing b-hadron decays are identified from the selected jet objects with a flat efficiency of 75% and the mis-tagging rates for non-b jets with ϵ(c → b) = 0.1 and ϵ(j → b) = 0.01 [137].The reconstructed τ had objects are seeded by jets, with an efficiency fixed at 60% [133].Light flavor lepton objects (electrons and muons) are reconstructed in the region |η| < 2.47 to have p T > 15 GeV.The missing transverse momentum p miss T , with magnitude E miss T , is defined as the negative vector sum of the transverse momenta of all identified physics objects, including electrons, photons, muons, jets and an additional soft term.
We apply a sequential overlap removal procedure to resolve ambiguities in which multiple electron, muon or jet candidates would otherwise be reconstructed from the same detector signature.This procedure uses the definition of the angular distance, ∆R = (∆y) 2 + (∆ϕ) 2 , where y is the rapidity of the object, and the following steps are applied in order: 1.If a b-jet is within ∆R = 0.2 of an electron or muon candidate, this electron or muon is rejected, as it is likely to be from a semi-leptonic b-hadron decay; 2. If the jet within ∆R = 0.2 of the electron is not b-tagged, the jet itself is discarded, as it likely originates from an electron-induced shower; 3. Any τ had object with p T ≤ 50 GeV within ∆R = 0.2 of an electron or muon is removed.Otherwise, the muon or electron is removed; 4. Any electron with ∆R = 0.2 of a muon is removed; 5. Any jet within ∆R = 0.2 of an electron or muon is removed; 6. Any jet with ∆R = 0.2 of a τ had is removed; 7. Muons with ∆R = 0.4 of a remaining jet are discarded to suppress muons from semi-leptonic decays of b-hardons.
The detection and discovery significance Z A is evaluated by Poisson formula where s and b are respectively the signal and background event numbers.

Searching strategy for HH channel
The HH channel targets events with two oppositely charged τ had candidates and two b-jets.The main sources of backgrounds in HH channel are t-quark processes, Z+jets, W +jets, di-bosons and multi-jet productions.The t-quark processes include the single-t process, t t and t tX, where X stands for a Z or H boson.The leptonic t t background event is complicated by an additional layer of decays, with t → W b in the first step followed by W → ℓν in the second.This results in an event topology with two identical branches, t → bℓν, each with a visible (bℓ) and invisible (ν) component.In order to suppress the background from the t t processes, the 'stransverse mass' m T2 [138,139] is an event variable used to bound the masses of an unseen pair of particles, defined as where the transverse mass m T is defined as with p T and q T being respectively the transverse momenta of the visible and invisible components, and p T and q T are their magnitudes correspondingly.In our analysis, for the pairing of two b-jets and two leptons, we select the pair that yields the smaller value for m T2 , marked as m min T2 .Several other variables used in this channel are the invariant mass of two leptons m ℓℓ , the invariant mass m bb and the angular distance ∆R bb of two b-jets.The distributions of the variables E miss T , ∆R bb , m bb and m min T2 are shown in Fig. 6.For the signal events, the following cut sequences are summarized in Table 1: • The basic cut requires two reconstructed b-jets plus two hadronic τ had candidates.
The event containing any other jet or lepton candidate will be vetoed.
• • The m min T2 > 175 GeV cut.The distribution of m T2 variable has an upper endpoint at the parent particle.The cut at 175 GeV here is to further suppress the background.
After the cut-flows, we find that the significance Z A of the signal point of m H = 600 GeV can reach about 2.4σ.

Searching strategy for ZH channel
The ZH channel targets the event containing exactly two b-jets and two leptons (e or µ) that can form one OSSF pair.The dominant backgrounds are t t, t tX, bbℓℓ, VH process, Z+jets and di-bosons.
The search strategy in the ZH channel is to reconstruct Z and h with large E miss T .As summarized in Table 2 and based on the distributions of E miss T , m ℓℓ , ∆R ℓℓ and m min T2 in Fig. 7, both the dominant background and the benchmark signal µ = 600 GeV are selected with the following cuts: • The basic cut requires the signal event to contain two b-jets and an OSSF lepton pair.Events containing another jet or an additional lepton are rejected.• The E miss T trigger.As shown in Fig. 7, the background basically follows an exponentially decreasing distribution.Here E miss T is required to be larger than 100 GeV.
• The Z and H mass peak cut.The invariant masses of the OSSF pair m ℓℓ and two b-jets of the signal event are required to located in the Z peak window 80 GeV ≤ m ℓℓ ≤ 100 GeV and in the Higgs window 90 GeV ≤ m bb ≤ 150 GeV, respectively.
• The OSSF pair cone ∆R ℓℓ cut.The Z → ℓℓ candidate of the signal event appears as a smaller cone compared with the SM backgrounds, so the cut ∆R ℓℓ < 1.5 greatly improves the statistical significance Z A .
• The m T2 cut.In order to further suppress the t-quark background, a cut of m min T2 ≥ 175 GeV is also applied in the ZH channel.
After the cuts, the signal for m H = 600 GeV and BR( χ0  1, but for ZH channel.The label '-' represent that no event are retained after applying the corresponding cut.For the signal we choose a benchmark point m H = 600 GeV, BR χ0 1 → Z G = 0.5, marked as S(600,0.5).

Searching strategy for ZZ channel
The ZZ channel targets the final states containing two OSSF pairs of leptons in Z peak plus large E miss T .The dominant irreducible backgrounds are the processes that can produce four charged leptons, including ZZ (off-shell Z/γ contributions are also included), t tZ, V V Z (V = W/Z), and the Higgs production via gluon-gluon fusion.
As shown in Fig. 8 and summarized in Table 3, for this channel we apply the following cuts sequentially: • The basic cut requires the signal event to have exactly two light flavor (e or µ) OSSF lepton pairs.Events with additional signal lepton candidate and any b-jet are discarded.
• The OSSF pair with mass closer to the Z-boson mass is labeled as the first Z can-  8, the narrower peak of the first Z candidate is due to the ordering of the Z candidates, so that widening the range of the m 2nd OSSF window increases the signal acceptance and the statistical significance.
• The final cut requires E miss T ≥ 100 GeV.Note that the distribution of E miss T depends on the value of µ.Although for parameter points with a relatively large µ value, e.g. the benchmark point in Fig. 8 with µ = 600 GeV, it seems to be possible to take a stronger cut on E miss T (e.g.E miss T ≥ 150 GeV) to suppress backgrounds more efficiently, a stronger cut on E miss T would hurt signal seriously for a lower value of µ.So we take a universal soft cut here for simplicity.
After the selection criterion, the detection sensitivity to BR χ0 1 → Z G of m H = 600 GeV can reach to about 0.3.

Searching strategy for WH channel
In order to search for the chargino-neutralino productions shown in Fig. 5(d  signal event requires a single isolated electron or muon, two b-jets and large E miss T .The major backgrounds arise from SM processes containing top quarks and W -bosons.As shown in Fig. 9, the kinematic variables E miss T , the invariant mass of the two b-jets m bb , the cotransverse mass variable m CT and the m T variable of lepton p T show great difference between the backgrounds and signals.The SR-WH is defined to suppress the SM backgrounds via the following cut strategies: • The basic cut requires exactly one light flavor lepton and two b-jets.The primary SM processes that contribute to the preselection region are t t, single top quark (mostly in the tW channel), W +jets and W Z production.Events with additional leptons or jets are rejected.
• The cut E miss T ≥ 100 GeV can suppress all SM backgrounds by about one order of magnitude after the basic cut.
• The invariant mass of the two b-tagged jets is required in the range of 90 GeV ≤ m bb ≤ 140 GeV, consistent with the mass of the SM Higgs boson.
• The angular distance of the two b-jets ∆R bb is required to be smaller than 1.5, in order to tag a SM Higgs fat cone.
• The cotransverse mass variable m CT of the b-jet pair is defined as where p b1 T and p b2 T are the magnitudes of the transverse momenta of the two b-jets and ∆ϕ bb is the azimuth angle between the two b-jets [140].For the t t background, the m CT variable has a kinematic endpoint of about 150 GeV.SR-WH requires a harder cut m CT > 200 GeV to effectively suppress the t t and tW backgrounds.
• The SM processes with one W -boson decaying to leptons, originating primarily from the semi-leptonic t t and W +jets background, are suppressed by requiring the transverse mass m T to be greater than 150 GeV in SR-WH.
As shown in Table 4, for the signal with µ = 600 GeV, the statistical significance Z A can reach to about 3.5σ in SR-WH.

Searching strategy for WZ channel
For this channel the final state with three leptons is known as the golden channel for electroweakinos.Since in our study we tend to make a conservative evaluation for the exploration of the simplified model depicted as in Fig. 5 (e), we use a common approach throughout most steps of our analyses as summarized in Table 5: • The basic cut requires exactly three leptons that can form an OSSF pair.In the cases of three leptons with same flavor, among the two combinations of the OSSF lepton pair, the combination with a smaller m T corresponding to the third lepton is selected to reconstruct the event.Events containing any b-jet are discarded.
• E miss T is required to be greater than 100 GeV.
• In order to suppress the non-resonance background, such as t t, the invariant mass of the OSSF lepton m ℓℓ is required to be in the Z-boson mass window [80,100] GeV.
• The angular distance between the two leptons in the OSSF lepton pair ∆R ℓℓ is required to be smaller than 1.5 to search for a boosted Z-boson from a heavier higgsino decay.
• As shown in Fig. 10, the m T variable of the third lepton for the SM background has a Jacobian peak with a sharp cut-off at the W -boson mass, while the signals have relatively flat distributions.Therefore, the cut m 3rd−ℓ T being greater than 110 GeV can greatly improve the statistical significance.
For m H = 600 GeV, the statistical significance via this WZ channel can reach to about 10σ.All the SRs defined in this work and the corresponding SM backgrounds are summarized in Table 6.The first three SRs, i.e., SR-HH, SR-ZH and SR-ZZ, target the first Table 6.Selection criteria and the event numbers of all SM backgrounds simulated in signal regions (SRs).The variables have mass-dimension in GeV.The "-" means that the corresponding background is negligible.From Fig. 11 (a), one can find that the SR-HH, SR-ZH and SR-ZZ channels are complementary to each other for probing a higgsino pair, and the parameter space with the higgsino mass m H less than 575 GeV can be fully covered regardless of the branching ratio.The exclusion potential of the higgsino mass m H can reach about 960 GeV using the ZZ channel with the assumption that BR( χ0 1 → Z G) = 1.And in the situation that BR( χ0 1 → h G) = 1, the HH channel can push m H to about 620 GeV.Searching di-Higgs in bbτ τ plus E miss T final state can distinguish the signal from t t background well.Meanwhile, the ZH and ZZ channel can reach the branching ratio of a higgsino into h G to about 5% and 1% correspondingly.

Results and discussions
From Fig. 11 (b) we see that the exclusion limits from the WH and WZ channel can reach to about 550 and 650 GeV, respectively.The statistical sensitivity of the WZ channel is better than the WH channel.
Roughly speaking, the detection limits of the channels containing Higgs bosons are a little weaker than the channels containing Z-bosons.The Higgs tagging techniques, as proposed in recent works [141][142][143][144][145][146][147][148][149] using the machine learning methods, can improve the detection ability of the WH, HH and ZH channels.
Finally, for the GAMBIT global fit result with the current LHC data for the G-EWMSSM, we reanalyse the samples to figure out their accessibility at the HL-LHC.As shown in Fig. 12, our search strategies can greatly probe the parameter space which survives the current LHC constraints.Comparing to Fig. 1, one can see that almost the whole parameter space with |µ| < 500 GeV will be covered by the HL-LHC.The fine-tuning measure ∆ EW , defined in Eq. (1.1), is required to be larger than 100, which indicates a fine-tuning worse than 1%.There are still a small group of samples with −µ ≲ M 1 ≃ 125 GeV that are not accessible at the HL-LHC.As shown in Eq. (2.9) and Eq.(2.10), the decay widths of χ0 .Other surviving samples are featured with a bino-dominated or wino-dominated χ0 1 , and its dominant decay mode is χ0 1 → γ + G, which corresponds to a decay branch ratio greater than 80%.Therefore, the remaining samples can be probed via the process in final states containing 2γ + E miss T at the HL-LHC.

Summary
In this paper, we studied the HL-LHC sensitivity to the production of higgsinos in the framework of the EWMSSM assuming the LSP to be the ultra-light gravitino ( G).We considered two simplified models: • In the first one, both the charged higgsinos χ± 1 and the second neutral higgsino χ0 2 decay to the lightest higgsino χ0 1 plus some invisibly soft particles.In this case we have an enhanced number of lightest higgsino pairs with each pair giving three decay channels χ0 1 → G G + ZZ/Zh/hh.
• In the second model, the charged higgsino decays directly to a gravitino plus a Wboson and the second neutral higgsino decays directly to a gravitino plus a Z/h boson.So the production of χ± In each model the different channels are complementary and they were all jointly considered in our study.Through detailed Monte Carlo simulations for the signals and the standard model backgrounds, we found that the HL-LHC can probe the higgsinos up to about 600 GeV.Therefore, for such a MSSM with gravitino LSP, the HL-LHC can either discover a light higgsino or exclude the natural SUSY scenario in case of gravitino LSP.However, we should note that for the natural SUSY scenario with higgsino LSP, a hadron collider is rather weak for exploring higgsinos, and, instead, we need a lepton collider which allows for a sensitive probe for compressed higgsinos [150,151].
We also reinterpreted our result with G-EWMSSM global-fit data set.We found that the HL-LHC can probe almost the whole parameter space with |µ| ≤ 500 GeV, and push up the electroweak fine-tuning factor ∆ EW to 100 in case of unobservation.

. 16 )
while all other soft SUSY breaking parameters are fixed at 3 TeV, except m 3/2 which is fixed at 1 eV.Fig.1shows the samples in light of current experimental data, plotted in the M 1 − µ and M 2 − µ planes with color coded in ln L. From Fig.1one can find that the best-fit region (the points marked with red or yellow colors) favors |µ| to be smaller than about 500 GeV and smaller than |M 1 | and |M 2 |.

Figure 1 .
Figure 1.The G-EWMSSM samples in light of current experimental data, plotted in the M 1 − µ plane (left) and on M 2 − µ plane (right).The color indicates the log-likelihood ratio.The data is taken from Ref. [88, 115].

Figure 2 .
Figure 2. Triangular presentations of the compositions of χ01 (left) and χ0 2 (right) with color coded with their masses.The data is taken from Ref.[88,115]

Figure 5 .
Figure 5.The typical processes predicted by the simplified models for the productions of higgsinodominated electroweakinos with a light gravitino LSP.The neutralinos and charginos are pair produced, resulting in final states with rich leptons, b-jets and E miss T .

( 4 )
WH channel.As shown in Fig. 5 (d), this channel is targeting the final state of the W -boson leptonic decay and the h-boson b b decay.

Figure 6 .
Figure 6.The E miss T , ∆R bb , m bb and m min T2 distributions of the dominated SM backgrounds and the signal with m H = 600 GeV after the HH channel preselection criteria.

•
Higgs jet cut.The signal event contains a Higgs boson H → b b.Therefore, the invariant mass m bb and the angular distance ∆R bb of the two b-jets requires 90 GeV ≤ m bb ≤ 140 GeV and ∆R bb < 2 to suppress the non-Higgs resonance b-jets background like t t.

Figure 7 .
Figure 7.The E miss T , m ℓℓ , m bb , ∆R ℓℓ and m min T2 distributions of the dominated SM backgrounds and the signal with m H = 600 GeV after the ZH channel preselection criteria.

Figure 8 .
Figure 8.The m 1st OSSF , m 2nd OSSF and E miss T distributions of the dominated SM backgrounds and the signal with m H = 600 GeV after the ZZ channel preselection criteria.

Figure 9 .
Figure 9.The E miss T , m bb , m CT and m T distributions of the dominated SM backgrounds and the signal with m H = 600 GeV after the WH channel preselection criteria.
), the searching strategy targets the W -boson decay to lepton and the Higgs boson decay to b b.Then the

Figure 10 .
Figure 10.The E miss T , m ℓℓ , ∆R ℓℓ and m 3rd−ℓ T distributions of the dominated SM backgrounds and the signal with m H = 600 GeV after the WZ channel preselection criteria.

Figure 11 .
Figure 11.Model dependent exclusion contour at 95% confidence level on the higgsino with light gravitino scenario at HL-LHC.(a) Exclusion limit for the first simplified model in Sec.2.3.1, as a function of the lightest higgsino mass m H and the branching ratio BR χ0 1 → Z G = 1 − BR χ0 1 → H G .(b) The exclusion limit for the second simplified model in Sec.2.3.2, as a function of the signal cross sections that reproduce the expected excesses in each SR and the higgsino mass m H .

Figure 12 .
Figure 12.Same as Fig.1, but the gray samples (containing almost all red samples with large likelihood in Fig.12) indicate they can be covered at the HL-LHC.

1 into
Z/h are greatly suppressed by a factor of 1 − m 2 h/Z /m 2 χ0 i 4

Table 1 .
The cut flows of HH channel at the 14 TeV HL-LHC.For the signal we choose a benchmark point m H = 600 GeV, BR χ01 → Z G = 0, marked as S(600, 0).The event number are normalized to 3000 fb −1 .

Table 2 .
Similar to Table

Table 3 .
Similar to Table.1, but for the ZZ channel.For the signal we choose a benchmark point m H = 600 GeV, BR χ0 1 → Z G = 1, marked as S(600,1).

Table 4 .
Similar to Table1, but for WH channel.For the signal we choose a benchmark point m H = 600 GeV, marked as S(600).

Table 5 .
Similar to Table1, but for WZ channel.The label '-' represent that no event are retained after applying the corresponding cut.For the signal we choose a benchmark point m H = 600 GeV, marked as S(600).