SPECTRAL AND MORPHOLOGICAL ANALYSIS OF THE REMNANT OF SUPERNOVA 1987A WITH ALMA AND ATCA

The evolution of the remnant of supernova (SN) 1987A in the Large Magellanic Cloud has been closely monitored since the collapse of its progenitor star, Sanduleak (Sk) −69°202, on 1987 February 23. Models of Sk −69°202 indicated that it had an initial mass of ∼20 M_☉ (Hillebrandt et al. 1987). The mass range of the progenitor is consistent with the formation of a neutron star (NS; Thielemann & Arnett 1985), and thus with the neutrino events reported by the Kamioka NDE (Hirata et al. 1987) and IMB (Bionta et al. 1987; Haines et al. 1988) detectors. Models by Crotts & Heathcote (2000) suggest a transition from red supergiant into blue supergiant to explain the hourglass nebula structure, which envelopes the SN with three nearly stationary rings (Chevalier & Dwarkadas 1995; Blondin & Lundqvist, 1993; Martin & Arnett 1995; Morris & Podsiadlowski 2007). The two outer rings, imaged with the Hubble Space Telescope (HST; Jakobsen et al. 1991; Plait et al. 1995), are located on either side of the central ring in the equatorial plane (equatorial ring, ER), and likely formed at the same time as the ER (Crotts & Heathcote 2000; Tziamtzis et al. 2011). The synchrotron emission, generated by the shock propagating into the clumpy circumstellar medium (CSM) close to the equatorial plane, was detected in the mid-1990s (Turtle et al. 1990; Staveley-Smith et al. 1992), and has become brighter over time (Manchester et al. 2002; Zanardo et al. 2010). Radio observations have stretched from flux monitoring at 843 MHz with the Molonglo Observatory Synthesis Telescope (Staveley-Smith et al. 1993; Ball et al. 2001) to images of sub-arcsec resolution with the Australia Telescope Compact Array (ATCA; Gaensler et al. 1997; Manchester et al. 2005a; Ng et al. 2008; Potter et al. 2009; Zanardo et al. 2013). ATCA observations at 94 GHz (Lakićević et al. 2012b) have been followed by observations from 100 GHz up to 680 GHz with the Atacama Large Millimeter/submillimeter Array (ALMA; Kamenetzky et al. 2013; Indebetouw et al. 2014).

The ongoing shock expansion has been monitored at 9 GHz since 1992 (Gaensler et al. 1997, 2007). Shock velocities of ∼4000 km s⁻¹ have been measured between day 4000 and 7000 (Ng et al. 2008), while signs of a slower expansion have been tentatively detected after day ∼7000, as the shock has likely propagated past the high-density CSM in the ER (Ng et al. 2013). Similar evidence of slower shock expansion since day ∼6000 has been found in X-ray data (Park et al. 2005, 2006; Racusin et al. 2009) as well as in infrared (IR) data (Bouchet et al. 2006).

Since the early super-resolved images at 9 GHz (Gaensler et al. 1997), a limb-brightened shell morphology has been characteristic of the remnant. The radio emission, over the years, has become more similar to an elliptical ring rather than the original truncated-shell torus (Ng et al. 2013). The radio remnant has shown a consistent east–west asymmetry peaking on the eastern lobe, which has been associated with higher expansion velocities of the eastbound shocks (Zanardo et al. 2013). The asymmetry degree appears to have changed with the shock expansion, as images at 9 GHz exhibit a decreasing trend in the east–west asymmetry since day ∼7000 (Ng et al. 2013). High-resolution observations at 1.4–1.6 GHz (Ng et al. 2011) via Very Long Baseline Interferometry (VLBI) with the Australian Large Baseline Array (LBA) have highlighted the presence of small-scale structures in the brightest regions in both lobes.

The relation between the radio emission and the synchrotron spectral indices, α (S_ν∝ν^α), has been investigated via both flux monitoring (Manchester et al. 2002; Zanardo et al. 2010) and imaging observations (Potter et al. 2009; Lakićević et al. 2012b; Zanardo et al. 2013) with the ATCA. The progressive flattening of the radio spectrum derived from 843 MHz to 8.6 GHz at least since day 5000, coupled with the e-folding rate of the radio emission, has pointed to an increasing production of non-thermal electrons and cosmic rays (CR) by the shock front (Zanardo et al. 2010). On the other hand, the association of steeper spectral indices with the brightest eastern sites implies a higher injection efficiency on the eastern side of the supernova remnant (SNR; Zanardo et al. 2013). Flatter spectral indices in the center of the remnant have been tentatively identified from low-resolution two-frequency spectral maps (Potter et al. 2009; Lakićević et al. 2012b), while at 18–44 GHz the central and center–north regions have −0.5 ≲ α ≲ −0.3 (Zanardo et al. 2013). With ALMA, the spectral energy distribution (SED) of the remnant has been mapped where the non-thermal and thermal components of the emission overlap, identifying cold dust in the SNR interior (Indebetouw et al. 2014), which accommodates 0.4–0.7 M_☉ of the dust mass discovered with the Herschel Space Observatory (Herschel; Matsuura et al. 2011) in the ejecta.

This paper combines the results presented by Indebetouw et al. (2014) with a comprehensive morphological and spectral analysis of SNR 1987A based on both ATCA and ALMA data. In Section 2, we present the ALMA Cycle 0 super-resolved images before and after subtraction, in the Fourier plane, of the synchrotron and dust components (Section 3). In Section 4, we assess the remnant asymmetry from 44 to 345 GHz. In Section 5, we update the SED derived by Indebetouw et al. (2014), while, in Section 6, we investigate the spectral index variations in the SNR across the transition from radio to far infrared (FIR). In Section 7, we discuss possible particle flux injection by a pulsar situated in the inner regions of the remnant.

The ATCA and ALMA observations used in this study were performed in 2011 and 2012. ATCA observations at 44 and 94 GHz are detailed in Zanardo et al. (2013) and Lakićević et al. (2012b), respectively. ALMA observations were made in 2012 (Cycle 0) from April to November, over four frequency bands: Band 3 (B3, 84–116 GHz, λ 3 mm), Band 6 (B6, 211–275 GHz, λ 1.3 mm), Band 7 (B7, 275–373 GHz, λ 850 μm), and Band 9 (B9, 602–720 GHz, λ 450 μm). Each band was split over dual 2 GHz wide sidebands, with minimum baselines of 17 m (B9) to maximum baselines of 400 m (B3). All observations used quasars J0538-440 and J0637-752 as bandpass and phase calibrators, respectively. Callisto was observed as an absolute flux calibrator in B3 and B6, while Ceres was used in B7 and B9 (see also Kamenetzky et al. 2013). It is noted that, while ALMA is designed to yield data with flux density calibration uncertainty as low as ∼1%, in Cycle 0 this uncertainty is estimated at ∼5% at all frequencies. Relevant observational parameters are listed in Table 1 (see also Table 1 in Indebetouw et al. 2014).

Each data set was calibrated with the casa²¹ package, then exported in miriad²² for imaging. After clean-ing (Högbom 1974), both phase and amplitude self-calibration were applied in B3 over a 2 minute solution interval, while only phase calibration was applied in B6 and B7. No self-calibration was performed in B9. As in Zanardo et al. (2013), we note that since the self-calibration technique removes position information, each image was compared with that prior to self-calibration and, in case of positional changes, the self-calibrated images were shifted. Further adjustments were made in the comparison with the ATCA observations at 44 GHz, based on prominent features on the eastern lobe and location of the remnant center. As from Zanardo et al. (2013), the 44 GHz image was aligned with VLBI observations of the SNR (G. Zanardo et al., in preparation). Adding in quadrature these positional uncertainties and the accuracy of the LBA VLBI frame, the errors in the final image position are estimated at ∼60 mas.

