Finding Candidate TeV Halos among Very-high-energy Sources

We search for possible pulsar TeV halos among the very-high-energy (VHE) sources reported in different VHE surveys, among which in particular we use the results from the first Large High Altitude Air Shower Observatory catalog of γ-ray sources. Six candidates are found. They share similar properties of containing a middle-aged, γ-ray–bright pulsar in their positional error circles (the respective pulsars are J0248+6021, J0359+5414, J0622+3749, J0633+0632, J2006+3102, and J2238+5903), being in a rather clean field without any common Galactic VHE-emitting supernova remnants or (bright) pulsar wind nebulae (PWNe), and showing an absence of any γ-ray emissions in 0.1–500 GeV after removing the pulsars’ emissions. Combining these candidates with several reported (candidate) TeV halos, we obtain relationships between their luminosity at 50 TeV (L 50TeV) and the corresponding pulsars’ spin-down energy ( Ė ), which are L50TeV∼Ė0.9 and L50TeV/Ė∼6.4×10−4 , respectively. The relationships are nearly identical to previously reported ones. We probe possible connections between the extent/sizes of the VHE sources and the pulsars’ ages, and find a weak older-and-smaller trend. By comparing to the VHE detection results for PWNe, it is clear that the (candidate) TeV halos have hard emissions by either having power-law indices smaller than 2 in 1–25 TeV or by only being detected in 25–100 TeV. In addition, we also consider seven other VHE sources as possible TeV halos based on the results from different studies of them, but they do not fit cleanly with the properties listed above, indicating their potentially complex nature.


INTRODUCTION
TeV halos are extended Very-High-Energy (VHE; ≥100 GeV) γ-ray emissions around middle-aged pulsars (∼100 kyr).Their existence has been firmly established due to the detection of extended emissions around the nearby pulsars Geminga and Monogem, as observed by the High-Altitude Water Cherenkov (HAWC) Observatory (Abeysekara et al. 2017).Since the detection revelation, various studies focusing on their general existence and possible properties have been carried out (see, e.g., Linden et al. 2017, Sudoh et al. 2019, Fang 2022, and Mukhopadhyay & Linden 2022 and references therein).Importantly, TeV halos can be a significant contributor of cosmic electrons and positrons in our Galaxy (Giacinti et al. 2020;López-Coto et al. 2022a;Yan et al. 2024).
In our recent studies of Galactic VHE sources for those whose nature is not clear (i.e., unidentified), we have focused on finding and studying their possible lower energy counterparts by analyzing available multi-energy data (e.g., Xing et al. 2022;Zheng et al. 2023b).The primary ones used are the GeV γ-ray data (in energy range of 0.1-500 GeV) obtained with the Large Area Telescope (LAT) onboard the Fermi Gamma-ray Space Telescope (Fermi).It has been ascertained that for a significant fraction of VHE sources, a known pulsar is often found located in the field, within the error circle of such a VHE source (e.g., Albert et al. 2020).Some of these pulsars with the positional coincidence are γ-raybright, which can cause difficulties in analyses because of the low-spatial resolution of the LAT data (e.g., ∼1 deg at 1 GeV).The strategy we have applied to overcome such difficulties is to remove the 'contamination' of a pulsar, the pulsed emission, by timing this pulsar at γrays.This helps reveal the residual emission in a source Note-References for X-ray fluxes and distances: (1) Marelli et al. (2011), (2) Theureau et al. (2011), (3) Zyuzin et al. (2018), (4) Prinz & Becker (2015), (5) Danilenko et al. (2020), (7) Nice et al. (2013), (7) Rigoselli et al. (2022), (8) Kargaltsev et al. (2008), (9) Klingler et al. (2020), (10) Marelli et al. (2014), (11) Karpova et al. (2019), ( 12) Pavlov et al. (2008), (13) Kargaltsev et al. (2012) field, which allows us to conduct clean studies of it as a candidate counterpart (e.g., Xing et al. 2022).However in the studies of the VHE sources 3HWC J0631+107 (Albert et al. 2020; or 1LHAASO J0631+1040), 1LHAASO J1959+2846u, and 1LHAASO J2028+3352, where 1LHAASO stands for the first Large High Altitude Air Shower Observatory (LHAASO; Cao et al. 2019) catalog of Gamma-ray sources (Cao et al. 2023), no significant residual emissions were found after we removed the pulsed emissions of PSRs J0631+1036, J1958+2846, and J2028+3332, respectively.The non-detections, combined with the pulsars' many similarities to Geminga and the fact that no primary Galactic VHE-emitting sources (e.g., H. E. S. S. Collaboration et al. 2018a), such as supernova remnants (SNRs) or pulsar wind nebulae (PWNe), are known in the fields, led to our identification of the three VHE sources as being TeV halos powered by their respective pulsars (Zheng et al. 2023a;Zheng & Wang 2023).
Our series of work and obtained results suggest that there could be more TeV halos among the unidentified VHE sources.Particularly for those reported in 1LHAASO, the energy range spreads from 1 TeV to ap- ).In addition, the Third Fermi Large Area Telescope Catalog of Gamma-ray Pulsars (3PC) has very recently been released (Smith et al. 2023).It provides the timing solutions for Fermi-LATdetected γ-ray pulsars.These solutions allow us to easily carry out our studies of VHE sources when pulsed γ-ray emissions need to be removed.We have thus conducted a further search for candidate TeV halos.We have found six good candidates (first six in Table 1) and report the results in this paper.
In this work, we mainly used the detection results in 1LHAASO, but also included results reported in the High Energy Spectroscopy System (HESS) Galactic plane survey (HGPS; H. E. S. S. Collaboration et al. 2018a) and in the third HAWC catalog (3HWC; Albert et al. 2020).In some cases, reported results from the observations conducted with the Very Energetic Radiation Imaging Telescope Array System (VERITAS) and the Milagro Gamma-ray Observatory (MGRO; Abdo et al. 2009b) were also used.To be as complete as possible, we essentially went through all the Galactic VHE sources that are likely associated with a pulsar and show some aspects of a TeV halo.As a result, we found another seven sources and list them in the lower part of Table 1.Some of these sources are in a complex region, such as being potentially associated with an SNR/PWN in the field, and some contain a pulsar that does not have γ-ray emission or clear off-pulse phases in the case of being γ-ray bright.We included these sources in our discussion (Section 4), and a brief introduction for each of them is provided in Appendix A.
In the following Section 2, we describe the Fermi-LAT data we used and our data analyses, which include how we obtained the off-pulse data of the six pulsar targets through pulsar timing.In Section 3, in conjunction with our analysis results, we provide the properties of each pulsar target and its associated VHE source, which helps identify the latter as a candidate TeV halo.In Section 4, we discuss these sources' general properties by considering the VHE sources as being TeV halos powered by the corresponding pulsars.

