EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 668 (December 2022) A&A, 668 (2022) A23 Full HTML
Open Access
Issue
A&A
Volume 668, December 2022
Article Number A23
Number of page(s) 9
Section Astrophysical processes
DOI https://doi.org/10.1051/0004-6361/202244841
Published online 29 November 2022

© De Sarkar et al. 2022

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe-to-Open model. Subscribe to A&A to support open access publication.

1. Introduction

Recent observations by state-of-the-art observatories, such as the Large High Altitude Air Shower Observatory (LHAASO), Tibet ASγ, the High Altitude Water Cherenkov (HAWC), among others, have paved the way for the detection of multiple Galactic ultra-high energy (UHE; Eγ ≥ 100 TeV) gamma-ray sources (Abeysekara et al. 2020; Amenomori et al. 2019, 2021; Cao et al. 2021). Upcoming observatories such as the Cherenkov Telescope Array (CTA; Cherenkov Telescope Array Consortium 2019) and the Southern Wide-field Gamma-ray Observatory (SWGO; Albert et al. 2019) will be of importance to identify and characterize these PeVatrons: Galactic CR sources that accelerate particles up to PeV energies.

LHAASO is a state-of-the-art dual-task facility designed for CRs and gamma-ray studies at a few hundred GeV to a few PeV, located at 4410 m above sea level in China (Cao 2010). The recent data reported by the LHAASO observatory show the existence of 12 significantly detected sources (> 7σ) that emit gamma rays with energies above several hundred TeVs (Cao et al. 2021). Most of the sources reported by LHAASO have significantly extended gamma-ray emission regions up to ∼1°. The very high energy (VHE; 100 GeV ≤ Eγ ≤ 100 TeV) counterparts of these sources residing in the Galactic plane have been associated with pulsar wind nebulae (PWNe), based on the spatial proximity with highly energetic pulsars and typically extended morphological features (H. E. S. S. Collaboration 2018). It has already been posited that UHE gamma-ray emission spatially coincident or in very close proximity of energetic pulsars with a high spin-down luminosity (̇E >  1036 erg s−1) may be a universal feature; for more details, readers can refer to Albert et al. (2021), for example. Moreover, the Crab nebula, associated with pulsar PSR B0531+21, was confirmed to be a PeVatron by recent LHAASO observations (Cao et al. 2021). Bearing all of this in mind, it is natural to consider PWNe as possible Galactic PeVatrons, from which UHE gamma rays are detected.

PWNe, considered to be one of the most efficient lepton accelerators in the Galaxy, are powered by highly energetic pulsars. Pulsars dissipate most of their rotational energy via the injection of ultrarelativistic electron-positron pairs, which form a cold, ultrarelativistic wind of particles. Since the bulk velocity of this ultrarelativistic wind is supersonic with respect to the ambient medium, this wind creates a termination shock. Injected particles can be accelerated to very high energies at this termination shock. The accelerated leptons can then interact with the ambient matter, photon fields, and the magnetic field through Bremsstrahlung, inverse-Compton (IC), and synchrotron processes. The cooling of the accelerated leptons results in a MWL spectrum ranging from radio to gamma-ray energies. Here, we consider LHAASO J2226+6057 as our source of interest, as this is the only UHE gamma-ray source for which data have been observed across radio, X-ray, and gamma-ray energy ranges. It has been observed that LHAASO J2226+6057 is situated in a complex morphological region. Due to the close spatial proximity of Boomerang PWN, as well as supernova remnant (SNR) G106.3+2.7 and associated molecular clouds (MCs), it is hard to confirm the exact source responsible for the UHE gamma-ray emission observed. From the observations of Imaging Air Cherenkov Telescopes (IACTs), such as MAGIC and VERITAS, it was found that the emission region was divided into two morphological regions: the head and the tail. The head region contains Boomerang PWN and PSR J2229+6114, and the tail region contains VER J2227+608, which is likely to be associated with SNR G106.3+2.7 and the MCs in the region. Faint and diffuse radio and X-ray emission were also observed from the tail region. It is possible that the emission from Boomerang PWN illuminates the head region, whereas hadronic interaction occurring between the SNR and the MCs are responsible for the gamma-ray emission in the tail. Although it is difficult to confirm which source is actually responsible for the UHE gamma-ray emission, comprehensive studies are needed to explore both scenarios. New data from upcoming observations with high angular resolution by MAGIC and VERITAS will be crucial to shed new light on this source. We discuss the features of LHAASO J2226+6057, as well as the previous works done for this source in Sect. 2. Then we comment on the PWN model used in this work in Sect. 3, and present the results in Sect. 4. We finally conclude in Sect. 5.

2. LHAASO J2226+6057 features

