Issue 
A&A
Volume 633, January 2020



Article Number  A138  
Number of page(s)  12  
Section  Galactic structure, stellar clusters and populations  
DOI  https://doi.org/10.1051/00046361/201936632  
Published online  22 January 2020 
Search for Galactic Pevatron candidates in a population of unidentified γray sources
Institut für Physik, HumboldtUniversität zu Berlin, Newtonstr. 15, 12489 Berlin, Germany
email: spengler@physik.huberlin.de
Received:
4
September
2019
Accepted:
10
December
2019
Aims. A list of Pevatron candidates is presented to enable deeper observations and dedicated analyses.
Methods. Lower limits on the energy cutoff for unidentified γray sources detected in the High Energy Stereoscopic System (HESS) Galactic plane survey were derived. Additional public data from the Very Energetic Radiation Imaging Telescope Array System, HESS, and Milagro experiments were used for MGRO J1908+06 to confirm the limit derived from the HESS Galactic plane survey data and to enable further conclusions on the presence of spectral breaks.
Results. Five Pevatron candidates are identified in the HESS Galactic plane survey. The cutoff of the γray spectrum for these sources is larger than 20 TeV at 90% confidence level. The γray sources MGRO J1908+06 and HESS J1641−463, found to be Pevatron candidates in the analysis of the HESS Galactic plane survey catalog, have already been discussed as Pevatron candidates. For MGRO J1908+06, the lower limit on the γray energy cutoff is 30 TeV at 90% confidence level. This is a factor of almost two larger than previous results. Additionally, a break in the γray spectrum at energies between 1 TeV and 10 TeV with an index change ΔΓ > 0.5 can be excluded at 90% confidence level for MGRO J1908+06. The energy cutoff of accelerated particles is larger than 100 TeV at 90% confidence level in a hadronic scenario for all five Pevatron candidates. A hadronic scenario is plausible for at least three of the Pevatron candidates, based on the presence of nearby molecular clouds and supernova remnants.
Key words: astroparticle physics / gammarays: general / methods: statistical
© ESO 2020
1. Introduction
Theoretical models for the origin of cosmic rays (CRs) consider young supernova remnants (SNRs) as Pevatrons, that is, as hadron accelerators for energies up to the “knee” of the cosmicray spectrum measured at around 3 PeV (Baade & Zwicky 1934; Hillas 2005). The detection of a neutral pion decay signature in γray data from SNR observations confirmed that SNRs accelerate hadronic particles (see e.g., Ackermann et al. 2013; Jogler & Funk 2016). However, the maximal energy to which SNRs can accelerate hadrons is not determined and, in particular, whether or not SNRs are sources of PeV CRs is unclear.
The association of a γray source with a Pevatron can be tested with a measurement of the γray spectrum above TeV energies. An exponential cutoff E_{cut, γ} in the γray spectrum at a few hundred TeV is expected when CRs are accelerated up to PeV energies and produce γrays in interactions with ambient material (Gabici & Aharonian 2018). The Pevatron nature of a γray source is in turn constrained when an exponential cutoff at an energy well below 100 TeV is detected in the observed γray spectrum of the source. For example, the recent detection of an exponential cutoff at E_{cut, γ} = 3.5 TeV in the γray spectrum of Cassiopeia A rules out simple Pevatron models for this young Galactic SNR (Ahnen et al. 2017).
The most constraining lower limits on the energy cutoff E_{cut, h} of Galactic hadron accelerators are currently derived from the diffuse γray emission in the vicinity of the radio source Sagittarius A* (SgrA*) and for the γray source HESS J1641−463. Assuming that hadronic particles generate the γray emission from these sources, E_{cut, h} > 600 TeV and E_{cut, h} > 100 TeV are inferred at 90% and 99% confidence level (CL) respectively (Abramowski et al. 2016, 2014). More than 100 isolated sources of TeV γrays are currently known. However, no systematic search for Pevatron candidates among the known γray sources has yet been presented.
This paper is organized as follows. Section 2 outlines the general Pevatron search method followed in this work. The later analysis relies in large parts on the public catalog of the High Energy Stereoscopic System (HESS) Galactic plane survey (HGPS). This catalog is introduced in Sect. 3. A method to constrain energy cutoffs in measurements of powerlaw energy spectra is discussed in Sect. 4. The method is argued to be applicable to the HGPS catalog introduced above. A search for Pevatron candidates in the HGPS catalog is summarized in Sect. 5. For one Pevatron candidate identified in the previous analysis of HGPS data, additional public data are available from the HESS, Very Energetic Radiation Imaging Telescope Array System (VERITAS), and Milagro experiments. A combined analysis of these data is presented in Sect. 6. Finally, Sects. 7 and 8 discuss the results and summarize the conclusions.
2. Pevatron candidate search method
The detection of TeV γrays from an isolated region of the sky can indicate cosmic sites where particles are accelerated to energies of greater than 1 TeV. However, the identification of the primary particle type is nontrivial when only TeV γray data are available. Both leptonic and hadronic scenarios are typically conceivable. Electrons are accelerated in the leptonic scenario and generate γrays via the inverse Compton upscattering of lowenergy photons. This scenario is not relevant in the following search for hadronic Pevatrons as origins of Galactic CRs. However, the unambiguous exclusion of the leptonic scenario for a given γray source is very challenging. For example, as recently discussed, even the production of γrays with an energy exceeding 100 TeV is possible in a leptonic scenario (Amenomori et al. 2019). The nondetection of a γray energy cutoff at energies below few 100 TeV is therefore not considered as a sufficient criterium for the identification of a hadronic Pevatron. Instead of an unambiguous identification of hadronic Pevatrons, the following search aims to find Pevatron candidates, that is, γray sources where the hadronic Pevatron scenario cannot be excluded. The resulting list of Pevatron candidates allows dedicated multimessenger analyses with reduced statistical trial factors, involving, for example, searches for highenergy neutrinos as fingerprints of hadronic interactions.
Many known Galactic γray sources are identified as pulsar wind nebulae (PWNe) for which the γray emission is modeled in a leptonic scenario (Gaensler & Slane 2006). Another large class of γray sources has yet to be identified, with sources of radiation in other energy ranges. The search for new Pevatron candidates presented here is restricted to these unidentified γray sources in the Galactic disk. The motivation for the restriction to unidentified γray sources is at least twofold. The exclusion of PWNe removes a large class of γray sources for which the leptonic scenario is very likely. Additionally, Gabici & Aharonian (2007) discuss that the time period where SNRs act as Pevatrons might be short, that is, less than 1000 years. Easier than the direct detection of Pevatrons might therefore be the indirect identification via the detection of delayed γray emission from molecular clouds in the vicinity of SNRs. Gabici & Aharonian (2007) discuss this scenario and suggest that unidentified γray sources can be associated with molecular clouds illuminated by nearby SNRs.
Following Kelner et al. (2006), the spectrum of γrays generated in the vicinity of a Pevatron in interactions of accelerated CRs with ambient material is modeled by
Here, is the flux normalization at energy E_{0} and the parent CR spectrum has a cutoff energy of E_{cut, h} = 1/λ_{h} (Kelner et al. 2006). Models for diffuse shock acceleration predict a powerlaw index Γ ≈ 2 (see e.g., Hinton & Hofmann 2009). Following this model, TeV γray sources that are selected as Pevatron candidates must have a spectrum at TeV energies which is compatible with a powerlaw with index Γ = 2. Additionally, the inferred lower limit on the energy cutoff E_{cut, h} in the hadronic model given by Eq. (1) must be at least O(100 TeV), that is, of the same order of magnitude as the current best constraints derived for HESS J1641−463 and the γray emission in the vicinity of SgrA*.
To further constrain a leptonic scenario for the generation of γrays, the presence of a break in the TeV γray spectrum where the powerlaw index changes can be excluded. For example, in a leptonic scenario in which inverse Compton losses dominate over synchrotron cooling, a Klein–Nishina (KN) break is expected in the γray spectrum at an energy E_{KN} which depends on the target photon field. For cosmic microwave background (CMB) target photons, the break is expected at a few hundred TeV. For infrared target photons or optical starlight, the energy break is expected around 10 TeV and 30 GeV respectively (Hinton & Hofmann 2009). The γray spectrum at energies far away from E_{KN} can be described by a broken powerlaw with E_{break} = E_{KN}
It is expected that ΔΓ = Γ for the index change at the KN break (see e.g., Hinton & Hofmann 2009). The selection of sources with a powerlaw index around 2 at TeV energies excludes leptonic scenarios where the inverse Compton target photons are of a much higher energy than CMB photons, for example optical starlight, and the KN break must occur below TeV energies. However, leptonic scenarios with CMB target photons cannot typically be constrained with TeV γray data.
The selection of a γray source as a Pevatron candidate in the following search is not sufficient to identify the source with a Pevatron. Furthermore, leptonic emission scenarios cannot be ruled out and the lower limits on E_{cut, h} are not sufficient to argue for the acceleration of hadrons to PeV energies. Additionally, the selection of a γray source as a Pevatron is not necessary for the association with a Pevatron. For example, the γray source Cassiopeia A would not pass the selection due to the measured TeV cutoff discussed in Sect. 1. However, Cassiopeia A is still being discussed as a Pevatron candidate in multizone models (Zhang & Liu 2019).
Despite the Pevatron candidate selection being neither necessary nor sufficient for the association of a γray source with a Pevatron, selected sources can be a particularly interesting target for further observations. For example, targeted signal searches with instruments sensitive to highenergy neutrinos towards a few selected Pevatron candidates can be performed.
3. HESS Galactic plane survey data and quality selection
The currently largest catalog of Galactic TeV γray sources is based on the HGPS. This dataset, as far as is relevant for the data analysis below, and additional quality selection criteria are discussed in this section.
3.1. HESS Galactic plane survey data
The HGPS comprises data acquired from 2004 to 2013 with the HESS I (Aharonian et al. 2006b) array of Cherenkov telescopes in Namibia. Large parts of the Milky Way plane were observed to detect γray sources in a nominal energy range of 0.2−100 TeV with an angular resolution better than 0.1°. Observations of the Milky Way plane were performed in a scanning mode with individual observation runs of typically 28 min towards different sky directions. In total, 78 γray sources were detected in the plane survey. The majority (47) of the detected sources are unidentified, that is, no firm identification is currently possible with objects detected in other energy ranges.
The HGPS catalog is based on a generalpurpose analysis without special adoption for each individual source. Systematic effects of this analysis cannot be ruled out, in particular for sky regions with multiple γray sources.
An adaptive ring background algorithm is used in the analysis of data from the HGPS. This algorithm estimates the number of signal (N_{ON}) and background (N_{OFF}) events and enables a calculation of the excess events N_{γ} = N_{ON} − αN_{OFF}. The acceptance normalization factor α is obtained from independent observations of different sky regions void of γray sources. Energydependent information on N_{γ}, N_{ON}, N_{OFF} and α is not part of the public HGPS catalog. However, spectral flux data for each γray source detected in the HGPS are typically available in six logarithmic energy bins. The bin centers, E_{i}, i = 0, …, N − 1, of spectral points are related via E_{i} = ξ^{i}E_{min} with a binning factor ξ where E_{min} is the center of the lowest energy bin. Statistical errors on the γray flux are given as asymmetric 68% confidence level (CL) intervals and originate from the Poisson counting error on N_{ON} and N_{OFF}. Other error sources, for example an error on the acceptance normalization, are not included in the statistical error on the γray flux. A complete description of the HGPS dataset can be found in Abdalla et al. (2018).
3.2. Data selection
A special quality selection is applied in this analysis of the public HGPS data. The energy binwidth for spectral points is required to be much larger than the instrumental energy resolution to avoid a statistical correlation between spectral points. Following Aharonian et al. (2006b), the energy resolution of the HGPS analysis is assumed to be 15% or better. Only sources in the HGPS with a spectral energy binning factor of ξ > 1.5 are selected.
In the HGPS catalog, data on the error of spectral points are given as a lower (σ_{low}) and an upper (σ_{high}) error. The average error σ = 1/2 (σ_{low} + σ_{high}) is used in the analysis. Spectral points with an error asymmetry η = max(σ − σ_{low}, σ − σ_{high})/σ > 10% are discarded from the analysis. Also discarded are spectral points ϕ which are detected at low significance level S = ϕ/σ < 1.5. Sources in the HGPS catalog are not analyzed when less than five spectral points remain after the quality selection. A total of 25 out of 78 sources were selected for further analysis. A summary of the selected HGPS data is given in Table 1.
Characterization of the HGPS spectral data for the 25 selected γray sources.
4. Constraints on the energy cutoff of a powerlaw spectrum
A statistical hypothesis test can be used to search for an exponential cutoff in a powerlaw spectrum. A confidence interval for the energy cutoff can be constructed via the inversion of the test. However, special care must be taken in the analysis of public data from the HGPS because the number of available spectral points per γray source is typically small. The small sample size questions the application of frequently applied asymptotic results for the distribution of a test statistic when the null hypothesis is true. Additionally, systematic effects on the acceptance normalization cannot be ruled out and the catalog data are based on a generalpurpose analysis. Differences between the estimated and the true acceptance normalization due to instrumental effects could, for some observation runs, lead to outliers in the number of excess events, which are not modeled with the published statistical error. Additional systematic effects can originate from the confusion of γray sources in regions with complex morphology.
This section contains a description of statistical methods to constrain energy cutoffs in powerlaw spectra. The application of the likelihood ratio (LR) test and the Ftest for the presence of cutoffs in powerlaw spectra are discussed in Sect. 4.1. In Sect. 4.2, it is argued that the Ftest is more robust against potential violations of the normal error model of the HGPS catalog data than the LRtest. A method to derive a lower limit on the energy cutoff of a powerlaw γray energy spectrum is discussed in Sect. 4.3.
4.1. Models and hypotheses tests
Consider a powerlaw model for the γray flux ϕ(E) at energy E with an exponential cutoff at E_{cut, γ} = 1/λ_{γ}. The model can be parameterized by
where again Γ is the spectral index and ϕ_{0} the flux normalization at energy E_{0}. The parameter E_{0} is assumed to be fixed to 1 TeV in the following. Also fixed is the parameter β which controls the shape of the exponential cutoff. The case β = 1, for example, is frequently discussed in the literature (see e.g., Ahnen et al. 2017). The hadronic model given by Eq. (1) reduces to Eq. (3) with β = 0.5 and λ_{γ} = 16λ_{h}. The cutoff parameter λ_{γ} is constrained to λ_{γ} > 0 because a γray flux suppression for energies above E_{cut, γ} is expected for Pevatrons. The model given by Eq. (3) has three free parameters θ_{1} = (ϕ_{0}, Γ, λ_{γ}). The question remains as to whether the powerlaw model with exponential cutoff (Eq. (3)) gives a better description of a dataset than a powerlaw model parameterized by
This model has two free parameters θ_{0} = (ϕ_{0}, Γ) while E_{0} is again fixed to 1 TeV. To select the bestfitting model, the null hypothesis H_{0}, that is, the absence of an exponential cutoff (λ_{γ} = 0), can be tested against the physically constrained alternative hypothesis H_{1}, λ_{γ} > 0, by means of a hypothesis test at a given confidence level (CL).
In the following ϕ_{i}, i = 0, …, N − 1 is used to denote binned γray flux measurements with bincenters at energies E_{i}. All N measurements of ϕ_{i} are assumed to be independent. In practice, the independence of the ϕ_{i} can be expected when the energy binwidth for the calculation of ϕ_{i} is much larger than the instrumental energy resolution. Further, σ_{i} is used to represent the measurement errors corresponding to the γray fluxes ϕ_{i}. The measurements ϕ_{i} are, in a first step, assumed to be normally distributed around a hypothetical true value . A model ϕ(E_{i}θ) with parameters θ predicts the true flux . Two frequently used tests for a powerlaw hypothesis (Eq. (4)) against a model with exponential cutoff (Eq. (3)) are presented in the following. Subsequently, the robustness of the tests against deviations from the assumed normal error model is discussed.
4.1.1. Likelihood ratio test
Because the flux measurements are assumed to be independent, the likelihood function for the parameters θ of a given model factorizes to . Let be the maximum likelihood estimates for the model parameters under the hypotheses H_{k} with k ∈ {0, 1}. Consider the LRtest with the statistic where
for the hypotheses H_{k}. The null hypothesis is discarded at confidence level CL if −2lnΛ > Λ_{crit}(CL). The critical value Λ_{crit} can be calculated using the results from Chernoff (1954) because the two models given by Eqs. (4) and (3) are nested, that is, Eq. (3) reduces to Eq. (4) when λ_{γ} → 0. When denotes the inverse cumulative density function of a χ^{2}distributed random variable with one degree of freedom, the critical value is given by ; see Appendix A. The test will have the expected false positive error rate 1 − CL when the normal error model holds and the sample size N is large such that when H_{0} is true.
4.1.2. Ftest
The Ftest is, in the case of N independent measurements and two nested models with two and three parameters, based on the test statistic
The null hypothesis λ_{γ} = 0 is discarded in favor of the alternative hypothesis λ_{γ} > 0 at confidence level CL when F > F_{crit}(CL). Similar to the LRtest, the critical value is given by
where is the inverse cumulative density function of a Fdistributed random variable with 1 (numerator) and N − 3 (denominator) degrees of freedom. The Ftest is exact for small samples N when linear models are being compared. For nonlinear models, such as in the case considered here, the falsepositive error rate can only be expected to be asymptotically equal to 1 − CL (Gallant 1975). Like the special LRtest constructed above, the Ftest relies on the assumption of a normal error model.
4.2. Robustness of the LR and Ftests
The upper plot in Fig. 1 shows spectral data points from the HGPS for the γray source HESS J1800−240. A fit of the data to a powerlaw model (Eq. (4)) results in the best fit parameters θ_{s} = (ϕ_{0}, Γ) = ((4.3 ± 0.1) × 10^{−13} TeV^{−1} cm^{−2} s^{−1}, 2.44 ± 0.02). The robustness of the LR and Ftests constructed above is tested in a Monte Carlo simulation of true null models where random spectra are created based on the bestfit model ϕ(Eθ_{s}). For each simulated spectrum, all five spectral points are sampled from a Student tdistribution t_{5} (location, scale) with location parameter ϕ(E_{i}θ_{s}), scale parameter σ_{i}, and 5 degrees of freedom. For comparison, the same simulation is repeated by sampling from a normal distribution with mean ϕ(E_{i}θ_{s}) and standard deviation σ_{i}. The Student error model predicts more flux outliers than expected in the normal error model on which the constructed LR and Ftests rely.
Fig. 1. Spectral fits to the data from the HGPS source for which the most (upper panel) and least (lower panel) constraining lower limit on E_{cut, h} is derived. Shown in green are the best fits to the data with a powerlaw spectrum. Red lines are best fits to the data with a γray powerlaw spectrum with exponential cutoff at the lower limit E_{cut, γ} derived from the HGPS data. Yellow lines are the corresponding best fits to the data in a hadronic model with cutoff at the lower limit E_{cut, h}. Limits for the lower plot are at 90% CL. In the upper panel, best fits with lower cutoff limit at 90% CL are hardly distinguishable from the best powerlaw fit. Therefore, 99% CL best fits are shown in this case. 