LAT Data and Source Model
Photon data files with timing analysis results for each of the Fermi-LAT-detected pulsars are provided in 3PC.There are two types of data with different sizes, one containing photons within 3 • with an energy band ranging from 50 MeV to 300 GeV 1 , and the other containing photons of energies from 20 MeV to 1 TeV within 15 •2 .Both types of data files are centered at the position of each pulsar.The data were selected from the latest Fermi Pass 8 database with Event Class of 128.Events with a zenith angle larger than 105 • and bad quality flags were excluded.We used both of the data files in the following analyses.
In our analyses, we chose photons in the energy band of 0.1-500 GeV with a zenith angle of < 90  different and are given in Table 2 for each of the pulsar targets.We set the regions of interests (RoIs) with a size of 15 • × 15 • , centered at each of the pulsar targets.The latest Fermi-LAT Fourth Source Catalog (4FGL-DR4, Ballet et al. 2023) was used to construct source models.
For each target, the sources within a range of 15 • radius were included, and their 4FGL-DR4 spectral forms were used.In addition, two background models, the Galactic and extragalactic diffuse emission, were included in the source models, which were the files gll iem v07.fits and iso P8R3 SOURCE V3 v1.txt, respectively.

Timing Analysis
In both the 3 • and 15 • photon data files, a photon's probability (to be from a pulsar) and spin phase (of the pulsar) were included.We used the 3 • photon files to construct the pulse profiles and define the on-and offpulse phase ranges.The pulse profiles, with the on-and off-pulse phase ranges marked, are shown in Figure 1.The values of the phase ranges, as well as the timing solutions, are given in Table 2. 2.3.Likelihood Analysis of the on-and off-pulse Data 2.3.1.On-pulse Data The standard binned likelihood analysis was performed to the on-pulse data of each pulsar in 0.1-500 GeV.The spectral parameters of sources in a source model within 5 • from a pulsar target were set as free, while those of the other sources were fixed at the values given in 4FGL-DR4.In addition, the ) , was used for PSRs J0359+5414, J0633+0632, J2006+3102 and J2238+5903 (hereafter J0359, J0633, J2006, and J2238, respectively), and the other, ) b ]}, was used for PSRs J0248+6021 and J0622+3749 (hereafter J0248 and J0622, respectively).In the two forms, Γ and d (or ExpfactorS) are the photon index and the local curvature at energy E 0 , respectively, and b is a measure of the shape of the exponential cutoff.Following 4FGL-DR4, we fixed the value of b at 2/3 for our analysis.The likelihood analysis results for each pulsar are given in Table 3.
We also obtained the spectral data points of the pulsars from their on-pulse data.The energy range from 0.1 to 500 GeV was evenly divided logarithmically into 10 bins.Binned likelihood analysis was performed to each bin's data to obtain the fluxes.In this analysis, the normalizations of the sources within 5 • of a pulsar and the two background components were set as free parameters, while the other parameters were fixed at the values obtained above from the binned likelihood analysis of the data in the whole energy range.When the test statistic (TS) value of a bin was <4, we replaced the flux with the 95% upper limit as derived from the data of the bin.The obtained spectra are shown in Figure 2.

