An alternative interpretation for cosmic ray peaks

We propose an alternative mechanism based upon dark matter (DM) interpretation for anomalous peak signatures in cosmic ray measurements, assuming an extended dark sector with two DM species. This is contrasted with previous effort to explain various line-like cosmic-ray excesses in the context of DM models where the relevant DM candidate directly annihilates into Standard Model (SM) particles. The heavier DM is assumed to annihilate to an on-shell intermediate state. As the simplest choice, it decays directly into the lighter DM along with an unstable particle which in turn decays to a pair of SM states corresponding to the interesting cosmic anomaly. We show that a sharp continuum energy peak can be readily generated under the proposed DM scenario, depending on dark sector particle mass spectra. Remarkably, such a peak is robustly identified as half the mass of the unstable particle. Furthermore, other underlying mass parameters are analytically related to the shape of energy spectrum. We apply this idea to the two well-known line excesses in the cosmic photon spectrum: 130 GeV gamma-ray line and 3.5 keV X-ray line. Each observed peak spectrum is well-reproduced by theoretical expectation predicated upon our suggested mechanism, and moreover, our resulting best fits provide rather improved chi-square values.


I. INTRODUCTION
Dark matter (DM) is a necessary component in the time evolution of our Universe, and its existence has been supported by a tremendous amount of astrophysical and cosmological evidence for which the relevant observations are made mostly based upon the gravitational interaction of DM (see Ref. [1] for a general review). In fact, none of the Standard Model (SM) particles can explain various DM-related phenomena, so that the detection of any DM candidates can be not only exciting per se. but a strong sign of new physics framework beyond the Standard Model (BSM). A great deal of experimental effort to detect DM signals has been made in three different directions such as i) direct detection experiments by observing recoil energy of target nuclei scattered off by DM, ii) indirect detection experiments by observing cosmic rays originating from DM annihilation or decay, and iii) collider searches (for example, at the Large Hardon Collider at CERN) by actively producing DM candidates and exploiting their collider signatures. These three avenues to DM detection are complementary to one another, and have set the bounds of the viable DM mass and the associated cross section (see Ref. [2] and references therein for review). 1 Of those experimental attempts, indirect detection ex- * email: immworry@ufl.edu † email: log1079@gmail.com 1 Recently, Ref. [3] proposed a general scenario, dubbed "Inflatable DM models", within the context of which many well-motivated DM models having too large production of DM can be remedied, hence evade the bounds without tuning of underlying parameters.
periments have received particular attention as many of them have reported anomalous observations potentially signalling the presence of DM candidates at the locus of cosmic ray sources. For instance, PAMELA [4], Fermi-LAT [5], and AMS-02 [6] found quite a marked rise of the positron fraction in the energy range from roughly 10 to 200 GeV, and similarly ATIC [7], FERMI-LAT [8], and HESS [9] reported an excess in the positron-electron combined energy spectrum between 100 and 1000 GeV. Several photon channels showed intriguing excesses such as 3.5 keV line [10,11], 511 keV line [12], GeV bump [13], and 130 GeV line [14,15]. The positron excesses and the Galactic Center (GC) GeV γ-ray excess are featured by a continuum bump, while the other three X/γ-ray excesses showed a sharp peak within a very narrow energy range.
The latter class of excesses are particularly interesting because they can be readily connected to the DM interpretation. As typical DM candidates behave nonrelativistically, the photon energy from a DM pair annihilation (or 2-body decay) is monochromatic, being the same as (half) the DM mass. 2 In this context, many DM models to address those excesses have been introduced and studied in literature: for example, Ref. [16] for 511 keV line, Refs. [17,18] for 130 GeV line, and Ref. [19] for 3.5 keV line. In reality, the relevant signal spectrum does not appear as a δ-function-like peak but is smeared to some extent because of imperfection in cosmic ray detectors. With the assumption of a Gaussian smearing, the resultant γ-ray energy spectrum (typically) becomes symmetric with respect to the nominal peak.
This broadening effect has motivated the possibility of non-minimal DM scenarios for interpreting the narrow width of peaks as the physical, not the one induced from instrumental uncertainties. The next-to-minimal DM models hypothesize the situation where DM particles annihilate or decay into on-shell intermediate particles that decay into photons [20,21]. As such an on-shell intermediate particle comes along with a fixed boost, the photon energy spectrum is characterized by a rectangular shape [20,21]. If the mass gap between the DM and the intermediate particle is sufficiently small, hence so is the boost factor, then the photon energy spectrum becomes narrow enough, potentially being indistinguishable from the signal spectrum by the minimal scenario.
Having similar philosophy and positing the DM interpretation, we here propose a new mechanism to develop a narrow continuum energy spectrum which would fake a sharp spike. The research program to explain the excesses in cosmic ray energy spectra with the "energypeak" emerging under non-minimal DM frameworks has been initiated by Ref. [22], in which various observations of the energy-peak made in the context of collider physics [23][24][25][26][27][28] have been applied to the GC γ-ray GeV excess. As in [22], we begin the discussion with noting that multiple DM species could exist in the Universe, and the relevant DM models constructed upon such a DM framework can bring about not only nontrivial cosmological implications, e.g., "assisted freeze-out" [29], but interesting phenomenology, e.g., "boosted DM" [30][31][32]. In this context, the proposed mechanism involves a nonminimal dark sector containing multiple DM particles.
For the purpose of simplicity, we introduce two DM species one of which is assumed heavier than the other, denoting henceforth the former and the latter as χ h and χ l , respectively. The heavier DM communicates to the SM sector not directly but through the lighter DM. In addition, the heavier DM pair-annihilates into a pair of on-shell intermediate states (denoted as A) each of which subsequently decays into the lighter DM together with a dark pion or an axion-like particle (ALP) (denoted as a) emitting a couple of photons in the final state. 3 FIG. 1 schematizes the DM scenario that we consider throughout this paper. 4 We point out that although we employ the photon final state as a concrete example for elaborating our mechanism, it is straightforwardly extensible to other visible particle final states, e.g., e + e − .
In this DM scenario, the on-shell intermediate particle 3 In general, A and χ l can be either dark or SM sector particles (and may be even unstable) unless they are stringently constrained by other observations. However, for the sake of simplicity, we assume that they are dark sector particles. 4 The main idea in this paper can be readily applied to the decaying DM scenario, but we keep the annihilating one as a concrete example. A comes with a fixed boost factor, leading to a rectangular energy spectrum, hence a rectangular boost distribution for particle a. Due to variations in the boost factor of a the emitted photons develop a continuum energy spectrum whose width is determined by the mass parameters involved in the entire process demonstrated in FIG. 1. We remark that the relevant scenario invokes from the (relatively) narrow energy spectrum to the broadly-distributed energy spectrum, depending on details of the DM models of interest. We shall briefly discuss the aspect of the wide continuum bump signature in Sec. III, while primarily focusing on the "peak-faking" interpretation and the associated mass spectrum. Not surprisingly, the mass spectrum of m χ h > ∼ m A > ∼ m a +m χ l renders a relatively narrow photon energy spectrum. As a characteristic feature, the resulting differential energy spectrum is symmetric in the logarithmic scale, and remarkably, its center position is identified as half the mass of particle a, m a /2 [23,33]. Therefore, these structural properties enable us to not only distinguish this DM model scenario from other standard DM interpretations, but probe/measure the mass parameters of some dark sector particles.
To present our main idea, this paper is organized as follows. We first begin with the DM model under consideration in the next section. In Sec. III, we discuss the energy spectra of relevant visible particles arising from the DM model introduced in Sec. II, mainly focusing on their functional structure. We then apply the main idea to a couple of photon peaks in Sec. IV: i) Fermi-LAT 130 GeV excess and ii) 3.5 keV excess. Sec. V is reserved for our conclusions.

