Search for supersymmetry in final states with two oppositely charged same-flavor leptons and missing transverse momentum in proton-proton collisions at s\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \sqrt{s} $$\end{document} = 13 TeV

A search for phenomena beyond the standard model in final states with two oppositely charged same-flavor leptons and missing transverse momentum is presented. The search uses a data sample of proton-proton collisions at s\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \sqrt{s} $$\end{document} = 13 TeV, corresponding to an integrated luminosity of 137 fb−1, collected by the CMS experiment at the LHC. Three potential signatures of physics beyond the standard model are explored: an excess of events with a lepton pair, whose invariant mass is consistent with the Z boson mass; a kinematic edge in the invariant mass distribution of the lepton pair; and the nonresonant production of two leptons. The observed event yields are consistent with those expected from standard model backgrounds. The results of the first search allow the exclusion of gluino masses up to 1870 GeV, as well as chargino (neutralino) masses up to 750 (800) GeV, while those of the searches for the other two signatures allow the exclusion of light-flavor (bottom) squark masses up to 1800 (1600) GeV and slepton masses up to 700 GeV, respectively, at 95% confidence level within certain supersymmetry scenarios.

is assumed to decay exclusively into a Z boson, while in the latter, the ZH final state is accompanied by the ZZ final state with 50% branching fractions of the χ 0 1 decaying into an H or a Z boson. Only ZH and ZZ final states are taken into account in the analysis, since the contribution of the HH topology to our signal regions is expected to be negligible. Such models predict the SUSY particles to be produced via EW interactions, with limited if any production of accompanying quarks in the final state.  Figure 2. Diagram for GMSB gluino ( g ) pair production (left), where each g decays into a pair of quarks and a neutralino. The neutralino then decays to a Z boson and an LSP. Diagrams for sbottom b (center) and squark q (right) pair production are also shown. Such models feature a mass edge from the decay of a χ 0 2 via an intermediate slepton, . In the central diagram, a pair of b quarks is present in the final state. In these models we assume a fixed χ 0 1 mass of 100 GeV, while the mass of the slepton is taken to be equidistant from the masses of the two neutralinos. Only the lightest b mass eigenstate, b 1 , is assumed to be involved in the models considered. All these models assume strong production of SUSY particles and predict an abundance of quarks in the final state. This paper is organized as follows. Section 2 provides a brief description of the CMS detector, while section 3 describes the datasets, triggers and object reconstruction in CMS. Section 4 describes the event selection criteria and the SRs used in the search, while the estimation of the SM background contribution is described in section 5. Section 6 describes the fit to the m distribution, used to extract a possible edge-like signal. The results of the search are described in section 7, and are interpreted in terms of constraints on the cross sections of the SMS models, as described in section 8. Finally, a summary of the analysis is given in section 9.

JHEP04(2021)123 2 The CMS detector
The central feature of the CMS apparatus is a superconducting solenoid of 6 m internal diameter, providing a magnetic field of 3.8 T. Within the solenoid volume are a silicon pixel and strip tracker, a lead tungstate crystal electromagnetic calorimeter (ECAL), and a brass and scintillator hadron calorimeter (HCAL), each composed of a barrel and two endcap sections. The tracker system measures charged particles within the pseudorapidity range |η| < 2.5. The ECAL is a fine-grained hermetic calorimeter with quasi-projective geometry, and is segmented into the barrel region of |η| < 1.48 and in two endcaps that extend up to |η| < 3.0. The HCAL barrel and endcaps similarly cover the region |η| < 3.0. Forward calorimeters extend the coverage up to |η| < 5.0. Muons are measured and identified in the range |η| < 2.4 by gas-ionization detectors embedded in the steel fluxreturn yoke outside the solenoid. A two-tier trigger system selects events of interest for physics analysis. The first level of the CMS trigger system, composed of custom hardware processors, uses information from the calorimeters and muon detectors to select the most interesting events in a fixed time interval of less than 4 µs. The high-level trigger processor farm further reduces the event rate from around 100 kHz to about 1 kHz, before data storage. A more detailed description of the CMS detector and trigger system, together with a definition of the coordinate system used and the relevant kinematic variables, can be found in refs. [29,30].

Data, triggers, and object reconstruction
We use events containing at least two OS leptons (e + e − , µ + µ − , or e ± µ ∓ ). Only SF leptons (e + e − or µ + µ − ) are used to define SRs, while e ± µ ∓ events are used in control regions (CRs). These events are preselected using dilepton triggers that require the leptons with the highest (next-to-highest) transverse momentum p T to pass respective thresholds ranging from 17-23 (8)(9)(10)(11)(12) GeV, depending on the data taking period and lepton flavor. In addition, these triggers require the leptons to pass isolation criteria. To retain high efficiency for highly boosted topologies that contain nearly collinear lepton pairs, we also use a second set of dilepton triggers with higher respective p T thresholds of 25-37  GeV but without any isolation requirement. Trigger efficiencies are measured in events selected using triggers based on the p miss T and found to be 85-95%. In addition, a γ+jets event sample is used as a CR to estimate the Drell-Yan (DY) background (as discussed in section 5). This sample is collected using a set of photon triggers with p T thresholds ranging between 50 and 200 GeV. A subset of these photon triggers with lower p T thresholds are prescaled to keep the rate under control. Events collected with prescaled triggers are reweighed accordingly.
The particle-flow algorithm [31] aims to reconstruct and identify each individual particle in an event, with an optimized combination of information from the various elements of the CMS detector. The energy of photons is obtained from the ECAL measurement. The energy of electrons is determined from a combination of the electron momentum at the primary interaction vertex as determined by the tracker, the energy of the correspond-

