Characterizing Dark Matter Signals with Missing Momentum Experiments

Fixed target missing-momentum experiments such as LDMX and M$^3$ are powerful probes of light dark matter and other light, weakly coupled particles beyond the Standard Model (SM). Such experiments involve $\sim$ 10 GeV beam particles whose energy and momentum are individually measured before and after passing through a suitably thin target. If new states are radiatively produced in the target, the recoiling beam particle loses a large fraction of its initial momentum, and no SM particles are observed in a downstream veto detector. We explore how such experiments can use kinematic variables and experimental parameters, such as beam energy and polarization, to measure properties of the radiated particles and discriminate between models if a signal is discovered. In particular, the transverse momentum of recoiling particles is shown to be a powerful tool to measure the masses of new radiated states, offering significantly better discriminating ability compared to the recoil energy alone. We further illustrate how variations in beam energy, polarization, and lepton flavor (i.e., electron or muon) can be used to disentangle the possible the Lorentz structure of the new interactions.

A particularly promising accelerator-based strategy involves the fixed-target missing-momentum (MM) concept [18]. In this setup, a low-current O(1 − 10) GeV lepton beam is passed through a thin target, which is surrounded on both sides by tracking material and positioned upstream of a veto detector; the energy and momentum of individual beam particles are measured on both sides of the target. If DM (or any other invisible or long-lived particle) is produced in the target, the beam particle loses a large fraction of its energy and momentum, and no other visible particles are observed in the veto detector.
This technique has several appealing features: • Experimental Control: Unlike direct detection searches, whose sensitivity is subject to to astrophysical uncertainties and environmental backgrounds, the MM signal strength is fully calculable and irreducible background events occur at a rate of ∼ 10 −15 per incident beam particle [14].
• Coverage Breadth: Since DM production at accelerators is relativistic, the signal strength is largely insensitive to the Lorentz structure of the underlying interaction. Consequently, the same MM search simultaneously covers a wide range of DM spin and coupling varieties.
• Parametric Enhancement: Unlike beam-dump DM searches whose signal involves DM production followed by its scattering in a downstream detector, the MM setup only requires DM production. Thus, the signal rate depends only on the production rate without the added penalty of a small scattering probability.
The electron beam MM strategy is currently being developed by the LDMX collaboration [14,34] and additional studies are underway to explore a future muon beam MM experiment (M 3 ) at Fermilab [16]. Collectively, these efforts have the potential to test nearly every model of sub-GeV freeze-out for which the MM signal strength is parametrically related to the DM annihilation rate in the early universe [12,17,23]. 1 In this work we explore the optimistic scenario in which one of these experiments convincingly discovers a new physics signal. The existing MM literature has largely focused on the kinematics of signals arising from the radiative production of massive spin-1 bosons with vector current interactions (e.g. kinetically mixed dark photons) [14,18,34]. In this paper we generalize these studies to explore the kinematics of recoiling beam leptons in fixed target interactions where = e or µ, and X represents a broad range of possible final state particles, including single boson emission in 2 → 3 scattering and direct DM pair production in 2 → 4 processes as depicted in Fig. 1 1. Schematic diagrams depicting the signal topologies outlined in Sec. II. Left: two-step DM production via 2 → 3 scattering in which a massive mediator particle X is produced on shell and subsequently decays invisibly via X →χχ. Middle: Direct DM production through virtual light mediator exchange corresponding to the (ēe)(χχ)/q 2 "millicharge" interaction in Eq. (4). Right: Direct DM production through a contact (ēe)(χχ)/Λ 2 interaction from Eq. (4) represented by the black circle. (2)-(4) for 4 (solid lines) and 8 (dashed lines) GeV electron beams on a tungsten target. The horizontal axis represents the mediator mass m A for the radiative production of A vector mediators through the interaction in Eq. (2) (top) and 2mχ for the DM pair production through the interactions in Eq. (4). The heavy mediator case corresponds to m A = 10 GeV and exhibits a systematic increase in the cross section as the beam energy is increased. This behavior reflects the higher-dimensionality of the operator mediating this interaction (first line of Eq. (4)).
we develop novel strategies and statistical tests to distinguish these different model hypotheses and quantify the signal samples required for model discrimination.
This paper is organized as follows: in Section II we introduce representative phenomenological interactions for MM signals and describe the details of our numerical simulations; in Section III we evaluate the model discrimination potential of final state kinematic variables; in Section IV we study how varying the initial state beam energy and polarization yields new observables; in Section V we explore the additional clues that muon beams can offer; we offer some concluding remarks in Section VI.

