Atmospheric resonant production for light dark sectors

Cosmic ray atmospheric showers provide an effective environment for the production of MeV-scale dark sector particles. We show that, when available, the resonant annihilation of positrons from the shower on atmospheric electrons is the dominant production mechanism by more than an order of magnitude. We provide a quantitative example based on dark photon production and update existing constraints on a corresponding light dark matter model from kilotons neutrino experiments and xenon-based direct detection experiments.


I. INTRODUCTION
The idea of using the abundant cosmic ray (CR) flux impinging on earth's atmosphere as an "atmospheric collider" has a rich history dating back to the discovery of muons. In recent years, this flux has been widely leveraged to search for Feebly interacting particles (FIPs) with mass ranging from the MeV to the GeV scale whose interaction with matter is suppressed enough that they may have escaped prior detection [1][2][3][4]. The case where this FIP acts as a mediator between the Standard Model (SM) particles and a light, sub-GeV, dark matter candidate is particularly compelling. Indeed, while such construction preserves most of the elegant features of the vanilla "WIMP" dark matter, the small dark matter mass coupled with the low velocity of the galactic dark matter halo dramatically weakens existing direct detection searches. If this light dark matter (LDM) is assumed to compose all or a large fraction of the total dark matter density, strong constrains can nonetheless be obtained by estimating the CR scattering on dark matter particles. This leads to a secondary, relativistic flux which can eventually leave a sizeable recoil signature, either in direct detection experiments [5][6][7][8][9][10][11][12][13][14] or in neutrino experiments [15][16][17][18]. Alternatively, one can specified a proper model for a stable light particle which may not constitute the actual relic dark matter or only a small fraction. Direct production of such light states from CR interactions in the atmosphere is however possible, with a subsequent detection via scattering signatures in direct detection or neutrino detectors [19][20][21][22][23][24].
In this work, we build on the second approach by pointing out a new production mechanism relying on the annihilation of CR shower-induced positrons on atmospheric electrons. We illustrate this mechanism for models of LDM based on the exchanged of a new massive vector mediator. For definiteness, we mostly focus on a massive vector sharing the same interaction pattern as a photon (thus called dark photon) and a dark matter model following an "inelastic dark matter" structure as advocated in [25,26]. Our result could however equally apply to scalar or Majorana dark matter up to order one factors. We will further translate this improved description of the production rates by estimating the limits and projections from electron scattering in SuperK [27] and HyperK [28], and from coherent nuclear scattering in xenon-based dark matter experiments [29][30][31].

A. Cosmic ray showers description
Cosmic ray showers originate primarily from high energy protons impinging on the high atmosphere. The resulting showers are therefore hadronic at first, with electromagnetic components appearing afterwards mostly as byproducts of the decay of π 0 mesons. We will focus on dark sector states that are either stable, or with a decay length significantly longer than the typical length of the shower. Since the cosmic ray flux is nearly isotropic, the longitudinal and transverse development of the shower can then be integrated out in favour of an effective energy-dependent "track-length flux" [32,33]: where dΦ CR dE0 is the differential CR flux , n e is the differential number density of e ± in the shower, function of the shower depth parameter which parametrises the longitudinal shower development, and ρ (resp. ρ 0 ) is the atmospheric density at shower depth (resp. at ground level). 1 The ratio in densities is used to normalised the track length flux to a ground level "target" atmosphere.
This quantity can be calculated once and for all and then used to estimate the dark sector differential fluxes dΦ dEχ for a variety of new physics models following where dσ dEχ is the relevant production differential crosssection from the interaction between the shower particles and the atmosphere (it hence typically includes a factor of the atomic number Z or Z 2 depending on the relevant process).
In practice, we obtain this track length flux from two different approaches. First a semi-analytical way, based on the method advocated in [33]. The differential number density of secondary neutral mesons n M0 (E) from a pN collision is obtained in this case from the QGSPJETII software [34] at energy below 10 GeV and by EPOS-LHC [35] above, as packaged in CRMC [36]. This energy-dependent spectrum is then convoluted with the cosmic ray energy flux as parametrised in [37]. The subsequent electromagnetic showers are then described analytically following the approach of Rossi and Griesen [38] as reviewed and expanded in [39]. Second, we use a purely numerical approach by relying on the software CORSIKA [40] to simulate cosmic ray showers for incoming protons with kinetic energy between 2.5 GeV and 15 TeV (using EPOS for the hadronic interaction). The track length is obtained directly from the simulated tracks of photons and e ± , and the meson differential energy distributions collected for both the π 0 and η mesons. 2 We average over azimuthal angles for the shower proton progenitors.
We use the first approach for cross-checking purposes and rely on the fully numerical second approach in our final limits and projections. We present the resulting meson distribution in Fig. 1 and the track length fluxes in Fig. 2, with the corresponding datasets available at a Zenodo repository. 3 The two approaches present a good agreement, with as expected larger fluxes at low energies for the fully numerical approach. This reflects the fact that it further includes the secondary mesons neglected in the semi-analytical method. We additionally make a rough estimate of the error of the semi-analytical procedure by using either the nuclear collision length or the nuclear interaction length when estimating the distance travelled by CR proton before its energy is significantly reduced (see [32]). Our results are in good agreement with the recent literature for mesons productions in the atmosphere [22,41]. Note that the production rates have a ∼ 50% theoretical uncertainty due to the dependence on the hadronic interaction models used to model mesons production in the shower [41].