LHAASO J2226+6057 was detected at RA = 336.75° and Dec = 60.95° with a significance of 13.6σ above 100 TeV. Its gamma-ray spectrum reaches up to a maximum energy of 0.57 ± 0.19 PeV (Cao et al. 2021). This source is spatially associated with SNR G106.3+2.7, as well as the pulsar J2229+6114 and its wind nebula, known as the “Boomerang” nebula (Kothes et al. 2001). PSR J2229+6114 is a bright gamma-ray pulsar with a spin period of 51.6 ms, a characteristic age of 10 460 yr, and a spin-down luminosity of 2.2 × 1037 erg s−1 (Halpern et al. 2001a). Pulsed GeV gamma-ray emission from this pulsar was detected by Fermi-LAT (Abdo et al. 2009a). In the VHE range, SNR G106.3+2.7 was observed by VERITAS as VER J2227+608, located 0.4° away from PSR J2229+6114 (Acciari et al. 2009). Gamma rays with GeV energies (Xin et al. 2019), diffuse nonthermal X-rays (Fujita et al. 2021), as well as radio data (Pineault & Joncas 2000) have also been observed from the source region. The distance of the source suffers from great uncertainty. It was estimated to be 7.5 kpc (based on the pulsar dispersion measure; Abdo et al. 2009a), 3 kpc (based on X-ray absorption measurements; Halpern et al. 2001a), and 0.8 kpc (based on measurements of radial velocities of atomic hydrogen and molecular material; Kothes et al. 2001). In this work we consider the distance to be 3 kpc, similar to Joshi et al. (2022) and Yu et al. (2022).

The possible connection between LHAASO J2226+6057 and PSR J2229+6114 has previously been studied by Breuhaus et al. (2022), Joshi et al. (2022), and Yu et al. (2022). Breuhaus et al. (2022) performed a time-independent one-zone treatment, which is a steady-state leptonic scenario from the PWN, to explain only the highest energy gamma-ray data. The evolution of the source was neglected, and also MWL data were not used. Joshi et al. (2022) used a time-dependent one-zone leptonic scenario from the PWN. However, the authors did not consider the effect of escape for LHAASO J2226+6057, which they only accounted for LHAASO J1908+0621, nor did they explore the effect of age and the braking index on the evolution of the injected leptonic population. The impact of the SNR reverse shock and its effects on the PWN radius evolution were not taken into account, assuming that such effects are important only if the age of the PWN is greater than 10 kyr, which is not necessarily the case; for more details, readers can refer to Martín et al. (2016), for example. Yu et al. (2022) also performed a similar study. The authors argue that a distorted nebula, which was created due to the impact of a SNR reverse shock, is responsible for the GeV gamma-ray emission observed from the source region by Xin et al. (2019). Their results are similar to those of Joshi et al. (2022). Our PWN model here intends to test their conclusions after relaxing assumptions or adding additional physical details.

A recent paper by Tibet ASgamma Collaboration (2021) has proposed a hadronic origin of the LHAASO source based on the spatial proximity of a molecular cloud (MC) with the gamma-ray centroid of the source. However, as pointed out by Breuhaus et al. (2022), the associated SNR is quite old to produce the observed hard gamma-ray spectrum at the highest energies. Consequently, a hadronic scenario from a SNR+MC association would need a peculiar modeling to explain the UHE gamma-ray spectrum observed from the source. A novel approach was explored in De Sarkar & Gupta (2022) to explain the hadronic origin of LHAASO J1908+0621 and, perhaps, a similar approach may participate here as well.

3. Brief description of the model

For this work, we have used the code TIDE, for which earlier applications can be found in Torres et al. (2014), Martín et al. (2016), and Torres & Lin (2018), for example. The numerical code solves the evolution of leptonic pair distribution in the PWN, as a function of Lorentz factor γ at time t described by the following equation:

(1)

The left-hand side of Eq. (1) describes the variation of lepton distribution in time. The first term on the right-hand side is the lepton injection function Q(γ, t), which is usually assumed as a broken power law:

(2)

The second and third terms on the right-hand side take radiative losses into account – such as synchrotron, IC, Bremsstrahlung, adiabatic losses, or heating – as well as the escape of the particles (we assume Bohm diffusion), respectively (see Martín et al. 2012 for the incorporated formulae). The normalization factor Q0(t) is calculated using the spin-down luminosity L(t) of the pulsar through the equation, is the fraction of spin down power that goes on to power the magnetic field of the PWN. The magnetic field varies in time as a result of the balance between the magnetic field power and adiabatic gains or losses of the field due to the contraction or the expansion of the PWN (Torres et al. 2012). The reverberation phase of the PWN included in the model has been considered to be the same as that in Torres (2017). The variation of the PWN radius is calculated by taking into account the age, SNR explosion energy, ambient medium density, expansion velocity of the nebula, and pressure profiles of the SNR at the position of the PWN. The model takes into account the change in pressure profiles depending on whether the PWN shell is surrounded by unshocked ejecta (i.e., R <  RRS), or the shocked ejecta (i.e., RRS <  R <  RSNR), where RRS and RSNR are the radii of the SNR reverse shock and the SNR, respectively. After reverberation, an assumed Sedov expansion follows when the PWN reaches the pressure of the SNR. For more details, readers can refer to Martín et al. (2016), Torres (2017), and Torres & Lin (2018), for example. For a more thorough description of the reverberation period, we need to account for the fact that the ejecta pressure is not constant (Bandiera et al. 2020, 2022). This is beyond the scope of this work since, at the moment, a prescription to deal with this fact in the context of the radiative PWN is unavailable. Given the likely young age of the source, a PWN would be before the time of largest compression (see below), and thus we expect the current approach to be acceptable.

4. Results

4.1. Braking index and true age exploration