Off-pulse Data
We also performed standard binned likelihood analysis to the off-pulse data of the pulsars.The parameter setup was the same as that in the above analysis of the on-pulse data (Section 2.3.1).We assumed a power law (PL) for any emission at the position of each pulsar, From the analysis, no significant emissions were detected during the off-pulse phase ranges of the pulsars.In Table 3, we provided the TS values when we assumed Γ = 2 as the exemplary result.
To show the non-detection results in the off-pulse data and provide a clear view of the source fields, we obtained 0.1-500 GeV TS maps for the pulsars' regions.As indicated by the TS maps (Figure 3), no significant residual emissions are seen in any of the pulsar regions.

PULSARS AND THEIR ASSOCIATED VHE SOURCES
Below, based on the analysis results we obtained for each pulsar target (cf., Figures 2 & 3), we briefly describe the properties of the pulsars and their likely associated VHE sources in the following sections.We also searched for X-ray observational results for the pulsars, which helps us learn about the properties of their PWNe.When needed, we analyzed the archival X-ray data by ourselves.

PSR J0248+6021
PSR J0248 is a γ-ray pulsar, with its radio pulsations first detected by the Nançay radio telescope (Foster et al. 1997).At the pulsar's position, no X-ray emission was detected, and no evidence showed any ex-  Figure 2. Spectra and spectral upper limits of the six pulsars during their on-and off-pulse phase ranges, which are shown as black data points (and black dashed curves, the best-fit PLSEC models) and red lines (assuming Γ = 2), respectively.In addition, we overplot the spectra of the LHAASO sources, and available HAWC flux measurements and/or VERITAS (Archer et al. 2019) and MGRO (Abdo et al. 2009b) flux upper limits of the VHE sources.For details, see Section 3 and Figure 3.
LHAASO detected an extended source, 1LHAASO J0249+6022, that is in positional coincidence with J0248.The extension of the source is ∼ 0. • 38.In this region, no excess γ-ray emission was detected in off-pulse phase data analysis (see Figure 3).We noted that the region is rather clean, within which no SNRs are listed in the SNR catalog SNRcat 3 .Therefore we suggest that the extended TeV emission, 1LHAASO J0249+6022, is a TeV halo candidate powered by PSR J0248.

