Geoneutrinos in Large Direct Detection Experiments

TL;DR

This paper explores the potential of large direct detection experiments, especially using germanium, to detect geoneutrinos, providing new insights into Earth's interior and formation history.

Contribution

It demonstrates the feasibility of detecting geoneutrinos in large direct detection experiments, highlighting germanium as the most promising target material.

Findings

Germanium offers the lowest recoil energy threshold for detection.

Detection of $^{40}$K geoneutrinos is possible with realistic thresholds.

Minimum required exposure is around 10 tonne-years.

Abstract

Geoneutrinos can provide a unique insight into Earth's interior, its central engine and its formation history. We study the detection of geoneutrinos in large direct detection experiments, which has been considered non-feasible. We compute the geoneutrino-induced electron and nuclear recoil spectra in different materials, under several optimistic assumptions. We identify germanium as the most promising target element due to the low nuclear recoil energy threshold that could be achieved. The minimum exposure required for detection would be $\mathcal{O}(10)$ tonne-years. The realistic low thresholds achievable in germanium and silicon permit the detection of $^{40}$K geoneutrinos. These are particularly important to determine Earth's formation history but they are below the kinematic threshold of inverse beta decay, the detection process used in scintillator-based experiments.

Figures

Tables

Table 1. Table 1 : Considered experimental configurations. The last column shows the maximum energy of our region of interest.

Target	$A\left(Z\right)$	Isotope	Energy	Max
Material		Fraction	Threshold	Energy
			[eVnr]	[eVnr]
xenon(Xe)	128(54)	0.019	100	267
	129(54)	0.264
	130(54)	0.041
	131(54)	0.212
	132(54)	0.269
	134(54)	0.104
	136(54)	0.089
argon(Ar)	40(18)	0.996	100	855
germanium(Ge)	70(32)	0.208	40	489
	72(32)	0.275
	73(32)	0.077
	74(32)	0.363
	76(32)	0.076
silicon(Si)	28(14)	0.922	78	1221
	29(14)	0.047
	30(14)	0.031

Table 2. Table 2 : Geoneutrino signal flux.

Neutrino	Flux	Max Energy
Component	$\left[{\mathrm{cm}}^{-2}{s}^{-1}\right]$	$E_{\nu }\left[\mathrm{MeV}\right]$
$\mathrm{Geo}({\overline{\nu }}_e,{}_{}{}^{40}K)$	$5.47\left(1\pm 0.006\right)\times {10}^7$	$1.32$
$\mathrm{Geo}({\overline{\nu }}_e,{}_{}{}^{238}U)$	$1.22\left(1\pm 0.012\right)\times {10}^7$	$3.99$
$\mathrm{Geo}({\overline{\nu }}_e,{}_{}{}^{232}\mathrm{Th})$	$1.14\left(1\pm 0.300\right)\times {10}^7$	$2.26$

Table 3. Table 3 : Solar and reactor neutrino background fluxes.

Neutrino	Flux	Max Energy
Component	$\left[{\mathrm{cm}}^{-2}{s}^{-1}\right]$	$E_{\nu }\left[\mathrm{MeV}\right]$
$\mathrm{Solar}({\nu }_e,\mathrm{pp})$	$6.03\left(1\pm 0.006\right)\times {10}^{10}$	$0.42$
$\mathrm{Solar}({\nu }_e,\mathrm{pep}\left[\mathrm{line}\right])$	$1.47\left(1\pm 0.012\right)\times {10}^8$	$1.45$
$\mathrm{Solar}({\nu }_e,\mathrm{hep})$	$8.31\left(1\pm 0.300\right)\times {10}^3$	$18.77$
$\mathrm{Solar}({\nu }_e,{}_{}{}^7\mathrm{Be}\left[\mathrm{line}1\right])$	$4.56\left(1\pm 0.070\right)\times {10}^8$	$0.39$
$\mathrm{Solar}({\nu }_e,{}_{}{}^7\mathrm{Be}\left[\mathrm{line}2\right])$	$4.10\left(1\pm 0.070\right)\times {10}^9$	$0.87$
$\mathrm{Solar}({\nu }_e,{}_{}{}^{8}B)$	$4.59\left(1\pm 0.140\right)\times {10}^6$	$16.80$
$\mathrm{Solar}({\nu }_e,{}_{}{}^{13}N)$	$2.17\left(1\pm 0.140\right)\times {10}^8$	$1.20$
$\mathrm{Solar}({\nu }_e,{}_{}{}^{15}O)$	$1.56\left(1\pm 0.150\right)\times {10}^8$	$1.73$
$\mathrm{Solar}({\nu }_e,{}_{}{}^{17}F)$	$3.40\left(1\pm 0.170\right)\times {10}^6$	$1.74$
$\mathrm{Reactor}\left({\overline{\nu }}_e\right)$	$4.95\left(1\pm 0.100\right)\times {10}^4$	$10.00$

Table 4. Table 4 : Fission fraction and average released energy for each nuclear reactor isotope.

Isotope	$f_k$	$E_k$ [MeV/fission]
${}_{}{}^{235}U$	0.58	$202.36\pm 0.26$
${}_{}{}^{238}U$	0.07	$205.99\pm 0.52$
${}_{}{}^{239}\mathrm{Pu}$	0.30	$211.12\pm 0.34$
${}_{}{}^{241}\mathrm{Pu}$	0.05	$214.26\pm 0.33$

Table 5. Table 5 : Fitted values of the reactor neutrino spectrum α i , k subscript 𝛼 𝑖 𝑘 \alpha_{i,k} coefficients, for i = 1 𝑖 1 i=1 to 6, for the dominant nuclear isotopes (denoted by k 𝑘 k ).