In the works of Joshi et al. (2022) and Yu et al. (2022), the true age and the braking index were chosen to be 7000 yr and 3, respectively. Since the choice of the true age is essential for a time evolutionary model of PWNe, we explore other values here as well. On the other hand, if the pulsar is assumed to be a dipole in vacuum, and considered to be spinning down by emitting only magnetic dipole radiation (MDR), then the corresponding braking index associated with the spin down is calculated to be 3 (Manchester & Taylor 1977). Alternatively, the spin down of a pulsar driven entirely by a particle wind would result in a braking index of 1 (Michel 1969; Manchester et al. 1985). A combination of magnetic dipole radiation and wind braking would result in a braking index with a value in between 1 and 3 (Archibald et al. 2016). Most of the observed pulsar braking indices fall within this range (Espinoza et al. 2011b, Pons et al. 2012, and references therein). Although, there are exceptions: for example, PSR J1640–4631, which has a braking index of n = 3.15 ± 0.03 (Archibald et al. 2016), and PSR J1734–3333, which has a braking index of n = 0.9 ± 0.2 (Espinoza et al. 2011a). For this work, we only consider the braking index range 1 <  n <  3. Since the choice of the braking index affects the characteristics of the pulsar spin down, it is also important to explore whether the variation of n affects the SED.

To explore the effect of n and tage on the MWL SED, we have chosen tage = 1000, 4000, and 7000 yr, and n = 1.5, 2.0, 2.5, and 3.0. Given the characteristic age τc ≈ 10.5 kyr and the present-day spin-down luminosity L(tage)≈2.25 × 1037 erg s−1, the initial spin-down age (τ0) and initial spin-down luminosity (L0) were calculated using the relations,

(3)

and,

(4)

for a specific choice of n and tage. We also assume the values of the following parameters are the same throughout the work: minimum Lorentz factor γmin = 1; SN explosion energy ESN = 1051 erg; interstellar medium (ISM) ambient density ρISM = 0.1 cm−3; SNR core density index = 0; SNR envelope density index = 9 (we assume a type II SN as the progenitor, as in Chevalier & Fransson 1992; Blondin et al. 2001; Gelfand et al. 2009); PWN adiabatic index = 1.333 and SNR adiabatic index = 1.667; containment factor ϵ (ratio of the Larmor radius of particles to the radius of the termination shock) = 0.5: and magnetic compression ratio κ = 3 (strong shock condition). The ejected mass of the progenitor SNR was also fixed at Mej = 8 M⊙ for the n − tage exploration. For the target radiation fields for the IC interaction, we considered the cosmic microwave background, far-infrared (FIR) and near-infrared (NIR) radiation fields at the source position with temperatures of TFIR = 25 K and TNIR = 5000 K, respectively. The energy densities associated with the FIR and NIR field were taken from Porter et al. (2006) and fixed at those values for n − tage exploration. It has almost always been observed (see Torres et al. 2014) that the energy densities of the FIR and NIR fields required to explain the IC emission from PWNe, differ from the interstellar radiation fields reported by Porter et al. (2006). Breuhaus et al. (2021) also argue that an enhancement of the radiation field is needed to explain the hard gamma-ray spectrum observed at highest energies. However, since the aim of this subsection is to compare the effect of n and tage on the computed SED, we have fixed the energy densities of FIR and NIR fields at the values given by Porter et al. (2006), that is ωFIR = 0.29 eV cm−3 and ωNIR = 0.45 eV cm−3, so there are fewer free parameters left to vary. The maximum Lorentz factor of the lepton distribution was fixed by the most restrictive condition between the synchrotron limit (de Jager et al. 1996) or gyroradius limit (de Jager & Djannati-Ataï 2009) during the evolution. For a specific choice of n and tage, we varied the injection function parameters, that is to say low energy index α1, high energy index α2, energy break γb, as well as magnetic fraction η, to describe the MWL SED of the source.

First, we fixed the braking index at the value of n = 2.5, and considered three different cases of tage. The MWL SED along with the computed model flux, injected lepton spectrum, and magnetic field at the considered tage are shown in Fig. 1. From the figure, it can be seen that the fit corresponding to tage = 7000 yr is the best out of the three ages considered. For tage = 1000 yr and tage = 4000 yr, the X-ray and radio data can be adequately explained. However, for tage = 1000 yr, the IC emission is not significant enough to explain the high-energy data. Also, for tage = 4000 yr, the IC emission only partially explains it. From this exploration, it is apparent that the age of the PWN lies within 4000 and 7000 yr. In this particular work, we consider tage = 7000 yr to compare our results with those obtained by Joshi et al. (2022) and Yu et al. (2022), albeit these authors did not explore other options. We also discuss the results of considering the true age as a free parameter in a later subsection.

thumbnail Fig. 1.

LHAASO J2226+6057 MWL SEDs along with the calculated model flux (top row), present lepton spectrum (middle row) and magnetic field evolution with time (bottom row) are given for tage = 1000 yr (left column), 4000 yr (middle column) and 7000 yr (right column), for a fixed braking index n = 2.5. The radio data (green) and X-ray data (royalblue) are taken from Pineault & Joncas (2000) and Fujita et al. (2021) respectively. Fermi-LAT (brown), VERITAS (teal), Tibet ASγ (turquoise), MILAGRO (orange) and LHAASO (crimson) data are taken from Xin et al. (2019), Acciari et al. (2009), Tibet ASgamma Collaboration (2021), Abdo et al. (2009b) and Cao et al. (2021), respectively.