Deconvolution was carried out via the maximum entropy method (MEM; Gull & Daniell 1978) in B3, B6, and B7. A weighting parameter of robust = 0.5 (Briggs 1995) was used in all bands. The resultant diffraction-limited images, which have central frequency at 102 GHz in B3, 213 GHz in B6, 345 GHz in B7, and 672 GHz in B9, were then super-resolved with a circular beam of 08 in B3, 06 in B6, and 03 in B7 and B9. The diffraction limited and super-resolved images are shown in the first column of Figures 1 and 2, below the ATCA image at 94 GHz (Lakićević et al. 2012b). Integrated Stokes I flux densities, dynamic range and related rms are given in Table 2.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

To decouple the non-thermal emission from that originating from dust, the synchrotron component, as resolved with ATCA at 44 GHz (Zanardo et al. 2013), and the dust component, as imaged with ALMA at 672 GHz (B9; Indebetouw et al. 2014), were separately subtracted from the data sets at 94, 102, 213, 345, and 672 GHz. All subtractions were performed in the Fourier plane, via miriad task uvmodel, where the model flux density at 44 GHz was scaled to fit the SNR emission over the ER (I₄₄), while the B9 model flux density was scaled to fit the emission localized in the central region of the remnant (I_B9). Scaling of the 44 GHz model was tuned by minimizing the flux density difference on the brighter eastern lobe, without over-subtracting in other regions of the remnant. To separate the emission in the SNR center, the central flux was firstly estimated by fitting a Gaussian model via miriad task uvfit. The image model at 672 GHz was then scaled to match the flux of the Gaussian model. The scaling factor was further tuned to minimize over-subtraction.

Super-Nyquist sampling was applied in all images, using a pixel size of 8 mas to avoid artifacts when sources are not at pixel centers. Deconvolution via MEM was carried out on the residual images obtained from the subtraction of I_B9, while standard CLEANing was applied to the residuals obtained from the I₄₄ subtraction. All diffraction-limited-subtracted images are shown in the central and right columns of Figure 1, while Figure 2 shows the residuals after super-resolution with the circular beam used for the original images. All image parameters are given in Table 2.

The flux densities were determined by integrating within polygons enclosing the SNR emission. Uncertainties in the flux densities include uncertainties in the image fitting/scaling process combined with the uncertainty in the flux density calibration. We note that the residual images from the subtraction of both models at 44 and 672 GHz were not considered, since the error attached to the double subtraction exceeds the total integrated flux density. Only at 213 GHz the error–flux margin is minor, as the diffraction-limited image obtained after the dual subtraction has integrated flux density of 2.1 ± 1.9 mJy (Figure 3).

Zoom In Zoom Out Reset image size

ALMA observations of SNR 1987A capture both the remnant emission from the ER and that from the SNR interior, where the dense ejecta sit (Kamenetzky et al. 2013; Indebetouw et al. 2014).

While the image at 102 GHz barely resolves the two-lobe distribution of the emission, the images at 213 and 345 GHz clearly show the ringlike emission morphology, localized around the ER (see Figures 1 and 2). It is understood that the radio emission over the ER is primarily synchrotron emission, generated by the interaction of the SN shock with the dense CSM near the equatorial plane, which results in a magnetic-field discontinuity where particles are accelerated (e.g., Zanardo et al. 2010, 2013). As shown in early models of Type II SNRs (Chevalier 1982), the region of interaction between the SN blast and the CSM consists of a double-shock structure, with a forward shock, where ambient gas is compressed and heated, and a reverse shock, where the ejecta are decelerated. Between the two shocks, the reflected shocks, due to the forward shock colliding with the dense ER (Borkowski, Blondin & McCray 1997), propagate inward (Zhekov et al. 2009, 2010). The reverse shock at first expanded outward, behind the forward shock, but might have been inward-moving since day ∼7000. As discussed by Ng et al. (2013), the ring synchrotron emission is currently localized between the forward and reverse shocks, and likely has components from both the ER and high-latitude material above the equatorial plane. Truncated-shell torus models of the remnant geometry at 9 GHz indicate that the half-opening angle has been decreasing since day ∼7000, and is estimated at ∼27° at day 9568 (Ng et al. 2013).

At 345 GHz the SNR interior is brighter, while at 672 GHz the emission is predominantly localized in the central region of the remnant. Since the emission from this region rises steeply with frequency, it has been identified with thermal dust emission (Indebetouw et al. 2014), as dust grains, probably heated by ⁴⁴Ti decay and X-ray emission from the reverse shock (Larsson et al. 2013), emit strongly in the FIR regime. The central emission, visible both in B7 and B9, appears to extend over the inner optical ejecta (see Figure 4). In particular, as noted by Indebetouw et al. (2014), this inner emission shows a north–south elongation, which in B9 can be identified between P.A. ∼20° and P.A. ∼30°, similar to that seen with HST (Larsson et al. 2011, 2013). From both Figures 1 and 2, it can also be noticed that the SNR emission at 672 GHz includes possible emission located to the NW (see contour overlays in Figure 5). This NW feature has signal-to-noise ratio (S/N) of 3.7σ, and integrated flux density of 6.8 ± 0.5 mJy, i.e., ∼10% of the total integrated flux density at 672 GHz.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

3.1. Subtracted Images

The east–west asymmetry of the synchrotron emission, primarily associated with the emission morphology over the ER, is investigated from 44 to 345 GHz, as shown in Figure 6, where all images are restored with a 07 circular beam. At 213 and 345 GHz the central dust emission is subtracted, via a scaled model of the flux density as resolved in B9 (see Section 2).

Zoom In Zoom Out Reset image size

Considering radial sections crossing the SN site (Reynolds et al. 1995), the SNR asymmetry is first estimated as the ratio between the eastern peak of the radial slice crossing the maximum emission on the eastern lobe (see black/white profile in Figure 6), and the western peak of the radial slice crossing the maximum emission on the western lobe (see red profiles in Figure 6). This ratio, A_p, emphasizes the hot spots on each side of the remnant. Alternatively, the ratio A_i between the total flux densities integrated over the eastern and western halves of the image, is derived by splitting the image at the R.A. associated with the SN site. For comparison with previous asymmetry estimates (e.g., Zanardo et al. 2013), the ratio A₉₀ between the eastern and western peaks of the radial slice at P.A. 90° (see dark green profile in Figure 6) is also derived. All values for A_p, A_i, and A₉₀, both for the images restored with a 07 circular beam (as in Figure 6) and the diffraction-limited images, are plotted in Figure 7 and listed in Table 3.

Zoom In Zoom Out Reset image size

Table 3. Asymmetry Ratios

