Dark matter signals from timing spectra at neutrino experiments

TL;DR

This paper introduces a new timing-based method to detect dark matter signals at neutrino experiments, demonstrating its potential with the COHERENT data and exploring future sensitivities.

Contribution

It proposes a novel resonance-like timing signature for dark matter detection and applies it to existing data, providing new insights and future experimental prospects.

Findings

Detected an excess in timing distribution consistent with dark matter hypothesis

Demonstrated the method's ability to distinguish signal from background

Compared sensitivity of current and future experiments for dark matter detection

Abstract

We propose a novel strategy to search for new physics in timing spectra, envisioning the situation in which a new particle comes from the decay of its heavier partner with a finite particle width. The timing distribution of events induced by the dark matter particle scattering at the detector may populate in a relatively narrow range, forming a "resonance-like" shape. Due to this structural feature, the signal may be isolated from the backgrounds, in particular when the backgrounds are uniformly distributed in energy and time. For proof of the principle, we investigate the discovery potential for dark matter from the decay of a dark photon in the COHERENT experiment, and show the exciting prospects for exploring the associated parameter space with this experiment. We analyze the existing CsI detector data with a timing cut and an energy cut, and find, for the first time, an excess in the timing distribution which can be explained by such dark matter. We compare the sensitivity to the kinetic mixing parameter ($\epsilon$) for current and future COHERENT experiments with the projected limits from LDMX and DUNE.

Figures

Tables

Table 1. Table 1: Best-fit ϵ q superscript italic-ϵ 𝑞 \epsilon^{q} for a few M ′ superscript 𝑀 ′ M^{\prime} values (in MeV) for the single-mediator scenario.

$M^{\prime }$	50	75	100	1000
${ϵ}^q$	$3.5\times {10}^{-4}$	$4.4\times {10}^{-4}$	$5.5\times {10}^{-4}$	$4.6\times {10}^{-3}$

Equations

$${\mathcal{L}}_{\mathrm{int}}\supset g_{\chi }{A}_{\mu }^{\prime }\bar{\chi }{\gamma }^{\mu }\chi +e_q{ϵ}_1^q{A}_{\mu }^{\prime }\bar{q}{\gamma }^{\mu }q,$$

$$\frac{d^2{N}_{A^{\prime }}}{dtd\cos \theta }=g(\cos \theta )\cdot \frac{1}{{\tau }_{A^{\prime }}}{e}^{-\frac{t-t_F}{{\tau }_{A^{\prime }}}}\Theta (t-t_F),$$

$$f(T)\propto \int d\cos \theta {\left|\frac{dT}{dt}\right|}^{-1}\frac{d^2{N}_{A^{\prime }}}{dtd\cos \theta }.$$

$$T=t+\frac{1}{v_{\chi }}\sqrt{x_0^2+v_{A^{\prime }}^2(t-t_F{)}^2-2x_0{v}_{A^{\prime }}(t-t_F)\cos \theta }$$

$$\frac{d\sigma }{d{E}_r}$$

$$Excess=\frac{\mathrm{signal}-\mathrm{SS}-\mathrm{BRN}-\mathrm{SM}}{\sqrt{2\mathrm{SS}+\mathrm{BRN}+\mathrm{SM}}}.$$

Full text

††preprint:

MI-TH-1925

Numerous theoretical and experimental ideas have been put forth to identify the mass and associated interactions of DM candidate particles. Since traditional WIMP-based searches have not yet detected DM Aprile et al. (2018), expanding the search of parameter space is well-justified Battaglieri et al. (2017). Many models of light DM ($\lesssim$ GeV) emerge from a hidden/visible sector where light mediators (e.g., a dark photon) interact with DM Huh et al. (2008); Pospelov et al. (2008); Hooper and Zurek (2008); Cheung et al. (2009); Essig et al. (2010, 2013a); Dutta et al. (2019a). Because the DM mass is light in these models, it is difficult to detect such DM in traditional WIMP-based direct detection experiments.

