Can the ANITA anomalous events be due to new physics?

TL;DR

This paper critically examines the possibility that new physics explains ANITA's anomalous events, highlighting challenges and proposing a specific dark matter decay model as a potential explanation.

Contribution

It provides general arguments against most new physics explanations and identifies a unique dark matter decay scenario that could account for the observations.

Findings

Most new physics models face significant challenges explaining the events.

A metastable dark matter decay model can produce the observed signals.

Standard Model neutrinos cannot account for the events.

Abstract

The ANITA collaboration has observed two ultra-high-energy upgoing air shower events that cannot originate from Standard Model neutrinos that have traversed the Earth. Several beyond-the-standard-model physics scenarios have been proposed as explanations for these events. In this paper we present some general arguments making it challenging for new physics to explain the events. One exceptional class of models that could work is pointed out, in which metastable dark matter decays to a highly boosted lighter dark matter particle, that can interact in the Earth to produce the observed events.

Tables

Table 1. Table 1: Maximum probability is computed for the AAE with θ em = 27 ∘ subscript 𝜃 em superscript 27 \theta_{\rm em}=27^{\circ} ( ℓ c = 5785 ​ km subscript ℓ 𝑐 5785 km \ell_{c}=5785\leavevmode\nobreak\ \mathrm{km} ). The p 𝑝 p -value and significance σ 𝜎 \sigma take into account the two AAEs. The values in the process ⟨ 2 ⟩  and ⟨ 3 ⟩  with different numbers of steps are given.

# steps	max. prob.  ⟨2⟩	p-value ⟨2⟩	significance ⟨2⟩	max. prob. ⟨3⟩	p-value ⟨3⟩	significance ⟨3⟩
1	$3.2\times {10}^{-4}$	$3.7\times {10}^{-8}$	5.5 $\sigma$	$1.4\times {10}^{-3}$	$5.1\times {10}^{-6}$	4.6 $\sigma$
2	$5\times {10}^{-4}$	$4.8\times {10}^{-8}$	5.5 $\sigma$	$2.1\times {10}^{-3}$	$1.5\times {10}^{-6}$	4.8 $\sigma$
3	$6\times {10}^{-4}$	$6.5\times {10}^{-8}$	5.4 $\sigma$	$2.5\times {10}^{-3}$	$1.3\times {10}^{-6}$	4.8 $\sigma$
4	$7\times {10}^{-4}$	$8.1\times {10}^{-8}$	5.4 $\sigma$	$3.0\times {10}^{-3}$	$1.4\times {10}^{-6}$	4.8 $\sigma$
5	$7\times {10}^{-4}$	$9.3\times {10}^{-8}$	5.3 $\sigma$	$3.3\times {10}^{-3}$	$1.5\times {10}^{-6}$	4.8 $\sigma$

Equations

$$d{H}_{\text{AN}}=\mathrm{d}A\mathrm{d}{\Omega }_d\sin {\theta }_{\mathrm{em}}\mathrm{d}t{P}_{\mathrm{obs}},$$

$$\mu =\int \mathrm{d}{H}_{\text{AN}}\mathrm{d}{E}_{\nu }\varphi (E_{\nu }).$$

$$\big{\langle}0\big{\rangle}\qquad\nu\rightarrow\tau\to{\rm hadrons}\,,\qquad\qquad\qquad\qquad\qquad\qquad$$

$$\big{\langle}1\big{\rangle}\qquad\nu\rightarrow X\rightarrow\tau\to{\rm hadrons}\,.\qquad\qquad\qquad\qquad\quad$$

$$\big{\langle}2\big{\rangle}\qquad\nu\rightarrow X_{1}\rightarrow X_{2}\to\dots\to X_{n}\ \to\tau\to{\rm hadrons}\,.$$

$$\big{\langle}3\big{\rangle}\qquad\nu\rightarrow X_{1}\rightarrow X_{2}\to\dots\to X_{n}\ \to{\rm hadrons}\ .\ \ \ \quad$$

$${\ell }_{BSM}=({\sigma }_{BSM}{n}_{\oplus }{)}^{-1}$$

$${\ell }_{SM}=({\sigma }_{SM}{n}_{\oplus }{)}^{-1}.$$

$$P_{\mathrm{SM}}$$

$$\frac{\mathrm{d}{P}_{\nu \to \tau }}{d{z}_1}$$

$$\frac{\mathrm{d}{P}_{\mathrm{exit}}}{\mathrm{d}{\bar{E}}_{\tau }}$$

$$\frac{\mathrm{d}{P}_{\tau }}{\mathrm{d}{z}_3}$$

$$P^{(1A)}$$

$$\frac{d{P}_{\nu \to X}}{\mathrm{d}{z}_1}$$

$$\frac{d{P}_{X\to \tau }}{\mathrm{d}{z}_2}$$

$$P_{SM}^{(B)}$$

$$\frac{\mathrm{d}{P}_{\nu \to \tau }^{\prime }}{\mathrm{d}{z}_1}={\ell }_{SM}^{-1}\exp \left(-z_1({\ell }_{SM}^{-1}+{\ell }_{BSM}^{-1})\right).$$

$$x\equiv {\ell }_{BSM}/{\ell }_{SM},y\equiv {\ell }_X/R_{\text{varoplus}}.$$

$$H_{\text{AN}}^{\text{ref}}\simeq 2.7{\text{km}}^2\text{yr}\text{sr}.$$

$$\mu \simeq \frac{\int \mathrm{d}A{P}_{\mathrm{obs}}({\theta }_{\text{em}})\sin {\theta }_{\text{em}}}{\int \mathrm{d}A\sin {\theta }_{\text{em}}}{H}_{\text{AN}}^{\text{ref}}{E}_{\nu }\varphi (E_{\nu }).$$

$$L(n,{\theta }_{\mathrm{em}})=\frac{e^{-\mu }{\mu }^n}{n!}\prod _{i=1}^{n}f({\theta }_{\text{em}}^i),$$

$$f({\theta }_{\text{em}})=\frac{P_{\mathrm{obs}}({\theta }_{\text{em}})\sin {\theta }_{\text{em}}\sin {\theta }_E\mathrm{d}{\theta }_E/\mathrm{d}{\theta }_{\text{em}}}{\int \mathrm{d}{\theta }_E\sin {\theta }_E{P}_{\mathrm{obs}}({\theta }_{\text{em}})\sin {\theta }_{\text{em}}},$$

$$p=\sum _{n=0}^{\infty }\int L(n,{\theta }_{\mathrm{em}})\prod _{i=1}^{n}d{\theta }_{\mathrm{em}}^i$$

