A Comprehensive Study of the Spectral Variation and the Brightness Profile of Young Pulsar Wind Nebulae

An isolated pulsar releases most of its rotational energy through high-energy relativistic particle outflow called pulsar wind. When the pulsar wind interacts with the surrounding medium and is overpressured, a pulsar wind nebula (PWN) is formed. This drives a wind-termination shock, at which the particles are heated and reaccelerated to the ultrarelativistic regime (see reviews by Gaensler & Slane 2006; Slane 2017, p. 2159). The shocked particles move through the magnetic field in a PWN. They cool down and emit synchrotron radiation from radio to X-rays and inverse-Compton radiation in γ-rays.

Observationally, the injected particle distribution could be inferred from the broadband spectral energy distribution of PWNe. A one-zone model without spatial information can reproduce the observed properties (see, e.g., Pacini & Salvati 1973; Tanaka & Takahara 2010; Bucciantini et al. 2011; Vorster et al. 2013). The improvement of spatial resolution of instruments makes it possible to investigate the detailed structure and the evolution of PWNe. When transporting outward, the high-energy particles lose their energies on a fast timescale and result in a softening of the X-ray spectrum. At the same time, the synchrotron flux decreases. The radial profiles in the radio band have been first modeled as a diffusion nebula with a particle source injected at the center (Gratton 1972; Wilson 1972). A one-dimensional magnetohydrodynamic (MHD) model with a time-independent steady flow has been established from this idea and tested on Crab (Kennel & Coroniti 1984a, 1984b; Porth et al. 2014, 2016).

However, this model is not applicable for a few other PWNe under the unprecedented spatial resolution of Chandra X-ray Observatory (see, e.g., Slane et al. 2004). Combining the photon index (Γ) profile and the surface brightness (Σ) profile, Ishizaki et al. (2017) found that the simple MHD model with steady stream fails to reproduce observed properties in 3C 58 and G21.5−0.9, indicating that spatial diffusion of particles, which has been modeled by Tang & Chevalier (2012), should be taken into account. The pure diffusion model was first proposed to explain the flux and the spectral index of Crab in the optical band (Wilson 1972). The morphology of a young PWN is approximated by a sphere that results from the diffusion of an input particle source into the surrounding space. Later X-ray observations on 3C 58 and G21.5−0.9 show similar Γ distribution and suggest the importance of diffusion in these PWNe (Slane et al. 2000, 2004; Bocchino et al. 2005; Guest et al. 2019). Moreover, Porth et al. (2016) and Ishizaki et al. (2018) conclude that the particle transfer in these two objects is dominated by diffusion through modeling the surface brightness. In addition, the γ-ray halo observed in Geminga also indicates that diffusion plays an important role in PWNe (Abeysekara et al. 2017; Linden et al. 2017; Bao et al. 2021).

Although these works suggest the importance of diffusion, a systematic test of the pure diffusion model is lacking. We aim to test whether this model can be applied to the Γ and Σ profiles in young PWNe that are observed with Chandra. The sample selection is described in Section 2, in which nine PWNe are selected. In Section 3, we describe the X-ray data reduction and the spectral analysis. We perform joint fits of multiple Chandra observations and calculate the observational Γ and Σ profiles. The detailed theoretical modeling is described in Section 4. We follow the model derived by Tang & Chevalier (2012) with transmitting boundary, i.e., particles are not reflected by the boundary and the emission from particles beyond the boundary is not taken into account, and use the full synchrotron spectrum to calculate the Γ and Σ profiles (Ishizaki et al. 2017, 2018). The fitting results are shown in Section 5 and discussed in Section 6. Finally, we summarize our work in Section 7.

Our PWN sample consists of young and bright sources obtained from Table 1 in Kargaltsev & Pavlov (2010) that are well resolved by Chandra with angular size over 20''. Nine PWNe, 3C 58, G21.5−0.9, G54.1+0.3, G291.0−0.1, G310.6−1.6, Kes 75, MSH 15−52, Vela, and Crab, are selected in this study. Their basic properties are listed in Table 1. The distances, B-field strengths, PWN radius, and age are adopted from the literature (see, e.g., Tanaka & Takahara 2011, 2013, and references therein). We obtain all the pulsar properties from the ATNF pulsar catalog ⁶ (Manchester et al. 2005).