In this paper we develop a novel strategy to search for light DM which couples with light mediators and apply it to the data from the ongoing COHERENT experiment Akimov et al. (2018a). COHERENT makes use of a proton beam which impinges on a Hg target at the Spallation Neutron Source (SNS). Among the produced pions, $\pi^{+}$ decays create prompt muon neutrinos and delayed anti-muon and electron neutrinos. The measured energy spectra have been used to investigate new physics associated with neutrino non-standard interactions (NSI) Ohlsson (2013); Miranda and Nunokawa (2015) due to heavy or light mediators Coloma et al. (2017a, b); Liao and Marfatia (2017); Dent et al. (2018); Billard et al. (2018); Lindner et al. (2017); Farzan et al. (2018); Brdar et al. (2018), generalized scalar and vector neutrino interactions Aristizabal Sierra et al. (2018), hidden sector models Datta et al. (2019), and sterile neutrinos Kosmas et al. (2017); Blanco et al. (2019). It also sets independent constraints on the effective neutron size distribution of CsI Ciuffoli et al. (2018); Aristizabal Sierra et al. (2019); Papoulias et al. (2019). Since the proton beam is pulsed, the measured timing spectra may be used to distinguish between prompt and delayed events. The combined timing and energy spectra have been utilized to understand new physics models with neutrino flavor-dependent NSI Dutta et al. (2019b).

We show how both the timing and energy data from the COHERENT experiment can be used to search for light, $\lesssim 1$ GeV, DM. The DM event under consideration is initiated by the production of a dark photon decaying into a pair of DM particles (e.g., Refs. deNiverville et al. (2015); Ge and Shoemaker (2018)). A DM particle would then induce a nuclear recoil event at the detector. The dark photon production can occur from both $\pi^{-}$ and $\pi^{0}$. Most of the $\pi^{-}$ are stopped inside the Hg target and can create a dark photon via the absorption process, $\pi^{-}+p\rightarrow n+A^{\prime}$, followed by the decay of the dark photon $A^{\prime}$ to a DM pair deNiverville et al. (2015). The dark photon is emitted isotropically in this $\pi^{-}$ absorption process. The $\pi^{0}$ may produce an ordinary photon and a $A^{\prime}$ Ge and Shoemaker (2018). Since the $\pi^{0}$ move somewhat relativistically, the resulting DM lies relatively in the forward direction. Nevertheless, we find that the DM flux reaching to the COHERENT detectors, which are located $\sim 90^{\circ}$ from the beam direction Akimov et al. (2018a) (see Appendix), is comparable to that from the $\pi^{-}$ absorption. Further, there are additional contributions from $\pi^{\pm}+p/n\rightarrow n/p+\pi^{0}$.

The method that we develop to search for DM utilizes both the energy and timing spectra of the DM-initiated nuclear recoil events. We focus on the timing and energy spectra for the DM produced from the $\pi^{-}$ absorption and the $\pi^{0}$ decay. In the COHERENT experiment, the $\pi^{-}$ ($\pi^{+}$) and $\pi^{0}$ abundances per proton on target are 0.05 (0.11) and 0.1, respectively deNiverville et al. (2015); Akimov et al. (2018b); Coh (Rebeca Rapp). The produced $A^{\prime}$ is mostly relativistic unless its mass is $\sim$ 138 MeV.

We first derive the timing spectrum of DM-induced nuclear recoil events along with their energy distribution, and then compare the DM case to that of Standard Model (SM) neutrinos. The signal under consideration is initiated by production of a dark photon $A^{\prime}$ from the decay of the $\pi^{-}$-$p$ mesic state and $\pi^{0}$ decay through kinetic mixing. $A^{\prime}$ production and its subsequent decay to DM $\chi$ are governed by the following interaction Lagrangian:

[TABLE]

where $e_{q}=eQ_{q}$, $g_{\chi}$ and $\epsilon_{1}^{q}$ are dark-sector gauge coupling and kinetic mixing parameter (associated with the mixing between the $\gamma$ and new gauge Boson ${\epsilon\over{2}}F^{\mu\nu\prime}F_{\mu\nu}$ Holdom (1986); del Aguila et al. (1988); Babu et al. (1998)), respectively. This generic-looking Lagrangian can be accommodated in the context of a model, e.g.,  deNiverville et al. (2015); Dutta et al. (2019a).