Open with DEXTER 
Figure 2 shows the distribution of LR and Ftest pvalues for Monte Carlo simulations with the Student (upper panel) and normal (lower panel) error models. The shown pvalues refer to tests of a powerlaw model against a model with exponential cutoff in the γray spectrum (Eq. (3) with β = 1). Similar simulations with unchanged conclusions were performed when a powerlaw model is tested against a model where the γray emission is modeled in a hadronic scenario (Eq. (1)). The lower panel of Fig. 2 shows that the distribution of pvalues for a true null hypothesis is compatible with the uniform expectation when the normal error model is used. However, severe discrepancies between the uniform expectation and the LRtest are observed in the upper panel of Fig. 2. In contrast, the Ftest shows no indication for a deviation from the uniform expectation when the normal error model is not true. The bias towards low pvalues of the LRtest in the Student error model is equivalent to a bias towards larger significances. Similar Monte Carlo simulations for all other γray sources that were selected for the analysis (see Sect. 3) confirm the observation that the LRtest is not robust against violations of the normal error model in the HGPS catalog.
Fig. 2. Distribution of pvalues obtained in LR and Ftests in a Monte Carlo simulation of true null hypotheses. The expectation for the pvalue distribution is uniform in the interval 0 ≤ p ≤ 0.5 where 50% of the pvalues are predicted. The other half of the pvalues are expected at p = 1, which is due to the constraint on positive energy cutoffs (see also Appendix A). Upper panel: resulting pvalue distribution in the range 0 ≤ p ≤ 0.5 when the simulated spectral points disperse around the true data according to a Student tdistribution with 5 degrees of freedom. It is clear that the LRtest (red) deviates from the uniform expectation for small pvalues. This means that the LRtest severely overestimates the significance in this case. The distribution of pvalues for the Ftest (blue) is compatible with the uniform expectation. The distribution of pvalues in the lower panel is generated when the simulated spectral points are normally distributed around their true value. In this case, the pvalue distribution for both tests is compatible with the uniform expectation. The parameters of the Monte Carlo simulation are discussed in Sect. 4.2. 