II. MODELS AND SIMULATIONS
Throughout this paper our goal is to distinguish between representative signal models categorized according to how they couple to leptons or the mass of new particles involved. In this section, we define these interactions by their fixed-target production mode as radiative 2 → 3 reaction N → N X, where X is an invisibly decaying particle, or a 2 → 4 reaction that pair-produces the DM χ directly through an off-shell mediator. 2

A. Mediator Production: 2 → 3 Processes
The nominal LDMX/M 3 signal process involves renormalizable interactions through which a single "mediator" particle X is emitted in N → N X reactions inside the target, where = e or µ is the lepton beam particle and N is a target nucleus.
The renormalizable (mass dimension ≤ 4) interactions with new spin-1 states are where we have adopted an arbitrary overall normalization e to match the interaction of a kinetically mixed dark photon; e = √ 4πα EM is the electric charge and parametrizes the strength of this interaction relative to electromagnetism. Here P L,R ≡ (1 ∓ γ 5 )/2 is a left (right) projector and we consider various mass choices for each possible scenario m A , m Z , m V . Constraints on light new vector boson interactions with anomaly-free couplings to Standard Model (SM) particles are summarized in Refs. [36,37]; bounds on chiral and axial vector interactions are more model dependent, but are typically very strong because such interactions require additional SM-charged field content for anomaly cancellation [38,39]. We also consider renormalizable spin-0 interactions for scalar and pseudoscalar interactions with corresponding masses m s and m a GeV. Note that unlike the vector interactions in Eq. (2), which can be gauge-invariant at high energies, the Yukawa couplings in Eq. (3) are not invariant under the SU (2) L × U (1) Y gauge group, so they must arise from higher dimension operators proportional to the source of electroweak symmetry breaking. We therefore expect ∝ v/F , where v = 246 GeV is the Higgs vacuum expectation value and F is the mass scale of the heavy particle whose quantum numbers restore SM gauge gauge invariance. As such, the bounds on these interactions depend on the details of the ultraviolet completion; they are summarized in Refs. [40][41][42][43][44][45] for many representative examples.

B. DM Pair Production: 2 → 4 Processes
If the mediators in Sec II A are too heavy for direct production or if their decays to DM pairs are kinematically forbidden, the DM can still be produced directly via N → Nχχ reactions with virtual mediator exchange. Here we consider two representative limiting cases in which the amplitudes for DM pair production are proportional to where Λ GeV represents the mass scale of a vector mediator that has been integrated out to generate a contact operator and q is the momentum imparted to theχχ system in the limit Λ < 2m χ where the mediator cannot decay to DM particles. Note that in our numerical studies the heavy and light mediator cases are modeled using renormalizable interactions from Sec. II A by taking the mediator mass m A E beam or m A 2m χ , respectively. These two classes of processes are represented schematically by the middle and right diagrams in Fig. 1. We note that the list of operators in Eq. (4) is not exhaustive and can include different Lorentz struc-   tures (e.g., additional γ 5 insertions from integrating out the axial vector in Eq. (2)).

C. Simulation Details
For the numerical studies in the remainder of this paper, we generate signal samples for electron and muon beams based on Eqs. (2) -(4) using version 2.6.4 of MadGraph 5 aMC@NLO [46] including elastic and inelastic atomic and nuclear form factors for the target [47,48]. The inclusive signal production cross sections for these scenarios are presented in Fig. 2 as functions of either the mass of the radiated particle for the vector mediator in Eq. (2) or as a function of 2m χ in Eq. (4) where appropriate.
Although our qualitative conclusions below are largely independent of any particular experimental setup, for electron beam studies, our simulation is designed to match the anticipated LDMX phase 1 design with a 4 GeV electron beam impinging on a thin tungsten target [34]. Similarly, our muon beam simulation is motivated by the M 3 concept, which has been studied for a 15 GeV beam energy and also with a tungsten target [16]. For our signal samples we select events with E e ≤ 1.2 GeV and E µ ≤ 9 GeV for LDMX and M 3 , respectively.