$$P(n;\mu )\prod _{i=1}^{n}f({\theta }_{\mathrm{em}}^i)<P(2;\mu )f(27^{\circ })f(35^{\circ })$$

$${\mu }_{\mathrm{ps}}$$

$$M_{\chi }=ϵ{\rho }_{\text{DM}}^{\text{local}}{v}_{\text{rel}}\pi R_{\text{varoplus}}^2{t}_{\text{varoplus}}\lesssim 6\times 10^{12}\mathrm{g},$$

$$M_{\chi }=\frac{m_{\chi }{\tau }_{\chi }{N}_{\text{A}}}{t_{\text{A}}{\Omega }_{\text{A}}{P}_{\text{A}}}\simeq 2\times 10^{17}\mathrm{g}\left(\frac{{\tau }_{\chi }}{2\times 10^{27}\mathrm{s}}\right)\left(\frac{10^{-3}}{P_{\text{A}}}\right),$$

$$\frac{1}{{\Lambda }^2}(\bar{\chi }{\chi }^c)(\bar{\chi }\Psi )+y_{\tau }{\bar{L}}_{\tau }\varphi \chi +h.c.$$

$$\Lambda =2\times 10^{21}\mathrm{GeV}{\left(\frac{m_{\Psi }}{3\mathrm{EeV}}\right)}^{5/4}{\left(\frac{{\tau }_{\Psi }}{5\times 10^{17}\mathrm{s}}\right)}^{1/4}.$$

$$\Psi$$

$$\chi$$

Full text

††preprint: CERN-TH-2019-057

I Introduction

High energy neutrino astrophysics could provide a window to new physics at energies beyond those accessible on Earth, thanks to the long propagation distances possible for ultra-high-energy (UHE) neutrinos, compared to charged particles. The Auger and IceCube experiments are sensitive to UHE neutrinos, with several events at PeV energies now observed by IceCube. Upper limits on the flux at higher energies are established, up to energies of $10^{11}$ GeV.

The cross section for neutrinos to interact with nucleons through deep inelastic scattering (DIS) Connolly et al. (2011) grows with energy, and at ultra high energies above $\sim 10^{8}$ GeV the Earth becomes opaque to neutrinos. A promising strategy then is to search for Earth-skimming neutrinos that produce charged leptons with high probability, that do not lose too much energy before they exit Feng et al. (2002). The vast majority of such events would be initiated by $\nu_{\tau}$ converting to $\tau$ leptons via DIS since electrons or muons would be absorbed by the ice over much shorter distances and have a negligible probability to exit compared to $\tau$’s.

ANITA is a balloon-borne instrument that detects polarized radio emission from the electromagnetic component of cosmic ray showers in the atmosphere, reflected from ice at the South Pole, or arriving without reflection from events near the horizon. Because of phase inversion of the signal upon reflection, the difference between such events is distinguishable, and the zenith angle of the $\tau$ whose decay produced the shower can be determined to within $0.3^{\circ}$.

ANITA has made four flights since 2007, and observed two anomalously steep events, henceforth called ANITA anomalous events (AAEs), apparently emerging from the Earth, of energy $E_{\nu}\sim 0.6$ EeV on the first Gorham et al. (2016) and third Gorham et al. (2018) flights, corresponding to particles that traversed chord lengths of 5700 and 7200 km, respectively.111Recently ref. de Vries and Prohira (2019) argued that the observed inverted polarity of the AAEs might in fact be explained from down-going air showers. The mean free path of neutrinos at this energy is a few hundred km, making it exceedingly unlikely that standard model processes could explain these events Fox et al. (2018).

Several models beyond the Standard Model (SM) have been suggested as possible explanations. These fall into three broad categories: (1) SM neutrinos of astrophysical origin could convert to beyond-the-standard-model (BSM) particles through interactions in the Earth, followed by propagation of the BSM particle until it reconverts to SM particles that initiate the observed hadronic air shower Fox et al. (2018); Collins et al. (2019); Chauhan and Mohanty (2018); Anchordoqui and Antoniadis (2019); (2) dark matter (DM) that has accumulated within the Earth and decays to BSM particles that reconvert to SM particles shortly below the Antarctic surface could induce the air showers; or (3) an exotic flux of BSM particles, such as sterile neutrinos Cherry and Shoemaker (2018); Huang (2018) incident upon the Earth interact to produce the observed particles near the Earth’s surface. A possible source of such flux is the decay of long-lived DM particles Dudas et al. (2018); Anchordoqui et al. (2018); Heurtier et al. (2019).

A generic challenge for models in the first category is that, despite the increased probability for a BSM particle to traverse the Earth compared to EeV-energy neutrinos, there is still a large reduction in efficiency because of the need for this particle to reconvert within a thin layer near the surface of the Earth before it exits into the atmosphere. A major goal of the present work is to carefully quantify this efficiency and to show that it puts such explanations of the ANITA anomalous events at odds with null searches by IceCube and other experiments such as Auger. Models in the second category face the problem that the amount of DM particles that may have accumulated within the Earth cannot be large enough to explain the two AAEs. On the other hand, third category models have the advantage that the small efficiency factor can be overcome by assuming a larger flux of sterile BSM particles, that is relatively less constrained than the neutrino flux. Nevertheless such a large flux of sterile BSM particles can be problematic if at the same time active neutrinos or other nonsterile particles are produced. We will argue that previously proposed models of this kind are unlikely to be viable, but demonstrate an example of a model that does work. The essential ingredient is an ultraheavy DM particle that decays exclusively to another, much lighter, DM particle, whose flux is unconstrained. This boosted decay product then converts to $\tau$’s within the Earth which in turn induce the hadronic air showers.

The paper is structured as follows: In Section II, we perform a model-independent analysis of scenarios where a neutrino scatters into weakly-interacting states that can traverse the Earth. These states will convert to $\tau$ or directly decay to hadrons seen by ANITA. Single conversion and cascade decay processes are considered here, but neither diffuse nor pointlike $\nu_{\tau}$ fluxes can explain the AAEs. Constraints on models that attempt to explain the AAEs from dark matter decays either within the Earth or in the galactic halo are discussed in Section III. In either case, we find that it is not possible to achieve a large enough flux to explain the ANITA events, while remaining consistent with limits from IceCube and other experiments. In Section IV we give an example of a class of models that can however be consistent, involving decays of a heavy subdominant component of metastable DM into the lighter dominant DM particle, that can interact with matter in the Earth to produce the anomalous events. We conclude in Section V.

II ANITA events from neutrino flux

In this section, we analyze the possibility that the ANITA events are consistent with BSM explanations, assuming that the BSM particles come from ultra-high-energy (UHE) neutrinos interacting with matter in the Earth. Other sources of BSM particles will be considered in following sections.