We analyzed the Chandra observations taken with the Advanced CCD Imaging Spectrometer (ACIS) and operated in the timed-exposure (TE) mode for each PWN. For G21.5−0.9, we excluded those observations for calculation purposes, with the focal plane offset or the CCD not at optimal temperature. Moreover, we only analyze recent Chandra observations (after 2012) of Crab to reduce the computing time. All the data sets used in this analysis are summarized in Table A1. The radial profiles of Crab, 3C 58, and G21.5−0.9 have been investigated before (see, e.g., Slane et al. 2004; Bamba et al. 2010a, 2010b; Tang & Chevalier 2012; Guest et al. 2019). We, however, reanalyzed the data on our own to ensure that the latest calibration is applied and consistent with other PWNe in our sample.

Note that G10.65+2.96 fits our selection criterion. However, we cannot obtain the difference in spectral properties of inner and outer regions. G11.2−0.3 is one of the youngest PWNe (Borkowski et al. 2016). It has an elongated shape, and most of the PWN emission overlapped the jet viewing from Earth. After removing the jet (see, e.g., Figure 5 in Borkowski et al. 2016), the remaining photons are not statistically significant to measure the spectral difference between the inner and outer PWNe. We excluded them from our samples.

The central pulsar of G291.0−0.1 remains controversial. A young pulsar J1105−6107 was suggested as the central pulsar of this PWN (Kaspi et al. 1997). However, a later broadband study suggested that J1105−6107 is well outside the PWN and unlikely related (Slane et al. 2012). A spin-down luminosity of ∼5.7 × 10³⁷ erg s⁻¹ of the central object CXOU J111148.6−603926 is suggested, although this is estimated from the broadband spectral fit to the PWN (Slane et al. 2012). We note that this is not a precise measurement since a large uncertainty is possible.

All the data reduction and basic analysis were carried out using the Chandra Interactive Analysis of Observations (CIAO) version 4.11 with the calibration database of version 4.8.3. We reprocessed the data using the task chandra_repro. We use reproject_obs to reproject all data sets of a PWN to create a merged event file, and we use flux_obs to create an exposure-corrected image. These images were used to detect the position of pulsars and other point sources. The innermost circular region within 25 radius of the pulsar was removed in the subsequent spectral analysis. We then extracted the PWN spectra from circular annulus regions centered on the pulsar. We adjusted the width of annuli to obtain a sufficient signal-to-noise ratio in each bin (Figure A1). Several targets have jets, which are known to have different spectral behaviors compared to other parts of PWNe and may have different mechanisms of particle transfer (see, e.g., Slane et al. 2004; Temim et al. 2010). Moreover, a few objects like MSH 15−52 have knots or blobs at the opposite direction of jets. We followed previous analysis to remove those structures and background point sources from our spectral analysis (see, e.g., Yatsu et al. 2009).

The source and background spectra, corresponding weighted response matrix files, and ancillary response files (ARFs) of each annulus region are extracted using the specextract command. The background is selected from a nearby source-free region on the same chip. Following the standard analysis procedure for extended sources, we generated weighted ARFs by setting weight=yes and applied no point-source aperture correction by setting correctpsf=no.

We performed the entire spectral analysis using the Sherpa environment. The energy range of the spectral fitting is restricted to 0.5–8 keV to optimize the signal-to-noise ratio except for 3C 58, where we fit its high-energy tail in the 1–8 keV range to avoid the contribution of a soft thermal component below ∼1 keV (Slane et al. 2004; Gotthelf et al. 2007). All the spectra are grouped to have at least 30 X-ray photons in each energy bin. We use an absorbed power-law (PL) model to fit the spectra jointly across all observations for each PWN. We chose absorption model tbabs and set the interstellar abundance according to Wilms et al. (2000), who use the cross section presented in Verner et al. (1993).

We first fit the entire PWN with a single absorbed PL model to constrain the hydrogen column density N_H. This is justified if we assume that the N_H is generally the same among the entire PWN and the spectral behavior of each PWN is steady. Then, we fixed N_H at the best-fit value and fit the spectra of each annulus by linking Γ and the normalization across all observations. The purpose is to reduce the computation time. A few time-variable softening/hardening areas of some PWNe have been reported, but they do not significantly change the overall spectrum of each annulus (see, e.g., Guest & Safi-Harb 2020). We calculated the 0.5–8 keV unabsorbed flux using the sample_flux function with 10⁴ iterations to obtain the Σ profile of each PWN. For each iteration, this command simulates parameter values from a normal distribution around best-fit parameters.