$i$	235	238U	239Pu	241Pu
1	3.217	0.4833	6.413	3.251
2	-3.111	0.1927	-7.432	-3.204
3	1.395	-0.1283	3.535	1.428
4	-0.3690	-0.006762	-0.8820	-0.3675
5	0.04445	0.002233	0.1025	0.04254
6	-0.002053	-0.0001536	-0.004550	-0.001896

Table 6. Table 6 : List of selected most relevant nuclear reactors for Jinping. Columns contain, from left to right, the name of the reactor, the GPS location, the distance to the laboratory in kilometers and the thermal output of the reactor in MW.

Reactor	Location	Distance	Output
		[km]	[MWt]
Fangchengang	21∘40′15′′N 108∘33′30′′ E	998	12110
Changjiang	19∘27′39′′ N 108∘53′60′′ E	1206	3806
Yangjiang	21∘42′30′′ N 112∘15′40′′ E	1239	17430
Taishan	21∘55′04′′ N 112∘58′55′′ E	1318	1980
Ling Ao	22∘36′17′′ N 114∘33′05′′ E	1422	11620
Fuqing	25∘26′45′′ N 119∘26′50′′ E	1763	17740
Hanbit	35∘24′54′′ N 126∘25′26′′ E	2463	16874

Equations

$${\varphi }_k(E)=\frac{R_{{\bar{\nu }}_e}}{4\pi d^2}\sum _k{f}_k{S}_k(E),$$

$$R_{\bar{\nu}_{e}}=N_{\nu,{\rm f}}\,\frac{P_{th}}{\sum_{k}f_{k}E_{k}}e\simeq 6\times 10^{20}\Big{(}\dfrac{P_{th}}{\rm 1~{}GW}\Big{)}e~{}{\rm s}^{-1}~{},$$

$$S_k(E_{\nu })=\frac{d{N}_{\nu }}{d{E}_{\nu }}=\exp \left(\sum _{i=1}^6{\alpha }_{i,k}{E}_{\nu }^{i-1}\right),$$

$$\frac{d{R}_{\nu }^I}{d{E}_R}=\frac{C_I}{m_I}{\int }_{E_{\nu }^{\min }}{\varphi }_{\nu }(E_{\nu })\frac{d{\sigma }^I(E_{\nu },E_R)}{d{E}_R}d{E}_{\nu },$$

$$\frac{d{R}_{\nu }}{d{E}_R}=\mathrm{MT}\sum _I\frac{d{R}_{\nu }^I}{d{E}_R},$$

$$\frac{d{R}_{\nu }}{d{E}_R}=\mathrm{MT}\sum _I\frac{C_I}{m_I}{Z}_I{\int }_{E_{\nu }^{\min }}{\varphi }_{\nu }(E_{\nu })\frac{d{\sigma }^e(E_{\nu },E_R)}{d{E}_R}d{E}_{\nu },$$

$$E_{\nu }^{\min }=\sqrt{\frac{m{E}_R}{2}}.$$

$$E_R^{\max }=\frac{2E_{\nu }^2}{m+2E_{\nu }}.$$

$$\frac{d{\sigma }^I(E_{\nu },E_R)}{d{E}_R}=\frac{G_f^2{m}_I}{4\pi }{Q}_w^2\left(1-\frac{m_I{E}_R}{2E_{\nu }^2}\right)F_{\mathrm{SI},I}^2(E_R),$$

$$\frac{d{\sigma }^e(E_{\nu },E_R)}{d{E}_R}=$$

$$\displaystyle+g_{\nu}g_{\nu}^{\prime}\dfrac{m_{e}E_{R}}{E_{\nu}^{2}}\Big{]}~{}.$$

$$g_{\nu (\mu ,\tau )}=2{\sin }^2{\theta }_{\mathrm{W}}-1,g_{\nu (e)}=2{\sin }^2{\theta }_{\mathrm{W}}+1$$

$$g_{\nu }^{\prime }(e,\mu ,\tau )=2{\sin }^2{\theta }_{\mathrm{W}}.$$

$$\langle P_{ee}\rangle=\cos^{4}\theta_{\rm 13}\Big{(}1-\dfrac{\sin^{2}(2\theta_{\rm 12})}{2}\Big{)}+\sin^{4}\theta_{\rm 13}\simeq 0.55~{},$$

$$\mathcal{L}=$$

$$\displaystyle\times\prod_{k^{\prime}=1}^{n_{\nu}-n_{g}}\dfrac{1}{\sqrt{2\pi}\sigma_{\nu_{k^{\prime}}}}{\rm exp}\biggr{[}-\left(\frac{\phi_{\nu_{k^{\prime}}}-\bar{\phi}_{\nu_{k^{\prime}}}}{\sqrt{2}\sigma_{\nu_{k^{\prime}}}}\right)^{2}\biggr{]}~{}.$$

$$\frac{d{R}_{\mathrm{tot}}}{d{E}_R}=\sum _{k_g}^{n_g}\frac{d{R}_{{\nu }_{k_g}}}{d{E}_{\mathrm{R}}}({\varphi }_{{\nu }_{k_g}})+\sum _{k^{\prime }}^{n_{\nu }-n_g}\frac{d{R}_{{\nu }_{k^{\prime }}}}{d{E}_{\mathrm{R}}}({\varphi }_{{\nu }_{k^{\prime }}}),$$

$$P_k=\frac{{\mu }_k^{n_k}{e}^{-{\mu }_k}}{n_k!},$$

$$\Big{(}\text{PDF}\Big{)}_{k}=\dfrac{1}{R_{k}}\dfrac{dR_{k}}{dE_{R}}~{}.$$