We assume that the ANITA anomalous events are air showers initiated by the hadronic decay of energetic $\tau$ leptons, or alternatively by the decay of BSM particles directly into hadrons.222The energetic $\tau$’s which induce air showers above the Antarctic surface can also pass through the IceCube detector, where they could be misidentified as highly energetic muon tracks. At most three such events have been seen Fox et al. (2018), while a rough estimate of the number of events is approximately one order of magnitude larger than the number of AAEs Huang (2018); Fox et al. (2018). It is however possible that a detailed analysis, including in particular instrumental effects, could lead to a lower number of expected events at IceCube, thus resolving this tension and making $\tau$-induced air showers viable (cf. ref. Fox et al. (2018)). Hence we discuss both the case of $\tau$-induced air showers and air showers directly induced by a BSM particle in this paper. The predicted number of such events is the convolution of the differential exposure of ANITA to up-going air showers,

[TABLE]

and the differential diffuse neutrino flux $\phi(E_{\nu})\,\mathrm{d}E_{\nu}$:

[TABLE]

Here

• •

$\mathrm{d}A=R_{\varoplus}^{2}\sin\theta_{E}\mathrm{d}\theta_{E}\mathrm{d}\phi_{E}$ is the differential area of the Earth observed by ANITA, with $\theta_{E}$ being the angle between the South Pole and a point on the surface as seen from the center of the Earth, $\phi_{E}$ being the longitude, and $R_{\varoplus}=6371$ km the Earth’s radius;

• •

$\theta_{\rm em}$ is the emergence angle of the $\tau$ (respectively the BSM particle that causes the hadronic air shower), measured with respect to the horizon—see fig. 1;

• •

$P_{\rm obs}$ is the probability that an UHE $\nu$ that enters the Earth with energy $E_{\nu}$ leads to a $\tau$ (or possibly some BSM particle) that exits the Earth near the south pole at an emergence angle $\theta_{\rm em}$ with energy $\bar{E}_{\tau}$ and initiates a hadronic air shower within $D=10\leavevmode\nobreak\ {\rm km}$;333For altitudes larger than $6\,$km the typical shower is already starting to exit the atmosphere before reaching its maximum Gorham et al. (2016). We choose $D=10\leavevmode\nobreak\ {\rm km}$ in order to be conservative in constraining BSM scenarios.

• •

${\rm d}\Omega_{d}$ is the detection solid angle of ANITA. It depends on $\bar{E}_{\tau}$ and the $\tau$’s decay position, which requires detailed simulations of Extensive Air Showers (EAS) Romero-Wolf et al. (2018). Here we neglect this dependence, making the approximation that ${\rm d}\Omega_{d}$ can be integrated independently.

The live time of the two experiments was $17.25$ days for ANITA I Schoorlemmer et al. (2016) and $7$ days for ANITA III Gorham et al. (2018). ANITA II was not sensitive to up-going air showers, and data from ANITA IV have not been released so far.

II.1 Conversion probability

In the following we will compute upper limits on $P_{\rm obs}$ for several “process topologies,” independent of a specific BSM realization. In the SM the process topology is

[TABLE]

which is understood to include $\tau$ regeneration and energy loss. The BSM process will occur in parallel with any new physics processes, and so its effects should be included alongside the latter for quantitative results, even if the SM contribution is quite small.

Assuming the presence of BSM physics, several topologies are conceivable. In the first, the UHE $\nu$ converts to a BSM particle $X$, which subsequently converts or decays to a $\tau$:

[TABLE]

A second possibility is to have cascade decays of several BSM particles,

[TABLE]

This includes ⟨1⟩ as the special case $n=1$. As we will show below, for $n>1$ one can increase $P_{\rm obs}$ somewhat. It is also possible that $X$ decays directly into hadrons that initiate the observed air showers. We denote this by

[TABLE]

The initial scattering process for BSM physics that converts $\nu$ to $X$ is parametrized by the interaction length

[TABLE]

where $\sigma_{\scriptscriptstyle BSM}$ is the BSM cross section for scattering of $\nu$ on nucleons (or possibly electrons) and $n_{\oplus}$ is the number density of target particles in the Earth. In case of more than one BSM conversion, such as cascade decays, the $i$th interaction length is denoted as $\ell_{i}$. Similarly, the mean free path for the standard model DIS interaction is denoted by

[TABLE]

These interaction lengths generally depend upon $E_{\nu}$.

⟨0⟩** The SM process $\nu\rightarrow\tau\rightarrow{\rm\bf hadrons}$.** The total conversion probability $P_{\rm SM}$ can be expressed as

[TABLE]

where $B_{h}\simeq 0.65$ is the branching ratio for a $\tau$ to decay into hadrons (we count only the hadronic decays for initiating the kind of air showers observed by ANITA) and the chord length between the entry point of the UHE $\nu$ and the exit point of the $\tau$ is $\ell_{c}=2R_{\oplus}\sin\theta_{\rm em}$. For illustration see the line labeled 0 in fig. 1. The probability distribution for the incoming neutrino to scatter into a $\tau$ at position $z_{1}$ is

[TABLE]

and the probability for $\tau$ to exit with energy $\bar{E}_{\tau}$, starting from initial energy $E_{\nu}$, is given by Feng et al. (2002)

[TABLE]

in terms of the $\tau$ decay length $\ell_{\tau}(E_{\tau})=c\tau_{\tau}E_{\tau}/m_{\tau}$. (Here we neglected the fact that the $\tau$ on average receives only 80% of $E_{\nu}$ in the interaction Gandhi et al. (1996).) The Earth density is taken as a constant for simplicity, $\rho_{\oplus}=2.7\leavevmode\nobreak\ {\rm g/cm^{3}}$, and $\beta_{\tau}\sim(0.4-0.8)\times 10^{-6}$ cm2/g in the energy range of interest. We take $\beta_{\tau}=0.6\times 10^{-6}$ cm2/g for our estimates. The probability for the $\tau$ to decay in the atmosphere at position $z_{3}$ is

[TABLE]

The probability distribution (5) is integrated over $\bar{E}_{\tau}$ between the initial value $E_{\nu}$ and a minimum value $E_{\rm min}\sim 0.1\,$EeV that could still be associated with the observed AAEs.