After limiting the age of the system, we subsequently explored the effect of the braking index on the evolution of the source. Assuming a fixed tage = 7000 yr, we changed the braking index and tried to describe the MWL SED. We calculated L0 and τ0 using Eqs. (3) and (4). We again changed α1, α2, γb, and η, similarly to the previous case. The calculated spectra for n = 1.5, 2.0, 2.5, and 3.0 are given in Fig. 2. We note that n does not significantly affect the computed SEDs and we have selected n = 2.5 for further study.

thumbnail Fig. 2.

MWL SEDs of LHAASO J2226+6057 for n = 1.5, 2.0, 2.5, and 3.0. The color scheme of the data points is same as given in Fig. 1.

4.2. χ2 fitting of the MWL SED

We have considered typical PWN parameters as initial input, and, using the TIDEFIT code Martín & Torres (2022), subsequently varied them to solve Eq. (1) and computed the best χ2-fitted model spectrum. To find the best-fit spectra, we used tage = 7000 yr and n = 2.5. Similar to the previous discussion, L0 and τ0 have been fixed at the value calculated by Eqs. (3) and (4). The PWN parameters that had already been fixed in the above discussion have also been fixed in this calculation as well, except for the cases of the ejected mass and FIR and NIR energy densities (between 0.01 eV cm−3 and 5 eV cm−3). Since the ejected mass of the progenitor SNR directly affects the size, as well as the magnetic field of the PWN, it was left free within the typical ejected mass range of 7 M⊙–15 M⊙. Apart from these changes and similarly to the previous discussion, α1, α2, γb, and η were left free to vary so as to find the best-fit MWL spectrum. The parameters used in the model are reported in Table 1, in which we have divided the parameters among measured or assumed, derived, and fitted values. The resulting plots are given in Fig. 3. The time evolution of the calculated MWL spectrum and injection spectrum are given in Fig. 4.

Table 1.

Physical parameters used by and resulting from the fit.

thumbnail Fig. 3.

Calculated best-fit MWL spectrum is shown at the top row, along with the MWL data points. The color scheme of the data points are the same as that in Fig. 1. In the bottom panel of the figure, the residuals are also plotted. The color scheme of the residuals are the same as the data points. The middle row shows the timescales of radiative losses, adiabatic losses, and the escape of particles considered in the model (left) and the injected lepton spectrum at the present age (right). Also the time evolution of the magnetic field (left), as well as the SNR forward shock, SNR reverse shock, and PWN radius (right) are given in the bottom row.

thumbnail Fig. 4.

Time evolution of calculated MWL spectrum (left), as well as the injected lepton spectrum (right) are shown in the figure, assuming the parameters given in Table 1.

From the top row of Fig. 3, it can be seen that the computed MWL spectrum using the model along with the 1σ confidence interval, matches with the observed MWL data well. The goodness of fit can also be seen from the bottom residual plot. The systematic uncertainty associated with the model is 0.32, and the χ2/d.o.f. for the given fit is 35.65/30. Also from the figure, as well as from Table 1, it can be seen that FIR and NIR radiation fields do not contribute to the IC emission needed to explain the VHE-UHE data, rather it was found that the CMB is most likely solely responsible as the target photon field required. A recent work by de Oña Wilhelmi et al. (2022) has also considered CMB photons as the most relevant target for IC scattering. However, such a result is somewhat different from what was assumed by Joshi et al. (2022) who considered radiation fields 1.5–3.0 times that of CMB to fit the data.

From the bottom row of Fig. 3 and Table 1, it can be seen that the PWN would be very extended at the present age according to the model, and concurrently the associated magnetic field would be very low and close to the Galactic average value. The large radius of the PWN may agree with the large extension measured by LHAASO, but it contradicts the radio and X-ray sizes observed for the Boomerang PWN (Halpern et al. 2001a,b), as these sizes are much smaller compared to the calculated PWN radius. The required magnetic field to fit the data is also uncomfortably low (we discuss this below), and begs the question of how the particles are confined in such a diluted PWN.

4.3. Possible impact of reverberation

From Fig. 3, and from the resulting values given in Table 1, it can be seen that the radius of the SNR reverse shock is smaller than the PWN radius, which means that the reverse shock has just reached the position of the PWN shell at the onset of the reverberation phase of the PWN. If the age of the system is increased further from the considered age, the effect of reverberation will be apparent on the MWL spectrum. Since the PWN has just started to contract at the considered age of 7000 yr, it is unlikely that it is heavily distorted.

We considered tage = 7000, 8000, 9000, and 10 500 yr to study the effect of reverberation on the MWL spectrum. The braking index was chosen to be 2.5, and L0 and τ0 were calculated using Eqs. (3) and (4) for a specific choice of tage. Since we are only interested in comparing the effect of reverberation on the MWL spectrum based on the assumption of different ages of the PWN, we assumed the values of Porter et al. (2006) for the FIR and NIR energy densities. Similarly to Sect. 4.1, α1, α2, γb, and η were varied to find the most adequately fitted MWL spectrum for each cases of tage. The resulting plot is given in the left panel of Fig. 5. The MWL spectrum gradually deviates away from the observed MWL SED, and consequently the fit worsens. Since the magnetic field increases due to the compression of the PWN during the reverberation phase, the chosen input of the magnetic fraction η was decreased with an increasing considered age, so as to control the fraction of spin-down luminosity that goes on to power the magnetic field of the nebula. Nevertheless, due to the compression of the nebula and high synchrotron burning, the power law index at high energy gradually softens with age, indicating the efficient cooling of the injected high energy electrons to lower energies during reverberation. This in turn affects the spectrum of the resultant synchrotron and IC photons produced. This fact is apparent from the figure, as the calculated MWL spectrum is unable to explain the observed X-ray and the VHE-UHE data present in the SED.