Image	$A_{p}\equiv \frac{S_{\nu }(E_{{\rm Max}})}{S_{\nu }(W_{{\rm Max}})}$ a		$A_{i}\equiv \left. \int \limits _E S_{\nu }/\int \limits _W S_{\nu }\right.$ b		$A_{90}\equiv \left. S_{\nu } \left(\frac{E_{{\rm Max}}}{W_{{\rm Max}}} \right) \right|_{90^{\circ }}$c
GHz	07d	DL e	07	DL	07	DL
44	1.61 ± 0.08	1.42 ± 0.05	1.36 ± 0.02	1.53 ± 0.07	1.69 ± 0.08	1.54 ± 0.05
94	1.55 ± 0.05	⋅⋅⋅	1.18 ± 0.04	⋅⋅⋅	1.88 ± 0.07	⋅⋅⋅
102	1.44 ± 0.04	⋅⋅⋅	1.10 ± 0.04	⋅⋅⋅	1.40 ± 0.07	⋅⋅⋅
213	1.28 ± 0.03	1.43 ± 0.05	1.24 ± 0.06	1.23 ± 0.02	1.43 ± 0.02	1.49 ± 0.05
213−IB9f	1.44 ± 0.03	⋅⋅⋅	1.30 ± 0.02	⋅⋅⋅	1.19 ± 0.03	⋅⋅⋅
345	⋅⋅⋅	1.28 ± 0.03	⋅⋅⋅	0.78 ± 0.04	⋅⋅⋅	1.16 ± 0.03
345−IB9f	1.30 ± 0.02	⋅⋅⋅	0.80 ± 0.02	⋅⋅⋅	1.48 ± 0.04	⋅⋅⋅

Notes. ^aRatio between the eastern peak of the radial slice crossing the maximum emission on the eastern lobe, and the western peak of the radial slice crossing the maximum emission on the western lobe (see black/white and red profiles in Figure 6). ^bRatio between the total flux densities integrated over the eastern and western halves of the image. The image is split at the R.A. associated with the SN site (Reynolds et al. 1995). ^cRatio between the eastern and western peaks of the radial slice at P.A. 90° (see dark green profile in Figure 6). ^dImages resolved with a 07 circular beam (see Figure 6). ^eDiffraction-limited (DL) images (see Figure 1). ^fImages derived after subtraction of the model flux density at 672 GHz (Band 9, B9), scaled to fit the central emission (see Figure 6).

Download table as:  ASCIITypeset image

From Figure 6 it can be seen that the A₉₀ ratio does not fully capture the asymmetry changes with frequency, since the hot spot in the western lobe, which becomes brighter from 102 to 345 GHz, is located southward of the 90° profile. In Figure 7, the linear fits derived for both A_p and A_i ratios show a consistent decrease as frequencies reach the FIR. At 345 GHz, A_i values indicate that the east–west asymmetry is reversed, thus matching the asymmetry trend seen in recent HST images (Larsson et al. 2011), where the western side of the ring is markedly brighter. As discussed Section 3, the morphology similarities between the optical image of the SNR and the super-resolved image at 345 GHz are evident (see Figure 4).

The change of the remnant's east–west asymmetry over time has been discussed by Ng et al. (2013), as the result of a progressive flattening of the shock structure in the equatorial plane, due to the shock becoming engulfed in the dense UV–optical knots in the ER, coupled with faster shocks in the east side of the remnant. While X-ray observations do not show significant difference between the NE and SW reverse shock velocities, although the NE sector is brighter (Fransson et al. 2013), faster eastbound outer shocks have been measured in the radio (Zanardo et al. 2013) and point to an asymmetric explosion of a binary merger as SN progenitor (Morris & Podsiadlowski 2007, 2009). As the SN blast is gradually overtaking the ER, faster expanding shocks in the east would exit the ER earlier than in the west.

The effects of the asymmetric shock propagation are likely to emerge in the transition from radio to FIR rather than at lower frequencies, due to the shorter synchrotron lifetime at higher frequencies. To estimate the synchrotron lifetime in the FIR range, we use the approximation that, in a magnetic field of strength B, all the radiation of an electron of energy E is emitted only at the critical frequency ν_c (Rybicki & Lightman 1979). Considering the electrons orbit is inclined at a pitch angle θ to the magnetic field, the synchrotron lifetime, τ_e, can be derived as a function of ν_c (e.g., Condon 1992)

$$\begin{equation} \tau _{e}\equiv {E \over |dE/dt| }\sim 1.06\times 10^{9} (B\,{\rm sin}\theta)^{-3/2} \nu _{c}^{-1/2}, \end{equation} \tag{ 1 }$$

where τ_e is expressed in years, B in μG, and ν in GHz. As the electrons are expected to have an isotropic distribution of pitch angles, 〈sin²θ〉 = 2/3; for 200 ≲ ν ≲ 400 GHz, and assuming a magnetic field strength at the shock front of 10 ≲ B ≲ 20 mG (Berezhko & Ksenofontov 2006; Berezhko et al. 2011), in the radio/FIR transition we can estimate 20 ≲ τ_e ≲ 80 yr. However, since it is likely that sub-diffusive particle transport (Kirk et al. 1996) is taking place in regions of the SNR, in conjunction with efficient CR acceleration (Glushkov et al. 2013), local magnetic-field amplifications could exceed the above limits by at least an order of magnitude (Bell & Lucek 2001). In this scenario, the local τ_e would be of the order of months. Therefore, with regions in the remnant where electrons might be unable to cross the emission sites within their radiative life-time, the synchrotron emission at FIR frequencies would require the presence of relatively fresh injected and/or re-accelerated electrons to match the emission distribution at lower frequencies.

The SED for ATCA data from 1.4 to 94 GHz and ALMA data is shown in Figure 8. To match the average epoch of the ALMA data, the ATCA data are scaled to day 9280, via exponential fitting parameters derived for ATCA flux densities from day 8000, as measured at 1.4, 8.6, and 9 GHz (G. Zanardo et al., in preparation; Staveley-Smith et al. 2014). Across the transition from radio to FIR frequencies, the observed spectrum consists of the sum of thermal and non-thermal components.

Zoom In Zoom Out Reset image size

As described in Section 2, to identify the dust component of the emission from the inner regions of the remnant, the B9 model flux density, S_B9, has been scaled to fit the emission measured in the SNR central region at 94–345 GHz (S_B9 fit, hollow red circles in Figure 8). The subtraction of S_B9 fit from the visibilities at 94–345 GHz yields the residual flux densities indicated as S_ν − S_B9 fit (purple bars in Figure 8), as for the images shown in the central column of Figure 1. The flux densities S_B9 fit derived in B6 and B7, together with the total integrated flux density measured in B9 (S_B9), although obtained via a different reduction technique, have been associated by Indebetouw et al. (2014) with dust grains, in conjunction with data from Herschel (Matsuura et al. 2011) and the Atacama Pathfinder Experiment (APEX; Lakićević et al. 2012a; Figure 8).

Similarly (see Section 2), to separate the non-thermal emission from that thermal, the 44 GHz model flux density, S₄₄, has been scaled to fit the toroidal component of the emission at 94–672 GHz (blue/cyan diamonds in Figure 8). By fitting the resulting S_44 fit components and the ATCA flux densities at 1.4–44 GHz, we obtain, at day 9280, the synchrotron spectral index α_s = −0.727 ± 0.020, with $S_{\nu }(\alpha _{s})\propto \nu ^{\alpha _{s}}$. The spectral index measured from 1.4 to 94 GHz, i.e., for ATCA data only, is α_ATCA = −0.735 ± 0.028. While α_s is slightly flatter than α_ATCA, both values are consistent with the progressive flattening of the radio spectrum measured since day ∼5000 (Zanardo et al. 2010). The subtraction of S_44 fit from the visibilities at 94–672 GHz, gives the residual flux densities S_ν − S_44 fit (red/orange diamonds in Figure 8), as for the images shown in the right column of Figure 1. In B3, B6, and B7, S_ν − S_44 fit > S_B9 fit, i.e., the residuals exceed the emission expected from the dust.