PSR J0359+5414
The region of this pulsar is clean with no residual emissions detected in the off-pulse data (Figure 3).In X-rays, a weak PWN was detected with a luminosity 3 http://snrcat.physics.umanitoba.ca of ≃ 2.8 × 10 31 erg s −1 at a pseudo distance of 3.45 kpc (Zyuzin et al. 2018).Extended TeV emission at the region was reported by HAWC, with a size of 0. • 2±0.• 1, and this detection was matched by the LHAASO detection results.
The VHE source was already posited to be a TeV halo candidate powered by J0359 (Albert et al. 2023b).To distinguish TeV halos and PWNe, it has been discussed that the VHE γ-ray emissions are from a larger region than those of PWNe (Linden et al. 2017;López-Coto et al. 2022b;Albert et al. 2023b).However, this case is complicated by the existence of a nearby radio pulsar B0355+54 (Figure 3), which has a spin-down luminosity of Ė = 4.5 × 10 34 erg s −1 and a characteristic age of τ = 564 kyr, and as such, this radio pulsar's possible association with the VHE source could not be excluded (Albert et al. 2023b).
In the region, LHAASO detected extended VHE γ-ray emission named LHAASO J0621+3755, and it is likely a TeV halo (Aharonian et al. 2021).In the 1LHAASO catalog, 1LHAASO J0622+3754 was assigned to be associated with LHAASO J0621+3755, with the separation between them being only 0. • 03.Our analysis of the offpulse data verified the emptiness of the field at GeV γ-rays.
The distance of the pulsar was estimated to be 1.6 kpc by Pletsch et al. (2012), where the pulsar's γray luminosity L γ in 0.1-100 GeV was estimated from a L γ -Ė relationship that was derived based on the γ-ray pulsars with distance measures (for details see Saz Parkinson et al. 2010;Pletsch et al. 2012).We reestimated the distance by using the flux value given in the recent 4FGL-DR4, and found a value of 1.4 kpc.However this value can be highly different from the actual one.Another method to estimate the distance is to require L γ ≤ Ė, which sets an upper limit of 3.47 kpc for the distance, and if considering L γ ∼ 0.1 Ė, the distance would be ∼1.1 kpc.We adopted 1.1 kpc for J0622 but with an upper limit of 3.47 kpc.

PSR J0633+0632
This pulsar is also radio quiet, first detected at γ-rays by Fermi-LAT (Abdo et al. 2009a).In this source region, diffuse X-ray emission was detected and identified as a PWN (Ray et al. 2011;Danilenko et al. 2020).
The LHAASO detection indicated that the VHE source has a hard emission, as the WCDA observations only provided a flux upper limit (Cao et al. 2023).Our analysis of the off-pulse data provided a flux upper limit of ∼ 10 −13 erg cm −2 s −1 in the GeV energies.Considering that the source 1LHAASO J0635+0619 (as well as 3HWC J0634+067) is a TeV halo candidate powered by PSR J0633, the VHE emission is likely from a larger region than that of the X-ray PWN.Khokhriakova et al. (2023) searched for X-ray counterparts of TeV halos (socalled X-ray halos) around 5 pulsars that include J0633.However, no such extended emission was found.

PSR J2006+3102
This radio pulsar was reported with a distance of 4.7 kpc in Nice et al. (2013), but the updated value is 6.035 kpc in the Australia Telescope National Facility (ATNF) pulsar catalog (Manchester et al. 2005).Very limited information is available for this pulsar.We searched the Chandra and XMM-Newton archival data, but no observations were found.Using a set of data (Obsid : 03103085001, exposure time = 1.1 ks) obtained with the X-ray Telescope (XRT) onboard the Neil Gehrels Swift Observatory (Swift), we derived a 3σ upper limit of 0.01 cts s −1 in 0.3-10 keV at the pulsar's position.The corresponding energy-flux upper limit was 9.0 × 10 −13 erg cm −2 s −1 , where we assumed a PL source spectrum with an index of 2 and hydrogen column density N H = 8.27 × 10 21 cm −2 (towards the source direction, from HI4PI Collaboration et al. 2016).
Close to the edge of the extension region given by the LHAASO KM2A, there is an SNR, G68.6−1.2 (Figure 3), which, however, is faint and poorly defined according to the SNRcat.Given its poorly known properties and relatively large separation (∼ 0. • 68) from the VHE source, it is not clear if the SNR can be connected to the latter.We noted that 3HWC J2005+311 is also located in this region (Albert et al. 2020), and its spectrum is similar to that of 1LHAASO J2005+3050.However, the positions of the two sources do not overlap.The relation between them remains to be resolved from further observational results.

J2238+5903
J2238 also has very limited information available.An X-ray flux upper limit was reported by Prinz & Becker (2015) to be 4.4×10 −14 erg cm −2 s −1 in 0.1-2 keV (where a PL with index = 1.7 was assumed).The LHAASO WCDA observations were influenced by the Galactic Diffuse Emission (Cao et al. 2023), and we did not consider the WCDA measurements.

