Evidence that Ultra-High-Energy Gamma Rays are a Universal Feature Near Powerful Pulsars
HAWC Collaboration, A. Albert, R. Alfaro, C. Alvarez, J.D. ?lvarez, J.R. Angeles Camacho, J.C. Arteaga-Velázquez, K.P. Arunbabu, D. Avila Rojas, H.A. Ayala Solares, V. Baghmanyan, E. Belmont-Moreno, S.Y. BenZvi, C. Brisbois, K.S. Caballero-Mora, T. Capistrán, A. Carramiñana, S. Casanova, U. Cotti, J. Cotzomi, S. Coutiño de León, E. De la Fuente, C. de León, R. Diaz Hernandez, B.L. Dingus, M.A. DuVernois, M. Durocher, J.C. Díaz-Vélez, R.W. Ellsworth, K. Engel, C. Espinoza, K.L. Fan, M. Fernández Alonso, N. Fraija, A. Galván-Gámez, J.A. García-González, F. Garfias, G. Giacinti, M.M. González, J.A. Goodman, J.P. Harding, S. Hernandez, B. Hona, D. Huang, F. Hueyotl-Zahuantitla, P. Hüntemeyer, A. Iriarte, A. Jardin-Blicq, V. Joshi, D. Kieda, A. Lara, W.H. Lee, J. Lee, H. León Vargas, J.T. Linnemann, A.L. Longinotti, G. Luis-Raya, J. Lundeen, K. Malone, V. Marandon, O. Martinez, J. Martínez-Castro, J.A. Matthews, P. Miranda-Romagnoli, J.A. Morales-Soto, E. Moreno, M. Mostafá, A. Nayerhoda, L. Nellen, M. Newbold, M.U. Nisa, R. Noriega-Papaqui, L. Olivera-Nieto, N. Omodei, A. Peisker, Y. Pérez Araujo, E.G. Pérez-Pérez, C.D. Rho, Y.J. Roh, D. Rosa-González, E. Ruiz-Velasco, H. Salazar, F. Salesa Greus, A. Sandoval, M. Schneider, H. Schoorlemmer, J. Serna-Franco, A.J. Smith, R.W. Springer, P. Surajbali, M. Tanner, K. Tollefson, I. Torres, R. Torres-Escobedo, R. Turner, F. Ureña-Mena, L. Villaseñor, T. Weisgarber, E. Willox, H. Zhou
DDraft version January 21, 2021
Typeset using L A TEX twocolumn style in AASTeX62
Evidence that Ultra-High-Energy Gamma Rays are a Universal Feature Near Powerful Pulsars
A. Albert, R. Alfaro, C. Alvarez, J.D. ´Alvarez, J.R. Angeles Camacho, J.C. Arteaga-Vel´azquez, K.P. Arunbabu, D. Avila Rojas, H.A. Ayala Solares, V. Baghmanyan, E. Belmont-Moreno, S.Y. BenZvi, C. Brisbois, K.S. Caballero-Mora, T. Capistr´an, A. Carrami˜nana, S. Casanova, U. Cotti, J. Cotzomi, S. Couti˜no de Le´on, E. De la Fuente, C. de Le´on, R. Diaz Hernandez, B.L. Dingus, M.A. DuVernois, M. Durocher, J.C. D´ıaz-V´elez, R.W. Ellsworth, K. Engel, C. Espinoza, K.L. Fan, M. Fern´andez Alonso, N. Fraija, A. Galv´an-G´amez, J.A. Garc´ıa-Gonz´alez, F. Garfias, G. Giacinti, M.M. Gonz´alez, J.A. Goodman, J.P. Harding, S. Hernandez, B. Hona, D. Huang, F. Hueyotl-Zahuantitla, P. H¨untemeyer, A. Iriarte, A. Jardin-Blicq,
16, 19, 20
V. Joshi, D. Kieda, A. Lara, W.H. Lee, J. Lee, H. Le´on Vargas, J.T. Linnemann, A.L. Longinotti,
11, 10
G. Luis-Raya, J. Lundeen, K. Malone, V. Marandon, O. Martinez, J. Mart´ınez-Castro, J.A. Matthews, P. Miranda-Romagnoli, J.A. Morales-Soto, E. Moreno, M. Mostaf´a, A. Nayerhoda, L. Nellen, M. Newbold, M.U. Nisa, R. Noriega-Papaqui, L. Olivera-Nieto, N. Omodei, A. Peisker, Y. P´erez Araujo, E.G. P´erez-P´erez, C.D. Rho, Y.J. Roh, D. Rosa-Gonz´alez, E. Ruiz-Velasco, H. Salazar, F. Salesa Greus, A. Sandoval, M. Schneider, H. Schoorlemmer, J. Serna-Franco, A.J. Smith, R.W. Springer, P. Surajbali, M. Tanner, K. Tollefson, I. Torres, R. Torres-Escobedo, R. Turner, F. Ure˜na-Mena, L. Villase˜nor, T. Weisgarber, E. Willox, and H. Zhou HAWC Collaboration Physics Division, Los Alamos National Laboratory, Los Alamos, NM, USA Instituto de F´ısica, Universidad Nacional Aut´onoma de M´exico, Ciudad de Mexico, Mexico Universidad Aut´onoma de Chiapas, Tuxtla Guti´errez, Chiapas, M´exico Universidad Michoacana de San Nicol´as de Hidalgo, Morelia, Mexico Instituto de Geof´ısica, Universidad Nacional Aut´onoma de M´exico, Ciudad de Mexico, Mexico Department of Physics, Pennsylvania State University, University Park, PA, USA Institute of Nuclear Physics Polish Academy of Sciences, PL-31342 IFJ-PAN, Krakow, Poland Department of Physics & Astronomy, University of Rochester, Rochester, NY , USA Department of Physics, University of Maryland, College Park, MD, USA Instituto de Astronom´ıa, Universidad Nacional Aut´onoma de M´exico, Ciudad de Mexico, Mexico Instituto Nacional de Astrof´ısica, ´Optica y Electr´onica, Puebla, Mexico Facultad de Ciencias F´ısico Matem´aticas, Benem´erita Universidad Aut´onoma de Puebla, Puebla, Mexico Departamento de F´ısica, Centro Universitario de Ciencias Exactase Ingenierias, Universidad de Guadalajara, Guadalajara, Mexico Department of Physics, University of Wisconsin-Madison, Madison, WI, USA Tecnologico de Monterrey, Escuela de Ingenier´ıa y Ciencias, Ave. Eugenio Garza Sada 2501, Monterrey, N.L., Mexico, 64849 Max-Planck Institute for Nuclear Physics, 69117 Heidelberg, Germany Department of Physics and Astronomy, University of Utah, Salt Lake City, UT, USA Department of Physics, Michigan Technological University, Houghton, MI, USA Department of Physics, Faculty of Science, Chulalongkorn University, 254 Phayathai Road, Pathumwan, Bangkok 10330, Thailand National Astronomical Research Institute of Thailand (Public Organization), Don Kaeo, MaeRim, Chiang Mai 50180, Thailand Erlangen Centre for Astroparticle Physics, Friedrich-Alexander-Universit¨at Erlangen-N¨urnberg, Erlangen, Germany University of Seoul, Seoul, Rep. of Korea Department of Physics and Astronomy, Michigan State University, East Lansing, MI, USA Universidad Politecnica de Pachuca, Pachuca, Hgo, Mexico Centro de Investigaci´on en Computaci´on, Instituto Polit´ecnico Nacional, M´exico City, M´exico. Dept of Physics and Astronomy, University of New Mexico, Albuquerque, NM, USA Universidad Aut´onoma del Estado de Hidalgo, Pachuca, Mexico Instituto de Ciencias Nucleares, Universidad Nacional Aut´onoma de Mexico, Ciudad de Mexico, Mexico Department of Physics, Stanford University: Stanford, CA 94305–4060, USA Tsung-Dao Lee Institute and School of Physics and Astronomy, Shanghai Jiao Tong University, Shanghai, China
Corresponding author: Kelly [email protected] a r X i v : . [ a s t r o - ph . H E ] J a n HAWC Collaboration
Submitted to The Astrophysical Journal LettersABSTRACTThe highest-energy known gamma-ray sources are all located within 0.5 degrees of extremely powerfulpulsars. This raises the question of whether ultra-high-energy (UHE; >
56 TeV) gamma-ray emissionis a universal feature expected near pulsars with a high spin-down power. Using four years of data fromthe High Altitude Water Cherenkov (HAWC) Gamma-Ray Observatory, we present a joint-likelihoodanalysis of ten extremely powerful pulsars to search for UHE gamma-ray emission correlated withthese locations. We report a significant detection ( > σ ), indicating that UHE gamma-ray emission isa generic feature of powerful pulsars. We discuss the emission mechanisms of the gamma rays and theimplications of this result. The individual environment that each pulsar is found in appears to play arole in the amount of emission. INTRODUCTIONUltra-high-energy (UHE; >
56 TeV) gamma-ray emis-sion can be created via hadronic or leptonic processes.In the hadronic mechanism, a neutral pion decays intotwo gamma rays. In the leptonic mechanism, a lower-energy photon scatters off a relativistic electron via in-verse Compton scattering. The electron transfers part ofits energy to the gamma ray, resulting in a higher-energyphoton.Traditionally, it was thought that UHE gamma-raysources would be hadronic in nature, as leptonic emis-sion is suppressed in this energy regime due to the Klein-Nishina (KN) effect. However, present-day gamma-rayobservatories have the sensitivity required to detect lep-tonic sources above 56 TeV. See, for example, the de-tections of the Crab Nebula as well as several knownpulsar wind nebulae (PWN) (Abeysekara et al. 2019;Amenomori et al. 2019; Abeysekara et al. 2020; Abdallaet al. 2019). The astrophysical spectrum of a leptonicsource typically exhibits significant curvature due to theKN effect (Moderski et al. 2005). Conversely, hadronicsources follow the spectrum of their parent cosmic-raypopulation, which may or may not include a cutoff orcurvature.The High Altitude Water Cherenkov (HAWC) Obser-vatory is a gamma-ray observatory with an wide instan-taneous field-of-view ( ∼ E > erg/s). This is much higher than the number of high- ˙ E pulsars that would be expected to befound near UHE gamma-ray sources (Abeysekara et al.2020). Emission near pulsars is expected to be domi-nantly leptonic.The proximity of these gamma-ray sources to themost powerful pulsars, along with the curvature in theirspectra, invites the question of whether UHE gamma-ray emission is a generic feature expected near thesesources. This is investigated in this paper through ajoint-likelihood analysis of ten pulsars with ˙ E > erg/s .The paper is organized as follows: Section 2 describesthe details of the joint-likelihood method. Section 3 con-tains the results of the analysis. Section 4 discusses im-plications of the results. ANALYSIS METHOD2.1.
Joint-Likelihood Method
In this paper, we perform a joint analysis of sev-eral high- ˙ E pulsars that do not have significanct detec-tions by HAWC. The analysis uses a binned maximum-likelihood method that searches for a gamma-ray excessabove the background, using a simple model that de-scribes the UHE emission from pulsars based on variousobservables such as the distance and pulsar age. Thebackground in this case comes from the cosmic raysthat are detected by the HAWC detector and survivegamma/hadron separation cuts. This is performed us-ing the HAWC Accelerated Likelihood (HAL) pluginto the Multi-Mission Maximum Likelihood Framework(3ML) (Vianello et al. 2015). Details of the backgroundrejection and likelihood method are described in Abey-sekara et al. (2019). Data collected over 1343 days be-tween June 2015 and June 2019 are used. The data arebinned in quarter-decade bins of reconstructed gamma-ray energy using the “ground parameter” energy esti- https://github.com/threeML/hawc hal eV emission and powerful pulsars dNdE = A i K (cid:18) EE (cid:19) − α , (1)where A i accounts for the model-dependent relative fluxof each source (see Section 2.3), K is the normalization,and E is the pivot energy, which is fixed at the cen-ter of each energy bin. The three values used for E in this analysis are 74.99 TeV, 133.35 TeV, and 237.14TeV. Setting the pivot energy at the center of each binmakes the analysis relatively insensitive to the choice of α , which is fixed at 2.5.Each source (see Table 1) is fit individually in eachenergy bin. A step function is placed at the bound-aries of the energy bin to ensure that events that aremis-reconstructed in energy are excluded. Sources areassumed to be spatially extended with a Gaussian mor-phology. The width is fixed at σ =0.34 ◦ . The values of α and σ are approximately the average values for thehighest-energy gamma-ray sources from the eHWC cat-alog (Abeysekara et al. 2020).To obtain a joint-likelihood result, the log-likelihoodprofiles for each individual source, in each energy bin,are added and the value of K (hereafter called ˆ K ) thatoptimizes this log-likelihood profile is found. The totalflux normalization, κ , for the ten sources combined isthen just: κ = n sources (cid:88) i =1 ˆ KA i . (2)Using this flux normalization, the total differentialflux from all ten sources is: dNdE = κ (cid:18) EE (cid:19) − α . (3)A test statistic (TS) is also calculated to show howsignificant the value of κ is:TS = 2ln L S + B ( κ ) L B , (4)where L S + B is the maximum likelihood for the signal-plus-background hypothesis and L B is the likelihood forthe background-only hypothesis. After an overall best-fit value of κ is determined, it isused as an input to a Markov chain Monte Carlo and adistribution of the value is determined. If the TS in agiven energy bin is significant (TS > σ .Otherwise, 95% credible interval quasi-differential upperlimit is plotted. A uniform prior is used in the Bayesiananalysis.In the bins where an upper limit is determined, therange of expected upper limits are obtained from simu-lations of Poisson-fluctuated background.2.2. Source Selection and Dataset
For this analysis, the source list is defined by select-ing all pulsars from the ATNF pulsar database, version1.62 (Manchester et al. 2005) in the inner Galacticplane that are within HAWC’s field-of-view ( | b | < ◦ ,5 ◦ < l < ◦ ) and have ˙ E > erg/s. There are 24pulsars that meet this criteria.This list is then downselected. First, all pulsars thatare located within a degree of sources from the eHWCcatalog (Abeysekara et al. 2020) are removed. This re-moves pulsars that already have significant UHE emis-sion detected in their vicinity. Since the gamma-rayemission is modeled as extended in nature, this also re-moves pulsars whose associated emission may overlapwith those known sources, which would require moredetailed modeling.Three additional pulsars are removed from the sourceselection. PSR J1813-1749 is removed because HAWChas added more data since the publication of Abey-sekara et al. (2020) and has detected a new UHE source(eHWC-J1813-176) ∼ ◦ away from the pulsar. PSRJ1913+1011 and PSR J1930+1852 have been removedbecause the known TeV emission in those regions islikely from a SNR, not a PWN (Abdalla et al. 2018a,b).The final source list is composed of ten pulsars. Fig-ure 1 shows HAWC’s >
56 TeV map with these sourceslabeled, while Table 1 contains relevant pulsar parame-ters. 2.3.
Models
The parameter A i in Equation 1 describes the rela-tive contribution each pulsar receives in the analysis.This parameter can be used to test different models ofgamma-ray emission near pulsars. For the models thatrely on pulsar parameters, the relevant quantities aretaken from the ATNF pulsar catalog (Manchester et al. HAWC Collaboration
Figure 1. √ T S for ˆ
E >
56 TeV emission from HAWC, for the field-of-view used in this analysis. A point source is assumed asthe morphology. The locations of the ten pulsars used in this analysis are labeled.PSR name RA ( ◦ ) Dec ( ◦ ) Age ( P P ) (kyr) ˙ E ( × erg/s) Distance (kpc) P (s) ˙ P (ss − )B1800-21 270.96 -21.62 15.8 2.2 4.40 0.134 1.35 × − J1828-1101 277.08 -11.03 77.1 1.6 4.77 0.072 1.48 × − J1831-0952 277.89 -9.87 128 1.1 3.68 0.067 8.32 × − J1833-1034 278.39 -10.57 4.85 34 4.10 0.062 2.02 × − J1838-0655 279.51 -6.93 22.7 5.5 6.60 0.070 4.93 × − J1856+0245 284.21 2.76 20.6 4.6 6.32 0.081 6.21 × − J1928+1746 292.18 17.77 82.6 1.6 4.34 0.069 1.32 × − J1935+2025 293.92 20.43 20.9 4.7 4.60 0.080 6.08 × − B1937+21 294.91 21.58 2.35e5 1.1 3.50 0.002 1.05 × − J2022+3842 305.59 38.70 8.94 30 10.00 0.049 8.61 × − Table 1.
Information on the ten pulsars. All information comes from the ATNF database, version 1.62 (Manchester et al.2005). Here, ˙ E is the spin-down energy loss rate, P is the barycentric period, and ˙ P is the time derivative of the period. • No model : Here, A i for all sources is set equal to1. All sources are treated equally and the emissionis expected to be uncorrelated with pulsar param-eters such as distance. • /d : In this model, A i is set to 1/ d , where d is the distance to the pulsar. This model assumesthat closer sources will produce observable emis-sion. • ˙ E/d : Here, the 1/ d model discussed above ismultiplied by the spin-down power of the pul-sar. Therefore, closer, more-energetic pulsars havemore gamma-ray emission. • Inverse age : In this model, A i is the inverse ofthe spin-down age. This is defined as P/ (2 ˙ P ),where P and ˙ P are the period and the time deriva-tive of the period, respectively. In this model,younger sources are more likely to have detectableemission. • Flux at 7 TeV : Here, A i is the HAWC fluxat 7 TeV. This model assumes that sources thatare bright at multi-TeV energies should also giveoff detectable emission above 56 TeV. The 0.5 ◦ extended source map from the third catalog ofHAWC sources (3HWC) (Albert et al. 2020b) isused to extract these values. A i is computed byaveraging the flux from all pixels within a 0.5 de-gree radius of the pulsar. RESULTSThe results of the joint-likelihood analysis (as de-scribed in Section 2) are given in Table 2. Each of themodels give a total TS value above 9, corresponding toa detection of more than 3 σ . For some models, the TSis much higher. The 1/ d model gives the highest TS:41.3 (6.4 σ ).Figure 2 shows the total flux for this best-case sce-nario. The figures for the other models can be seen inthe Appendix. Table 4, also located in the Appendix,gives the total flux normalizations for each model.3.1. Predicted Gamma-ray Flux from an ArbitraryPulsar eV emission and powerful pulsars Model TS (56 < E <
100 TeV) TS (100 < E <
177 TeV) TS (177 < E <
316 TeV) Total TSNo model 27.93 8.33 1.59 37.891/ d E / d Table 2.
The test statistic for the joint-likelihood analysis for each model. Note that κ is fit individually in each energy bin sothe last column is not a joint TS. Figure 2.
The total combined flux in each of the energybins above 56 TeV for the scenario where the pulsars have arelative contribution defined by 1/ d , where d is the distancebetween the pulsar and the Earth. For the first two energybins (spanning 56 < E <
100 TeV and 100 < E <
177 TeV,respectively), the TS > < E <
316 TeV), there is no significant detection (TS=1.29) so a95% upper limit is plotted. The dark and light grey bandsare 68% and 90% containment for expected upper limits fromPoisson fluctuated background. The dotted lines are system-atic uncertainties on the central value (i.e., the solid red linemay be as low or as high as the dotted line once systematicuncertainties are included). See Section 3.2 for a descriptionof systematic uncertainties.
For each model, we can now calculate the expectedgamma-ray flux above 56 TeV from an arbitrary pul-sar using Equation 1. The values of K , derived fromEquation 2, are given in Table 3.3.2. Systematic Uncertainties
Systematic uncertainties are broken down into twocategories: detector systematics and modeling system-atics. Detector systematics, described in Abeysekaraet al. (2019), stem from mis-modeling of detector quan-tities such as the photomultiplier tube threshold andcharge in simulated Monte Carlo events. Each system-atic is treated independently; the results are added inquadrature to get a total uncertainty. This process isrepeated in each energy bin. Depending on the energy bin and model assumed, the detector systematic rangesfrom 10% to 25%.Modeling systematic uncertainties investigate how theanalysis would change if the emission near the pulsarsis different from what is assumed in the main analysis.Several modeling systematics are considered: • Spectral indices in the power law ( α in Equation1) of 2.0 and 3.0 • Replacing the power-law spectral model (Equation1) with a power-law with an exponential cutoff: dNdE = A i K (cid:18) EE (cid:19) − α e − E/E cut . (5)Here, α is fixed at 2.5 and E cut is fixed at 60 TeV,which are average values from Abeysekara et al.(2020) • Changing the Gaussian width to 0.23 ◦ and 0.45 ◦ .These are ± ◦ at the edges of the field-of-view. This meansthat the TS values presented in this paper may be un-derestimated; the size of this effect is ∼ HAWC Collaboration
Model K (56 < E <
100 TeV) K (100 < E <
177 TeV) K (177 < E <
316 TeV) Units for A i No model 8.47 × − × − - dimensionless1/ d × − × − - kpc − ˙ E/d × − × − - erg m − s − Inverse age 7.19 × − - - yr − Flux at 7 TeV 2.28 × − × − - TeV − cm − s − Table 3.
The proportionality constants needed to calculate the expected gamma-ray flux near a given pulsar. The dimensionsfor A i K are energy − distance − time − . With the value of K known and A i in the units given by the last column, the expectedgamma-ray flux can easily be calculated using Equation 1. K is reported at the pivot energy in each bin. Proportionalityconstants are not given for bins where only an upper limit is reported. Testing of Random Backgrounds
We investigate how often randomly chosen non-sourcelocations give TS values as high as in the study presentedabove. Sets of ten randomly chosen points from the anal-ysis region (4 ◦ < l < ◦ ; | b | < ◦ ) are run through thejoint-likelihood analysis. Points that are within a degreeof either a known high-energy source or any of the se-lected ATNF pulsars are excluded. Because the gamma-ray emission is assumed to be spatially extended, this isnecessary to avoid contributions from the known sourcesand the pulsars of interest.This analysis is performed for 2,877 sets of ten ran-dom source locations. Only five of these randomly cho-sen sets have a TS greater than 37.89. This is the TSfrom the “no model” case and is used as the comparisonbecause the other models require known pulsar informa-tion. This means that the results that we have obtainedare significant at the 3.12 σ level.We also investigate sets of ten randomly chosensources from HAWC’s third catalog (3HWC) (Albertet al. 2020b) to see if any of HAWC’s previously de-tected TeV sources emit at UHE. Once again, sourcesthat are known to emit above 56 TeV and sources withina degree of the ATNF pulsars used in the nominal anal-ysis are removed.Due to the relatively small number of remaining3HWC sources after this downselection (48 of the orig-inal 65 3HWC sources), it is not possible to run thou-sands of sets of ten sources without repeating a largesubset of the sources. We instead run 100 sets of tenrandomly chosen 3HWC sources.None of these 100 trials have a TS above the valuefrom the joint-likelihood unmodeled case. The highestTS from this study is 26.9 and the mean is 12.2 (stan-dard deviation 5.5). This implies that known gamma-ray sources that are not located near high- ˙ E pulsars areunlikely to emit at high energies, or if they do, theirfluxes are very low and combining ten sources is notenough for a detection. We also explore combining 20 3HWC sources at atime, instead of the ten that were combined in the pre-ceding paragraph. The entire TS distribution is shiftedto higher values. The highest TS is 45.9 (higher thanin the “no model” case); the mean is 24.6 and the stan-dard deviation is 7.2. The higher TS values when more3HWC sources are combined implies that UHE emissionmay be a generic feature of known gamma-ray sources,but with a very low flux. This should be investigated,but may be hard to do with current-generation experi-ments. Proposed experiments such as SWGO (Huente-meyer et al. 2019), will be important here. DISCUSSIONRegardless of which model is used, it is clear that theareas around high- ˙ E pulsars show hints of emitting atultra-high energies ( >
56 TeV). Interestingly, some mod-els based on the pulsar parameters, such as ˙ E or inverseage, have much lower TS values than the “no model”case, implying that they do not describe the emissionwell. Two models that perform better than the “nomodel” case are the gamma-ray flux at 7 TeV and 1/ d .It is unclear at this time why some pulsar parametersdo not seem to be good predictors of emission. Some ofthese parameters have fairly large uncertainties, whichcould be a contributing factor. For example, the char-acteristic age of the pulsar can differ from the true age.However, the uncertainties on the pulsar distances arerelatively small. For the ATNF pulsar database, 95% ofpulsars will have a relative error of less than a factor of0.9 in their distance estimate (Yao et al. 2017).Alternatively, the individual environment that eachpulsar is in could play a large factor in the amount ofUHE emission. While the pulsar itself is the particle ac-celerator, diffusion of electrons and positrons away fromthe pulsar is strongly dependent on quantities such asthe density and magnetic field of the environment. Also,the emission mechanisms are still uncertain. While it iscommonly assumed that the bulk of emission from thevicinity of a pulsar is from a PWN and therefore pre- eV emission and powerful pulsars a priori discarded. Some have raised the possibility ofhadronic emission mechanisms in or near PWN (Am-ato et al. 2003; Di Palma et al. 2017). Several of theHAWC sources known to emit above 56 TeV, most no-tably eHWC J1908+063 and eHWC J1825-134, havemolecular clouds nearby (Duvidovich et al. 2020; Voisinet al. 2016). These molecular clouds may be serving asa target for a portion of the gamma-ray production.Multi-wavelength and multi-messenger campaigns canhelp disentangle emission mechanisms. A neutrino de-tection coincident with one of these pulsars would bea smoking gun for hadronic emission mechanisms in ornear PWN. However, a recent stacked analysis lookingfor neutrino emission from PWN by IceCube did notyield a detection (Aartsen et al. 2020).Electromagnetic multi-wavelength studies could alsobe helpful. A leptonic source emitting above 56 TeVwill have a vastly different signature at lower energiesthan a hadronic one. For example, 100 TeV gamma raysimply a synchrotron peak in the keV regime, assuming a3 µ G field (Hinton & Hofmann 2009), with the emissionextending up to the MeV energy range. If the emissionis instead predominantly hadronic, there will be no suchpeak at these energies. Proposed experiments such asAMEGO (McEnery et al. 2019) will be important indistinguishing emission mechanisms. CONCLUSIONSIn this study, we have searched for UHE gamma-rayemission in the vicinity of pulsars with an ˙
E > erg/s. We find, with high significance ( > σ ) , thatUHE gamma-ray emission is a generic feature in thevicinity of this class of pulsars. 1/ d is the model thatgives the highest TS (i.e., a source that is closer to us is more likely to have observed UHE emission). Otherpulsar parameters such as the characteristic age and the˙ E do not seem to be good predictors of emission. Thisimplies that the environments the pulsars are located inmay play a role in the amount of emission.The TS values obtained are higher than would be ex-pected from combining random points in the Galacticplane. There is the possibility that all known gamma-ray sources emit above this energy threshold, albeit atan extremely low flux.Multimessenger and multiwavelength studies areneeded to disentagle the origin of the UHE emissionfrom these pulsars.We acknowledge the support from: the US Na-tional Science Foundation (NSF); the US Departmentof Energy Office of High-Energy Physics; the Labo-ratory Directed Research and Development (LDRD)program of Los Alamos National Laboratory; Con-sejo Nacional de Ciencia y Tecnolog´ıa (CONACyT),M´exico, grants 271051, 232656, 260378, 179588, 254964,258865, 243290, 132197, A1-S-46288, A1-S-22784,c´atedras 873, 1563, 341, 323, Red HAWC, M´exico;DGAPA-UNAM grants IG101320, IN111315, IN111716-3, IN111419, IA102019, IN112218; VIEP-BUAP; PIFI2012, 2013, PROFOCIE 2014, 2015; the Universityof Wisconsin Alumni Research Foundation; the Insti-tute of Geophysics, Planetary Physics, and Signaturesat Los Alamos National Laboratory; Polish ScienceCentre grant, DEC-2017/27/B/ST9/02272; Coordi-naci´on de la Investigaci´on Cient´ıfica de la UniversidadMichoacana; Royal Society - Newton Advanced Fel-lowship 180385; Generalitat Valenciana, grant CIDE-GENT/2018/034; Chulalongkorn University’s CUni-verse (CUAASC) grant. Thanks to Scott Delay, LucianoD´ıaz and Eduardo Murrieta for technical support.REFERENCES Aartsen, M. G., Ackermann, M., Adams, J., et al. 2020,The Astrophysical Journal, 898, 117Abdalla, H., Abramowski, A., Aharonian, F., et al. 2018a,Astronomy and Astrophysics, 612, A8—. 2018b, Astronomy and Astrophysics, 612, A1Abdalla, H., Aharonian, F., Benkhali, F. A., Ang¨uner,E. O., et al. 2019, Astrononomy and Astrophysics, 621Abeysekara, A., Albert, A., Alfaro, R., et al. 2020, Phys.Rev. Lett., 124, 021102Abeysekara, A. U., Albert, A., Alfaro, R., et al. 2019, TheAstrophysical Journal, 881, 134—. 2017, The Astrophysical Journal, 843, 39 Albert, A., Alfaro, R., Alvarez, C., et al. 2020a, PhysicalReview Letters, 124, 131101—. 2020b, arXiv e-prints, arXiv:2007.08582Amato, E., Guetta, D., & Blasi, P. 2003, Astronomy andAstrophysics, 402, 827Amenomori, M., Bao, Y. W., Bi, X. J., et al. 2019, PhysicalReview Letters, 051101,doi:10.1103/PhysRevLett.123.051101Di Palma, I., Guetta, D., & Amato, E. 2017, TheAstrophysical Journal, 836, 159Duvidovich, L., Petriella, A., & Giacani, E. 2020, MNRAS,491, 5732
HAWC Collaboration
Hinton, J., & Hofmann, W. 2009, The Annual Review ofAstronomy and Astrophysics, 47, 523Huentemeyer, P., Abreu, P., Albert, A., et al. 2019, inBulletin of the American Astronomical Society, Vol. 51,109Manchester, R. N., Hobbs, G. B., Teoh, A., & Hobbs, M.2005, The Astronomical Journal, 129, 1993McEnery, J., van der Horst, A., Dominguez, A., et al. 2019,in Bulletin of the American Astronomical Society,Vol. 51, 245Moderski, R., Sikora, M., Coppi, P. S., & Aharonian, F.2005, Monthly Notices of the Royal AstronomicalSociety, 363, 954 Smith, A. J. 2015, Proceedings of Science (34thInternational Cosmic Ray Conference), 966Vianello, G., Lauer, R. J., Younk, P., et al. 2015,Proceedings of Science (34th International Cosmic RayConference), 1042Voisin, F., Rowell, G., Burton, M. G., et al. 2016, MNRAS,458, 2813Wilks, S. S. 1938, The Annals of Mathematical Statistics,9, 60Yao, J. M., Manchester, R. N., & Wang, N. 2017, TheAstrophysical Journal, 835, 29 eV emission and powerful pulsars
Figure 3.
Identical to Figure 2 from the main text, but for the “no model” case.
Figure 4.
Identical to Figure 2 from the main text, but for the ˙
E/d model. Figure 5.
Identical to Figure 2 from the main text, but for the inverse age model. HAWC Collaboration
Figure 6.
Identical to Figure 2 from the main text, but for the model defined by the gamma-ray flux at 7 TeV.Model κ (56 < E <
100 TeV) κ (100 < E <
177 TeV) κ (177 < E <
316 TeV)No model 8.47 +1 . − . × − +0 . − . × − × − d +1 . − . × − +0 . − . × − × − ˙ E / d +1 . − . × − +4 . − . × − × − Inverse age 3.99 +1 . − . × − × − × − Flux at 7 TeV 7.80 +1 . − . × − +0 . − . × − × − Table 4.
The total combined flux normalization for each of the models ( κ from Equation 3). The units are TeV − cm − s −1