A New Probe of Relic Neutrino Clustering using Cosmogenic Neutrinos

We propose a new probe of cosmic relic neutrinos (C$\nu$B) using their resonant scattering against cosmogenic neutrinos. Depending on the lightest neutrino mass and the energy spectrum of the cosmogenic neutrino flux, a Standard Model vector meson (such as a hadronic $\rho$) resonance can be produced via $\nu\bar{\nu}$ annihilation. This leads to a distinct absorption feature in the cosmogenic neutrino flux at an energy solely determined by the meson mass and the neutrino mass, apart from redshift. By numerical coincidence, the position of the $\rho$-resonance overlaps with the originally predicted peak of the Greisen-Zatsepin-Kuzmin (GZK) neutrino flux, which offers an enhanced absorption effect at higher redshifts. We show that this absorption feature in the GZK neutrino flux may be observable in future radio-based neutrino observatories, such as IceCube-Gen2 radio, provided there exists a large overdensity in the C$\nu$B distribution. This therefore provides a new probe of C$\nu$B clustering at large redshifts, complementary to the laboratory probes (such as KATRIN) at zero redshift.


INTRODUCTION
The detection of the cosmic neutrino background (CνB) is an extremely important problem in fundamental physics [1]. So far, there has been only indirect evidence for CνB from precise measurements of the primordial elemental abundances in big bang nucleosynthesis (BBN) [2], cosmic microwave background (CMB) [3], and large-scale structure (LSS) [4,5]. However, the direct detection of CνB remains an open challenge, which is often dubbed as the "Holy Grail" of neutrino physics. A direct detection of CνB will provide a strong validation of the Hot Big Bang cosmological model, and moreover, will provide a window into the first second of creation before the pre-recombination age. There exist several proposals for the direct detection of CνB [6][7][8][9], but none of them are expected to be experimentally feasible in the foreseeable future.
Recall that in the standard cosmological picture, the effective temperature of the CνB today is inherently connected to the CMB temperature: From neutrino oscillation data [10], we know two masssquared differences, while the absolute value of the neutrino mass is still unknown. Assuming that the lightest neutrino mass m lightest ≳ T ν,0 , the CνB today can be thought of as a non-relativistic gas of fermions with the number density of n ν,0 = 3 4 ζ(3) π 2 g ν T 3 ν,0 ≃ 56 cm −3 (2) per neutrino flavor (and similarly for antineutrinos).
Here g ν = 1 is the number of degrees of freedom for each neutrino. Although the CνB has the largest flux among all natural and man-made neutrino sources [11], it is their small kinetic energy that makes any direct probe of CνB extremely challenging. An interesting idea for CνB detection is via its scattering off ultra high-energy (UHE) neutrinos or cosmic rays. In particular, it was pointed out long ago [12] (see also Refs. [13][14][15][16]) that resonant annihilation of CνB with UHE neutrinos can produce a Standard Model (SM) Z boson on-shell (the so-called 'Z-burst') for UHE neutrino energies of the order of where z is the redshift at which the νν annihilation occurs in the sky. The cross-section for this process is pretty large: ⟨σ ann.
However, given the cosmological [3] and laboratory [17] constraints on the neutrino mass which is required to be ≲ O(eV), the Zburst energy (3) is clearly beyond the so-called GZK cutoff [18,19]. As a reminder, although the origin of UHE cosmic rays (UHECR) are unknown [20], it is widely believed that the highest-energy neutrinos originate in the scattering of protons from the UHECR sources (with energy E p ≥ 5.5 × 10 19 eV) off the CMB photons, which resonantly produces ∆ + baryons (with m ∆ ≈ 1232 MeV) that subsequently decay to pions and nucleons, with neutrinos being one of the final decay products -these are the so-called cosmogenic or GZK neutrinos [21][22][23]. Therefore, it is unlikely for any absorption features in the UHE neutrino spectrum due to the Z-burst to be observed, unless some super-GZK UHECR sources are established [24][25][26][27][28]. It is worth mentioning here that the resonance energy can be lowered to sub-GZK scale FIG. 1. Feynman diagram for the resonant production of neutral vector mesons (e.g., ρ 0 ) in νν annihilation, followed by their hadronic decay. and even down to the PeV scale currently being probed at neutrino telescopes such as IceCube, if the resonant state is sufficiently light. This has been discussed in the literature in the context of non-standard neutrino interactions with light mediators [29][30][31][32][33][34][35][36][37][38][39][40][41][42][43], however the Standard Model (SM) itself provides light vector mediators below the GeV scale that can mediate a resonance and therefore supplies a mechanism for sub-GZK absorption.