DISCUSSION
Following our previous work on identifying candidate pulsar TeV halos (Zheng et al. 2023a;Zheng & Wang 2023) by mainly analyzing the off-pulse GeV data of γray pulsars in the fields of VHE sources, from which any residual emissions may help reveal their nature as possibly being primary Galactic sources, such as SNRs or PWNe, we further found six candidates because of the non-detection of any significant residual emissions.The pulsars' properties, including information for their X-ray emissions, are summarized in Table 1.As discussed in Zheng & Wang (2023), there may contain a relationship between the TeV halos' luminosity at 50 TeV, L 50TeV , and corresponding pulsars' spin-down energy Ė.This relationship helps indicate the fraction of the total energy spent on powering the TeV halos.Since most of the sources (including those presented in the Appendix) in this work have been detected by LHAASO KM2A in 25-100 TeV, we thus also estimated their L 50TeV from the differential fluxes at 50 TeV given in the LHAASO results (Cao et al. 2023).The L 50TeV values are provided in Table 1.Fitting the data points that include four sources in Zheng & Wang (2023) and five sources in this work (excluding J0622 whose distance is highly uncertain), we obtained L 50TeV = 2.27 +1.82 −1.72 Ė0.90 +0.02 −0.01 , with a reduced χ 2 value of ≃0.8 for 7 degree of freedom (DoF), where we assumed a 30% uncertainty for distances (this uncertainty was dominant).We used the Markov Chain Monte Carlo (MCMC) code emcee (Foreman-Mackey et al. 2013) for the fitting, since it conveniently provides error ranges.It can be noted that the L 50TeV ∼ Ė0.9 relationship (see Figure 4), reported in Zheng & Wang (2023), still holds.
Another relationship we tested was L 50TeV / Ė being either a function of the pulsars' characteristic ages τ or a constant.Fitting the data points, we obtained kyr (where τ is in units of kyr) with reduced χ 2 ≃0.8 for DoF=7, or L 50TeV / Ė = 6.4 ± 0.8 × 10 −4 with reduced χ 2 ≃ 0.6 for DoF=8.Both results are also very similar to those previously obtained in Zheng & Wang (2023).For the first result, the large uncertainty for the τ index indicates its value close to zero, and thus the second result, L 50TeV / Ė being a constant (as in Zheng & Wang 2023), is preferred.
In addition, we also tested the physical sizes S of the VHE sources as a function of τ .The sizes were derived from the extension sizes in degrees, as summarized in Table 1 from the LHAASO KM2A measurements.We obtained S = 64.51+21.54  −21.04 × τ −0.25 +0.09 −0.07 kyr pc, with reduced χ 2 ≃ 2.1 for DoF=6.The uncertainties are large, and there is a source, J2028+3332 (Zheng & Wang 2023), significantly deviating from the relationship (although the source's distance is uncertain).In any case, there is a possible older-and-smaller trend, which could be an interesting feature that may reveal the evolutionary processes of electron/positron ejection of pulsars and halos.Further observational results obtained from more collected data with LHAASO may verify this trend.
As we also searched for other potential TeV halo candidates from among mainly 1LHAASO sources, there are seven of them whose properties may provide hints as to their possible TeV-halo nature based on different studies (See Appendix Section A).We showed their corresponding properties (Table 1) in Figure 4.As can be seen, they generally have large scattering around the relationships we obtained above.In particular, five of them were compact sources (see Appendix Section A and Figure A2) in the LHAASO KM2A measurements.Because different VHE observational facilities have different sensitive energy bands and spatial resolutions, we did not try to replace the KM2A results with those from other facilities.Thus, most of these sources currently do not fit in the S ∼ τ relationship at all.From this comparison, we may conclude that either they are not TeV halos or their emissions may contain significant contributions from other sources, which would be in agreement with the various results from many multi-energy studies about them (Appendix Section A).Finally, H. E. S. S. Collaboration et al. (2018b) studied all PWNe and candidates in 1-10 TeV.On the basis of their results, one conclusion that may be drawn is PWNe tend to have a soft spectrum with the PL index Γ s > 2. By comparison, as pointed by Zheng & Wang (2023), candidate TeV halos often show hard spectra with Γ s < 2. We further explored this possible feature by constructing Figure 5, in which the PL indices of the candidate TeV halos (as well as the sources described in the Appendix) and the HESS-confirmed PWNe are shown, where the hard PL indices Γ h are from the LHAASO KM2A 25-100 TeV measurements.Some of the sources, in particular those HESS PWNe, were only detected in one energy band (such as the soft 1-10 TeV band), and we put these sources at the Γ h = Γ s line; note that the error bars indicate the measurements at which energy band are known.It is clear to see that most candidate TeV halos either show emissions with Γ s < 2, or simply have detectable hard TeV emissions (only with known Γ h in 25-100 TeV).By comparison, PWNe have soft emissions with Γ s > 2 or do not have any detectable hard TeV emissions (those at the Γ h = Γ s line).One notable source of the PWNe is that of the Vela pulsar, which has Γ s < 2 (the data point at the low left along the Γ h = Γ s line in Figure 5).On the other hand, one exception among the candidate TeV halos is 1LHAASO J0249+6022 (associated with PSR J0248), which has Γ s > 2. Detailed studies of this source may help understand the cause of the deviation.In any case, the comparison strengthens our previous suggestion in Zheng & Wang (2023) that TeV halos are different from PWNe by having hard emissions.
We thank the anonymous referee for helpful comments.This research is supported by the Basic Research Program of Yunnan Province (No. 202201AS070005), the National Natural Science Foundation of China (12273033), and the Original Innovation Program of the Chinese Academy of Sciences (E085021002).D.Z. acknowledges the support of the science research program for graduate students of Yunnan University (KC-23234629).We also collected information for seven additional VHE sources, given their particular features revealed from observational and theoretical studies.Three of the corresponding pulsars have γ-ray emissions, but two of them (J1740+1000 and J1813−1246) do not have clear off-pulse phase ranges.The timing solutions we used are given in Table A1, and the pulse profiles are shown in Figure A1.The other four pulsars do not have detectable γ-ray emissions in Fermi-LAT data.We analyzed the Fermi-LAT data using the same processes described in Section 2, and calculated the TS maps for each of the sources.For the field of PSR J1826−1256, the TS map was from the pulsar's off-pulse data, and for all the other sources, the TS maps were from the whole data but with all known sources given in 4FGL-DR4 removed.As the information in the following sections indicates, most of the sources are in complex regions (see Figure A2) and their nature is still under different investigations.

A.1. PSR B0540+23
This pulsar was first detected by Jodrell Bank Mark IA radio telescope (Davies et al. 1972).Analyzing the XMM-Newton archival data, Prinz & Becker (2015) reported a faint X-ray point source associated with the radio position.No X-ray diffuse emission was detected by the extended ROentgen Survey with an Imaging Telescope Array (Khokhriakova et al. 2023).HAWC detected ∼ 0. • 5 extended TeV emission, named HAWC J0543+233, around the pulsar, which was suggested to be a potential TeV halo (Riviere et al. 2017).In 3HWC, the source was resolved to two sources, 3HWC J0540+228 and J0543+231 (Albert et al. 2020).LHAASO detected a source, 1LHAASO J0542+2311u, positionally coincident with 3HWC J0543+231 in the energy range of from 25 TeV to above 100 TeV, and the detection showed a large extension of ∼ 0. • 98.The pulsar is away from the reported positions of 1LHAASO J0542+2311u and 3HWC J0543+231 by ∼ 0. • 29 and ∼ 0. • 37, respectively.We analyzed the Fermi-LAT data, but no GeV emissions were detected at the pulsar's position.In the extension region given by LHAASO, only one Fermi-LAT source, 4FGL J0544.4+2238, is known but at the region's edge.This GeV source is an unidentified source with soft PL emission (photon index was 3.33±0.17).It is rather hard to make connections between the soft emission with the VHE sources.Thus in the rather clean region, we considered the VHE emission possibly originates from a candidate TeV halo.

A.2. PSR J1740+1000
The Arecibo Telescope observations discovered this young pulsar, but it is located away from the Galactic Plane (McLaughlin et al. 2002).An X-ray tail structure was found behind the pulsar (Kargaltsev et al. 2008).VERITAS searched for TeV γ-ray emission from the tail, but no emissions were detected (Benbow et al. 2021).HWAC detected a source in 14.8-274.0TeV in the region, which is likely also the LHAASO source reportedly detected in the energy range from 25 TeV to above 100 TeV.The pulsar is listed as a Fermi-LAT source but with a detection significance of only 6.7σ, and its pulsed γ-ray emission is likely too faint to be clearly identified (Figure A1).We tested to remove the GeV source and checked if there were any residual emissions (the Fermi-LAT data analyzed were in 0.1-500 GeV during the time period of from 2008-08-04 15:43:36 (UTC) to 2023-02-16 00:00:00 (UTC)).No such emission was found.Given the non-detection of the PWN tail in VHE energies, we considered the HAWC/LHASSO VHE source is a candidate TeV halo possibly powered by the pulsar.It should be noted that if this scenario is the case, it suggests that TeV halos could emit ultra-highenergy (UHE, > 100 TeV) γ-rays, since the source had a TS value of 37.2 above 100 TeV in the LHAASO's detection.

A.3. PSR J1809−1917
The region is rather complex because of different VHE emission detections plus additional sources revealed by related multi-energy studies.A VHE source, HESS J1809−193, was first reported in the region and was considered as the emission originating from a candidate PWN associated with the pulsar (Aharonian et al. 2007).Indeed, the PWN-like extended emission and its variations were observed in X-rays (Kargaltsev & Pavlov 2007;Anada et al. 2010;Klingler et al. 2018Klingler et al. , 2020)).In addition, there are at least three SNRs, G11.0−0.0,G11.1+0.1, and G11.4−0.1, and several molecular clouds (MCs) located in this region (Castelletti et al. 2016).The interaction between G11.0−0.0 and a nearby MC was proposed to be the process producing the VHE emission (Castelletti et al. 2016).
HGPS reported the source with a size of ∼ 0. Correspondingly, HAWC and LHAASO each detected a source in the region, while detailed properties, such as the extension size and flux, are different.
Given the complexness of the region, any clear identification for the connections between different sources at multi-energies is not straightforward.
In any case, given the work in H. E. S. S. Collaboration et al. ( 2023), we considered the extended component (matches the LHAASO WCDA detection) as being a possible TeV halo.
No obvious off-pulse phase range could be determined in the GeV energy band (Figure A1), and no PWN was found in X-ray observations (Marelli et al. 2014).VHE emissions were detected by HESS, HAWC, and LHAASO (H.E. S. S. Collaboration et al. 2018a;Albert et al. 2020;Cao et al. 2023) with a high positional coincidence.In H. E. S. S. Collaboration et al. (2018a), the VHE source was suggested to be associated with a relic PWN that is only detected in the TeV band.With limited information for the sources and the region, if we consider a relic PWN, an extended halo forming from escaping particles from the PWN could be a possible scenario.We thus included this source in the studies.
We were able to define an off-pulse phase range for the γ-ray emission of the pulsar (Figure A1), but no significant residual emission was detected in the off-pulse data.In 3PC, the pulsar's distance was estimated to be 1.55 kpc.However, Karpova et al. (2019) used interstellar reddening relationships towards the source direction and estimated a value of ∼3.5 kpc.We included this pulsar in our studies, although it should be noted that its characteristic age is rather young, 14.4 kyr, and its PWN was clearly detected in X-rays (Karpova et al. 2019).