thumbnail Fig. 5.

Impact of reverberation on the PWN MWL spectrum with increasing age is shown (left). The time evolution of the PWN magnetic field (red) and radius (blue) is shown (right) for tage = 6500 yr (dashed), 6800 yr (dotted), 7000 yr (solid), 7200 yr (dot-dashed), and 7500 yr (long-dashed).

We have also considered tage = 6500, 6800, 7200, and 7500 yr and fitted the SED to study how the magnetic field and radius change. The corresponding plot is given in the right panel of Fig. 5. The final PWN radius increases with increasing age, and it also corresponds to the onset of the PWN contraction due to the reverberation in all cases. The PWN must be at the beginning stage of the contraction if the observed SED is to be explained. In all cases, the magnetic field is still low and close to, or lower than, the Galactic average value.

4.4. tage as a free parameter

From the above discussion, it is apparent that the true age of the PWN cannot be much greater than 7000 yr due to reverberation beginning afterwards. As discussed in Sect. 4.1, the true age cannot also be lower than 4000 yr. Thus, we have used TIDEFIT, leaving the true age parameter free. Apart from tage, the free parameters considered in this case are the same as that discussed in Sect. 4.2. The best-fit tage from the fitting of the MWL SED is 4880.3 (3871.4, 6060.4) yr. This value is consistent with the range of the said parameter discussed above. The best fitted model spectrum is similar to that shown in Fig. 3. Although, the best-fit values of the free parameters obtained from the χ2-fitting of the MWL SED, in this case, are not exactly the same as those obtained assuming tage = 7000 yr (see Sect. 4.2 and Table 1), the 1σ confidence intervals of these free parameters for the two cases overlap with each other. Both solutions are essentially the same within uncertainties. So we decided not to report the best-fit results obtained in this case. We have also obtained a comparatively large PWN radius (RPWN = 8.26 pc) and low magnetic field (BPWN = 1.87 μG) associated with the PWN at the best fitted tage, similar to the case discussed in Sect. 4.2 for tage = 7000 yr, which we discuss next.

5. Concluding remarks: Large radius and low magnetic field issue

In this paper, we have provided a detailed, time-dependent, one-zone model to explain the UHE gamma-ray emission observed from the direction of LHAASO J2226+6057 using the emission from PWN associated with PSR J2229+6114. We summarize the main points obtained in this work, and compare them with previous studies.

  1. The effects due to the variation of the true age and braking index were not explored in previous studies performed for this source. From our study, we found that a true age between 4000 and 7000 yr is most suitable for the fitting of the SED of the source. Considering tage as a free parameter during the fitting of the MWL SED also revealed a compatible result. Our study found that no significant effects can be seen on the MWL spectrum if we consider different braking indices.

  2. We used the TIDEFIT code to fit the observed MWL SED by computing the best χ2 fit model spectrum. From the fit, we found that CMB is the target photon field responsible for the IC cooling of the injected leptons from the PWN, whereas the effects of FIR and NIR radiation fields are negligible in this case. Additionally, we have found that the PWN at its current age must be extended (RPWN  ∼  10 pc), as are the extended source regions observed by VERITAS (∼14 pc (0.27°); Acciari et al. 2009) and LHAASO (∼25.6 pc (0.49°); Cao et al. 2021). The obtained PWN radius is larger when compared to, for example, X-ray radii observed for the Boomerang PWN, which is normal in one-zone models.

  3. We have taken the effect of reverberation into account in our modeling. It was found that the PWN is at the onset of compression due to the impact of a reverse shock hitting the shell of the PWN. The effect of reverberation on the MWL spectrum was also explored by gradually increasing the true age of the PWN to the characteristic age. The MWL description worsened with increasing age. The true age cannot be much larger than the considered age of 7000 yr if a PWN is responsible for the gamma-ray emission detected. Since reverberation does not provide a better fit, its inclusion does not solve the large radius and low magnetic field issue.

  4. We have also found that a very low magnetic field (∼2 μG) is needed to explain the MWL SED of the source, which is comparable with the Galactic average magnetic field value. A similarly low magnetic field (a few μG) was found in previous studies as the one needed to describe this source (Joshi et al. 2022; Yu et al. 2022). Moreover, a low magnetic field was found for other LHAASO detected PeVatron cadidates as well (De Sarkar & Gupta 2022; Joshi et al. 2022; Crestan et al. 2021; Li et al. 2021; Burgess et al. 2022), which is an a priori obvious outcome of requesting a leptonically generated high energy emission. Further complications of the model, as we have discussed, do not significantly alleviate this.

As discussed earlier, an uncertainty regarding the distance of the source remains. We have explored this uncertainty by considering two very different values D = 800 pc and D = 7.5 kpc. For D = 7.5 kpc, we have found the fit is comparatively worse compared to that discussed in Sect. 4.2, with χ2/d.o.f. = 36.78/30 and a systematic uncertainty of 0.41. On the other hand, the fit in the case of D = 800 pc is comparable to that discussed in Sect. 4.2, with χ2/d.o.f. = 35.88/30 and a systematic uncertainty of 0.32. There is no clear evidence to overrule one distance over the other, so we only report the case of D = 3 kpc in detail to directly compare our results with those of Joshi et al. (2022) and Yu et al. (2022). It is to be noted that the power law lepton injection spectrum is more favored in the case of D = 800 pc, rather than a broken power law spectrum. However, the issue of a large radius and low magnetic field remains in this case as well.