While the subtraction of the flux densities is inevitably affected by errors (see Table 2), given that the I₄₄-subtracted images have a S/N ≳ 7σ and the residual emission appears primarily located westward of the optical ejecta (see Section 3.1), we investigate the nature of this emission excess as: (1) free–free emission from an ionized fraction of the inner ejecta, (2) synchrotron emission from a compact source located in the inner regions of the remnant, and (3) emission from grains of very cold dust.

5.1. Free–Free Emission

To estimate the free–free radiation in the SNR as imaged with ALMA, we hypothesize an ionized portion of the ejecta as an approximately spherical region, located inside the ER. Using the beam size of the super-resolved I₄₄-subtracted image at 102 GHz as an upper limit, we consider the radius of the spherical region up to R_s ∼ 040 (≈3.08 × 10¹⁷ cm). Such radius covers the extent of the inner ejecta as imaged in the optical (Larsson et al. 2013), and stretches to the likely radius of the reverse shock, qualitatively identified with the inner edge of the emission over the ER (see Figure 4). Given that the pre-SN mass has been estimated between 14 M_☉ and 20 M_☉ (Smartt 2009), if we assume that the Hii region in the ejecta has a uniform density ρ_ej(r) = 3 M_H ii(4πr³)⁻¹, we set 0.7 M_☉ ≲ M_H ii ≲ 2.5 M_☉. The lower limit, M_H ii ∼ 0.7 M_☉, i.e., ρ_ej(r = R_s) = 1.1 × 10⁻²⁰, represents a partial ionization of the ejecta by X-ray flux, either within the inner region or on the outer layer. The upper limit, M_H ii ∼ 2.5 M_☉, i.e., ρ_ej(r = R_s) = 4.1 × 10⁻²⁰ g cm⁻³ matches the density model by Blinnikov et al. (2000) scaled to the current epoch (see Figure 21 in Fransson et al. 2013), and corresponds to complete ionization of the H and He mass within R_s.

The optical depth associated with the Hii region along the line of sight (los) can be estimated as

$$\begin{equation} \tau _{{\rm ff}}\approx 3.28\times 10^{-7} \,T_{4}^{-1.35} \nu ^{-2.1} N_{e}^{2} R_{s}, \end{equation} \tag{ 2 }$$

where T₄ = T_e/(10⁴ K), ν is in GHz, N_e is in cm⁻³ and R_s ≈ ∫_losdl is in pc. For T₄ ∼ 1, given the emission measure 4.6 × 10⁶ ≲ EM ≲ 5.9 × 10⁷ cm⁻⁶ pc, where ${\rm EM}=N^2_eR_s$, at frequencies 102 ≲ ν ≲ 672 GHz the emission becomes nearly transparent as 2.28 × 10⁻⁶ ≲ τ_ff ≲ 1.2 × 10⁻³. The flux associated with the ionized component of the ejecta, can then be derived as

$$\begin{equation} S _{\nu _{{\rm ff}}} \propto \tau _{{\rm ff}}(\nu) {2 k T \over \lambda ^2} \Omega \propto \nu ^{-0.1}. \end{equation} \tag{ 3 }$$

Considering the solid angle Ω subtended by the same radius, R_s, at all frequencies, for the lower limit M_H ii ≈ 0.7 M_☉, i.e., N_e = 6.8 × 10³ cm⁻³, $S _{\nu _{{\rm ff}}} \sim (S_{\nu }-S_{44\,{{\rm fit}}})$ at 102 and 213 GHz, as Equation (3) yields $S_{102_{{\rm ff}}}\sim 3.9$ mJy, $S_{213_{{\rm ff}}}\sim 3.2$ mJy, $S_{345_{{\rm ff}}}\sim 3.1$ mJy, and $S_{672_{{\rm ff}}}\sim 2.9$ mJy. These $S _{\nu _{{\rm ff}}}$ values (blue stars in in Figure 8) would well fit the SED (see blue fit in Figure 8). If the Hii density is considerably higher, as given by the upper limit M_H ii ≈ 2.5 M_☉, the derived $S _{\nu _{{\rm ff}}}$ fluxes would exceed the emission residuals by an order of magnitude.

As a constraint to the free–free emission component, the hypothesized Hii region would also produce optical and near-IR Hi recombination lines. The resultant Hα flux can be estimated as

$$\begin{equation} S_{\rm H\alpha } \approx \epsilon \, N_{e} N_{+} \, \left ({{4}\over {3}} \pi R_{s}^3 \right) {{f} \over {4 \pi d^{2}}}, \end{equation} \tag{ 4 }$$

where ∝α_A is the emissivity of the Hα line per unit volume, with α_A = α_A(T) = ∑_nα_n(T) the total recombination coefficient, f is the volume filling factor and d is the distance to the source. If one adopts $\epsilon \approx 3.53\times 10^{-25}\,T_{4}^{-0.92}$ erg cm⁻³ s⁻¹ for T₄ ≃ 0.5 (e.g., Storey & Hummer 1995), the number of ions N₊ = N_e = 6.8 × 10³ cm⁻³ as for the lower limit assumed for the density of the Hii region, and R_s ≈ 3.08 × 10¹⁷ cm, for f ∼ 1 one obtains S_Hα ≈ 1.26 × 10⁻¹¹ erg cm² s⁻¹. To match the Hα flux from the core measured by Fransson et al. (2013) (see Figure 8 therein) at ∼1.4 × 10⁻¹⁴ erg cm² s⁻¹ on day 9000, the estimated S_Hα would have to undergo 7.4 mags of extinction. Such extinction is possible but improbable. We note that the above S_Hα also exceeds the Hα emission by the reverse shock, measured at ∼2.0 × 10⁻¹³ erg cm² s⁻¹ on day 9000, and associated with a density 250 ≲ N₊ ≲ 750 cm⁻³ (Fransson et al. 2013).

As regards the magnitude of the flux from the Br-γ line, the model estimate by Kjær et al. (2010) on day 6840 is a factor of ∼10⁻² smaller than the S_Hα derived from Equation (4), while this is mainly associated with emission from the hot spots in the ER.

In the absence of external X-ray heating, both the heating and ionization would be powered by radioactive decays in the core. The resultant flux, estimated via models as in Kozma & Fransson (1998), would be several orders of magnitude smaller than the lower limit for $S_{\nu _{{\rm ff}}}$.

Another possible source of ionizing emission is a pulsar wind nebula (PWN) in the SNR interior. The properties of any PWN are very uncertain (see Section 7), but there is the expectation that line emission would accompany free–free emission also in this case.

5.2. Flat-spectrum Synchrotron Emission

5.3. Dust Emission

If the excess emission is due to an additional synchrotron component, S_ν(α_p)∝ν^−0.2, this would provide a constraint to the net dust emission, $S_{\nu _{d}}$, in ALMA data. The subtraction $S_{\nu }-S_{44_{\rm fit}}-S_{\nu }(\alpha _{p})$, where $S_{\nu _{p}}=S_{0_{p}}(\nu / \nu _{0})^{\alpha _{p}}$ with $S_{0_{p}}=3.5\pm 3.1$ mJy as from Table 2, leads to $S_{213_{d}}\sim 0.7$ mJy, $S_{345_{d}}\sim 4.9$ mJy and $S_{672_{d}}\sim 45.0$ mJy (black triangles in Figure 8). We note that $S_{672_{d}}$ coincides with the integrated flux density of the central feature of the related image, which extends over the inner ejecta as seen with HST (see Figure 4). The net dust can be fitted via a modified Planck curve of thermal emission, as