$$q_0=\{\begin{array}{ll}-2\text{ln}\left(\frac{\mathcal{L}({\vec{\varphi }}_{{\nu }_{k_g}}=0,{\widehat{\hat{\vec{\varphi }}}}_{{\nu }_{k^{\prime }}})}{\mathcal{L}({\hat{\vec{\varphi }}}_{{\nu }_{k_g}},{\hat{\vec{\varphi }}}_{{\nu }_{k^{\prime }}})}\right),& {\hat{\vec{\varphi }}}_{{\nu }_{k_g}}\ge 0\\ 0,& {\hat{\vec{\varphi }}}_{{\nu }_{k_g}}<0,\end{array}$$

$$p_0={\int }_{q_0^{\mathrm{obs}}}^{\infty }d{q}_{0}f(q_0|H_0),$$

$$F({\varphi }_{{\nu }_{k_g}},\mathrm{MT})\equiv \sum _{i=1}^{N_{\mathrm{sim}}}\frac{1}{N_{\mathrm{sim}}}\{\begin{array}{ll}1& \text{mbox}if{p}_0^i\le \alpha (Z)\\ 0& \text{mbox}if{p}_0^i>\alpha (Z)\end{array},$$

Full text

I Introduction

The nature of dark matter (DM) remains one of the most pressing issues in modern physics. While weakly interacting massive particles (WIMPs) remain a theoretically well-motivated DM candidate, despite significant efforts no convincing detection signatures have been observed and a sizable portion of the allowed parameter space has been constrained. Planned next generation large scale direct detection experiments will further explore the uncharted corners of WIMP interaction type, candidate mass and cross section. The experiments will eventually encounter an irreducible background due to coherent neutrino-nuclear interactions of solar as well as other neutrinos, i.e. the “neutrino floor” (e.g. Billard et al. (2014); Monroe and Fisher (2007); Gelmini et al. (2018)). The neutrino floor can obfuscate the potential DM signal. Hence, this could constitute a major obstacle in making definitive statements regarding DM observations. However, the neutrino signal, which experiments will surely observe, as well as possible other sources, constitute in themselves interesting targets of exploration. It is thus imperative and timely to fully explore the scientific scope of large direct detection experiments beyond their main scope as DM detectors.

The next generation of large scale direct detection experiments constitute a promising venue as neutrino detectors. This has been bolstered by the recent detection of the coherent neutrino-nucleus scattering by the COHERENT experiment Akimov et al. (2017). Both coherent neutrino-nucleus as well as elastic electron-neutrino scatterings in direct detection experiments have already been employed in a range of studies to assess future detection capabilities within the context of neutrino-related physics, including sterile neutrinos (e.g. Pospelov (2011); Billard et al. (2015)), non-standard neutrino interactions (e.g. Harnik et al. (2012); Dutta et al. (2017)), supernovae explosions Chakraborty et al. (2014), solar neutrinos (e.g. Billard et al. (2015); Cerdeno et al. (2018); Essig et al. (2018); Newstead et al. (2019)) as well as DM annihilations McKeen and Raj (2018).

Geoneutrinos (see e.g. Bellini et al. (2013) for review) are anti-neutrinos $\overline{\nu}_{e}$ that originate from radioactive decays of uranium 238U, thorium 232Th and potassium 40K within the Earth’s interior. They provide a unique insight into the Earth’s composition, heat production mechanisms as well as its thermal evolution and formation. In particular, geoneutrinos can reveal crucial information about the Earth’s energy budget, a fundamental quantity in geology. They allow one to establish which fraction of energy originates from ongoing radiogenic decays of heat producing elements and which is primordial, associated with Earth’s formation history. Furthermore, geoneutrinos can test the hypothetical Earth’s central reactor, proposed to source the Earth’s magnetic field Hollenbach and Herndon (2001).

The study of geoneutrinos finally became a possibility with recent developments of neutrino detectors. Geoneutrino observations have been reported by the KamLAND Araki et al. (2005) and Borexino experiments Bellini et al. (2010), through inverse beta decay (IBD) $\overline{\nu}_{e}+p\rightarrow e^{+}+n$ interactions within scintillator materials. Large uncertainties, originating from low statistics and high backgrounds, do not allow experiments to make definitive statements regarding geophysics at the current stage. Further progress is expected with future scintillator-based experiments such as SNO+ Andringa et al. (2016), JUNO An et al. (2016); Han et al. (2016), HANOHANO Cicenas and Solomey (2012) and Jinping Beacom et al. (2017); Wan et al. (2017); Sramek et al. (2016). Directional detection experiments that utilize electron-neutrino scattering have also been put forward Leyton et al. (2017); Wang and Chen (2017).

Following earlier work in the field (e.g. Monroe and Fisher (2007)), a potential geoneutrino signal in direct detection experiments has been largely neglected, assumed to be too small for detection. Given the status of experimental searches, the plans for larger detectors as well as recent interest in geoneutrinos, it is timely to revisit this issue and reassess the geoneutrino observation potential of future large-scale direct detection experiments. This is the main purpose of this work.

II Large Direct Detection Experiments

A variety of proposals for the next generation multi-ton scale direct detection experiments have been put forth Cushman et al. (2013). Due to the inherent uncertainty in the final designs, making definitive statements about their scientific capabilities requires making assumptions about their composition, energy threshold and backgrounds. We explore several of the most likely detector configurations, based on xenon(Xe), germanium(Ge), argon(Ar) and silicon(Si), assuming optimistic design realizations.

In Table 1 we list the experimental setups we consider, each with a different target element. We do not include isotopic abundances of elements that contribute less than $1\%$. We optimistically assume that experiments have perfect detection efficiency and energy resolution. Furthermore, when considering neutrino coherent interactions with the nuclei the expected background is assumed to arise exclusively from neutrinos. This allows to treat all the analyses on the same footing and provide results independent of specific configurations that could change in the future. The same assumption is not realistic for neutrino-electron scattering, since the backgrounds there are expected to be significant. Thus, we do not analyze this channel.

The reactor neutrino background as well as the geoneutrino signal depend sensitively on the location of the experiment. We assume that the experiment is located in the Jinping Underground Laboratory in China, where the expected geoneutrino signal is the largest and the reactor background is small. Its depth (6720 m.w.e.) ensures that the background due to cosmogenic muons will be highly suppressed.

Recent analyses of xenon-based detectors, such as XENON1T Aprile et al. (2017, 2018) and LUX Akerib et al. (2017), set nuclear recoil thresholds of a few keV, when employing both scintillation and ionization signals. These experiments can go to even lower recoil energies when using only their ionization signal, although their background rejection in this operation mode is less efficient. In this way the XENON100 collaboration lowered its threshold to 0.7 keV for nuclear recoils Aprile et al. (2016). Thresholds as low as 0.1 keV have been used in some projections for future detectors, such as LZ Akerib et al. (2018). Similarly, the argon-based DarkSide-50 was able to lower its threshold for nuclear recoils to 0.6 keV with an analysis of only the ionization signal Agnes et al. (2018a) (although their official analysis has some-what higher thresholds Agnes et al. (2018b)). Assuming optimistic realizations of future noble-gas based experiments, we use a threshold of 0.1 keV for xenon as well as argon.

Germanium and silicon will both be used in the future SuperCDMS SNOLAB experiment. With their HV (high voltage) detectors SuperCDMS has the capacity to reach the nuclear recoil thresholds of 40 eV and 78 eV (see Table VII of Ref. Agnese et al. (2017)). We adopt these values here.

In order to maximize the geoneutrino signal when studying nuclear recoils, we will use data within a region of interest (ROI) defined by the recoil energies between the assumed experimental threshold energy and the maximum energy that can be imparted by any geoneutrino (see Table 1).

III Geoneutrino Signal and Background

The geoneutrino signal as well as reactor and solar neutrino background components are displayed in Figure 1. We describe them in more detail in the section below.

III.1 Geoneutrino signal

III.1.1 Spectrum

The radioactive decays of 238U, 232Th and 40K are responsible for more than 99% of the Earth’s radiogenic heat. These elements decay via a series of $\alpha$ and $\beta^{-}$ processes, which terminate at stable isotopes. Each $\beta^{-}$ decay results in the emission of $\overline{\nu}_{e}$. Comprehensive computations are needed to obtain the overall resulting $\overline{\nu}_{e}$ spectrum for each decay chain, including the spectral shapes and rates of individual decays for more than $80$ different branches in each chain Ludhova and Zavatarelli (2013). In Fig. 2 we display the resulting geoneutrino spectra for each element used in our study, as shown in Ref. Araki et al. (2005). Note that the kinematic threshold of 1.8 MeV for inverse-beta decay (IBD) interactions makes detectors based on this reaction insensitive to geoneutrinos generated from the 40K decays.

III.1.2 Flux normalization

The geoneutrino flux expected in different locations strongly depends on the amount and distribution of heat producing elements within Earth. The Earth’s surface heat flow of $47\pm 2$ TW Davies and Davies (2010) is composed partly of heat generated by radioactive decays in the crust and mantle and partly of heat from the secular planet cooling (i.e. primordial heat that is byproduct of the Earth’s formation). While the crust contributions are relatively well known, the contributions from the mantle remain highly uncertain. The predicted radiogenic heat, originating from the heat producing elements, varies greatly between different bulk silicate Earth models Sramek et al. (2013), which do not include the Earth’s metallic core. The models can be generally classified as cosmochemical, geochemical and geodynamical, which make radiogenic heat predictions of $11\pm 2$ TW, $17\pm 3$ and $33\pm 3$ TW, respectively.

To maximize the geoneutrino signal, we consider the most optimistic scenario of $100\%$ of Earth’s heat being radiogenic, i.e. Earth’s radiogenic heat comprising $\sim 47$ TW Bellini et al. (2013). This choice is compatible with current geoneutrino measurements within few standard deviations Agostini et al. (2015). The location-dependent geoneutrino flux is particularly sensitive to the amount of crust beneath the laboratory site, which contains the largest portion of heat producing elements. Due to this reason, as well as substantial distance to the closest nuclear reactors, Jinping Underground Laboratory, situated near the Himalayan mountains in China where the amount of crust is maximal, has been previously identified as a favorable location for geoneutrino detection Sramek et al. (2016). In principle, the geoneutrino flux is calculated by summing contributions from each Earth grid component, using a detailed geological map specifying heat producing element distribution and abundance. Assuming a nearly fully radiogenic Earth, we take the geoneutrino flux at Jinping (see Fig. 6 of Ref. Wan et al. (2017)) corresponding to $\gtrsim 90$ TNU (terrestrial natural units), which is $\sim 2.5$ larger (see Fig. 2 of Ref. Sramek et al. (2016)) than the computed flux for each radiogenic component at the SNOLab (Sudbury, Canada) using the Reference Earth Model (see Table 2 of Ref. Huang et al. (2013)). The geoneutrino fluxes we adopt are provided in Table 2. The corresponding flux uncertainties are taken from Ref. Huang et al. (2013).

III.2 Neutrino background

Since the geoneutrino signal ceases at neutrino energies above $\sim 5$ MeV, we shall only consider the backgrounds due to solar and reactor neutrinos (see Fig. 1), which dominate at energies $E_{\nu}\lesssim 20$ MeV. At higher neutrino energies, contributions from atmospheric and diffuse supernova neutrinos also become important.

III.2.1 Solar neutrinos

Solar electron neutrinos $\nu_{e}$ are produced as a byproduct of nuclear fusion reactions in the Sun (see Haxton et al. (2013) for review). Multiple reaction chains contribute, varying in resulting neutrino flux and energy. More than 98% of the Sun’s energy is produced through the proton-proton cycle. Starting with the initial reaction of $p+p\rightarrow~{}^{2}{\rm H}+e^{+}+\nu_{e}$, this cycle results in $pp$, $hep$, $pep$, 7Be and 8B neutrino components. The remaining solar energy is released through the Carbon-Nitrogen-Oxygen (CNO) cycle, which results in 13N, 15O and 17F neutrinos. The solar neutrino background is independent of the laboratory location.

Specific solar neutrino fluxes can be inferred from a solar model. The Standard Solar Model (GS98) Grevesse and Sauval (1998) showed a good agreement with helioseismological studies. On the other hand, more recent constructions (AGSS09) Asplund et al. (2009) that posses a higher degree of internal consistency, appear to disagree with results from helioseismology. This discrepancy is the so-called “solar metallicity” problem. The CNO neutrinos are particularly promising in this context, given their ability to shed light on the Sun’s internal composition, and their observation potential in future direct detection experiments has been recently explored Billard et al. (2015); Cerdeno et al. (2018); Newstead et al. (2019). For our analysis we assume the solar flux normalizations and uncertainties as predicted by AGSS09 (see Table 2 of Haxton et al. (2013)), which we specify in Table 3. The background in our study is dominated by 8B neutrinos. Since the 8B flux difference between the two solar models is only $\sim 20$%, we do not expect a significant change in our sensitivity results if the alternative model is employed.

III.2.2 Reactor neutrinos

Reactor electron anti-neutrinos $\overline{\nu}_{e}$ (see Hayes and Vogel (2016) for a review) originate from fission $\beta$-decays of uranium (235U and 238U) and plutonium (239Pu and 241Pu) within the reactor fuel. Since the considered isotopes are short lived, the corresponding neutrino flux follows the reactor operation. The fraction of isotope contributions as well as the average energy release per fission are given in Table 4, following Ref. Ma et al. (2013). As a reactor operates, its fuel composition and relative isotope fractions change with time, in a manner specific for each individual reactor. We neglect these effects (an analysis including these contributions can be found in Ref. Murayama and Pierce (2002)).

Unlike the flux of solar neutrinos, the reactor neutrino flux depends on the laboratory location and in particular on the characteristics and distances of the main contributing reactors Baldoncini et al. (2015). The flux of reactor anti-electron neutrinos from an isotope $k$ is given by

[TABLE]

where $R_{\bar{\nu}_{e}}$ is the emission rate, $d$ is the distance of a given reactor to the laboratory, $S_{k}(E)$ is the corresponding neutrino spectrum and $f_{k}$ is the fraction of the isotope $k$ in the reactor fuel. The neutrino emission rate is

[TABLE]

where $N_{\nu,{\rm f}}\simeq 6$ is the average number of anti-neutrinos produced per fission, $P_{th}$ is the thermal power output of the reactor, $E_{k}$ is the fission energy and $e=0.75\,(1\pm 0.1)$ is the average reactor operation efficiency that includes shut-downs use . The uncertainty of the efficiency reflects an uncertainty of about $8\%$ in the reactor anti-neutrino spectrum (see Table 3 of Ref. Kopeikin (2012)) as well as an uncertainty of few percent on the fission emission. This is the source of the reactor flux uncertainty in Table 3. Approximate analytic expressions for neutrino spectra $S_{k}(E_{\nu})$ have been constructed Vogel and Engel (1989); Mueller et al. (2011); Huber (2011); Murayama and Pierce (2002). We adopt the spectra from the model of Ref. Mueller et al. (2011), which is based on a phenomenological fit to data with an exponentiated polynomial. The resulting expression is

[TABLE]

where $\alpha_{i,k}$ is the respective fit coefficient of order $i$. The best values of the fit coefficients, as obtained in Ref. Mueller et al. (2011), are displayed in Table 5. Strictly, these fits are only valid for $E_{\nu}\gtrsim 1.8$ MeV. For uranium, Ref. Vogel and Engel (1989) does not find substantial deviations in the fit results for energies above $0.5$ MeV. For other elements, the flux at lower energies appears enhanced. Since this could be a model artifact and not of physical significance, we conservatively set $S_{k}(E_{\nu}<1.8~{}\text{MeV})=S_{k}(E_{\nu}=1.8~{}\text{MeV})$ for all elements. This effectively appears as a “kink” feature visible in the resulting flux figures. Recent measurements of the reactor anti-neutrino spectra from Daya Bay Adey et al. (2018) indicate agreement with the predictions of the leading theoretical models of Huber-Muller Huber (2011); Mueller et al. (2011) and ILL-Vogel Vogel et al. (1981); Schreckenbach et al. (1985); Von Feilitzsch et al. (1982); Hahn et al. (1989) at the level of $\sim 10$% difference on the predicted anti-neutrino yield. Since reactor background constitutes a sub-dominant contribution in our study, we do not expect reactor neutrino modeling uncertainty to significantly affect our results.

We list the reactors considered in our study for the Jinping location in Table 6.

III.3 Recoil rates

We consider two types of interactions in direct detection experiments, elastic electron neutrino scattering and coherent neutrino-nucleus scattering.

For a neutrino flux $\phi_{\nu}(E_{\nu})$ the resulting differential event rate per unit time and detector mass as a function of the recoil energy $E_{R}$, per unit time and mass $m_{I}$ of a target nuclide $I$ in a detector is given by

[TABLE]

where $d\sigma^{I}(E_{\nu},E_{R})/dE_{R}$ is the coherent neutrino-nucleus scattering differential cross-section and $C_{I}$ is the fraction of nuclide $I$ in the material. When several nuclides are present, for each type neutrino flux contribution $\phi_{\nu}(E_{\nu})$ the differential event rate is given by summing over them

[TABLE]

where MT is the experimental exposure.

A similar expression holds for electrons if we neglect the effects of their atomic and molecular binding. In the approximation of considering electrons as free in the material, we need to take into account the charge number $Z_{I}$ in every nuclide $I$. Hence, the interaction rate in a detector becomes

[TABLE]

where $d\sigma^{e}(E_{\nu},E_{R})/dE_{R}$ is the neutrino-electron elastic scattering differential cross-section.

For a target mass $m$ at rest ($m_{I}$ or the electron mass $m_{e}$), the minimum neutrino energy required to produce a recoil of energy $E_{R}$ is

[TABLE]

The maximum recoil energy due to a collision with neutrino of energy $E_{\nu}$ is

[TABLE]

III.3.1 Coherent nuclear scattering cross-section

The Standard Model coherent-scattering neutrino-nucleus cross-section is given by Freedman et al. (1977)

[TABLE]

where $m_{I}$ is target nuclide mass, $G_{f}$ is Fermi coupling constant, $F_{{\rm SI},I}(E_{R})$ is the form factor, which we take this to be the Helm form factor Helm (1956), $Q_{w}=[(1-4\sin^{2}\theta_{\rm W})Z_{I}-N_{I}]$ is the weak nuclear charge, $N_{I}$ is the number of neutrons, $Z_{I}$ is the number of protons and $\theta_{\rm W}$ is the Weinberg angle. Since $\sin^{2}\theta_{\rm W}=0.223$ Patrignani et al. (2016), the coherent neutrino-nucleus scattering cross-section follows an approximate $N_{I}^{2}$ scaling.

III.3.2 Elastic electron scattering cross-section

The Standard Model neutrino-electron scattering cross-section is Vogel and Engel (1989)

[TABLE]

Here, the flavor-dependent $g_{\nu}$ weak couplings

[TABLE]

are enhanced for the electron flavor, while

[TABLE]

The corresponding anti-neutrino cross-sections are obtained via $g_{\overline{\nu}(\nu)}\leftrightarrow g_{\nu(\overline{\nu})}^{\prime}$ interchange.

The electron-neutrino cross-section of Eq. (III.3.2) is valid only for a neutrino scattering with a free electron. However, the target electrons of interest in this study are bound to nuclei. For larger recoil energies and smaller nuclei the atomic effects are expected to become negligible. On the other hand, for recoil energies of few keV in germanium Chen et al. (2014a, b, 2015) and $\mathcal{O}(10)$ keV in xenon Chen et al. (2017) the cross sections can become noticeably suppressed. Since we do not perform a statistical analysis for electron recoils, this effect is neglected.

III.3.3 Neutrino oscillations

The coherent neutrino-nucleus scattering is mediated by the $Z$ boson. Thus, this interaction is independent of the neutrino flavor and the corresponding oscillation effects. On the other hand, electron-neutrino interactions are sensitive to oscillation effects.

Due to the long range of production distances, the survival probability of $\overline{\nu}_{e}$ geoneutrinos due to averaged out oscillation effects is approximately

[TABLE]

assuming $\theta_{\rm 13}=8.5^{\circ}$ and $\theta_{\rm 12}=34.5^{\circ}$ de Salas et al. (2018). The Earth’s Mikheyev-Smirnov-Wolfenstein (MSW) oscillation effects for geoneutrinos are negligible Wan et al. (2017). On average, $\langle P_{ee}\rangle\simeq 0.55$ is also valid for reactor neutrinos. For solar neutrinos we also take an electron neutrino survival probability of $\langle P_{ee}\rangle\simeq 0.55$ for all fluxes aside 8B, since they lie in the region of vacuum oscillations of $E_{\nu}\lesssim 1$ MeV. For 8B neutrinos we take $\langle P_{ee}\rangle\simeq 0.35$, since their energies are in the matter-enhanced oscillation region of $E_{\nu}\gtrsim 5$ MeV. This choice is consistent with the large mixing angle (LMA-)MSW solution, as indicated by measurements Agostini et al. (2018).

IV Predicted recoil spectra

IV.1 Electron recoils

The left panel of Fig. 3 shows the recoil spectrum of the three types of geoneutrinos in germanium (dashed lines). The electron recoil spectrum is almost the same for all target elements. This can be seen in the right panel of the same figure, where we show the total electron recoil spectrum (the sum of all three geoneutrino contributions) for several elements. The reason is that the ratio $Z_{I}/m_{I}\simeq Z_{I}/A_{I}m_{p}$, where $m_{p}$ is the proton mass, which enters in the rate is approximately 1/2$m_{p}$ for all elements and the differential cross section Eq. (III.3.2) is the same, neglecting binding effects.

The backgrounds for a geoneutrino signal in electron recoils due to contributions other than neutrinos are expected to be large. Thus, we do not examine detectability of geoneutrinos through electron scattering.

IV.2 Nuclear recoils

Fig. 4 shows the three geoneutrino nuclear recoil spectra (dashed lines) for Xe, Ar, Ge and Si. The total (sum of all three contributions) geoneutrino nuclear recoil spectra for these and other target elements are shown in Fig. 5. These spectra, and the corresponding differential cross sections, differ considerably between various elements. With the threshold we are using (see Table 1), Fig. 4 shows that 40K geoneutrino signal cannot be detected either in Xe or Ar, and that 232Th signal cannot be observed in Xe. The lowest threshold of 0.6 keV currently employed in Xe-based experiments would not allow to detect any geoneutrino signal. The most promising targets to detect geoneutrinos are therefore Ge and Si, because of the achievable low threshold as demonstrated by the SuperCDMS SNOLAB HV (high voltage) detectors.

V Detection Sensitivity for Nuclear Recoils

V.1 Statistical analysis

To determine the discovery sensitivity of a direct detection experiment we use a frequentist analysis based on the profile likelihood ratio test Cowan et al. (2011a); Rolke et al. (2005); Cowan et al. (2011b) used recently in Ref. Gelmini et al. (2018) and other studies Billard et al. (2012, 2014); Aprile et al. (2011); Ruppin et al. (2014). This test is performed by generating simulated data sets for the experiment, assuming a normalization (i.e. an energy integrated neutrino flux) for each of the considered neutrino fluxes ${\phi}_{\nu_{k}}$, which we take here to be the average theoretical predicted value $\overline{\phi}_{\nu_{k}}$. Here $k=1,2,\dots,n_{\nu}$, with $n_{\nu}$ being the total number of neutrino contributions in Table 3 as well as the geoneutrino fluxes in Table 2 that produce recoils above the particular experimental threshold. More specifically, for each of the ${\phi}_{\nu_{k}}$ parameters we generate fake data consisting of a total number of “observed” events $N_{\rm o}$ at particular recoil energies $E_{j}$ with $j=1,2,\dots,N_{\rm o}$ that we use to define the following likelihood function

[TABLE]

Here, $n_{g}$ is the number of geoneutrino species that produce recoils above threshold and $k^{\prime}$ numbers the fluxes of neutrinos other than geoneutrinos, with uncertainties $\sigma_{\nu_{k^{\prime}}}$ taken into account. The total differential rate is the sum of the contributions of the signal due to geoneutrinos, denoted by the index $k_{g}$, and the background, assumed here to be due only to all other types of neutrinos,

[TABLE]

defined in Eq. (5), respectively. The total number of predicted events $N_{\rm E}$ is obtained by integrating Eq. (15) over our recoil energy ROI and multiplying by the assumed exposure MT. In Eq. (V.1) the extended likelihood (i.e. the term in the curly brackets) is multiplied by a Gaussian product of likelihoods, centered around the mean predicted flux normalization $\overline{\phi}_{\nu_{k^{\prime}}}$, for each background neutrino species $\nu_{k^{\prime}}$ ($k^{\prime}=1,2,\dots,n_{\nu}-n_{g}$) to take into account the uncertainty in the $\phi_{\nu_{k^{\prime}}}$ fluxes. In the Gaussian likelihoods the $\sigma_{\nu_{k^{\prime}}}$ adopted is the 1$\sigma$ uncertainty in the particular flux given in Table 3 (see Section III.2).

As explained in Ref. Gelmini et al. (2018), the procedure we follow to obtain each set of simulated data involves two steps: 1) find the total number of events of each type $k$ (where $k$ numbers all neutrino types above threshold), and 2) find the recoil energy for each of the fake data events. For each realization of fake data, the number of events of a specific type, $n_{k}$, is found from a Poisson distribution $P_{k}$,