A.6. PSR J1826−1334
This pulsar is relatively young, with a characteristic age of ∼ 21 kyr.HESS detected largely extended VHE emission with a size of ∼ 1 deg around this pulsar.Together with 4FGL J1824.4−1350e, the VHE emission has been identified as arising from the PWN of the pulsar (Aharonian et al. 2006;Khangulyan et al. 2018;Duvidovich et al. 2019;H. E. S. S. Collaboration et al. 2020).However, from the theoretical point of view, the large extension is not easily explained with the typical PWN modeling (see details in Khangulyan et al. 2018 andCollins et al. 2024).We included this source as a potential TeV halo case.

A.7. PSR J1928+1746
HAWC detected VHE emission associated with this pulsar (Albert et al. 2020), named 3HWC J1928+178.This radio pulsar (Cordes et al. 2006) was not found with an X-ray counterpart in Chandra and NuSTAR observations (Kargaltsev et al. 2012;Mori et al. 2020).Albert et al. (2023a) re-analyzed the latest HAWC data, and two components were revealed in the region of 3HWC J1928+178, one with a size of ∼ 0. • 18 and another with an extremely large size of ∼ 1. • 43.We ana-lyzed Fermi-LAT data, but no significant GeV emission was detected at the position of the pulsar or in the surrounding region (Figure A2).Considering non-detection of X-ray emission, the age of the pulsar, the extendedness of 3HWC J1928+178, and the low energy density compared to the local interstellar medium, the HAWC source was suggested to be in a transitional phase from a PWN to a TeV halo (Albert et al. 2023a).
LHAASO also detected a source close to PSR J1928+1746, 1LHAASO J1928+1746u, but the extension was only found in WCDA observations, not in KM2A ones.In addition, there was a nearby source 1LHAASO J1928+1813u that has a large extension of ∼ 0. • 63, only detected by KM2A (Table 1).Given the high positional coincidence of this source with that of the largely extended component identified by HAWC, we considered that these two could be associated with each other.We included 1LHAASO J1928+1813u instead as a candidate TeV halo that is associated with the pulsar.