VECTOR-MESON RESONANCES
Within the SM, even at much lower energies, there is the possibility of a vector (or axial-vector) meson resonance in νν annihilation. Recall that in the SM, neutrinos and left-handed components of charged leptons are doublets of SU (2) L : ψ i = (ν i , ℓ − i ) T , where i = 1, 2, 3 is the family index. After electroweak symmetry breaking, the relevant piece of the SM Lagrangian reads [44] where A µ and Z µ are the photon and the Z-boson fields respectively, g is the SU (2) L gauge coupling, θ w is the weak mixing angle, e = g sin θ w is the positron electric charge, Q i is the fermion charge relative to the positron, and g V,A are the vector and axial-vector couplings for the Z-boson. Let us first recap the familiar case of e + e − collisions. For the center-of-mass energy √ s ≪ m Z , the e + e − → ff (f being any SM fermion) cross-section falls as 1/s for dimensional reasons (as long as √ s ≫ m f ). Therefore, the ratio R of the cross-sections into quark-antiquark and µ + µ − pairs is a constant that depends on the underlying number of quark degrees of freedom (see e.g. Fig.  52.2 in Ref. [44]). One can readily see that for s ≪ m 2 Z , the e + e − cross-section is dominated by the vector meson resonances with J P C = 1 −− (as predicted long ago in Refs. [45,46]), such as ρ, ω , ϕ, J/Ψ, Υ etc. -an experimentally well-established fact.
In much the same way, when it comes to νν collisions, the weak current contains an admixture of both vector and axial couplings [cf. Eq. (4)]. The weak vector current will produce the J P C = 1 −− resonances, such as the ρ 0 meson as shown in Fig. 1, while the axial-vector current will produce J P C = 1 ++ resonances, such as a 1 (1260), f 1 (1285), f 1 (1420), f 1 (1510), χ c1 (3872), χ b1 , etc. [47][48][49]. The corresponding cross-sections for both vector and axial-vector meson resonances can be calculated using the Breit-Wigner resonance formula [44]: where BR stands for the branching ratio of the resonance particle X to neutrinos, m X and Γ X being its mass and total width respectively. The squared center-of-mass energy needed for a given resonance is s = 2m ν E ν , where E ν is the incident UHE neutrino energy: It is a numerical coincidence that the cosmogenic neutrino flux typically peaks around 10 18 eV (depending on the exact value of the galactic-extragalactic crossover between the "second knee" and the "ankle" in the CR spectrum) [50], whereas the lightest vector meson (ρ 0 ) resonance (with m ρ ≈ 775 MeV) requires a neutrino energy Therefore, we will mostly focus on the ρ-resonance in this work. Note that although the next lightest vector meson (ω) mass (m ω = 782.7 MeV) is close to the ρ-mass, its narrow width (Γ ω = 8.7 MeV) makes a small difference to the broad resonance feature induced by the large ρ-width (Γ ρ ≈ 150 MeV). The ρ production cross-section can then be estimated from Eq. (5), where the corresponding BR can be obtained from the partial width where f ρ ≈ 216 MeV is the ρ-meson decay constant [51] (see also Refs. [52,53]). We find that the cross-section at resonance is O(10 −38 cm 2 ) [49], much smaller than the Z-burst cross section, primarily due to the hugely suppressed BR(ρ → ν i ν i ) ∼ 10 −13 . Nevertheless, for charged current resonances e.g.ν e e − → ρ − with a similar cross-section, O(100) events from the hadronic decay of the ρ-meson can occur within the volume of IceCube for a 10-year exposure; whether these events can be detected is an open question [54]. For the neutral current resonances like νν → ρ 0 under discussion, the corresponding event rates from ρ-decay are entirely negligible on Earth due to the low density of the CνB, i.e. n ν,0 ≪ n e . Instead, we propose to look for absorption features in the GZK neutrino spectrum itself due to the resonant production of the ρ-meson. As we show below, an observable effect requires an extremely large CνB overdensities at higher redshifts (i.e. n ν,z /n ν,0 ∼ 10 11 ) which we will assume to be the case. Possible origins of this overdensity will be discussed later.

