1. Introduction
More than a century after evidence for the existence of Dark Matter (DM) in our universe started to emerge [
1], we still do not know what DM is made out of. Weakly Interacting Massive Particles (WIMPs) remain one of the best-motivated DM candidates to date [
2,
3,
4]. Among the most sensitive probes of the WIMP DM paradigm are
direct detection experiments [
5,
6,
7,
8,
9], where one searches for interactions of DM with atomic nuclei. Over the past few decades, direct detection experiments have made impressive progress: liquid noble gas detectors have reached
(tonne−yr) exposures with nuclear recoil energy thresholds as small as a few keV [
10,
11,
12,
13], while cryogenic bolometric detectors have pushed the nuclear recoil energy threshold below 100 eV, albeit as of yet only with detectors accumulating
(kg days) of exposure [
14,
15,
16,
17,
18]. A large experimental program is underway to extend the sensitivity of direct detection experiments both to lower WIMP masses and smaller WIMP-nucleon cross sections [
19,
20,
21,
22,
23], as well as to develop directional detectors [
24,
25,
26,
27,
28,
29,
30,
31,
32,
33,
34]. However, these detectors are becoming ever more costly and more challenging to operate.
Paleo-detectors [
35,
36,
37] are a proposed alternative technique to search for interactions of DM with nuclei: instead of operating a nuclear recoil detector in real time, as in a conventional direct detection experiment, one would search for the persistent traces of WIMP–nucleon interactions recorded in natural minerals over timescales as long as
Gyr. Many minerals are excellent solid state nuclear track detectors [
38,
39,
40,
41] and retain these traces (i.e., damage tracks) caused by a recoiling atomic nucleus for timescales that can exceed the age of the Earth by many orders of magnitude (i.e., ≫4.5 Gyr). Natural minerals commonly found on Earth are up to
Gyr old and can be dated to a few-percent accuracy [
42,
43,
44]. Modern microscopy technologies, including electron microscopy [
40,
45], Helium Ion Beam Microscopy (HIBM) [
46,
47], Small Angle X-ray scattering (SAXs) [
48,
49,
50], and optical microscopy [
41,
51,
52], potentially allow for the readout of such nuclear damage tracks in macroscopic samples of target material. For example, HIBM may allow one to read out
mg of material with
nm precision when combined with pulsed-laser and fast-ion-beam ablation techniques [
53,
54,
55,
56]. On the other hand, SAXs promises spatial resolutions of
nm in
g of target material [
49,
50] (see the discussion in References [
35,
36]). These advanced microscopy techniques enable paleo-detectors to reach much larger exposures with better track length resolutions than previous ideas to probe rare events using natural minerals; see References [
57,
58,
59,
60,
61,
62,
63,
64,
65,
66,
67,
68,
69,
70,
71] and, in particular, the work by Snowden-Ifft et al. searching for damage tracks from WIMP–nucleus interactions in ancient phlogopite mica samples [
72,
73,
74,
75]. While not the focus of this article, paleo-detectors have also been proposed as a probe of astrophysical neutrino fluxes from Galactic supernovae [
76], cosmic ray interactions in Earth’s atmosphere [
77], and our Sun [
78]. There has been recent related work on using crystal-defects created in certain minerals by WIMP–nucleus interactions to search for light DM [
79,
80,
81] or to monitor nuclear reactors [
82], and to use very long damage tracks in ancient quartz samples to search for ultra-heavy DM [
83].
Using natural minerals found on Earth to record DM–nucleus interactions over
Gyr timescales provides paleo-detectors with two crucial advantages over conventional direct detection experiments. First, the enormous integration times promise exposures many orders of magnitude larger than what is feasible in laboratory experiments. Reading out 100 g of material that has been recording nuclear damage tracks for 1 Gyr yields an exposure of
g Gyr, or, in units more appropriate for a real-time experiment,
kilotonne yr. The
nm spatial resolution achievable with SAXs readout corresponds to a nuclear recoil energy threshold of
keV. Thus, paleo-detectors would combine the recoil energy threshold achievable in direct detection experiments utilizing liquid noble gases as target materials with the exposures achievable in large (planned) neutrino detectors such as Super-/Hyper-Kamiokande or DUNE. Second, while conventional direct detection experiments measure events in real time, paleo-detectors integrate over timescales up to
Gyr, comparable to the age of the Solar System,
∼
Gyr, and the orbital period of the Solar System around the Milky Way,
∼
Myr. Using a series of samples with different ages, paleo-detectors offer a unique pathway to gaining information about a signal’s time-evolution over Gyr timescales [
76,
77,
78,
84].
Of course, paleo-detectors must pay a price for these advantages: compared to a conventional detector operated in a controlled laboratory environment, background mitigation might be much more challenging. The background sources for DM searches in paleo-detectors are qualitatively the same as in conventional direct detection experiments since both approaches are, at their heart, nuclear recoil searches. Extensive background studies for paleo-detectors have been performed in previous work [
35,
36,
76], and we will discuss the different sources further below. As we will see, the most important background sources are cosmic rays, solar neutrinos, and radioactivity. Regarding the first, the fundamental strategy for paleo-detectors is to use natural minerals extracted from deep within the Earth where the rock-overburden has shielded them from cosmogenic backgrounds. Conventional direct detection experiments are operated in deep underground laboratories for the same reason. However, paleo-detectors require only (at most) a few kg of target material compared to actively operating a tonne-scale experiment; samples of such modest sizes can be extracted from much greater depths than those of existing underground laboratories (for example, from existing boreholes drilled for geological R&D and oil exploration). For the backgrounds from solar neutrinos and radioactivity, it is important to note that, due to their large exposures, paleo-detectors would contain a sizable number of background events. A rare event search, such as for DM-induced nuclear recoils, would thus not proceed as in conventional direct detection experiments, where one typically attempts to find a background-free signal region. Instead, a paleo-detector analysis is more similar to a collider physics search: one would model the distribution of background events (for example in track-length space) and search for deviations of the observed data from the background-only prediction using a signal template. One important consequence of the large number of background (and possible signal) events is that one can characterize many aspects of the backgrounds from the data themselves. Just as in collider searches, paleo-detectors can use
control regions, where the backgrounds dominate, to better characterize the background contribution in the
signal region of the search.
The remainder of this article is organized as follows: in
Section 2, we discuss the calculation of DM signals in paleo-detectors followed by the most important background sources in
Section 3. In
Section 4, we show projections for the reach of paleo-detectors in the conventional WIMP mass–cross section plane. Compared to previous estimates of the sensitivity of paleo-detectors to DM signals [
35,
36,
37], in this work, we make use of both improved background modeling and new statistical techniques. In this sense, the results presented here supersede the sensitivity projections shown in References [
35,
36,
37]. In
Section 4, we also explore the robustness of the sensitivity projections to changes in the assumptions we make in our analysis; for example, we remove external constraints on the age of the minerals and the parameters controlling the background normalizations and explore the effect of systematic mismodeling of the spectral shape of the backgrounds. Furthermore, we show how the projected sensitivity changes with the radiopurity of the samples, and demonstrate that, even if the uranium concentration in paleo-detector samples is
(per weight, denoted as “g/g” throughout), many orders of magnitude larger than what we expect in the most radiopure samples obtained from ultra basic rock or marine evaporite deposits, paleo-detectors could still probe DM–nucleon cross sections below current limits. Specifically, for WIMP masses
GeV/c
, the sensitivity of paleo-detectors could still reach down all the way to the conventional neutrino floor in a Xe-based direct detection experiment. We summarize and conclude in
Section 5. We make the code used in this work available: paleoSpec [
85] for the computation of the signal and background spectra, and paleoSens [
86] for the sensitivity forecasts.
2. Dark Matter Signals in Paleo-Detectors
The signature of WIMP DM in a paleo-detector is the spectrum of damage tracks left by nuclei in the mineral sample after receiving a “kick” from WIMP-nucleus scattering. In this section, we discuss how this signal is computed. The differential rate of nuclear recoils per unit target mass for a WIMP with mass
scattering off target nuclei
N with recoil energy
is given by [
87,
88]
where
is the atomic mass number of
N,
is the reduced mass of the WIMP-proton system with the proton mass
, and we have assumed isospin-conserving spin-independent WIMP–nucleon interactions with the (zero momentum transfer) WIMP-proton cross section
. Expressions analogous to Equation (
1) for isospin-violating SI interaction and spin-dependent WIMP–nucleon interactions can be found in References [
87,
88,
89,
90,
91], and for more general WIMP–nucleus interactions in References [
92,
93]. Note that while we focus on spin-independent interactions in this work, paleo-detectors could also probe spin-dependent WIMP–nucleus interactions; see Reference [
36].
In Equation (
1), the nuclear form factor
accounts for the finite size of the nucleus; we will use the Helm parameterization [
94,
95,
96] in our numerical calculations. (Note that more refined calculations of the form factors are available, although only for a few isotopes; see e.g., References [
97,
98,
99,
100]). The distribution of WIMPs is described by the local WIMP (mass) density,
, and the
mean inverse speed,
where
is the velocity distribution of WIMPs with
(
) the velocity (speed) of the WIMP relative to the detector. In Equation (
2),
accounts for the scattering kinematics;
is the minimal speed, with respect to a nucleus of mass
, that a WIMP must have to give rise to a recoil with energy
(where
is the reduced mass of the WIMP–nucleus system). In our numerical calculations, we set
GeV/cm
, while for
, we use a Maxwell–Boltzmann velocity distribution in the Galactic rest frame with velocity dispersion
km/s [
101], truncated at the escape speed
km/s [
102], and boosted to the Solar System frame by
km/s [
103] as in the Standard Halo Model (SHM) [
7,
95,
104].
In order to obtain the differential track length spectrum (per unit target mass),
, we use the range of a nucleus
N with energy
,
where
is the stop** power of
N in the target material. The true length of a damage track in a given target may differ from the range if, for example, the nucleus creates a lasting damage track only along some portion of the length it travels through the material, or if the nucleus’ trajectory deviates significantly from a straight line. Our previous studies suggest that such effects should be small (see the discussion in Reference [
36]), making Equation (
3) a good proxy for the track length. Furthermore, the effects of thermal annealing could potentially be significant over geological timescales; fortunately, any associated modifications to the track lengths would be similar for both the signal and background recoils. While detailed experimental studies for each target material will be required to account for these possible corrections, for the purposes of this work we will use Equation (
3) for the track length and
SRIM [
105,
106] to calculate stop** powers. The left panel of
Figure 1 shows the range of the different constituent nuclei in gypsum [Ca(SO
(H
O)] as an example. (Note that we do not show the range of H in
Figure 1 since we do not expect ions with charge
(in particular, H and He) to produce lasting damage tracks in most materials; see the discussion in [
36]). We see that a 1 keV recoil gives rise to a track a few nm long, while a 100 keV recoil produces a track with a length of a few hundred nm.
In a composite target material composed of different isotopes with mass fraction
, the differential rate at which tracks are produced (per unit target mass) is then given by
with
given in Equation (
1). In the right panel of
Figure 1 we show
for a 5 GeV/c
(with cross section
) and a 500 GeV/c
(with
) WIMP in gypsum as an example together with various background spectra, which we will discuss in
Section 3. Note that for the
GeV/c
case, the largest contribution to the signal is at track lengths of a few nm, while heavier WIMPs, such as the 500 GeV/c
case shown in
Figure 1, give rise to longer tracks with typical lengths of order
nm.
The nuclear recoil tracks in a paleo-detector sample can be read out with a variety of readout techniques as briefly discussed in the introduction (see Reference [
36] for a more detailed discussion). Here, we will consider two scenarios:
High-resolution readout scenario: We assume that 10 mg of material can be read out with a track length resolution of
nm. This scenario may be achievable with Helium Ion Beam Microscopy (HIBM) combined with pulsed-laser and fast-ion-beam ablation techniques [
46,
47,
53,
54,
55,
56]. As we will see, this scenario is advantageous for low-mass (
GeV
) WIMP searches.
High-exposure readout scenario: We assume that 100 g of material can be read out with
nm. This scenario could be realized via Small Angle X-ray scattering (SAXs) tomography at a synchrotron facility [
48,
49,
50]; it is better suited for heavier (
GeV
) WIMP searches.
These are the same readout scenarios as what we have used in previous paleo-detector studies [
35,
36,
37,
76,
77]. We note that these scenarios represent challenging applications of the respective techniques; demonstrating that their feasibility and especially their scalability to the proposed sample masses would be an important step towards the realization of a paleo-detector experiment. There are ongoing efforts to demonstrate the readout of few-keV recoils using a variety of microscopy techniques at JAMSTEC [
108] and at Rensselaer Polytechnic Institute [
109].
In order to assess the sensitivity of paleo-detectors, we will consider the binned track length spectrum as the experimental observable. Accounting for the finite resolution of the readout method, the number of tracks in the
i-th bin with reconstructed length
in a sample of mass
M, which has been recording tracks for a time
, is
where
is given by Equation (
4) and we have assumed that
is independent of time. (For time-dependent signals, replace
in Equation (
5)). The
window function,
W, accounts for the finite resolution effects of the readout. In our numerical calculations, we will assume that the probability of measuring a track length
x from a track with true length
is Gaussian-distributed with variance
. The corresponding window function is
Note that in order to avoid artificial sensitivity to tracks which are much shorter than the spatial resolution of the readout, we remove all tracks with true track length
from the spectra (
) before applying Equation (
5). For all numerical results in this work, we use 100 log-spaced bins between
and
nm.
Figure 2 shows the same 5 GeV/c
and 500 GeV/c
WIMP DM spectra as the right panel of
Figure 1 together with the various background sources after taking into account finite resolution effects. The left panel is for the high-resolution scenario, while the right panel is for the high-exposure scenario. We will discuss this figure further at the end of
Section 3 after introducing the relevant background sources. Let us here only note that the high-resolution scenario is better suited for relatively light WIMPs with masses
GeV/c
, while the high-exposure scenario is aimed at heavier WIMPs,
GeV/c
. This is for two reasons. First, lighter WIMPs give rise to softer nuclear recoils, resulting in (on average) shorter signal tracks, making excellent track length resolution crucial to resolving the details of the signal and background spectra in the region where the ratio of signal to background events is largest for such light WIMPs. Second, current bounds from direct detection experiments rule out much smaller cross sections for heavier WIMPs with
GeV/c
; hence, searches for WIMPs with
GeV/c
and cross sections smaller than current upper limits require very large exposures.
3. Backgrounds
The background sources for DM searches with paleo-detectors are the same as for conventional direct detection experiments. However, the relative importance of the respective sources is different because, compared to conventional direct detection experiments, paleo-detectors use much smaller target masses (≲1 kg) and integrate over much longer timescales (
0.1–1 Gyr) without information about the timing of individual signal or background events. The experimental observables in paleo-detectors are nuclear damage tracks with, for all practical purposes, perfect rejection of electronic recoils. Furthermore, natural defects in minerals are either single-site or span across the entire (mono)-crystalline volume and do not resemble nuclear damage tracks. Hence, the only relevant backgrounds are nuclear recoils. In the remainder of this section, we will give a brief overview of the most important background sources; we refer the reader to References [
36,
76] for more detailed discussions.
Cosmogenics: In order to mitigate backgrounds from cosmic rays, minerals to be used as paleo-detectors must have been shielded by a large overburden for as long as they have been recording nuclear damage tracks (just as conventional direct detection experiments are operated in deep underground laboratories). Paleo-detectors require only (at most) a few kg of target material rather than actively operating a tonne-scale experiment; such modest amounts of materials can be sourced from even greater depths than those of existing underground laboratories. One promising source of samples is (existing) boreholes drilled for geological R&D or oil exploration. For example, for an overburden of 5 km rock, the cosmogenic-muon-induced neutron flux is
[
110], leading to negligible backgrounds for the purposes of a paleo-detector. We note that at depths
km rock, backgrounds from atmospheric neutrinos producing neutrons in the vicinity of the target become comparable to the backgrounds from cosmogenic muons [
111]. We also note that paleo-detector samples can be stored near the surface after extraction from deep in the Earth before readout. For example, the cosmogenic-muon-induced neutron flux in a ∼50 m deep storage facility is ≲
.
Astrophysical Neutrinos: Neutrinos scattering off nuclei in a paleo-detector sample give rise to nuclear damage tracks. We include neutrinos from the Sun, supernovae, and the interactions of cosmic rays with the Earth’s atmosphere in our background modeling. We take the solar and atmospheric neutrino fluxes from Reference [
112]. Because of the long integration times, paleo-detectors are sensitive not only to neutrinos from supernovae in far-away galaxies throughout our Universe (the Diffuse Supernova Neutrino Background, DSNB), but also to neutrinos from local supernovae—the supernova rate in the Milky Way is estimated to be 2–3 per 100 years [
113,
114,
115,
116,
117,
118]. We compute the DSNB spectrum and the contribution from Galactic supernovae as in Reference [
76]; see also References [
118,
119,
120,
121,
122]. The neutrino-induced background is dominated by solar neutrinos at track lengths below ∼100 nm, with supernova neutrinos giving the largest contributions for tracks of a few 100 nm, before atmospheric neutrinos dominate the neutrino-induced backgrounds at even longer track lengths. Note that while neutrino-induced nuclear recoils are a background for DM searches, they can also be an interesting signal for paleo-detectors; see References [
76,
77,
78].
Radiogenics: Any natural mineral used as a paleo-detector will contain trace amounts of radioactive materials. In order to mitigate radiogenic backgrounds, it is crucial to use minerals with as low a concentration of radioactive elements as possible. The most important radioactive isotope for paleo-detectors is
. Typical minerals formed in the Earth’s crust have
concentrations of order
∼
g/g, leading to prohibitively large backgrounds for DM searches. However, minerals that constitute so-called Ultra Basic Rocks (UBRs), formed from the material of the Earth’s mantle, and Marine Evaporites (MEs), formed from evaporated sea water, have typical
concentrations orders of magnitude lower [
123,
124,
125,
126,
127,
128,
129], (see also Reference [
36] and especially the discussion in the appendix of Reference [
76]). As in previous works on paleo-detectors [
35,
36,
37,
76,
77], we will assume benchmark
concentrations of
g/g for UBRs and
g/g for MEs.
Since paleo-detectors record only nuclear recoils, the most relevant radiogenic background arises from
-decays and spontaneous fission. Note that it is not the
-particles themselves or the fission fragments that lead to backgrounds for DM searches:
-particles (i.e., He nuclei) do not give rise to lasting damage tracks in most materials, and even if they do, the range of the
’s is of order
μm in natural minerals, orders of magnitude larger than the tracks from DM-induced nuclear recoils. Similarly, the fission fragments from spontaneous fission of heavy elements have typical energies of tens of MeV, producing much longer tracks than DM-induced recoils. However, in an
-decay, the decaying nucleus recoils with typical energies of some tens of keV against the
-particle. Due to the long integration time of paleo-detectors, almost all of the
nuclei in the minerals that undergo the first (
) decay in the
decay chain will have continued on to decay to the stable
Pb. There are eight
-decays in the
decay chain, leading to a characteristic pattern of connected nuclear recoil tracks. We expect such patterns of connected tracks to be easily distinguishable from the isolated recoil tracks DM–nucleus interactions would cause and, hence, do not consider them a relevant background for DM searches. However, the second
-decay in the
decay chain (
) has a relatively long half-life,
Myr. This leads to a population of isolated tracks from
atoms that have undergone the initial (
) decay, but not the second
-decay (
) in the uranium series [
130,
131]. The
nucleus receives
keV of recoil energy in the (
) decay, leading to a population of monochromatic isolated tracks with approximate number density
[
36].
The second important source of radiogenic backgrounds is neutrons: spontaneous fission of heavy isotopes as well as (
)-reactions between the
-particles produced in
-decays and nuclei comprising the target material give rise to neutrons with MeV energies. These neutrons in turn scatter off the nuclei in the paleo-detector sample, giving rise to a broad spectrum of nuclear recoil tracks. In order to model this background, we use
SOURCES-4A [
132] to compute the neutron spectra from spontaneous fission and
-reactions of
’s from all
-decaying nuclei in the
decay chain; note that the contributions from spontaneous fission and (
)-reactions to the neutron flux are typically of the same order of magnitude, and the precise balance depends on the chemical composition of the target. From the neutron spectra, we calculate the induced nuclear recoil spectra using
TENDL-2017 [
133,
134,
135,
136] neutron-nucleus cross sections tabulated in the
JANIS4.0 [
137] database. (We take only elastic neutron-nucleus scattering into account; this yields a conservative estimate of the background because neutrons typically lose a larger fraction of their energies through inelastic processes than in elastic scattering. Note also that our Monte Carlo simulation of the nuclear recoils induced by radiogenic neutrons has been validated with a calculation of the nuclear recoils induced by the same neutron spectra with
FLUKA [
138,
139,
140] in Reference [
77] for the particular case of a halite target.) The neutron-induced backgrounds are strongly suppressed in target materials containing hydrogen: since neutrons and H nuclei have approximately the same mass, neutrons lose a large fraction of their energy in a single interaction with hydrogen, leading to efficient moderation of the radiogenic neutrons even if hydrogen comprises only a small fraction of the atoms in the target material [
35,
36].
In
Figure 1 and
Figure 2, we show the track-length spectra produced by these various background sources together with two benchmark DM signals in gypsum [Ca(SO
(H
O)]. At short track lengths, the background budget is dominated by (solar) neutrino-induced nuclear recoils. At intermediate ranges between a few tens and a hundred nm (the precise values depend on the particular mineral), the isolated
tracks from
decays contribute, before the background budget becomes dominated by radiogenic-neutron-induced backgrounds at larger track lengths. Comparing to the benchmark WIMP spectra, we see that solar neutrinos are the most relevant background source for WIMPs with masses in the few GeV/c
range, while radiogenic neutrons are the limiting background in searches for heavier WIMPs with masses larger than ∼
GeV/c
.
In
Figure 2, we show the (binned) number of tracks one would observe in the high-resolution (left panel) and high-exposure (right panel) scenarios after accounting for the effects of finite readout resolution. Because of the large exposure of paleo-detectors, a sizable number of background events will be present. A DM search with paleo-detectors would thus proceed quite differently from a traditional direct detection experiment, in which one typically attempts to construct a signal region with very few (or, ideally, zero) background events, and then searches for an
number of DM-induced events. Given the large number of background events, in a paleo-detector analysis one would model the distribution of background events (for example, in track-length space) and search for deviations of the observed data from the background-only prediction. Thus, the
uncertainty in the background prediction is much more important than the absolute number of background events.
For illustration, we include error bands around the respective background components in
Figure 2. These bands include both the (statistical) Poisson noise and a
systematic error in the background prediction. For the latter, we make well-motivated assumptions in order to illustrate how the sensitivity of paleo-detectors to signals from DM of various masses would depend on different background components in a cut-and-count analysis; see previous studies of paleo-detectors [
35,
36,
37,
76,
77]. However, we note that in the spectral analysis discussed in
Section 4, the sensitivity of paleo-detectors to a WIMP DM signal is virtually independent of any prior knowledge of the background normalizations. Let us briefly describe our choices for the systematic errors in
Figure 2. Neutrino-induced backgrounds are controlled by sources that are not associated with the environment of the paleo-detector target sample: the fluxes of solar, supernova, or atmospheric neutrinos. These fluxes may vary by an
factor over the ∼Gyr integration time of paleo-detectors [
76,
77,
78]; hence, we assign a (relative) systematic uncertainty of
to the neutrino-induced backgrounds. (In
Figure 2, we plot all neutrino-induced backgrounds together to avoid excessive clutter. In our sensitivity projections, we treat the solar, DSNB, Galactic supernova, and atmospheric-neutrino-induced backgrounds as four independent background contributions.) The radiogenic backgrounds, on the other hand, are controlled by the
concentrations, which do not change over time. Furthermore, the shape of the track length spectra induced by radiogenic backgrounds can be calibrated by studying the track length distributions in mineral samples with large concentrations of
. Finally, the
concentration in a given sample can be measured using, for example, mass spectrometry or gamma ray spectroscopy [
141,
142]. Hence, we expect the uncertainty on the radiogenic background predictions to be much smaller, and we assign them a (relative) systematic error of
.
From the bottom panels in
Figure 2, we can directly read off the
signal-to-noise ratio (per bin), i.e., the ratio of the number of signal events to the uncertainty of the background prediction: in these panels, the lines show the ratio of the number of signal to the (total) number of background events (
S/
B), while the sand-colored shaded region shows the relative uncertainty of the background prediction. For a given signal, the ratio of
to the (relative) background uncertainty directly corresponds to the significance one would obtain in a simple cut-and-count analysis. As discussed in References [
35,
36], any signal for which
S/
B is larger than the shaded region, indicating the background uncertainty in the bottom panels of
Figure 2 in at least one bin, could be probed with a simple cut-and-count analysis.
Given the large number of events expected in a paleo-detector, one can go beyond the simple cut-and-count analysis and perform a
spectral analysis in track-length space: Armed with predictions for the spectral shape of the various background components and of a possible signal, one can fit all of these spectral templates to the data and ask if the data is better described by either the backgrounds-only or the backgrounds+signal hypothesis. Heuristically, such a search uses
control regions, where some background dominates, to measure that background component, and uses this information to reduce the error of the background prediction in the
signal region, where the signal-to-background ratio is largest. For example, short track lengths can serve as a control region for the solar-neutrino-induced backgrounds (as well as the
tracks in the high-exposure scenario), while track lengths longer than ∼
nm can provide an excellent control region for the neutron-induced backgrounds. We have first used such a spectral analysis in Reference [
37] and will discuss an alternative approach for the spectral analysis in the next section.
4. Sensitivity Forecasts
In
Section 2 we described how to compute the DM signal in paleo-detectors. In particular, we showed how to obtain the track length spectrum one would observe in a paleo-detector from the differential recoil rate per recoil energy and how to account for finite track-length resolution effects. In
Section 3, we discussed the most relevant background sources and briefly sketched how one would perform a DM search with a paleo-detector. Due to the potentially large exposure, a sizable number of background (and, possibly, signal) events will be present in any paleo-detector sample. In order to search for DM, one would perform a spectral analysis of the observed recoil spectrum in track-length space and ask if that spectrum is better fit by a background-only model or by a background+signal hypothesis. In References [
35,
36] we used a much simpler cut-and-count analysis to project the sensitivity of paleo-detectors to DM, while in Reference [
37], we estimated the sensitivity using a spectral analysis based on the Fisher-matrix formalism using the
swordfish package [
143,
144,
145]. In the remainder of this section, we describe a profile likelihood ratio-based approach to calculate the projected sensitivity of paleo-detectors to a DM signal with a spectral analysis. As we will see, one advantage of using this statistical approach is that we can easily incorporate the effect of background mismodeling into our sensitivity forecast. We show results for the sensitivity forecasts in
Figure 3,
Figure 4 and
Figure 5.
The predicted track length spectrum in a paleo-detector is
where the
s are the track length spectra described in Equation (
5). Each
is computed from the appropriate recoil energy distribution, given by Equation (
1), for the DM signal and the corresponding expressions for the different background components discussed in
Section 3. In Equation (
7), the terms in the first line describe the contributions induced by solar neutrinos (
), the DSNB (
), the Galactic Supernova Neutrino Background (
), and atmospheric neutrinos (
). The terms in the second line describe the radiogenic backgrounds, separated into: the isolated
decays from
nuclei, which have undergone only the inital
-decay during the time the sample has been recording nuclear damage tracks (
), and the background induced by radiogenic neutrons (
). Finally, the term in the third line is the DM signal. For each
, we indicate in parantheses the parameter controlling its normalization: the flux (
) for each neutrino contribution, the
concentration in the sample (
) for the radiogenic backgrounds, and the WIMP-nucleon cross section (
) for the DM signal. Note that all contributions are also proportional to the mass of the sample,
M, and all contributions except for
are proportional to the length of time the sample has been recording damage tracks,
; we have suppressed these dependencies in Equation (
7) for compactness. Finally, the DM signal also depends on the mass of the WIMP,
.
To easily compare to existing upper limits and the projected sensitivities of conventional direct detection experiments, we calculate the projected upper (exclusion) limit on
at 90% confidence level as a function of
for a paleo-detector. The set of nuisance parameters in Equation (
7) is then
Since paleo-detectors are ultimately a counting experiment, we expect the data to be Poisson-distributed. The appropriate log-likelihood to observe the data
given the parameters
is thus (here and in the following, we drop constant factors in the expression of the likelihood, which cancel in the likelihood ratio we are ultimately interested in):
The nuisance parameters will be constrained by other measurements. For example, the neutrino fluxes will be constrained by existing neutrino experiments and theoretical modeling. The mass, age, and uranium concentration of the sample will be constrained by direct measurements of the sample mineral. In a frequentist analysis, these constraints can be taken into account by jointly analyzing these ancillary measurements with the track length spectrum data. To mimick such a joint analysis, we include a set of external (Gaussian) constraints on the nuisance parameters:
Here,
is the central value for the
j-th nuisance parameter supplied by the ancillary measurement(s), and
is the
relative uncertainty of the measurement(s). The likelihood is then
To calculate the projected sensitivity of paleo-detectors, we consider the maximum log-likelihood ratio [
146,
147,
148]:
In the numerator,
is the set of nuisance parameters that maximizes
for given values of
(and
), while in the denominator,
is maximized over both
and
(for a given value of
), with best-fit values
and
. In order to calculate the projected exclusion limit, we use the
Asimov data set [
146] under the assumption of no signal:
where
is a set of fiducial choices for the nuisance parameters. The 90% confidence level exclusion limit for a given value of
is obtained by finding the smallest value of
for which
. (By Wilks’ theorem [
149], the distribution of
q is (asymptotically)
distributed, and for a one-dimensional
-distribution,
for a
p-value of 0.1.)
We describe the the fiducial values for the nuisance parameters listed in Equation (
8) and the relative uncertainty of the external constraints on those parameters,
, defined by Equation (
10). For the high-resolution (high-exposure) scenario, we choose a fiducial value of
mg (
g) for the sample mass, and for both scenarios,
Gyr for the time the sample has been recording damage tracks. We assume that these parameters can be constrained by ancillary measurements to
and
; while measuring the mass of a target sample to high precision is trivial, its age can be inferred with few-percent precision using standard geological dating techniques [
42,
43,
44]. For the normalization of the neutrino fluxes,
, and the
concentrations in our target minerals,
, we take fiducial values described in
Section 3 (the latter are also summarized in
Table 1). As a benchmark, we assume that the associated nuisance parameters can be constrained to
and
; these choices are equivalent to the respective systematic uncertainties of the backgrounds assumed in previous studies [
35,
36,
37,
76,
77]. However, as we will discuss below, the projected sensitivity of paleo-detectors to a DM signal is largely independent of external constraints on the nuisance parameters.
In
Figure 3, we show the projected 90% confidence level upper limit under these fiducial assumptions in both the high-resolution (sample mass
mg, track-length resolution
nm; left panel) and the high-exposure (
g,
nm; right panel) readout scenario. We show results for a selection of minerals (see
Table 1): three examples of Ultra-Basic Rocks (UBRs), the very common olivine and two less common UBRs, which contain hydrogen, phlogopite mica, and nchwaningite. We also show the projected sensitivity for three examples of Marine Evaporites (MEs): halite, gypsum, and sinjarite, where the latter two contain hydrogen. As discussed above, the high-resolution scenario has better sensitivity for WIMP masses
GeV/c
where the dominant background stems from solar neutrinos; here, excellent track length resolution is crucial to resolving the differences between the shapes of the signal and background track length spectra. For
GeV/c
, the high-exposure scenario has better sensitivity, since for such WIMP masses, radiogenic backgrounds are dominant and maximizing the exposure is crucial in order to achieve the best sensitivity. Comparing the sensitivity for the different target materials, we can note that materials containing hydrogen generally have better sensitivity than those without hydrogen. This is due to the effective suppression of the radiogenic neutron backgrounds by hydrogen; see the discussion in
Section 3. Furthermore, we note that MEs typically have better sensitivity than UBRs; this is because we assume a lower
concentration (
g/g) for MEs than for UBRs (
g/g). Out of the minerals considered in
Figure 3, throughout the shown WIMP mass range and for both the high-resolution and the high-exposure readout scenarios, sinjarite promises the best sensitivity. Primarily, this is due to sinjarite being an ME containing hydrogen, and hence having relatively small
-induced backgrounds. However, the sensitivity is also enhanced by its particular chemical composition: the other example of an ME containing hydrogen shown in
Figure 3, gypsum, is very oxygen-rich, and hence, both the DM-induced and the background spectra are dominated by the contribution of oxygen. Sinjarite, on the other hand, has one calcium, two oxygen, and two chlorine atoms; each of these elements can therefore contribute comparably to the recoil spectra. This adds additional structure to the track length spectra, which is advantageous in a spectral analysis.
In
Figure 3, we also show existing limits from conventional direct detection experiments [
12,
14,
17,
107,
150] as well as the neutrino floor for a Xe-based experiment [
151] with the shaded regions. At small WIMP masses,
GeV/c
, paleo-detectors could probe WIMPs with WIMP-nucleon cross sections many orders of magnitude smaller than those currently probed by conventional experiments; for example, at
GeV/c
the best existing upper limit is
[
14,
17]. Moreover, paleo-detectors offer sensitivity deep into the conventional (solar) neutrino floor in the
GeV/c
mass range. Due to the large exposures, a paleo-detector would observe
(100–1000) recoil tracks from solar neutrinos in the high-resolution readout scenario; see
Figure 2. In a spectral analysis, this large number of neutrino-induced events combined with the excellent track-length (and, in turn, recoil energy) resolution would allow one to characterize the solar neutrino background in the “signal region” and hence probe WIMP-nucleon cross sections deep into the conventional neutrino floor. For heavier WIMPs (
GeV/c
), we see that paleo-detectors could probe WIMP-nucleon cross sections roughly one and a half orders of magnitude below current upper limits, but somewhat larger than the conventional neutrino floor in a Xe-based direct detection experiment. In this WIMP mass region, the sensitivity of paleo-detectors is limited by the radiogenic neutron backgrounds, as we will discuss further below.
In the remainder of this section, we will explore the robustness of the sensitivity projections shown in
Figure 3 to the assumptions we have made. Let us first note that the projected exclusion limits are virtually independent of the ancillary measurements. For fixed normalization of the backgrounds in the Asimov data, removing the Gaussian constraints (defined by Equation (
10)) on the parameters controlling the normalization of the neutrino-induced backgrounds (
) and the radiogenic backgrounds (
). as well as removing the external constraint on the age of the samples (
) leads to a depreciation of the projected upper limits by less than a factor two across the WIMP mass range shown in
Figure 3. We find the most pronounced depreciation of the upper limits in the
GeV/c
WIMP range when the constraints on the
are removed. Hence, the data from the track-length spectrum alone suffices to search for DM with paleo-detectors, even in the absence of any ancillary information about the normalization of the backgrounds.
In order to explore how robust the sensitivity forecasts are with respect to mismodeling of the
shape of each background component, we replace the Poisson log-likelihood, Equation (
9), with a
distribution,
where
is the variance in the
i-th bin. If we include only the Poisson error of the data,
, the Poisson log-likelihood and the
likelihood become asymptotically equivalent in the large-
limit. Given the number of signal events one would observe in a paleo-detector (see
Figure 2) we thus expect to find the same sensitivity for both forms of the likelihood. The
distribution allows us to include additional uncertainties, which we model as
where
parameterizes the (relative) systematic modeling error. We note that it is straightforward to extend this to
, allowing for the assignment of a systematic error to the
j-th contribution to the data in the
i-th bin. In this work, we will only employ a simple overall systematic error
to explore the robustness of the sensitivity projections.
In
Figure 4, we show the projected upper limit for a gypsum paleo-detector under different assumptions for the background-shape mismodeling. First, we note that for
, i.e., accounting only for the Poisson error in the variance,
, both the Poisson log-likelihood and the
distribution yield the same projected limit. In the high-resolution scenario (left panel of
Figure 4) we find that if we include an error of
, corresponding to 10% bin-to-bin mismodeling of the shape of the backgrounds, the sensitivity does not depreciate significantly. For
, the projected upper limit is up to an order of magnitude worse than in the
case for
GeV/c
, with smaller changes at larger WIMP masses. In the high-exposure scenario (right panel of
Figure 4), the sensitivity is not quite as robust against background-shape mismodeling. For small systematic errors,
, we find that the sensitivity becomes only slightly worse than for
, while for
, the projected upper limit is a factor of ∼5 larger than for
. The high-exposure-readout scenario is most relevant for the mass range
GeV/c
, where the dominant background stems from radiogenic neutrons. The corresponding track length spectrum is straightforward to calibrate in paleo-detectors by measuring the tracks in samples of the same mineral but with larger
concentration. Thus, it seems unlikely that the systematic (bin-to-bin) shape errors would be as large as
in an actual paleo-detector experiment. However, we note from
Figure 4 that even for errors as large as
, a gypsum paleo-detector could still probe WIMP-nucleon cross sections almost an order of magnitude below current experimental limits for
GeV/c
.
In
Figure 5, we show the effect of assuming different
concentrations on the projected upper limits in a gypsum paleo-detector. Trivially, we find that paleo-detectors could probe smaller WIMP cross sections if the samples are more radiopure. For most of the range of
concentrations,
g/g shown in
Figure 5 and most of the shown WIMP mass range, the projected upper limit scales like
. This scaling arises from the uncertainty of the radiogenic backgrounds in the “signal region” being given by the associated Poisson error, with the normalization of the radiogenic backgrounds inferred from the relevant “control regions”. The
scaling of the sensitivity is broken if other backgrounds provide significant contributions in the “signal region.” For example, the radiogenic background contribution can be subdominant to that of solar neutrinos for
GeV/c
or, for very small
concentrations, atmospheric neutrinos for
GeV/c
. Note that in the high-resolution scenario, paleo-detectors could probe WIMP-nucleon cross sections down to the conventional (solar) neutrino floor even for
concentrations as large as
g/g in the
GeV/c
mass range. For heavier WIMPs, we find that a gypsum paleo-detector in the high-exposure readout scenario could still probe WIMP-nucleus cross sections below current experimental limits even if the
concentration is orders of magnitude larger than our fiducial assumption of
g/g.
Before moving on, let us note that the statistical framework used here to project the sensitivity of paleo-detectors can easily be adapted to answer different statistical questions or to incorporate additional uncertainties. For example, by using the Asimov data as defined above (backgrounds only) but also including a DM signal contribution, one can calculate the projected
discovery reach rather than the exclusion limit. Using the
likelihood, it is straightforward to include arbitrarily complicated uncertainties; another approach to including (modeling) uncertainties would be to use random realizations of noisy mock data instead of the Asimov data set. It is also straightforward to extend the framework to simultaneously analyze a series of paleo-detectors. Samples of different ages can be used to explore the sensitivity of paleo-detectors to time-dependent signals. For instance, the Solar System could have traversed some denser region of the DM halo during its past few revolutions around our Galaxy [
84]. Time-dependent neutrino signals may also arise from the evolution of our Sun [
78], the Galactic supernova rate [
76], or the flux of cosmic rays hitting Earth’s atmosphere [
77]. The code used in this work is available here [
86]. Some of the options mentioned in this paragraph are already implemented in the code, and we invite the reader to add the capabilities in order to answer their own analysis questions.