The large estimated PWN size appears to be a caveat of the PWN interpretation of the LHAASO source. It is likely that both the head and the tail region contribute to the total observed emission from the source. In light of such a complicated source morphology, a one-zone treatment of the source region considered in the model proves to be a simplistic take on the same. The large PWN radius obtained from the calculation could be a consequence of such a simplistic assumption. Nevertheless, our model tends to be the most definitive PWN approach considered thus far to explain the MWL emission from the LHAASO source, complete with a reverberation consideration and PWN true age estimation. It is to be noted that the PWN radius obtained by Joshi et al. (2022) (∼3.1 pc) is different from that obtained in this work. Although one zone model has been considered in both cases, the differences in the formalism adopted to compute the radius evolution in the presence of the background SNR may be the reason behind this. In any case, further investigation will be important to properly address the complexity of the source morphology.

The Galactic magnetic field (GMF) plays an important role in the cosmic ray propagation. The intensity and orientation of GMF are constrained by several methods such as Zeeman-splitting observations (Crutcher 1999); infrared, synchrotron, and starlight polarization studies (Nishiyama et al. 2010; Jaffe et al. 2010; Heiles 1996); and Faraday rotation measures (Han et al. 2006; Pshirkov et al. 2011). The GMF model typically has three components, namely the disc (BDisc), halo (BHalo), and turbulent (BTurb) contributions. Typical values for BDisc and BHalo lie within the range of 2–11 μG (De Sarkar et al. 2021; Di Bernardo et al. 2013); although, the value of the turbulent component depends on the halo height of the Galaxy (Di Bernardo et al. 2013), for example, for a typical halo height of 8 kpc, the value of BTurb comes out to be ≈6 μG (De Sarkar et al. 2021). Although the exact structure of the small-scale GMF is not known yet, from the models given by Pshirkov et al. (2011) and Jansson & Farrar (2012), and also from the observed secondary-to-primary ratio (De Sarkar et al. 2021, and references therein), the large-scale, average GMF can be estimated to be in the range 2–6 μG. The PWN magnetic field (∼2 μG) would be marginally close to, or even lower than, any of the estimates of the average GMF value.

We see two possible ways out of this. On the one side, the local environment of the PWN could have been evacuated by the explosion of the supernova, and/or by earlier explosions, so that the local field in the vicinity is actually much lower than a few μG, allowing for a magnetic field contrast to appear between the PWN and its environment. On the other hand, there is still a chance that the approximate representation of reverberation that we have adopted here is still misleading. We know that assuming the ejecta pressure as constant (as all models of PWNe have so far done) is in fact an oversimplification (Bandiera et al. 2020, 2022). A better treatment of the ejecta pressure could plausibly change the reverberation results, and we shall explore this in the future. Of course, it is also possible that the PWN explanation of the source is not realized at all, and an additional hadronic component provides the dominant contribution. For the moment, a conflicting low value for the PWN magnetic field (not solved by age, braking indices, or the currently assumed behavior of the reverberation process) leaves the PWN origin of LHASSO J2226+6057, and other similar sources, in search for further observational tests.

Acknowledgments

This work has been supported by grant PID2021-124581OB-I00. D.F.T. acknowledges as well USTC and the Chinese Academy of Sciences Presidential International Fellowship Initiative 2021VMA0001. This work was also supported by the Spanish program Unidad de Excelencia María de Maeztu CEX2020-001058-M. A.D.S. acknowledges hospitality provided by the Institute of Space Sciences, where majority of this work was carried out. A.D.S. also acknowledges RRI HPC2020 cluster for providing facilities to perform numerical computations reported in this work. A.D.S. thanks Nayantara Gupta for helpful comments.