$$\begin{equation} S_{\nu _{d}}(\beta,T)= {M_{\rm dust} \over d^{2}} \, \kappa _{\nu } \, B_{\nu }(T)\propto \nu ^{\,\beta }\,B_{\nu }(T), \end{equation} \tag{ 5 }$$

where M_dust is the dust mass, κ_ν = 3 Q_ν/4ρ a_ν is the dust mass absorption coefficient, with Q_ν the absorption efficiency for spherical grains of density ρ and radius a_ν, 1 ≲ β ≲ 2 for insterstellar dust (e.g., Cortese et al. 2014 and references therein) and B_ν(T) is the Planck function. As shown in Figure 8, the best fit of both ALMA fluxes, as reported in this paper, and Herschel fluxes, as from Matsuura et al. (2011), yields β = 1.23 and dust temperature T = 22.2 K. The thermal peak of the SED has been previously fitted (Matsuura et al. 2011; Lakićević et al. 2011, 2012b, 2012a) with temperatures estimated between 17 and 26 K, while Indebetouw et al. (2014) fit amorphous carbon dust at T = 26 ± 3 K.

The sum of the main synchrotron component, S_ν∝ν^−0.73, and the emission component from dust grains at T ∼ 22 K is lower than the measured ALMA flux densities at 213 and 345 GHz (see dashed gray fit in Figure 8). To match the emission excess of ∼3 mJy in this frequency range, we could also postulate a second dust component. This would require very cold dust at T ≲ 5 K, i.e., at temperatures where the assumption of either amorphous carbon or silicates leads to dust masses implausibly large for physically realistic grains. In particular, as from Equation (5), at 345 GHz a flux density of ∼3 mJy would require dust at T ∼ 3 K with M_dust ∼ 50 M_☉, as obtained by using k_ν ∼ 2.5 cm² g⁻¹ for amorphous carbon (Zubko et al. 1996, 2004). We note that warmer dust at T ∼ 180 K has been identified by Dwek E. et al. (2010) in the ER, where the dust grains are likely collisionally heated by the expanding radiative shocks (Bouchet et al. 2006).

The spectral index distribution across the remnant is investigated via multi-frequency spectral index maps and the T–T plot method (Costain 1960; Turtle et al. 1962).

Spectral index maps are derived from images at 44, 94, 102, 213, and 345 GHz. In Figure 9, the maps resulting from the combination of images at two, three, four and five frequencies are shown. All images are reduced with identical procedure and restored with a circular beam of 08. ATCA data at 44 and 94 GHz are scaled to day 9280, via exponential fitting parameters derived for ATCA flux densities from day 8000, as measured at 8.6 and 9 GHz (G. Zanardo et al., in preparation). While the two-frequency maps are derived from direct division of the flux densities, the maps derived from three to five frequencies are more accurate, as the solution minimizes the function

$$\begin{equation} {\boldsymbol{\varepsilon }}_{pq} ={\sum _{k=1}^{n}}\,{\left |{ \,S_{\nu _{k}}^{^{\,pq}}\left(\alpha _{M}^{^{pq}}\right) - S_{\nu _{0}}^{^{\,pq}} \left ({\nu _{k} \over \nu _{0} } \right)^{\alpha _{M}^{^{pq}}} } \right|}^{2} \end{equation} \tag{ 6 }$$

at each pixel of coordinates (p, q), for images at ν_k frequencies, with n ⩾ 3.

Zoom In Zoom Out Reset image size

In Figure 9, it can be seen that the spectral indices overall become flatter as the frequencies reach the FIR. For all maps, the spectral index, α_M, varies between −1.25 and 0.25, while in the five-frequency map the index range narrows to −0.96 ⩽ α_M ⩽ 0.18. In most of the multi-frequency maps, α_M is steeper on the eastern half of the SNR, while flat spectral indices surround the center, where the bulk of dust sits (see Section 3), and extend onto the NW and SW regions of the remnant. With reference to the five-frequency map, −0.9 ≲ α_M ≲ −0.6 on the eastern lobe and −0.4 ≲ α_M ≲ 0 on the western side of the SNR, with α_M ∼ 0 around the central region and −0.4 ≲ α_M ≲ −0.3 predominantly in the NW quadrant (P.A. ∼330°) and in the SW quadrant (P.A. ∼210°).

The T–T plot method is applied to two images at different frequencies, each obtained with identical reduction process and with the same angular resolution. The spectral variations are assessed over image regions not smaller than the beam size, where the spectral index, α_TT, is determined from the flux density slope m, where S₂ = m S₁ + q. The regions used for the T–T plots are shown in Figure 10. Six frequency pairs are considered: 44−102 GHz, 102−213 GHz, and 213−345 GHz; and 44−94 GHz, 44−213 GHz, and 44−345 GHz. Different box sizes are used to suit the controlling beam of each frequency pair. All derived α_TT values are listed in Table 4.

Zoom In Zoom Out Reset image size

Table 4. Regional Spectral Indices

Regionc	Fitd	αTT(S1, S2)a e						$\alpha _{M_\mu }(S_{1,2,3,4})$b e
		(44, 102)	(102, 213)	(213, 345)	(44, 94)	(44, 213)	(44, 345)	(44, 94, 213, 345)
N	y	⋅⋅⋅	⋅⋅⋅	−0.33 ± 0.02	⋅⋅⋅	−0.38 ± 0.03	⋅⋅⋅	⋅⋅⋅
	ODR	⋅⋅⋅	⋅⋅⋅	−0.41 ± 0.02	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅
	Abs	⋅⋅⋅	⋅⋅⋅	−0.40 ± 0.02	⋅⋅⋅	−0.37 ± 0.03	⋅⋅⋅	⋅⋅⋅
	(0, 0)	⋅⋅⋅	⋅⋅⋅	−0.05 ± 0.04	⋅⋅⋅	−0.43 ± 0.05	⋅⋅⋅	⋅⋅⋅
N1	y	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.28 ± 0.08	⋅⋅⋅	−0.39 ± 0.02	−0.39 ± 0.09
	ODR	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.31 ± 0.09	⋅⋅⋅	−0.42 ± 0.02	⋅⋅⋅
	Abs	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.27 ± 0.08	⋅⋅⋅	−0.35 ± 0.03	⋅⋅⋅
	(0, 0)	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.36 ± 0.14	⋅⋅⋅	−0.40 ± 0.04	⋅⋅⋅
N2	y	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.73 ± 0.03	⋅⋅⋅	−0.81 ± 0.08	−0.85 ± 0.14
	ODR	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.76 ± 0.03	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅
	Abs	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.77 ± 0.03	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅
	(0, 0)	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.84 ± 0.04	⋅⋅⋅	−0.89 ± 0.07	⋅⋅⋅
E	y	−0.95 ± 0.02	−0.07 ± 0.02	−1.05 ± 0.04	⋅⋅⋅	−0.48 ± 0.02	⋅⋅⋅	⋅⋅⋅
	ODR	−0.95 ± 0.01	−0.07 ± 0.03	⋅⋅⋅	⋅⋅⋅	−0.60 ± 0.05	⋅⋅⋅	⋅⋅⋅
	Abs	−0.95 ± 0.02	−0.06 ± 0.03	−1.12 ± 0.04	⋅⋅⋅	−0.47 ± 0.02	⋅⋅⋅	⋅⋅⋅
	(0, 0)	−0.97 ± 0.02	−0.03 ± 0.03	0.11 ± 0.07	⋅⋅⋅	−0.46 ± 0.02	⋅⋅⋅	⋅⋅⋅