[TABLE]

where the mean $\mu_{k}$ is the number of events predicted by the model being tested, defined by the particular average total neutrino flux $\bar{\phi}_{\nu_{k}}$, the target element, energy range and assumed exposure MT. Choosing a random number for the cumulative probability distribution of each $P_{k}$ one value of $n_{k}$ is randomly generated, according to an inverse transform sampling. The total number $N_{0}$ of “observed” events is then $N_{0}=\sum_{k=1}^{n_{\nu}}n_{k}$.

To determine the energy of each one of the $n_{k}$ events of each neutrino type we use as probability density function (PDF) the corresponding differential recoil rate $dR_{k}/dE_{R}$ normalized by the total rate (i.e. the rate integrated over the specified energy ROI for the given experiment $R_{k}$) as

[TABLE]

The corresponding $n_{k}$ recoil energies are again obtained with inverse transform sampling. With the above procedure, we simulate 600 data sets for each exposure and experimental configuration.

For each simulated data set we define a test statistic $q_{0}$ that allows one to reject the background only hypothesis $H_{0}$, assuming it is true, with a probability not larger than some value $\alpha$ that denotes the significance level of the test. In our case, $H_{0}$ is the assumption that the flux of geoneutrinos is zero, ${\phi}_{\nu_{k_{g}}}=0$, For our analysis we chose $\alpha$ to correspond to either 2$\sigma$, or 3$\sigma$ (i.e. $\alpha=0.0228$ and $\alpha=0.00135$, respectively). The test statistic $q_{0}$ for each simulated data set is the profile likelihood ratio, defined as