II. DARK MATTER MODEL
We here discuss a class of DM models for which our main idea is applied, and offer a viable DM model that realizes the relevant scenario. As briefly mentioned earlier, we imagine that the dark sector is non-minimal, meaning that there exist more than one DM candidate. Although the arbitrary number of DM species could be introduced, we employ only two types of DM particles, the heavier (χ h ) and the lighter (χ l ), for simplicity.
The lighter DM is assumed to directly communicate to the SM sector, whereas the heavier DM is set out to have interactions with the SM sector via the lighter DM. In this sense, the relic abundance of χ h can be computed by the scheme of "assisted freeze-out" [29]. We further assume that χ h has a contact with χ l not directly, but via an on-shell intermediate state A, i.e., a pair of heavier DM particles annihilate into a pair of A's each of which subsequently decays into χ l together with a dark pion or an ALP as depicted in FIG. 1. Finally, the dark pion or the ALP decays into a photon pair whose energy spectrum is the major interest of this paper. We remark that it is not a necessary condition for particle A to get pair-produced, i.e., A could be produced in association with another particle A as the detailed dynamics of A is irrelevant to the later argument and formalism. Similarly, A does not need to directly decay into the lighter DM, i.e., χ l could be replaced by other heavy particles whose detailed dynamics does not affect the later argument and formalism. In this context, the model set-up demonstrated in FIG. 1