Open with DEXTER 
The improved robustness of the Ftest compared to the LRtest comes at the price of a reduced test power. The test power can be estimated with a toy Monte Carlo simulation of true alternative hypotheses. For this, five spectral points are calculated at energies E_{i} = ξ^{i}E_{min} where ξ = 2.4 and E_{min} = 0.4 TeV are the median values of the binning factor and the minimal energy bin for the selected HGPS data (see Table 1). At each energy, the flux is sampled from ϕ_{0}(E/E_{0})^{−Γ} exp(−E/E_{cut, γ}) where ϕ_{0} = 3.51 × 10^{−12} cm^{−2} s^{−1} TeV^{−1} and Γ = 2.11 are the median flux normalization and powerlaw index of the selected HGPS data (see Table 1). The flux normalization energy E_{0} = 1 TeV is fixed and the cutoff energy E_{cut, γ} is varied. The error of each sampled flux point is set to 20% of the flux in this toy simulation. Figure 3 compares the energy cutoff detection power between the LR and Ftests. It is clearly seen that, compared to the LRtest, the Ftest has less or at most equal power for all considered true energy cutoffs. For example, for a true energy cutoff of 100 TeV, the LRtest has 10% more power than the Ftest.
Fig. 3. Simulated coverage of the interval I given by Eq. (8) at 90% CL (blue). The coverage is compatible with the nominal confidence level for true energy cutoffs below 50 TeV where the Ftest power (magenta) is large. The interval I overcovers for large true energy cutoffs where the Ftest power is small. The LRtest power is shown in yellow for comparison. 