[TABLE]

where the $\phi_{\nu_{k^{\prime}}}$ are treated as nuisance parameters. The hats refer to values that maximize the likelihood, with double-hats referring to values maximizing the likelihood subject to the constraint ${\phi}_{\nu_{k_{g}}}=0$. Note that by definition $q_{0}\geq 0$, with larger values of $q_{0}$ indicating greater incompatibility of the simulated data with the background only hypothesis $H_{0}$. Thus, we require that the probability $p_{0}$ of having a $q_{0}$ value larger (i.e. more incompatible with the data if due to background only) than the “observed” (i.e. than the $q_{0}$ of the simulated dataset) $q_{0}^{\rm obs}$

[TABLE]

is not larger than $\alpha$, $p_{0}\leq\alpha$. Here, $f(q_{0}|H_{0})$ is the PDF of obtaining $q_{0}$ under the background-only hypothesis $H_{0}$. In the large sample limit, Wilks’ theorem ensures that $f(q_{0}|H_{0})$ is given by one half times a delta function at $q_{0}=0$ plus one half times a $\chi^{2}$ distribution with $n_{g}$ degrees of freedom Cowan et al. (2011a), with $n_{g}$ corresponding to the number of geoneutrino species contributing to the signal. For xenon $n_{g}=1$, argon $n_{g}=2$ and for germanium and silicon $n_{g}=3$. In the analyses detailed below, we require a value of $\alpha$ corresponding to a $Z\sigma$ significance with $Z=2$ or 3, namely $\alpha(2)=0.0228$ and $\alpha(3)=0.00135$. The value of $\alpha$ corresponding to Z$\sigma$ for $n_{g}=1$ (namely for 1 degree of freedom) can be obtained by requiring, approximately, that ${q_{0}^{\rm obs}}\geq Z^{2}$. No similar approximate condition exists for $n_{g}=2$ and $n_{g}=2$ and the lower limit of ${q_{0}^{\rm obs}}$ is found numerically.