is a simplified version.
A possible realization of the DM scenario at hand can be summarized as follows. We consider two fermionic DM species, χ h and χ l , with an intermediate fermionic state ψ A and a singlet pseudo-scalar a (e.g., a dark pion or an ALP as mentioned earlier). Then, the effective Lagrangian required for the DM scenario exhibited in FIG. 1 is simply described by the following operators: where F µν denotes the dual field strength tensor as usual, and Λ and f a describe the associated suppression scales whose details can be revealed by appropriate UV completion. The first term ensures an s-wave annihilation of the heavier DM, i.e., χ h χ h → ψ A ψ A , the second induces the decay of ψ A into χ l and a, and the last corresponds to two photon decay of a as in FIG. 1, respectively. The stability of χ h and χ l can be easily achieved with separate symmetries, e.g., U(1) ⊗ U(1) [29] or Z 2 ⊗ Z 2 [30]. We again stress that the Lagrangian in Eq. (1) is a simple realization and there exist a host of other possibilities to accommodate the event topology in FIG. 1. Exhausting all of them is, however, beyond the scope of this paper.

III. ENERGY SPECTRUM
In this section, we derive the analytic expression for the gamma-ray energy spectrum and discuss its properties which could be distinguished from other (standard) scenarios, given the scenario in FIG. 1.

A. Derivation of the analytic expression
Assuming that the heavier DM particles are nonrelativistic, their pair annihilation into two A's leads to a fixed boost of particle A (denoted as γ A ) relating the two mass parameters by γ A = m χ h /m A with m i symbolizing the mass of particle species i. Since A obtains a non-zero boost factor, the energy of particle a is not monochromatic, but given by a broad spectrum. The a energy, E a , measured in the laboratory frame is parameterized as where θ * a is the emission angle of a in the A rest frame with respect to the boost direction of A and E * a is the fixed a energy measured in the rest frame of particle A: If A is either a scalar or produced in an unpolarized way, then cos θ * a becomes a flat variable, resulting in a rectangular distribution in E a by a simple chain rule whose range is given by Similarly to Eq. (2), the observed photon energy for a fixed γ a is expressed as where θ * γ denotes the intersecting angle between its emission direction and the boost direction of particle a in the a rest frame and E * γ is the fixed photon energy measured in the a rest frame, that is, half the mass of particle a: Unlike γ A in the case of particle a, γ a is not single-valued but distributed. Denoting its distribution by g(γ a ), from Eq. (4) we find the unit-normalized expression (or equivalently, probability distribution function) for g(γ a ) as where Θ(x) is the usual Heaviside step function, and γ ± a are defined by Here we used the fact that cos θ * γ is a flat variable so that g(γ a ) develops a rectangular distribution as well.
We then find that for any fixed γ a , the unit-normalized differential energy distribution is where E ± γ can be obtained by setting cos θ * γ to be ±1 in Eq. (5): Denoted by f (E γ ), the expression for the unit-normalized total energy spectrum can be obtained by summing Eq. (9) over all relevant γ a 's, that is, where γ min a and γ max a are defined as Since g(γ a ) is upper-bounded, γ + a determines the spanning range of f (E γ ) as follows: We finally remark that in the actual data analysis with concrete examples, all prefactors in Eq. (12) are eventually absorbed into the overall normalization parameter N , and as a consequence the shape of f (E γ ) is completely determined by γ + a , γ − a , and E * γ , i.e., there are four independent fit parameters.