A model based on diffusion was proposed to improve the MHD model, in which the cross-field scattering is small and the particle diffusion can be ignored (de Jager & Djannati-Ataï 2009, p. 451). In the pure diffusion model, the diffusion coefficient (D) is energy independent. We further assume that D and the B-field strength (B_PWN) are constants in the PWN and set the outer boundary to be transparent for particles. Previous studies have suggested that additional advection of particles could improve the fit, but this only dominates the innermost region of Crab (Tang & Chevalier 2012; Ishizaki et al. 2017). On the other hand, diffusion dominates most of the volume of the PWN.

We assume power-law distributed particles $\dot{N}(E,r=0)={{KE}}^{-\alpha }$ injected from the center of a PWN, where E is the particle energy and α is the power-law index. The particles diffuse from the center and distribute spherically symmetry. As the particles extend to a radius R, which corresponds to the outer bound of a PWN, we adopt a simple case that the particles transmit out and leave the system. The B-field strength and the PWN radius R are predetermined parameters using the broadband spectral energy distribution fitting results obtained from Tanaka & Takahara (2011, 2013). (see Table 1). Following Gratton (1972), the particle number density distribution is

$$\begin{eqnarray}&&N(E,r,t)=\displaystyle \frac{K}{4\pi {rD}}{E}^{-\alpha }{f}_{\alpha }(u,v),\end{eqnarray} \tag{ 1 }$$

where f_α (u, v) is an integral

$$\begin{eqnarray}&&{f}_{\alpha }(u,v)=\displaystyle \frac{1}{\sqrt{\pi }}{\int }_{v}^{\infty }{\left(1-\displaystyle \frac{u}{x}\right)}^{\alpha -2}\displaystyle \frac{{e}^{-x}}{\sqrt{x}}{dx},\end{eqnarray} \tag{ 2 }$$

and u is defined as

$$\begin{eqnarray*}&&u=\displaystyle \frac{{r}^{2}{QE}}{4D},\end{eqnarray*}$$

where

$$\begin{eqnarray*}&&Q=\displaystyle \frac{4{\sigma }_{T}}{3{m}^{2}{c}^{3}}\displaystyle \frac{{B}^{2}}{8\pi }\end{eqnarray*}$$

is the factor of synchrotron energy loss. The lower limit of the integral in Equation (2) is

$$\begin{eqnarray*}v=\left\{\begin{array}{ll}\tfrac{{r}^{2}}{4{Dt}} & t\lt \tfrac{1}{{QE}}\\ \tfrac{{r}^{2}{QE}}{4D} & t\gt \tfrac{1}{{QE}}\end{array}\right.,\end{eqnarray*}$$

where t is the age of the PWN. Once we obtain the energy distribution of particles, we calculate the synchrotron spectrum at any radius with an assumed magnetic field. The normalization constant K can be derived from the luminosity of the PWN. We assume a fraction η of the spin-down energy that is converted to the energy of particles above 1 TeV. This can be written as

$$\begin{eqnarray}&&{\int }_{{E}_{0}}^{{E}_{\max }}E\dot{N}(E,r=0,t){dE}=\eta {L}_{\mathrm{sd}},\end{eqnarray} \tag{ 3 }$$

where E₀ = 1 TeV is the lower bound of the integration and ${E}_{\max }$ is the maximum value of the electrostatic energy in the polar cap of the central pulsar. By performing an integration to $\dot{N}(E,r=0,t)={{KE}}^{-\alpha }$, we obtain

$$\begin{eqnarray}&&K=\displaystyle \frac{\eta {L}_{\mathrm{sd}}\left(\alpha -2\right)}{1-{\left({E}_{\max }/{E}_{0}\right)}^{2-\alpha }}{E}_{0}^{\alpha -2}.\end{eqnarray} \tag{ 4 }$$

In our analysis, we use the full synchrotron spectrum to normalize the energy distribution of particles for the emission calculation. We followed the procedure presented in Tanaka & Takahara (2010, 2011) and Ishizaki et al. (2017) to calculate the synchrotron power. The spectral flux density j_ν (r) can be written as

$$\begin{eqnarray}&&{j}_{\nu }(r)=\int N(E,r,t){p}_{\nu }(E){dE},\end{eqnarray} \tag{ 5 }$$

where

$$\begin{eqnarray}&&{p}_{\nu }(E)=\displaystyle \frac{\sqrt{3}{e}^{3}B}{{{mc}}^{2}}F\left(\displaystyle \frac{\nu }{{\nu }_{c}}\right),\end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray}&&{\nu }_{c}=\displaystyle \frac{3}{4\pi }\displaystyle \frac{{eB}}{{mc}}{\gamma }^{2},\end{eqnarray} \tag{ 7 }$$