Let us suppose that $A^{\prime}$ is produced at $t_{F}$ where $t_{F}$ is the timing of $\pi^{-,0}$ production induced by the 1 GeV SNS beam which is $0.6~{}\mu$s wide and pulsed at 60 Hz. We then assume that $A^{\prime}$ flies for $v_{A^{\prime}}(t-t_{F})$ along the $\theta$ direction with respect to the line joining the Hg target and the detector (see Appendix), and decays to a $\chi$ pair. One of the $\chi$’s then may travel towards the detector for $v_{\chi}t^{\prime}$. Denoting the timing measured at the detector by $T$, we see that $T$ is the sum of $t$ and $t^{\prime}$, i.e., $T=t+t^{\prime}(v_{A^{\prime}}(t-t_{F}),t-t_{F},\cos\theta)$, where we explicitly express $t^{\prime}$ as a function of $t-t_{F}$ and $\cos\theta$. We are interested in the differential number of events in $T$, or equivalently the DM flux at the detector of interest, $f(T)=dN_{\chi}/dT$. Parameterizing the angular distribution of dark photons by $g(\cos\theta)$, we find

[TABLE]

where $\Theta(x)$ is the step function. $g(\cos\theta)$ is 1/2 for dark photons from the $\pi^{-}$ absorption. We then obtain

[TABLE]

A simple geometry consideration gives

[TABLE]

with $x_{0}$ being the distance between the Hg target and the detector. We consider both $\pi^{-}$ and $\pi^{0}$ contributions. In our calculation, we use the GEANT4 Agostinelli et al. (2003) simulations for the COHERENT geometry to determine the angular and energy spectra of photons from pion absorption and decays Coh (Rebeca Rapp). We find that the $\pi^{0}$ contribution is bigger than the $\pi^{-}$ absorption.

The top panel of figure 1 demonstrates example timing spectra for a CsI detector, with three different choices for the rest-frame mean lifetime of $A^{\prime}$. The solid and dashed histograms are for a relativistic dark photon ($m_{A^{\prime}}=75$ MeV) and a non-relativistic dark photon ($m_{A^{\prime}}=138$ MeV), respectively, with $m_{\chi}$ fixed to 5 MeV. Here the $\pi^{-}$ flux – which is approximated by a Gaussian distribution with a mean value of 0.7 $\mu$s and a width of 0.15 $\mu$s to model the arriving time of the proton on target which reproduces the timing spectrum in Akimov et al. (2018b) – is convoluted. For the non-relativistic case, most of the $\chi$’s can reach the detector (modulo a factor of $(4\pi x_{0}^{2})^{-1}$). Not surprisingly, as $A^{\prime}$ is shorter-lived, the spectrum width gets narrower, manifesting in a resonance-like bump feature more visibly. By contrast, for the relativistic case, if $A^{\prime}$ is long-lived, it decays far away from the detector so that only a small fraction of the $\chi$’s can reach the detector, contributing to the upper tail of the spectrum. Therefore, relatively short-lived $A^{\prime}$ would give more statistics. Indeed, we see that most of DM events populate within $\sim 1.5~{}\mu$s which roughly corresponds to the mean value plus the width of the beam pulse. Note that prompt neutrinos leave events within $\sim 1.5~{}\mu$s whereas delayed neutrinos spread out over a broad range Akimov et al. (2018a, b). So, requiring $T\lesssim 1.5~{}\mu$s essentially rejects most of delayed neutrino events while a large portion of prompt neutrino events and relativistic (non-relativistic) $A^{\prime}$-induced DM events irrespective of $\tau_{A^{\prime}}$ (with $\tau_{A^{\prime}}\lesssim 0.1~{}\mu$s) are kept.

Regarding DM-nucleus scattering, we remark that in principle DM scattering can be governed by physics different from that for dark photon production encoded in Eq. (1). Introducing a generic mediator of mass $M^{\prime}$, DM-mediator coupling $g_{D}$, and quark-mediator coupling $e_{q}\epsilon_{2}^{q}$, we find that the differential spectrum in recoil energy $E_{r}$ of the target nucleus can be expressed as

[TABLE]

where $F$ denotes the form factor and where $Z$ and $m_{N}$ are the atomic number and the mass of the target nucleus. The underlying interaction is of dark-photon type for illustration. We neglected $m_{\chi}$ in the curly brackets as $m_{N}\gg m_{\chi}$. Clearly, the spectral behavior is (nearly) independent of $m_{\chi}$. The bottom panel of figure 1 displays the expected nuclear recoil spectrum for $(M^{\prime},m_{\chi})=(75,5)$ MeV (green). For comparison we show the $E_{r}$ distributions of the prompt neutrinos (blue) and the delayed neutrinos (orange) with (solid) and without (dashed) experimental efficiencies. We see that prompt neutrino events occur almost entirely in the region $E_{r}\lesssim 14$ keV, so employing a lower cut at this energy removes the remaining prompt neutrinos, while retaining a large portion of the DM candidate events.