B. Functional properties and discussions
To discuss functional properties of the photon energy spectrum, we first revisit the expression of E γ for a fixed γ a shown in Eq. (5). Since cos θ * γ spans −1 to +1, the range of E γ is trivially given by as also expressed in Eq. (10). One remarkable feature from the above range is the fact that the lower (upper) end in the right-hand side is smaller (greater) than 1, implying that E * γ is the only commonly-included energy value for any γ a . Moreover, we observe that E * γ is the geometric mean of minimum and maximum energy values, i.e., (E * γ ) 2 = E min γ E max γ . This again implies that the E * γ value is located at the center of the E γ distribution for a given γ a in the logarithmic scale. As mentioned before, the shape of E γ distribution for a fixed γ a is rectangular due to the flatness of cos θ * γ in Eq. (5). In order to obtain the overall energy distribution, one should "stack up" all of such rectangles contributing to the energy distribution. This statement is already formulated in Eq. (11) as a Lebesque-type integral representation. Hence, one would expect that the resultant photon energy distribution contains a unique peak at E γ = E * γ in conjunction with the observation made in the previous paragraph. Indeed, there arises a subtlety here: the validity of this expectation depends whether or not the smallest boost factor of particle a approaches 1.
The condition for γ a = 1 (i.e., solving Eq. (8) with γ − a being 1) is E * a = m a γ A , representing a hypersurface formed by m A , m χ h , and m χ l for any fixed m a , i.e., only delicately chosen mass spectra can attain this condition. Once γ a = 1 is available, clearly, E * γ appears as a unique and cusp-structured peak [23,26]. The left panel in FIG. 2 demonstrates an example spectrum in this category in the logarithmic scale. The chosen mass spectrum is (m χ h , m A , m χ l , m a ) = (237.5, 200, 50, 100) GeV, and events were generated by pure phase space. The simulated events are binned to the blue histogram, and the associated theory prediction is shown by the red line. We clearly see that the theory expectation can reproduce the data well enough, and the spectrum is symmetric with respect to E * γ = m a /2 = 50 GeV (indicated by a black dashed line) in this scale. We also remark that both sides of the distribution look like straight lines, which can be easily seen from Eq. (12) with γ min in logarithmic E γ , so that the whole spectrum appears as an isosceles triangle.
On the other hand, g(γ a ) starts from the value away from γ a = 1, no peak is developed in the middle of the energy distribution. Instead, a plateau region emerges because even the narrowest rectangle corresponding to the smallest γ a has a finite-sized width. Nevertheless, it is straightforward that the center of the relevant photon energy spectrum can be identified as E * γ in the logarithmic base. These expectations are manifestly shown in the right panel of FIG. 2. The relevant mass spectrum is the same as the previous case with m A replaced by 170 GeV. Also, the existence of a plateau region makes the whole energy spectrum appears as an isosceles trapezoid in the logarithmic base. This plateau structure is a distinguished feature from other energy spectra. However, its presence may be invisible in actual indirect detection experiments due to the issue of energy resolution. Clearly, if the plateau is smaller than the relevant resolving power, its existence is rarely identifiable so that the relevant spectrum easily fakes a unimodal distribution like the previous case. Even for the energy spectrum with a sufficiently sizable plateau, its identification may be unavailable with small statistics, i.e., more data accumulation may be required. As the proposed mechanism, which is denoted as Scenario iii), is used for explaining narrow peaks, it is interesting to compare and contrast it with the other two conventional scenarios enumerated below.
• Scenario i): Photons directly come from DM annihilation/decay. • Scenario ii): Photons are emitted as a decay product of the on-shell intermediate particle into which DM annihilates/decays.
In Scenario i), the width of the peak is typically caused by the intrinsic energy resolution of detectors, and thus the final energy spectrum is symmetric about the nominal peak. Even in the logarithmic scale, the spectrum is described by a smooth curve, while it is no more symmetriclooking about the peak. In Scenario ii), the width of the energy spectrum is physical. However, identifying the peak position is ambiguous due to the box-like spectral behavior in both the linear and the logarithmic scales. Therefore, defining the symmetry property of the shape is not available. The comparisons thus far are summarized in Table I. One can easily see that their respective morphological features differ from one another, and therefore, one is able to pin down the underlying DM scenario with a reasonable amount of signal statistics.
Before closing the current section, we briefly discuss the case of wide continuum bump spectra. Certainly, such types of spectra arise within a broad realm of relevant parameter space so that the theory prediction in Eq. (12) can be employed to explain continuum cosmic ray excesses. Due to the unique morphological features discussed so far, the relevant signal spectrum could be easily distinguished from other continuum bumps, taken as strong evidence of a non-trivial dark sector.