Figure 1 .
Figure 1.Pulse profiles (top) and two-dimensional phaseograms (bottom) of six pulsar targets.The on-and off-pulse phase ranges we defined based on the pulse profiles are marked by the dashed lines.

Figure 3 .
Figure3.TS maps of the regions of the six pulsar targets in 0.1-500 GeV, calculated from the off-pulse data of the pulsars.Each panel has a size of 3 • × 3 • centered at a pulsar target.Green diamonds and crosses mark the positions of the pulsars and nearby Fermi -LAT sources, respectively.The positional error circles and extension regions of LHAASO, HAWC, and HESS sources are marked by solid and dash circles respectively, and the name of the corresponding 1LHAASO source is given at top-left of each panel.In the region of 1LHAASO J0359+5406, a non-γ-ray pulsar is indicated by a yellow diamond, and in that of 1LHAASO J2005+3050, an SNR is shown as the magenta dashed circle.

Figure 5 .
Figure 5. Power law indices of the (candidate) TeV halos and HESS-confirmed PWNe.Values of the first group are from LHAASO WCDA (Γs in 1-25 TeV) and KM2A (Γ h in 25-100 TeV), and those of the latter are mostly from HESS in 1-10 TeV (H.E. S. S. Collaboration et al. 2018b).When there is only one measurement, either Γs or Γ h , the source is put at the Γ h = Γs line (dash-dotted).Because two WCDA measurements suffered GDE, the sources are also put at the Γ h = Γs line as well, indicated by the dotted lines.HESS 1825−137 (or 1LHAASO J1825−1337u, associated with PSR J1826−1334) and CTA 1(Aliu et al. 2013; or  1LHAASO J0007+7303u)  have different reported Γs values, and both values of each of them are shown (connected with a grey line and pointed with an arrow).The PWN of the Vela pulsar is the lowest grey data point along the Γ h = Γs line (with Γs < 2).
are obtained from 3PC APPENDIX A. TIMING ANALYSIS AND BRIEF INTRODUCTION FOR SEVEN CANDIDATES