In order to analyze the COHERENT data using both the energy and timing spectra Akimov et al. (2018b) from neutrinos and DM, we adopt the statistical method described in Ref. Dutta et al. (2019b). We allow for Poisson fluctuations of the background in each energy and time bin (model (c) of Ref. Dutta et al. (2019b)), and fix the size of the neutron distribution to $R_{n}=4.7$ fm. We also quote our results for $R_{n}=5.5$ fm which is the model independent central value obtained from the fit to the COHERENT data Cadeddu et al. (2018). We examine two limiting cases: i) the specific part of the energy and timing data in which the DM signal is predicted to appear, after removing as many neutrino-induced events as possible, and ii) the full energy and timing data.

As discussed in the previous section, we apply cuts $E_{r}>14$ keV (16 photoelectrons) and $T<1.5~{}\mu$s to substantially suppress both prompt and delayed neutrino events, but keep the DM candidate events 111Unless $m_{A^{\prime}}\approx 138$ MeV and $\tau_{A^{\prime}}\gtrsim 0.03~{}\mu$s. for the published COHERENT data Akimov et al. (2018b). We also apply an upper-cut $E_{r}<26$ keV since the background is well understood for COHERENT in the range 5 to 26 keV Akimov et al. (2018b). The experimental efficiency is also given in Akimov et al. (2018b). After these cuts, we find 97 total events. Out of them 49 events have been classified as the steady-state (SS) background, while 19 may be identified as delayed neutrino events forming the SM (i.e., neutrino) background. There are also 3 events in the cut window arising from beam related neutron (BRN) backgrounds. There is then an “excess” of 26 events which corresponds to a 2.4$\sigma$ statistical uncertainty. For $R_{n}=5.5$ fm, the significance becomes $\sim$3$\sigma$. For calculating the significance we apply Scholberg (2018)

[TABLE]

We also calculate the significance from the likelihood ratio test for the DM fit to the excess and find the significance to be 1.98. The SS and BRN backgrounds emerge from measuring beam-on anti-coincident events and the simulation from GEANT4, respectively. We use the same systematic uncertainty $\sim 28\%$ Akimov et al. (2017) which incorporates flux, form factor, quenching factor and the signal acceptance uncertainties.

We first attempt to explain the excess with a DM hypothesis, again assuming that the DM scattering is governed by a different mediator. We fit the selected events, varying the associated (effective) coupling constant $\epsilon$ and mediator mass $M^{\prime}$ which is responsible for the interaction between the DM and the nucleus. The left panel of figure 2 shows $1\sigma$-best fits to the data set with the cuts implemented (blue band). The scattered “islands” come from our MCMC sampling methods, with limited live points which cause some “islands” in the result. For comparison, the orange band shows the parameter space when performing a fit to the full energy and timing data at $1\sigma$. We see that there exists an overlapping region, and further find that both “before-cut” and “after-cut” data sets are well accommodated by the parameter points with $M^{\prime}\gtrsim 100$ MeV. For comparison to the DM case, we determine whether a NSI neutrino hypothesis is able to fit both the before-cut and after-cut data. For the neutrino case we consider a non-zero coupling $g_{e}$, the NSI in the $\nu_{e}$ neutral-current interaction. As shown in the right panel of figure 2, it is not possible to simultaneously fit both the before-cut and after-cut data sets with this neutrino hypothesis. In fact, this NSI model does not show a good fit for the excess in the prompt timing bin (i.e., $T<1.5~{}\mu$s). The fuzzy region at low $g_{e}$ shows at there is some statistical consistency with the SM in this region, in particular for the before-cut data. The situation becomes even worse with $g_{\mu}\neq 0$, since it affects not only the delayed but also the prompt spectrum.

In the DM case, the parameter $\epsilon$ is defined as $\epsilon=\epsilon^{q}_{1}\epsilon^{q}_{2}\epsilon_{D}\sqrt{{\rm BR}_{A^{\prime}\rightarrow\chi\chi}}$, where $\epsilon^{q}_{1}$ is the $q$-$A^{\prime}$ kinetic mixing which describes the dark photon production from the $\pi^{-}$ absorption, $\epsilon^{q}_{2}$ is the quark-mediator kinetic mixing for the DM-nucleus scattering cross-section, and $g_{D}=e\epsilon_{D}$ is the DM-mediator coupling. This is the most general description, since in a realistic model there can be more than one mediator, e.g., scalar and gauge boson mediators commonly occur in models with spontaneous symmetry breaking. Of course, the best-fit contour can also be interpreted in the case where there exists only a single mediator, i.e., $M^{\prime}=m_{A^{\prime}}$.