In the Fig. 6 we show the fraction $F$ of simulated data sets generated with a particular exposure MT that fulfill this requirement, for $Z=2$ or 3. Namely, we determine the fraction of simulated data sets that produce p-values $p_{0}^{i}\leq\alpha(Z)$ (i.e. those resulting in $\geq Z\sigma$ detection of geoneutrinos), by computing

[TABLE]

where $N_{\rm sim}$ is the number of simulated data sets with fixed parameters $({\phi}_{\nu_{k_{g}}},{\rm MT})$.

Notice that the fraction $F$ is the “power” of the test of $H_{0}$ with respect to the alternative hypothesis $H_{\sigma}$, in which ${\phi}_{\nu_{k_{g}}}\neq 0$ (e.g. Patrignani et al. (2016), Sec. 40), i.e. $F=(1-\beta)$. This means that the probability of not rejecting $H_{0}$ when the “geoneutrino detection” alternative hypothesis $H_{\sigma}$ is true is less than the value $\beta=1-F$.

V.2 Results

Fig. 6 shows the two bands of fractions $F$ of simulated data sets generated with a particular exposure MT in which there is a detection of geoneutrinos at an x$\sigma$ level with x $\geq Z$, for $Z$=2 (in blue) and $Z$=3 (in purple), as a function of the exposure MT, for Ar, Ge and Si compositions. The higher exposure boundary (solid line) in each band assumes the systematic errors shown in Table 3. The lower exposure boundary (dashed line) of each band assumes instead a very optimistic scenario in which the errors have been divided by a factor of 100. The gray shaded area highlights the region of $F>0.8$, where the chance of a type-2 error is small (i.e $\beta<0.2$).