JHEP04(2021)123
ing ECAL cluster, and the energy sum of all bremsstrahlung photons spatially compatible with originating from the electron track. The momentum of muons is obtained from the curvature of the corresponding track. The energy of charged hadrons is determined from a combination of their momentum measured in the tracker and the matching ECAL and HCAL energy deposits, corrected for the response function of the calorimeters to hadronic showers. Finally, the energy of neutral hadrons is obtained from the corresponding corrected ECAL and HCAL energies. The particles reconstructed with this algorithm are referred to as PF candidates. Events selected for further study require the presence of at least one reconstructed vertex. Due to the presence of additional pp interactions within the same or nearby bunch crossings (pileup) several vertices are reconstructed. The candidate vertex with the largest value of summed object p 2 T is taken to be the primary pp interaction vertex (PV). The physics objects considered for the construction of the candidate vertex are the jets, clustered using the anti-k T jet finding algorithm [32,33] with the PF tracks assigned as inputs, and the associated missing transverse momentum, taken as the negative of the vector p T sum of those jets.
Electrons and muons are identified among the PF candidates by exploiting specific signatures in the CMS subdetectors [34,35]. Leptons reconstructed in the transition region between the barrel and endcap of the ECAL (1.4 < |η| < 1.6) are rejected to reduce efficiency differences between electrons and muons. Muons are required to pass the medium identification criteria described in ref. [34], while electrons are selected according to a multivariate discriminant based on the shower shape and track quality variables [35]. These criteria maintain approximately 99 (90)% efficiency for muons (electrons) produced in the decay of W or Z bosons [34,36]. For both lepton flavors, the impact parameter relative to the PV is required to be <0.5 mm in the transverse plane and <1 mm along the beam direction. To reject lepton candidates within jets, leptons are required to be isolated from other particles in the event. The lepton isolation variable is defined as the scalar p T sum of all PF candidates in a cone around the lepton. The cone size, defined as ∆R = (∆η) 2 + (∆φ) 2 , where φ is the azimuthal angle in radians, changes as a function of the lepton p T : ∆R = 0.2 when p T < 50 GeV, ∆R = 10 GeV/p T when 50 < p T < 200 GeV, and ∆R = 0.05 otherwise. This choice prevents efficiency loss due to the overlap of leptons and jets in events with high jet multiplicity. In order to mitigate the effect of pileup, charged particles that originate only from the primary vertex are taken into account in the calculation of the isolation variable. In addition, residual contributions from pileup to the neutral component of the isolation are subtracted using the method described in ref. [35]. The isolation variable is required to be <10 (20)% of the electron (muon) p T . The electron and muon selections are optimized to maximize the corresponding selection efficiencies, in addition to retaining similar selection efficiencies for the two flavors, in order to enhance the statistical power of some of the CRs described in section 5 that are employed to estimate SM backgrounds.
Photons are required to pass identification criteria based on the cluster energy distribution in the ECAL and on the fraction of their energy deposited in the HCAL [37]. Photons must have p T > 50 GeV, and be within |η| < 2.4, excluding the "transition region" of 1.4 < |η| < 1.6 between the ECAL barrel and endcap. Photons are required to -5 -