The parameter choices that we use to obtain the best-fit points are $\tau_{A^{\prime}}=1$ ns, $m_{A^{\prime}}=75$ MeV and $m_{\chi}=5$ MeV. However, we find that the best-fit points do not change in the $\epsilon-M^{\prime}$ plane for the following variations: i) for $\tau_{A^{\prime}}\lesssim 4$ ns, since the DM flux maximizes for $\tau_{A^{\prime}}\lesssim 4$ ns with $m_{A^{\prime}}<138$ MeV, ii) for the non-relativistic case, i.e., $m_{A^{\prime}}=138$ MeV, with $\tau_{A^{\prime}}\lesssim 30$ ns and iii) for any $m_{\chi}$ smaller than $m_{A^{\prime}}/2$. For relativistic scenarios with large $\tau_{A^{\prime}}~{}(\geq 4~{}{\rm ns})$, the best-fit regions get scaled by the appropriate associated DM flux (see the top panel of figure 1). For non-relativistic scenarios with large $\tau_{A^{\prime}}~{}(\geq 30~{}{\rm ns})$, it is not possible to fit before-cut and after-cut data sets simultaneously because DM will contribute to both before $1\mu$s and after $1\mu$s events. Figure 2 is shown for $R_{n}=4.7$ fm. However, the best fit contours do not change for $R_{n}=5.5$ fm.

Based on the above discussions, we describe the best-fit parameters for the following two scenarios.

• •

Single-mediator scenario: In this case, $\epsilon^{q}\equiv\epsilon_{1}^{q}=\epsilon_{2}^{q}$ and the dark photon $A^{\prime}$ should decay fast. Otherwise, $\epsilon$ is small, meaning that the DM-nucleus scattering is so small that a very small number of events would occur. Here $\epsilon=\epsilon^{q}_{1}\epsilon^{q}_{2}\epsilon_{D}\sqrt{{\rm BR}_{A^{\prime}\to\chi\chi}}\to({\epsilon^{q}})^{2}\epsilon_{D}\sqrt{{\rm BR}_{A^{\prime}\to\chi\chi}}$. We can choose $\epsilon_{D}=1/e$ to make $g_{D}=1$ which makes $\tau_{A^{\prime}}$ small, and we can still make use of the left panel of figure 2 (where $\tau_{A^{\prime}}$ is set to be $\leq 1$ ns). Table 1 shows the best $\epsilon^{q}$ for a few $M^{\prime}(=m_{A^{\prime}})$ values for which the resulting branching ratios for $\pi^{0}\to\gamma A^{\prime}$ and $\pi^{\pm}\to e^{\pm}\nu A^{\prime}$ agree with current precision data of $\pi^{0}$ and $\pi^{\pm}$ Tanabashi et al. (2018). We do not report any numbers below $M^{\prime}=50$ MeV as we find that the best-fit region with the before-cut data does not overlap with that with the after-cut data.

• •

Multi-mediator scenario: Unlike the previous scenario, $\tau_{A^{\prime}}$ is not necessarily small, since $\chi$ scatters off the target nucleus via a new mediator with large coupling while the dark photon can decay to a pair of DM particles with a longer lifetime. Table 1 for a single mediator scenario still holds with $\epsilon^{q}$ identified as $\sqrt{\epsilon^{q}_{1}\epsilon^{q}_{2}\epsilon_{D}e}$.