Open with DEXTER 
Despite the reduced power of the Ftest when compared to the LRtest, the possible influence of systematic effects in HGPS data motivates the use of the Ftest in the following discussion.
4.3. Lower limit on the energy cutoff
In practice, a hypothesis test for the presence of a cutoff is performed at high confidence level, for example at a significance level of 5σ. If the powerlaw hypothesis is not discarded, a lower limit on the energy cutoff E_{cut} is required at reduced confidence level (e.g., 90%) to constrain theoretical models. Depending on the alternative hypothesis, E_{cut} can refer to E_{cut, γ} (Eq. (3)) or to E_{cut, h} (Eq. (1)). A lower limit on the energy cutoff E_{cut} corresponds to a onesided confidence interval
is in the following used to refer to the maximum likelihood estimate of the energy cutoff as determined in a leastsquare fit. The lower limit is the solution to the equation
in the interval where
and RSS is given by Eq. (5). The parameters are the maximum likelihood estimators of ϕ_{0} and Γ for a powerlaw model with fixed energy cutoff at E. Similarly, is the maximum likelihood estimator for the parameters of a powerlaw model with variable exponential energy cutoff, again with the constraint that E_{cut} > 0. The critical value F_{crit} = F_{crit}(CL) is given by Eq. (7).
Equation (10) is the result of an Ftest of the hypothesis H_{0} : E_{cut} = E against the alternative hypothesis H_{1} : E_{cut} ≠ E and E_{cut} > 0. The existence of , defined in Eq. (8), is shown in Appendix B. In the following, it is assumed that is the unique solution of Eq. (10) in the energy interval (see also the discussion in Appendix B).
It is shown in Appendix C that the interval I is expected to overcover, that is, the probability to contain the true value is larger than the desired confidence level, when the true energy cutoff is large. However, excellent frequentist coverage properties are expected when the true energy cutoff is small such that the Ftest power at the given CL is large. Figure 3 shows the result of a Monte Carlo simulation where γray energy spectra with energy cutoff shape β = 1 were simulated following the respective discussion in Sect. 4.2. The result of this simulation shows that the frequentist coverage of the interval I agrees with the nominal 90% confidence level when the true energy cutoff is smaller than 50 TeV. For energy cutoffs larger than 50 TeV, the test power at the confidence level of 90% decreases fast and the interval overcovers.
It is concluded that the lower limit defined in Eq. (8) has very good coverage properties when the true cutoff of the γray spectrum is small. Otherwise it is conservative. For typical source and instrumental parameters representative of the HESS Galactic plane survey, good coverage can be expected when the true cutoff of the γray energy spectrum is below 50 TeV.
4.4. Systematic errors
The lower limit on the energy cutoff is the lower end of the confidence interval I given by Eq. (8). This is a confidence interval for the maximum likelihood estimator of the energy cutoff. Systematic effects on the energy cutoff become dominant when the systematic variation of is not covered by the confidence interval I. For data from Cherenkov telescope arrays like HESS and VERITAS, possible reasons for systematic variations of can be atmospheric and analysis effects or instrumental problems like broken camera elements. In case of the HGPS, the systematic errors on the fit parameters of a powerlaw model are estimated to be 30% on the flux normalization and 0.2 on the powerlaw index (Abdalla et al. 2018). These results were obtained from the comparison of the fit results obtained with different analysis methods. However, no public information is available on the systematic variation of in the HGPS.
5. Analysis of HGPS data
The analysis of the HGPS data, after the quality selection discussed in Sect. 3.2, starts with a consistency check. It is ensured that the HGPS spectral points can be fit with results that are in reasonable agreement with the HGPS analysis. Afterwards, lower limits on the spectral cutoffs are derived for the selected γray sources.
The HGPS catalog contains bestfitting powerlaw parameters for each source. The fit in the HGPS analysis is performed with an unbinned forward folding method (Piron et al. 2001). In contrast, binned leastsquare fits of the public spectral HGPS data are used in this analysis because unbinned data are not publicly available. However, the best leastsquare fit parameters are in reasonable agreement with the cataloged bestfit parameters. The maximum absolute difference between the powerlaw index obtained in this analysis and the powerlaw index given in the HGPS catalog is 1.8σ for HESS J1708−443 where σ denotes the HGPS catalog index error. The median absolute difference between the fitted powerlaw indices is 0.6σ. The corresponding largest absolute difference between the flux normalizations is 1.2σ for HESS J1026−582. The median absolute flux difference for the dataset is 0.1σ.
Table 2 shows the results of a search for spectral cutoffs for the 25 selected γray sources in the HGPS. The bestfit exponential cutoff models (Eq. (3) with β = 1 and Eq. (1)) have a pvalue larger than 1% in a χ^{2} test for all considered spectra. All calculated solutions to Eq. (10) are unique. No energy cutoff is detected with a statistical significance exceeding 5σ for any of the analyzed sources. The lower limit on the energy cutoff in a hadronic scenario is larger than 100 TeV for five sources. For these five sources, the best fit powerlaw index is in the range between 2.0 and 2.5. Figure 1 shows the spectral data together with the bestfit powerlaw and the excluded cutoffs for the sources HESS J1800−240 and HESS J1746−285. For these sources, the most and least constraining lower limit on the energy cutoff is derived in this analysis.
Results of the energy cutoff analysis for 25 unidentified γray sources in the HGPS catalog.
No systematic error on the lower limit on the energy cutoff can be derived based on HGPS data alone since the systematic variation of is not public. However, in the case of the source HESS J1908+063, an at least partially independent dataset is available. This dataset is discussed in the following section.
6. Analysis of MGRO J1908+06 data
The γray source MGRO J1908+06 is associated with HESS J1908+063. This γray source is one of the unidentified HGPS sources for which the lower limit on the energy cutoff in a hadronic scenario is found to be larger than 100 TeV in Sect. 5. Compared to the HGPS catalog, additional public data for this γray source are available. These data, acquired with different experiments, are described and analyzed in this section.
6.1. Data for MGRO J1908+06
Milagro detected the γray source at a median energy of 20 TeV with a flux of (8.8 ± 2.4) × 10^{−15} TeV^{−1} cm^{−2} s^{−1} (Abdo et al. 2007). HESS spectral data are published in Aharonian et al. (2009b), based on the analysis of 27 h of data. This dataset will be denoted by HESSd1 in the following while HESSd2 will be used to identify the spectral HGPS data on MGRO J1908+06. The 27 h HESSd1 dataset must be at least partially independent from the 12 h HESSd2 because of the greater livetime. VERITAS also observed MGRO J1908+063 and collected 62 h of data between the years 2007 and 2012. Spectral datapoints derived from the VERITAS dataset are published in Aliu et al. (2014). The spectral properties of the HESSd1, VERITAS, and Milagro datasets are compatible with each other (Aliu et al. 2014). Two lower limits on the energy cutoff E_{cut, γ} of the γray spectrum were derived before. A 90% CL lower limit on E_{cut, γ} is derived from a combination of the VERITAS and Milagro data at 17.7 TeV (Aliu et al. 2014). Also at 90% CL, the combination of the HESSd1 and Milagro data leads to a lower limit of E_{cut, γ} = 19.1 TeV (Aharonian et al. 2009b).
The additional Cherenkov telescope data for MGRO J1908+063 from HESS and VERITAS must pass the same quality selection as discussed for the HGPS spectral data in Sect. 3.2. The HESSd1 dataset passes this selection. However, the public spectral data for MGRO J1908+063 from VERITAS is binned with a factor ξ = 1.4 while the energy resolution is similar to that of HESS (Holder 2011). Only every second spectral point from the public VERITAS data is therefore used in this analysis. The following analysis results are checked with the full VERITAS dataset to ensure that no bias towards better results is introduced with this spectral point selection.
6.2. Data analysis
The analysis of the selected spectral data from HESSd1, VERITAS, and Milagro results in a bestfit powerlaw spectrum with ϕ_{0} = (4.2 ± 0.2) × 10^{−12} TeV^{−1} cm^{−2} s^{−1} and Γ = 2.17 ± 0.04. The powerlaw index is compatible with the index Γ = 2.26 ± 0.06 cataloged in the HGPS while the flux normalization in the HGPS, (1.05 ± 0.09) × 10^{−11} TeV^{−1} cm^{−2} s^{−1}, is incompatible. The fit to a powerlaw model has a pvalue of 0.91 in a χ^{2} goodnessoffit test. The lower 90% CL limits on the spectral energy cutoffs derived from this dataset are E_{cut, γ} = 29.5 TeV and E_{cut, h} = 141.5 TeV. These limits confirm the limits derived from the HGPS data. Overall, the agreement between this analysis and the HGPS analysis result is very good with regards to the derived lower limits on spectral energy cutoffs. Data and bestfitting models are shown in Fig. 4.
Fig. 4. Spectral γray data for the source MGRO J1908+06. Shown as green circles are the HESS data from Aharonian et al. (2009b). Shown in cyan open circles are the spectral data from the HGPS. Blue spectral points are from VERITAS (Aliu et al. 2014). VERITAS data marked with a square are not selected for the analysis to avoid correlations between spectral points. The Milagro data point is from Abdo et al. (2007). Shown in black, yellow, and red are best fits to all datapoints indicated by filled circles. The black line is a fit to a powerlaw model. Yellow and red lines are best fits to a powerlaw model with cutoff at the 90% CL lower limits E_{cut, γ} = 29.5 TeV (yellow) and E_{cut, h} = 141.5 TeV (red). 