E0	y	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.55 ± 0.03	⋅⋅⋅	−0.68 ± 0.03	−0.60 ± 0.03
	ODR	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.59 ± 0.03	⋅⋅⋅	−0.72 ± 0.04	⋅⋅⋅
	Abs	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.54 ± 0.03	⋅⋅⋅	−0.66 ± 0.03	⋅⋅⋅
	(0, 0)	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.69 ± 0.07	⋅⋅⋅	−0.71 ± 0.06	⋅⋅⋅
E1	y	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.65 ± 0.04	⋅⋅⋅	−0.66 ± 0.04	−0.52 ± 0.13
	ODR	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.68 ± 0.04	⋅⋅⋅	−0.68 ± 0.04	⋅⋅⋅
	Abs	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.70 ± 0.05	⋅⋅⋅	−0.65 ± 0.04	⋅⋅⋅
	(0, 0)	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.78 ± 0.08	⋅⋅⋅	−0.68 ± 0.04	⋅⋅⋅
E2	y	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.60 ± 0.03	⋅⋅⋅	−0.73 ± 0.03	−0.55 ± 0.12
	ODR	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.63 ± 0.03	⋅⋅⋅	−0.74 ± 0.02	⋅⋅⋅
	Abs	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.65 ± 0.03	⋅⋅⋅	−0.72 ± 0.03	⋅⋅⋅
	(0, 0)	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.74 ± 0.04	⋅⋅⋅	−0.73 ± 0.03	⋅⋅⋅
S	y	⋅⋅⋅	⋅⋅⋅	−0.05 ± 0.05	⋅⋅⋅	−0.45 ± 0.02	⋅⋅⋅	⋅⋅⋅
	ODR	⋅⋅⋅	⋅⋅⋅	−0.10 ± 0.04	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅
	Abs	⋅⋅⋅	⋅⋅⋅	−0.23 ± 0.03	⋅⋅⋅	−0.44 ± 0.02	⋅⋅⋅	⋅⋅⋅
	(0, 0)	⋅⋅⋅	⋅⋅⋅	−0.22 ± 0.03	⋅⋅⋅	−0.53 ± 0.07	⋅⋅⋅	⋅⋅⋅
S1	y	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.61 ± 0.02	⋅⋅⋅	−0.52 ± 0.07	−0.63 ± 0.07
	ODR	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅
	Abs	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.60 ± 0.02	⋅⋅⋅	−0.56 ± 0.07	⋅⋅⋅
	(0, 0)	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.95 ± 0.18	⋅⋅⋅	−0.67 ± 0.09	⋅⋅⋅
W	y	−0.77 ± 0.02	−0.11 ± 0.03	0.25 ± 0.04	⋅⋅⋅	−0.40 ± 0.06	⋅⋅⋅	⋅⋅⋅
	ODR	−0.78 ± 0.02	−0.11 ± 0.02	0.17 ± 0.02	⋅⋅⋅	−0.43 ± 0.06	⋅⋅⋅	⋅⋅⋅
	Abs	−0.77 ± 0.02	−0.10 ± 0.03	0.08 ± 0.04	⋅⋅⋅	−0.43 ± 0.06	⋅⋅⋅	⋅⋅⋅
	(0, 0)	−0.76 ± 0.02	−0.18 ± 0.03	0.16 ± 0.03	⋅⋅⋅	−0.49 ± 0.05	⋅⋅⋅	⋅⋅⋅
W0	y	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.69 ± 0.03	⋅⋅⋅	−0.49 ± 0.09	−0.53 ± 0.05
	ODR	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.77 ± 0.04	⋅⋅⋅	−0.57 ± 0.10	⋅⋅⋅
	Abs	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.54 ± 0.10	⋅⋅⋅
	(0, 0)	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.66 ± 0.03	⋅⋅⋅	−0.49 ± 0.11	⋅⋅⋅
W1	y	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.93 ± 0.09	⋅⋅⋅	−0.61 ± 0.09	−0.37 ± 0.11
	ODR	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.97 ± 0.08	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅
	a	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.99 ± 0.08	⋅⋅⋅	−0.60 ± 0.09	⋅⋅⋅
	(0, 0)	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.75 ± 0.10	⋅⋅⋅	−0.51 ± 0.15	⋅⋅⋅
W2	y	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.40 ± 0.11	⋅⋅⋅	−0.50 ± 0.02	−0.47 ± 0.11
	ODR	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.44 ± 0.11	⋅⋅⋅	−0.52 ± 0.02	⋅⋅⋅
	Abs	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.35 ± 0.12	⋅⋅⋅	−0.49 ± 0.02	⋅⋅⋅
	(0, 0)	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.52 ± 0.14	⋅⋅⋅	−0.51 ± 0.02	⋅⋅⋅
C	y	−0.86 ± 0.08	−0.00 ± 0.04	1.20 ± 0.04	⋅⋅⋅	−0.43 ± 0.02	⋅⋅⋅	⋅⋅⋅
	ODR	−0.90 ± 0.08	−0.00 ± 0.03	1.97 ± 0.09	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅
	Abs	−0.86 ± 0.08	−0.00 ± 0.03	2.25 ± 0.11	⋅⋅⋅	−0.44 ± 0.01	⋅⋅⋅	⋅⋅⋅
	(0, 0)	−0.84 ± 0.08	−0.09 ± 0.04	⋅⋅⋅	⋅⋅⋅	−0.43 ± 0.02	⋅⋅⋅	⋅⋅⋅
C1	y	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.40 ± 0.08	⋅⋅⋅	⋅⋅⋅	−0.01 ± 0.06
	ODR	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅
	Abs	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.47 ± 0.08	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅
	(0, 0)	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅	−0.73 ± 0.15	⋅⋅⋅	⋅⋅⋅	⋅⋅⋅

Notes. ^aThe spectral indices α_TT are derived from five images at 44, 94, 102, 213, and 345 GHz. The images are analyzed in six frequency pairs, as indicated in Figure 10. ^bThe median spectral index, $\alpha _{M_\mu }$, is derived by fitting a Gaussian to the histogram of the α_M values resulting in the spectral map obtained from images at 44, 94, 213, and 345 GHz in the corresponding T–T region (see the top map in Figure 10). ^cThe regions selected for the T–T plots are designated in Figure 10. ^dThe given spectral indices are derived from different linear interpolations: least-square fit from vertical squared errors, α_y; orthogonal distance regression, α_ODR; robust fitting from the square root of absolute residuals, α_Abs; linear regression with forced zero interception (S₁ = 0, S₂ = 0), α_(0, 0). ^eThe errors on α_TT and $\alpha _{M_\mu }$ are the 1σ error on the slope of the linear fit combined with the uncertainty in the flux calibration.

Download table as:  ASCIITypeset image

Similarly to the trend of α_M, α_TT values become flatter at higher frequencies (see Figure 11). From 102 to 213 GHz, α_TT ≈ −0.1 across the whole remnant (see Figure 10), while from 213 to 345 GHz the spectral distribution appears markedly split in two larger regions (see also the two-frequency map in Figure 9), with very flat indices on the western half of the SNR and steep indices on the eastern side. For the frequency pairs 102–213 and 213–345 GHz, the T–T plots for the eastern region (see E in Table 4) give Δα_TT = α_(102, 213) − α_(213, 345) = 1.0 ± 0.07. This could be indication of a local spectral break at 213 GHz. Using ν_c = 213 GHz in Equation (1), for B = 20 mG, τ_e ∼ 35 yr. However, given the high Mach number ($\mathcal {M}\sim 10^{4}$) of the eastbound shocks (Zanardo et al. 2013), likely sub-diffusive particle transport (Kirk et al. 1996) by the shock front and, consequently, local magnetic-field amplifications (Bell 2004), it is possible that the CR in the eastern lobe, radiating at ν_c, are already past their synchrotron lifetime.