and $F(x)=x{\int }_{x}^{\infty }{K}_{5/3}(y){dy}$, with x = ν/ν_c and K_5/3 the modified Bessel function of order 5/3 (Rybicki & Lightman 1986). The specific surface brightness M_ν can be integrated as

$$\begin{eqnarray}&&{M}_{\nu }({r}_{\perp })={\int }_{C}\displaystyle \frac{{j}_{\nu }(r)}{4\pi }{ds},\end{eqnarray} \tag{ 8 }$$

where ds is a line element along the line of sight C and r_⊥ is the distance perpendicular to the line of sight from the central pulsar, i.e., angular distance. In a spherically symmetric system, this line integral can be written explicitly as

$$\begin{eqnarray}&&{M}_{\nu }({r}_{\perp })=2{\int }_{{r}_{\perp }}^{{r}_{N}}\displaystyle \frac{{j}_{\nu }(r)}{4\pi }\displaystyle \frac{{rdr}}{\sqrt{{r}^{2}-{r}_{\perp }^{2}}},\end{eqnarray} \tag{ 9 }$$

where r_N is the radius of the nebula. The surface brightness Σ(r_⊥) can be obtained by integrating over the energy range from ν₁ to ν₂ as follows:

$$\begin{eqnarray}&&{\rm{\Sigma }}({r}_{\perp })={\int }_{{\nu }_{1}}^{{\nu }_{2}}{M}_{\nu }({r}_{\perp })d\nu .\end{eqnarray} \tag{ 10 }$$

Finally, the spectral index Γ (from ν₁ to ν₂) can be calculated as

$$\begin{eqnarray}&&{\rm{\Gamma }}=\displaystyle \frac{\mathrm{log}[{M}_{{\nu }_{1}}({r}_{\perp })/{M}_{{\nu }_{2}}({r}_{\perp })]}{\mathrm{log}({\nu }_{1}/{\nu }_{2})}.\end{eqnarray} \tag{ 11 }$$

We fit the radial profiles of Γ and Σ with the pure diffusion model described in Section 4. Three parameters, α, D, and η, are tuned to fit the profiles, while other parameters are fixed. Our objective is to determine whether the simple diffusion could reproduce the characteristic of the profile. Therefore, we did not perform a formal χ² fit because it is difficult to discuss the goodness of fit statistically owing to the large systematic uncertainty of this simple model. The best-fit results are listed in Table 2. We performed two fits for each source: one fits the Γ profile (Γ fit), and another fits the Σ profile (Σ fit). Both profiles for each PWN are shown in Figure 1. Note that these two fits are not independent. The shapes of both profiles are determined once α and D are chosen. The parameter η only tunes the normalization of the Σ profile. We discuss the parameter dependence in Section 6.

Zoom In Zoom Out Reset image size

The pure diffusion model generally well describes the sample PWNe. Both Γ and Σ fits yield consistent results in G21.5−0.9, G291.0−0.1, G310.6−1.6, and Kes 75. These PWNe have no complex structure compared to Crab, Vela, and MSH 15−52. However, some PWNe are relatively faint with large uncertainties in Γ that could smear detailed structures in the profile. On the other hand, large discrepancies between the parameters yielded from the Γ and Σ fits are seen in 3C 58, G54.1+0.3, and possibly MSH 15−52. There are, however, a couple of poor fits, in particular Crab and Vela. They have significant changes in slope in both profiles, and their brightness profiles are concave down. Moreover, these are especially observed in the innermost region, where the advection may play an important role (Tang & Chevalier 2012; Ishizaki et al. 2017).

In Figure 2 we plot the Γ and Σ profiles versus the physical size. All the PWNe have a similar shape in both profiles except for Crab and Vela. The Γ monotonically and smoothly increases with distance from the pulsar in nearly all sources. The overall slope resides in a narrow range. All the sampled PWNe have similar Σ profiles except for Crab and Vela. We discuss these apparent features and the physical interpretation using the diffusion model in Section 6.

Zoom In Zoom Out Reset image size

6.1. Applicability of the Pure Diffusion Model

6.2. Degeneracy of Model Parameters