References

  1. Abdo, A. A., Ackermann, M., Ajello, M., et al. 2009a, ApJ, 706, 1331 [NASA ADS] [CrossRef] [Google Scholar]
  2. Abdo, A. A., Allen, B. T., Aune, T., et al. 2009b, ApJ, 700, L127 [CrossRef] [Google Scholar]
  3. Abeysekara, A. U., Albert, A., Alfaro, R., et al. 2020, Phys. Rev. Lett., 124, 021102 [NASA ADS] [CrossRef] [Google Scholar]
  4. Acciari, V. A., Aliu, E., Arlen, T., et al. 2009, ApJ, 703, L6 [NASA ADS] [CrossRef] [Google Scholar]
  5. Albert, A., Alfaro, R., Ashkar, H., et al. 2019, arXiv e-prints [arXiv:1902.08429] [Google Scholar]
  6. Albert, A., Alfaro, R., Alvarez, C., et al. 2021, ApJ, 911, L27 [NASA ADS] [CrossRef] [Google Scholar]
  7. Amenomori, M., Bao, Y. W., Bi, X. J., et al. 2019, Phys. Rev. Lett., 123, 051101 [NASA ADS] [CrossRef] [Google Scholar]
  8. Amenomori, M., Bao, Y. W., Bi, X. J., et al. 2021, Phys. Rev. Lett., 127, 031102 [NASA ADS] [CrossRef] [Google Scholar]
  9. Archibald, R. F., Gotthelf, E. V., Ferdman, R. D., et al. 2016, ApJ, 819, L16 [NASA ADS] [CrossRef] [Google Scholar]
  10. Bandiera, R., Bucciantini, N., Martín, J., Olmi, B., & Torres, D. F. 2020, MNRAS, 499, 2051 [NASA ADS] [CrossRef] [Google Scholar]
  11. Bandiera, R., Bucciantini, N., Martín, J., Olmi, B., & Torres, D.F. 2022, MNRAS, submitted [Google Scholar]
  12. Blondin, J. M., Chevalier, R. A., & Frierson, D. M. 2001, ApJ, 563, 806 [NASA ADS] [CrossRef] [Google Scholar]
  13. Breuhaus, M., Hahn, J., Romoli, C., et al. 2021, ApJ, 908, L49 [NASA ADS] [CrossRef] [Google Scholar]
  14. Breuhaus, M., Reville, B., & Hinton, J. A. 2022, A&A, 660, A8 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Burgess, D. A., Mori, K., Gelfand, J. D., et al. 2022, ApJ, 930, 148 [NASA ADS] [CrossRef] [Google Scholar]
  16. Cao, Z. 2010, Chin. Phys. C, 34, 249 [Google Scholar]
  17. Cao, Z., Aharonian, F. A., An, Q., et al. 2021, Nature, 594, 33 [NASA ADS] [CrossRef] [Google Scholar]
  18. Cherenkov Telescope Array Consortium (Acharya, B. S., et al.) 2019, AAS, 240, 32 [Google Scholar]
  19. Chevalier, R. A., & Fransson, C. 1992, ApJ, 395, 540 [NASA ADS] [CrossRef] [Google Scholar]
  20. Crestan, S., Giuliani, A., Mereghetti, S., et al. 2021, MNRAS, 505, 2309 [NASA ADS] [CrossRef] [Google Scholar]
  21. Crutcher, R. M. 1999, ApJ, 520, 706 [NASA ADS] [CrossRef] [Google Scholar]
  22. de Jager, O. C., & Djannati-Ata\"{i}, A. 2009, Astrophys. Space Sci. Lib., 357, 451 [NASA ADS] [CrossRef] [Google Scholar]
  23. de Jager, O. C., Harding, A. K., Michelson, P. F., et al. 1996, ApJ, 457, 253 [Google Scholar]
  24. de Oña Wilhelmi, E., López-Coto, R., Amato, E., & Aharonian, F. 2022, ApJ, 930, L2 [CrossRef] [Google Scholar]
  25. De Sarkar, A., & Gupta, N. 2022, ApJ, 934, 118 [NASA ADS] [CrossRef] [Google Scholar]
  26. De Sarkar, A., Biswas, S., & Gupta, N. 2021, J. High Energy Astrophys., 29, 1 [NASA ADS] [CrossRef] [Google Scholar]
  27. Di Bernardo, G., Evoli, C., Gaggero, D., Grasso, D., & Maccione, L. 2013, J. Cosmol. Astropart. Phys., 2013, 036 [CrossRef] [Google Scholar]
  28. Espinoza, C. M., Lyne, A. G., Kramer, M., Manchester, R. N., & Kaspi, V. M. 2011a, ApJ, 741, L13 [NASA ADS] [CrossRef] [Google Scholar]
  29. Espinoza, C. M., Lyne, A. G., Stappers, B. W., & Kramer, M. 2011b, MNRAS, 414, 1679 [NASA ADS] [CrossRef] [Google Scholar]
  30. Fujita, Y., Bamba, A., Nobukawa, K. K., & Matsumoto, H. 2021, ApJ, 912, 133 [NASA ADS] [CrossRef] [Google Scholar]
  31. Gelfand, J. D., Slane, P. O., & Zhang, W. 2009, ApJ, 703, 2051 [NASA ADS] [CrossRef] [Google Scholar]
  32. H. E. S. S. Collaboration (Abdalla, H., et al.) 2018, A&A, 612, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  33. Halpern, J. P., Camilo, F., Gotthelf, E. V., et al. 2001a, ApJ, 552, L125 [NASA ADS] [CrossRef] [Google Scholar]
  34. Halpern, J. P., Gotthelf, E. V., Leighly, K. M., & Helfand, D. J. 2001b, ApJ, 547, 323 [NASA ADS] [CrossRef] [Google Scholar]
  35. Han, J. L., Manchester, R. N., Lyne, A. G., Qiao, G. J., & van Straten, W. 2006, ApJ, 642, 868 [Google Scholar]
  36. Heiles, C. 1996, ApJ, 462, 316 [NASA ADS] [CrossRef] [Google Scholar]
  37. Jaffe, T. R., Leahy, J. P., Banday, A. J., et al. 2010, MNRAS, 401, 1013 [NASA ADS] [CrossRef] [Google Scholar]
  38. Jansson, R., & Farrar, G. R. 2012, ApJ, 757, 14 [NASA ADS] [CrossRef] [Google Scholar]
  39. Joshi, J. C., Tanaka, S. J., Miranda, L. S., & Razzaque, S. 2022, arXiv e-prints [arXiv:2205.00521] [Google Scholar]
  40. Kothes, R., Uyaniker, B., & Pineault, S. 2001, ApJ, 560, 236 [NASA ADS] [CrossRef] [Google Scholar]
  41. Li, J., Liu, R.-Y., de Oña Wilhelmi, E., et al. 2021, ApJ, 913, L33 [NASA ADS] [CrossRef] [Google Scholar]
  42. Manchester, R. N., & Taylor, J. H. 1977, Irish. Astron. J., 13, 142 [NASA ADS] [Google Scholar]
  43. Manchester, R. N., Durdin, J. M., & Newton, L. M. 1985, Nature, 313, 374 [NASA ADS] [CrossRef] [Google Scholar]
  44. Mart{\’{\i}}n, J., & Torres, D. F. 2022, J. High Energy Astrophys., 36, 128 [CrossRef] [Google Scholar]
  45. Martín, J., Torres, D. F., & Rea, N. 2012, MNRAS, 427, 415 [NASA ADS] [CrossRef] [Google Scholar]
  46. Martín, J., Torres, D. F., & Pedaletti, G. 2016, MNRAS, 459, 3868 [CrossRef] [Google Scholar]
  47. Michel, F. C. 1969, ApJ, 158, 727 [NASA ADS] [CrossRef] [Google Scholar]
  48. Nishiyama, S., Hatano, H., Tamura, M., et al. 2010, ApJ, 722, L23 [NASA ADS] [CrossRef] [Google Scholar]
  49. Pineault, S., & Joncas, G. 2000, AJ, 120, 3218 [NASA ADS] [CrossRef] [Google Scholar]
  50. Pons, J. A., Viganò, D., & Geppert, U. 2012, A&A, 547, A9 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  51. Porter, T. A., Moskalenko, I. V., & Strong, A. W. 2006, ApJ, 648, L29 [CrossRef] [Google Scholar]
  52. Pshirkov, M. S., Tinyakov, P. G., Kronberg, P. P., & Newton-McGee, K. J. 2011, ApJ, 738, 192 [Google Scholar]
  53. Tibet ASgamma Collaboration (Amenomori, M., et al.) 2021, Nat. Astron., 5, 460 [NASA ADS] [CrossRef] [Google Scholar]
  54. Torres, D. F. 2017, ApJ, 835, 54 [NASA ADS] [CrossRef] [Google Scholar]
  55. Torres, D. F., & Lin, T. 2018, ApJ, 864, L2 [NASA ADS] [CrossRef] [Google Scholar]
  56. Torres, D. F., Cillis, A. N., & Rodriguez, J. M. 2012, ApJ, 763, L4 [Google Scholar]
  57. Torres, D. F., Cillis, A., Martín, J., & de Oña Wilhelmi, E. 2014, J. High Energy Astrophys., 1, 31 [NASA ADS] [Google Scholar]
  58. Xin, Y., Zeng, H., Liu, S., Fan, Y., & Wei, D. 2019, ApJ, 885, 162 [NASA ADS] [CrossRef] [Google Scholar]
  59. Yu, H., Wu, K., Wen, L., & Fang, J. 2022, New Astron., 90, 101669 [NASA ADS] [CrossRef] [Google Scholar]