IV. APPLICATIONS
Armed with the argument in the previous section, we apply the basic idea to a couple of existent cosmic ray peaks: i) 130 GeV line [14,15] and ii) 3.5 keV line [10,11]. Although claiming that these examples could be understood by the proposed mechanism, we admit that the underlying DM models for them might not fall into the scenario depicted in FIG. 1. Nevertheless, we emphasize that the applicability of the relevant technique is restricted to neither these two examples nor gamma-ray spectra. In other words, nothing precludes us from applying the DM interpretation at hand for any of future excessive cosmic ray signals with a narrow peak.

A. 130 GeV line
Our first example is the famous 130 GeV line whose original data was collected by the FERMI-LAT collaboration. We basically conduct the fit to the spectrum of the observed 130 GeV gamma-ray excess with the ex-pected shape in Eq. (12) with the overall normalization parameter N S being added: The relevant data points and errors are taken from Reg4 with the ULTRACLEAN event class in Ref. [15]. Since the measured bin counts contain the contributions from backgrounds, the relevant fit is performed simultaneously with the background template. We assume that the backgrounds can be parameterized by a simple (graduallyfalling) power law such as where N B is the normalization parameter for backgrounds and p encodes the background shape. Here E * γ is the same E * γ contained in γ min a of Eq. (16). Therefore, the entire data set is fitted with the combination of f S (E γ ) in Eq. (16) and f B (E γ ) in Eq. (17), i.e., Our fit result is shown in the upper-left panel of FIG. 3. The data points and associated error bars are represented by black dots and black lines, correspondingly, while the red solid curve shows the best-fitted model. We clearly see that our fit is in a rather good agreement with the gamma-ray energy spectrum, i.e., the template given in Eq. (18) reproduces the data sufficiently well. More quantitatively speaking, we find that the χ 2 value is 19.7 for 19 degrees of freedom (i.e., 25 data points subtracted by 6 fitting parameters such as N S , N B , γ + a , γ − a , E * γ , and p) between 80 and 210 GeV. This number is quite comparable to that from conventional scenarios, suggesting that our DM scenario be considered (equally) plausible in explaining the observed data.
The fit also tells us useful information. First of all, we extract the mass of particle a from the measurement of E * γ (= m a /2), that is, where the errors here are evaluated by 1σ statistical uncertainty. The other mass parameters can be also estimated by the measured γ + a and γ − a , for which the bestfitted numbers are reported as (γ − a ) ext = 1.0000 +0.0012 , respectively. Again, the errors are reported by 1σ statistical uncertainty. Note that the lower error for (γ − a ) ext is not provided because there is no sensitivity to the values of γ − a < 1 according to the definition of γ min a in Eq. (13).
From Eq. (8), we obtain the expressions for E * a and p * a in and the difference between the former squared and the latter squared leads the following mass relation: with which we perform a scan of allowed mass space. As γ − a shows the smallest error, we simply fix it to be unity for phenomenological purpose. Instead of doing three-dimensional scanning, we choose three different The red, green, and blue regions are for m χ l = 10, 100, and 1000 GeV, respectively. Solid curves are the contours evaluated with the best-fitted m a and γ + a . On the other hand, dashed (dot-dashed) curves are the ones evaluated with m a + δm a and γ + a + δγ + a (m a − δm a and γ + a − δγ + a ), so that one can get some intuition on the mass spectrum allowed by 1σ variations of the relevant parameters. We observe that the viable mass spectra are more or less compact. This is not surprising because the narrow peak enforces a degenerate mass spectrum not to obtain too large boost in any of the steps in FIG. 1. Dark sector scenarios featuring such a mass spectrum could be achieved by a symmetry. Building realistic DM models is, however, beyond the scope of this paper, so we do not further pursue this direction here.