ATTENUATION OF GZK NEUTRINOS
As we discuss below, for an observable GZK absorption feature to appear in neutrino telescope data, it is necessary that GZK neutrinos have a high probability of encountering an overdense neutrino cloud. The most natural way for this to occur is if neutrino clustering is spatially correlated with UHECR production such that the GZK neutrinos are produced directly within the ultradense neutrino cloud. For an observable signature, it turns out to be necessary for UHECR production to be dominated by only a few sites in the Universe, and for a substantial fraction of the total neutrino population in one Hubble patch to be concentrated into clouds around these production sites. One may ask whether or not such a model is consistent with existing CR data. Due to the loss of directional information for charged CRs, due to e.g. deflection by galactic [55] or intergalactic [56] magnetic fields, we argue that existing UHECR observations cannot rule out such a scenario. In fact, UHECR anisotropies have already been observed [57,58], and future measurements by large, ground-based experiments, like LHAASO [59], SWGO [60] and IceCube-Gen2 [61] can in principle identify the origin of these anisotropies. The scenario we propose can therefore be tested both by searching for the resonant absorption we discuss below and by the necessary anisotropy in the GZK neutrino flux.
Let us imagine, at some redshift z, the source of GZK neutrinos is surrounded by a cluster of CνB with density n ν = ξn ν,0 (1 + z) 3 ≡ n ν,z (1 + z) 3 , where n ν,0 is given by Eq. (2) and ξ quantifies the overdensity. The GZK neutrino flux will then be attenuated due to the resonant production of ρ-mesons, which promptly decay into pions and subsequently into charged leptons and neutrinos but typically with much lower energy than the original neutrino energy. The attenuation factor is given by where λ = 1/σn ν is the mean free path, σ is the νν cross section given by Eq. (5), L = (1/H 0 )(qξ) −1/3 is the average cluster length traversed by the GZK neutrinos (assuming a spherical cluster), q is the total number of clusters in the Universe (assuming that each cluster is of the same size, and that the total number of relic neutrinos in the Universe is conserved), and H 0 the Hubble constant at z = 0. The maximum attenuation at the resonant energy [cf. Eq. (7)] as a function of the redshift and overdensity is shown in Fig. 2. It is clear that the attenuation effect is more pronounced at larger redshifts and/or larger overdensities simply because those two parameters directly impact the total CνB density at a given location.
Arbitrarily large redshifts are not realistic because the GZK neutrino flux, which is believed to originate from the progenitor cosmic-ray flux, should be correlated with the star-formation rate which peaks between z = 1 and 2, and is essentially non-existent after z = 4 [62]. Similarly, large overdensities are difficult to realize in concrete models, as we will discuss below, however they are not beyond the realm of possibilities for certain beyond the SM (BSM) scenarios. As a concrete example, we choose z = 2, ξ = 10 11 and q = 1 in Eq. (9) in order to illustrate the effect of resonant absorption on the GZK neutrino flux, as shown in Fig. 3. Here the thick solid black curve is the unattenuated diffuse flux prediction from Ref. [50] which is marginally consistent with present IceCube [63] and Pierre Auger [64] upper limits; the green band corresponds to the CR models considered in Ref. [50] with minimal and maximal energy density at the 99% confidence level (CL). Note that this flux assumes primarily proton-dominated sources. Heavier elements like iron in the source composition could in principle increase the flux uncertainty by an order of magnitude or even more [65][66][67][68][69][70][71][72][73]. However, it is not clear which source composition is preferred by data, because the Pierre Auger data [64] with lower UHECR flux prefers heavier composition, whereas the TA data [74] prefers a higher flux. For a review, see e.g., Ref. [75].
In Fig. 3, the thick solid blue, purple and orange curves show the attenuated flux due to resonant ρ-production for three benchmark values of the lightest neutrino mass m 1 = 0.01 eV, 0.05 eV and 0.5 eV, respectively. Normal mass ordering is assumed, as seems to be moderately preferred by the latest oscillation data from T2K [78,79] and NOvA [80,81]. For m 1 = 0.01 eV, we can see two dips corresponding to the lighter mass eigenstates (m 1 , m 2 ) combined and the heavy one (m 3 ), respectively. This is because the gauge couplings are flavor diagonal, and therefore, each mass eigenstate contributes separately to a single dip (unlike the case of secret neutrino interactions which can be flavor off-diagonal, and all mass eigenstates undergo absorption in all dips [42,82]). But as the m 1 value increases, the three neutrino mass eigenstates become quasi-degenerate and the two peaks merge into one, as shown for the other benchmark points.