With a Xe target even a 2$\sigma$ detection of geoneutrinos is not reachable below a 500 tonne-year exposure, the largest we consider. The reason is that for the 0.1 keV threshold that we assume only a small portion of 238U geoneutrino signal is observable (see Fig. 4).

Fig. 6 shows that the most promising detector configuration for geoneutrino detection is based on Ge, because of the lowest potential energy threshold.

VI Summary

A measurement of the level of geoneutrino emission is of utmost importance to determine the fraction of the heat emitted by Earth that remains from its formation. The major uncertainty is due to the unknown mantle contribution to the radiogenic heat.

Direct detection experiments have a clear advantage over scintillator-based experiments, which have been used to study geoneutrinos thus far, in that their threshold could be low enough (at least in Ge- and Si-based detectors) to observe geoneutrinos with energies below the kinematic threshold of 1.8 MeV, the lowest neutrino energy required for inverse beta decay. This means that 40K geoneutrinos can also contribute to the signal in direct detection experiments. A measurement of the present level of 40K with respect to the levels of other geoneutrino types could provide insight into Earth’s “volatility curve”, a central quantity for Earth’s formation history Bellini et al. (2013).

In this paper we computed the electron and nuclear recoil spectra in different target materials. In Fig. 1 one can see that the geoneutrino electron recoil rate is about one order of magnitude lower than the rate of CNO neutrinos, whose detection sensitivity has been recently explored Billard et al. (2015); Cerdeno et al. (2018); Newstead et al. (2019).