Figure A1 .
Figure A1.Pulse profiles (top) and two-dimensional phaseograms (bottom) of three γ-ray pulsars.Only for PSR J1826−1256, an off-pulse phase range can be defined, which are marked by dashed lines.
• 4 (H.E. S. S. Collaboration et al. 2018a), but H. E. S. S. Collaboration et al. (2023) re-analyzed the data and revealed 2 components, an extended plus a compact one.While the compact component could be either the VHE emission of the PWN or due to the SNR-MC interaction, the extended component was suggested to be a halo around the PWN of the pulsar (H.E. S. S. Collaboration et al. 2023).

Figure A2 .
Figure A2.TS maps of a size of 3 • × 3 • in 0.1-500 GeV for the seven VHE sources and corresponding pulsars.Green diamonds indicate the positions of γ-ray pulsars and yellow ones the positions of the other non-γ-ray pulsars.The positional error circles and extension region of LHAASO, HAWC, and HESS sources are marked with solid and dash circles respectively.All the catalog Fermi -LAT sources are also indicated with their catalog names.

Table 1 .
Properties of pulsars and likely associated VHE sources.

Table 2 .
Timing solutions and phase ranges for six pulsar targets.
Zheng et al. 2023a andwide eneWang 2023 coverage allows us to select sources with spectra harder than those of the PWNe (H.E.S. S. Collaboration et al. 2018b), which is a possible feature that may be used to differentiate the TeV halos from the PWNe (seeZheng et al. 2023a and Zheng &Wang 2023

Table 3 .
Binned likelihood analysis results from the on-and off-pulse phase data.