III. KINEMATIC VARIABLES
This section generalizes earlier studies of dark photon production in missing momentum experiments [14,18,34] by using kinematic variables -E beam (beam energy), E e (electron recoil energy), and p T,e (recoil electron transverse momentum) -as tools for distinguishing various new physics scenarios. The differential distributions of E e and p T,e are shown in Figs. 3 and 4 for the on-shell emission of spin-1 and spin-0 particles (Sec. II A), and in Fig. 5 for direct DM production (Sec. II B) via heavy and light mediators; in Fig. 6 we directly compare the distributions of these models. We discuss these results in more detail below.

A. On-Shell Mediators: Mass Measurement
On shell mediator emission arises in 2 → 3 processes as shown in the left panel of Fig. 1. The corresponding kinematic distributions shown in Fig. 3 and 4 are sensitive to the mass of the emitted particle. It is therefore interesting to investigate the ability of missing momentum experiments to distinguish different mass hypotheses. We quantify this discriminating power using a simple likelihood ratio test as described below.
We start by representing a histogram of experimen- tal events (either in E e , p T,e or both) as a vector d = (d 1 , . . . , d N ). We treat each bin as an independent counting experiment with Poisson statistics; the d i are then drawn from the Poisson distribution f P (d i , ν i ) where ν are the predicted means in each bin which define a model hypothesis. The joint probability distribution or likelihood for d to be observed given the means ν is [49,50] Models that describe the data better have a larger likelihood. For computational simplicity we work with the logarithm of the likelihood and define a test statistic (TS) λ to compare two different models A and B for a given data set This test statistic is negative when model A is preferred over B and positive otherwise (in the limit of large statistics in all bins λ simply becomes a difference of χ 2 values for the two models). We use our Monte Carlo event samples to generate many mock experiments and study the resulting distribution of λ for various combinations of hypotheses. If there are enough signal events, one can reject the A hypothesis with a given confidence level CL compared to B if the probability of obtaining λ < 0 is 1 − CL. In other words, we can find the number of events N sig such that where f (λ; N sig ) is the distribution of the test statistic for a given number of signal events. Note that due to the finite size of any Monte Carlo sample, a sufficiently fine binning of the events (or binning in multiple variables) will result in small statistics and significant fluctuations in individual bins. We address this issue by using a modification of Eq. (7) based on Ref. [51] as described in Appendix A. We also do not account for detailed experimental effects such as energy and momentum smearing, so our results should be viewed as an idealized best-case scenario. However, our choices of p T, and E binning are motivated by the detailed simulation-based detector studies of Refs. [14,34]. In particular, we use 5 MeV p T bins and 1% E bins. While the latter is almost certainly an overly optimistic choice, we will show that p T spectra still offer superior model discrimination ability compared to recoil energy.
We show an application of this TS in Fig. 7. The upper left panel shows a histogram of a mock dataset with 100 events drawn from a high-statistics MC p T,e distribution for a 100 MeV vector mediator emitted on shell as in Fig 1 (left). This distribution is compared to expected distributions for a set of different masses. It is clear that some hypotheses fall outside of the gray 1σ statistical uncertainty band of the mock data and are therefore disfavored, while others are indistinguishable within errors. Equations (8) and (9) allow us to quantify this observation. The result is shown in the upper right panel of Fig. 7, which illustrates the number of signal events needed to distinguish a test mass hypothesis on the vertical axis from the "true" model (i.e., the one that was used to generate the mock data set) on the horizontal axis. We see that very different masses can already be distinguished with the background-free "discovery threshold" number of signal events N sig = 3, while comparable masses will require O(10−100) events to disentangle.