A PROBE OF EARLY CνB OVERDENSITY
We now translate the result in Fig. 3 into an experimental sensitivity plot for the CνB overdensity as a function of the lightest neutrino mass, as illustrated in Fig. 4 for IceCube-Gen2 radio. Here we calculate the number of expected GZK neutrino events at IceCube-Gen2 radio with the unattenuated GZK flux (black curve in Fig. 3) and compare it with the attenuated number of events for different values of the lightest neutrino mass assuming normal mass ordering. We have further assumed that forthcoming measurements of the GZK neutrino flux will converge to the presently unconstrained theoretical prediction shown in Fig. 3 and this is a realistic goal, provided the flux is primarily from proton-dominated CR sources. We perform a single-bin analysis by constructing an optimal bin size around the resonance energy for given values of lightest neutrino mass and redshift. If we take the full available energy range E ν ∈ [10 7 , 10 11 ] GeV (for which the effective area is known) as our bin size, then we are simply performing a counting experiment in which a precise determination of the local flux attenuation in a narrow energy band is hindered. In the opposite limit, if we choose a bin width too narrow to capture the whole energy range corresponding to the width of the ρ resonance, the sensitivity would also be poor. Therefore, we construct a bin width that matches the energy range corresponding to four widths of the ρ meson (two on either side of the resonance energy where the bin is centered). We have checked that our results are stable against choosing slightly different bin widths, as well as against experimental energy resolution and smearing effects, which are expected to be at the level of O(10%) (or better) for shower events -so at least a factor of two better than the width-to-mass ratio of the resonance Γ ρ /m ρ ≃ 19%.
The number of unattenuated (N wo ) and attenuated (N w ) events is given by where T is the exposure time (taken to be 10 years here), Ω = 4π is the solid angle of coverage, A eff is the IceCube-Gen2 radio effective area, Φ is the unattenuated GZK neutrino flux (cf. the solid black curve in Fig. 3), and R is the attenuation factor given by Eq. (9) (for the unattenuated case, R = 1). Then we compute the χ 2 using the log likelihood method: The 90% CL sensitivities (χ 2 = 2.71 for one degree of freedom) for the IceCube-Gen2 radio with T = 10 years of exposure are shown in Fig. 4 for two representative cases: (i) A single cluster of overdense neutrinos at z = 2 (solid red curve); and (ii) A stochastic distribution of 10 neutrino clusters at z = 1 (dashed red curve). For scenario (i), the sensitivity is clearly stronger because the cluster size is bigger [for q = 1, the cluster length L given below Eq. (9) is of order O(Mpc) for ξ = 10 11 ]; therefore, the likelihood for resonant absorption increases. The bigger the number of clusters, the smaller their size (as we are restricted by the total number of CνB neutrinos in the Universe); therefore, we need a larger overdensity to see the absorption effect, as illustrated by the scenario (ii). The sensitivity in both cases is strongest for the lightest neutrino mass around 0.05 eV because the corresponding resonant energy [cf. Eq. (7)] coincides with the peak of the predicted unattenuated GZK flux (see the purple benchmark in Fig. 3). The flux goes down rapidly at higher energies which correspond to lower neutrino mass for resonance; therefore, the sensitivity curve goes up to the left. On the other hand, the IceCube-Gen2 radio effective area drops significantly for smaller energies, in particular below E ν ∼ 10 17 eV; for these energies, resonance is achieved for larger neutrino masses and this is why the sensitivity in Fig. 4 worsens for increasing neutrino mass. The solid purple curve in Fig. 4 is the current KA-TRIN upper limit on the overdensity: ξ < 1.1 × 10 11 at 95% CL, derived using the possibility of CνB capture on tritium nuclei [93]. It is important to point out that the KATRIN limit is on the local overdensity, i.e. at z = 0, whereas the IceCube-Gen2 radio limit derived here provides a complementary probe at higher redshifts. The vertical brown line in Fig. 4 is the KATRIN upper limit on the electron antineutrino mass: mν e < 0.8 eV at 90% CL [17], which is equivalent to an upper limit on each of the mass eigenvalues in the quasi-degenerate limit using the relation m 2 νe = i |U ei | 2 m 2 i (U being the PMNS mixing matrix [44], assumed here to be unitary). Again, this bound is strictly applicable only at z = 0. Note that neutrinoless double beta decay experiments have also imposed an upper limit on the effective neutrino mass m ββ = | i U 2 ei m i | which translates into a bound on m 1 ≲ 0.2 − 0.6 eV (depending on the nuclear matrix element used) [94,95]; however, if neutrinos are Dirac particles, this bound does not apply, and therefore, is not included in Fig. 4.
Similarly, the vertical gray line in Fig. 4 is the 95% CL Planck upper limit derived from the most stringent cosmological constraint on the sum of neutrino masses [3]. However, it is important to keep in mind that the this bound on m i strongly depends on the combination of the cosmological datasets used, and varies from 0.12 eV to 0.60 eV (95% CL) [3]. More importantly, the CMB bound is valid at a very high redshift of z ≈ 1100, and there IceCube-Gen2 radio sensitivity at 90% CL to CνB overdensity as a function of the lightest neutrino mass. The red dashed and solid curves correspond to (z = 1, q = 10) and (z = 2, q = 1) respectively. The purple line is the current 95% CL KATRIN upper limit on local overdensity [93]. The vertical brown line is the 90% CL KATRIN upper limit on the neutrino mass [17]. The vertical gray line is the 95% CL Planck upper limit [3].
Since neutrinos are fermions, they cannot be clustered to arbitrarily high densities in a stable configuration due to the Pauli exclusion principle [113]. Under the assumption that the relic neutrinos behave like an ideal Fermi gas, the maximum possible clustering depends on their Fermi energy E F : n max ν,0 = (2m ν E F ) 3/2 /(3π 2 ).
Using semi-analytic arguments which suggest an upper bound on the Fermi momentum p F ≲ 0.9 m ν , Ref. [114] obtained a maximum CνB density of n max ν,0 ≃ 1.5 × 10 9 cm −3 (m ν /0.1 eV) 3 for neutrino bound states (in presence of Yukawa interactions of neutrinos with a new light scalar). This theoretical bound is more stringent than the limits/sensitivities shown in Fig. 4. Let us note, however, that since neutrino clustering at the scales we consider does require BSM physics the phase space considerations are somewhat model dependent. Redshift dependent BSM effects, as well as the possible impact of non-standard cosmology on the evolution of the total neutrino energy density, or other exotic non-standard neutrino interactions could all conceivably modify the naive Fermi-gas phase space limit discussed above; we therefore, do not include phase space considerations in Fig. 4.

