Investigating the nature of MGRO J1908+06 with multiwavelength observations
Jian Li, Ruo-Yu Liu, Emma de Ona Wilhelmi, Diego F. Torres, Qian-Cheng Liu, Matthew Kerr, Rolf Buehler, Yang Su, Hao-Ning He, Meng-Yuan Xiao
IInvestigating the nature of MGRO J1908 +
06 with multiwavelengthobservations
Jian Li † , , Ruo-Yu Liu ‡ , Emma de O˜na Wilhelmi § , Diego F. Torres , , , Qian-Cheng Liu ,Matthew Kerr , Rolf B ¨ u hler , Yang Su , Hao-Ning He , Meng-Yuan Xiao ABSTRACT
The unidentified TeV source MGRO J1908 +
06, with emission extending fromhundreds of GeV to beyond 100 TeV, is one of the most intriguing sources in theGalactic plane. MGRO J1908 +
06 spatially associates with an IceCube hotspot ofneutrino emission, though not significant yet, indicating a possible hadronic originof the observed gamma-ray radiation. Here we describe a multiwavelength analysison MGRO J1908 +
06 to determine its nature. We identify, for the first time, anextended GeV source as the counterpart of MGRO J1908 +
06, discovering possiblyassociated molecular clouds (MCs). The GeV spectrum shows two well-di ff erentiatedcomponents: a soft spectral component below ∼
10 GeV, and a hard one ( Γ ∼ . − CAS Key Laboratory for Research in Galaxies and Cosmology, Department of Astronomy, University of Scienceand Technology of China, Hefei 230026, China † [email protected] School of Astronomy and Space Science, University of Science and Technology of China, Hefei 230026, China Department of Astronomy, Nanjing University, 163 Xianlin Avenue, Nanjing 210023, China ‡ [email protected] Deutsches Elektronen-Synchrotron DESY, D-15738 Zeuthen, Germany § [email protected] Institute of Space Sciences (ICE, CSIC), Campus UAB, Carrer de Can Magrans, 08193, Barcelona, Spain Institut d’Estudis Espacials de Catalunya (IEEC), 08034 Barcelona, Spain Instituci´o Catalana de Recerca i Estudis Avanc¸ats (ICREA), E-08010, Barcelona, Spain Space Science Division, Naval Research Laboratory, Washington, DC 20375, USA Purple Mountain Observatory and Key Laboratory of Radio Astronomy, Chinese Academy of Sciences, Nanjing210034, China Key Laboratory of Dark Matter and Space Astronomy, Purple Mountain Observatory, Chinese Academy ofSciences, 210023 Nanjing, Jiangsu, China a r X i v : . [ a s t r o - ph . H E ] F e b Subject headings:
Gamma-rays: pulsars: individual (PSR J1907 +
1. Introduction
MGRO J1908 +
06 is an extended bright TeV source of unknown nature. It was first discoveredby the Milagro water Cherenkov telescope in its sky survey results after seven years of operation(Abdo et al. 2007). MGRO J1908 +
06 was subsequently detected in the TeV range by the HighEnergy Stereoscopic System (H.E.S.S., Aharonian et al. 2009), the Very Energetic RadiationImaging Telescope Array System (VERITAS, Ward 2008; Aliu et al. 2014), the AstrophysicalRadiation with Ground-based Observatory at YangBaJing (ARGO-YBJ, Bartoli et al. 2012), andrecently by the High-Altitude Water Cherenkov Observatory (HAWC, Abeysekara et al. 2017).The TeV luminosity of MGRO J1908 +
06 is comparable to the Crab Nebula (Bartoli et al. 2012),making it one of the most luminous Galactic gamma-ray sources in the TeV range. MGROJ1908 +
06 is among the four sources detected above 100 TeV in a recent study by HAWC,indicating its possible ability to accelerate particles to PeV energy (Abeysekara et al. 2020).Besides its possible PeVatron nature, MGRO J1908 +
06 is spatially associated with an IceCubeneutrino emission hotspot (Aartsen et al. 2019; Aartsen et al. 2020), though the post-trialsignificance is low.The nature of MGRO J1908 +
06 remains unrevealed. Searches for its multiwavelengthcounterparts in radio, X-ray and GeV gamma rays have been unsuccessful (Abdo et al. 2010;Kong et al. 2007; Pandel 2015; Duvidovich et al. 2020). MGRO J1908 +
06 is spatially associatedwith a middle-aged (20-40 kyr, Downes et al. 1980) supernova remnant (SNR) G40.5 − + ∼ . × erg s − and a characteristic age of 19.5 kyr. The distance toPSR J1907 + ± − Σ -D relation, 5.5 to 8.5kpc, Downes et al. 1980; 6.1 kpc, Case & Bhattacharya 1998). PSR J1907 + − + − . The nature of MGRO J1908 +
06 is under debate. A leptonic pulsar wind nebula using the electron-density model of Yao et al. (2017) ( http://119.78.162.254/dmodel/index.php ), the + +
06 and modelled theemission processes, and we report our results below.
2. Multiwavelength observations
Multiwavelength observations were investigated in this paper, including CO observationsfrom the Milky Way Imaging Scroll Painting (MWISP) project (Su et al. 2019),
Fermi
LargeArea Telescope (LAT) observations of GeV gamma rays,
XMM-Newton observations in X-ray, andthe Very Large Array Galactic Plane Survey (VGPS, Stil et al. 2006) in radio. The details of theseobservations are reported in Appendices.
3. Results3.1. Molecular clouds in the region of MGRO J1908 + Molecular clouds (MCs) towards the SNR G40.5 − CO ( J = CO( J = O ( J = +
06. MCs are discovered to be spatially associated withMGRO J1908 +
06 in the CO ( J = CO ( J = / s (Figure1; also see Appendix A for details), which is consistent with the possible distance, 3 . − + − CO ( J = CO ( J = CO ( J = CO ( J = − to 66 km s − with a coverage of 4 km s − in Appendix A , Figure 5. Di ff erent regionswith MCs are studied in detail and their astrophysical properties are estimated in Appendix A . Themean density is estimated to be ∼
45 cm − . In the 46 km s − to 66 km s − velocity slices, thereare apparent filaments of MCs positionally located next to PSR J1907 + CO (J = −
0) and CO (J = −
0) line profiles of the MCs toward the SNR G40.5 − distances of PSR J1907 + + + CO line (i.e., a wing part deviating from the main Gaussianin a CO profile where the signal-noise-ratio of the CO is less than 3) in the grid of the COspectra. For a further search of interaction signals, we studied the data from the INT PhotometricH α Survey of the Northern Galactic Plane (IPHAS) covering the MGRO J1908 +
06 region. No H α ascribable to interaction was detected. + Though being very bright at TeV energies, no GeV counterpart for MGRO J1908 + Fermi -LAT studies (Abdo et al. 2010a; Ackermann et al. 2011).Two gamma-ray sources PSR J1907 + + +
06, as listed in the 4FGL catalog (Abdollahi et al. 2020). For a deeper searchof its GeV counterpart, we analyzed more than 11 years of
Fermi -LAT data in the 0.1–300 GeVband. To minimize the contamination from the bright gamma-ray pulsar PSR J1907 + ff -peak phases of this pulsar. A full description of our dataanalysis is presented in Appendix B . We discovered a previously-undetected, extended GeV sourcespatially associated with MGRO J1908 +
06, hereafter referred as Fermi J1906 + + +
06 strongly suggests a common origin.To investigate this, we derived the morphological and spectral properties of Fermi J1906 + + Fermipy , and the disk morphology yields more significant extensiondetection. The disk morphology resulted in a center position of R.A. = ◦ ± ◦ , Decl. = ◦ ± ◦ , and a radius of 0.83 ◦ ± ◦ , yielding a TS ext =
44 (see
Appendix B ). Comparingthe likelihood values, the disk morphology is significantly preferred over a two-point-sourcemodel (PSR J1907 + + ∆ TS =
23 (see AppendixB). We further considered a two-component spatial model, constituted by a template usingVERITAS counts map (Aliu et al. 2014) plus the point source 4FGL J1906.2 + ∆ TS =
18. Additionally, in this two-component spatial model4FGL J1906.2 + =
25, 4FGL). 5 –We carried out the spectral analysis in 0.1-300 GeV. Adopting the disk morphology, wetested a power-law ( dN / dE = N ( E / E ) − Γ cm − s − MeV − ) and a log-parabola spectral model( dN / dE = N ( E / E b ) − ( α + β log ( E / E b )) cm − s − MeV − ). The log-parabola spectral model is preferredwith a ∆ TS =
89, indicating a significant spectral curvature. With the the disk morphology andlog-parabola spectral model, Fermi J1906 + ± × − erg cm − s − between 0.1 −
300 GeV . Adopting thedisk morphology, the best-fitted parameters are shown in Table 1 and the gamma-ray spectralenergy distribution (SED) of Fermi J1906 + ∼ Fermi -LAT data above 30 GeV further without pulsar gating.With a spectral cuto ff at 2.9 GeV, the magnetospheric emission from PSR J1907 + + = ◦ ± ◦ , Decl. = ◦ ± ◦ with aradius of 0.51 ◦ ± ◦ , yielding a TS of 47 assuming a power-law spectral model and a TS ext = + + ± ± × − erg cm − s − (Table 1), confirming the existence of an additional spectralcomponent at higher energies beyond the log-parabola component (Figure 2, top right panel).We added this component to the analysis of Fermi J1906 + ff -peakdata, fixed it and constructed a two-component spectral model (log-parabola plus power law). Theresults are reported in Table 1 and shown in Figure 2, top right panel. The multiwavelength SEDof MGRO J1908 +
06 is shown in Figure 3. The GeV SED shows a natural continuity to the TeVrange. The TS map in the 30 GeV to 1 TeV energy range is shown in Fig. 2, bottom right panel.No clear CO ( J = + = ◦ ± = ◦ ± ◦ ± CO ( J = To explore the influence of di ff erent Galactic di ff use emission models, we tested the previous Galactic di ff useemission component (“gll iem v06.fits”; Aecro et al. 2015) in our background spectral-spatial model. FermiJ1906 + ± × − erg cm − s − between 0.1 −
300 GeV. − range, suggesting a hadronic origin. We produced a template from the 62 – 66 km s − CO ( J = +
06 counterparts in X-ray and radio wavelengths
The
XMM-Newton
X-ray satellite has covered MGRO J1908 +
06 with 5 observations (ObsID 0553640101, 0553640201, 0553640701, 0553640801 and 0605700201), providing a totalexposure of 109 ks. Combining all available MOS data, we produced a particle background-subtracted and exposure-corrected count rate map for the MGRO J1908 +
06 region (Figure 4, leftpanel). No X-ray emission coincident with the gamma-ray emission region could be detectedin the 0.2 to 10 keV energy range. The morphology measured by H.E.S.S. (a Gaussian profilewith σ of 0.34 ◦ ; Aharonian et al. 2009) is the only one in TeV range which is fully covered bythese XMM-Newton observations. Adopting the H.E.S.S. morphology , excluding point sourceswithin, assuming a spectral index of 2 and a H i column density of N H = × cm − (Abdo etal. 2010a), we calculated a 95% unabsorbed upper limit for MGRO J1908 +
06 as 1.2 × − erg cm − s − (0.2 −
10 keV). From our CO observation, assuming a distance of 3.2 kpc, the local H i column density to the MGRO J1908 +
06 region is N H = × N H ∼ (1 − × cm − ( Appendix A ,Table 2), which is an order of magnitude lower than the total Galactic H i column density in thisdirection estimated using the HEASARC tool and with Chandra (N H = × cm − ; Abdo etal. 2010a), ruling out local absorption as the cause of any lack of X-ray counterpart.With a 19 ks Chandra observation (ObsID 7049 on 2009-08-19), Abdo et al. (2010a) reportedthe detection of PSR J1907 + − + XMM-Newton data, Obs. ID 0605700201 (on 2010-04-26) provided the longest exposure (52.5ks) and it is the only observation in which PSR J1907 + Chandra (Abdo et al. 2010a), we produced counts map from 2 − + ∼ (Figure 4, top right). The detected emission is thus consistent with apoint source and no associated PWN could be identified. https://heasarc.gsfc.nasa.gov/cgi-bin/Tools/w3nh/w3nh.pl https://xmm-tools.cosmos.esa.int/external/xmm_user_support/documentation/uhb/ .800 287.400 287.000 286.600 286.200 285.8007.2007.0006.8006.6006.4006.2006.0005.800 PSR J1907+0602 4FGL J1906.2+0631
SNR G40.5-0.5 .800 287.400 287.000 286.600 286.200 285.8007.2007.0006.8006.6006.4006.2006.0005.800
PSR J1907+0602 4FGL J1906.2+0631
SNR G40.5-0.5
Fig. 1.— CO ( J = CO ( J = − . 4FGL J1906.2 + + − σ to 6 σ by 1 σ steps. The x and y axes are R.A. and decl. (J2000) in degrees.Table 1: Best-fit Spectra Parameters of Fermi J1906 + O ff -peak Analysis from 0.1-300 GeVModel Γ α β Energy Flux TS(Power law) (Log-parabola, E b = b = − erg cm − s − )Log-parabola – 2.90 ± ± ± ◦ disk)Log-parabola + Power law (fixed) 1.61 (fixed) 3.16 ± ± ± ◦ disk) (0.51 ◦ disk) Analysis from 30 GeV–1 TeVModel Γ α β Energy Flux TS(Power law) (Log-parabola, E b = b = − erg cm − s − )Power law 1.61 ± ± ◦ disk) SNR G40.5-0.5
PSR J1907+0602 4FGL J1906.2+0631
PS 1
Energy (MeV) ) - s - d N / d E ( e r g c m E − − − − PSR J1907+0602 4FGL J1906.2+0631
PS 1
PSR J1907+0602 4FGL J1906.2+0631
PS 1
Fig. 2.—
Top left : Fermi -LAT TS map of MGRO J1908 +
06 region in the 0.1–300 GeVrange. Background 4FGL sources are shown with white crosses. The disk morphology of FermiJ1906 + + + Top right :GeV SED of Fermi J1906 + + ff -peak and phase averaged analyses, respectively (see text for detail). The dashed red and bluecurves on two instances indicate the two-component spectral modelling (Table 1) with the sumshown with a black line. Bottom left : Fermi -LAT TS map of MGRO J1908 +
06 region in the0.1–2 GeV range. The black contours correspond to TS values staring from 25 with a step sizeof 5. The green contours correspond to the CO ( J = − intensity map of thesurrounding region, starting from 15K with a step of 5K. The labels are as in the top left panel. Bottom right : Fermi -LAT TS map of MGRO J1908 +
06 region from 30 GeV–1 TeV smoothedwith 0.3 degree Gaussian. The labels are as in the top left panel. 9 –Fig. 3.— Multi-wavelength SED of MGRO J1908 +
06 with hadronic and leptonic hybrid modeling.Besides the GeV and X-ray measurements in this paper, other data are taken from Abdo et al.(2007) (Milagro), Aharonian et al. (2009) (HESS), and Abeysekara et al. (2020) (HAWC). TheVERITAS SED data are consistent with HESS and are not shown for concision. The solid greycurve shows the LHAASO point-source sensitivity of one-year exposure (Bai et al. 2019). Thedashed grey curve represents the IceCube 8-yr point-source sensivity for muon neutrino (Aarstenet al. 2019). 10 –We checked for variability with spectral analysis of the eight-month separated
XMM-Newton (Obs. ID 0605700201 on 2010-04-26) and
Chandra (ObsID 7049 on 2009-08-19) observations.We adopted a simple power-law model plus absorption with H i column density fixed at N H = × cm − (Abdo et al. 2010a) because of the low overall counts. The 1 −
10 keV
XMM-Newton
MOS 1&2 data yield a best-fit spectral index of 1.79 + . − . with unabsorbed total energyflux of 6.29 + . − . × − erg cm − s − , whereas the Chandra
ACIS spectrum is well fitted withcompatible spectral parameters, a best-fit spectral index of 1.28 + . − . and unabsorbed total energyflux of 5.14 + . − . × − erg cm − s − . No spectral or flux variability can be claimed. The XMM-Newton and
Chandra spectra are shown in
Appendix C , Figure 8.During this
XMM-Newton observation, the PN camera was operating in small window mode,providing su ffi cient time resolution (5.7 ms) to search for X-ray pulsations. We extract photonsfrom PN data using a radius of 20 arcsec in 0.2 −
10 keV and 3 −
10 keV with barycenter correctionapplied. Using
Tempo2 (Hobbs et al. 2006) and the photons plugin and contemporaneous Fermi -LAT gamma-ray ephemeris, we have assigned pulsar rotational phase to each extracted photon. Nosignificant X-ray pulsation is detected.The VGPS at 1420 MHz (Stil et al. 2006) has covered MGRO J1908 +
06 region and theimage is shown in Figure 4, bottom panel. PSR J1907 + − ff use radiolarge-scale emission associated with MGRO J1908 +
06 is seen.
4. Discussion
Detailed analysis of the
Fermi -LAT data revealed that the gamma rays from the direction ofMGRO J1908 +
06 follow a two-component spectrum. Based on the TS maps and multiwavelengthobservations, we propose two accelerators on the field: one region related to SNR G40.5 − + Fermi -LAT data shows a hard spectrum of α he ∼ . +
06 with a PWN powered byPSR J1907 + >
30 GeV) component measured by
11 –
SNR G40.5-0.5PSR J1907+0602
PSR J1907+0602
10 12 14 16 18 20
SNR G40.5-0.5PSR J1907+0602 4FGL J1906.2+0631
PS 1
Fig. 4.—
Top left : Gaussian-smoothed ( σ = +
06 region from 0.2 −
10 keV,combining all available
XMM-Newton
MOS 1 & 2 data. White contours are the VERITASsignificance map, as in Figures 1 and 2. The labels are as in Figure 2.
Top Right : Gaussian-smoothed ( σ = XMM-Newton
MOS 1 & 2 combined counts map of PSR J1907 + − ∼
90% fractional encircled energyof MOS 1 and MOS 2, which should contain ∼
90% of the counts from a point source.
Bottom :VGPS 1420 MHz image of the MGRO J1908 +
06 region. Legends are the same as the first panel. 12 –
Fermi -LAT, combined with the very-high energy spectrum measured by Cherenkov instruments, doindeed resemble the spectral signature associated with inverse Compton emission from GeV / TeVPWNe (e.g. Crab, Abdo et al. 2010b; MSH 15-52, Abdo et al. 2010c; HESS J1825-137, Grondinet al. 2011). The low-energy component (below a few GeV) is described by a soft spectrum,similar to the ones observed on evolved SNRs (Acero et al. 2016) , and it shows a significant peakcoincident with an enhancement of molecular material (see Fig. 2 Bottom Left and results shownin the Appendix ), implying a tentative hadronic origin. Within 50 pc region of SNR G40.5 − + + ∼ + ∼ pp collision, the gamma-ray flux is calculatedusing the parameterized formulae provided by Kamae et al. (2006), for a proton spectrum of theform of a power-law function in energy space, with a slope s p . Assuming a distance of 8 kpc,the steep gamma-ray spectrum at low energy can be well-reproduced by employing s p = . ff usion-modifiled spectrum, i.e., the injection spectrum is softened by theenergy-dependent di ff usion of CR protons with a di ff usion coe ffi cient D ( E ) ∝ E δ . We do notfurther specify the value of the injection spectral slope of CR protons and the index of the di ff usioncoe ffi cient δ , but just note that an injection spectral slope of 2 . − . δ = . − . n =
45 cm − in the surrounding MC, resulting in a pp collisioncooling timescale of t pp (cid:39) (cid:16) . σ pp nc (cid:17) − = × ( n / − ) − yr. The total proton energy neededto account for the gamma-ray emission is W p (cid:39) . × ergs, which is well consistent with thereach of the usual 10% of the kinetic energy released in SNRs (Aharonian et al. 2004). Note that ifthe distance of the SNR and MCs is 3.2 kpc, the inferred gas density of the MC would be 2.5 timeshigher (see Table 2), and this change reduces the requirement for the proton energy to 10 ergs. SNR analysis in Acero et al. (2016) starts from 1GeV. Fermi J1906 + ±
13 –For the leptonic component, we used an electron / positron broken-power law distribution, i.e., dN / dE e ∝ E − s e , e for E e < E b and dN / dE e ∝ E − s e , e for E e ≥ E b , which is usually chosen todescribe the SED of PWNe (Tanaka & Takahara 2010; Bucciantini et al. 2011; Martin et al. 2012;Torres et al. 2013; Torres et al. 2014). The time-independent spectra of synchrotron radiationand IC radiation are calculated following Blumenthal (1970). The IC emissivity is calculatedin the optically thin case, which is appropriate for these objects, using the general Klein-Nishinadi ff erential cross-section. We adopt the interstellar radiation field modelled in Popescu et al. (2017)as well as the cosmic microwave background as the target photon field for the IC scattering. The Fermi -LAT data above 10 GeV together with the HAWC data can be reproduced with s e , = . s e , = . E b = . W e (cid:39) × erg in theemitting electrons / positrons . Assuming the age of the system equal to the characteristic age ofPSR J1907 + t cool (cid:39) E / . − (cid:104) (U ph + U B ) − / − (cid:105) − kyr. Given U ph = − and U B = B / π , the inferred magnetic field strength is B (cid:39) µ G . Such amagnetic field is consistent with the X-ray upper limit posed by XMM-Newton . The low magneticfield strength is similar to some other relic nebulae in the TeV regime (Aharonian et al. 2006;H.E.S.S. Collaboration 2012; Liu et al. 2019) associated to intermediate-aged pulsars. Some ofthose PWNe display energy-dependent morphology in the TeV regime (Aharonian et al. 2006;H.E.S.S. Collaboration 2012; H.E.S.S. Collaboration 2019). Above 30 GeV, the best fit to the LATdata is the template adopted from VERITAS data, consistent with the morphology obtained in TeVrange. If the emission is indeed related to the pulsar, we expect the nebula to be more compactand closer to the pulsar at TeV energies. The energy-dependent morphology could be the key tounderstand the transport mechanism of particles within the PWN and the evolution of the PWN(H.E.S.S. Collaboration 2019; Liu & Yan 2020). Deep observations with TeV observatories suchH.E.S.S., HAWC or LHAASO will provide crucial input to disentangle the origin of the gamma-ray emission observed.We also note that there are actually many relevant physical processes which can influence themodelling, such as the particle injection history, particle spectral evolution and particle transport.We leave such a more realistic modelling to the future study and here we simply test the feasibilityof the hybrid interpretation. In the considered hybrid scenario, the neutrino emission fromMGRO J1908 +
06 would not be detectable by current instruments which are operating above 100GeV (e.g. IceCube). This is because the neutrino spectrum arising from pp collisions generallyresemble that of the pionic gamma-ray spectrum, which is important only below ∼
10 GeV and Note that the required total energy is not sensivite to the chosen maximum energy and minimum energy of thespectrum given s e , < s e , >
14 –drops quickly with energy. However, we should also note that it is not clear for the time beingthat whether the gamma-ray spectrum above 100 TeV would decline as our model expectation.If a hardening of the gamma-ray spectrum beyond 100 TeV presents in the future observation byHAWC or LHAASO, a hadronic gamma-ray component as well as neutrino detection may then beexpected.The
Fermi
LAT Collaboration acknowledges generous ongoing support from a number ofagencies and institutes that have supported both the development and the operation of the LAT aswell as scientific data analysis. These include the National Aeronautics and Space Administrationand the Department of Energy in the United States, the Commissariat `a l’Energie Atomique andthe Centre National de la Recherche Scientifique / Institut National de Physique Nucl´eaire et dePhysique des Particules in France, the Agenzia Spaziale Italiana and the Istituto Nazionale di FisicaNucleare in Italy, the Ministry of Education, Culture, Sports, Science and Technology (MEXT),High Energy Accelerator Research Organization (KEK) and Japan Aerospace Exploration Agency(JAXA) in Japan, and the K. A. Wallenberg Foundation, the Swedish Research Council and theSwedish National Space Board in Sweden. Additional support for science analysis during theoperations phase is gratefully acknowledged from the Istituto Nazionale di Astrofisica in Italyand the Centre National d’ ´Etudes Spatiales in France. This work performed in part under DOEContract DE-AC02-76SF00515.J. L. acknowledges the support from the Alexander von Humboldt Foundation and theNational Natural Science Foundation of China via NSFC-11733009. R.-Y. L. acknowledges thesupport from the National Natural Science Foundation of China via NSFC-U2031105. E. O.W. acknowledges the support from the Alexander von Humboldt Foundation. The work of D.F. T. has been supported by the grants PGC2018-095512-B-I00, SGR2017-1383, and AYA2017-92402-EXP. Q. C. L. acknowledges support from the program A for Outstanding PhD candidateof Nanjing University. Work at NRL is supported by NASA.
REFERENCES
Aartsen M. G. et al. 2019, Eur. Phys. J. C, 79, 234Aartsen M. G. et al. 2020, PRL, 124, 051103Abeysekara A. U. et al. 2017, ApJ, 843, 40 15 –Abeysekara A. U. et al. 2020, ApJ, PRL, 124, 021102Abdo, A. A. et al. 2009, ApJ 700, 1059Abdo, A. A. et al. 2007, ApJ 664, L91Abdo, A. A. et al. 2010a, ApJ 711, 64Abdo, A. A. et al. 2010b, ApJ, 708, 1254Abdo, A. A. et al. 2010c, ApJ, 714, 927Abdo, A. A. et al. 2010d, ApJ, 718, 348Abdo, A. A. et al. 2013, ApJS, 208, 17Abdollahi, S . et al. 2020, ApJS, 247, 33Acero, F. et al. 2015, ApJS, 218, 23Acero, F. et al. 2016, ApJS, 224, 8Ackermann, M., Ajello, M., Baldini, L. et al. 2011, ApJ, 726, 35Ackermann, M., Ajello, M., Allafort, A., et al. 2013, Science, 339, 807Aguilar, M., Ali Cavasonza, L., Ambrosi, G., et al. 2016, Phys. Rev. Lett., 117, 231102Aharonian, F. A. & Atoyan, A. M. 1996, A&A, 309, 917Aharonian, F. A., Akhperjanian, A. G., Aye, K.-M., et al. 2004, Nature, 432, 75Aharonian, F., Akhperjanian, A. G., Bazer-Bachi, A. R., et al. 2006, A&A, 460, 365Aharonian F. et al. 2009, A&A, 499, 723Aliu E. et al. 2014, ApJ, 711, 64Bai, X., Bi, B. Y., Bi, X. J., et al. 2019, arXiv e-prints, arXiv:1905.02773Bartoli B. et al. 2012, ApJ, 760, 110Blumenthal, G. R., & Gould, R. J. 1970, Reviews of Modern Physics, 42, 237Brand, J., & Blitz, L. 1993, A&A, 275, 67Bucciantini, N., Arons, J., & Amato, E. 2011, MNRAS, 410, 381 16 –Case, G. L., & Bhattacharya, D. 1998, ApJ, 504, 761Celli, S., Morlino, G., Gabici, S., et al. 2019, MNRAS, 490, 4317Duvidovich, L., Petriella, A. & Giacani, E. 2020, MNRAS, 491, 5732Downes A. J. B., Pauls T., Salter C. J., 1980, A&A, 92, 47Frerking, M. A., Langer, W. D., & Wilson, R. W. 1982, ApJ, 262, 590G´enolini, Y., Boudaud, M., Batista, P.-I., et al. 2019, Phys. Rev. D, 99, 123028Grondin, M.-H., Funk, S., Lemoine-Goumard, M., et al. 2011, ApJ, 738, 42H. E. S. S. Collaboration, Abramowski, A., Acero, F., et al. 2012, A&A, 548, A46H. E. S. S. Collaboration, Abdalla, H., Aharonian, F., et al. 2019, A&A, 621, A116Hobbs, G., Edwards, R., & Manchester, R. 2006, ChJAS, 6, 189Huang, Z.-Q., Liu, R.-Y., Joshi, J. C., et al. 2020, arXiv e-prints, arXiv:2001.02973Kamae, T., Karlsson, N., Mizuno, T., et al. 2006, ApJ, 647, 692Kong A. K. H. 2007, ATel.1251Larson, R. B. 1981, MNRAS, 194, 809Li J., Torres D., Liu R.-Y et al. 2020, Nature Astronomy, 4, 1177Liu, Q.-C., Chen, Y., Zhou, P., et al. 2020, ApJ, 892, 143Liu, R.-Y., & Yan, H. 2020, MNRAS, 494, 2618Lyne A. G. et al. 2017, ApJ, 834, 137Mart´ın, J., Torres, D. F., & Rea, N. 2012, MNRAS, 427, 415Pandel D. 2015, arXiv:1512.08140Popescu, C. C., Yang, R., Tu ff s, R. J., et al. 2017, MNRAS, 470, 2539Ray, P. S., Kerr, M., Parent, D., et al. 2011, ApJS, 194, 17Shan, W., Yang, J., Shi, S., et al. 2012, IEEE Transactions on Terahertz Science and Technology,2, 593 17 –Stil J.M. et al. 2006, AJ, 132, 1158Su Y. et al. 2019, ApJS, 240, 9Tanaka, S. J., & Takahara, F. 2010, ApJ, 715, 1248Torres, D. F., Cillis, A. N., & Mart´ın Rodriguez, J. 2013, ApJ, 763, L4Torres, D. F., Cillis, A., Mart´ın, J., et al. 2014, Journal of High Energy Astrophysics, 1, 31Ward, J. for the Veritas collaboration, 2008, arXiv:0810.0664Wood, M. et al., 35th International Cosmic Ray Conference. 10-20 July, 2017. Bexco,Busan, Korea, Proceedings of Science, Vol. 301. Online at https: // pos.sissa.it / cgi-bin / reader / conf.cgi?confid = Appendix A: CO data analysis
The CO data used in this work are part of the Milky Way Imaging Scroll Painting (MWISP)project (Su et al. 2019), including CO ( J = CO ( J = O ( J = × (cid:48)(cid:48) at 115 GHz. The spectral resolution is 61 KHz,corresponding to velocity resolutions of 0 .
16 km s − for CO and 0 .
17 km s − for CO and C O.The typical rms noise level is about 0.5 K for CO ( J = CO ( J = O( J = CO ( J = CO ( J = − to 66 km s − , with a coverage of 4 km s − .We made an estimation of the astrophysical properties of the MCs in 4 regions and haveparameterized the distance as d = . d . kpc. We estimated the kinematic distance to the MCs This preprint was prepared with the AAS L A TEX macros v5.2.
18 –using the Milky Way’s rotation curve suggested by Brand & Blitz (1993), assuming the Sun’sGalactocentric distance to be 8.5 kpc and orbital speed to be 220 km s − . Therefore, the velocityof each MC could indicate two candidate kinematic distances, the near side one and the far sideone. CO( J = CO( J = CO lines, under the assumption of local-thermal-equilibrium (LTE) condition,optically thin conditions for CO ( J = CO ( J = T ex =
12 K, 12 K, 21 K, and 15 K, respectively, thevalue estimated from the maximum CO ( J = N (H ) ≈ × N ( CO) (Frerking et al. 1982) has been used.The mean density of each region is calculated by dividing the column density toward the COemission peaks by the cloud size along the line-of-sight, which is estimated from the FWHM ofthe CO line with Larson’s law (Larson 1981). Furthermore, a mean density of the four regions,when taking the weight of the volume into consideration, is estimated to be ∼
45 cm − . Appendix B:
Fermi -LAT data analysis
The analysis shown in this paper uses more than 11 years of
Fermi -LAT P8R3 data, from2008 August 4 (MJD 54682) to 2019 November 09 (MJD 58796). All gamma-ray photons withina circular region of interest (ROI) of 15 ◦ radius centered on PSR J1907 + Fermi
Science Tools 11-07-00 release. In the data reduction, a zenith angle threshold of 90 ◦ isadopted to reject contamination from gamma rays from the Earth’s limb. The selected Fermi -LAT instrument response functions (IRFs) is “P8R3 V2 Source”. Known gamma-ray sourcesfrom the
Fermi
Large Area Telescope Fourth Source Catalog (4FGL, Abdollahi et al. 2020)within 20 ◦ of PSR J1907 + ff use emission components (“iso P8R3 SOURCE V2 v1.txt”).The spectral parameters of the sources within 4 ◦ of PSR J1907 + ff use emission components were all left free. PSR J1907 + + +
06, thus not included in the spectral-spatial model. Thespectral parameters of sources with larger angular separations were fixed at the 4FGL values. Thespectral analysis was performed using a binned maximum likelihood fit (spatial bin size 0.1 ◦ and30 logarithmically spaced bins in the 0.1–300 GeV range) For the analysis in 30 GeV–1 TeV, thespectral parameters of the sources within 4 ◦ of PSR J1907 + ff use emissioncomponent are fixed to 4FGL values except for the prefactor (spectral normalization) because oflow statistics.The significance of the sources was evaluated by the Test Statistic (TS). This statistic is 19 – PSR J1907+0602
PS 1
SNR G40.5-0.5
PSR J1907+0602
PS 1
SNR G40.5-0.5
PSR J1907+0602
PS 1
SNR G40.5-0.5
PSR J1907+0602
PS 1
SNR G40.5-0.5
Fig. 5.— From top to bottom: CO ( J = CO ( J = − , 50 – 54 km s − , 54 – 58 km s − , 58 – 62 km s − ,62 – 66 km s − . The color denotes the intensity. White contours correspond to the VERITASsignificance map starting from 3 σ with a step of 1 σ . The labels are as in Figure 2. The x and yaxes are R.A. and decl. (J2000) in degrees. 20 – PSR J1907+0602
PS 1
SNR G40.5-0.5
PSR J1907+0602
PS 1
SNR G40.5-0.5
PSR J1907+0602
PS 1
SNR G40.5-0.5
PSR J1907+0602
PS 1
SNR G40.5-0.5
PSR J1907+0602
PS 1
SNR G40.5-0.5
PSR J1907+0602
PS 1
SNR G40.5-0.5 -Figure 5 continues 21 –
43 2 1
PS 1
Velocity (km/s)40 45 50 55 60 65 70 M B ) → K ( T Region 1
Velocity (km/s)40 45 50 55 60 65 70 M B ) → K ( T Region 2
Velocity (km/s)40 45 50 55 60 65 70 M B ) → K ( T Region 3
Velocity (km/s)40 45 50 55 60 65 70 M B ) → K ( T Region 4
Fig. 6.—
Top panel : CO ( J = − ,shown twice. White contours correspond to VERITAS significance map (Aliu et al. 2014) startingfrom 3 σ to 6 σ by 1 σ steps. The green circle shows the disk morphology of Fermi J1908 +
06 (seeSection 3.2). The x and y axes are R.A. and decl. (J2000) in degrees. Four regions delineatedin green and labelled with roman numerals “1” to “4” are used to estimate the astrophysicalparameters for the molecular gases (see Table 2).
Lower panels : CO ( J = CO( J = = − L max , / L max , ), where L max , is the maximum likelihood value for a modelin which the source studied is removed (the “null hypothesis”), and L max , is the correspondingmaximum likelihood value with this source being incorporated. The square root of the TS isapproximately equal to the detection significance of a given source. The significance of sourceextension was defined as TS ext = − L point / L ext ), where L ext and L point are the gtlike globallikelihood of the extended source hypotheses and the point source, respectively. The threshold forclaiming the source to be spatially extended is set as TS ext >
16, which corresponds to a significanceof ∼ σ . The source localization, extension fitting and TS maps production were carried out usingthe Fermipy analysis package (version 0.17.4; Wood et al. 2017). Energy dispersion correction hasbeen applied in the analysis. The SEDs are computed assuming a power-law shape with spectralindex fixed at 2.PSR J1907 + + ff -peakphases of PSR J1907 + Tempo2 (Hobbs et al. 2006) with the
Fermi plug-in (Ray etal. 2011), we have assigned pulsar rotational phases for each gamma-ray photon that passed theselection criteria, adopting the most updated ephemeris for PSR J1907 + + ff -peak definition in Li et al. (2020),which is φ = − − ff -peak phase selection, the prefactorparameter of all sources were scaled by 0.439, the width of the o ff -peak interval.In the o ff -peak analysis of the 0.1–300 GeV band, we searched for significant TS excessbeyond Fermi J1906 + + = ◦ ± ◦ , decl. = ◦ ± ◦ (Figure 2). Assuminga power-law spectral shape, the likelihood analysis of PS 1 resulted in a TS value of 31, spectralindex of 2.27 ± ± × − erg cm − s − .Since VERITAS observations have the deepest exposure on MGRO J1908 +
06 in TeV rangeand provided most detailed TeV morphology, we included the VERITAS counts map as a templatein our GeV morphology analysis. HAWC observations provided the only SED data points onMGRO J1908 +
06 above 50 TeV. Thus HAWC data are adopted in the multi-wavelength SEDmodelling. We noticed that in ∼ + Appendix C: X-ray data analysis
XMM-Newton data sets were reduced with the Science Analysis System (SAS, version 23 –Table 2: Fitted and Derived Parameters for the MCs around 50 km s − in 4 regions as indicated inFig. 6Region N (H ) n (H ) M (H ) FWHM Line Center Near / Far Distance(10 cm − ) (cm − ) (10 M (cid:12) ) ( km s − ) ( km s − ) (kpc)1 1.2 16 d − . d . / d − . d . / d − . d . / d − . d . / Pulse Phase C oun t s off peak no r m a li z e d c oun t s Fig. 7.— Pulse profile of PSR J1907 + ◦ . Two rotationalpulse periods are shown, with a resolution of 100 phase bins per period. The Bayesian blockdecomposition from Li et al. (2020) is shown by red lines. The o ff -peak intervals ( φ = − − emproc for MOS data were used to process the raw observationdata files (ODFs). XMM-Newton data were also filtered to avoid the periods of hard X-raybackground flares.The
Chandra data were reduced using CIAO version 4.7 and CALDB version 4.7.7. Wereprocessed the
Chandra data to level = Energy (keV) −4 −3 n o r m a li z e d c o u n t s s − k e V − −101 ( d a t a − m o d e l ) / e rr o r −5 −4 −3 n o r m a li z e d c o u n t s s − k e V − ( d a t a − m o d e l ) / e rr o r Energy (keV) Fig. 8.— PSR J1907 + XMM-Newton
MOS 1 (black) and MOS 2 (red) spectra (left) and