JHEP04(2021)123
be isolated from other PF candidates within a cone of ∆R = 0.3. To distinguish photons from electrons, we reject photons that can be associated to a pattern of hits in the pixel detector indicating the presence of a charged-particle track. To reduce the contamination due to mismeasurements of the photon energy that can create a significant p miss T , events with ∆φ( p γ T , p miss T ) < 0.4 are rejected. The vector p miss T is defined as the negative vector p T sum of all the PF candidates in the event.
To further identify additional leptons and isolated charged hadrons in the final state, isolated charged particle tracks that are identified by the PF algorithm as leptons (charged hadrons) and having p T > 5 (10) GeV are used.
Jets are clustered from PF candidates using the anti-k T clustering algorithm [32] with a distance parameter of 0.4, unless specified otherwise, implemented in the FastJet package [33,38]. Jet momentum is determined as the vectorial sum of all particle momenta in the jet, and is found from simulation to be, on average, within 5 to 10% of the true momentum over the whole p T spectrum and detector acceptance. Pileup interactions can contribute additional tracks and calorimetric energy depositions, increasing the apparent jet momentum. To mitigate this effect, tracks identified to be originating from pileup vertices are discarded and an offset correction is applied to correct for remaining contributions. Jet energy corrections are derived from simulation studies so that the average measured energy of jets becomes identical to that of particle level jets. In situ measurements of the momentum balance in dijet, photon+jet, Z+jet, and multijet events are used to determine any residual differences between the jet energy scale in data and in simulation, and appropriate corrections are made [39]. Additional selection criteria are applied to each jet to remove jets potentially dominated by instrumental effects or reconstruction failures. Jets are required to satisfy |η| < 2.4 and p T > 25 or 35 GeV, where the 25 GeV threshold is considered in regions in which the presence of jets is vetoed, in order to efficiently reject SM processes with jets, while the 35 GeV threshold is used to construct regions aiming for topologies with jets. Corrections to the jet energy are propagated to p miss T using the procedure developed in ref. [40]. As isolated prompt leptons or photons may be included in the jet definition, jets are removed from the event if they point within ∆R < 0.4 of any of the selected leptons or the highest p T photon. The DeepCSV algorithm [41] is used to identify jets produced by the hadronization of b quarks, using a working point that yields an identification efficiency of about 70% and misidentification probabilities of 1 and 12% for light-flavor or gluon jets and charm jets, respectively. These efficiencies are measured in data samples enriched in tt and multijet events as a function of jet p T and η [41] and are used to correct the prediction from simulated events. Jets passing the b-tagging criteria are required to have |η| < 2.4 and p T > 25 or 35 GeV, depending on the SR, as described in section 4.
Jets reconstructed using the anti-k T clustering algorithm with a distance parameter of 0.8 are used to identify energetic W and Z bosons that decay to qq , since their decay products are collimated into a single large radius jet. The V (V = W or Z) boson candidates have p T > 200 GeV and soft-drop masses between 65 and 105 GeV; the soft-drop mass is a groomed jet mass calculated using the mass drop algorithm [42,43] with the angular exponent β = 0 and a soft cutoff threshold z cut < 0.1. Additional selection criteria are -6 -
Samples of simulated events are used to model signal and background processes. The BSM signal events are generated using the MadGraph5_amc@nlo program 2.3.3 [47] at leading order (LO) precision, with up to two additional partons in the matrix element calculation. Samples of DY processes and photons produced in association with jets (γ+jets) are generated using the MadGraph5_amc@nlo event generator at LO precision, with up to four additional partons in the matrix element. The ttV and VVV events are produced with the same generator at next-to-LO (NLO) precision. Other SM processes, such as WW, qq → ZZ, tt, and single top quark production, are generated at NLO precision using powheg (v1.0, or v2.0) [48][49][50]. A generator-level p T -dependent next-to-NLO (NNLO)/NLO k-factor [51][52][53], ranging from 1.1 to 1.3, is applied to simulated qq → ZZ events to account for the missing higher-order matrix element contributions. Finally, the gg → ZZ process is generated at LO using mcfm 7.0 [54][55][56].
For modeling fragmentation and parton showering, generators described above are interfaced to pythia [57]  For all SM processes, the detector response is simulated through a Geant4 model [64] of the CMS detector, while BSM samples are processed using the CMS fast simulation framework [65,66]. The simulation programs account for different detector conditions in the three years of data taking. Multiple pp interactions are superimposed on the hard collision, and the simulated events are reweighed in a way that the number of collisions per bunch crossing accurately reflects the observed distribution.
Cross sections at NLO and NNLO [47,50,[67][68][69][70] are used to normalize the simulated background samples, while signal cross sections are implemented at NLO using next-toleading-logarithmic (NLL) order in α S [71][72][73][74][75][76][77][78] soft-gluon for the EW processes, or at approximately NNLO + next-to-NLL (NNLL) order in α S [79][80][81][82][83][84][85][86][87][88][89][90] for the strong production. The production cross sections for the EW GMSB model are computed in a limit of massdegenerate higgsino states, the lightest chargino ( χ ± Table 1. Summary of search category selections. In regions with the additional lepton veto selection, events containing additional leptons or charged isolated tracks are rejected. All the regions besides the edge search samples implement a veto to an additional lepton. The numbers in the rightmost column represent the edges of the bins that define the signal regions. Events in the edge search sample are further categorized as tt -like and non-tt -like as described in section 4.2.1. interpretation of the results is performed separately in each sample, we do not require the samples to be disjoint. The first (second) sample targets strong (EW)-production SUSY processes with an on-shell Z boson in the decay chain. Another sample, referred to as the "edge" sample, targets strong SUSY production with an off-shell Z boson or a slepton in the decay chain. The requirements for the fourth sample are designed to be sensitive to the direct production of a slepton pair. The selections used to define all samples are summarized in table 1. In addition to the SRs, we also define a set of CRs to be used in the estimation of the main SM backgrounds.
The baseline selection requires the presence of two OS leptons within |η| < 2.4 and with p T > 25 (20) GeV for the highest (next-to-highest) p T lepton. Each event must contain lepton flavors consistent with the corresponding requirement at the trigger level; e.g., if an event is preselected using a dilepton e + e − trigger, both leptons are required to be electrons. To avoid differences in reconstruction and isolation efficiencies between electrons and muons in boosted topologies, the two highest p T leptons are required to be separated by a distance ∆R > 0.1. The m of the dilepton system, its transverse momentum p T , and p miss T are required to be greater than 20, 50 and 50 GeV, respectively. In the SRs, the -8 -

JHEP04(2021)123
two highest p T leptons are also required to have the same flavor, e + e − or µ + µ − , while for a number of CRs we require the presence of different-flavor (DF) leptons, e ± µ ∓ .
To suppress backgrounds where instrumental p miss T arises from mismeasurements of jet energies, the two highest p T jets in the event are required to have a separation in φ from p miss T of at least 0.4, or ∆φ( p In regions with only one jet, this criterion is only applied to the single jet. If the aforementioned jet is a V boson candidate, the selection is modified to ∆φ > 0.8.

The on-Z search sample
Events with a Z boson candidate define the "on-Z" SRs and must have an invariant mass of 86 < m < 96 GeV. Events containing additional leptons or isolated tracks, as described in section 3, are rejected.

Strong-production on-Z search samples
Six disjoint event categories are defined that are expected to be sensitive to strong production of SUSY particles. These are defined on the basis of the number of jets (SRA, SRB and SRC) reconstructed with a distance parameter of 0.4 having p T ≥ 35 GeV (henceforth called n j ) and the presence of b-tagged jets (b veto and b tag). This selection is made targeting the gluino ( g ) pair production mode considered in section 1, in cases where one of the Z boson decays leptonically and the remaining, hadronically. Further requirements are made on the M T2 variable defined below, as well as H T , the scalar sum of jet p T . Each category is divided into multiple bins of p miss T , as indicated in table 1. The M T2 variable [91,92] is used to reduce the tt background contribution. It is constructed from p miss

Electroweak-production on-Z search samples
The first EW on-Z event category (referred to as "VZ" category) targets final states with a diboson pair (ZZ or ZW), with one leptonically decaying Z boson, and with the second boson decaying into jets. Depending on its momentum, the decay products of the decaying boson can either be collimated and reconstructed within a large radius jet, or resolved into two jets. For this reason, we define two subcategories, "boosted" and "resolved" that are subdivided into several bins of p miss T .

JHEP04(2021)123
For the resolved subcategory we require the presence of at least two jets, and require the two that are closest in φ to have an invariant mass m jj < 110 GeV, consistently with being a V boson decaying into jets. To reduce the tt background contribution, we reject events that have b-tagged jets with p T > 25 GeV or M T2 ( ) < 80 GeV.
In the boosted subcategory we require the presence of a large-radius jet with p T > 200 GeV, consistent with a hadronically decaying V boson candidate (n boosted V > 0). In order to ensure that the boosted and resolved categories are disjoint, events with n j > 1 are not accepted.
The last EW-production on-Z category, referred to as "HZ", is designed to be sensitive to events with an H → bb decay. Events in this category must have exactly two b-tagged jets with p T > 25 GeV and an invariant mass m bb < 150 GeV. To reduce the tt background contribution, the M T2 variable is calculated using combinations consisting of one lepton and one b-tagged jet as visible objects. Each lepton is paired with a b-tagged jet, and M T2 is evaluated for all possible b-b combinations. The smallest value of M T2 is denoted by M T2 ( b b). We require M T2 ( b b) > 200 GeV, since in tt events this variable has an endpoint at the top mass. The events are finally subdivided in bins of p miss T .

The off-Z search samples
Additional samples ("edge" and "slepton") are defined targeting models without on-shell Z bosons in the final state. The edge SRs are designed for signals with several jets in the final state and with a kinematic edge in the dilepton invariant mass distribution. The slepton SRs do not require significant jet activity in the final state.

Edge search sample
The edge sample is constructed with events with at least two jets, p miss T > 150 or 200 GeV, and M T2 ( ) > 80 GeV to reject DY and tt events. Two approaches are used to search for a kinematic edge in the m spectrum. The first one is based on a fit to the m distribution in events with p miss T > 200 GeV as described in section 6. In the second approach, we count the number of events with p miss T > 150 GeV distributed across 28 disjoint regions as described below. A looser selection on p miss T is applied, since with this categorization we can define regions with improved signal purity. First, we define seven bins in m , excluding the region 86 < m < 96 GeV, to be able to probe different positions of a possible kinematic edge. For each m bin, events are further categorized according to the b-tagged jet multiplicity, counting b-tagged jets with p T > 25 GeV. Events are also categorized as tt-like or nontt-like based on a likelihood discriminant that exploits different kinematic properties of tt events relative to a range of possible BSM contributions. We construct this discriminant as a product of probability density functions in the observables p miss T , p T , the ∆φ between the two leptons |∆φ |, and m b .
The m b variable is defined as the sum of the invariant masses of two lepton-jet pairs. Priority is given to pairs consisting of a lepton and a b-tagged jet. However, if there are no b-tagged jets in the event, we use jets without b tags. The first lepton-jet pair is selected as the one with the minimum invariant mass. The second pair is obtained by repeating the same procedure, after the exclusion of the already selected lepton and jet.

JHEP04(2021)123
The likelihood is constructed from probability density functions for each observable obtained from DF CR enriched in tt events. We use a sum of two exponential distributions for p miss T , a third-order polynomial for |∆φ |, and Crystal Ball (CB) [93] functions for both p T and m b . These distributions are found to model well those observed. The negative logarithm of the likelihood is then taken as the discriminator value used to categorize the event as either tt-like or non-tt-like.

Slepton search sample
The slepton SRs seek BSM signatures with two leptons, with p miss T (> 100 GeV), no btagged jets, and moderate jet activity. The threshold on the highest p T lepton is raised from the baseline requirement of 25 to 50 GeV, in order to further suppress the DY+jets contribution. In addition, m is required to be < 65 or >120 GeV and M T2 ( ) must be >100 GeV. Events are categorized on the basis of the jet multiplicity (n j = 0 or n j > 0), but events with one jet or more are kept only if p 2 T /p The n j > 0 category serves to recover possible BSM events characterized by moderate initial-state radiation (ISR). Events are then further split into bins of p miss T , as shown in table 1.

Standard model background
Three independent sources of SM backgrounds contribute to the SRs. The first consists of flavor-symmetric backgrounds from SM processes where SF and DF lepton pairs are produced at the same rate. The dominant process contributing to such a category is tt production. Additional contributions arise from WW, Z/γ * → τ + τ − and tW production as well as events with leptons from hadron decays. Flavor-symmetric backgrounds are estimated by constructing DF control samples in data.
The second source of backgrounds results from DY+jets events with significantly mismeasured p miss T (referred to as instrumental p miss T in what follows). This background is estimated from photon data samples in combination with CRs enriched in DY+jets events.
The third type of SM backgrounds consists of processes yielding final states with an SF lepton pair produced in the decay of a Z boson or a virtual photon accompanied by neutrinos (ν) produced in the decay of a W or Z bosons. The main process contributing here is VZ production. Rarer processes, such as ttZ production, also contribute to certain SRs. These backgrounds are referred to as Z+ν backgrounds and are estimated from simulation. The prediction is validated in dedicated data control regions.

Flavor-symmetric backgrounds
As already mentioned, the estimation of flavor-symmetric backgrounds exploits the fact that in such processes, the DF and SF events are produced at the same rate. The CRs are defined in data with the same selections as the corresponding SRs, but requiring the presence of a DF lepton pair instead of an SF pair. The background contribution in the SR is then predicted by means of a transfer factor, denoted by R SF/DF , that accounts for the differences in reconstruction, identification and trigger efficiencies between DF and SF events. These are caused by the residual differences in the efficiencies between electrons and -11 -

JHEP04(2021)123
Year  Table 2. Summary of the r µ/e parameters obtained by fitting the lepton p T and η, in a DYenriched control region, for different data taking years. Only statistical uncertainties are shown.
muons. The transfer factor consists of the product of two correction factors, determined from CR data. The first correction factor, r µ/e , is the ratio of muon and electron reconstruction and identification efficiencies measured in a region enriched in DY events, requiring two SF leptons, at least two jets, p miss T < 50 GeV, and 60 < m < 120 GeV. Assuming that the efficiency for each of the two leptons in the event is independent of the other lepton, r µ/e can be defined as is the number of µ + µ − (e + e − ) events. The r µ/e factor is parametrized as a function of the lepton p T and η by the following empirical form: and The constants a 1 , a 2 , b 1 , c 1 , c 2 , and r 0 µ/e are extracted in a fit to the r µ/e computed in data in bins of the η and p T of the positive lepton in the DY-enriched sample. The fit is performed iteratively, in which the p T and η dependencies, and the normalizations, are extracted in separate steps. These values, shown in table 2, are obtained separately for each data taking year and found to be statistically consistent with those predicted from simulation. A greater dependency on η is observed in the r µ/e factor in data collected in 2016 and 2017 that is caused by a loss in the transparency of the ECAL endcap crystals, which affected trigger performance and was corrected in the 2018 data. The transparency loss and its effects are stronger in data collected in 2017. We assign systematic uncertainties of 5% to the measured r µ/e value and an additional 5% for each of its p T and η parametrizations that cover possible residual kinematic dependence.
Neglecting differences in trigger efficiencies, r µ/e can be used to estimate the number of SF (e + e − and µ + µ − ) events from the observed number of DF events (N DF ) in the DF CR as follows: N est.
where T µµ , T ee and T eµ are the trigger efficiency of the di-muon, di-electron and muon-electron triggers respectively, is JHEP04(2021)123 then used to account for differences between SF and DF dilepton trigger efficiencies T . These efficiencies are measured in data and found to be 85-95%, depending on the lepton flavor and the data taking period. The resulting R T coefficient is measured to be 1.03-1.05, with an uncertainty of 4-5%.
The transfer factor R SF/DF , used to predict SF events from DF ones, is finally defined as: The background estimation method is validated in data in a CR enriched in flavorsymmetric tt events. This region is defined by requiring an SF lepton pair, exactly two jets, and 100 < p miss T < 150 GeV. Events with 70 < m < 110 GeV are rejected to reduce the contribution from DY events. Figure 3 compares the prediction from a DF selection in SF data in this region, as a function of different kinematic variables. An agreement within the uncertainties is observed, thus validating the background estimation method.
The statistical uncertainty arising from the limited size of the DF control sample represents the dominant contribution to the total uncertainty in the flavor-symmetric background prediction. For the estimation of this background in the on-Z SRs, where 86 < m < 96 GeV, the m requirement in the DF control sample is relaxed to m > 20 GeV, and an additional multiplicative factor, κ = N DF (86 < m < 96 GeV)/N DF (m > 20 GeV), is used to account for the different m selection in CRs and SRs. This factor is determined from dedicated DF CRs in data, defined by relaxing or merging a subset of selection requirements described in section 4. The regions of interest (SRA, SRB and SRC strong-production SRs, and the HZ and resolved VZ SRs) are defined in table 1. The boosted VZ SR is also considered, relaxing the veto of additional jets. In these regions, κ is measured to be in the range 0.045-0.067. We also determine κ as a function of several kinematic variables to assess the possible dependencies. Based on these measurements, we assign a systematic uncertainty of 20% to the value of κ to cover such effects.

Drell-Yan+jets backgrounds
The contribution from DY+jets events to the SRs mainly arises from mismeasurements of momenta of reconstructed objects affecting p miss T . In regions where jets in the final state are required, instrumental p miss T arises mainly from jet energy mismeasurement, and the p miss T "templates" method [23-26] is used to estimate the resulting background contribution. In the slepton SRs, since only jets with low p T are present, we use a different method exploiting a CR enriched in DY+jets events.
The p miss T "templates" method relies on the fact that instrumental p miss T in DY+jets events is caused by limited detector resolution in measuring the p T of the jets recoiling against the leptonically decaying Z boson. Since the p T resolution of leptons and photons is much better than that of jets, the p miss T distribution in DY+jets events can be estimated directly from γ+jets data.
The γ+jets events are selected with jet requirements identical as those used in defining the SRs in section 4. We assume that the γ+jets events are not affected by potential contamination from any of the BSM physics considered in this search.  The M T2 variable used to select events in several SRs requires the presence of two visible objects and therefore cannot be defined in the γ+jets sample. Instead, its behavior is emulated by mimicking the decay of the photon into two leptons. The decay is modeled assuming that the leptons arise from a particle that has the mass of a Z boson and the momentum of the selected photon, with the angular distributions in the decay as expected at LO in perturbation theory. The simulated leptons are used to calculate the M T2 ( ) variable in the γ+jets data sample.
Events with genuine p miss T may, in fact, be present in the γ+jets sample, originating from EW processes such as Wγ + jets production, where the W boson decays to ν.

JHEP04(2021)123
However, such contributions can be suppressed by rejecting events that contain additional leptons. The residual EW contamination in the γ+jets sample, which is larger at large p miss T , is subtracted using simulation.
The photon p T distribution in γ+jets events is expected to differ from that of the Z boson in DY+jets, mainly because of the different boson masses. Thus, simulation is used in each SR to obtain a set of weights that match the photon p T distribution to the expected Z boson p T distribution. These weights are then used to reweigh the p miss T templates in γ+jets data in the SRs. After this correction, the corrected p miss T template in each SR is normalized based on the observed yield in dilepton data in the range 50 < p miss T < 100 GeV, where DY+jets events dominate the data sample. We note that to account for potential contamination from BSM physics the 50 < p miss T < 100 GeV bin in each SR is included in the signal extraction fit described in section 8.
Several sources of uncertainty are considered for the DY+jets background prediction: the statistical uncertainty arising from the limited size of the γ+jets sample in each p miss T bin, the systematic uncertainty in the EW subtraction, and the statistical uncertainty in the template normalization arising from the dilepton data yield in the range 50 < p miss T < 100 GeV. An additional systematic uncertainty is assessed through a closure test of the method in simulation, where the p miss T distribution in simulated DY+jets events is compared to the distribution obtained by applying the background prediction method to a γ+jets simulated sample. In each p miss T bin, we assign an uncertainty equal to the largest of the differences between the predicted and simulated yields, and the statistical uncertainty reflecting the size of the samples. The resulting uncertainty ranges between 20 and 100% across the search bins with the largest values obtained in bins affected by the limited number of simulated events.
The validity of the method is further tested in data CRs enriched in events containing instrumental p miss T . These samples are defined by inverting the ∆φ( p j 1,2 T , p miss T ) selection (or, in the boosted VZ region, ∆φ(V boson candidate, p miss T )). In addition, the b-tagged jet multiplicity categorization is removed from the on-Z strong-production regions yielding a total of six validation regions (VRs) with the same p miss T binning as used in the corresponding SRs. The observed p miss T distribution is compared to the prediction in the VRs in figure 4 showing agreement within the uncertainties.
The method described above is also used to predict the DY+jets background in the edge SRs, where events with 86 < m < 96 GeV are rejected, and therefore the contribution from DY+jets events is expected to be small. In this case, the prediction is obtained from a CR with inverted m selection, by means of a transfer factor r in/out defined as the ratio of the DY+jets yield in a given m bin over the yield in the range 86 < m < 96 GeV. The r in/out ratio is measured in a data control sample enriched in DY+jets events, obtained by requiring at least two jets, with p miss T < 50 GeV and M T2 ( ) > 80 GeV, after subtracting the flavor-symmetric contribution estimated as described in section 5.1. The r in/out value is measured to be in the range 0.003-0.06, depending on the m bin. We assign a systematic uncertainty in r in/out to cover its possible dependence on p miss T and n j , of 50 (100)% in m bins below (above) 150 GeV.    Table 3. Summary of the uncertainties in background estimations performed on data.
In the slepton SRs, the DY+jets background is estimated in each p miss T bin using a CR in data enriched with DY+jets events, obtained by applying the same selection criteria as used in the SRs, but inverting the selection on m (65 < m < 120 GeV). The prediction is then obtained by means of a transfer factor, r in/out , which is measured in data, after relaxing the p miss T and n j selections applied in the SRs. The r in/out value is measured to be 0.07, with a 50% uncertainty obtained from a closure test performed using simulated DY+jets events. To account for possible contamination from BSM physics in the 65 < m < 120 GeV region, that region is included in the signal extraction fit described in section 8.
The systematic uncertainties associated with the flavor-symmetric and DY+jets background estimation are summarized in table 3.

T
Backgrounds from events with Z/γ * bosons and genuine p miss T such as WZ, ZZ, and ttZ can be important in SRs of large p miss T , and are estimated directly from simulation. Dedicated data CRs of trileptons and two pairs of OSSF leptons are used to determine the overall normalization and to check the modeling of such events in simulation. Systematic uncertainties as large as 50% are assessed for each process to cover differences between data and simulation. In predicting the ZZ yield we also assign an additional uncertainty given by the difference between the nominal NLO simulation and the NNLO prediction achieved applying the k-factor as described in section 3. Finally, we include statistical uncertainties

Edge fit to the dilepton invariant mass distribution
We perform a simultaneous unbinned maximum likelihood fit as a function of m in e + e − , µ + µ − , and e ± µ ∓ data to search for a kinematic edge. The fit is performed in the "edge fit" SR defined in section 4. The functional forms used to model the signal and the two main SM background components (flavor-symmetric background and backgrounds arising from other SM processes containing a Z boson) are described below.
The flavor-symmetric background component is modeled using the CB function:

JHEP04(2021)123
This model has five free parameters: the overall normalization, the mean µ CB and the full width σ GB at half maximum of the Gaussian core component, the transition point α between the Gaussian core and the power-law tail, and the power-law parameter n. Backgrounds containing a leptonically decaying Z boson (Z/γ * +X) are modeled through a sum of an exponential function, which describes the rise at small mass, and a Breit-Wigner function with the mean and the width set to the nominal Z boson values [94] convolved with a double-sided CB function, P DSCB (m ) to account for the experimental resolution: where µ DSCB and Γ DSCB are the mean and width, respectively, of the CB function, and α 1 and α 2 are the transition points. The model for the Z/γ * +X background line shape is thus: where P BW and P exp are the Breit-Wigner and exponential functions, respectively. The complete DY+jets background model has therefore nine free parameters each for the e + e − and µ + µ − final states. The signal component is described by a triangular form, inspired by the slepton edge models [12], convolved with a Gaussian function to account for the experimental resolution: The signal model has two free parameters: the fitted signal yield and the position of the kinematic endpoint, m edge , as the experimental resolution Γ for each leptonic final state is obtained from the CB function of the DY+jets model. In an initial step, a fit to data is performed in a DY+jets-enriched CR with at least two jets, M T2 ( ) > 80 GeV, and p miss T < 50 GeV, separately for e + e − and µ + µ − events, to determine the parameters for backgrounds containing a Z boson. The final fit is then performed simultaneously to the invariant mass distributions in the e + e − , µ + µ − , and e ± µ ∓ data samples. The model for the flavor-symmetric background is varied consistently in the SF and DF samples. The relative normalization of SF and DF events is given by the R SF/DF factor, which is treated as a nuisance parameter, constrained by a Gaussian prior with the mean value and standard deviation (s.d.), as determined in section 5.1. In total, the final fit has ten parameters: a normalization parameter for each of the three fit components, four parameters for the distribution of the flavor-symmetric background, R SF/DF , the relative fraction of dielectron and dimuon events in the flavor-symmetric prediction, and the position of the signal edge. Out of these, only R SF/DF is constrained, while the others are treated as free parameters of the fit.

JHEP04(2021)123 7 Results
The observed yields in each SR are compared to the SM predictions for the on-Z, edge, and slepton SRs. In the search for an edge, a fit is also performed to the m distribution in data to find a kinematic edge in the m spectrum as discussed in section 6.

Results for the on-Z samples
The results for the strong production on-Z SRs are summarized in table 5. The corresponding p miss T distributions are shown in figure 5. No significant deviations are observed relative to the SM background. The largest disagreement corresponds to one of the SRA b tag categories in which 42 events are observed and 31.4 ± 3.8 background events are expected, corresponding to a local significance of 1.4 s.d.
The results for the EW-production on-Z SRs are summarized in table 6. The corresponding p miss T distributions are shown in figure 6. The observed data yields are consistent with the SM background predictions. The largest discrepancy between data and prediction occurs in the highest p miss T bin of the resolved VZ regions, where 2 events are observed while 6.3 ± 2.2 are predicted, corresponding to a deficit with a local significance of 1.2 s.d.

Results for the edge search samples
Comparisons between the SM predictions and the observed data in the 28 edge SRs are summarized in table 7. A graphical representation of the same results is displayed in figure 7.
We find an agreement between the observed data and SM predictions in all SRs. The largest deviation is observed in the tt-like region for 300 < m < 400 GeV and n b > 0, in which 49 events are observed and 76 +10 −9 were expected, corresponding to a deficit in data with a local significance of 2.4 s.d. We also observe a slightly larger number of events than the background prediction in the high m non-tt-like category, but the predictions agree within one s.d.
The dilepton mass distributions and the results of the kinematic edge fit are shown in figure 8 while Table 8 presents a summary of the fit results. A best fit signal yield of 27±22 events is obtained when evaluating the signal hypothesis in the edge fit SR with a fitted edge position of m = 294 +12 −20 GeV, assuming the signal normalization to be nonnegative. To test the compatibility of this result with the background-only hypothesis, we estimate the global p-value [95] of the result using the test statistic −2 ln Q, where Q denotes the ratio of the fitted likelihood value for the signal+background hypothesis to that for the background-only hypothesis. The test statistic −2 ln Q is evaluated in data and compared to the corresponding quantity computed using a large sample of background-only pseudoexperiments where the edge position is not fixed to any particular value. The resulting p-value is interpreted as the one-sided tail probability of a Gaussian distribution, and corresponds to an excess in the observed yields relative to the SM background prediction at a global significance of 0.7 s.d. If unphysical negative signal yields are permitted, the best fit corresponds to a negative signal yield with an edge position of 34.4 GeV and a global significance of 1.8 s.d.
-20 -     Table 6. Predicted and observed event yields in the EW-production on-Z search regions, for each p miss T bin as defined in table 1 before the fits to data described in section 8. Uncertainties include both statistical and systematic sources. The p miss T template prediction in each SR is normalized to the first p miss T bin of each distribution in data.         Table 9. Predicted and observed event yields in the slepton search and control regions. A background-only fit to observation in the CR is performed to determine the DY+jets contribution as described in section 8. Uncertainties include both statistical and systematic sources.

Results in the slepton search regions
The p miss T distribution of data events in the slepton SRs is shown together with the SM background predictions in figure 9. Results are also summarized in table 9. The observed data yields are consistent with the SM predictions. The largest discrepancy between data and SM prediction is observed in the highest p miss T bin of the SR without jets where 17 events are observed and 9.3 ± 2.3 are predicted, corresponding to a local significance of 1.6 s.d.

Interpretation of the results
The results are interpreted in the context of models of BSM physics presented in section 1. hypotheses to the data in the SRs and some of the CRs: the event yield with 50 < p miss T < 100 GeV in the on-Z SRs and the on-Z category in the slepton SRs. The uncertainties in the modeling of the backgrounds, summarized below, are inputs to the fitting procedure. The likelihoods are constructed as the product of Poisson probability density functions, one for each SR, with additional nuisance parameters constrained by log-normal terms that account for uncertainties in the background predictions, and in signal yields when the signal is included in the hypothesis. When a CR for a given background is included in the fit, the normalizations of both signal and that background are treated as free parameters. This accounts for the possible presence of signal events in the CRs (signal contamination). The signal+background fits are used to set 95% confidence level (CL) upper limits on the production cross sections of the modeled signals. We employ a modified frequentist approach, using the CL s criterion and relying on asymptotic approximations, to calculate the distribution of the profile likelihood test statistic [96][97][98][99]. The limits are then used, together with the theoretical cross section, to exclude ranges of masses for the BSM particles involved in each model.

Systematic uncertainties in the signal
We include uncertainties in the expected signal for all of the SMS models, as summarized in table 10. The integrated luminosities of the 2016, 2017, and 2018 data-taking periods are individually known with uncertainties in the 2.3-2.5% range [100][101][102], while the total Run 2 (2016-2018) integrated luminosity has an uncertainty of 1.8%, the improvement in precision reflecting the (uncorrelated) time evolution of some systematic effects. We also include uncertainties in the lepton trigger, identification, and isolation efficiencies, in the b tagging efficiencies and the mistag probability. To check the modeling of ISR in the EW (strong) signal simulations we obtain distributions in the number of ISR jets in data samples enriched in DY+jets (tt). We derive weights as ratios of these distributions to simulation. The systematic uncertainty in the ISR modeling is then given by the difference between weighted and unweighted simulations of the signals.
Additional uncertainties arise from the potential mismodeling of pileup, from JES, from the choice of the µ R and µ F scales used in the event generator [103][104][105], and from the limited number of simulated events. Finally, any further possible mismodeling of lepton efficiencies, jet energy response, b tag efficiency, and the p miss T distributions associated with the CMS fast simulation framework is accounted for with an additional uncertainty. The assumed correlations in the signal uncertainties across SRs are the same as those described in section 5 for the background processes. Uncertainties in ISR modeling and fast simulation are treated as correlated across the SRs and across the data samples.

Interpretations of the results using simplified SUSY models
The results in the strong-production on-Z search regions are interpreted using an SMS model of gluino pair production, discussed in section 1. This signal is characterized by final states with substantial activity (energy in jets). All strong-production on-Z SRs are included in the maximum likelihood fit in order to maximize the acceptance in models in which the g and χ 0 1 masses are close, where less jet activity is expected. The upper limit at 95% CL on the signal production cross section is shown in figure 10, as a function of the g and χ 0 1 masses, together with the expected and observed exclusion contours. We exclude gluino masses up to 1600-1870 GeV, depending on the mass of χ 0 1 , extending the reach of previous CMS results [21] by approximately 100 GeV.
Upper limits at 95% CL on the cross section of models for χ ± 1 χ 0 2 and for χ 0 1 pair production are set using the results of the EW-production on-Z SRs. For χ ± 1 χ 0 2 production with decays to WZ, the VZ search regions provide most of the sensitivity. While the resolved VZ SR is sensitive to a wide range of χ ± 1 / χ 0 2 and χ 0 1 mass hypotheses, the use of the boosted VZ region improves the sensitivity for masses of the χ ± 1 / χ 0 2 much larger than the mass of χ 0 1 , where the bosons produced in the decay chain receive a large Lorentz boost. Figure 11 shows the cross section upper limits and the exclusion contours at 95% CL obtained for this model as a function of the χ ± 1 / χ 0 2 and χ 0 1 masses. We exclude χ ± 1 / χ 0 2 masses up to 750 GeV, extending the reach of ref.
[21] by approximately 100 GeV. Two models are considered for χ 0 1 pair production. One assumes that both χ 0 1 decay into a Z boson with a 100% branching fraction. The other assumes that each χ 0 1 can decay to either Z or H with equal probability. The first model leads to signatures with a pair of Z bosons, with most of the signal events expected to populate the VZ SRs. On the other hand, signal events where an H decays to bb are expected to populate the HZ region. The observed and expected upper limits at 95% CL on the production cross section times branching fraction product for both models are shown in figure 12. In these two scenarios, we are able to exclude χ 0 1 masses up to 800 and 650 GeV respectively, extending the reach of ref.
[21] by approximately 150 GeV. Figures 11 and 12 show observed exclusion limits that are more stringent than the expected ones. In both cases this arises from the downward fluctuation in the observed data yields in the two highest p miss T SRs of the resolved VZ search region, as discussed in section 7.1.
The edge search regions used in the counting experiment serve to constrain the two slepton edge models presented in section 1 (figure 2 left and middle diagrams). Figure 13 shows the upper limits at 95% CL on the production cross section for both of these models. We exclude bottom (light-flavor) squark ( q ) masses up to 1300-1600 (1600-1800) GeV, depending on the assumed χ 0 2 mass. For the case of the bottom squark pair production, we improve the results from ref.
[21] by up to 300 GeV. The observed exclusion limits are more stringent than expected for models with small χ 0 2 mass. This is caused by a mild deficit of observed events in the non-tt-like edge regions for 20 < m < 60 GeV.

Summary
A search is presented for phenomena beyond the standard model in events with two oppositely charged same-flavor leptons and missing transverse momentum in the final state. The search is performed in a sample of proton-proton collisions at √ s = 13 TeV collected with the CMS detector corresponding to an integrated luminosity of 137 fb −1 . Search regions are defined to be sensitive to a wide range of signatures. The observed event yields and distributions are found to be consistent with the expectations from the standard model. The results are used to set upper limits on the production cross sections of simplified models of supersymmetry. Gluino masses are excluded up to 1870 GeV, light-flavor (bottom) squark masses up to 1800 (1600) GeV, chargino (neutralino) masses up to 750 (800) GeV, and slepton masses up to 700 GeV, typically extending the reach over previous CMS results by a few hundred GeV.

JHEP04(2021)123
we gratefully acknowledge the computing centers and personnel of the Worldwide LHC Computing Grid for delivering so effectively the computing infrastructure essential to our analyses. Finally, we acknowledge the enduring support for the construction and operation of the LHC and the CMS detector provided by the following funding agencies: