Indirect dark-matter detection with MadDM v3.2 – Lines and Loops

Automated tools for the computation of particle physics’ processes have become the backbone of phenomenological studies beyond the standard model. Here, we present MadDM v3.2. This release enables the fully automated computation of loop-induced dark-matter annihilation processes, relevant for indirect detection observables. Special emphasis lies on the annihilation into γX\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma X$$\end{document}, where X=γ,Z,h\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$X=\gamma , Z, h$$\end{document} or any new particle even under the dark symmetry. These processes lead to the sharp spectral feature of monochromatic gamma lines – a smoking-gun signature of dark matter in our Galaxy. MadDM provides the predictions for the respective fluxes near-Earth and derives constraints from the gamma-ray line searches by Fermi-LAT and HESS. As an application, we discuss the implications for the viable parameter space of a top-philic t-channel mediator model and the inert doublet model.


Introduction
Indirect detection is one of the three main strategies to search for dark-matter interactions beyond the gravitational one.It explores traces of dark-matter annihilation or decay in cosmic messengers such as photons, neutrinos, positrons, anti-protons or heavier anti-nuclei.While primary annihilation into unstable standard-model particles typically leads to a continuum spectrum of messenger particles, direct annihilation into neutral messengers can give rise to prominent spectral features, such as monochromatic lines [1][2][3].This is a smokinggun signal for dark matter as the pronounced peak is hardly mimicked by any astrophysical background.While neutrino monochromatic lines typically arise at tree level, see e.g.[4][5][6][7], annihilation into photons is a loop-induced process as dark matter is commonly considered to be electrically neutral.Even though the signal is suppressed by a loop factor, the experimental sensitivity to such sharp energy spectra is high, such that current gamma-ray telescopes, like Fermi-LAT [8] and HESS [9,10], set strong constraints on the dark-matter parameter space.
Corresponding to the lowest order in perturbation theory, in renormalizable models, the leading diagrams for annihilating into pairs of photons (or of a photon and a neutral particle such as Z and h) are nondivergent.Nevertheless, their computation can be involved.For instance, the number of diagrams to be computed increases considerably for intricate models that contain many (charged) particles.This calls for the need for an automatized procedure to obtain predictions for gamma-ray lines within generic dark matter models.
In this paper, we present MadDM v3.2.It features the automatized computation of loop-induced processes for indirect detection within any dark-matter model for which a UFO [11] model at next-to-leading order (NLO) can be generated.This new module builds upon the indirect detection module, released with MadDM v3.0 [12] (see also Ref. [13] for a recent update including a short user guide).At the time of writing, it is the only numerical tool with such capability.The currently available public packages provide the computation of loop-induced annihilation cross section for specific models and operators, partly based on analytic expressions obtained in the literature.For instance, mi-crOMEGAs [14,15] and DarkSUSY [16,17] include expressions for the supersymmetric neutralino dark matter in the minimal supersymmetric standard model (MSSM) [18,19], and the inert doublet model (IDM) [20].In addition, micrOMEGAs provides automatized com-arXiv:2107.04598v2[hep-ph] 23 Mar 2023 putations in several extensions of the MSSM by being linked to the package SloopS [21,22].Further results have been obtained for simplified models [23][24][25][26][27], the next-to-MSSM [28][29][30], and Kaluza-Klein dark matter [31][32][33].Alternatively, there have been efforts in obtaining analytic expression in rather model-independent approaches, see e.g.[34,35] or utilizing an effective theory approach, where the heavy scale physics which gives rise to the interaction between dark matter and the photon is integrated out and the loop processes are reduced to an effective vertex, see e.g.[36][37][38][39][40][41][42][43][44][45][46].
Within MadDM v3.2 we utilize MadLoop [47] for the computation of loop-induced annihilation processes.MadLoop has been embedded in MG5 aMC [48] and has been extensively used to compute QCD and EW corrections in the standard model, and QCD corrections e.g.[49][50][51][52], and loop-induced processes [53,54] in the framework of collider searches for dark matter.The automated computation of cross sections for gamma-lines within MadDM v3.2 has already been introduced in [55] for t-channel simplified models, a category of models where loop-induced processes are highly relevant for indirect dark-matter searches, see e.g.[26,27,52].Together with this paper, we release MadDM v3.2 (hereafter simply MadDM).It includes the automated computation of the respective fluxes for a variety of dark-matter density profiles and exclusion limits from utilizing results from the Fermi-LAT satellite and the HESS telescope.
While the focus of this new MadDM release is on sharp spectral features, its capability is more general and accounts for any loop-induced dark-matter annihilation channel: an example is annihilation into pair of gluons, which arises at one loop in case of top-philic dark matter, see e.g.[50,56,57].These processes lead to a continuum gamma-ray energy spectrum because of the showering and hadronization of the gluons and can, for instance, be constrained by the Fermi-LAT dwarf spheroidal analysis [58,59].The new MadDM is designed to recognize the type of energy spectrum produced by the loop-induced process and automatically assign it to the correct analysis pipeline.Note that the automatized loop process computation is performed exclusively for indirect detection.The automated inclusion of loop-induced processes and higher-order corrections in the relic density computation 1 goes beyond the scope of this article and is left for future work.
The remainder of the paper is organised as follows.In section 2, we discuss the generation of one-loop processes within MadDM, provide details on the required UFO model files and outline the main functionalities of the program.Relevant astrophysical aspects of gamma- 1 See refs.[60][61][62][63] for related work.
ray line searches are detailed in 3, while section 4 provides a phenomenological study of gamma-ray line signatures within two simple dark-matter models: a topphilic t-channel mediator model and the IDM.We conclude in section 5. Further details and technicalities can be found in the appendices.In appendix Appendix A, we provide a short user guide focusing on the new features.In appendix Appendix B, we describe the implementation of automatic J-factor computation and provide more information about the experimental data.In appendix Appendix C we describe the treatment of multiple line signatures within MadDM.

MadLoop interface
One-loop computation is now a very mature topic with a large number of tools and libraries designed for such computations (see the review [64] and references therein).In the context of this work, we use MadLoop [47], more precisely, the loop-induced implementation of it [53].The task of MadLoop is to generate all the Feynman diagrams and to numerically evaluate the numerator of such diagram either as a complex number or as a polynomial in the loop-momenta.Such information can be used by Collier [65], Ninja [66][67][68], CutTools [69] or IREGI [70], to decompose the loop in a sum of scalar integrals either employing tensor integral reduction (TIR, introduced by Passarino & Veltman [71]) or by performing the reduction at the integrand level (OPP method [72]).Finally, the library OneLoop [73] or QCDLoop [74] are used to return the finite part and the pole of the associated loop.
A key characteristic of MadLoop is its flexibility to use all the above-mentioned alternative tools to find the most suitable for the decomposition of a given loop integral.Such a decomposition can be numerically unstable (for example in TIR when the Gram determinant vanishes) and can lead to unreliable results.For a given computational tool, MadLoop assesses the numerical stability of the result either using the estimator of the associated method or by simply re-evaluating the same quantity after a boost or rotation.If the error is higher than a user-specified threshold, another tool is used.If none of the methods encoded in MadLoop returns a result precise enough, the code redoes the computation in quadruple precision (available for Ninja and Cut-Tools) which is two orders of magnitude slower.The capacity of MadLoop to test various reduction algorithms allows minimizing the number of times quadruple precision is needed.
In general, loop computations require special attention to properly handle the finite parts coming from dimensional regularization (where all quantities are defined in d dimension with d = 4 − 2 ).Such terms have two origins.They arise from the epsilon part of the denominator (called R 1 ) and of the numerator (called R 2 ).Thanks to the fact that the denominators have a simple and well-defined analytical structure, R 1 can always be systematically reconstructed in the loop evaluation procedure.The computation of R 2 is more complex and is typically done analytically before the computation of the loop, i.e. at the stage of generating the UFO model files, as described below.
The computation of loop-induced processes is often more challenging than the one related to NLO processes associated with a Born amplitude B. This has multiple reasons.First, in standard computations, one can do the reduction of the loop directly on the amplitude summed/averaged over the helicities of the initial/final states, while in the case of loop-induced processes, the reduction needs to be performed for each helicity combination independently.To mitigate the slowing down, the strategy used in MG5 aMC and MadDM is to Monte-Carlo over helicity contributions.Second, in computations at NLO one can avoid evaluating the loop contribution for each phase-space point.The trick is to use the Born amplitude (up to a normalization factor) as an approximation of the loop amplitude L, which corresponds at the integral level at the trivial formula L = αB+ (L−αB).The application of importance sampling on the latter equality reduces considerably the number of times the loop is computed.This trick is obviously not possible for loop-induced processes rendering their evaluation computationally intensive.A dedicated effort on the phase-space parametrization (and on the parallelization) has been achieved in [53].
To be able to perform loop computations within MG5 aMC and MadDM it is necessary to import NLO UFO model files.This can be achieved by using FeynRules [75], FeynArts [76] and NLOCT [77].NLOCT is a specific software enabling the analytical computation of the R 2 terms and ultra-violet (UV) counter-terms for a given model.The latter are required for NLO computation and need to be included in the NLO UFO format, although they are not used for the computation of loop-induced processes.To be more specific, the user has to implement the model in FeynRules by writing the Lagrangian in the appropriate form.Afterwards, it is renormalized within Feyn-Rules.To enable electroweak loops the flag QCDonly needs to be set to False.Feynman gauge should be used.The renormalized Lagrangian is then passed to FeynArts, which expresses the various interaction ver-tices in terms of their couplings and Lorentz structures.Subsequently, NLOCT is called to solve the renormalization conditions and compute the UV counter-terms and the R 2 terms.The lists of the R 2 terms and UV counter-terms are written on an external file by NLOCT, which must be imported in FeynRules to obtain the NLO UFO format of the model.Depending on the model, the file .nloproduced by NLOCT could be sizeable, making the exporting of the model slow or unstable.It is advisable to make use of the Assumptions list when running NLOCT to specify the possible relations between the parameters of the model (for example by writing explicitly the relations between the particles of the standard model).Those relations will be used during the computation to shorten the analytical expressions of the UV and R 2 terms.It is possible to export a model using either the MS or MS renormalization schemes; the latter is recommended.Note that in the presence of unstable internal particles in the loop the use of the complex mass scheme can improve numerical stability.It can be enabled by choosing ComplexMass->True in the WriteCT command, see [77] for details.

Main functionality
For dark-matter phenomenology, loop-induced annihilation processes with one or two photons in the final state are of particular interest as it provides a smokinggun signature of monochromatic gamma-lines.MadDM automatically generates all contributing diagrams for a given dark-matter model by running the command:

MadDM> generate indirect_spectral_features
More precisely, it generates all diagrams for the final states γγ and γX, where X includes the Z and h particles of the standard model as well as all additional beyond standard model (BSM) particles, that are lighter than twice the dark-matter mass and transform even under the dark symmetry that stabilizes dark matter.Individual channels can be generated by explicitly specifying the final state, e.g.

MadDM> generate indirect_spectral_features a z
Besides the annihilation cross section, MadDM computes the fluxes and the experimental constraints as discussed in section 3.More details on the commands introduced with this release and the corresponding output are provided in appendix Appendix A.
Loop-induced annihilation into other final states may contribute to the continuum flux of cosmic messengers.A prominent example is annihilation into a pair of gluons, gg, that subsequently shower and hadronize.This channel can, for instance, become relevant in top-philic dark-matter models, where dark matter couples only to the top quark at tree level, see e.g.[50,56,57].In the absence of a tree-level diagram for the channel, MadDM automatically switches to the loop-induced mode.Hence, this annihilation channel is considered by executing the command:

MadDM> generate indirect_detection g g
It can also be computed together with tree-level diagrams:

MadDM> generate indirect_detection MadDM> add indirect_detection g g
Here, the first line leads to the computation of all 2 → 2 tree-level annihilation processes.For more details about the existing syntax of MadDM see [12,13].)After the computation of the annihilation cross section and generation of events, MadDM proceeds with the indirect detection analysis pipeline as introduced in [12]: The energy spectra are either obtained using Pythia 8 [78] or the PPPC4DMID tables [79].Subsequently, the energy spectra of all contributions (including tree-level contributions if present) are summed up to obtain the total fluxes at source and (if requested) near Earth.The corresponding gamma-ray flux is confronted with constraints from Fermi-LAT observations of dwarf spheroidal galaxies [58].
3 Gamma-ray line phenomenology

Gamma-ray flux
Photons travel straight (i.e. on geodesics) from the production to the detection point and, hence, they trace their sources.In general, the differential gamma-ray flux integrated over the region of interest (ROI) in the sky is given by: where σv i is the velocity averaged cross-section of dark-matter particles with a mass m DM ) into final states labeled by i. E γ and N γ is the photon energy and the number of photons per annihilation, respectively.Accordingly, dN i γ /dE γ is the differential gamma-ray energy spectrum per annihilation.For dark-matter candidates that are not self-conjugated, eq.(1) has to be multiplied by an additional factor of 1/2.MadDM provides both the differential flux as well as the total integrated flux.The second part of the equation defines the J factor: where ρ(r) denotes the dark-matter density distribution.The second integral integrated is performed over the line of sight (los) l.The most commonly assumed dark-matter density profiles are spherically symmetric and given by: -Generalised Navarro-Frenk-White (gNFW) [80] ρ gNFW (r) = ρ s r s r -Einasto [81] ρ -Burkert [82] ρ -Isothermal [83] ρ Iso (r In all density profiles, the parameters r s and ρ s are the scale radius and the scale density, respectively.For a given r s , ρ s is normalized to match the specified energy density measured at the Sun position (by default R = 8.5 kpc and ρ = 0.4 GeV/cm 3 is chosen).From the gNFW, eq. ( 3), the original NFW and the contracted NFW (NFWc) density profiles are obtained for the choice γ = 1 and γ = 1.3, respectively.For the Einasto profile, α defines the curvature of the density profile, and it is usually fixed at the value α = 0.17.
MadDM allows for the computation of the J-factors for the above-listed dark-matter density profile with general parameters and a generic parametrization of the ROI, centered around the galactic center.Details about this parametrization and the J-factor computation are provided in appendix Appendix B.1.

Spectral feature
Dark matter is assumed to be non-relativistic in galactic halos: typical values within the Milky Way [84] and for similar galaxies are of the order of 10 −3 c, while dwarf spheroidal galaxies are expected to host even colder dark matter, v ∼ 10 −5 c.Hence, for annihilation into γγ or γX, the photon energy spectrum, dN line γ /dE γ in eq. ( 1), is a sharp spectral line: as a consequence of the kinematics of the two-body annihilation process [85][86][87].Accordingly, for γγ, and for γX, where X is a neutral standard-model or BSM particle with mass m X .This characteristic shape allows for peak searches within the experimental data.The implementation of the latter is described in section 3.3 together with the specific experimental ROIs for the Fermi-LAT satellite and HESS telescope, optimized for gamma-ray line searches.

Experimental constraints
Gamma-ray line signals from dark-matter annihilation are hardly mimicked by any standard astrophysical background.As a consequence, the sensitivity of peak searches pursued by the astrophysical experiments is high, typically compensating the loop-suppression of the annihilation rate.As the dark-matter density profile peaks at the Galactic center, it is expected to provide the largest signal.However, depending on the morphology of the continuum gamma-ray background and the cuspiness of the dark-matter density profile, different choices for the ROI maximize the expected sensitivity.Here we consider current upper limits at 95% confidence level (CL) on the dark-matter annihilation cross section into a pair of photons from the Fermi-LAT satellite [8] and the HESS telescope [9,10].Note that gamma-ray line searches have also been proposed in dwarf spheroidal galaxies, see e.g.[88].However, we do not consider this possibility here.Fermi-LAT data from the Galactic center, collected over 5.8 years of observation, provide an exclusion bound reaching σv UL γγ ∼ 2 × 10 −30 cm 3 /s for m DM = 1 GeV and σv UL γγ ∼ 3 × 10 −28 cm 3 /s for m DM = 500 GeV for NFWc, the cuspiest dark-matter density profile.The analysis takes into account four dark-matter density profiles and associates to each of them an optimized ROI.The latter is defined as a circular region centered on the Galactic center, and the galactic plane is masked except for a 12 • ×10 • box centered on the Galactic center.The ROIs are R3, R16, R41, and R90, which are optimized for the NFWc (eq.( 3) with γ = 1.3),Einasto (eq.( 4)), NFW (eq.( 3) with γ = 1.0), and isothermal (eq.( 6)) density profiles, respectively.In the galactic plane, longitudes further than 6 • from the Galactic center are removed from all ROIs larger than R3, because it is not expected a large dark-matter signal in that region and the photon emission is dominated by standard astrophysical sources.Note that the dependence on the dark-matter density profile is sizeable although the use of optimized ROIs mitigates the effect.In the case of isothermal profile (cored profile) the upper bound excludes annihilation cross-section σv UL γγ ∼ 8 × 10 −29 cm 3 /s for m DM = 1 GeV and σv UL γγ ∼ 3 × 10 −27 cm 3 /s for m DM = 500 GeV.Fermi-LAT data constrain dark-matter masses in the range 200 MeV to 500 GeV.
The data from the center of the Milky Way collected by the HESS telescope concern heavier darkmatter masses, starting from 300 GeV and going up to a maximum of 70 TeV and are based on 254 h of live time observation.The analysis is performed for the Einasto density profile and gives an upper limit of σv UL γγ ∼ 4 × 10 −28 cm 3 /s at m χ = 500 GeV and σv UL γγ ∼ 2 × 10 −25 cm 3 /s at m χ = 70 TeV.The search has one optimized ROI, called R1, defined as a circular region centered on the Galactic center, with the galactic plane masked at latitudes lower than 0.3 • .The collaboration also quotes the upper limits for a NFW density profile, obtained for R1 simply by rescaling σv UL γγ obtained with the Einasto dark-matter density profile by the appropriate J-factor.Also for these exclusion limits, the choice of the dark-matter density profile introduces large astrophysical uncertainties.This is known and common to all searches pointing towards the Galactic center, a region where the determination of the spatial distribution of the dark matter is rather uncertain.
The experimental exclusion limits are provided as a function of m DM in two forms.First, in terms of the flux Φ(E γ ) and, secondly, translated into annihilation cross section σv γγ .The latter is computed assuming the γγ annihilation channel only.They depend on the chosen dark-matter density profile, i.e. they are antiproportional to the J-factor.We make use of both information which is stored within the ExpClass module of MadDM, for further details see appendix Appendix C.
The MadDM indirect spectral features command computes all 2 → 2 processes at one-loop order characterized by at least one photon in the final state, i.e. γγ and γX, and confronts it with the abovementioned experimental constraints.Two types of results are shown.On the one hand, we display the cross sections for all individual channels together with the respective experimental limits. 2Note that we neither combine different channels (if their E line γ is equal within 2 The limits for γX are obtained from the experimental ones for γγ by the rescaling σv UL γX = the experimental resolution) nor check the applicability of the experimental analysis (in case of multiple distinguishable lines; see below) for this type of results.
On the other hand, we list all spectral lines as they would occur in the experimental analysis and display the respective fluxes together with the corresponding 95% CL limits for both experiments.To this end, MadDM combines all channels whose peaks are sufficiently close in energy such that they are indistinguishable given the experimental resolution.This is, in particular, relevant for heavy dark matter, i.e. for 4m 2  DM /m 2 X 1.Note that the energy resolution is experiment-specific.Furthermore, (combined) spectral lines lacking a sufficient separation from each other question the applicability of the experimental limitsetting procedure and are flagged accordingly.The minimal separation and combination conditions (in terms of the experimental resolution) can be adjusted by the user.The merging procedure is detailed in appendix Appendix C.
While MadDM uses the Einasto profile by default, the user can freely choose the dark-matter density profile (among the ones listed in Sec. 3) and its parameters in the file maddm card.dat as well as the ROI for the Fermi-LAT analysis, see appendix Appendix B.1 for details.The program automatically computes the corresponding J-factors.However, note that the cross section upper limits from Fermi-LAT are only displayed for the combination of ROI and profile chosen in the analysis.Similarly, the display of limits is omitted outside the energy range for which the analyses have been performed.Finally, we would like to mention that the user can easily implement other experiments e.g. for the computation of projected limits.A template experiment is provided with the code, see appendix Appendix B.2 for details.
On top of the above results, MadDM computes an approximate likelihood from the predicted flux for two of the ROIs of the Fermi-LAT analysis, namely R3 and R16.In Ref. [89], a likelihood profile has been derived from the Fermi-LAT data as a function of the integrated flux in these ROIs.We utilize these results within MadDM.The likelihood functions extend the MadDM experimental module ExpConstraints.

Physics applications
To demonstrate the physics impact, in the following we apply the new feature of MadDM to two BSM scenarios and derive gamma-ray line constraints on their , where E line γ is the position of the line according to eq. ( 9).model parameter space.We consider a simplified topphilic dark-matter model as well as the IDM as illustrative examples.

Simplified top-philic dark-matter model
This scenario is a simplified dark-matter model that supplements the standard model by a top-philic scalar mediator, t, and a Majorana fermion, χ.The former has gauge quantum numbers identical to the right-handed top quark, while the latter is a singlet under the standardmodel gauge interactions.Imposing a Z 2 symmetry under which only χ and t transform oddly, we can stabilize χ for m χ < m t, rendering it a viable dark-matter candidate.
The interactions with the standard model are described by the Lagrangian where D µ denotes the gauge covariant derivative, t the top quark Dirac field and λ χ is the dark-matter coupling.The two masses and λ χ are considered free parameters of the model. 3uch t-channel mediator or charged parent models have been widely studied in the literature.In particular, the top-philic case considered here has been studied in detail in a large range of λ χ [56].Through thermal freeze-out, the model can explain the relic density in a wide range of masses with perturbative couplings.For the slice in parameter space with ∆m = m t − m χ = m t , shown in figure 1, the smallest coupling is required around m χ m t where the annihilation into a pair of top-quarks opens up.For smaller dark-matter masses, this channel becomes Boltzmann suppressed during freeze-out and the three-body and loop-induced annihilation processes tW b and gg, respectively, are dominant [56].
We generate the NLO UFO model with FeynRules [75], FeynArts [76] and NLOCT [77] and compute the loop-induced annihilation into γγ, γZ and γh that occur via triangle and box diagrams involving the scalar mediator and the top-quark in the loop.They involve a total of 14, 14, and 6 diagrams, respectively.The computation of the γγ contribution has been performed in the literature for the first time in [2,3] (see also [56,90]).Our results agree with the ones from [2,56].
The resulting constraints on the dark-matter coupling, λ χ , are shown in figure 1.The green and blue lines denote the limits from Fermi-LAT and HESS, respectively, assuming an Einasto profile.The corresponding shaded bands around them indicate the uncertainty from the density profile.Its lower and upper boundaries mark the limits assuming the NFWc and the isothermal profile, respectively.Notably, the band for the case of Fermi-LAT is considerably smaller as a result of optimized ROIs for the different profiles, see section 3.3.The dot-dashed curve denotes the projection for CTA [91] taken from Ref. [92].
For m χ 240 GeV, the energy of the spectral lines arising from annihilation into γh and γZ are equal within the detector resolution of Fermi-LAT and hence combined for the limit setting.They are combined with γγ above m χ 260 GeV.However, the γh annihilation cross section is p-wave suppressed and hence negligible (the corresponding limit for the individual channel is far outside the displayed window).The limits arising from combinations of channels are drawn as solid lines.
For m χ 260 GeV we display the limits for the individual channels; the dotted and dot-dashed lines correspond to γγ and γZ, respectively.However, in the range 110 GeV m χ 260 GeV, the peaks are not sufficiently separated while their integrated fluxes are of similar size questioning the validity of the limit-setting procedure in the Fermi-LAT analysis (see appendix Ap-pendix C for further details).This is indicated by the missing green band which is only shown in the ranges providing a reliable limit.Note that for m χ 110 GeV it corresponds to the stronger limit (γZ).
Current limits only constrain the region below 110 GeV for the NFWc profile.Note that this region is also challenged by direct detection experiments (see [56]) and -independent on astrophysical parameters -by LHC searches for supersymmetric top partners [93][94][95][96].However, for the considered slice ∆m m t the signal acceptance is subject to large uncertainties and, hence, this region is often blanked out in the experimental limits, see e.g.[93].

Inert doublet model
The IDM supplements the standard model by an additional (inert) Higgs doublet, Φ, that is odd under an exact Z 2 symmetry, stabilizing the lightest of its state.The interactions with the standard model arise from the gauge kinetic terms for Φ and the scalar potential, which reads After electroweak symmetry breaking, the model contains a total of five physical scalar states with masses given by where Imposing m h 125 GeV, five free model parameters remain.We express them as {m H 0 , m A 0 , m H ± , λ L , λ 2 }.For details about the parameters and the IDM, we refer to [97][98][99].
Despite its simplicity, the IDM leads to a rich phenomenology.It provides a viable dark-matter candidate in various regions in parameter space involving all "exceptional" regimes [100] of dark-matter freeze-out -annihilation near a resonance, close to thresholds and coannihilation, see, e.g.[101][102][103][104].
Here, we assume H 0 to be the dark-matter candidate and focus on the region around m H 0 72 GeV, which has been previously studied by some of the authors in [104].In this region, the measured relic density can be explained by annihilation into W W * , ZZ * via the gauge kinetic interaction alone [99,105], where V * denotes an off-shell vector boson.Accordingly, the H 0 -Higgs coupling, λ L , can be arbitrarily small.Interestingly, it is among the two regions that provide a good fit to the Fermi-LAT Galactic center gammaray excess when interpreted as a signal of dark matter (the other one being the resonant region which, however, requires highly tuned parameters) [104].At the same time, it is unchallenged by current limits from the LHC [106] and direct detection -due to the small λ L , diagrams with massive gauge bosons in the loop [107] dominate the direct detection cross section which are around an order of magnitude smaller than the current limit from LZ [108].While upcoming observations of continuum gamma-ray emissions will provide the sensitivity to probe this region [104], in the presence of the yet unexplained excesses in the continuum emissions, an independent probe of the model via gammaline searches is highly desirable.Hence, loop-induced annihilation into γγ and γZ is a smoking-gun signature that allows us to test its dark-matter origin.Again, we generate the NLO UFO model with Feyn-Rules [75], FeynArts [76] and NLOCT [77].The processes H 0 H 0 → γγ and H 0 H 0 → γZ involve 140 and 172 diagrams, respectively.Note that H 0 H 0 → γh is forbidden due to charge-conjugation invariance (i.e. a generalization of Furry's theorem).We have validated our numerical setup using existing results for σv γγ in the literature [20,35] and found agreement within the numerical precision. 4e consider points within the 2σ region from the parameter scan and fit performed in [104], which takes into account constraints from the relic density [109], electroweak precision observables [110,111], new physics searches at LEP-II [112,113], indirect detection searches for continuum γ-ray spectra from dwarf spheroidal galaxies [58] and theoretical requirements of unitarity, perturbativity and vacuum stability computed with 2HDMC [111].The scan was performed with Multinest [114,115] for efficient parameter sampling.
The resulting cross sections for σv γγ and σv γZ are shown in figure 2 together with the corresponding limits from Fermi-LAT [8] and future projections for GAMMA-400 (2 years) as well as a combination of both observations (with 12 and 4 years observational time, respectively) [116].All lines assume the Einasto profile while the shaded band around the Fermi-LAT limit illustrates the uncertainty due to the choice of the profile.The upper and lower boundary of the shaded band corresponds to the limit assuming the isothermal and NFWc profile, respectively.In the γγ channel, the current limits only constrain the considered region for the case of NFWc profile.For the Einasto profile, future observations are expected to provide sensitivity.In the case of an excess in γγ, the confirmation of a linesignal around E γ 43 GeV in γZ would be important to establish the model.However, the sensitivity of the combination of Fermi-LAT and GAMMA-400 (purple, dot-dashed curve in right plot) is still too low by about a factor of two to reach the expected signal originating from σv γZ .

Discussion and conclusions
Gamma-ray lines are smoking-gun signatures of direct dark-matter annihilation into photons in galactic halos.Additionally, they can complement searches for continuum gamma-ray emission from dark matter.This is, in particular, relevant given the presence of unexplained excesses in the continuum gamma-ray flux, like the long-standing Fermi-LAT Galactic center excess.
For electrically neutral dark matter, annihilation into γγ or γX is loop-induced.While the respective cross sections have been computed for a variety of models, with this paper, we present MadDM v3.2 that enables their automated computation for arbitrary darkmatter models implemented in the NLO UFO format.It utilizes an interface to MadLoop, a MG5 aMC tool for computing loop processes, optimized for the involved kinematics of non-relativistic dark-matter annihilation.Furthermore, MadDM computes the integrated fluxes and applies experimental constraints from Fermi-LAT and HESS to the model.Allowing for the computation of the J-factor for a variety of dark-matter density profiles as well as regions of interest, it provides the user with the necessary flexibility for comprehensive dark-matter studies.
We demonstrated the capabilities of the program by applying it to two relevant physics scenarios.First, we studied the constraints from line-signatures on a simplified top-philic dark-matter model.It provides a viable dark-matter candidate over a large range of masses.For dark matter lighter than the top-quark, 2 → 2 annihilation processes at tree-level are suppressed, while loopinduced processes become important and consequently require large dark-matter couplings to match the relic density.This also enhances the strength of the gammaline signals relative to the continuum emission for in- direct detection.However, current line searches from Fermi-LAT only constrain the region of dark-matter masses below roughly 110 GeV.
According to the given experimental energy resolution, MadDM merges the spectral feature of various channels improving the constraining power.For instance, for the channels γγ and γZ, this happens above m DM = 260 GeV in the case of the Fermi-LAT constraints, while their signals could be merged for the entire reported range for the current HESS and projected CTA limits.Nevertheless, the model remains unchallenged by these limits.
As a second scenario, we considered the inert doublet model.While it provides several viable regions in the parameter space, we concentrate on the region 71 GeV m DM 74 GeV that yields the measured relic density without any particular tuning of the involved parameters, providing the required annihilation rate via annihilation into W W * , ZZ * due to the standard-model gauge interactions.Furthermore, this region fits the Fermi-LAT Galactic center excess and is currently not challenged by direct detection or collider searches.We found that current limits from gamma-line observation constrain the scenario for the very cuspy NFWc profile while future observations by GAMMA-400 and Fermi-LAT are expected to probe it for a slightly broader range of dark-matter profiles.An intriguing possibility is the observation of two peaks, from γγ and γZ which are well separated in the considered parame-ter region.However, the γZ signal still appears out-ofreach according to the above-mentioned projections.
While the main focus of this paper is on gammaline signatures, we would like to stress that the current release of MadDM allows for the computation of arbitrary (one-)loop-induced dark-matter annihilation processes, i.e. it is not restricted to photon final states.Annihilation into a pair of gluons is a relevant example that falls into this class.Unlike annihilation into photons, it leads to a continuum spectrum which is confronted with experimental constraints that we introduced in an earlier release of MadDM.
This release can be seen as the first step towards incorporating automated loop-level computations in MadDM.In the light of current data from the cosmic microwave background enabling a relic density measurement below percent level precision, higher-order corrections in the corresponding theoretical predictions become relevant and call for automated computational tools.While MG5 aMC utilizes a suitable framework for this quest, an efficient computation of thermally averaged cross sections at NLO requires further research that is left for future work.Another planned development of MadDM concerns generic one-loop computations of direct detection cross sections.They are of phenomenological relevance for a variety of dark-matter models.In fact, in the two physics scenarios considered in this paper, the leading contribution to dark-matter scattering off nuclei is loop-induced.
As described in section 2, MadLoop requires several additional tools, utilized for the evaluation of loop integrals, namely Ninja [66][67][68], Collier [65], Cut-Tools [69] and IREGI [70].When computing a loopinduced process for the first time, MadDM asks the user to confirm their automatic installation within the launch interface (see below). 6However, the installation can also be performed independently at any time in the MadDM command line: MG5_aMC_v2_9_4$ python2.7 bin/maddm.pyMadDM> install looptools MadDM> quit MG5_aMC_v2_9_4$ The computation of loop processes for indirect detection of dark matter using MadDM calls all the above external packages.Hence they should be appropriately cited.Now everything is set for the user to start computing loops and lines for her/his favorite dark-matter models.

Appendix A.2: Generate command
Once the user has entered the MadDM command line by executing "python2.7 bin/maddm.py",he/she can navigate through the program via the following commands.First, the user has to load a NLO UFO model.In the following, we will consider the model topphilic NLO EW CM UFO as our reference benchmark model.This model is already available in the model repository of MG5 aMC and is automatically downloaded, when typing the appropriate command: With the current version, MadDM has the new feature enabling the automated recognition of all BSM particles that are odd under the Z 2 -symmetry that stabilizes dark matter.For instance, for the model under consideration, after performing the define darkmatter command, MadDM displays: display z2-odd INFO: z2-odd: chi, t1 ...

MadDM> generate indirect_spectral_features MadDM> output my_process_dir
The indirect spectral features command performs the automatic generation of all diagrams for the final states γγ and γX, where X includes the standardmodel Z and h as well as all additional BSM particles, that are kinematically accessible and even under the dark symmetry.Hence, for the model under consideration the generated final states are γγ, γZ and γh only, as t1 is Z 2 -odd.
The command output creates all the necessary code to perform the loop-induced computations.The above commands create the folder structure illustrated in figure 3. The indirect detection directory called Indirect LI line is created, where LI and line stand for loop-induced process and line signal, respectively.Inside this directory, in the Cards folder, all cards for settings are stored including the MadLoopParams.datcard, which serves to set the MadLoop run parameters.The output folder inside my process dir contains a sub-directory for each run performed for the generated process (run 01, run 02, . . . ) which, in turn, contains all outputs of the computation (stored in MadDM results.txtand maddm.out)as well as a copy of the maddm card.datused.Further output is linked to the directory Output Indirect, see the next section for more details.
Instead of computing all possible channels that provide gamma-lines signatures at once the user may want to compute individual channels.This is possible by specifying the final state as an argument of the command.For instance, for the γγ final state, the syntax is MadDM> ... MadDM> generate indirect_spectral_features a a MadDM> output my_process_dir_aa (Note that the photon is denoted by "a" here.)According to the general MadDM syntax (see e.g.[13]) a specific subset of channels could be defined, for instance:

MadDM> ... MadDM> generate indirect_spectral_features a a MadDM> add indirect_spectral_features a z MadDM> output my_process_dir_aa_az
As the computation of loop-induced process can be expensive (in particular for sophisticated models) we recommend specifying the final state if the computation of all channels is not needed.

Appendix A.3: Launch command
The loop-induced computations for a given parameter point are performed via the launch command:

MadDM> launch my_process_dir
This opens the launch interface that allows the user to change settings and model parameters as it is shown in figure 4. As in previous releases of MadDM, these changes can be made in two ways.First, (repeatedly) entering a number 1-5 allows to alternate between the options displayed: here we focus on item 4 ( indirect detection of spectral lines) which is the only one turned on in figure 4. Furthermore, entering 6 or 7 opens the files param card.dat or maddm card.dat,respectively, with a command-line editor (vim by default) to modify the settings.The former file contains all model parameters, while the latter allows the user to change most of the MadDM settings.Many new parameters related to the gamma-line analysis have been introduced with this release.They are described below.Alternatively, the user can enter the command set <parameter> <value> in the launch interface to control the various settings.For our considered benchmark model, we choose the following model parameters as an example: > set mchi 450.> set mst 622.> set lambdachi 0.85 These changes are stored in param card.dat.Additionally, the astrophysical parameters for the gamma-ray line analysis can be specified, for instance: > set vave_indirect_line 7.5e-4 > set profile einasto > set alpha 0.17 > set r_s 20.0 > set nevents 10000 These settings are stored in maddm card.dat.In the first line, the mean average velocity (in units of c) relevant for the considered gamma-ray line observations is set.For the applicability of the Fermi-LAT and HESS constraints implemented in MadDM, a value for the center of the Galaxy has to be entered.Accordingly, if vave indirect line lies outside [5 × 10 −4 ; 10  warning is displayed.Here we choose the corresponding value for the center of the Galaxy from [84] which is chosen by default.Note that vave indirect line is independent of vave indirect cont the latter of which being the mean average velocity relevant for the com-putation of continuum spectra (e.g. for gamma-ray observations in dwarf spheroidal galaxies). 7n the second to fourth lines, the dark-matter density profile is set to Einasto, with α = 0.17 and scale radius of 20 kpc.These are the default parameters.These parameters are in accordance with the one chosen for the cross-section limit-setting in the Fermi-LAT analysis for the R16 ROI.In the last line, the number of generated events is set to 10000, which is the recommended ball-park to obtain a precise determination of the line annihilation cross-sections.
Note that the content of the file maddm card.datcan be adjusted according to the level of detail to which the user needs to control the settings.While the most common settings are contained in maddm card.datby default (in particular the settings mentioned above) the user can switch to the full mode in which an extended set of parameters is displayed that are of interest for specific applications only.This can be achieved by typing MadDM> update to_full directly after entering the launch interface.Note that the maddm card.datfile can be opened in a commandline editor by typing 7. Examples of parameters that can only be changed in the full mode appear in appendix Appendix B. An example for the full content of the maddm card.datassociated with gamma-line signatures can be found in figure 10.
Once the user has completed the settings, the run is launched by pressing enter.The launch interface is closed and the results are displayed on screen.The screen output is shown in figure 5.It contains the following information: -The annihilation cross-section has been computed with the reshuffling method.Note that the MadDM fast mode is disabled for loop-induced processes unless the loop process is described by an effective vertex, in which case the process arises at tree level; -MadDM displays the velocity at which the processes is computed (as set above); -MadDM prints the value of the annihilation crosssection for the aa, az and ah separately, for the scenario called All DM.The latter assumes that the dark-matter particles make up 100% of the required relic abundance, being agnostic on how this is achieved (see [13] for more details); -MadDM displays the corresponding experimental 95% CL upper limits for the selected ROI and density profile (here R16 and Einasto, respectively).
This information is only provided if the user chooses a density profile and ROI that match the experimental ones.Note that this is the case for the default settings.Here the model point is flagged ALLOWED as the theoretical predictions are below experimental bounds.It would be flagged EXCLUDED otherwise.Note that for model points that can be both constrained by Fermi-LAT and HESS, MadDM shows the most constraining limit among the two.-Below the line INFO: *** Line limits from Fermi-LAT 2015, HESS 2018, MadDM displays the results associated with the flux computation.
First, MadDM provides information about the darkmatter density profile.Below the results for the Fermi-LAT and HESS analyses, respectively, are shown.Note that results are only shown if E γ of at least one peak lies within the experimental range of the corresponding experiment (± half of full width at half maximum).This is the case for both experiments for the chosen benchmark model.For each experiment, MadDM prints the J-factor for the chosen density profile and the corresponding ROI (see appendix Appendix B.1 for details) and shows the energy range of the experiment.Finally, it displays the predicted flux and the corresponding 95% CL flux upper limits together with the flag ALLOWED (or EXCLUDED).For the considered parameter point, the spectral lines for the three final states appear at the same energy (within the resolution of both experiments).Hence, they have been merged to a single one (see appendix Appendix C for details).In the case of several distinguishable peaks, the corresponding results for each peak would be listed.Note that due to the combination of channels the constraints from the flux computation can be stronger.
The same output information is stored in the file my process dir/output/run 01/MadDM results.txt,see figure

MadDM> add indirect_detection g g MadDM> add indirect_spectral_features MadDM> output my_process_dir_complete
Here, we perform the relic density computation without specifying a potential coannihilator, because we consider a benchmark point in which chi and t1 have a relatively large mass splitting rendering coannihilation to be irrelevant.In general, coannihilators should be considered, see [13] for details.The first command add indirect detection generates the code necessary to perform tree-level indirect detection computations.The second command add indirect detection g g generates the code to compute the process into pairs of gluons, which is loop-induced but produces a continuum energy spectrum of cosmic messengers.All continuum spectra are added and the ones for gamma rays are automatically confronted with limits from Fermi-LAT observations of dwarf spheroidal galaxies.The third command add indirect spectral features adds the computation of spectral lines as discussed above.The combination of these three commands represents all possible ways of generating gamma rays from the annihilation of dark matter for this specific top-philic model.The execution of these commands creates the folder structure displayed in figure 6, inside the directory my process dir complete.In the example above, there are the following three folders corresponding to the three classes of processes generated: 1. Indirect tree cont: this folder contains the treelevel indirect detection processes leading to continuum spectra of cosmic messengers.These processes are run with the dark-matter halo velocity vave indirect cont; 2. Indirect LI line: this folder contains the gammaray line loop-induced processes.This run is pursued with the velocity vave indirect line; 3. Indirect LI cont contains loop-induced processes leading to continuum spectra (here the gg final state).This run is pursued with vave indirect cont.The photon energy spectrum generated from gg is merged together with the gamma-ray spectrum used for the continuum photon analysis based on the Fermi-LAT data.As in previous versions the spectra can be generated using Pythia 8 or be taken from PPPC4DMID (in the case of standard-model final states).
Each indirect detection folder is run separately by MadDM, one after the other.The different relevant settings are contained in the Cards folder of each indirect detection directory.The user can change settings in the launch interface, as described below.Note that there is also a fourth possibility if the annihilation into photons proceeds via an effective vertex in the model.In this case, MadDM generates a directory Indirect tree line that indicates that the process occurs at tree level but provides a line signature processed by the corresponding analysis pipeline.We will not comment on this case further as the behavior is similar to the example being discussed.Note, however, that in this case, the MadDM fast mode is available.
After performing the command Other options are flux source and sigmav.For more information about these options we refer the reader to [12,13].Settings for the run are done by typing > set mchi 450.> set mst 622.> set lambdachi 0.85 > set vave_indirect_cont 2.e-5 > set vave_indirect_line 7.5e-4 > set profile einasto > set alpha 0.17 > set r_s 20.0 > set nevents 10000 These are almost the same settings described in appendix Appendix A.3, the one addition concerns the dark-matter velocity for the computation of continuum fluxes (vave indirect cont), which is set to a value compatible with the one measured in dwarf spheroidal galaxies.The run is performed in the MadDM precise mode.The number of generated events, nevents, is common to all indirect detection observables, i.e. for each class of processes, 10000 events are generated.Note that this choice should be adjusted to the user's needs.For the computation of smooth energy spectra for the continuum emission (in particular in the tails) a significantly larger number of events may be required.
The final output screen of MadDM is displayed in figure 8, and is constituted by three main blocks: -The first block concerns the relic density.The chosen benchmark point predicts the measured relic density within experimental errors, Ωh 2 = 0.12.Furthermore, MadDM displays the temperature parameter at freeze-out x f , the total thermally averaged annihilation cross section evaluated at the point x f and the dark-matter fraction of the considered candidate ξ = Ωh 2 theo /Ωh 2 Planck used for the rescaling of (in)direct detection limits.See [13] for details.With the current version, MadDM also lists all annihilation channels and their contribution to the thermally averaged annihilation cross-section in percent.
-The second block shows the output for the (treelevel and loop-induced) processes leading to continuum fluxes.Since the relic density has been computed and the point is (slightly) under-abundant, the results are shown for both the Thermal and All DM scenario.In the former case, the limits are rescaled according to the dark-matter fraction ξ, while the latter assumes the measured relic density and hence no rescaling.For the benchmark point under consideration ξ = 0.99, and hence no significant difference between the two cases.The corresponding exclusion limits at 95% CL are given, for each final state separately.These are obtained by a simple comparison with the tabulated exclusion limits in the MadDM experimental module.The line DM DM > all takes into account all computed channels (here t t and gg) and performs a likelihood analysis to provide the 95% CL exclusion limit from the Fermi-LAT data from dwarf spheroidal galaxies.Note that to perform the likelihood analysis, the gamma-ray energy spectra inside Output Indirect LI cont are combined with those inside Output Indirect tree cont.Their sum is given directly in the folder output/run 01/, as illustrated in figure 6. -In the third block MadDM displays the annihilation cross sections and fluxes for the gamma-line signals, as already discussed in appendix Appendix A.3.Here, however, since the relic density variable has been computed, both dark-matter scenarios, Thermal and All DM, appear, see the previous bullet point for details.
Finally, we recall that auto-completion is available to conveniently find commands and parameter names.Note further that all commands and settings described in the previous sections can be performed via a script whose path is passed as an argument with the execution of maddm.py,see [13] for details.follows: A two-dimensional projection of the integration region A at the Galactic center is shown in figure 9.It visualizes the parameters of the two possible masks considered in MadDM.On the one hand, it allows one to mask the Galactic plane up to a latitude β and down to a Galactic longitude λ.On the other hand, it allows for a circular mask around the Galactic center with an opening angle α 1 ∈ (0, α 2 ).The outer opening angle of the ROI is denoted by α 2 ∈ (0, π).We use galactic coordinates for the computation: longitude l ∈ (−π, π) and latitude b ∈ (−π/2, π/2).For definiteness the galactic coordinates are: Here, we assume that the integration cone is centered around the x axis.The mask parametrization introduced above allows one to describe a variety of ROIs, in particular, those considered by the Fermi-LAT and HESS collaborations.For instance, the R1 ROI used by the latter is obtained for α 1 = 0 • , α 2 = 1 • , λ = 0 • and β = 0.3 • .We have validated the computation of J-factors using the ones reported by the collaborations [8,9] and elsewhere in the literature [79] and find agreement below the percent level. 8e have included four dark-matter density profiles, specified in eqs.(3) (gNFW), (4) (Einasto), (5) (Burkert) and ( 6) (isothermal).These dark-matter profiles are all normalized at the Sun position, the normalization is given by the value of the scale radius r s and by the scale density ρ s .The latter is not a free parameter and is automatically computed once the user has fixed ρ , R and r s .The dark-matter density profile and the scale radius can be changed by the user by editing the file maddm card.dat or via the launch interface, for example:

MadDM> set profile nfw MadDM> set r_s 25
These commands change the dark-matter density profile from its default Einasto to NFW with (which corresponds to gNFW with γ = 1) and the scale radius to 25 kpc (the default value is r s = 20 kpc).The default values for the measured energy density at the Sun position and the Sun's distance from the Galctic center are ρ = 0.4 GeV/cm 3 and R = 8.5 kpc, respectively, following the choice of the Fermi-LAT collaboration.Note that these parameters are not contained in the standard maddm card.dat.The change of these parameters requires enabling the full mode via typing update to full directly after entering the launch interface, see appendix Appendix A.3.An example of the file maddm card.dat in full mode is given in figure 10.However, in the given example, MadDM computes the value of ρ s that corresponds to the default ρ , R and to r s = 25 kpc.
For NFWg, the inner slope parameter γ is free.Users can change it by typing:

MadDM> set profile gnfw MadDM> set gamma 1.3
This choice corresponds to NFWc density profile.Considering the Einasto, the slope α is a free parameter and can be set by MadDM> set profile einasto MadDM> set alpha 0.2 which change α from its default value 0.17 to 0.2.

Appendix B.2: Details about experimental data
In general, the user can set the ROI for a given experiment.Concretely, the experimental bounds are derived for a ROI which is optimized for a certain dark-matter density profile.In the case of Fermi-LAT there are four ROIs: 1. R3 associated with NFWc 2. R16 associated with Einasto

R41 associated with NFW 4. R90 associated with isothermal
For HESS, only one ROI is defined, namely R1.While the user is free to choose any of the four ROIs of the Fermi-LAT analysis by the variable roi_fermi_2015 in the file maddm_card.dat,we recommend choosing the ROI associated with the dark-matter density profile considered.This can automatically be achieved by setting the above parameter to default. 9Note that if a ROI is chosen that is not optimized for the considered profile, MadDM raises a warning.If the user only chooses the ROI without specifying the dark-matter density profile, MadDM automatically chooses the profile for which it was optimized.Note further that the approximate Fermi-LAT likelihood and p-value computations are only performed if the user chooses the R3 or R16 ROI.
The experimental data is stored in the directory MG5 aMC v2 9 4/PLUGIN/maddm/ExpData.It contains the respective data files: Fermi_2015_lines_R3_NFWcontracted.dat Fermi_2015_lines_R16_Einasto.dat Fermi_2015_lines_R41_NFW.dat Fermi_2015_lines_R90_Isothermal.dat for the Fermi-LAT analyses and HESS_2018_lines_R1_Einasto.dat for the HESS analysis.Each file contains the parameters of the density profile, the corresponding J-factor, the ROI and the 95% CL upper limits on σv γγ and the flux in the 120 (60) energy bins considered by Fermi-LAT (HESS).
In addition to the implemented constraints from Fermi-LAT and HESS, we allow the user to consider further experimental data in the analysis pipeline.To this end, we have included a template experiment.It can, for instance, be used to compute the projected limits for future observations.The associated parameters 9 For the Burkert profile, the R90 is chosen by default.are among the extended set of parameters accessible in maddm_card.dat in the full mode (see appendix Appendix A.3).The respective block is shown in figure 10 below the line Template line experiment.The relevant parameters are: name, arxiv: specify the name and (if applicable) the arXiv number of the analysis (optional); roi, profile: it is possible to specify the ROI (i.e.α 2 in degrees) and the dark-matter density profile associated with it (the user can choose between the same profiles as for the main profile parameter discussed above); r s, gamma, alpha: these are the profile parameters to be specified for the default dark-matter density profile chosen as in the point above; mask parameters: allow the user to specify a mask, according to figure 9, where mask latitude is λ, mask longitude is β and mask inner angle is α 1 ; they are all expressed in degrees; energy resolution: the relative energy resolution of the experiment in per-cent; taken to be constant over the detection range; detection range parameters: allow the user to specify the minimum (min) and maximum (max) energy range of the experiment (in GeV); constraints file: additionally, the user can provide the (projected) constraints of the experiment via a data file in the same format as e.g.Fermi 2015 lines R90 Isothermal.dat.The file should be included in the directory MG5 aMC v2 9 4/PLUGIN/maddm/ExpData of MadDM.It should contain three columns: the darkmatter mass (in GeV), σv γγ upper limits (in cm 3 s −1 ) and Φ γγ upper limits (in cm −2 s −1 ).
To enable the template experiment, the parameter toggle template line experiment has be to set to on.In this case, MadDM considers it in the same way as the constraints from Fermi-LAT and HESS.Note that the dark-matter density profile associated with the template experiment is supposed to be the one for which the cross section upper limit has been derived (if provided by the user).The density profile considered for the J-factor computation within MadDM remains the one specified through the parameter profile in maddm card.dat.However, if the two profiles differ, a warning is raised.Note that it is also possible to switch the Fermi-LAT and HESS analysis on/off.When switched off the respective experiment is not considered in the analysis pipeline.

Appendix C: Merging algorithm for peaks
In the presence of multiple spectral lines, the application of experimental constraints (performed under the assumption of a single gamma-ray line signal) requires particular attention.In the following, we provide details on the implemented merging algorithm for multiple spectral lines that are close in energy.Furthermore, we discuss the considered criteria for the applicability of the experimental analyses in this case.Due to the non-relativistic nature of cold dark matter, the width of the gamma-ray line signal is expected to be very small.In particular, it is assumed to be small compared to the experimental energy resolution.In the experimental measurement, the sharp peak, eq. ( 7), is hence smeared.We approximate it by a Gaussian, where σ corresponds to the experimental resolution, related to the full width at half maximum by FWHM = 2 √ 2 ln 2 σ.The resolution of Fermi-LAT is reported in [8].We consider the EDISP3 energy resolution, associated with the best energy reconstruction estimator.It ranges from ∼ 10% at 300 MeV down to 5.5% at 500 GeV.Note that the width of each bin corresponds to the energy resolution (68% containment) for that bin following EDISP3.The HESS telescope has a resolution in the energy of 10% above 300 GeV [9].
The merging algorithm follows an iterative process.In each iteration, it considers all possible pairs of spectral lines (called peaks in the following) and selects the pair that is closest in energy, i.e. for which is minimal, where a, b index the peaks and FWHM a is the FWHM of peak a.If s ab < 1, the Gaussian signals, eq.(C.5), of the peaks are summed and considered as a single (merged) peak.The iterative process is pursued with the set of merged and remaining un-merged peaks until no merging is possible anymore (under the above criterion).An illustrative example of a merging of two peaks is shown in figure 11.
The set of peaks obtained after merging represents the signal as potentially being observed by the experiment.However, if the signal still contains multiple peaks the application of the experimental constraints to each of them is, in general, questionable as the limit-setting procedure has been performed under the hypothesis of a single peak.Similar issues can arise if a merged peak has become too broad.We address these concerns by introducing the following criteria the peaks have to satisfy.First, we require the FWHM of any merged peak to be smaller than 1.5 times the FWHM corresponding to the experimental resolution at its peak energy.If this criterion is not met, the peak is flagged with an error (error type 1).Secondly, we require a minimal separation between the (merged) peaks.The minimal separation in units of the FWHM is controlled by the parameter n_fwhm_separation whose default value is set to 5. Peaks that do not satisfy the separation criterion are flagged with the error type 2. Note that the separation criterion is bypassed for a large hierarchy between the corresponding fluxes of the peaks, i.e. if the ratio of fluxes (each normalized to the upper limit at the corresponding energy) is larger than the variable peak_height_factor (with default value 10).In this case, only the peak with the smaller flux is flagged with the error 2.

Fig. 1
Fig.1Parameter space for the top-philic model, including the constraints for loop-induced dark-matter annihilation into monochromatic photons.The black solid curve (and shaded band around it) shows the dark-matter coupling that yields the measured relic density, assuming Ωh 2 = 0.1199 (with a 10% uncertainty on the theory prediction), taken from[56].The green (blue) curve denotes the Fermi-LAT (HESS) limits for the Einasto profile.The lower and upper boundaries of the shaded bands indicate the corresponding limits for the NFWc and Isothermal profile, respectively.The dotted and dot-dashed curves show the upper limits for the individual final states γγ and γZ, respectively, in the range where the signals are not combined (see text for details).The purple curve shows projections for CTA assuming the Einasto profile.

Fig. 2
Fig.2Cross section for the loop-induced dark-matter annihilation into γγ (left panel) and γZ (right panel) for the allowed IDM parameter points explaining all of dark matter (Ωh 2 = 0.12) in the range 70 ∼ 75 GeV.Parameter points in the 1 and 2σ regions are drawn in dark brown and red, respectively.The green solid curve denotes the Fermi-LAT limit for the Einasto profile.The lower and upper boundary of the green shaded band shows the respective limits for NFWc and isothermal profile, respectively.The blue dotted and purple dot-dashed curves show projections for GAMMA-400 (2 years) and the combination of Fermi-LAT (12 years) and GAMMA-400 (4 years), respectively, both assuming the Einasto profile.
7 ./bin/maddm.pyMadDM> import model topphilic_NLO_EW_CM_UFO MadDM> define darkmatter chi where we have specified the dark-matter candidate (chi) in the last line.The phenomenology of this model has been considered in Sec.4.1.It has only two BSM particles.Besides dark matter, it contains a coloured scalar mediator called t1.This model is a simplified version of the more general DMSimpt framework [52, 55] which describes very generic t-channel models coupling to all standard-model quarks.

Appendix A. 4 :
A complete example While we concentrated on the gamma-line features in the previous examples we finally display a complete example involving the computation of the freeze-out relic density, of continuum fluxes (at tree level and loopinduced) and gamma-line fluxes.The corresponding commands are: MG5_aMC_v2_9_4$ python2.7 ./bin/maddm.pyMadDM> import model topphilic_NLO_EW_CM_UFO MadDM> define darkmatter chi MadDM> generate relic_density MadDM> add indirect_detection

Fig. 5
Fig. 5 Example of output provided by MadDM in gamma-ray line computations after the successful run of MadLoop.The numbers quoted correspond to the benchmark model point taken as an example in appendix Appendix A.3.

Fig. 6
Fig.6Example of the indirect detection folders and files of MadDM.Left: Structure of the user folder my process dir complete that contains all relevant setting cards inside the individual indirect folders (Indirect LI line, Indirect LI cont and Indirect tree cont), and the output files in output/.Right: Zoomed view of the output directory, where the main results are stored.The file MadDM results.txtrecaps the value of all observables computed by the user.The Output Indirect LI cont, Output Indirect LI line and Output Indirect tree cont contain the indirect detection lhe event files for loop-induced processes gg, aX and tree level, respectively.

the
MadDM interface is displayed, as shown in figure 7. The relic density, indirect detection observables for continuum emission and line signatures are turned on, while the other modules are off.For the contin-uum emission module the mode flux earth is chosen.

Fig. 7
Fig. 7 Example of the launch interface after performing the launch command in MadDM for both tree level and loop-induced computations.The relic density, standard indirect detection and indirect detection spectral feature calculations are turned on, while direct detection and the scan capability with MultiNest are switched off.

Fig. 8
Fig. 8 Example of output provided by MadDM when running relic density and all indirect detection computations.The numbers quoted correspond to the benchmark model point taken as an example.

↵ 1 ↵ 2 2 2 Fig. 9
Fig.92D projection at the Galactic center of the integration region A defined in the J-factor integral.

# 3 #Fig. 10
Fig. 10 Example of the content of maddm card.datassociated with gamma-line signatures including the extended set of parameters (full mode).

Fig. 11
Fig. 11 Graphical example of the merging algorithm: on the left the γγ and γZ peaks are shown separately, while on the right the obtained merged peak is shown.