B. Dark sector productions
We will focus for definiteness on simple LDM scenarios with a dark photon mediator V µ (with mass m V ) interacting with the SM electromagnetic current J µ em via a small kinetic mixing parameter ε and with a dark sector 2 In practice we sum the path lengths for each track, obtained from the travelled length and the mass overburden at its start and end point, see [40].
Since direct on-shell production of the dark photon dominates the production rates, the precise nature of the dark matter does not impact significantly the result of this work. In order to compare with existing limits, we consider the case of an inelastic dark matter structure: with very small splitting between both states m χ1 ∼ m χ2 ≡ m χ . Assuming m V 2m χ , the dark photon decays mostly to dark matter states. The first source of dark photons is the decay of neutral mesons, π 0 → γV , η, η → γV . The typical branching ratio is given by and similarly for η → γV . Thanks to the long life-time of these mesons, this process is only suppressed by two powers of the small kinetic mixing parameter.
On the other hand, cosmic ray showers develop a large electromagnetic component from the radiative decays of those same light neutral mesons π 0 , η. These showers convert the large initial energy of the primary meson into a high number of low energy electrons and positrons. The latter are particularly interesting as they can annihilate with the electrons present in the air with a cross-section given by: We show in Fig. 3 the resulting dark photon fluxes for both the resonant and meson decay production channels for two typical dark matter masses. At low mass the resonant production dominates the π 0 , η → V γ process by more than an order of magnitude, and generate a large flux of "monochromatic" dark photons with energies m 2 V /2m e . 4 4 In the limit of a very light mediator, one recovers the case of millicharge particles. In particular, the relevant processes become e + e − → χχ and π 0 → γχχ. We expect a similar enhancement of the production fluxes for MeV-scale millicharge particles and leave to future work a thorough study of this scenario.

C. LDM detection
We consider first the χ − e scattering process. The differential cross section for χe → χe scattering with respect to the outgoing electron energy E f in the laboratory frame is [42]: where E is the incoming LDM energy and f and s stand for the Dirac fermion and scalar χ respectively; The total signal yield can then be obtained analytically, convoluting the differential cross section with the incoming LDM distribution and the cut efficiency for electron recoil detection. 5 Given the relatively low kinetic energy of the dark matter particles produced in the late stage of atmospheric showers, we also study the coherent nuclear scattering on a nucleus of mass M A . We follow the treatment of [44][45][46] (see also [47,48]) and use: where Q = √ 2M A E r is the momentum exchanged with the nucleus, F Helm is the Helm form factor andf f (E r ) = (E r − E) 2 + E 2 − M A E r (resp.f s (E r ) = (E r − 2E) 2 /2) for a dirac fermion DM (resp. scalar DM). This signal is broadly similar to the recoil from a heavy non-relativist dark matter and can therefore be searched for directly in the various dark matter experiment.
We show in Fig. 4 the differential number of interactions for a LDM with mass m χ = 5 MeV and m χ = 25 MeV as function of the squared exchanged momentum −q 2 , which is linked to the recoil energy E r by −q 2 ≡ 2m Xe E r in the coherent scattering case, and by −q 2 ≡ 2m e (E f − m e ) for electron scattering processes. For the LDM curves with m χ = 5 MeV the rates are dominated by LDM particles from the decay of resonantly-produced dark photon. These dark photons are "mono-chromatic" in that they have all the same energy m 2 V /2m e . Thus the dark matter energy spectrum, and consequently the distribution of recoil energies is bounded by this energy, leading in particular to the threshold seen in the χe curves in Fig. 4.