B. Missing Momentum vs. Missing Energy
Armed with the statistical test introduced above, we can now compare the model discrimination power of missing momentum experiments like LDMX and M 3 against other fixed-target lepton beam techniques that only measure the missing energy of the beam (e.g., the NA64 experiment [52,53]).
As in the previous section, we use simulated data to estimate the number of signal events required to distinguish a test hypothesis from the "true" model using recoil energy alone at a given confidence level, yielding a twodimensional histogram similar to the upper right panel of Fig. 7. We compare this to p T,e -only discrimination by subtracting the two histograms. The resulting histogram is shown in the bottom right panel of Fig. 7. This difference is strictly positive, implying that transverse momentum enables superior model discrimination with a much smaller sample of signal events.
For illustration, in the left column of Fig. 7 we also show various p T,e (top) and E e (bottom) test hypothesis templates (colored curves) plotted alongside N sig = 100 events of simulated data drawn from a vector mediator sample with m A = 100 MeV (black curve). The gray band surrounding the black curve is the statistical uncertainty of the mock data. Visually it is clear that, relative to the E e distributions, the p T,e distributions span a greater variety of shapes for the same model parameters and thereby offer more discriminating power. For nearly the entire mass range shown in this figure, recoil energy distributions require O(100) events to discriminate between widely separated mass hypotheses (e.g., MeV and GeV), whereas transverse momentum can already exclude closely spaced choices (e.g., 10 and 100 MeV) requiring only a few events, near the discovery threshold for a zero-background environment.