With these analytic formulae, we can reproduce the conversion probability given by the Monte Carlo simulations in ref. Alvarez-Muñiz et al. (2018) with vanishing ice thickness for \theta_{\rm em}\lesssim$$$. Unlike the simulations, we neglect \tauregeneration Halzen and Saltzberg ([1998](#bib.bib19)), which is important for larger\theta_{\rm em}$and low$\bar{E}{\tau}$. After exiting the Earth, the$\tau$does not suffer further energy loss, and it will decay in the atmosphere. The resulting probabilities are shown as a function of$\theta{\rm em}in fig. [2](#S2.F2), in whichP_{\rm SM}$drops steeply at large$\theta_{\rm em}$. Even though$\tau$regeneration increases the probabilities for large$\theta_{\rm em}$, this falls far short of being able to explain the AAEs within the SM Fox et al. (2018).

⟨1⟩** The BSM process $\bm{\nu\rightarrow X\rightarrow\tau\rightarrow{\rm\bf hadrons}}$.** To obtain a higher conversion probability at large $\theta_{\rm em}$, we consider the possibility that the UHE $\nu$ can convert to a weakly interacting particle $X$ with a larger mean free path than neutrinos, followed by $X$ decaying or rescattering into a $\tau$. The conversion probability can be written as

[TABLE]

with superscript $A$ denoting that we neglect the SM contribution to the total conversion probability here (see eq. (12) below). The probability distributions are given by

[TABLE]

$dP_{\nu\to X}$ is the differential probability that $\nu$ converts to $X$, determined by the mean free path for BSM interactions $\ell_{\scriptscriptstyle BSM}$ (in addition to the damping by the SM contribution to the $\nu$ conversions). $\mathrm{d}P_{X\to\tau}$ quantifies the decay or rescattering of $X$ (produced at $z_{1}$) into $\tau$ at $z_{2}$, with mean free path $\ell_{X}$. In the rescattering scenario, there are two possible subcases: the interaction may simply be the inverse of the one that produced $X$ (if the original interaction involved the full third generation SU(2)L doublet $(\nu_{\tau},\tau)$ and not just the right-handed component of $\tau$), or it may be a different interaction. In the first case, one has $\ell_{X}\cong\ell_{\scriptscriptstyle BSM}$, so this can be considered as a special case of the more general situation. The general situation is like the second case, where $\ell_{X}$ and $\ell_{\scriptscriptstyle BSM}$ are independent parameters. $\mathrm{d}P_{\rm exit}/\mathrm{d}E_{\tau}$ is to account the energy loss of $\tau$ as in eq. (7) with $z_{1}\to z_{2}$.

We must also include the SM component, similar to \big{\langle}0\big{\rangle},

[TABLE]

which is however modified by the new physics because the mean free path for $\nu\to\tau$ is decreased by the BSM contribution,

[TABLE]

We omit the number of the process in the superscript of $P^{(B)}_{SM}$, since this same modification will apply to all three BSM scenarios.

The total probability $P^{(1A)}+P^{(B)}_{SM}$ can be expressed in terms of dimensionless variables $x$ and $y$,

[TABLE]

For an incident $3\leavevmode\nobreak\ \mathrm{EeV}$ neutrino to generate a $\tau$ decaying in the lower atmosphere, at a fixed value of $\theta_{\rm em}=27^{\circ}$ representative of an AAE, we find the largest possible probability $\simeq 3.2\times 10^{-4}$, around $x\sim 0$, $y\sim 0.906$. The value of $y$ is close to the chord length $\ell_{c}(27^{\circ})\simeq 0.906R_{\varoplus}$, which makes $X$-decays close to the surface of the Earth likely. At $27^{\circ}$, the BSM probability $P^{(1A)}$ (which is an increasing function of $1/x$) dominates over $P^{(B)}_{SM}$ by far.

Fig. 2 illustrates the dependence of $P_{\textrm{\rm obs}}$ on the emergence angle $\theta_{\rm em}$ for the two choices $(x,y)=(0,0.9)$ and $(2,0.9)$, respectively, and $E_{\nu}=3\,$EeV. The green (solid) line is the total probability, while the other curves show the individual contributions $P^{(1,A)}$ and $P^{(B)}_{SM}$, as well as the SM contribution (red dashed curve). (Since the latter neglects neutrino regeneration effects, the true SM probability is higher for $\theta_{\rm em}\gtrsim 3^{\circ}$.)

⟨2⟩** The BSM process $\bm{\nu\rightarrow X_{1}\rightarrow\cdots\to X_{n}\to\tau\rightarrow{\rm\bf hadrons}}$.** Rather than a single BSM particle $X$ in the conversion, it is conceivable to have a cascade of decays or several consecutive conversions. One expects that increasing the number of steps increases the conversion probability, since the final step can then have a higher likelihood of producing a $\tau$ near the surface of the Earth. Going to a very large number of steps is theoretically unlikely, and part of the initial $E_{\nu}$ is lost in each decay. Therefore we give the numerical results up to five steps, fixing $\theta_{\rm em}=27^{\circ}$ and $\ell_{c}=5785\leavevmode\nobreak\ \mathrm{km}$, by neglecting the energy loss in every conversion, except for the energy loss of $\tau$, taking the initial $\nu$ energy to be $3\leavevmode\nobreak\ \mathrm{EeV}$. In table 1, second column, we show the maximum probabilities as a function of the number of steps in the cascade. They are achieved when the conversion rate of $\nu\to X_{1}$ is high (small $x$) and all the decay lengths of the $X_{i}$ are roughly $\ell_{i}\sim\ell_{c}/n$. The probabilities include the hadronic branching ratio $B_{h}$.

⟨3⟩** The BSM process $\bm{\nu\rightarrow X_{1}\rightarrow\cdots\to X_{n}\rightarrow{\rm\bf hadrons}}$**. The AAEs do not necessarily originate from $\tau$ decays. A BSM particle decaying in the atmosphere to hadrons could also give rise to the observed events, for example a dark photon decay $A^{\prime}\to h$. As in case 2, we take $\theta_{\rm em}=27^{\circ}$ and $E_{\nu}=3\leavevmode\nobreak\ \mathrm{EeV}$, and neglect the energy loss of BSM particles while they traverse the Earth. The maximum probabilities are shown in column 3 of table 1, where we take the branching ratio $X_{n}\to{\rm hadrons}$ as $1$. These can be approximated by an analytic formula, $P_{n}\simeq\sqrt{n/2\pi}\,D/(\ell_{c}\sin\theta_{\rm em})$, which works well if $\ell_{c}/n\gg D/\sin\theta_{\rm em}$.

In the following, we will refer to the model-dependent probability defined in this section generically as $P_{\rm obs}(\theta_{\rm em})$, considering it as a function of the emergence angle.

II.2 Limits on BSM assuming a diffuse $\nu$-flux

The predicted number of ANITA events is in general given by eq. 2, where conversion probabilities for different models are given in table 1. We use the exposure estimate obtained from the reflective events observed by the ANITA-I and III flights Fox et al. (2018),

[TABLE]

For the incoming neutrino flux, we take the limits from IceCube Aartsen et al. (2018a) and Auger Aab et al. (2015): the energy times diffuse isotropic flux is $\lesssim 6$ km*-2* sr*-1* yr*-1* at $E_{\nu}=3\,$EeV. Assuming that $\nu_{\tau}$ constitutes $1/3$ of the total diffuse flux, due to oscillations, this implies $E_{\nu_{\tau}}\,\phi(E_{\nu_{\tau}})\lesssim 2$ km*-2* sr*-1* yr*-1*. Our conclusions below do not change significantly on varying the energy in the range of [$0.1\,$EeV, $10\,$EeV].

$H_{\textrm{\tiny AN}}^{\textrm{\tiny ref}}$ does not include the conversion probability $P_{\textrm{\rm obs}}$. The latter is a function of $\theta_{\textrm{em}}$, so we average over it in order to obtain the estimated total exposure. As mentioned above we assume that $\mathrm{d}\Omega_{d}$ does not depend on $\theta_{\textrm{em}}$. Then the predicted number of events from a diffuse flux in the region visible to ANITA can be estimated as

[TABLE]

To check whether the AAEs can be compatible with the BSM scenarios discussed above we perform an extended likelihood analysis. The extended likelihood takes into account the total number of observed AAEs $n$ as well as the measured emergence angles $\theta_{\textrm{em}}^{i}$ of these events, and is defined as

[TABLE]

where

[TABLE]

is the probability distribution function to observe the emergence angle $\theta_{\textrm{em}}$, and $L(0,\theta_{\rm em})\equiv 1$. From $L$ one can compute the $p$-value

[TABLE]

for ANITA to observe at least two events with the respective emergence angles $27^{\circ}$, $35^{\circ}$. The angular integration region is defined to be such that

[TABLE]

in terms of the Poisson distribution $P(k,\mu)=\mu^{k}e^{-\mu}/k!$. In practice, the important contributions to (19) come from $n\geq 2$, and for $n\geq 3$ one can neglect the angular dependence and take $L(n)=P(n;\mu)$. For BSM model ⟨1⟩, we find that the parameters $x$ and $y$ that maximize the extended likelihood are $x\sim 0$, $y\sim 1.13$. As expected $y$ is slightly different from the one given in the previous section since we are considering two events with emergence angles $27^{\circ}$ and $35^{\circ}$, respectively. The corresponding $p$-value is $3.7\times 10^{-8}$, which implies that the BSM topology ⟨1⟩ is inconsistent with the data at the level of $5.5\,\sigma$. As shown in table 1, for the model ⟨2⟩, the discrepancies are roughly the same, $\sim 5.5\,\sigma$; for the model ⟨3⟩, the discrepancy can be reduced to the level of $4.6\,\sigma$. It turns out that increasing the number of steps does not correspond to an decrease of the $p$-value, which depends upon the total acceptance and not just on $P_{\rm obs}$.

We conclude that if a diffuse neutrino flux is assumed to be the source, BSM physics is insufficient to explain the ANITA anomalous events. The minimum incompatibility, tension at the level of 4.6 $\sigma$, is obtained for a cascade of several BSM particles $X_{i}$ in which the last one directly decays to hadrons with a branching fraction close to unity.

II.3 Limits on BSM from point sources

The IceCube diffuse flux upper limit for $\nu_{\mu}$ (that we presume to also hold for $\nu_{\tau}$), $E_{\nu}\,\phi(E_{\nu})\lesssim 2$ km*-2* sr*-1* yr*-1*, rules out both SM and BSM explanations of the AAEs. A loophole could in principle be bright neutrino point sources. If there exists a large number of point sources with an isotropic distribution, their total flux cannot exceed the diffuse flux limit. Thus for an isotropic distribution of many point sources, one reaches the same conclusion as for the diffuse flux. On the other hand, when the number of point sources $N_{\rm ps}$ is low, the single source flux limit is bounded by total diffuse flux divided by $N_{\rm ps}$; see fig. $8$ of Aartsen et al. (2019). The optimal situation for explaining the AAEs is then to have only two points sources, both in the Northern hemisphere, that give rise to the two AAEs.

For a single point source of muon neutrinos, assuming a spectrum of the form $\phi_{\rm ps}\sim E^{-2}$, the sensitivity of IceCube is $\phi_{\rm ps}\,E_{\nu}=0.13\,\textrm{km}^{-2}\,\textrm{yr}^{-1}$ at $E_{\nu}=1\,\rm{EeV}$ Aartsen et al. (2019). The IceCube sensitivity to $\nu_{\tau}$ is lower because of poorer angular resolution, but we will assume that the two fluxes are equalized by oscillations. Since the angle is fixed for point sources, the unit of flux does not include $\textrm{sr}^{-1}$, and this changes the estimated number of observed events from eq. (16) to the form

[TABLE]

where $\mathrm{d}A_{\rm ps}=(h/\sin\theta_{\textrm{em}})^{2}\mathrm{d}\Omega_{d}$ is the differential area over which ANITA observes an air shower originating from a point source at an angle $\theta_{\textrm{em}}$ and $h\cong 35\,$km is its altitude. Taking $\int\mathrm{d}A_{\rm ps}\sin\theta_{\textrm{em}}\simeq 9\times 10^{4}\leavevmode\nobreak\ \mathrm{km}^{2}$, we can estimate in the BSM model ⟨1⟩, the maximum $\mu_{ps}=1\times 10^{-5}$ for two point sources in the northern hemisphere, excluding the model at the $6.5\,\sigma$ confidence level. The reason for this stronger exclusion is the assumption of the spectrum $\phi(E_{\nu})\sim E_{\nu}^{-2}$ Aartsen et al. (2019). It leads to fewer events at $E_{\nu}\sim 1\,$EeV from anisotropic point sources than does the diffuse flux.

This suggests another possible loophole (apparently taken in ref. Collins et al. (2019)), namely to have a non-monotonic spectrum with a peak at energies $\sim 1$ EeV. Ref. Aartsen et al. (2019) computes the sensitivity and discovery potential of IceCube to such sources, as a function of $E_{\nu}$ and declination $\delta$ (see fig. 4 in that reference). For $\delta=30^{\circ}$, similar to the angles of the AAEs, the differential flux limits are $\phi_{ps}\,E_{\nu}\sim 6\times 10^{3}\,{\rm km}^{-2}\,{\rm yr}^{-1}$ at $E_{\nu}=0.5\,{\rm EeV}$ and $\phi_{ps}\,E_{\nu}\sim 1.6\times 10^{5}\,{\rm km}^{-2}\,{\rm yr}^{-1}$ at $E_{\nu}=1\,{\rm EeV}$. Taking the former flux we find that BSM model ⟨1⟩ is in tension with the ANITA data at only $1.6\,\sigma$, while with the latter, all the models including the SM are consistent with the data.

However at declination $\delta=30^{\circ}$ and $E_{\nu}\sim 1$ EeV, the Auger experiment is more sensitive than IceCube, since events from this direction are down-going or Earth-skimming Pierre Auger Collaboration (2012). A limit of $E_{\nu}^{2}\,\phi_{ps}<140\,$EeV km*-2* yr*-1* is reported for an $E_{\nu}^{-2}$ flux from a point source. For a flux concentrated at $E_{\nu}\sim 1\,$EeV, the limit is presumably weaker, but probably not greatly so since Auger’s sensitivity to differential spectra is peaked around 1 EeV. It is reasonable to assume the differential limit at $\sim 1\,$EeV is a factor of $5$ weaker than the $E_{\nu}^{2}$ flux. Taking the differential limit of $E_{\nu}\phi_{\nu}<700$ km*-2* s*-1* at $E_{\nu}\sim 1\,$EeV, this would be in $3.1\,\sigma$ tension for being able to explain the AAEs. We conclude that nonpower-law flux point sources at high declination, while perhaps not ruled out, can only be marginally consistent with the AAEs.

Finally, intense transient sources lasting $\ll$ 24 hours, coming from the northern hemisphere and peaked in the EeV range, might have been missed by both Auger and IceCube, so that we are unable exclude this as a possible explanation, either within the SM or coming from new physics. A study of this possibility is beyond the scope of the present work.

III ANITA events from decaying Dark Matter

Since very large fluxes are required to overcome the small transmission probability through the Earth, one is motivated to look for an exotic source of BSM particles that could produce AAEs, whose flux could be much larger than the experimental limits on the neutrino diffuse or pointlike fluxes. Decaying dark matter has been suggested. In one scenario, the signal is dominated by decays of DM that has been trapped in the Earth Anchordoqui et al. (2018), while in another it is DM decaying in the galactic halo that dominates Heurtier et al. (2019). We will argue that both of these are in varying degrees of tension with complementary constraints. However it is possible that decaying DM gives rise only to highly boosted BSM particles, that interact in the Earth to produce the AAEs. We point out in the following Section that this third scenario, heretofore unconsidered, can be viable.

III.1 DM decays in the Earth

DM particles $\chi$ can accumulate inside the Earth if they have a sufficiently large cross section $\sigma_{\chi N}$ for scattering on nucleons. The current upper limit from XENON1T, extrapolated to DM mass $m_{\chi}=1$ EeV, is $\sigma_{\chi N}<10^{-39}\,$cm2 Aprile et al. (2018), giving a mean free path in the Earth of $\lambda>0.3\times 10^{10}$ km, which is $5\times 10^{5}$ times greater than the radius of the Earth. The collection efficiency is expected to be $\epsilon=\pi R_{\oplus}/(2\lambda)\lesssim 3\times 10^{-6}$ (accounting for the average distance through the Earth traversed by a DM particle). The total mass of DM which the Earth can collect is

[TABLE]

where $\rho_{\textrm{\tiny DM}}^{\textrm{local}}\simeq 0.3\,\textrm{GeV}/\textrm{cm}^{3}$ is the local DM density, $v_{\textrm{rel}}\simeq 220\,\textrm{km/s}$ is the DM velocity and $t_{\varoplus}\simeq 4.5$ Gyr is the age of the Earth.

On the other hand, in order to explain the $N_{\textrm{\tiny A}}$=2 anomalous events within $t_{\textrm{\tiny A}}=$ 34 days, in a solid angle $\Omega_{\textrm{\tiny A}}\simeq 4\textrm{\,km}^{2}/(4\pi R_{\varoplus}^{2})$, and with DM lifetime $\tau_{\chi}$, one needs a total DM mass accumulation of

[TABLE]

where $P_{\textrm{\tiny A}}$ is the probability for a DM decay product to propagate through the Earth and produce a $\tau$ that initiates a hadronic shower in the lower atmosphere. We take an optimistic value of $P_{\textrm{\tiny A}}\sim 10^{-3}$ based on our previous estimates. The lifetime is constrained by IceCube searches for neutrinos from decaying DM Aartsen et al. (2018b) to be $\tau_{\chi}\gtrsim 2\times 10^{27}$ s. This is the weakest limit (corresponding to $\chi\to\nu\nu$) which in general depends upon the final state particles in the decay. The required mass is five orders of magnitude greater than the amount of DM that could collect inside the Earth, so even with an unrealistically large value for the conversion probability $P_{\textrm{\tiny A}}\sim 1$, the scenario is still ruled out.

One might ask whether this negative conclusion could be evaded by considering strongly interacting dark matter, that is stopped in the Earth before reaching underground detectors. However the larger cross sections needed have been excluded by other experiments and astrophysical constraints, as summarized in ref. Mack et al. (2007). That reference further closes any remaining loophole by showing that too much heat would be deposited in the Earth relative to the measured heat flow, if the direct detection bound is significantly violated; hence our upper bound on $\epsilon$ is robust.

III.2 DM decays in the galactic halo

Even without collecting in the Earth, decays in the DM halo of our galaxy can provide a flux of sterile BSM particles that is far bigger than the experimental limit on the diffuse neutrino flux. This could potentially overcome the small survival probability for traversing the Earth, if these BSM particles can convert to $\nu_{\tau}$, $\tau$ or some other particle that produces a hadronic shower in the atmosphere. This is the strategy of ref. Heurtier et al. (2019), in which a bosonic DM particle $\chi$ is presumed to decay into sterile right-handed neutrinos $\nu_{R}$, that have a small mixing angle $\theta\sim 0.01$ with active neutrinos. Through this mixing, $\nu_{R}$ can interact with nucleons in the Earth to produce $\nu_{\tau}$, leading to the observed AAEs.

However because of the mixing, there are also direct decays $\chi\to\nu_{R}\nu_{\tau}$, with a branching ratio of $\theta^{2}$, rendering the decays themselves detectable by IceCube, regardless of $\nu_{R}$ interactions in the Earth. Ref. Heurtier et al. (2019) explains the observed AAEs with $m_{\chi}=20\,$EeV and $\tau_{\chi}=\theta^{2}\times 10^{27}$ s. The effective lifetime for the subdominant decays $\chi\to\nu_{R}\nu_{\tau}$ is therefore simply $10^{27}$ s, independent of $\theta$. This is in tension with IceCube limits on decaying DM in the $\chi\to\nu\nu$ channel Aartsen et al. (2018b). Although the latter reference stops slightly short of $m_{\chi}=1\,$EeV, the limits have been extended to higher masses Kachelriess et al. (2018), giving $\tau_{\chi}>2\times 10^{28}$ s, at $m_{\chi}=20\,$EeV. As a consequence, the expected number of AAEs is at least an order of magnitude too low. Ref. Heurtier et al. (2019) also identifies a region of parameter space around $m_{\chi}=5\times 10^{4}\,$EeV, $\tau_{\chi}=\theta^{2}\times 10^{25}\,$s, corresponding to an effective lifetime of $10^{25}\,$s for the $\chi\to\nu_{R}\nu_{\tau}$ channel. This too is disfavored, since the limits on the lifetime from Kachelriess et al. (2018) are two orders of magnitude stronger than what one would need to obtain an order one number of AAEs.

IV A viable model: DM decays to DM

The previous example suggests that one might achieve a viable explanation if DM decays solely into a BSM particle whose interactions in the Earth can produce an exiting $\tau$, with sufficiently high probability. Here we present a phenomenological model with a heavy metastable DM particle $\Psi$ decaying to a highly boosted, light DM particle $\chi$, that can interact in the Earth to produce $\tau$ leptons.

The heavy component $\Psi$ must have a mass $m_{\Psi}\sim 1$ EeV to produce the AAEs. This is above the mass limit for thermal production of DM, but there exist several mechanisms for producing super-heavy DM, for example by gravitational particle production Kuzmin and Tkachev (1998); Chung et al. (1999a); Kannike et al. (2017) or reheating/preheating Kuzmin and Rubakov (1998); Chung et al. (1998); Greene et al. (1997); Chung et al. (1999b) at the end of inflation, freeze-in Kolb and Long (2017), dilution by entropy production Bramante and Unwin (2017), or bubble collisions during a phase transition Chung et al. (1998). We suggest a scenario where $\Psi$ may be a subdominant component of the total dark matter, while the lighter $\chi$ particle constitutes most of the DM and gets its relic density from thermal freezeout.

IV.1 Particle physics model

For a specific model, we take both $\Psi$ and $\chi$ to be fermionic, with $\chi$ coupling to lepton doublets, in particular $L_{\tau}$, and an inert Higgs doublet $\phi$. The interaction Lagrangian is

[TABLE]

There is a $Z_{4}$ charge carried by the new particles that guarantees the stability of $\chi$. The scale of the dimension-6 operator is related to the lifetime of $\Psi$ by

[TABLE]

(If desired, a generalization of the model that brings $\Lambda$ below the Planck scale is to introduce a $Z_{2n}$ discrete symmetry under which

[TABLE]

so that $\Psi\to(2n-1)\chi$ through an operator of dimension $3n$, $\bar{\psi}\chi(\bar{\chi}^{c}\chi)^{n-1}$, with $n\geq 3$.)

For superheavy DM with $m_{\Psi}\sim 3\,$EeV, the decay product $\chi$s are highly energetic, and can interact with nucleons in the Earth by the diagram shown in fig. 3. A novel feature is that the internal $L_{\tau}$ can go on shell, so the cross section is logarithmically sensitive to the width of $\nu_{\tau}$, $\Gamma_{\nu}$. Since the decay rate of $\nu_{\tau}$ is negligible, this width is dominated by the interaction rate of $\nu_{\tau}$ in the Earth,

[TABLE]

where $n_{N}\cong 3.3\times 10^{24}/$cm3 is the density of nucleons and the energy dependence of $\sigma_{\nu N}$, valid around $E_{\nu}\sim$EeV, is taken from ref. Alvarez-Muñiz et al. (2018).

At high energies, the scattering is dominated by low-virtuality $W$ or $Z$ exchange. We find that the charged-current cross section is (see appendix A for details)

[TABLE]

where $\hat{s}=xs$ is the parton-level invariant, with momentum fraction $x$. The function $\ln(\sqrt{\hat{s}}/\Gamma_{\nu})$ denotes the leading behavior at large values of the argument, but we use the more exact expression (41) for the following estimates.

Averaging (28) over the parton distribution functions of the nucleons gives a result that depends upon the inert Higgs doublet mass $m_{\phi}$, because the cross section is dominated by small $x$, whose minimum value is $x_{\rm min}=m_{\phi}^{2}/s$. The enhancement factor ${\cal E}$ takes the place of $\ln(\sqrt{\hat{s}}/\Gamma_{\nu})$,

[TABLE]

and ${\cal E}$ varies from $12,000$ for $m_{\phi}=100\,{\rm GeV}$ to $3,000$ for $m_{\phi}=1\,{\rm TeV}$, as shown in fig. 4. We find that it is accurately described by the analytic fit

[TABLE]

IV.2 ANITA anomalous events

Our model can be considered to be of the type ⟨1⟩ in our classification of section II.1, in the limit where $x\to 0$ (see eq. (14)), which indicates that the BSM particle $X=\Psi$ is already present upon entry of the Earth. The scattering process actually produces two sources of $\tau$: first from the primary vertex giving the doublet $L_{\tau}$ shown in fig. 3, and second through the fast decays of the produced $\phi\to\chi L_{\tau}$. In section II.1 we showed that the probability for these $\tau$s to be observed is maximized when $y=\ell_{x}/R_{\oplus}\cong 1$. In the remainder we will assume this condition is satisfied, to partly constrain the parameter space of the model.

The cross section (29) leads to a scattering length in the Earth of

[TABLE]

taking the density of nucleons to be $n_{N}=3.3\times 10^{24}$/cm3. Imposing $y=1$ fixes $y_{\tau}$ as a function of $m_{\phi}$ via the relation (30). We achieve lengths of order the Earth radius with reasonable values of $y_{\tau}\lesssim 1.2$, as shown in fig. 5.

The procedure of section II.1, with $x=0$ and $y=1$ in process ⟨1⟩ gives a probability of $P\sim 5\times 10^{-4}$ at large $\theta_{\rm em}$, like for the AAEs. Then using the likelihood method of section II.2, we find that the best-fit flux of $\chi$ particles is $\Phi_{\chi}=2230\,$km*-2* yr*-1* sr*-1*. For a rough estimate, we can compute the flux coming from DM decays by assuming that there is a constant mean density of $\chi$, $n_{\chi}=3n_{\Psi}(\tau_{u}/\tau_{\Psi})$ from the decays, where $n_{\Psi}=\rho_{\Psi}/m_{\Psi}=f_{\Psi}\Omega_{\rm cdm}\rho_{\rm crit}/m_{\Psi}$ and $\tau_{u}$ is the age of the universe. Here $f_{\Psi}\ll 1$ is the abundance relative to that if $\Psi$ constituted all of the DM.

The isotropic flux is then

[TABLE]

which determines the $\Psi$ lifetime when equated with the best-fit flux:

[TABLE]

A more quantitative estimate using the $J$-factor for decays in the galactic halo gives a similar estimate. One needs a very high scale $\Lambda\sim f_{\Psi}^{1/4}\times 10^{24}$ GeV in the dimension-6 operator, that could be lowered by taking a larger value of $n$ in the generalized version of the model. We do not concern ourselves here with trying to build a UV complete model, but instead emphasize that the general framework may be promising for understanding the AAEs via new physics.

IV.3 Relic density and other constraints

It is interesting that the same interaction that induces the scattering of $\chi$ in the Earth can also determine its thermal abundance through the annihilations $\chi\bar{\chi}\to L_{\tau}\bar{L}_{\tau}$. The cross section at threshold is

[TABLE]

As a rough estimate, we equate this to the nominal level $(\sigma v)_{0}\cong 2\times 10^{-26}\mathrm{cm}^{3}/\mathrm{s}$ that leads to the observed relic density. In conjunction with the determined value of $y_{\tau}$ for fitting the AAEs, this imposes a relation between $m_{\chi}$ and $m_{\phi}$ that is plotted as the solid curve in fig. 5. It is approximately fit by $m_{\chi}\cong 0.3\,m_{\phi}-22\,$GeV, consistent with the requirement that $m_{\phi}>m_{\chi}$ so that $\phi\to\chi L_{\tau}$ decays occur.

The inert doublet $\phi$ could be pair-produced though electroweak interactions at the LHC, leading to $\tau$ pairs with missing energy. This is the same signature as for stau pairs in the MSSM that decay to $\tau$s and neutralinos, $\tilde{\tau}\to\tilde{\chi}^{0}\tau$. ATLAS has searched for this signal in the 8 TeV data Aad et al. (2016) and more recently CMS has searched in the 13 TeV run Sirunyan et al. (2018). So far the sensitivity does not reach the expected production cross section, though in a future run at the HL-LHC it is projected that $m_{\phi}$ (in the guise of $\tilde{\tau}$) will be excluded below 650 GeV Collaboration (2018).

The annihilation channel $\chi\bar{\chi}\to\tau^{+}\tau^{-}$ is constrained by Fermi observations of gamma rays from dwarf spheroidal galaxies Ackermann et al. (2015). These exclude $m_{\chi}\lesssim 80$ GeV, as indicated on fig. 5.

V Conclusions

The origin of the two anomalous events observed by the ANITA balloon experiment remains unexplained. In this paper we explored if and how these AAEs might be explained by new physics, within the assumption that BSM states transverse the Earth and converts to $\tau$ or hadrons in the atmosphere. Our conclusions are:

a) A diffuse flux of $\nu_{\tau}$ cannot explain the AAEs, and it is disfavored by $\sim 5\sigma$ by considering single conversion and cascade decay models. For point sources, the tension is reduced to $\sim 3\sigma$, since the flux limit is weaker.

b) DM decaying inside the Earth cannot account for the AAEs because not enough DM could have been accumulated inside the Earth during its history.