Based on our pure diffusion model with transmitting boundary, three parameters are critical for the Γ profile: B-field in the PWN, D, and the injected power-law index α of the particle energy. Moreover, η plays an important role in the Σ profile. However, the B-field strength is a predefined value and varies dramatically among studies; see, e.g., G21 (Tanaka & Takahara 2011, 2013; Guest et al. 2019). This parameter is coupled with D, and therefore it is meaningless to directly compare D among PWNe. Taking 3C 58 as an example, we demonstrate how the model profile changes as we change these parameters in Figure 3. From Figures 3(a) and (b), the effects of B and D are coupled with each other (Tang & Chevalier 2012). A flatter Γ profile, which implies a lower energy dissipation, can be achieved by both decreasing B and increasing D, and vice versa. This is because the spatial extent of the X-ray-emitting particles is defined by the distance over which they spread before exhausting their energy via synchrotron radiation. The ratio of the diffusion length of a particle emitting a photon of a given frequency ν to a radius R can be written as $\sqrt{4{{Dt}}_{\mathrm{cool}}}/R\propto {\left({R}^{-2}{\nu }^{-1/2}{{DB}}^{-3/2}\right)}^{1/2}$, where t_cool = 1/QE is the cooling timescale of the synchrotron emission. The Σ profile is the same if DB^−3/2 remains constant. Therefore, DB^−3/2 is the key parameter for the shape of the Γ profile. Since the PWN B-field strength is difficult to estimate and it often depends on many assumptions, we will not discuss the effect of D. Instead, we focus on DB^−3/2 in the following discussions. The DB^−3/2 values of the sampled PWNe are listed in Table 2. The Γ fitting results of 3C 58 and G21.5−0.9 are in the same order as those derived in Tang & Chevalier (2012) with full synchrotron spectrum.

Zoom In Zoom Out Reset image size

We further test the sensitivity of the effect of parameters on the Σ profile (Figures 3(e) and (f)). The DB^−3/2 value changes not only the shape but also the displacement (or interception) of model profiles. Comparing to the Γ profile, we notice that the shape of the Σ profile becomes insensitive to DB^−3/2 as long as DB^−3/2 is large enough. From Figure 3(c), the offset of Γ, which implies the general hardness of the PWN, is dominated by α of the injected high-energy particles. On the other hand, the displacement of the Σ profile could be affected by DB^−3/2, α, and η. Finally, η only affects the displacement in the Σ profile because it represents the energy-transferring efficiency from the spin-down power of the pulsar to the PWN.

6.3. Profile Characteristics and Model Parameters

Based on the discussion in Section 6.2, we test whether the characteristics, i.e., slope and mean value of the Γ and Σ profiles, can be interpreted with the pure diffusion parameters. To achieve this goal, we simply fit the Γ profile of each PWN with a straight line. Crab and Vela are not included in this analysis since their profiles cannot be described with our model. The slope can roughly represent the Γ gradient. We take the mean of Γ as the value at the half of the PWN radius predicted from the best-fit straight line. For the Σ profile, we define the same slope and the mean Σ as we did for the Γ profile. The difference is that we use the logarithmic scale for the surface brightness.

Then, we plot the apparent features (slope and mean value of the Γ and the Σ profiles) versus the model parameters in Figure 4. The most significant correlation is, not surprisingly, between the mean of Γ and α. The linear correlation coefficient is ρ = 0.89 with a null hypothesis probability p = 0.006 for the Γ fit, where the Σ fitting result is ρ = 0.96 with p = 6 × 10⁻⁴. A relatively weak correlation can be seen between DB^−3/2 and the slope of the Γ profiles with ρ = −0.85 and p = 0.02 for the Γ fit. The correlation is insignificant for the Σ fit with ρ = −0.40 and p = 0.3. This is caused by the large discrepancy between the Γ and Σ fitting results of a few PWNe where the Σ fit cannot well describe the Γ profile.

Zoom In Zoom Out Reset image size

The slope of the Σ profile is correlated with DB^−3/2, where the correlation analysis suggests ρ = 0.75 and p = 0.05 for both the Σ and Γ fits. The relation between these two parameters might be nonlinear. Therefore, we also calculate Kendall's rank correlation and obtain a Kendall's τ = 0.81 with p = 0.01 for the Γ fit and τ = 0.52 with p = 0.1 for the Σ fit. This suggests a possible correlation and is consistent with the result in Section 6.2, in which we show that a larger D (or smaller B) results in a flatter Σ profile. Interestingly, the slope of the Σ profile anticorrelates with η with ρ = −0.84 and p = 0.01 for the Σ fit (ρ = −0.88 and p = 0.009 for the Γ fit), suggesting a similar significance to that between DB^−3/2 and the Σ profile. This is not expected, as we discussed in Section 6.2, and could indicate a physical connection underneath. Finally, the mean surface brightness likely anticorrelates with DB^−3/2 and positively correlates with η; these are expected in Section 6.2. However, the mean surface brightness does not correlate with α, suggesting that α is not the dominant factor for the Σ profile.