DISCUSSION AND CONCLUSION
In principle, the resonant absorption mechanism proposed here is a purely SM phenomenon and does not require any BSM physics per se. However, as shown in Fig. 4, the usefulness of this method in probing the allowed parameter space crucially depends on the existence of a large overdensity of CνB around or along the line of sight of the GZK neutrino source(s). As mentioned above, this is possible, which however most likely depends on new neutrino interactions with BSM fields. It is interesting to note that some of the proposals for evading the neutrino mass bounds also predict significant CνB clustering [105,106,115]. However, it still remains to be seen if O(10 10 ) overdensities as needed in Fig. 4 can be achieved in practice. One might wonder if gravitational clustering around the GZK source could be of any help; but by itself, it can only provide overdensities up to 10 3 [116] (gravitational clustering at z ≳ 0.5 may be much larger than gravitational clustering around Earth at z = 0). Another possibility is early neutrino decoupling; since the neutrino number density goes as T 3 [cf. Eq. (2)], an earlier decoupling temperature by even a factor two in a non-standard cosmology could gain us an order of magnitude in the number density. A combination of several of these mechanisms may be necessary to produce the large overdensities that are experimentally accessible at near-term experiments.
To conclude, we have proposed a new probe of relic neutrino overdensity which extends to higher redshifts. While the relic neutrino clustering will necessarily require BSM physics, the absorption mechanism relies exclusively on SM physics, namely, the scattering of relic neutrinos with cosmogenic neutrinos. Hadronic vector-resonances enhance the cross section and allow for neutrino telescopes to compete with bounds from KA-TRIN while simultaneously probing large redshifts. For large local neutrino densities that are correlated with the production sites of UHE cosmic rays, future radiotechnology-based neutrino observatories like IceCube-Gen2 radio can detect this "dip" feature. This provides a new probe of non-standard cosmologies beyond the current laboratory and cosmological bounds on neutrino masses.

ACKNOWLEDGMENTS
We thank Carlos Argüelles, Diego Aristizabal Sierra, Peter Denton, Abhish Dev, Francis Halzen, Shunsaku Horiuchi, Pedro Machado, Kohta Murase, Shmuel Nussinov, Jack Shergold and Ian Shoemaker for useful comments and discussions. In addition, we are grateful to Steffen Hallmann and Brian Clark for providing us the IceCube-Gen2 radio effective area. We are also grateful to Ibragim Alikhanov for pointing out a missing factor of 6 in Eq. (5) of v1, which now leads to slightly better results in Figs. 3 and 4. BD and RP acknowledge the Fermilab theory group for local hospitality during the completion of this work. RP also thanks the Weizmann Institute of Science for their hospitality. This work was partly performed at KITP, which is supported in part by the National Science Foundation under Grant No. NSF PHY-1748958, and at the Aspen Center for Physics, which is supported by the National Science Foundation