c) A large flux of “sterile” BSM particles whose interactions in the Earth can produce a $\tau$ near the Antarctic surface could in principle explain the AAEs. Such a flux of sterile particles could originate from the decay of DM in the galactic halo. One must ensure however that the DM decay does not also induce a flux of non-sterile particles, such as active neutrinos, in excess of the limits from IceCube. As an illustration of a viable model, we discuss the case of a metastable EeV-scale dark matter particle that decays exclusively to a lighter dark matter particle $\chi$. The latter can in turn scatter within the Earth to produce a $\tau$ air shower, consistently with all constraints. The relic density of the stable DM can be explained by the same interaction that induces the air showers, which can also lead to events at the LHC resembling supersymmetric $\tilde{\tau}$ production and decay.

Acknowledgements.

The authors would like to thank Markus Ahlers, Kfir Blum, Jack Collins, Derek Fox, Peter Gorham, Joachim Kopp, Kohta Murase, Matteo Puel, Sarah Shandera, Ian Shoemaker and Donglian Xu for useful discussions. CG and WX are supported by the European Research Council grant NEO-NAT. JC is supported by the Natural Sciences and Engineering Research Council of Canada.

Appendix A DM-nucleon cross section

Here we give some details of the computation of the charged-current scattering process shown in fig. 3. The spin-averaged squared matrix element takes the simple form