The values of $\epsilon^{q}$ shown in Table 1 are obtained assuming that the dark photon couplings to up and down quarks are proportional to their charges. If, however, we want to use the universal charge (e.g., 1), then we need to scale the $\epsilon^{q}$ by $\sqrt{2Z/(9A)}$ where $Z=54$ and $A=130$ for CsI. The best-fit values of $\epsilon^{q}$ are below any existing bounds Essig et al. (2013b); Dror et al. (2017) arising from meson decays, e.g., $K\rightarrow\pi+$invisibles Artamonov et al. (2009). The model details become important for this constraint, i.e., whether it contains fully conserved current, additional Higgs sector Davoudiasl et al. (2014) and the value of $g_{D}$, etc. The excess can be explained in the allowed regions of parameter space of the realistic models, e.g., $U(1)_{T3R}$, $U(1)_{B-L}$, etc. deNiverville et al. (2015); Dutta et al. (2019a). For example, for 50 MeV dark photon, the coupling ($e\epsilon$)$\sim 10^{-4}$, needed to explain the excess, is well allowed by all the existing data in a generic vector-portal DM model deNiverville et al. (2015). The COHERENT limit for NSI of neutrinos is better than any existing limit from various experiments using the timing plus energy data where the SM backgrounds cannot be sufficiently suppressed Dutta et al. (2019b). However, for the DM analysis, since we have vetoed the SM neutrino backgrounds using the energy and timing cuts, we can obtain an even better reach in terms of new physics coupling. Now assuming no excess above the backgrounds (which could be caused by an improper estimation of the time of creation of the pions at the target 222The COHERENT collaboration is investigating this possibility [Kate Scholberg, private communication].), the values of $\epsilon$ become smaller by a factor of 1.5 compared to those in Table 1, as shown in figure 3.

The future LDMX experiment Åkesson et al. (2018) will investigate the sub-GeV DM parameter space which arises from a dark photon decaying to DM, using an electron beam dump. We note that this parameter space is already being probed via nuclear recoils at COHERENT, therefore representing a complementary approach. In figure 3, we compare the reach of $(\epsilon^{X})^{2}$ as a function of mediator mass for the current COHERENT data and for a future argon detector with the LDMX reach assuming that $\epsilon^{q}=\epsilon^{e}$. We also show the existing limits from NA64 Banerjee et al. (2019) relevant to both single- and multi-mediator scenarios in the plot assuming $\epsilon^{q}=\epsilon^{e}$.

Our current and projected limits are derived using the formalism of Ref. Datta et al. (2019), and they are essentially governed by the $\pi^{0}$ contribution. We show two scenarios: (i) the dark photon coupling ($\epsilon^{q}_{1}$) is the same as the mediator-nucleus coupling ($\epsilon^{q}_{2}$) and (ii) $\epsilon^{q}_{1}$ is fixed at $10^{-2}$ (current experimental constraint Artamonov et al. (2009)) with $\alpha_{D}\equiv g_{D}^{2}/(4\pi)=0.5$. We use a dark photon mass $m_{A^{\prime}}=75~{}\mathrm{MeV}$ and a DM mass $m_{\chi}=5~{}\mathrm{MeV}$. The figure, however, is unchanged for $m_{A^{\prime}}\leq 138$ MeV, $m_{\chi}\leq m_{A^{\prime}}/2$ and $\tau_{A^{\prime}}\leq 4$ ns. We also note that the reach of the current COHERENT data in probing $(\epsilon^{X})^{2}$ in figure 3 is competitive with DUNE experiment reach De Romeri et al. (2019).

In conclusion, we have argued that the timing information available in neutrino experiments with pulsed beam such as the COHERENT data is a powerful probe of new physics. We have shown that the combination of energy and timing cuts can eliminate SM neutrino events very efficiently, thereby allowing the possibility of isolating DM-induced events. As applied to the published COHERENT data, we find a considerable number of excess events over the expected backgrounds. This excess of events may be explained by a dark matter hypothesis, and is unlikely to be explained by SM neutrino interactions. We note that this conclusion is distinct from the results presented in Ref. Dutta et al. (2019b), who showed that using the full energy and timing data without using the cuts, a neutrino model is able to explain the data. Even though we have presented a DM interpretation of the COHERENT data, it remains possible that the events may be explained by an unidentified background, by a systematic uncertainty on the observed steady-state background or by exotic beyond the SM scenarios. Distinguishing a background hypothesis from a DM hypothesis may be possible with timing and energy information on individual nuclear recoil events. Our analysis strategy can be used to understand dark photon decaying to DM in similar COHERENT type set-ups with timing measurements 333As far as the beam pulse duration is less than $\sim 1~{}\mu$s to populate DM signals in the prompt timing bins and the pulse interval is greater than $\sim 10~{}\mu$s to sufficiently attenuate the flux of delayed neutrinos., e.g., JSNS2 Ajimura et al. (2017), CCM Aguilar-Arevalo et al. (2018).

Acknowledgments