Open with DEXTER 
Compared to the HESSd2 dataset with 6 spectral points, the combined data from HESSd1, VERITAS, and Milagro with 15 spectral points considered in this section contain many more flux measurements. This enables the search for a spectral break in the TeV energy range. Multiple Ftests are performed where the fit to the data of a powerlaw model (H_{0} : ΔΓ = 0) and a broken powerlaw (H_{1} : ΔΓ > 0, see Eq. (2)) are compared. The comparison cannot be performed in a single test because the break energy E_{break} is not defined under the null hypothesis. Using E_{i} again to denote the energies of the N available spectral flux data points, the energy range E_{2} = 0.75 TeV to E_{N − 2} = 14.8 TeV is scanned for a break energy E_{break} with a logarithmic binning factor of 1.05. The Ftest for the comparison of the fit quality of the powerlaw model and the broken powerlaw model, assuming the respective energy break energy, is performed in each scanning step. The Ftest is inverted to derive an upper limit on the index change ΔΓ as a function of the assumed break energy E_{break}. The result is shown in Fig. 5. An index change ΔΓ > 0.5 is ruled out at 90% CL for energies between 1 TeV and 10 TeV.
Fig. 5. Upper limit on a break in the powerlaw spectrum by ΔΓ as a function of the break energy. Shown as red line is the 90% CL upper limit on ΔΓ derived from the combined data discussed in Sect. 6. 