A. Experiments considered
Given that the atmospheric dark matter flux is uniformly spread on earth, the best sensitivity arises from experiments with the largest or most sensitive detection volume. We will illustrate this by considering both the Super-K neutrino telescope program, and the XENON detectors program.
a. Super-K and upgrades The super-K detector has an ample dataset of neutrino interactions which can be mimicked by the scattering of a LDM particle. Given that the bulk of the events occurs at low recoil energy, we use the supernova neutrino search from [27]. 6 The analysis focused on electron recoils between ∼ 16 MeV to 88 MeV, and included 2853 days of data. Following [15], we consider a sensitivity of 23 events for the full run and project these results to the HyperK [28] design (with a 190 kt fiducial volume) and the Super-K upgrade with gadolinium doping [49] (SuperK-Gd), assuming an improved sensitivity of 0.84 (resp. 0.6 for SuperK-Gd) events per year.
b. XENON1T and future upgrades The coherent nuclear scattering signatures discussed previously can be observed in direct detection experiments. We focus on the XENON program, and in particular the data from the XENON1T [29] standard WIMP search, which shares most of the characteristic of the coherent scattering signal emphasised above. Including the efficiencies for the selection cuts in the fiducial volume [29], we focus on the intermediate background search corresponding to the 0.9t reference volume and put a 95% C.L limit at 5 signal events. 7 We further make projection for 5 years of data for XENONnT based on the recent projections for the WIMP case [30]: we use a 4t fiducial volume and estimated nuclear recoil background of 2.0 events for 20t · yr in the energy range 4 − 50 keV. We then scale this result to obtain a projection for five years of data-taking in the DARWIN experiment [31] with a 30t fiducial mass.
We stress that we have focused on these experimental programs primarily as an illustration of the relevance of resonant atmospheric LDM fluxes. Several other experiments could have also sensitivities to this flux. For instance, on the dark matter side, the liquid argon-based DarkSide20k and Argo projects [51] or for neutrinos detectors, the DUNE far detector [52], KM3-Net [53] and JUNO [54] programs.

B. Results
Resonant production plays a key role in enhancing the production at low masses. We show in Fig. 5 the resulting 95% C.L. limits and projections. The Super-K analysis constrains parameter space comparable to beam neutrino experiments MiniBooNE [55], COHERENT [56] and CCM [46] and slightly below the existing limit from the NA64 [57,58], and BaBar [59] analysis (see the review [60] for a recent summary).
The limits for 3m χ = m V 40 MeV are dominated by the resonant production in the secondary electromagnetic showers, while for 40 MeV m V pion, then η meson decays become the dominant production mode. The sensitivity of nuclear-recoil based analysis in XENON is about an order of magnitude lower than the electron-based in kilotons neutrinos detectors. However, when considering only mesons decays as the production mechanism, it relies only on the LDM interactions with hadronic states and can therefore be used to constrain models where the vector mediator may be lepton or electron-phobic.
We further present projections for the successor of these experiments, starting from the SuperK-Gd experiment with a two-year run with water enriched with Gadolinium [49] and the Hyper-K project, including a long 5-year run and following the background level considered in [20]. Both improvements push the parameter space in ε accessible by around a factor of two each. Regarding the case of coherent scattering signal, we also  = . , = FIG. 5. Limit at the 95% C.L. on our LDM scenario from both the Super-K [27] (orange area) and XENON1T [29] (purple area) experiments, using the ratio mV = 3mχ and αD. The purple dashed line represent the limit from XENON1T without resonant production. We also show projection for SuperK-Gd with Gadolinium (orange dashed line), and Hyper-K (orange dotted line) following [20], and for 5 years of data in the XENONnT and DARWIN projects [31] (blue dashed and dotted lines). The grey regions are the limits from NA64 [58], MiniBooNE [55], COHERENT [56], CCM [46] and BaBar [59]. The green dashed-dotted line represents the projection from DM scattering in the near detector of the DUNE experiment as derived in [32].
show the projection for both the XENONnT [30] and DARWIN [31] (assuming a 5-year run) projects. The above results can be readily recasted into different types of vector mediators. A particularly relevant case is a baryon-number gauge boson, which presents a natural lepton-phobia. The XENON1T limit ε X1t based on meson decays only can for instance be projected into a limit on the baryon gauge coupling g D as g D,X1t = ε X1t × 4πα em Z Xe /A Xe comparable to the recent result from the CCM [46] collaboration.

IV. SUMMARY
We have presented in this letter a new atmospheric production mechanism for light dark sector states based on the annihilation of positrons from CR showers on atmospheric electrons. Comparing with the full production from mesons decays (including secondary mesons) that we have obtained from the complete simulation of the CR showers, resonant production dominates at small masses by more than an order of magnitude. It provides an abundant, albeit low energetic, dark matter flux which can be subsequently searched for in detectors with low recoil thresholds. The SM-only distributions as derived from our full numerical simulation are available on a Zenodo database. 8 They can be used to estimate the production rates of a large range of bosonic Feebly Interacting Particles, from ALPs to millicharge particles, based on the abundant flux of low energy electrons/positrons and photons generated in the showers.
We have updated the present and projected constraints from various neutrinos telescope experiments on this class of new physics candidates, focusing in particular on the dark photon-mediated light dark matter scenario. The future reach of next generation neutrino telescope is remarkable and on par with accelerator-based experiments, providing a strong incentive for the experimental collaborations to consider this type of analysis in the future. Coherent scattering in next generation dark matter experiments was found to provide weaker limits, but can be used to probe models where the vector mediator may exhibit electron-phobia.