C. Off-Shell Mediators
Ref. [23] has shown that missing momentum experiments are sensitive to a wide range of mediator masses in direct DM production, including off-shell production through heavy mediators. In this section, we instead focus on the kinematic features of these signals and extend the discussion to the emission of invisible particles through a massless off-shell mediator (e.g., production of millicharged particles via the photon).
In Fig. 5 we show the distributions of electron recoil energies and transverse momenta for direct DM or millicharged particle production for a few different masses. These distributions look qualitatively similar to on-shell Two dimensional histogram of the number of signal events required to distinguish a given test mass hypothesis on the vertical axis from the true mass on the horizontal axis (the one used to generate the mock dataset) at 95% confidence level assuming only statistical errors. Darker colors correspond to more signal events needed. This result makes use of the pT,e spectral information only. Bottom left: same as top left, but showing the electron recoil energy (Ee) instead of transverse momentum (pT,e). Bottom right: comparison of recoil energy-only versus transverse momentum-only analyses to discriminate between different mass hypotheses. The density plot shows the difference between the number of signal events required to distinguish a test model from the "true" scenario using only Ee and pT,e information. This difference is strictly positive, implying that transverse momentum enables superior model discrimination with a much smaller number of signal events. This improvement is particularly pronounced for models with similar masses where the recoil energy spectra are nearly identical. mediator emission. However, in the off-shell case, the invariant mass of the emitted particles is not fixed, so one naively expects broader distributions in p T,e compared to the on-shell case. This is explicitly illustrated in right plot of Fig. 6 in the heavy mediator regime. A massless mediator leads to an enhancement of the amplitude at invariant mass close to the kinematic minimum of 2m χ , which compensates for the expected broadening. Consequently, the transverse momentum distributions in the millicharge and on-shell cases look very similar. Figure 6 illustrates another important point: the kinematic distributions are somewhat degenerate in "theory space." That is, a given distribution can be interpreted as coming from an on-shell mediator of a certain mass, or from off-shell mediators and DM of different masses. Thus, it is essential to have a complementary array of experiments that can probe these interactions at different √ s. For example, it is clear from Fig. 6 that a MM experiment with a few signal events will not be able to distinguish between the production of an on-shell mediator with m A = 350 MeV or direct DM or millicharge production with 2m χ = 350 MeV via an off-shell mediator. the photon energy spectrum if the mediator is massive and m A < √ s, or a broader excess if the mediator is massless [3,[54][55][56][57].
However, we note that if the beam energy in a MM experiment is varied over a sufficiently broad range of values, the energy dependence of the new-physics cross section can be extracted, in principle, as shown in Fig.  2. Indeed, the interaction for a sufficiently heavy mediator arises from a higher-dimension operator, so the signal rate grows more prominently with energy. Consequently, the signal's beam energy dependence can be used to distinguish this class of models from on-shell mediator production or millicharge-like production through a virtual light mediator.

IV. BEAM POLARIZATION
In this section we quantify the discriminating power of the incident beam energy and polarization. While the nominal LDMX setup is not designed to support polarization, electron beam polarimetry for O(few-10 GeV) beams is well established and there is no a priori impediment to the measurements we consider here [58,59]. Although detailed studies are ultimately needed to firmly establish the feasibility of a polarized source that can satisfy the other LDMX beam requirements (e.g., ∼ 100 pA currents required to avoid pileup-related backgrounds), such efforts are beyond the scope of the present work.
We define the differential left/right energy asymmetry in terms of a new observable where σ L/R are signal production cross sections for left/right polarized beam particles. For the vector and axial vector interactions it can been shown that the differential asymmetry vanishes since the amplitudes for leftand right-handed electron beam are independent of the polarization.
On the other hand, for chiral vector interactions involving a P L operator, the cross section for a right-handed polarized electron beam vanishes when m e E beam . This is can be seen in Fig. 8 where we plot the differential left/right polarization asymmetry in signal events as a function of electron recoil energy E e for different chiral vector mediator masses, and we only show statistical error bars. We see that there is noticeable asymmetry for all masses, making the chiral vector scenario easily distinguishable from the vector and axial vector scenarios, which do not exhibit an asymmetry.
Note that all potential backgrounds from QED processes (such as photonuclear hadron production) involve vector Lorentz structures, so no asymmetry is expected even if such events cannot be vetoed; electroweak "invisible" backgrounds from direct neutrino production via off-shell Z exchange eN → eNνν will contribute to the asymmetry, but the event rates for these processes are negligible for a Phase 1 LDMX run [14,34].
Independently of how the asymmetry in Eq. (10) is binned, we can also define an inclusive total asymmetry observable where N L/R is the number of signal evens observed using a left/right beam polarized beam, assuming equal luminosities for the two data samples. The uncertainty on the total asymmetry is where we have used standard error propagation and Poisson uncertainties for the signal events δN L/R = N L/R ; to claim nσ evidence of a polarization asymmetry, we require A tot > nδA tot . Note that, while the total asymmetry defined in Eq. (11) is related to A LR (E e ) from Eq. (10), the former is not the integral of the latter. Furthermore, in the limit of a purely chiral interaction (∝ P L or P R ), signal events will overwhelmingly be produced with only one beam polarization, so an asymmetry can be identified with a small number of signal events, as soon as the number of signal events exceeds the Poisson error on the total event count.

V. ELECTRON VS. MUON BEAMS
In the previous section we have seen that the observation of a polarization asymmetry in the electron energy distribution indicates that DM interacts with the SM via a chiral vector mediator, while non-observation of an asymmetry indicates that the mediator has either vector or axial vector interactions (similar conclusions hold for scalar mediators). Distinguishing between vector and axial vector interactions would require modifying the beam in a different way. In particular, the difference between the differential energy distributions of the vector and axial vector mediator scenarios is proportional to mass of the beam particle m 2 . This can be easily seen by considering the − γ → − A process for the vector and axial couplings of A ; 3 in the limit of m / √ s 1 the squared amplitudes differ by where s is the Mandelstam variable for the 2 → 2 process (rather than the N → N + A collision). This expression illustrates two important points: the difference between the axial and vector interactions is amplified for heavier beam particles, and the axial interaction tends to produce more events at lower − recoil energies (this can be confirmed by expressing t in terms of the recoil − energy and p T ). These differences disappear at larger m A or beam energies, both of which increase s. This intuition is borne out in the full MC simulation of a muon missing momentum experiment from Ref. [16] which features a 15 GeV muon beam colliding with a tungsten target. We show the recoiling muon energy distributions in Fig. 9 for various A masses in the vector (upper left panel) and axial models (upper right panel) . We see that for m A m µ , the distributions are visibly different at low recoil energies as expected. An analogous behavior is also shown in the recoil energy distributions for scalar and pseudo-scalar interactions in the lower row of Fig. 9. We therefore conclude that experiments utilizing muon beam and electron beams are complimentary for probing both the flavor and Lorentz structure of BSM interactions.
Although the muon beam missing momentum concept, as demonstrated in [16] for M 3 , is similar to an LDMX-style electron beam experiment, there are some important differences. For example, LDMX anticipates a monochromatic electron beam, utilizes a thin target ( radiation length), and requires ∼ 10 14 − 10 16 electrons on target for phases 1 and 2, respectively. By contrast, muon beams are typically prepared from boosted pion decays with a broader spread of beam energies and lower luminosities; as a result, M 3 specifically is designed to run with ∼ 10 13 total muons delivered to the target. The optimal recoil energy cut to suppress SM backgrounds was also found to be different. Using the M 3 set up as a benchmark we can estimate the number of signal events required to distinguish between axial and vector interactions using the likelihood ratio method described in Sec. III. As before, we generated mock data from our MC samples for A of a given mass and the vector interaction, and then studied the distribution of log likelihood ratios for the axial model (with A of the same mass) and vector model. The result is shown in Fig. 10. As the mediator mass becomes larger, the kinematic distributions of axial and vector interactions become more and more similar, requiring larger event samples to disentangle. This figure also illustrates the importance of experimental resolution on the kinematic quantities. The solid and dashed lines correspond to two different binnings in recoil energy and transverse momentum. Higher resolution (i.e., finer binning) enables model discrimination with a fewer number of observed events, as the spectral differences between models are spread over a larger number of bins.

VI. CONCLUSION
Fixed target missing-momentum experiments involving electron or muon beams are powerful dark matter (and dark sector) discovery tools. In this paper we have studied the model discriminating potential of these techniques using various simplified models categorized according to production topology and Lorentz structure.
In our numerical studies we have varied the masses of BSM particles in the models above and developed statistical tests to discriminate between different signal hypotheses using kinematic variables and beam parameters: lepton recoil energy E , transverse momentum p T, , and incident beam energy and polarization. We have also studied the discriminating potential of using electron vs. muon beams. Our main conclusions can be summarized as follows: 1. Kinematic Variables: For radiative 3-body processes, we find that mediator mass can be well determined using a combination of recoil energy and transverse momentum variables. Using a likelihood analysis, we compared E and p T, distributions against mock data to determine how many signal events are required to distinguish various mediator masses at given confidence levels.
We also find that, over the full MeV-GeV mediator mass range, p T, is a superior kinematic variable and greatly enhances signal discriminating power relative to E -only measurements (e.g., NA64 [52,53]). Indeed, for most test masses in this range, p T, enables clear discrimination between most hypotheses with only O(10) signal events near the nominal discovery threshold for a zerobackground experiment. In contrast, the same analysis using only E distributions typically requires several hundred events or more. These conclusions hold independently of the mediator's identity (e.g., scalar vs. vector) as the p T, distributions are similar across the different particle variations studied here. Consequently, E and p T, are good variables for determining mass assuming a given BSM scenario, but these alone are not effective at discriminating between different models. This conclusion also applies to efforts to distinguish between on-shell mediator production in 2 → 3 processes and 2 → 4 processes that produce DM pairs through either light or heavy mediators; for an equivalent invariant mass of new states, the kinematic distributions of these models are not easily distinguished. This highlights the complementarity of different accelerator experiments that can potentially resolve some of these ambiguities with higher center-of-mass energies.

Beam Polarization:
We also find that for different electron-vector Lorentz structures, there are new asymmetry observables that can discriminate between certain scenarios. In particular, subtracting the energy distributions for signal events with incident left and right-handed electrons yields residuals for chiral electron-mediator couplings of the form Z µē γ µ P L,R e, where P L,R is a left/right projector; if the interaction is a general linear combination of vector and axial-vector couplings, these observables only yield nonzero residuals for the chiral component proportional to P L or P R . Although the feasibility of beam polarization in missing momentum experiments has not yet been studied, the discriminating power of this method warrants future work to assess its compatibility with other accelerator requirements (e.g., beam current and structure).

Beam Flavor:
Finally, we find that a combination of electron and muon beam missingmomentum searches can be used to determine whether a BSM interaction is parametrically enhanced by the mass of the beam particle. As representative examples, we consider pseudoscalar a¯ γ 5 and axial-vector V µ¯ γ µ γ 5 interactions for which N → N a/V emission requires a lepton chirality flip. Consequently, the corresponding cross sections for these processes are proportional to m 2 and the difference between electron and muon beam signals can be used to distinguish these models from scenarios that do not require a chiral flip in order to produce the new BSM state.
Our motivation in this paper has been to study model discrimination power of fixed target missing-momentum experiments in a variety of well motivated scenarios. However, the studies should be interpreted with great care as we only consider simulated signal distributions with statistical uncertainties to assess the distinguishability of various BSM scenarios assuming negligible backgrounds from SM processes. Since such experiments are expected to have negligible irreducible backgrounds and low levels of reducible backgrounds from photonuclear reactions [14,16] this is a well motivated working assumption. However, detailed future studies should include systematic uncertainties, detector level smearing effects, and various levels of background contamination to more realistically quantify the model discrimination power of the techniques outlined in this paper.

ACKNOWLEDGMENTS
In this appendix we describe the likelihood analysis used in Sec III to distinguish between particle masses using final state beam recoil energy E e and transverse momentum p T,e . The simple likelihood given in Eq. (7) correctly accounts for statistical fluctuations in the observations. However, it assumed that each bin contains a large number of simulated events, i.e., the relative fluctuations in the bin counts are small. The predictions for missing recoil energy and transverse momentum of the electrons are obtained using finitely-sized MC samples. Even if the sample is large, if we bin it finely enough (and especially if we bin in multiple kinematic variables at the same time) each bin will contain a small number of events which are subject to Poisson fluctuations due to the finite size of the MC sample. If we are comparing two distributions that only differ in such bins where fluctuations are important, we miss the theoretical uncertainty due to MC statistics. This uncertainty can be incorporated following Ref. [51] (see also Ref. [60] for a concise description of the problem). The idea is to introduce nuisance parameters that encode the true expected values in each bin; the count of MC events is then a specific random Poisson realization of this expected value. Since we cannot evaluate the cross-section exactly in that kinematic bin, we do not know what this true value is, and therefore we must marginalize over it. The likelihood function that takes this effect into account is where a i are the "raw" (unnormalized) MC bin counts (i.e., the number of events in each kinematic bin as it comes out of the simulation -these numbers grow as the MC sample gets larger); A i are the true expected values in each bin (equal to a i in the limit of an infinitely large MC sample, and therefore not known); p is the signal strength, such that pA i corresponds to the theory prediction of the expected count in bin i. A i and p are nuisance parameters over which we must maximize the likelihood -luckily we will be able to do this analytically. The notation used here corresponds to that of Ref. [51]. In other words, Eq. (A2) treats the MC sample as another dataset; both the "real" data d and the MC data a constrain p and A. Next we maximize this with respect to p and A i : The solutions of these equations are Note that the signal strength p takes the natural value which ensures that the predicted number of events i pA i matches the number of observed data events N d . At the maximum of the likelihood, A i = a i , unless a i d i ; this encodes the fact that for a finite MC sample we do not quite know what the actual expected value of the bin counts is. Plugging these solutions into Eq. (A3) we find ln L ≡ max p,Ai ln L The last line contains factors that are commonly dropped when discussing real data (i.e., d); however since a i depends on the model, we must keep this factor such that we can meaningfully compare the likelihoods in different models (which generically have differently-sized MC samples). d i !, however, is the same in each model, so we can drop it. The likelihood in Eq. (A8) is not very intuitive so it is useful to consider some illustrative limits: • a i d i ∀i, and therefore N mc N d . Expanding in small quantities one finds which is just the standard Poisson likelihood of Eq. (7) with ν i = (N d /N mc )a i (there are corrections that go like ln a i from approximating ln a i !). Thus in the limit in which theory (MC statistical) uncertainties are not important we recover the naive result.
• Without taking theory uncertainties into account, if MC predicts a i = 0 events in a bin, but the data d i is not 0 there, the standard log likelihood for that bin is d i ln 0 − d i = −∞, i.e., the model is immediately ruled out. Using the above log likelihood instead, we find that such a bin would instead contribute to the log likelihood. The would-be infinity is regularized by the fact that N d /N mc is not 0. This is the desired behavior, since it prevents us from overstating the discriminating power of certain bins if we have insufficient MC statistics.
The full likelihood of Eq. (A8) is used to define the test statistic in Eq. (8), which is then studied in the various model discrimination examples in Sections III and V.