Open with DEXTER 
The largest local significance obtained in the 61 tests of the energy break scan is p_{local} = 0.14 with an Ftest statistic c = 1.23. The local significance must be transformed into a global significance p_{global} considering the number of performed tests. Based on results for the global pvalue in a multiple testing scenario (see Gross & Vitells 2010, and references therein), Algeri & van Dyk (2017) derive a global pvalue correction for an Fprocess. In the case considered here, this correction reads
where E[N_{c0}] is the expected number of upcrossings over the level c_{0} of the process of Ftest statistics during the energy cutoff scan when the null hypothesis is true. Equation (11) allows us to extrapolate the global pvalue correction from a low test statistic level c_{0}, where a small number of simulations is required for the estimation of E[N_{c0}], to a larger test statistic value c. Following Algeri et al. (2016b), E[N_{c0}] is estimated in a parametric bootstrap simulation at c_{0} = 0.3 to be E[N_{c0}] = 0.52 ± 0.02. With this correction, the global pvalue of the Ftest for a spectral break in the selected MGRO J1908+06 data from HESS, VERITAS, and Milagro can be estimated to be p_{global} = 0.5. The powerlaw hypothesis is therefore not discarded in favor of a broken powerlaw.
7. Discussion
The analysis of HGPS data resulted in five γray sources for which the energy spectrum is compatible with a powerlaw, and in a hadronic scenario the lower limit on the energy cutoff of the accelerated particle population is larger than 100 TeV. These five γray sources are discussed in detail in this section. An emphasis of the discussion is on the plausibility of the hadronic scenario for these sources.
7.1. HESS J1800−240
The region around HESS J1800−240 is found to have a complex γray morphology in Aharonian et al. (2008) based on a HESS dataset of 42 h livetime. Four γray sources are detected in this region. The SNR W28 is associated with the source HESS J1801−233. The region around HESS J1800−240 itself is subdivided into three hotspots that are spatially coincident with molecular clouds.
Based on a smaller dataset of 10 h, the spectral HGPS analysis detects only one γray source in the region. The lower limit on the energy cutoff for accelerated particles derived above is E_{cut, h} = 510 TeV at 90% CL. This constraint is comparable to the lower limit of E_{cut, h} = 600 TeV at 90% CL on the energy cutoff of particles in the vicinity of the Galactic center in a hadronic scenario (Abramowski et al. 2016).
A physics case for the association of this γray source with a Pevatron is the model discussed in Gabici & Aharonian (2007) where delayed CRs from a SNR illuminate a molecular cloud. The plausible physics case and the interesting constraint on the energy cutoff in a hadronic scenario may motivate a systematic investigation of differences between the HGPS analysis and the analysis in Aharonian et al. (2008) as well as a deeper exposure of the region with current Cherenkov telescopes.
7.2. HESS J1641−463
The γray source is discussed in Abramowski et al. (2014), where the TeV spectrum extracted from 72 h of HESS data is found to be compatible with a powerlaw model. Additionally, a lower limit on the energy cutoff at E_{cut, h} = 100 TeV is derived at 99% CL in a hadronic model. The lower limit of E_{cut, h} = 185 TeV at 90% CL on the energy cutoff derived in this work from 27 h of HGPS data is compatible with the result in Abramowski et al. (2014).
A detailed hadronic model for the γray emission of HESS J1651−463 is presented in Tang et al. (2015). Within this model, runaway CR particles accelerated inside the young and nearby SNR G338.3−0.0 interact with a molecular cloud where γrays are produced via the decay of neutral pions. The SNR G338.3−0.0 itself is associated with the γray source HESS J1640−465.
Given the presence of a nearby SNR and a molecular cloud, a hadronic scenario for this γray source is plausible. Observations with instruments which have an improved sensitivity, at a few tens of TeV to a few hundred TeV, or a neutrino detection are necessary to confirm a Pevatron scenario.
7.3. MGRO J1908+06
The γray source MGRO J1908+06 was discovered by Milagro (Abdo et al. 2007) and later confirmed by HESS (Aharonian et al. 2009b). Abdo et al. (2009) find a spatial coincidence between MGRO J1908+06 and the PWN of the pulsar PSR 1907.5+0602 detected by Fermi/LAT. However, three significant emission regions are detected in the field with data from VERITAS (Aliu et al. 2014). This raises the question as to whether the entire emission can originate from the PWN or an additional γray source is present. The SNR G40.5−0.5 is discussed in Aliu et al. (2014) as a candidate for an additional γray source.
The HESS spectrum is confirmed by the measurement with VERITAS (Aliu et al. 2014). However, the spectrum measured by the Astrophysical Radiation with Groundbased Observatory (ARGOYBJ) is incompatible with the spectra measured by HESS and VERITAS (Bartoli et al. 2012). The HighAltitude Water Cherenkov GammaRay Observatory (HAWC) detects the γray source at an energy of 7 TeV and derives a flux which varies within a factor of 2.5 depending on the assumed source extension (Abeysekara et al. 2017). The γray source is discussed as a possible source of highenergy neutrinos (Halzen et al. 2017). Aartsen et al. (2017b) find a pretrial pvalue of 2.5% when testing a backgroundonly hypothesis against the emission of astrophysical highenergy neutrinos from MGRO J1908+06 based on data from IceCube.
It is concluded that the spectral properties of the γray emission of MGRO J1908+06 are currently not unambiguously measured. Moreover, the source identification is not finally resolved, although an at least partial association with the PWN of PSR 1907.05+0602 is likely.
The lower limit on the exponential cutoff of the γray spectrum derived above, E_{cut, γ} ≈ 30 TeV, is a factor of almost two larger than previous results discussed in Aharonian et al. (2009b) and Aliu et al. (2014). The lower limit E_{cut, γ} ≈ 30 TeV is derived from HGPS data (E_{cut, γ} = 31.1 TeV at 90% CL) and confirmed in an analysis of combined data from HESS, VERITAS, and Milagro (E_{cut, γ} = 29.5 TeV at 90% CL). The compatibility of the lower limits on the energy cutoff derived from two different datasets shows that the systematic error on the energy cutoff is not likely to be dominant in the case of this source.
A sharp energy break is ruled out at 90% CL in the energy range between 1 TeV and 10 TeV (see Fig 5). A search for such an energy break could, for example, be motivated by the sharp energy break observed around 1 TeV in the electron spectrum measured on earth (see e.g., Aharonian et al. 2009a. The powerlaw index of the electron spectrum measured on earth changes within measured errors by ΔΓ_{e} = 1 at E_{break} = 1 TeV. This energy break can be interpreted as a cooling break (Recchia et al. 2019). When electrons generate the γray signal detected towards MGRO J1908+06 via inverse Compton scattering in the Thomson regime, a break in the spectrum of electrons by ΔΓ_{e} = 1 would translate into a change of the γray spectrum index by ΔΓ = 0.5 (see e.g., Hinton & Hofmann 2009). An energy break of this kind is ruled out at 90% CL in the energy range between 1 TeV and 10 TeV.
Although the derived lower limits on the γray energy cutoff are more constraining than previous lower limits and a cooling break can be ruled out in the energy range between 1 and 10 TeV, the results cannot rule out leptonic scenarios for the γray emission. Deeper exposures and observations with more sensitive instruments are necessary to reveal the nature of this source.
7.4. HESS J1852−000
It is discussed in Abdalla et al. (2018) that the spectral HGPS data for this source must be treated with caution due to deviations between the main HGPS result and an independent cross check data analysis. HESS data for this source have been analyzed before (Kosack et al. 2011), but an energy spectrum has not been derived. No other spectral data in the TeV γray energy range are publicly available. The analysis of the spectral cutoff for this source is therefore not conclusive.
Different scenarios for the production of γrays are discussed in Kosack et al. (2011) and Bamba et al. (2016). Among them is an association with the SNR Kes 78 and a nearby molecular cloud.
It is concluded that a hadronic source scenario in which runaway CRs from the SNR Kes 78 illuminate a molecular cloud can currently not be ruled out. The inconclusive spectral data may motivate a more detailed analysis and possibly further observations with currentgeneration Cherenkov telescopes.
7.5. HESS J1634−472
HESS J1634−472 is cataloged as an unidentified source in the HGPS without further discussion (Abdalla et al. 2018). The source was also detected in a previous survey of the inner Galaxy with HESS (Aharonian et al. 2006a). Here, the SNR G337.2+0.1 and a source of Xrays are found to be within 0.2° of HESS J1634−472 but no association is claimed. Acero et al. (2013) discuss a detection of the source with Fermi/LAT and find that there is no counterpart pulsar in this region which is energetic enough to power the γray source HESS J1634−472.
No further public TeV γray data are available for this source. The interesting constraint on E_{cut, h} > 108 TeV at 90% CL derived from the analysis of HGPS data above may motivate a dedicated reanalysis of the HESS data and possibly a deeper exposure of the source.
8. Conclusion
Five γray sources are found in the HGPS catalog for which the maximal energy of accelerated particles is at least 100 TeV in a hadronic model. For at least three of these sources, a hadronic scenario for the γray emission is plausible, that is, as a result of the interaction of runaway CRs with a molecular cloud. One of the Pevatron candidates found in the HGPS is MGRO J1908+06. The γray spectrum of this source extends to at least E_{cut, γ} = 30 TeV without indication of a cutoff. The lower limit on E_{cut, γ} for this source is a factor of almost two larger than previous results. A break ΔΓ > 0.5 for this γray source can be ruled out at 90% CL in the energy range between 1 and 10 TeV.
The search presented in this work can only find Pevatron candidates. A conclusive identification must be performed with data from neutrino telescopes (Aartsen et al. 2017a) or γray data around and above 100 TeV. The extended air shower array HAWC (DeYoung 2012), now operating for more than 3 years, can add spectral γray data at a few tens of TeV. An improved sensitivity and higher energies will be accessible with the upcoming γray detectors, for example the Cherenkov Telescope Array (CTA; The CTA Consortium 2011) and the Large High Altitude Air Shower Observatory (LHAASO; Sciascio 2016). Finally, the development of new experimental techniques, such as for the Tunka Advanced Instrument for Gammaray and Cosmic ray Astrophysics (TAIGA; Tluczykont et al. 2014), will help to identify Galactic Pevatrons, thereby aiding the quest for the origin of cosmic rays.
Acknowledgments
The author acknowledges financial support by the German Ministry for Education and Research (BMBF).
References
 Aartsen, M. G., Ackermann, M., Adams, J., et al. 2017a, J. Instrum., 12, P03012 [CrossRef] [Google Scholar]
 Aartsen, M. G., Abraham, K., Ackermann, M., et al. 2017b, ApJ, 835, 151 [NASA ADS] [CrossRef] [Google Scholar]
 Abdalla, H., Abramowski, A., Aharonian, F., et al. 2018, A&A, 612, A8 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Abdo, A. A., Allen, B., Berley, D., et al. 2007, ApJ, 664, L91 [NASA ADS] [CrossRef] [Google Scholar]
 Abdo, A. A., Allen, B. T., Aune, T., et al. 2009, ApJ, 700, L127 [NASA ADS] [CrossRef] [Google Scholar]
 Abeysekara, A. U., Albert, A., Alfaro, R., et al. 2017, ApJ, 843, 39 [NASA ADS] [CrossRef] [Google Scholar]
 Abramowski, A., Aharonian, F., Ait Benkhali, F., et al. 2014, ApJ, 794, L16 [NASA ADS] [CrossRef] [Google Scholar]
 Abramowski, A., Aharonian, F., Ait Benkhali, F., et al. 2016, Nature, 531, 476 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Acero, F., Ackermann, M., Ajello, M., et al. 2013, ApJ, 773, 77 [NASA ADS] [CrossRef] [Google Scholar]
 Ackermann, M., Ajello, M., Allafort, A., et al. 2013, Science, 339, 807 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Aharonian, F., Akhperjanian, A. G., BazerBachi, A. R., et al. 2006a, ApJ, 636, 777 [NASA ADS] [CrossRef] [Google Scholar]
 Aharonian, F., Akhperjanian, A. G., BazerBachi, A. R., et al. 2006b, A&A, 457, 899 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Aharonian, F., Akhperjanian, A. G., BazerBachi, A. R., et al. 2008, A&A, 481, 401 [NASA ADS] [CrossRef] [EDP Sciences] [MathSciNet] [PubMed] [Google Scholar]
 Aharonian, F., Akhperjanian, A., Anton, G., et al. 2009a, A&A, 508, 561 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Aharonian, F., Akhperjanian, A. G., Anton, G., et al. 2009b, A&A, 499, 723 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Ahnen, M. L., Ansoldi, S., Antonelli, L. A., et al. 2017, MNRAS, 472, 2956 [NASA ADS] [CrossRef] [Google Scholar]
 Algeri, S., & van Dyk, D. A. 2017, ArXiv eprints [arXiv:1701.06820v4] [Google Scholar]
 Algeri, S., Conrad, J., & van Dyk, D. A. 2016a, MNRAS, 458, L84 [NASA ADS] [CrossRef] [Google Scholar]
 Algeri, S., van Dyk, D., Conrad, J., & Anderson, B. 2016b, J. Instrum., 11, P12010 [CrossRef] [Google Scholar]
 Aliu, E., Archambault, S., Aune, T., et al. 2014, ApJ, 787, 166 [NASA ADS] [CrossRef] [Google Scholar]
 Amenomori, M., Bao, Y. W., Bi, X. J., et al. 2019, Phys. Rev. Lett., 123, 051101 [NASA ADS] [CrossRef] [Google Scholar]
 Baade, W., & Zwicky, F. 1934, Proc. Natl. Acad. Sci., 20, 259 [NASA ADS] [CrossRef] [Google Scholar]
 Bamba, A., Terada, Y., Hewitt, J., et al. 2016, ApJ, 818, 63 [NASA ADS] [CrossRef] [Google Scholar]
 Bartoli, B., Bernardini, P., Bi, X. J., et al. 2012, ApJ, 760, 110 [NASA ADS] [CrossRef] [Google Scholar]
 Chernoff, H. 1954, Ann. Math. Statist., 25, 573 [CrossRef] [Google Scholar]
 DeYoung, T. 2012, Nucl. Instrum. Methods Phys. Res. Sect. A, 692, 72 [NASA ADS] [CrossRef] [Google Scholar]
 Gabici, S., & Aharonian, F. 2018, A&A, 612, A8 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gabici, S., & Aharonian, F. A. 2007, ApJ, 665, L131 [NASA ADS] [CrossRef] [Google Scholar]
 Gaensler, B. M., & Slane, P. O. 2006, ARA&A, 44, 17 [NASA ADS] [CrossRef] [Google Scholar]
 Gallant, A. R. 1975, Am. Stat., 29, 73 [Google Scholar]
 Gross, E., & Vitells, O. 2010, Eur. Phys. J. C, 70, 525 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Halzen, F., Kheirandish, A., & Niro, V. 2017, Astropart. Phys., 86, 46 [NASA ADS] [CrossRef] [Google Scholar]
 Hillas, A. M. 2005, J. Phys. G: Nucl. Part. Phys., 31, R95 [NASA ADS] [CrossRef] [Google Scholar]
 Hinton, J. A., & Hofmann, W. 2009, ARA&A, 47, 523 [NASA ADS] [CrossRef] [Google Scholar]
 Holder, J. 2011, Int. Cosmic Ray Conf., 12 [Google Scholar]
 Jogler, T., & Funk, S. 2016, ApJ, 816, 100 [NASA ADS] [CrossRef] [Google Scholar]
 Kelner, S. R., Aharonian, F. A., & Bugayov, V. V. 2006, Phys. Rev. D, 74, 039901 [NASA ADS] [CrossRef] [Google Scholar]
 Kosack, K., Chaves, R. C. G., & Acero, F. 2011, Int. Cosmic Ray Conf., 7, 195 [Google Scholar]
 Piron, F., DjannatiAtai, A., Punch, M., et al. 2001, A&A, 374, 895 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Recchia, S., Gabici, S., Aharonian, F. A., & Vink, J. 2019, Phys. Rev. D, 99, 103022 [NASA ADS] [CrossRef] [Google Scholar]
 Sciascio, G. D. 2016, in Proceedings of the 9th Cosmic Ray International Seminar, Nuclear and Particle Physics Proceedings, 279, 166 [Google Scholar]
 Tang, Y., Yang, C., Zhang, L., & Wang, J. 2015, ApJ, 812, 32 [NASA ADS] [CrossRef] [Google Scholar]
 The CTA Consortium (Actis, M., et al.) 2011, Exp. Astron., 32, 193 [NASA ADS] [CrossRef] [Google Scholar]
 Tluczykont, M., Hampf, D., Horns, D., et al. 2014, Astropart. Phys., 56, 42 [NASA ADS] [CrossRef] [Google Scholar]
 Wilks, S. S. 1935, Ann. Math. Stat., 9, 60 [CrossRef] [Google Scholar]
 Zhang, X., & Liu, S. 2019, ApJ, 874, 24 [NASA ADS] [CrossRef] [Google Scholar]
Appendix A: Critical values
Consider, as in Sect. 4.1, the two nested hypotheses H_{0} : λ = 0 and H_{1} : λ > 0 for the single parameter λ. The parameter set for H_{0} is not in the interior of the parameter set of H_{1}. The asymptotic result for the LR test statistic −2lnΛ under H_{0} discussed in Wilks (1935) is therefore not applicable. However, a modification is discussed in Chernoff (1954) where it is predicted that −2lnΛ is asymptotically distributed like an equal mixture of two random variables when H_{0} is true. In the considered case, one of the random variables is χ^{2}distributed with one degree of freedom. The other random variable has a Dirac delta function probability density function δ(0) (see also Algeri et al. 2016a). With this asymptotic distribution f_{H0} of the test statistic −2lnΛ under H_{0}, the critical value Λ_{crit} of the LR test statistic can be inferred from . Let denote the cumulative density functions of a χ^{2}distribution with one degree of freedom. Subsequently, . This leads to the equality which can be used to calculate the critical value Λ_{crit} of the LR test statistic at a given CL.
A similar result follows for the critical value of the asymptotic F statistic in Sect. 4.1.2: the F statistic (Eq. (6)) and the LR statistic −2lnΛ are related via (see Sect. 4.1.1). Let now F_{1, N − 3} be the cumulative density function for a Fdistributed random variable with 1 and N − 3 degrees of freedom. Because the LR statistic is asymptotically expected to be zero in half of the tests when the null hypothesis is true, the critical value F_{crit} of the Ftest can again be inferred from to be .
Appendix B: Existence of the lower limit
Let F(E) and F_{crit} be given by Eqs. (7) and (10) and. The lower limit on the energy cutoff is defined as a solution to the equation
in the interval . Here, is the maximum likelihood estimator for the energy cutoff (see Sect. 4.3). It holds that F_{crit} > 0 and F(E) is a continuous function. A solution to Eq. (B.1) does always exist in the given interval because and F(E) → ∞ for E → 0 by definition of F(E).
However, the solution to Eq. (B.1) in the interval must not necessarily be unique when general data is considered. Multiple solutions can, for example, occur in case of a mismatch between the model with an exponential cutoff and the true model where the likelihood function has multiple local maxima. In this case, the meaning of a lower limit on the energy cutoff is questionable.
In any case, special care must be taken when multiple solutions to Eq. (B.1) are found in the interval . In the data analyzed in this work, all lower limits were unique when using a confidence level of 90%.
Appendix C: Frequentist coverage of the lower limit
Consider the interval I given by Eq. (8). It is argued that I has excellent frequentist coverage properties when the true energy cutoff is small such that the test power is large. When the true energy cutoff is large such that the test power is small, the interval I overcovers and the lower limit on the energy cutoff is conservative.
In the following, F, F(E), and F_{crit} are given by Eqs. (6), (7), and (10). It is assumed that only one solution for Eq. (B.1) in the interval exists. Similarly, it is also assumed that at most one solution to Eq. (B.1) exists in the interval . Moreover, is the solution to Eq. (B.1) in , if it exists, and if no solution exists. Now consider the interval . This interval is the acceptance interval A = {EF(E) ≤ F_{crit}} of an Ftest, and as such defines a confidence interval for . However, the confidence level of the Ftest to which A is the acceptance interval depends on E. This is a result of the fit being constrained to λ > 0.
When the true energy cutoff is large, F(E) will need to be evaluated at large energies E to obtain A. By definition of F(E), it holds that F(E)→F when E → ∞ and F(E) is the Ftest statistic for a test of a powerlaw hypothesis against a powerlaw model with exponential cutoff. This means that A is, in the limit E → ∞, the acceptance interval for an Ftest constructed at confidence level CL. In other words, it holds that P(E_{true} ∈ A) = CL. Additionally, by definition of the interval I, it holds that P(E_{true} ∈ I) ≥ P(E_{true} ∈ A). Together it follows that P(E_{true} ∈ I) ≥ CL, which is the overcoverage of I as a frequentist confidence interval when the true energy cutoff is large.
Now, let the true energy cutoff be small but positive such that the powerlaw hypothesis is discarded at the given CL. The fit constraint λ > 0 is in this case irrelevant because the best fit λ is always large. By definition of F_{crit}, it holds that P(E_{true} ∈ A) = 2CL − 1. This is because without the constraint λ > 0, A is the acceptance interval of an Ftest that is constructed at confidence level 2CL − 1. Close to , F(E) can be quadratically approximated; see also Fig. C.2. Using the nomenclature of Sect. 4.1, it holds that where σ^{2} is the inverse Fisher information of L(θ_{1}) evaluated at . This means that the confidence interval A will be symmetric around when E is small such that the solution to Eq. (B.1) is found in the range where the quadratic approximation holds. It follows that . For the coverage probability, it follows that . This states that very good coverage properties of I as a frequentist confidence interval are expected when the Ftest power is large at the given CL.
Fig. C.1. Example for the dependence of F(E) given by Eq. (10) on the cutoff energy E. The example is based on the HGPS data for the source HESS J1858+020. The maximum likelihood estimate for an exponential cutoff in the γray spectrum (Eq. (3) with β = 1) is indicated with a green dashed line. The blue and the magenta lines indicate the critical value F_{crit} and the inferred lower limit on E_{cut, γ} at 90% CL. The magenta region is the excluded range of cutoff energies. 

Open with DEXTER 
Fig. C.2. Same as in Fig. C.1 but based on the HGPS data for the source HESS J1507−622. The powerlaw hypotheses is discarded at 90% CL and the acceptance interval A = {EF(E) ≤ F_{crit}} is now a finite interval roughly symmetric around . Shown in yellow is the quadratic approximation of F(E) around . 

Open with DEXTER 
All Tables
Results of the energy cutoff analysis for 25 unidentified γray sources in the HGPS catalog.
All Figures
Fig. 1. Spectral fits to the data from the HGPS source for which the most (upper panel) and least (lower panel) constraining lower limit on E_{cut, h} is derived. Shown in green are the best fits to the data with a powerlaw spectrum. Red lines are best fits to the data with a γray powerlaw spectrum with exponential cutoff at the lower limit E_{cut, γ} derived from the HGPS data. Yellow lines are the corresponding best fits to the data in a hadronic model with cutoff at the lower limit E_{cut, h}. Limits for the lower plot are at 90% CL. In the upper panel, best fits with lower cutoff limit at 90% CL are hardly distinguishable from the best powerlaw fit. Therefore, 99% CL best fits are shown in this case. 

Open with DEXTER  
In the text 
Fig. 2. Distribution of pvalues obtained in LR and Ftests in a Monte Carlo simulation of true null hypotheses. The expectation for the pvalue distribution is uniform in the interval 0 ≤ p ≤ 0.5 where 50% of the pvalues are predicted. The other half of the pvalues are expected at p = 1, which is due to the constraint on positive energy cutoffs (see also Appendix A). Upper panel: resulting pvalue distribution in the range 0 ≤ p ≤ 0.5 when the simulated spectral points disperse around the true data according to a Student tdistribution with 5 degrees of freedom. It is clear that the LRtest (red) deviates from the uniform expectation for small pvalues. This means that the LRtest severely overestimates the significance in this case. The distribution of pvalues for the Ftest (blue) is compatible with the uniform expectation. The distribution of pvalues in the lower panel is generated when the simulated spectral points are normally distributed around their true value. In this case, the pvalue distribution for both tests is compatible with the uniform expectation. The parameters of the Monte Carlo simulation are discussed in Sect. 4.2. 

Open with DEXTER  
In the text 
Fig. 3. Simulated coverage of the interval I given by Eq. (8) at 90% CL (blue). The coverage is compatible with the nominal confidence level for true energy cutoffs below 50 TeV where the Ftest power (magenta) is large. The interval I overcovers for large true energy cutoffs where the Ftest power is small. The LRtest power is shown in yellow for comparison. 

Open with DEXTER  
In the text 
Fig. 4. Spectral γray data for the source MGRO J1908+06. Shown as green circles are the HESS data from Aharonian et al. (2009b). Shown in cyan open circles are the spectral data from the HGPS. Blue spectral points are from VERITAS (Aliu et al. 2014). VERITAS data marked with a square are not selected for the analysis to avoid correlations between spectral points. The Milagro data point is from Abdo et al. (2007). Shown in black, yellow, and red are best fits to all datapoints indicated by filled circles. The black line is a fit to a powerlaw model. Yellow and red lines are best fits to a powerlaw model with cutoff at the 90% CL lower limits E_{cut, γ} = 29.5 TeV (yellow) and E_{cut, h} = 141.5 TeV (red). 

Open with DEXTER  
In the text 
Fig. 5. Upper limit on a break in the powerlaw spectrum by ΔΓ as a function of the break energy. Shown as red line is the 90% CL upper limit on ΔΓ derived from the combined data discussed in Sect. 6. 

Open with DEXTER  
In the text 
Fig. C.1. Example for the dependence of F(E) given by Eq. (10) on the cutoff energy E. The example is based on the HGPS data for the source HESS J1858+020. The maximum likelihood estimate for an exponential cutoff in the γray spectrum (Eq. (3) with β = 1) is indicated with a green dashed line. The blue and the magenta lines indicate the critical value F_{crit} and the inferred lower limit on E_{cut, γ} at 90% CL. The magenta region is the excluded range of cutoff energies. 

Open with DEXTER  
In the text 
Fig. C.2. Same as in Fig. C.1 but based on the HGPS data for the source HESS J1507−622. The powerlaw hypotheses is discarded at 90% CL and the acceptance interval A = {EF(E) ≤ F_{crit}} is now a finite interval roughly symmetric around . Shown in yellow is the quadratic approximation of F(E) around . 

Open with DEXTER  
In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.