We thank Phil Barbeau, Pilar Coloma, Yuri Efremenko, Pedro Machado, Grayson Rich, and Kate Scholberg for useful discussions. BD, DK and LES acknowledge support from DOE Grant de-sc0010813. The work of DK was supported in part by the Department of Energy under Grant No. DE-FG02-13ER41976/DE-SC0009913. SL acknowledges support from COS-STRP (TAMU). The work of JCP is supported by the National Research Foundation of Korea (NRF-2019R1C1C1005073 and NRF-2018R1A4A1025334). The work of SS is supported by the National Research Foundation of Korea (NRF-2017R1D1A1B03032076 and in partial by NRF2018R1A4A1025334). This work was supported by IBS under the project code, IBS-R018-D1 (SS).

Appendix

Brief description of COHERENT experiment: The SNS in Oak Ridge National Laboratory produces neutrons using a 1.4 MW, 1 GeV proton beam pulsed at 60 Hz on a mercury (Hg) target. A typical proton beam trace has a full-width at half-maximum (full-width) of 360 ns ($\sim 600$ ns). Various detectors are placed in the neutrino alley (see Fig. 10 of Ref. Akimov et al. (2018a) for their exact locations). For example, in Fig. 4, we show the experimental configuration with the CsI detector at $90^{\circ}$ with respect to the beam direction. The dark photon flies at direction $\theta$ and decays to dark matter candidates which are then detected at the detector.

Data analysis: In Fig. 5, we show the significance of excess in the COHERENT data Akimov et al. (2018b) over the SM and other background (beam related neutron and steady state) varying energy and timing cuts. We find the significance of the excess is maximized for $T<1.5\,\mu$s and $E_{r}>12$-$14$ keV. We use this set of cuts for our analysis.

Now, the region of $T<1.5\,\mu$s is preferred by dark matter events as we see from Fig. 1 (top panel) of the main text and the neutrinos emerging from the pion decay ($\pi^{+}\rightarrow\mu^{+}+\nu_{\mu}$) in the SM (but not the delayed neutrinos emerging from the $\mu$ decay with a lifetime being 2.2 $\mu$s), while our selected region with $E_{r}>14$ keV and $T<1.5\,\mu$s is only preferred by the DM candidate particles since the DM candidates have more energy compared to the $\nu_{\mu}$ emerging from the pion decay ($<14$ keV).

Dark photon angular and energy distributions: In Fig. 6, we show the unit-normalized angular distributions of dark photons from $\pi^{0}$ decays and $\pi^{-}$ absorption. We use GEANT4 10.5 with the FTFP_BERT library Agostinelli et al. (2003) to simulate the case of $\pi^{0}$ decays. On the other hand, the produced monochromatic photon via the $\pi^{-}$ absorption process induces nuclear fragmentation processes, eventually resulting in a bunch of soft photons ($\lesssim 20$ MeV). We take into account such photons in our data analysis, although their contribution to our signal sensitivities is almost negligible since the nuclear $E_{r}$ produced by the DM emerging from the decays of dark photon produced by these soft photons is about 1 keV which is outside the range $14<E_{r}<26$ keV considered for the events for the DM analysis.

Further, once a dark photon is monochromatically created in the $\pi^{-}$ absorption process instead of a photon, it hardly participates in further nuclear disintegration processes because its interaction with nucleons inside the mercury atom is very weak/negligible due to the kinetic mixing suppression.

$\pi^{-}$ may get involved in the charge exchange process, $\pi^{-}+p\to\pi^{0}+n$, as mentioned in the main text. Our GEANT4 simulation suggests that most of $\pi^{-}$ do not result in $\pi^{0}$, so we assume that the absorption process happens to them. Therefore, we develop our own simulation code to take care of these single-energy-valued dark photons from the $\pi^{-}$ absorption case. However, $\pi^{0}$ decays provides the dominant contribution compared to the $\pi^{-}$ absorption for the dark photon parameter space limits since we have about $2$ times more $\pi^{0}$ compared to $\pi^{-}$ per proton on target. The dark photon is emitted isotropically in this $\pi^{-}$ absorption process whereas the dark photons emerging from the $\pi^{0}$ decay have a relatively large contribution in the forward direction.

In the top panel of Fig. 7, we show the unit-normalized dark photon energy spectra from $\pi^{0}$ decays and $\pi^{-}$ absorption. By contrast, in the bottom panel of Fig. 7, we show the number of DM events per second per nucleus from the $\pi^{0}$ decays and the $\pi^{-}$ absorption process where the former surpasses the latter by at least a factor of 2 at the detector.