B. 3.5 keV line
Our next example is the well-known 3.5 keV line. As in the previous case, the fit is performed to the spectrum of the observed 3.5 keV gamma-ray excess with the signal template given in Eq. (16). The relevant data points and associated errors are taken from the processed data for the MOS spectrum of the central region of the Andromeda galaxy (M31) found in Ref. [11]. The numbers to be used contain only the DM component, i.e., the background component is already subtracted. Therefore, we execute the fit procedure only with the signal template, not introducing any background template unlike the previous case.
The fit result is exhibited in the upper-right panel of FIG. 3. Again, the data points and the error bars are represented by black dots and black lines, respectively, while the red solid line describes the best-fitted model. Our fit is in a very good agreement with the measured energy spectrum. The reported χ 2 value is 4.61 for 13 degrees  [15]. The fit is performed with 25 data points (black dots), and the best-fitted is represented by the red curve. Upper-right panel: the 3.5 keV X-ray energy spectrum of the DM components taken from Ref. [11]. The fit is performed with 17 data points. Lower-left panel: the allowed mass space for the 130 GeV line along with 1σ variations of ma and γ + a for three different masses of particle χ l . Lower-right panel: the allowed mass space for the 3.5 keV line.
of freedom (i.e., 17 data points subtracted by 4 fitting parameters such as N S , γ + a , γ − a , and E * γ ) between 3 and 4 keV. This number is significantly improved, compared with that from standard interpretations, and therefore, our DM framework can be considered as a plausible scenario accommodating the observed spectrum.
Speaking of various best-fit parameters, we first find that the extracted mass parameter for particle a is m ext a = 7.09 +0.08 −0.06 keV, where errors are estimated by 1σ statistical uncertainty. We find the best-fit values for γ + a and γ − a are (γ + a ) ext = 1.0021 +0.0015 −0.0006 , (γ − a ) ext = 1.0000 +0.0014 , respectively, together with 1σ statistical uncertainty. In order to obtain the allowed mass space within 1σ variations of the relevant parameters, we follow the same procedure in the previous case for three different m χ l masses, 1, 10, and 100 keV. The lower-right panel of FIG. 3 demonstrates the allowed region again in the plane of (m A − m χ l − m a ) versus (m χ h − m A − m χ l − m a ). The red, green, and blue regions correspond to m χ l = 1, 10, and 100 keV, respectively.

V. CONCLUSIONS
Dark matter indirect detection offers an excellent opportunity to confirm the existence of DM. Several experimental collaborations have already reported anomalous phenomena in the relevant cosmic ray energy spectrum. Particular attention has been paid to sharply-peaked signals due to the easiness of DM interpretation. Typical models assume that (non-relativistic) DM particles directly annihilate or decay into visible ones. In this minimal DM scenario, the narrow width of the peak is typically understood as the one stemming from the imperfection of cosmic ray detectors.
We have rather taken the viewpoint that this width can be physically induced, and proposed a new mechanism to realize it with the assumption of a non-minimal dark sector. Two DM species were introduced, and the heavier one is assumed to annihilate to the on-shell intermediate state which subsequently decays into the lighter one and an unstable particle. This unstable particle further decays into a pair of visible particles which can be a source of anomalous peaks in the cosmic ray energy spectra. We showed that the signal spectrum can be narrow enough to fake a sharp spike with a suitable choice of the associated mass parameters.
The shape of the full signal spectrum was derived, and several interesting functional properties were discussed. We pointed out that the peak position, one of the fit parameters, is immediately identified as half the mass of the above-mentioned unstable particle. We also showed that other mass parameters can be estimated by other fit parameters. We then enumerated various morphological features to be utilized for distinguishing several DM sce-narios in which the sharp peak signature is available. The viability of the relevant strategy was assessed with two real observational data sets, 130 GeV line and 3.5 keV line. We found that both of them can be well-described by the theoretical expectation in our DM scenario, and that each allowed parameter space is fairly large.
Finally, we emphasize that our DM scenario and the associated data analysis are not restricted to photon energy peaks, i.e., any cosmic ray energy peaks can take advantage of them. We strongly encourage people to pursue the direction exploited in this paper as well whenever cosmic ray peaks are observed.