6.4. Correlation between Model Parameters and Physical Parameters

We tested whether the best-fit parameters and physical parameters correlate with each other. We performed correlation analysis between four parameters, including DB^−3/2, α, η, and the PWN age. The results are shown in Figure 5. It is worth mentioning that η of Kes 75 is the largest among the selected PWNe. Given that PSR J1846−0258 is a high B-field pulsar that shows a magnetar-like burst and outburst behavior, it is possible that the dissipation of the extremely strong magnetic field results in additional energy output from the pulsar (Gavriil et al. 2008; Kumar & Safi-Harb 2008; Ng et al. 2008).

Zoom In Zoom Out Reset image size

In general, we found no significant correlations between these parameters, indicating that these parameters are independent of each other. We found that DB^−3/2 may weakly anticorrelate with $\mathrm{log}(\eta )$ with ρ = −0.7 and p = 0.07 for the Γ fit. This possible correlation is not observed for the Σ fit with ρ = −0.5 and p = 0.3, suggesting no significance. We could not draw a strong conclusion with our current sample. However, this could be connected to the anticorrelation between the slope of the brightness profile and η in Section 6.3. A larger DB^−3/2 implies a steeper Γ gradient in the radial profile, i.e., stronger synchrotron energy loss. Therefore, a high energy conversion rate might not be needed.

We applied a pure diffusion model to fit the Γ and Σ profiles of nine young PWNe, including well-studied 3C 58 and G21.5−0.9; a few faint sources, G54.1+0.3, G291.0−0.1, G310.6−1.6, and Kes 75; and three bright PWNe, Crab, MSH 15−52, and Vela. We analyzed archival Chandra observations and fit the spectra of each PWN with a power-law model. Except for Crab and Vela, all the sampled PWNe show similar Γ and Σ profiles. Γ increases slowly with the distance from the central pulsar. This behavior cannot be described by the MHD model, which suggests an abrupt increase in Γ (Slane et al. 2004; Tang & Chevalier 2012; Ishizaki et al. 2017) and therefore suggests that particle diffusion dominates the particle transfer rather than a steady flow for most PWNe. We derived the brightness profile from the pure diffusion model and found that the observed profiles are consistent with the model prediction. However, the parameters that can describe the Γ profile and those for the Σ profiles are not consistent with each other in a few cases like 3C 58 and G54.1+0.3. The discrepancy implies that a simple diffusion model may not be enough to produce the entire distribution of PWNe and future improvement is necessary. The discrepancy could be attributed to the intrinsic spherical asymmetry, the particle reflection by the PWN boundaries, energy-dependent diffusion, and spatially dependent B-field. The last two could be specified with multiwavelength studies. A model with a spatially varying magnetic field could be constructed in the future. We further investigate the relationship between the model parameters and the apparent feature of both profiles. The effects of the diffusion coefficient and the magnetic field are coupled with each other, and the observed Γ and the Σ gradient highly depend on the DB^−3/2 value. The mean hardness of a PWN significantly correlates with the distribution of the injected particles, while the brightness of a PWN is dominated by the energy conversion efficiency η. Finally, we found no significant correlations between model parameters and the physical quantities of PWNe.

We thank the anonymous reviewer for valuable comments that improved this paper. This research is based on the data obtained from the Chandra Data Archive and has made use of the software provided by the Chandra X-ray Center (CXC) in the application packages CIAO, ChIPS, and Sherpa. C.-P.H. acknowledges support from the Ministry of Science and Technology in Taiwan through grant MOST 109-2112-M-018-009-MY3. W.I. is supported by JSPS KAKENHI grant No. 21J01450. C.-Y.N. is supported by a GRF grant of the Hong Kong Government under HKU 17305416. S.J.T. acknowledges support by the Aoyama Gakuin University-Supported Program "Early Eagle Program."

Facility: CXO (ACIS). -

Software: CIAO (Fruscione et al. 2006), Sherpa (Freeman et al. 2001).

We list all data sets used in this analysis in Table A1.

Figure A1 shows the Chandra 0.5–8 keV images of all the sampled PWNe. To study the radial profile, we use multiple annuli centered on the central pulsar to select X-ray photons and collect X-ray spectra.

Zoom In Zoom Out Reset image size