Zoom In Zoom Out Reset image size

Flat spectral indices in the western lobe extend north and south at both 213–345 GHz and 44–213 GHz (see Figure 10), while T–T plots from higher resolution images show a narrower north–south alignment of the flat regions (see N1, C1 and W2 in Figure 10 and Table 4), with −0.5 ≲ α_TT ≲ −0.3. The T–T plots for 44–94 GHz and 44–345 GHz also yield the steepest spectral indices in region W1, this might be due to the local emission drop in the NW sector of the ER, visible in the images at 94 and 345 GHz at P.A. ∼300°.

As discussed in Section 5, the flat-spectrum western regions could be linked to a PWN. We note that spectral maps of the Crab Nebula via observations centered at 150 GHz (Arendt et al. 2011) have identified spectral indices around −0.2 in the inner central regions of the PWN, while spectral indices ∼ − 0.3 have been associated with the PWN periphery.

As discussed in Section 3.1, the emission at 102 and 213 GHz in the I₄₄-subtracted images appears to peak west of the SN site (Reynolds et al. 1995; Figures 1 and 2), and to mainly extend west of the optical ejecta (see Figure 5). Besides, both the spectral maps and T–T plots (Figures 9 and 10) show that flat spectral indices can be associated with the center–west regions of the SNR (see Section 6). These results could be explained by a possible PWN, powered by a pulsar likely located at a westward offset from the SN position.

The pulsar-kick mechanism has been linked to asymmetries in the core collapse or in the subsequent SN explosion, presumably due to asymmetric mass ejection and/or asymmetric neutrino emission (Podsiadlowski et al. 2005; Wongwathanarat et al. 2013; Nordhaus et al. 2012). As evidence for the natal kick, NS mean three-dimensional (3D) speeds have been estimated at $\overline{v}_{\rm NS}\sim 400\pm 40$ km s⁻¹ (Hobbs et al. 2005), while a transverse velocity of ∼1083 km s⁻¹ has been detected by Chatterjee et al. (2005), which would imply a 3D NS birth velocity as high as 1120 km s⁻¹ (Chatterjee et al. 2005). In the context of SNR 1987A, by day 9280 the NS could have travelled westward of the SN site by ∼20–80 mas, while for an impulsive kick of the same order as $\overline{v}_{\rm NS}$ a distance of ∼42 ± 5 mas would have been covered. With a western offset of ∼005, the NS would be situated inside the beam of the I₄₄-subtracted images associated with the peak flux density, both at 102 and 213 GHz. If we take into account the ∼60 mas uncertainty intrinsic to image alignment (see Section 2), as well as the error of 30 mas in each coordinate of the SN VLBI position (Reynolds et al. 1995), the NS could be located near the emission peak as seen in the I₄₄-subtracted images at 102 and 213 GHz. We note that the inner feature of fainter emission detected in the SNR at 44 GHz, as aligned with VLBI observations of the remnant (Zanardo et al. 2013), is centered ∼60 mas west of the SN site.

If a pulsar is embedded in the unshocked ejecta, the PWN would be in its early stages of evolution, likely surrounded by uniformly expanding gas (Chevalier & Fransson 1992). Diffuse synchrotron emission from the PWN would be due to the relativistic particles, produced by the pulsar, accelerated at the wind termination shock (Kirk et al. 2009).

Assuming a power-law energy distribution of electrons, i.e., the particle density N_e is expressed as N_e(E)∝KE^−s, where s = 1–2α and K∝(m_ec²)^{s − 1}, the synchrotron emission of the PWN can be written as

$$\begin{equation} S_{\nu _{P}} \propto {K B_{{\rm {PWN}}}^{ { {1\over 2} (s+1) } }} {\nu \,^{{ {1\over 2} (1-s) } }} , \end{equation} \tag{ 7 }$$

where B_PWN is the nebular magnetic field strength. Noting that, in radio observations, the energy in electrons cannot be separated from that in the magnetic field (Reynolds et al. 2012), the equipartition magnetic-field strength could be derived as (e.g., Longair 2011; see revised formula by Beck & Krause 2005; Arburtina et al. 2012)

$$\begin{equation} B_{\rm PWN} \approx \left [ G_{0} \,G \,(\mathcal {K}+1) \, { S_{\nu } \over f d\, \theta _{\rm PWN}^{3}} \,\, \nu ^{{ {(1-s)\over 2} } } \right ] ^{{{{ { 2 \over (5+s)} }}}} , \end{equation} \tag{ 8 }$$

where G₀ is a constant, G = G(ν, s) is the product of different functions varying with the minimum and maximum frequencies associated with the spectral component and the synchrotron spectral index (Beck & Krause 2005; Longair 2011), $\mathcal {K}$ is the ion/electron energy ratio, f is the volume filling factor of radio emission, and θ_PWN = R_PWN/d is the angular radius. Considering 005 ≲ R_PWN ≲ 015, −0.4 ⩽ α ⩽ −0.1 (as from Section 5), 102 ⩽ ν ⩽ 672 GHz, S_ν ≈ 3 mJy, and taking $\mathcal {K}\approx 100$ (Beck & Krause 2005) while f ≈ 0.5, Equation (8) leads to 1 ≲ B_PWN ≲ 7 mG. For 0.2 ≲ f ≲ 0.9, 2 ≲ B_PWN ≲ 5 mG. Since the equipartition is a conjecture for young SNRs and no longer valid when the spectral index is flatter than −0.5, these B_PWN estimates might be inaccurate.

The energy inside the PWN, due to the PWN magnetic field, can be simplified as

$$\begin{equation} E_{{\rm PWN}, B}(t)\sim V_{\rm PWN}(t)\frac{B^{2}_{\rm PWN}}{8\pi }, \end{equation} \tag{ 9 }$$

where the magnetic field is considered uniform and isotropic inside the PWN volume, $V_{\rm PWN}(t)=4/3\, \pi \,R_{\rm PWN}^{3}(t)$. For 005 ≲ R_PWN ≲ 015, and, as for the parameters used in Equation (8), 1 ≲ B_PWN ≲ 7 mG, at t = t₂₅ = 8.0 × 10⁸ s (≈9280 days), we estimate 0.9 × 10⁴³ ≲ E_{PWN, B} ≲ 1.2 × 10⁴⁶ erg. According to models by Chevalier & Fransson (1992), about 30% of the total energy input into the PWN, E_PWN, goes to the internal magnetic pressure in the PWN, while most of the remaining pulsar spin-down energy would drive the PWN expansion into the ejecta.

In terms of integrated radio luminosity, calculated as

$$\begin{equation} L_{\rm rad}(\nu)=4 \pi d\,^{2} \int _{\nu _{\rm min}}^{\nu _{\rm max}}\,S_{\nu }(\nu)\, \mathrm{d}\nu, \end{equation} \tag{ 10 }$$