[TABLE]

in the limit of massless particles, where

[TABLE]

for the diagram with virtual $\nu$. The parton-level cross section in terms of the three-body phase space can be written as

[TABLE]

where we choose coordinates such that (in the center-of-mass frame)

[TABLE]

with $E_{1}=E_{2}=\sqrt{\hat{s}}/2$ and

[TABLE]

being fixed by energy-momentum conservation, since $E_{5}=\sqrt{\hat{s}}-E_{3}-E_{4}=|\vec{p}_{5}|=|\vec{p}_{3}+\vec{p}_{4}|$. Here $\hat{s}$ is the partonic center-of-mass energy squared,

[TABLE]

where $x$ is the quark momentum fraction and $m_{N}$ the nucleon mass. Only for $x$ close to $x_{\rm min}=m_{\phi}^{2}/s$ does it matter that the external particles are massive; otherwise the masses can be neglected since they are $\ll 1$ EeV, the energy scale of interest.

The momentum squared in the $\nu_{\tau}$ propagator is $2p_{1}\cdot p_{3}$ $=2E_{1}E_{3}(1-\cos\alpha)$, so the integral over $\cos\alpha$ is dominated by $\cos\alpha\cong 1$ at high energies. One finds that

[TABLE]

where $\delta=(\sqrt{\hat{s}}/2-E_{3})/\Gamma_{\nu}$. Since the most sensitive $\cos\alpha$ dependence is in this factor, we take $\cos\alpha\to 1$ in the remaining part of $|{\cal M}|^{2}$, which depends upon $E_{4}$. The integral over $E_{4}$ can be done analytically, and in the limit $E_{3}\gg m_{W}$ it gives

[TABLE]

The remaining integral over $E_{3}$ can also be done analytically,

[TABLE]

where $\hat{\Gamma}=\Gamma_{\nu}/\sqrt{\hat{s}}$.

Finally one must integrate over the parton distribution functions,

[TABLE]

where we evaluated them at scale $Q=m_{W}$ in accordance with the low virtuality of the internal $W$ boson indicated by (40).