All Tables

Table 1.

Physical parameters used by and resulting from the fit.

In the text

All Figures

thumbnail Fig. 1.

LHAASO J2226+6057 MWL SEDs along with the calculated model flux (top row), present lepton spectrum (middle row) and magnetic field evolution with time (bottom row) are given for tage = 1000 yr (left column), 4000 yr (middle column) and 7000 yr (right column), for a fixed braking index n = 2.5. The radio data (green) and X-ray data (royalblue) are taken from Pineault & Joncas (2000) and Fujita et al. (2021) respectively. Fermi-LAT (brown), VERITAS (teal), Tibet ASγ (turquoise), MILAGRO (orange) and LHAASO (crimson) data are taken from Xin et al. (2019), Acciari et al. (2009), Tibet ASgamma Collaboration (2021), Abdo et al. (2009b) and Cao et al. (2021), respectively.
In the text
thumbnail Fig. 2.

MWL SEDs of LHAASO J2226+6057 for n = 1.5, 2.0, 2.5, and 3.0. The color scheme of the data points is same as given in Fig. 1.
In the text
thumbnail Fig. 3.

Calculated best-fit MWL spectrum is shown at the top row, along with the MWL data points. The color scheme of the data points are the same as that in Fig. 1. In the bottom panel of the figure, the residuals are also plotted. The color scheme of the residuals are the same as the data points. The middle row shows the timescales of radiative losses, adiabatic losses, and the escape of particles considered in the model (left) and the injected lepton spectrum at the present age (right). Also the time evolution of the magnetic field (left), as well as the SNR forward shock, SNR reverse shock, and PWN radius (right) are given in the bottom row.
In the text
thumbnail Fig. 4.

Time evolution of calculated MWL spectrum (left), as well as the injected lepton spectrum (right) are shown in the figure, assuming the parameters given in Table 1.
In the text
thumbnail Fig. 5.

Impact of reverberation on the PWN MWL spectrum with increasing age is shown (left). The time evolution of the PWN magnetic field (red) and radius (blue) is shown (right) for tage = 6500 yr (dashed), 6800 yr (dotted), 7000 yr (solid), 7200 yr (dot-dashed), and 7500 yr (long-dashed).
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.