For the detection of nuclear recoils from geoneutrino signal we considered several potential experimental setups with optimistic thresholds and located in the Jinping Underground Laboratory in China, where the expected geoneutrino flux is significantly larger than in other laboratory sites, because of its proximity to the Himalayas, where the crust is the thickest. Assuming that radiogenic heat accounts for most of the heat emitted by our planet, a viable possibility, we study the minimal required exposures for a geoneutrino observation at the 2 and 3$\sigma$ level, if the background is due only to other neutrinos (solar and reactor neutrinos).

The results of our analysis are shown in Fig. 6. The most promising target material analyzed in this study for geoneutrino detection is germanium, which could have the lowest experimental threshold, if the exposure reaches $\mathcal{O}(10)$ tonne-years.

Acknowledgements

We thank Carolina Lithgow-Bertelloni, Lars Stixrude, Gerald Schubert and Linyan Wan for discussions and Sergio Palomares-Ruiz for comments. The work of GBG and VT was supported in part by the U.S. Department of Energy (DOE) under Grant No. DE-SC0009937. SJW acknowledges support from Spanish Ministry of Economy and Competitiveness (MINECO) national grants FPA2014-57816-P and FPA2017- 85985-P, and from the European projects H2020-MSCAITN-2015//674896-ELUSIVES and H2020-MSCA-RISE2015.