if ν_min = 102 GHz and ν_max = 672 GHz bracket the frequency range in which the PWN is detected, S_ν ≈ 3 mJy leads to L_rad ≈ 5.4 × 10³³ erg s⁻¹. The derived L_rad is comparable with the limit of L_opt ⩽ 5 × 10³³ erg s⁻¹ given by Graves et al. (2005) for a compact source in the optical band from 290 to 965 nm at t = 6110 days. A similar limit has been placed on the 2–10 keV X-ray luminosity, L_X ⩽ 5.5 × 10³³ erg s⁻¹, using Chandra images (Park et al. 2004). Since these L_ν estimates are upper limits and free–free emission may be a substantial component of the radio luminosity (see Section 5.1), we can take L_rad ≈ 10³³ erg s⁻¹ as a realistic upper limit.

In the free-expansion regime (Chevalier & Fransson 1992), ∼1% of the pulsar power is emitted by the shock wave in the SN and additional emission from the pulsar nebula is expected. Given this, we can set $\dot{E}_{\rm spin} \approx 10^{35}$ erg s⁻¹ as an upper limit for the spin-down power of the pulsar. For a typical pulsar surface dipole magnetic field B_s ∼ 10¹² G (Manchester et al. 2005b), this spin-down luminosity corresponds to a pulsar period P ∼ 0.15 s and characteristic age $\tau _{c}= P/(2 \dot{P}) \sim 10^{5}$ years. Lower luminosities would imply lower dipole magnetic fields and/or longer pulsar periods. Such parameter ranges are plausible for the putative pulsar at the center of SN 1987A, as there is good evidence that many pulsars are born with a spin period not much different to their present period (Popov & Turolla 2012; Gotthelf et al. 2013).

Following Chevalier (1977), the velocity at the outer edge of the nebula can be defined as

$$\begin{equation} v_{\rm PWN}={6 \over 5} \left ({125 \over 132} \,\,{\dot{E}_{\rm spin}\, t \over \pi \rho _{ej} \, t^{3}} \right)^{{ { 1 \over 5} }}, \end{equation} \tag{ 11 }$$

assuming the PWN is freely expanding and has constant density ρ_ej. At t = t₂₅, setting 10⁻¹⁸ ≲ ρ_ej ≲ 10⁻¹⁹ g cm⁻³ in the central region of the SNR, as from the density model by Blinnikov et al. (2000; see Figure 21 in Fransson et al. 2013), the swept-up shell velocity becomes 260 ≲ v_PWN ≲ 410 km s⁻¹, which leads to R_PWN not greater than 005. Since this is well below the resolution of the ATCA and ALMA images presented here, the emission from a possible PWN would appear as a point source.

As mentioned in Section 5, a pulsar embedded in the SNR interior would emit ionizing radiation within the inner layers of the ejecta. While an X-ray pulsar has yet to be detected (Helder et al. 2013), illumination of the inner ejecta by X-ray flux has been reported by Larsson et al. (2011) though attributed to the reverse and reflected shocks, as well as to shocks propagating into the ER.

Given the early stages of a possible PWN, the ionized ejecta would be mainly due to the radiation from the various shocks. If the H_II region is assumed to be spherical (see Section 5), the related dispersion measure (DM) can be derived as $\,{\rm DM} (M_{\rm H\, \scriptsize{II}}, R_{\rm H\, \scriptsize{II}})\approx 3 \, M_{\rm H\, \scriptsize{II}}/ (m_{p} \, 4\pi \, R_{\rm H\, \scriptsize{II}}^{2})$. Neglecting clumping in the ejecta, for 0.1 M_☉ ≲ M_H ii ≲ 2.5 M_☉ and 0.01 ≲ R_H ii ≲ 0.10 pc (i.e., 005 ≲  R_H ii ≲ 042), the resulting DM is shown in Figure 12.

Zoom In Zoom Out Reset image size

We have presented a comprehensive morphological and spectral analysis of SNR 1987A based on both ATCA and ALMA data (Cycle 0), from 1.4 to 672 GHz. We have investigated the components of the SNR emission across the transition from radio to FIR, as the combination of non-thermal and thermal emission. A summary of our findings is as follows.

• 1.  
  To decouple the non-thermal emission from that originating from dust, the synchrotron component, as resolved with ATCA at 44 GHz, and the dust component, as imaged with ALMA at 672 GHz, were subtracted from the data sets at 94, 102, 213, 345, and 672 GHz. The images derived from the subtraction of the scaled model flux density at 672 GHz, highlight the ring-like synchrotron emission that originates over the ER. The images similarly obtained from the subtraction of the 44 GHz model, show residual emission mainly localized west of the SN site.
• 2.  
  An analysis of the emission distribution over the ER, in images from 44 to 345 GHz, shows a gradual decrease of the east-to-west asymmetry ratio with frequency. We attribute this to the shorter synchrotron lifetime at high frequencies. We estimate that, at frequencies higher than 213 GHz, the electrons might be unable to cross the eastern emission region.
• 3.  
  Across the transition from radio to FIR, the SED suggests additional emission components beside the synchrotron main component (S_ν∝ν^−0.73) and the thermal emission originating from dust grains at T ∼ 22 K. We argue that this excess emission could be due to a second flat-spectrum synchrotron component with −0.4 ≲ α ≲ − 0.1. This could imply the presence of a PWN originating from an embedded pulsar.
• 4.  
  Spectral index measurements across the SNR, from 44 to 345 GHz, show predominantly flat spectral indices, −0.4 ≲ α ≲ −0.1, in the western half of the remnant at frequencies above 102 GHz, while α ∼ 0 around the central region. From 102 to 213 GHz, −0.1 ≲ α ≲ 0 across the whole remnant. From 213 to 345 GHz, the steepening of the spectral indices over the eastern lobe might be indication of a local spectral break.
• 5.  
  Results from both the morphological and spectral analysis might suggest the presence of a PWN in the SNR interior, powered by a pulsar likely located at a westward offset from the SN position. In this scenario, a NS embedded in the unshocked ejecta might be emitting diffuse synchrotron emission $S_{\nu _P} \sim 3$ mJy above 102 GHz. The possible PWN would have magnetic field strength 1 ≲ B_PWN ≲ 7 mG. The upper limit of the integrated radio luminosity is derived as L_rad ∼ 5.4 × 10³³ erg s⁻¹, which leads to an upper limit of the pulsar spin-down power $\dot{E}_{\rm spin} \approx 10^{35}$ erg s⁻¹. While effects such as luminosity, beaming, scattering, and absorption might have prevented the detection of pulsed emission in the SNR for over 20 years (Manchester 1988, 2007), a renewed pulsar search with the Parkes telescope is currently in progress. Future observations with ALMA, ATCA and VLBI are needed to probe the SNR for a possible PWN.

We thank Philipp Podsiadlowski and Bruno Leibundgut for useful discussions; Richard McCray, Craig Wheeler, and Eli Dwek for feedback on the drafts. We thank Robert Kirshner and Peter Challis for providing the HST images. Figures 9 and 10 utilize the cube helix color scheme introduced by Green (2011). This paper makes use of the following ALMA data: ADS/JAO.ALMA #2011.0.00273.S (PI: Indebetouw). ALMA is a partnership of ESO (representing its member states), NSF (USA), and NINS (Japan), together with NRC (Canada) and NSC and ASIAA (Taiwan), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO, and NAOJ. The Australia Telescope Compact Array is part of the Australia Telescope, which is funded by the Commonwealth of Australia for operation as a National Facility managed by CSIRO. L.S.S. and B.M.G. acknowledge the support of CAASTRO, through project number CE110001020. R.A.C. acknowledges the support of NASA grant NNX12AF90G. M.L. acknowledges an ESO/Keele University studentship.