EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 690 (October 2024) A&A, 690 (2024) A161 Full HTML
Open Access
Issue
A&A
Volume 690, October 2024
Article Number A161
Number of page(s) 31
Section Catalogs and data
DOI https://doi.org/10.1051/0004-6361/202450730
Published online 03 October 2024

© The Authors 2024

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

The interstellar media (ISM) of galaxies are dotted with supernova remnants (SNRs), which are emission nebulae produced by supernova (SN) shocks propagating through the ambient ISM. Over 300 SNRs have been discovered in the Milky Way (Green 2019), and more than a thousand SNRs and SNR candidates have been found in external galaxies from multiwavelength surveys (see Vučetić et al. 2015; Long 2017, and references therein). SNRs are vital to our understanding of collisionless shocks (e.g., Raymond 1979; Chevalier et al. 1980; McKee & Hollenbach 1980; Marcowith et al. 2016), cosmic ray acceleration (Blandford & Eichler 1987; Blasi 2013; Bell et al. 2013), progenitor models and explosion physics of supernovae (SNe) (e.g., Chevalier 2005; Schaefer & Pagnotta 2012; Krause et al. 2008; Patnaude & Badenes 2017), the formation and survival of dust grains in SN ejecta (e.g., Dwek & Arendt 1992; Williams & Temim 2016; Sarangi et al. 2018), the properties of central compact objects (e.g., Dubner et al. 1998; Holland-Ashford et al. 2017; Katsuda et al. 2018), and the effects of shocks on interstellar clouds (e.g., White & Long 1991; Reach et al. 2006; Slane et al. 2015; Koo et al. 2020).

Extragalactic SNR surveys provide a unique perspective into the evolution of SNR shocks, their stellar progenitors, and impact on the surrounding ISM. They are complementary to Galactic SNRs, which undoubtedly offer the most detailed multiwavelength views of shocks (e.g., Milisavljevic et al. 2024), but are also affected by distance uncertainties, variable line-of-sight extinction, and source confusion at low latitudes where most objects are located (e.g., Green 2005; Reach et al. 2006; Green et al. 2008). Examples of well-studied extragalactic SNR populations include those in the Magellanic Clouds (e.g., Mathewson & Clarke 1973; Long et al. 1981; Maggi et al. 2016; Maggi et al. 2019), M31 (Magnier et al. 1995; Lee & Lee 2014a; Galvin & Filipovic 2014), M33 (Gordon et al. 1999; Long et al. 2010; Lee & Lee 2014b; Long et al. 2018; White et al. 2019), M83 (e.g., Dopita et al. 2010; Blair et al. 2014; Winkler et al. 2017; Long et al. 2022), NGC 6946 (e.g., Matonick & Fesen 1997; Lacey et al. 1997; Long et al. 2019, 2020), M51 (Winkler et al. 2021), and the Sculptor Group galaxies (e.g., Blair & Long 1997; Pannuti et al. 2000, 2002; Kopsacheili et al. 2024).

With extragalactic surveys, one can observe the distribution of SNRs of different ages across entire galaxies for face-on galaxies at the known fixed distances of the galaxies. This leads to reliable estimates of SNR sizes and luminosities at different wavelengths, from which one can infer a variety of statistical information about the SNR population, such as the distribution of evolutionary stages, visibility times, ambient densities, and magnetic field properties (Gordon et al. 1999; Chomiuk & Wilcots 2009; Thompson et al. 2009; Badenes et al. 2010; Long et al. 2010; Asvarov 2014; Sarbadhicary et al. 2017; Elwood et al. 2019). The spatial correlation of these SNRs with ISM and stellar population surveys of nearby galaxies have led to unique constraints on the progenitor age distribution of SNe (e.g., Badenes et al. 2009; Maoz & Badenes 2010; Jennings et al. 2014; Williams et al. 2014; Díaz-Rodríguez et al. 2018; Auchettl et al. 2019; Williams et al. 2019; Koplitz et al. 2021), highly relevant to the many open questions about progenitor models of Type Ia and core-collapse SNe (see Maoz et al. 2014; Smartt 2015, for detailed discussion).

The discovery and identification of these extragalactic SNRs have been primarily at optical wavelengths from observations of forbidden emission lines. The most widely used diagnostic to distinguish SNRs from photoionized nebulae has been the [S II]/Hα>0.4 criterion (e.g., Mathewson & Clarke 1973; Parker et al. 1979; Dodorico et al. 1980), since low-ionization species like S+ are more abundant in the extended recombination zones of radiative shocks compared to H II regions (Raymond 1979; Dopita et al. 2005; Allen et al. 2008). This contrasts with the Milky Way, where SNRs are primarily identified with radio continuum observations that can pierce through the midplane dust, and the power-law index of their synchrotron spectra can be measured (Green 2014; Anderson et al. 2017; Dokara et al. 2021; Cotton et al. 2024). The Magellanic Clouds also have a significant number of SNRs identified at X-ray (for young SNe) and radio wavelengths (e.g., Maggi et al. 2016; Maggi et al. 2019; Kavanagh et al. 2022; Bozzetto et al. 2023; Zangrandi et al. 2024). Beyond the Local Group, however, radio and X-ray observations have yet to reach sufficient spatial resolution and sensitivity to independently discover large numbers of SNRs like optical surveys (Russell et al. 2020).

Most optical SNR surveys have relied on narrowband filter imaging, which, despite its success in discovering many SNRs, is not without its limitations. A major issue is the continuum subtraction, which is done by scaling the flux from a broadband image or an off-band narrowband image. This was often suboptimal in surveys, due to variable seeing conditions and color properties of stars in the fields (e.g., Dopita et al. 2010; Blair et al. 2014), which particularly affected the recovery of faint, low surface brightness SNRs (Long et al. 2022). Typically the filters blended adjoining lines, such as the [N II]λ6549,6583 and Hα, which compromised the [S II]/Hα ratio (Blair et al. 2014). As a result, narrowband imaging and band ratios alone were often not sufficient to confirm objects as SNRs, and additional medium-or high-resolution follow-up spectroscopy of individual objects was necessary to provide additional line ratios and kinematic indicators (Winkler et al. 2017; Long et al. 2018).

The proliferation of wide-field integral field unit (IFU) spectroscopy with instruments such as Multi Unit Spectroscopic Explorer (MUSE) (Bacon et al. 2010) and Spectromètre Imageur à Transformée de Fourier pour l’Etude en Long et en Large de raies d’Emission (SITELLE) (Drissen et al. 2019) has significantly improved our ability to survey extragalactic nebulae with comparable efficiency as narrowband imaging, while at the same time circumventing many of its shortcomings. The ability to obtain a spectrum per pixel at high sensitivity and medium spectral resolution (R > 103) has led to more accurate modeling and subtraction of the underlying stellar continuum, and separation of adjoining emission lines. Integral field spectrograph has also opened up access to additional line diagnostics for identifying SNRs, such as [O I] and [O II] (Fesen et al. 1985); linear combinations of S, N, and O forbidden lines (Kopsacheili et al. 2020); and the widths of forbidden emission lines, which are expected to be broader than in H II regions due to the presence of high-velocity shocked material (Points et al. 2019). The promise of IFU-based SNR searches have already been demonstrated in individual galaxies such as NGC 300 (Roth et al. 2018), NGC 3344 (Moumen et al. 2019), NGC 4030 (Cid Fernandes et al. 2021), M83 (Long et al. 2022), and NGC 4214 (Vicens-Mouret et al. 2023), with more recent work building toward cataloging thousands of shock-ionized nebulae across larger galaxy samples (Congiu et al. 2023, hereafter C23).

In this paper we push the boundaries of extragalactic SNR studies even further, to a sample of 19 galaxies in the range of 5–20 Mpc, as part of the Physics at High Angular resolution in Nearby GalaxieS (PHANGS)-MUSE survey (Emsellem et al. 2022). We have already characterized and constructed catalogs of likely H II regions and diffuse ionized gas (DIG) in these galaxies (see, e.g., Emsellem et al. 2022; Belfiore et al. 2022; Santoro et al. 2022; Groves et al. 2023; Congiu et al. 2023; Egorov et al. 2023), which is critical for identifying SNRs. Although at these distances the SNRs are unresolved or barely resolved (average spatial resolution of 67 pc), this is the largest single survey of extragalactic SNRs ever conducted, with a diverse sample of host galaxies that span a range of stellar masses and star formation rates. These galaxies have substantial multiwavelength data at ~1″resolution or better from facilities including, but not limited to, Hubble Space Telescope (HST), James Webb Space Telescope (JWST), Atacama Large Millimeter/submillimeter Array (ALMA), Karl G. Jansky Very Large Array (VLA), AstroSAT, and Chandra X-ray Observatory (CXO). This multiwavelength synthesis provides the most complete observational census of the multiphase ISM and the multigeneration stellar population that are associated with SN explosions.

The focus of the current paper is to introduce the catalog, the selection procedure adopted, and the basic properties of the PHANGS-MUSE SNRs. Subsequent papers will investigate the correlation of these SNRs with stellar population and ISM maps of the galaxies. The paper is organized as follows. In Sec. 2 we introduce the PHANGS-MUSE dataset for our sample. Then in Sec. 3 we explore methods for identifying SNRs in these galaxies and describe how we constructed our catalogs. The resulting catalogs and SNRs properties are presented in Sec. 4. Further discussion and comparison with literature works are included in Sec. 5. Finally, we summarize the main conclusions in Sec. 6.

Table 1

Properties of PHANGS-MUSE galaxies used in this work.

2 Data

We selected SNRs using data from the PHANGS-MUSE survey of 19 nearby (≲20 Mpc) galaxies, conducted using the Multi Unit Spectroscopic Explorer (MUSE) instrument mounted on the Very Large Telescope (VLT) in Chile (Bacon et al. 2010). MUSE provides ~arcsecond-resolution integral field spectroscopy in the wavelength range of 4800 to 9300 Aring; at a spectral resolution of R ~ 3000. Full details of the sample, observations, reductions, and data products are found in Emsellem et al. (2022). These 19 galaxies (see Table 1) populate the star-forming main sequence of galaxies, with stellar masses in the range of 109.4 − 1011 M⊙. The pixel scale is 0.2“/pixel and the typical seeing in R-band is 0.8” (Emsellem et al. 2022), corresponding to an average spatial resolution of 67 pc at the distance of our targets (Anand et al. 2021; Scheuermann et al. 2022). Our ~100 pc resolution is larger than the typical size of SNRs (Long et al. 2010; Tüllmann et al. 2011; Asvarov 2014), but is well suited to spatially distinguish SNRs from nearby unrelated H II regions and surrounding DIG. To limit the effects of spatial blending between SNRs and surrounding line emitting regions we used the native resolution data, which provide the best possible spatial resolution at the expense of a lightly varying point spread function (PSF) between fields. As many of the regions are found to be blended with luminous H II regions (see Sec. 5.3), our source identification is limited more by seeing than by sensitivity. Using smoothed data (e.g., the “copt” data products), we obtain a sample of objects that is 60% smaller.

The main PHANGS-MUSE data products we used are the total line intensity and kinematic (velocity and velocity dispersion) maps for the set of strong lines discussed below. A detailed description of the datasets, including the data analysis pipeline (DAP) used to extract these parameters, is given in Emsellem et al. (2022), and here we provide a summary. Spectral fitting is performed by using the penalized pixel-fitting code pPXF (Cappellari & Emsellem 2004; Cappellari 2017). The stellar continuum in the spectra was identified and removed using stellar population templates from the E-MILES simple stellar population models (Vazdekis et al. 2016). Emissions lines were fitted with Gaussian templates to recover line fluxes as well as line kinematics (velocity and velocity dispersion). We note that at the MUSE spectral resolution, we generally find that the line emission associated with SNRs is well described by single Gaussian fits. In pPXF, we performed fitting for 23 lines along with their corresponding gas kinematics. The kinematics of some emission lines were grouped together and fitted simultaneously, such as Balmer lines, and both low- and high-ionization lines. In this work we used [S II]λ6716, [S II]λ6730, Hα, [O I]λ6300, [O III]λ5006 fluxes and [S II]λ6716 velocity dispersion. We note that Hα and [S II], lines commonly used for SNR identification, are fit independently. Variations in seeing due to the differences in wavelength between these line maps are negligible (<0.01″; Emsellem et al. 2022). For fluxes, a signal-to-noise ratio (S/N) cut of S/N ≥ 5 is applied, while for the velocity dispersion, a higher cut S/N ≥ 20 is used (Egorov et al. 2023), hereafter E23. These S/N cuts are applied before identifying SNRs in all maps.

We also made use of the Groves et al. (2023) and Santoro et al. (2022) catalog of H II regions to compare with our SNR properties. Groves et al. (2023) identified 23 244 nebulae as H II regions, imposing on a S/N cut of five for all emission lines and using Baldwin, Phillips & Terlevich (BPT) diagnostics ([O III]λHβ versus [O I]/Hα and [O III]/Hβ versus [S II]/Hα) for classification (Baldwin et al. 1981; Kewley et al. 2001; Kauffmann et al. 2003).

Before selecting objects from the images, we also masked and excluded regions that suffer from contamination by foreground stars (provided as star masks by Emsellem et al. 2022), as well as regions where the ionized gas may have shocks from processes other than SNRs. For example, regions close to the galactic centers and bars can exhibit line ratios and physical conditions consistent with widespread shocks caused by largescale supersonic inflows or outflows of gas (Athanassoula 1992; López-Cobá et al. 2022). In addition, some of the galaxies in our sample (e.g., NGC 1365) contain active galactic nucleus (AGN), which also strongly affects the line ratios in the central regions. As a result, we used the environments defined in Querejeta et al. (2021) to mask the centers of all galaxies to avoid shocked emissions from potential AGN, and the bar region for galaxies showing increased gas turbulence over large areas. For bar regions that were without increased shocked regions, we did not mask them off (see Table A.1 in Appendix A for a complete list of which environmental mask is applied for each galaxy).

3 Methods

3.1 Parent sample identification

Distinguishing SNR from other nebulae is a long-standing problem in nearby galaxies (Sabbadin et al. 1977; Rhea et al. 2023). At ~100 pc resolution, most ionized nebulae appear to be luminous and compact sources. However, with multiple emission-line diagnostics, it is possible to find and classify SNRs.

The [S II]/Hα line ratio is historically the most commonly used diagnostic (Blair et al. 1981; Blair & Long 2004) to separate SNRs from H II regions. This diagnostic works because the shocked region is richer in warm [S +] with respect to H II regions, which translate into a much higher [S II]/Hα line ratio (Allen et al. 2008). A similar situation applies to [O I]. Indeed, the [O I]/Hα line ratio has been predicted to be a better diagnostic for distinguishing SNRs and H II regions than the more commonly used [S II]/Hα line ratio (Kopsacheili et al. 2020). However, this emission line is normally weaker than [S II], and in the very local universe, it can be affected by blending with sky line airglow emission. Thanks to the redshift of our sample, the intrinsic [O I] emission is shifted from sky emission, and the depth of our MUSE data allows us to obtain sufficient S/N in our [O I] fluxes to consider it in our selection of SNRs.

Another prominent feature of SNRs is their broad emissionline profiles arising from the presence of fast radiative shocks of ~200 km s−1 (Long 2017), which significantly exceeds the typical ~20 km s−1 velocity dispersion of H II regions (Winkler et al. 2015). This kinematic information provides a promising alternate identification method (Points et al. 2019; Long et al. 2022; Winkler et al. 2023; Egorov et al. 2023). Finally, by combining multiple line ratios (such as in BPT diagrams) it is possible to determine the dominant ionization mechanism, whether it is photoionization-dominated or shock-dominated (Congiu et al. 2023; Makarenko et al. 2023).

We explore the use of five different criteria in this paper. These are the [S II]/Hα line ratio, the [O I]/Hα line ratio, the [S II]λ6716 velocity dispersion (σ) and two BPT diagnostics: [O I]-[O III] and [S II]-[O III]. We identify objects using each of these five different criteria separately and combine these catalogs to construct our “parent sample” of objects, which we later classify as either SNRs or SNR candidates. All these criteria are explained below in more detail.

The Python package astrodendro1 was used to select objects from line ratio maps and line kinematic maps. astrodendro is a tool used for calculating dendrograms (Rosolowsky et al. 2008), which visually represent the hierarchical clustering of data points, particularly in astronomical data analysis. We use astrodendro to identify regions around local maxima (above a threshold min_value), selected based on their contrast with the surrounding pixels (min_delta) and required to have a minimum associated area that indicates the line ratio to be coherent over about the area of a PSF. We require regions used in further analysis to contain at least 10 pixels, which is the approximate area of the PSF. This ensures that all considered regions are detected over the whole area of the PSF and do not reflect an artifact or ‘hot pixel.’ We further require all regions to have a diameter of less than 200 pc, as otherwise they would be unphysically large compared to literature catalogs of SNRs (Asvarov 2014).

In this analysis, we do not aim to produce a complete catalog of SNRs in these galaxies but rather try to provide a sample that is as pure as possible. Many of the choices made in the following sections are designed to minimize the selection of spurious objects, such that all of the objects in our final catalog are credible SNRs and SNR candidates, as judged in a subsample of objects by visual inspection. The well-studied nearby galaxy M83 is used as a benchmark to validate our method in Sec. 3.3, and we compare with other approaches to identify ‘shocked’ nebulae in Sec. 5.4.

3.1.1 Selection via [SII]/Hα and [OI]/Hα

[S II] emission arises not only from shocks in SNRs but also through photoionization, such as in H II regions or the DIG. We initially tried to identify SNRs directly from the line ratio maps, and ran astrodendro, identifying peaks where the min_value was set to a threshold of 0.4 in [S II]/Hα and a threshold of 0.017 in [O I]/Hα, based on shock models (Kopsacheili et al. 2020). We found 36 826 objects in the [S II]/Hα maps and 15 891 objects in the [O I]/Hα maps. This implies nearly 1000–2000 SNRs per galaxy, which appears unphysically high (although see Sec. 5.1) considering that even galaxies like M83, NGC 6946, and M51, which have similar star formation rates as the PHANGS sample and are about a factor of 2–4 closer, only have on the order of a few hundred SNRs that have been spectroscopically identified (Winkler et al. 2021; Long et al. 2022).

Visual inspection reveals that the majority of these objects found in the PHANGS-MUSE sample are at the boundaries of H II regions and at low Hα surface brightness, implying they are either artifacts due to noise or associated with DIG. Figure 1 shows cutout images of four areas in NGC 628 and the regions identified within them. These images compare the region boundaries identified in the Groves et al. (2023) nebular catalog, where objects were identified via their Hα morphology, with the objects identified by astrodendro using the [S II]/Hα and [O I]/Hα ratio maps. The first column illustrates objects identified by Groves et al. (2023) and classified as H II regions (red lines) or unclassified regions2 (blue lines), overlaid on the Hα image. The second to fifth columns show peaks identified by astrodendro overlaid on the corresponding line-ratio or velocity dispersion map. The last column indicates the [S II] flux in log scale. It is apparent that the 0.4 threshold in [S II]/Hα selects not only SNRs but also spurious regions, which in shallower imaging would not have been detected (see also Winkler et al. 2023). After accounting for the contributions of photoionized emission to our maps (see below), our region selection becomes much cleaner, and can often be more directly associated with the ‘unclassified’ nebular objects.

At low surface brightness, most star-forming galaxies show widespread line emission arising from DIG, which makes up ~50% of the integrated Hα emission of galaxies (Zurita et al. 2000). DIG typically shows elevated levels of [S II]/Hα, [N II]/Hα and [O I]/Hα (Haffner et al. 2009), and can thus be confused for low surface brightness SNRs (e.g., Blair & Long 1997; Long et al. 2022). As shown in Belfiore et al. (2022), both the [S II]/Hα and [O I]/Hα line ratios in the DIG can be modeled by leaking radiation from H II regions, suggesting the DIG is principally powered by leakage of radiation from H II regions with a typical mean free path for the ionizing radiation of ~2 kpc.

At high surface brightness, SNRs might be blended in projection with giant H II regions (see the third row in Fig. 1). if a SNR is seen in projection against an H II region, there is the possibility that the SNR will be outshone by the H II region, and thus be challenging to distinguish.

To exclude the contribution from DIG and H II regions to the [S II]/Hα and [O I]/Hα line ratio maps, we build on the results of Belfiore et al. (2022), who carried out a detailed characterization of DIG in the PHANGS-MUSE galaxies. Belfiore et al. (2022) spatially binned the line emission to study the faint DIG, and found strong anti-correlations between the Hα surface brightness and common diagnostic line ratios, which they show can be well fitted by photoionization models as a function of ionization parameter (in this context, a proxy for Hα intensity).

We examine the correlations between line ratio and Hα surface brightness on a pixel-by-pixel level in each of the 19 galaxies (Figs. 2 and 3), in order to search for regions exhibiting anomalously high line ratios at fixed Hα surface brightness. We expect that most pixels are dominated by photoionization, as SNRs typically only contribute ~5% of the Hα flux of an entire galaxy (Vučetić et al. 2015). Overall, we find very good agreement with the binned relations found in Belfiore et al. (2022) for the [S II]/Hα line ratios, as our data is sufficiently deep that the [S II] line is detected in a large fraction of the pixels. For the fainter [O I] line, it is apparent that we are limited by sigma-clipping at low Hα surface brightness, where [O I] is less frequently detected, and the majority of our pixels reflect ~3σ outliers in the binned distribution. At high surface brightness, however, the relations agree. As is apparent in Figs. 2 and 3, we expect that in many ways our pixel-based fits are more conservative, as they generally result in steeper slopes that do not correctly model the DIG line ratios at low surface brightness. We emphasize here that as a result of our sensitivity limits, we expect that our sample is both incomplete and biased toward selecting SNR with high line ratios.

Given the observed trends, we derive linear scaling relations to characterize the expected line ratios for log [O I]/Hα (taking all galaxies jointly) and log [S II]/Hα (taking each galaxy separately) as a function of log Hα surface brightness. We used a linear regression model to fit these pixels to get a slope and an intercept, provided in Figs. 2 and 3. As we find that the [O I]/Hα versus Hα surface brightness shows uniform trends for all 19 galaxies, we fit all galaxies together, although we caution that these fits are not necessarily representative of the full DIG due to the S/N limitations on detecting [O I]. Pixels belonging to SNRs are apparent in some cases, as shown in the horizontal stripes above the fitted relation. We are generally searching for spatially coherent objects that are outliers with high line ratios, and we tailor our object selection to identify regions based on the pixel distribution and our own visual inspection.

Using these relations, we construct line-ratio residual maps as ([S II]Hα)residual=log⁡([S II]Hα)obs−log⁡([S II]Hα)fit,$\[\left(\frac{\left[\mathrm{S}~{\mathrm{\small{II}}}\right]}{H \alpha}\right)_{\text {residual}}=\log \left(\frac{\left[\mathrm{S}~{\mathrm{\small{II}}}\right]}{H \alpha}\right)_{\mathrm{obs}}-\log \left(\frac{\left[\mathrm{S}~{\mathrm{\small{II}}}\right]}{H \alpha}\right)_{\text {fit}},\]$(1)

where the subscripts obs denotes the observed line flux and fit is calculated as a function of the Hα surface brightness from the linear fit equations. ([O I]/Hα)residual is calculated in a similar way.

In our residual maps, we measure the standard deviation (σ) and use astrodendro to identify peaks assuming a min_value of 2σ in [S II]/Hαresidual and 1σ in [O I]/Hαresidual. min_delta is set to 0.1 dex to avoid selecting the second maximum peak that is close to the highest peak. We chose these thresholds to guarantee we did not miss objects (as judged by visual inspection). We identify 880 objects selected from the [S II]/Hα line ratio residual, and 1368 objects selected from the [O I]/Hα line ratio residual.

To judge the impact of any sigma-clipping bias at low Hα surface brightness, which principally affects the [O I]/Hα line ratios, we repeated this analysis using the fit and σ determined from the DIG binned data, and identifying 3σ pixel-scale outliers. This resulted in a very similar number of objects selected by their [S II]/Hα residuals, and only slightly fewer objects selected by their [O I]/Hα residuals, indicating that our fine-tuning of the astrodendro parameter min_value results in a roughly equivalent selection of objects.

thumbnail Fig. 1

Zoomed-in image on four regions, showing (from left to right) log Hα flux (10−20 cm−2ergs−1), [S II]/Hα, [S II]/Hα residual, [O I]/Hα residual, and [S II] velocity dispersion maps and log [S II] flux. Objects identified by our method as SNRs and SNR candidates are marked with black solid ellipses. Objects identified in the Groves et al. (2023) nebular catalog as H II regions are enclosed by red solid lines, while regions unclassified in that catalog are outlined with blue solid lines. Region contours in the first column are from Groves et al. (2023). The second column shows region contours identified by astrodendro without any photoionization correction. The remaining columns show regions identified by astrodendro after using our final selection methods. It is apparent that many of the regions in the second column are spurious.
thumbnail Fig. 2

Log([S II]/Hα) as a function of logarithm Hα surface brightness for all pixels in each galaxy. The red solid line is the fitted correlation for each galaxy, with slope and intercept values given in the legend of each panel. The black solid line is representative of the DIG relation, as measured between 0.40 Rmax (maximum radial coverage) and 0.60 Rmax for each galaxy from Belfiore et al. (2022). The gray area indicates the 3σ range of this relation. SNRs typically lie above the fitted lines and show a looping structure toward the upper right. The horizontal black dashed line indicates the typical value of [S II]/Hα = 0.4 used to select SNRs. Many low Hα surface brightness pixels associated with the DIG have [S II]/Hα ≥ 0.4.
thumbnail Fig. 3

log([O I]/Hα) as a function of log Hα surface brightness for all pixels in all 19 galaxies. The red solid line is the fitted global correlation for all galaxies, with slope and intercept provided in the lower right corner of this plot. The black solid line is representative of the DIG relation, as measured between 0.40 Rmax (maximum radial coverage) and 0.60 Rmax for each galaxy from Belfiore et al. (2022). The gray area indicates the 3σ range of this relation. The black dashed line indicates the theoretical value of [O I]/Hα = 0.017 to select SNRs (Kopsacheili et al. 2020). The horizontal blue dotted line indicates the empirical value of [O I]/Hα = 0.1 used select SNRs. For NGC 1672, the individual pixels lying above the majority come from an SNe that just happened several days before our observation.

3.1.2 Selection via velocity dispersion maps

As mentioned earlier, the spectral resolution of our MUSE data also enables the selection of SNRs from broadened emission lines in their spectra. In particular, we use the velocity dispersion of the [S II] line, corrected for instrumental broadening, to select our SNRs. We note that the kinematics of [S II]λ6716 has been tied in the fit with other low-ionization lines to further improve the reliability of the measurement; however, all of these (e.g., [S II]λ6730, [N II]λλ6548,84, [O I]λλ6300,64; see Emsellem et al. 2022) should be strongly associated with shocks. While the Hα velocity dispersion can, in principle, be utilized, we find that the [S II] line kinematics can be more directly and uniquely associated with SNR shocks. On average, the line kinematics typically agree between both tracers within 5 km s−1.

An S/N cut of 20 for [S II]λ6716 is applied to ensure errors in the measured velocity dispersion remain below 10% when deconvolved from the ~50 km s−1 instrumental velocity dispersion (Bacon et al. 2017; Egorov et al. 2023). With smoothing and binning, higher S/N is achieved across the map, but SNRs become blended with H II regions and DIG, and it is thus more challenging to identify by their kinematic signatures. At low S/N, DIG might be hard to distinguish kinematically from SNRs, as it also produces broadened line emission due to its turbulent nature and increased scale height (Haffner et al. 2009; Levy et al. 2019). After masking galaxy centers, we do not see evidence that we are biased toward selecting DIG emission with this method as our S/N cut on [S II] has effectively masked most of the DIG; however, with sufficiently deep data it remains challenging to distinguish these components (Winkler et al. 2023).

We take a conservative approach and set our astrodendro threshold to select peaks with velocity dispersion larger than 50 km s−1 (full width at half maximum (FWHM) = 118 km s−1) for our parent sample. We visually tested different thresholds and then chose this value that optimized the number of detected SNR candidates without including too many spurious detections. A min_delta of 10 km s−1 is used to make sure local peaks are identified properly. A total of 352 objects are identified using this technique.

3.1.3 Selection via BPT distances

In addition to the aforementioned criteria of [S II], [O I] and velocity dispersion that are traditionally used for selecting SNRs, we also explore the properties, and the potential for selecting SNRs using BPT diagrams, widely used to assess the line emission mechanism in galaxies (Baldwin et al. 1981; Veilleux & Osterbrock 1987). Here we did not consider the [O III]/Hβ versus [N II]/Hα diagram as [N II]/Hα varies more with metallicity (see Sec. 4). On the ~100 pc scale of an individual nebula (top row of Fig. 4), regions that are dominated by shock heating lie above the extreme starburst line (Kewley et al. 2001), while regions dominated by photoionization from OB stars lie below this line (Congiu et al. 2023). These diagrams make use of the same [S II]/Hα and [O I]/Hα line ratios as we describe in Sect. 3.1.1; however, from the shape of the extreme starburst line it is apparent that a single fixed threshold may not be appropriate in all physical conditions.

For this reason, we select pixels as shock-heated versus photoionized based on their distances from the extreme starburst line in the [O III]/Hβ versus [O I]/Hα and [O III]/Hβ versus [S II]/Hα BPT diagrams (Kewley et al. 2006), as shown in Fig. 4. Here “distance” refers to the Euclidean distance (square root of the sum of vertical and horizontal separation squared) to the closest point on the extreme starburst line in both BPT diagrams. While the DIG also populates a similar parameter space in such diagrams (Zhang et al. 2017), the [O III] emission is particularly faint and is effectively excluded with our S/N>5 cut.

We construct maps of these two BPT distances, and run astrodendro on each, with a min_delta of 0.05 dex to avoid selecting the local maxima. We use a minimum distance of ±0.2 dex to classify pixels as shocks (positive distance) versus photoionized (negative distance), where the chosen threshold of 0.2 is roughly the 3σ uncertainty in the ratios measured for all pixels with S/N≥5. Figure 4 (bottom row) shows the integrated line ratios for all objects in our final catalog. SNRs do preferentially populate the shock-heated regime, and those SNRs found to be populating the photionized regime have typically been selected via one of our other criteria, or are blended with H II regions(see Sec. 5.3). A total of 635 objects are identified using the BPT distance in the [O III]/Hβ versus [O I]/Hα diagram, and a total of 372 objects are identified using the BPT distance in the [O III]/Hβ versus [S II]/Hα diagram.

3.1.4 Final parent sample construction and integrated line flux measurements

We employ five criteria to select objects, and using the central coordinates determined through astrodendro, we perform a crossmatch within a 1” radius across all objects. This process yields a parent sample consisting of 2233 objects selected by at least one of our criteria. In cases where the astrodendro centers differ across the five maps, we adopt the center of the first criterion that identified the object as an SNR or SNR candidate (following the order in which we present the criteria in this paper). A summary of how these different criteria each contribute to the total sample is presented in Sec. 3.2 and revisited in Sec. 5.2.

Figure 5 shows a map of the parent sample identified within NGC 628 using each criterion, and Fig. 6 shows a zoomed-in SNR of NGC 628, demonstrating the distinct range of values of each criterion at the SNR location compared to the surrounding region. For all 19 PHANGS-MUSE galaxies, the locations of 2233 objects in the parent sample are shown in Fig. 7 and Appendix B. We find on average 118 objects per galaxy, which are located throughout the disks of each galaxy.

For this parent sample we measure line fluxes, velocities, and velocity dispersions from the integrated spectra of these objects (see Table C.1 in Appendix C) using the PHANGS-MUSE data analysis pipeline or DAP (described in Emsellem et al. 2022). We expect most SNRs will not be not resolved in our galaxies, so we use a fixed radius of 50 pc for all objects, comparable to or slightly larger than the PSF in all galaxies, and measure emission-line fluxes and kinematics integrated over this area (see Emsellem et al. 2022 and Groves et al. 2023 for more details on how emission-line fluxes and gas dynamics are measured using the DAP).

It is important to note a caveat in the interpretation of these emission-line fluxes. The integrated line flux ratios [S II]/Hα and [O I]/Hα that we measure from the integrated spectra here are different from the [S II]/Hα and [O I]/Hα residual maps that we use to identify objects in Sect. 3.1.1. While we use astrodendro to identify SNRs from line flux ratios residual maps, when running the DAP we are unable to apply such a correction. Corrections for background photoionized emission have not been extensively undertaken in previous research, so we believe the line measurements in our catalog are directly comparable to literature results. We emphasize that these line fluxes may be biased by the DIG background, or any H II region they overlap with. However, disentangling these background effects is challenging as they do not exhibit a smooth morphology. As can be seen in Fig. 1, when the identified object overlaps with a H II region, the [S II]/Hα ratio is less distinct from the environment; however, after our photoionization correction, the [S II]/Hα ratio becomes more prominent. Hereafter, line fluxes for our SNRs and SNR candidates all refer to integrated flux without background subtraction.

thumbnail Fig. 4

BPT diagrams with demarcation line from Kewley et al. (2006) with an indication of BPT distances as defined as the distance to the red solid extreme star-burst lines. The plus sign (+) indicates positive distance where shock-heating dominates, while the minus sign (-) indicates negative distance where photoionization dominates. The top left and top right plots show the distribution density of nebular regions in these 19 galaxies (Groves et al. 2023). In the bottom left and bottom right plots, SNRs and SNR candidates are indicated by red and blue dots, respectively. Nebular regions tend to occupy the photoionized region, while SNRs are concentrated in the shocked areas.

3.2 Parent sample classification

We have used five criteria in Sect. 3.1 to identify objects for our parent sample. While the [S II]/Hα line ratio has most commonly been used in the literature to identify SNRs (e.g., Long et al. 2022), all of the diagnostics we employ are valid indicators of SNR shocks. To increase our confidence in our identification of these objects as true SNRs, we consider how to combine our different criteria and the integrated line fluxes to identify a robust subsample of objects that we can classify as SNRs.

Table 2 presents a breakdown of how many objects are identified by each criterion. We identify the largest number of objects using the [O I]/Hα line ratio maps (1368), and the smallest number of objects using the [S II]λ6716 velocity dispersion (352). Among our five criteria, 1476 objects were only selected by only one criterion, 757 objects were selected by at least two criteria, 365 objects were selected by at least three criteria, 184 objects were selected by at least four criteria and 68 objects were selected by all five criteria. Given the range of characteristics apparent in our parent sample, we investigate using different pairs of criteria, as this provides a secondary confirmation of our SNR classification.

Figure 8 (left) presents a visualization of what fraction of objects are matched pairwise between the five criteria. Clearly, the [O I]/Hα residual and [S II]/Hα residuals produce the most overlap (~20% of the parent sample), but are still not perfectly in agreement, and also show non-negligible ~10% overlap with the remaining indicators. We also note that a large number of objects (414) do not match between any two criteria, but exhibit high line ratios in their integrated line fluxes ([O I]/Hα>0.1 and [S II]/Hα>0.4), such that in previous work, they would have resulted in clear SNR classifications.

Based on this analysis of the 2233 objects from our parent sample, we define two subsamples. The “SNR candidates” are those objects only selected by one criterion, and either [S II]/Hα<0.4 or [O I]/Hα<0.1. The rest of the objects we classify as “SNRs”, which means that at least two of the five methods identify the object, or both the integrated [S II]/Hα and [O I]/Hα are high enough that the line emission is best explained by shock excitation. Figure 8 (right) presents a visualization of how objects are matched pairwise for just the SNR sample. [O I]/Hα residual and [S II]/Hα residuals select about half of the objects, but other pairs of diagnostics are also contributing significantly (10–20%) to the SNR sample.

We present spectra for four representative SNRs and one typical H II region in Fig. 9 in order to show the typical sensitivity and spectral features that we observe. In all panels, the strong sky lines at 6300Å and 6363Å are masked. The top two panels show the spectra for two SNRs with high [S II]/Hα > 0.8, the middle two panels for two SNRs with high-velocity dispersion > 80 km s−1, and the bottom panel for a representative H II region. The broadened emission lines are apparent by eye for the two high-velocity dispersion SNRs compared to the H II region, as is the increased [S II]/Hα ratio in the top two SNR.

For the SNRs and SNR candidates, when we place them in the [O I]/Hα v.s [S II]/Hα parameter space, we can clearly see that (by design) the majority of SNRs lie above [O I]/Hα>0.1 and [S II]/Hα>0.4 (see Fig. 10a) while the SNR candidates distribution overlaps with the H II regions (see Fig. 10b) making it difficult to distinguish SNR candidates from H II regions. The distinct boundary at 0.1 for [O I]/Hα in H II regions arises from the definition of these objects in Groves et al. (2023) using BPT diagram cuts.

Using the region footprints identified by Groves et al. (2023), we measure how many of the objects in our parent sample fall in projection on top of nebulae. In total, ~35% of our objects overlap with nebulae classified as H II regions. ~30% overlap with unclassified nebulae (typically labeled ‘unclassified’ as they are found to be inconsistent with photoionization based on diagnostic line ratios), indicating that these objects are likely also identifiable by their Hα morphology. The remaining ~35% of the parent sample is not included in the nebular catalog at all, and half of these are classified as SNRs, indicating how many objects would be missed if only line flux intensity is used for object selection. Within our catalog, we label those objects that overlap with H II regions as ‘blended’, while those that do not overlap (or overlap with unclassified regions in the Nebular catalog) we identify as ‘isolated’.

We provide integrated properties for all objects in our parent sample as a machine readable catalog with the list of columns provided in Table 3. The catalog includes the host galaxy name and distance, object classification, object distance to the galaxy center in units of effective radius, environment mask that SNR lies on, coordinates, which criterion selected the object, emission-line fluxes, gas kinematics, isolated flag, [S II]/Hα value, and Hα luminosity.

thumbnail Fig. 5

Zoomed-in map (3 kpc × 3 kpc) of the SNR population recovered in one of the galaxies in our sample, NGC 628, using our five selection criteria (described in Sec. 3). Panels from top left to bottom right show the log Hα emission, [O I]/Hα residual map, [S II]/Hα residual map, [S II]λ6716 velocity dispersion, BPT distance (OI-OIII) map, and the BPT distance (SII-OIII) map. Objects selected by each individual criterion are marked with black ellipses in the corresponding subplot and the locations of objects classified as SNR are marked with red boxes. Pixels with low S/N (≤5), and those corresponding to stellar sources, and the central region of the galaxy, have been masked (see text for details).
Table 2

Number of identified SNRs and SNR candidates in 19 galaxies.

thumbnail Fig. 6

Example of what one SNR (in NGC 628) looks like when using the five criteria. This SNR is selected in the [O I]/Hα residual map, [S II]/Hα residual map, [S II]λ6716 velocity dispersion, BPT distance (OI-OIII), and BPT distance (SII-OIII) maps. The selection is marked with a black ellipse.

3.3 Validation of methods: Comparison with supernova remnants in M83

Our SNR identification method differs somewhat from past SNR searches in nearby galaxies, which were done using visual identification of Hα bright sources that exhibit high [S II]/Hα line ratios. In addition to our automated source detection, we also differ by using line-ratio “residual” maps to account for diffuse emission, and make use of multiple optical line diagnostics for selection, including the kinematics (from the velocity dispersion maps) of our sources.

To test how well our identification method can recover SNRs, we apply our technique to the well-studied SNR population in M83. This galaxy has an extensive MUSE mosaic (Della Bruna et al. 2022) and is well matched to our sample in terms of stellar mass, star formation rate, and morphology, but is closer in distance (4.61 Mpc; Saha et al. 2006) and thus reaches significantly higher physical resolution (~20 pc). M83 is one of the best-studied galaxies in terms of its SNR population, with a total of 366 SNRs identified using a combination of X-ray, IR and optical data (Dopita et al. 2010; Blair et al. 2012, 2014; Long et al. 2014; Williams et al. 2019), including recent work that used MUSE data to expand the SNR catalog (Long et al. 2022) (L22 hereafter).

The comparison between the objects we detect and the SNRs in the L22 catalog is visualized in Fig. 11. After masking the galaxy center, we identify a parent sample of 125 objects that satisfy at least one of our five criteria (see Appendix D for a complete table). Within the MUSE footprint (excluding the center), L22 identifies 188 SNRs and we recovered 96 of these (51%). If we adjust our threshold for selecting sources in the [S II]/Hα residual map then we can recover more of the known SNRs, but introduce a lot (>200) of additional SNR candidates that are not in the L22 catalog and may just be spurious detections or dynamically induced shock networks in the ISM. In all, 96 out of 125 (76.8%) of our parent sample correspond to known SNRs.

Statistically, the unrecovered 92 objects have lower median Hα surface brightness and [S II]/Hα than the recovered 96 ones, which makes their identification difficult using our method. In addition, the SNRs listed in L22 include 16 SNRs selected using methods other than optical lines (e.g., X-ray and [Fe ii]1.644 μm). We do not expect to detect these objects in the optical, and do not recover any of them.

If we construct a SNR sample in M83, consisting of objects that satisfy at least two criteria using our methods, then 58 SNRs are identified. All of them are also identified as SNRs in L22. This demonstrates that our method can recover SNRs with high confidence by leveraging multiple-line diagnostics.

The 67 objects in the parent sample that are identified using only one method make up our SNR candidate sample. 57% of these are contained in L22. This suggests that we can generally consider many of these SNR candidates to be true SNRs. Of the objects not in the L22 catalog, based on visual inspection some look like reasonable SNR candidates, while others do not. For example, we identify a cluster of regions to the north-east of the galaxy center, along the bar and close to the bar-end, none of which exist in the literature catalogs and may instead be due to shocks associated with the bar dynamics.

This comparison demonstrates that overall our technique produces a reliable SNR catalog. Even our SNR candidates are often correctly identifying SNRs. By applying an automated peak-finding algorithm, we save time by not having to carry out a visual inspection and have a reproducible, homogeneous selection. However, the choices we have made in how we identify objects and construct our SNR catalog clearly do not produce a complete catalog and are also biased toward selecting only ‘typical’ SNRs. We do not recover the very young SNR in M83 identified by Blair et al. (2014), although it could in principle be identified by its broadened [S II] lines. However, it is so compact (1 pc in diameter), that when blended at 100 pc scales with neighboring H II regions it is no longer selected by its kinematic signatures, and we expect such young objects are systematically missed in our catalog.

thumbnail Fig. 7

Parent sample of 2233 objects identified across the 19 PHANGS-MUSE galaxies, compared to the Hα emission. The galaxy name is in the upper left corner. SNRs are indicated by red dots while SNR candidates are in blue. The background shows the Hα map. Regions excluded from our search are outlined in solid black. Most galaxies have the center masked (see Appendix A for more details). We observe SNRs to be distributed across the full field of view, and good qualitative correspondence with Hα bright sites of star formation (see Appendix B for the rest of the galaxies).
thumbnail Fig. 8

Fraction of objects in the parent sample (left) and SNR sample (right) that are identified by specific pairwise combination of our five criteria ([O I]/Hα residual, [S II]/Hα residual, [S II]λ6716 velocity dispersion, BPT distance (OI-OIII), and BPT distance (SII-OIII)).
Table 3

Parent sample catalog generated in this work (see text for details).

thumbnail Fig. 9

Spectra for four SNRs and an H II region in NGC 628. Top two: SNRs with high [S II]/Hα value > 0.8; middle two: SNRs with high-velocity dispersion > 80 km s−1; bottom one: H II region. Characteristic lines, e.g., [O I], [S II], [N II], and Hα are labeled in the plot. All SNRs have broader emission lines than the H II region, though this is most apparent in the two SNRs (NGC0628_17 and NGC0628_48) as they were selected to have particularly broad lines. The sky emission at [O I]λ6300Å, λ6363Å has been masked.
thumbnail Fig. 10

Distribution of identified objects and H II regions density (in grayscale) in the [O I]/Hα versus [S II]/Hα plane. Panel (a): SNRs (red dots); Panel (b): SNR candidates (blue dots). For objects that overlap with a H II region, the symbol is an empty circle in light red or blue. Identified objects tend to occupy different spaces than H II regions in the [O I]/Hα versus [S II]/Hα plane. However, this tendency is more obvious for a SNR than a SNR candidate.

4 Results

4.1 Comparison with H II regions

We compare the integrated emission-line properties of our sample of SNRs with H II regions identified in Groves et al. (2023). As we can see in the comparison in Fig. 12, SNRs have higher [O I]/Hα and [S II]/Hα ratios than H II regions at the same Hα luminosity. These ratios decrease with increasing Hα luminosity. This reflects our selection biases, as a faint SNR blended with a brighter H II region will no longer exhibit the distinctive line ratios and line kinematics that we select for with our methods. This can be seen directly, as SNRs appearing in projection with H II regions (light red circles in the figure) have higher Hα luminosities and correspondingly lower [O I]/Hα and [S II]/Hα ratios.

In Fig. 13a, we see a concentration of both SNRs and H II regions at small reff. This is expected as both closely follow the star formation rate distribution (Cronin et al. 2021), which is more concentrated in the inner regions of galaxies (Leroy et al. 2021). We find hints of radial trends in the SNR line ratios, with higher line ratios found in SNRs compared to H II regions at all radii, and higher values of [S II]/Hα in SNRs at smaller radii. Generally, the [S II] line is brighter for higher shock velocities (Allen et al. 2008). Since SNR shocks typically slow down over time, this is indicative of a larger population of younger SNRs. Given the overall larger number of SNRs at small radii, this likely simply reflects the higher probability of detecting these younger SNRs, although higher ambient densities around SNRs at smaller galactic radii could also be playing a role.

We also find similar trends in the metallicity-sensitive [N II]λ6583/Hα line ratio with radius (Fig. 13b). In most previous SNR surveys that have relied on narrowband imaging, the [N II]λ6583 lines often ended up convolved within the Hα filter and not studied separately; however, with MUSE we can cleanly separate the [N II]λ6583 and Hαλ6562 lines. The higher line ratios at smaller radii could be a consequence of higher ISM metallicities (Kreckel et al. 2019; Pessa et al. 2021).

In Figs. 14a and 14b, we compare the velocity dispersion (corrected for instrumental broadening) of the Hα line and the [S II]λ6716/[S II]λ6730 emission-line ratio, which is a density-sensitive diagnostic3. We can clearly see that while H II regions are concentrated at low-velocity dispersions of 30–40 km s−1, SNRs can have velocity dispersions that reach more than 100 km s−1, owing to the presence of higher-velocity shocked material in the SNRs. Velocity dispersion is also one criterion that we used in Sect. 3.1.2 to identify SNRs. Regarding ambient densities, most of the SNRs have line ratios consistent at the 3σ level with the low-density limit (~10–100 cm−3), with only ~15% of regions showing ratios corresponding to significantly higher densities (>100 cm−3). A similar distribution of values has also been seen in the SNRs of nearby galaxies (e.g., Long et al. 2018, 2022). No visible trend exists between the [S II] doublet ratio and the velocity dispersion of the SNRs.

Figure 14b shows the power of combining velocity dispersion and [S II]/Hα ratio together to distinguish SNRs from H II regions. Our SNRs populate a region of parameter space that has both high-velocity dispersions and high [S II]/Hα values. The distinction is clearer in this regime; however, it becomes complicated when SNRs overlap with H II regions as indicated by light red circles in the figure.

We compare our observed SNRs with the MAPPINGS (Allen et al. 2008) model in Fig. 15. Different grids show different ISM metallicities. Here we adapted the Dopita et al. (2005) abundance model as opposed to the default solar models. Each grid samples a range of shock velocities and magnetic field strengths (0~10 μG). We see that our SNRs are consistent with grids sampling between Large Magellanic Cloud (LMC) (half solar metallicity) and twice solar metallicities (as defined in MAPPINGS model). This is an important validation of our selection process, that we are indeed primarily picking objects that are consistent with the expected parameter space of shocks (Long et al. 2022). This is consistent with the range of metallicities observed in the H II regions detected in these galaxies (Groves et al. 2023).

thumbnail Fig. 11

SNRs in M83. The locations are indicated of SNRs identified by Long et al. (2022) (black boxes) and those identified by applying methods used in this work (red circles for SNRs and blue circles for SNR candidates). All the objects we identify as SNRs are also identified as SNRs in Long et al. (2022). Overall, 77% of our parent sample were previously identified as SNRs in the literature.

4.2 Comparison with the literature: Supernova remnants in other galaxies

In order to further validate our selected SNR population in the 19 MUSE galaxies, we compare their optical line emission properties with the well-studied SNR populations in four nearby galaxies − M33 (Long et al. 2018), M51 (Winkler et al. 2021), M83 (Long et al. 2022), and NGC 6946 (Long et al. 2019). The SNR populations in these galaxies have been vetted not only with optical spectroscopy, but also multiwavelength information from radio, X-ray, and in some cases near-IR lines (e.g., [Fe II]1.64 μm in M83).

In the process of identifying SNRs, a selection bias may exist in more distant galaxies, as we are prone to selecting more luminous SNRs. These could be bright either because they are young or because the explosion happens in a dense environment. At larger distances (and lower physical resolution) we are also more likely to miss SNRs, as they can be blended with brighter H II regions.

We note that within our sample, there is a set of SNRs that exhibit very high Hα luminosities (≥1038 erg s−1) and relatively low [S II]/Hα ratios (≤0.4). In many respects, these characteristics more strongly resemble H II regions than SNRs; however, as they are selected by multiple criteria we remain confident in their classification as SNRs. We see that these SNRs are almost all identified as blended with H II regions in our sample. If we exclude these, we identify very few SNRs with Hα luminosities above 1037.5 erg s−1.

In Fig. 16, we compare the line ratios [S II]/Hα and [N II]/Hα of the SNRs in our PHANGS-MUSE galaxies and the four nearby galaxies. Similar to Fig. 15, we find that our SNRs show a range of values of different line diagnostics that are also observed in SNR surveys of more nearby well-studied galaxies. All four galaxies exhibit a similar range in [S II]/Hα, which likely implies that we are sampling SNRs within a similar range of shock velocities.

Interestingly, we find a bimodal distribution for the metallicity-sensitive [N II]/Hα in our PHANGS-MUSE SNRs. The top branch matches well with SNRs in M83, a high-mass and high-metallicity galaxy, while the bottom branch overlaps with SNRs in M33, a low-mass and low-metallicity galaxy. Some SNRs in M51 occupy the same location in this plot with the top branch of our SNRs but some SNRs have much higher [N II]/Hα ratios than any other SNRs. For NGC 6946, the SNRs lie in between the two branches of our SNRs.

We also attempted to determine whether this trend exists within our 19 galaxies, as shown in Fig. 17. In this figure, SNRs belonging to the same galaxy are represented by the same color, ranging from red to purple as the stellar mass of the galaxy decreases. We observe the same trend as identified from the four literature targets, with higher mass systems exhibiting higher [N II]/Hα ratios than lower mass galaxies, which likely relates to metallicity offsets between these targets. A more detailed investigation of these trends with metallicity is beyond the scope of this work, but also agrees qualitatively with the trends expected based on shock models (see Fig. 15).

In summary, the comparison with MAPPINGS models in Fig. 15 and with the nearby SNR population in Fig. 16 provides a crucial validation of our SNR selection process and the purity of our SNR catalog. Despite not being able to resolve their structure, it is clear that we are picking objects that are consistent with the parameter space of shocks, as well as the observed range of properties of well-studied SNRs in nearby galaxies.

thumbnail Fig. 12

[S II]/Hα (top) and [O I]/Hα (bottom) versus Hα luminosity for SNRs (red dots) compared with the H II regions density distribution (in grayscale). For SNRs that overlap with a H II region, the symbol is an empty circle (in light red). Thresholds of [O I]/Hα=0.1 and [S II]/Hα=0.4 (horizontal dash-dotted lines) are shown. SNRs overlapping with H II regions are selected by other criteria so they can still be identified as SNRs. Emission lines are not corrected for the photoionization.
thumbnail Fig. 13

Radial trends of SNRs and H II regions. r/reff is the ratio of SNR distance to the effective radius of the corresponding galaxy. Panel (a): [S II]/Hα changes with the distance to their corresponding galactic centers for SNRs (red dots) and H II regions (in grayscale). Panel (b): [N II]λ6583/Hα change with the distance to their corresponding galactic centers for SNRs (red dots) and H II regions (in grayscale). For SNRs that overlap with a H II region, the symbol is an empty circle in light red. The drop in the numbers at low r/reff is likely influenced by the mask in the centers of the galaxies.
thumbnail Fig. 14

Comparison of emission-line ratios and velocity dispersion for supernova remnants (SNRs) and H II regions. Panel (a): Line ratios of [S II] versus velocity dispersion of SNRs and H II regions. The top axis shows the FWHM. The Hα velocity dispersion changes with the [S II]λ6716/[S II]λ6730 ratio for SNRs (red dots) and H II regions (in grayscale). The ratios above 1.4 are not physical and are due to noise. The horizontal dash-dotted lines indicate different electron densities with the corresponding value next to it. Panel (b): Velocity dispersion of Hα changes with [S II]λ6716,6730/Hα ratio for SNRs (red dots) and H II regions (in grayscale). For SNRs that overlap with a H II region, the symbol is an empty circle in light red. These two criteria work effectively together in distinguishing between SNRs and H II regions.

5 Discussion

5.1 Expected number of supernova remnants and inferred supernova frequencies

We have identified 2233 objects that we classify as SNRs or SNR candidates in 19 nearby galaxies. To explore the completeness of our sample, we compare it with a theoretical estimation of the number of SNRs we expect in each galaxy. The number of core-collapse (CC) SNRs correlates with the star formation rate (SFR) surface density. To estimate the number of SNRs in each galaxy, we calculated the SN rate per unit area following the formalism presented in Bacchini et al. (2020) as Rcc=ΣSFRfcc,$\[R_{\mathrm{cc}}=\Sigma_{\mathrm{SFR}} f_{\mathrm{cc}},\]$(2)

where ∑SFR is the SFR surface density and fcc≈1.3×10−2M⊙−1$\[f_{\mathrm{cc}} \approx 1.3 \times 10^{-2} M_{\odot}^{-1}\]$ represents the number of core-collapse SNRs for a unit of stellar mass formed under a Kroupa (2002) initial mass function. We neglect SNe associated with type Ia here; however, as estimates of the ratio of CC to type Ia SNe ranges from about 3:1 to 1:1 (Mannucci et al. 2005; Li et al. 2011), accounting for these only decreases the value of fcc by a factor of two at most.

We estimate the SFR surface density from the dust-corrected Hα luminosity (Kennicutt Jr 1998) as SFR⁡(M⊙yr−1)=7.9×10−42L(Hα)(erg s−1).$\[\operatorname{SFR}\left(M_{\odot} y r^{-1}\right)=7.9 \times 10^{-42} L(\mathrm{H} \alpha)\left(\mathrm{erg} \mathrm{~s}^{-1}\right).\]$(3)

We account for dust attenuation using the Balmer decrement, assuming a fixed Hα/H/β = 2.86 (Storey & Hummer 1995) based on Case B recombination and the Calzetti attenuation curve and Rv=3.1 (Calzetti 2001). Combining SFR, the MUSE field of view (FoV) area, and the estimated time (10000–20000 yr; Badenes et al. 2010; Sarbadhicary et al. 2017) for which SNRs remain visible in the optical, we estimate the expected number of SNRs for each galaxy. The comparison of the estimated SNR number and identified SNR number is provided in Fig. 18. For the majority of galaxies, the number of identified SNRs is around half of the expected number of SNRs. When taking the parent sample as a whole, including both SNRs and SNR candidates, we find a better agreement with our simplistic estimation. This suggests that we are recovering about half of the expected SNRs, even more if a large number of our SNR candidates are true SNRs.

In almost all cases, we find fewer SNRs than would be predicted. This is not surprising, as our approach was not designed to be maximally complete. The few galaxies that show more SNRs than predicted do not show any consistent global properties, although one of these (NGC 5068) is the closest galaxy in our sample and another (NGC 4303) is one of the seven galaxies that hosts an AGN. Within the Milky Way, only ~30% of identified SNRs were chosen based on optical observations (Green 2014), suggesting we can expect to miss a large number of objects. On the other hand, we cannot distinguish SNRs due to core collapse from Type Ia ones, which might cause a (small) excess of observed SNRs when compared to our theoretical estimation (Eq. (2)). Overall, given the tantalizing correlation between our expected and observed SNR numbers, this suggests that our parent sample is surprisingly indicative of the total estimated SNR population in the surveyed areas.

If we believe we are achieving a fairly complete census of SNe, we can also directly compute an upper limit to the SN frequency for our sample of galaxies. If we assume SNRs remain visible in the optical for 10000 yr (Sarbadhicary et al. 2017), then on average we find a SN frequency of one per 85 yr per galaxy, with values of the average time between SNe within individual galaxies ranging from 30 to 476 yr (Fig. 19). The galaxies with less frequent SNe (lowest number of SNe per galaxy) are found to be our most distant (D > 18 Mpc) targets, which makes sense as we are more likely to have missed SNe in these cases where our physical resolution is worse. The average SN frequency increases to one per 65 yr per galaxy if we only consider the 13 galaxies with distances less than 18 Mpc. We find no significant secondary correlation with integrated galaxy properties like global SFR or stellar mass, but our sample is still small. These SN frequencies are only upper limits due to the uncertainties in the SNR visibility timescales but are similar to the frequency of one per 40±10 yr estimated for our own Milky Way (Tammann et al. 1994), which is reasonably well matched to the median stellar mass and SFR in our sample of galaxies. Scaling this to the stellar mass contained within the MUSE field of view, this corresponds to a median SN rate per unit mass (SNu; Mannucci et al. 2005) of 1 SNe per 1010 M⊙ per century.

thumbnail Fig. 15

Comparison of supernova remnant properties with MAPPINGS models. Left: Observed [O III]λ5006/Hβ versus [N II]λ6583/Hα for SNRs (red dots). Right: Observed [N II]λ6583/Hα versus [S II]/Hα for SNRs (red dots). For SNRs that overlap with a H II region, the symbol is an empty light red circle. Background grids are MAPPINGS models with different metallicities. Our identified SNRs mostly lie between models for LMC and twice solar metallicities.
thumbnail Fig. 16

Relation between [S II]/Hα ratio and [N II]/Hα ratio for SNRs identified in four nearby galaxies from the literature (Long et al. 2018; Winkler et al. 2021; Long et al. 2022; Long et al. 2019) and 19 MUSE galaxies. SNRs are marked with red dots; SNRs that overlap with a H II region are marked as an empty circle in light red.
thumbnail Fig. 17

Relation between [S II]/Hα ratio and [N II]/Hα ratio for SNRs in 19 galaxies. SNRs in the same galaxy are marked in the same color and with the symbol. In the legend, from top to bottom, the stellar mass of the galaxy decreases, from 1011 M⊙ for NGC 1365 to 10941 M⊙ for NGC 5068.
thumbnail Fig. 18

Comparison between the estimated SNR number from the SFR and objects identified in this paper. The lower limit of the identified object is the number of SNRs and the upper limit is the number in the parent sample including SNR candidates. The lower limit of the estimated SNR number is given by the optically visible time of SNRs for 10 000 yr and the upper limit corresponds to 20 000 yr. The dashed black line is the one-to-one relation. The shaded area provides coverage from two times the recovery rate to 50% of the recovery rate of SNRs.

5.2 Lessons learned when identifying supernova remnants with integral field units

By using multiple criteria to select objects for our parent sample, we find we are able to nearly double the number of objects identified compared to using any one diagnostic alone. The combination of [O I]/Hα residuals and [S II]/Hα residuals gives the largest number of sources, 1749, corresponding to ~78% of the objects in the parent sample. The combination of [O I]/Hα residuals and BPT (OI-OIII) gives the second largest number of sources, 1703. The BPT diagnostics select fewer objects than the line ratio residual approaches (Table 2); however, they are almost all (~90%) identified as robust (isolated) SNRs by our method. Taken together, all five methods have distinct diagnostic powers.

Our results when using these five different techniques to identify SNRs and SNR candidates reflect the fact that the optical properties of SNRs present a variety of emission-line characteristics, particularly when seen in projection against nearby H II regions and DIG. There is no singular feature that provides high confidence in our SNR classification, but instead, the emissionline characteristics provide a continuous spectrum, with increasing confidence of a correct classification the more independent features we observe. However, we lack an independent method to confirm our classification, making it challenging to robustly determine ‘best practices’ for SNR determination. Shock models like MAPPINGS provide some guidance (Fig. 15); however, a richer set of highly resolved emission-line maps of SNRs, either via observations (e.g., Drory et al. 2024) or simulations (e.g., Makarenko et al. 2023), would significantly improve our understanding of how to best identify and classify SNRs that are blended with background emission. By introducing an automated approach to source detection, we hope that future work can begin a more systematic exploration of the completeness limits probed by IFU datasets through the injection of simulated objects to test how well SNRs with different characteristics are recovered.

There are other line ratio combinations, not explored in this work, that could also potentially serve as useful diagnostics for the classification of SNRs in IFU datasets. For example, [Fe II]λ8616/Paschen11λ8862 should be higher in SNRs than in H II regions, as iron lines are strong shock indicators (Oliva et al. 1989; Reipurth et al. 2000). Unlike other optical [Fe II] lines that are produced by continuum pumping fluorescence (such as λ4287 or λ5158), λ8617 is essentially produced by collisional excitation at low excitation levels (Rodríguez 1999; Verner et al. 2000; Méndez-Delgado et al. 2021; Mendoza et al. 2023). This characteristic, combined with the low ionization potential of Fe+, makes it a promising shock indicator analogous to other [Fe II] lines at far-infrared wavelengths, such as the bright [Fe II] 1.64 μm line (Antoniucci et al. 2014). Additionally, while Fe is typically depleted into dust grains in HII regions (Rodríguez 2002; Izotov et al. 2006), shock waves are capable of efficiently destroying these dust grains and releasing Fe into its gaseous phase (Mesa-Delgado et al. 2009; Méndez-Delgado et al. 2022). Therefore, in the case of overlapping between SNRs and H II regions, this indicator could favor the conditions of the SNRs. We see clear evidence of this line in some of our spectra; however, the strong neighboring sky lines make a robust extraction of line fluxes challenging and beyond the scope of this work.

The [S III]/9069/[S II]λ6716,6730 ratio might also be used to diagnose the ionization mechanism (Long et al. 2022). In H II regions, [S III]/[S II] correlates strongly with ionization parameter (Kewley et al. 2001), but it can span a wide range of values depending on the local ISM density and the radiation field produced by the ionizing source (Groves et al. 2023). As shown in L22, high values of [S III]/[S II] provide strong indications of H II regions, but low values are more challenging to interpret as the separation with H II regions is not clear and [S III] can be too weak.

thumbnail Fig. 19

Number of SNe as a function of stellar mass (left) and SFR (right) considering the entire parent sample, colored by galaxy distance and considering galaxy integrated properties measured only within the MUSE field of view. These SN counts can also directly be converted to a SN frequency by assuming SNRs are visible in the optical for 10 000 yr. Lines of constant SN rate per unit mass (SNu) are also shown. The star markers represent the number for the parent sample while the circles are the number of SNRs. The line between them shows the range of possible true SNR numbers.

5.3 Classification of the candidate sample and blending with HII regions

In Sec. 4, we only focus on properties of 1166 SNRs, whereas there are still 1067 SNR candidates. Even though we are less certain these are SNRs, they are robust detections using diagnostics that give strong evidence for shock excitation. We also saw in Sec. 5.1 that including these objects provides better agreement with the expected number of SNRs given the SFR in each galaxy, and in Sec. 3.3 we showed that in M83 about half of our candidates are true SNRs in Long et al. (2022).

One reason it is difficult for us to confirm the nature of individual objects is due to the blending of emission arising from different ionizing sources within our ~100 pc apertures. As shown in Table 2, 202 SNRs (17%) and 569 SNR candidates (53%) are blended with H II regions. In Fig. 20, we see that more than half of the SNR candidates overlap with H II regions and that these objects exhibit systematically lower [S II]/Hα line ratios. This is expected if bright photoionized emission, which exhibits systematically lower ratios, is blended with an unresolved SNR. Clearly, we have more confidence in classifying objects as SNRs when they are isolated. At our resolution, we cannot determine if blended objects are physically associated with H II regions or they just appear in projection against the H II region along our line of sight. In total, 35% of objects in the parent sample are not isolated sources.

As is shown in Sec. 4.1, the majority of SNRs that exhibit large Hα luminosities (Fig. 12) are found to be blended with H II regions and have the lowest [S II]/Hα and [O I]/Hα line ratios. Interestingly, if we remove these blended SNRs from the sample then we find a relatively sharp upper Hα luminosity threshold of ~1037.5 erg s−1. At these distances, it remains a challenge to properly decompose SNRs from the surrounding photoionized emission, but from our combination of multiple line diagnostics and novel use of the [S II]/Hα residual to account for blended photoionized emission, we have succeeded in identifying a large sample of objects that would have been missed by previous searches.

This can be seen more directly in the four SNR examples shown in Fig. 1. In the third row, it is clear from the Hα intensity that this object sits in a complex star-forming environment, as reflected by the extended photoionization-dominated nebulae outlined in red (left panel). In the [S II]/Hα map this object is not pronounced, but in the residual map it becomes more distinct from the background. Detailed deblending of such objects is beyond the scope of this work but it could be improved in the future by the use of higher-resolution imaging of the ionized hydrogen line emission from HST (Barnes et al. 2022) and JWST (Barnes et al. 2023).

thumbnail Fig. 20

Probability distribution functions of [S II]/Hα for SNRs and SNR candidates that are isolated or overlap with H II region. More SNR candidates (bottom) overlap with H II regions (blended) and have lower [S II]/Hα values than SNRs (top).

5.4 Characterizing the shocked regions

Congiu et al. (2023) and Egorov et al. (2023), have both presented catalogs of objects classified as ‘shocked’ using the same MUSE data for the same sample of 19 galaxies, and here we compare our catalog with their results.

C23 used a contour-based algorithm to identify peaks in a weighted average of [O III]λ5007, Hα and [S II]λ6716, 6731 flux maps. Regions are defined by the surface brightness contours around the closest peak. In this way, the same set of 19 MUSE galaxies are segmented into 40 920 regions and by comparing with different model grids the regions are classified into five types: H II regions, shocked regions (6336), ambiguous, planetary nebula, and unknown. We cross-matched these regions with our parent sample and compare their classification.

Our sample is 35% smaller in total than the ‘shocked’ sample in C23. 1676 objects (75%) of our parent sample can be associated with the C23 segmentation mask. The 25% of objects missed is comparable to the ~35% of objects that had not been identified in Groves et al. (2023), suggesting that many SNRs are missed when selecting by emission-line flux and morphology alone. Of the objects that can be cross-matched to the C23 catalog, 676 (38%) of our regions are classified as shocked. This is only about 10% of the total number of shocked regions in C23. The SNRs sample shows a better overall agreement, with 512 (56%) of our SNRs identified as shocked in C23. The differences in classification likely reflect a difference between our empirical diagnostics and the model grids used in that work. Moreover, line ratios from our integrated spectra can appear consistent with photoionization (Fig. 10b), given the frequent blending with H II regions and the contamination from the DIG background, which would result in a C23 classification of ‘H II region’.

To further demonstrate that the usage of a combination of emission-line ratios and gas kinematics helps to improve the identification of SNRs, we compare our objects with about 1500 regions of locally elevated Hα velocity dispersion identified by E23 in the same sample of galaxies. Only 23% of objects in our catalog reside in the regions identified by E23. This indicates that the kinematic information alone is usually not sufficient for SNR identification. Meanwhile, about 320 unresolved regions with high-velocity dispersion from E23 are considered to be probable SNR candidates because of their high [S II]/Hα line ratio. 65% of these regions coincide with the objects in our catalog. Visual inspection of the remaining objects demonstrates that the majority of them reside in very low surface brightness regions making their classification uncertain, and leaving us with only 25 (about 7%) SNR candidates from E23 that are not present in our catalog.

6 Conclusions

Using optical IFU data from the PHANGS-MUSE survey (Emsellem et al. 2022), we constructed a new catalog identifying 1166 SNRs plus 1067 SNR candidates in 19 nearby galaxies, doubling the number of known SNRs in the local universe.

We summarize the key results below:

  • We used a total of five criteria to identify objects with astrodendro in our ~50 pc resolution observations. Four criteria are based on line-ratio diagnostics ([S II]/Hα, [O I]/Hα, and two BPT diagnostics), and the fifth criterion is based on the [S II] emission-line velocity dispersion. Hα morphology is not used as a criterion in our selection;

  • Instead of identifying sources directly in the line ratio maps, we constructed [S II]/Hα and [O I]/Hα residual maps, where we could identify objects that emit brightly at these line ratios given their Hα surface brightness. This allowed us to naturally account for background and foreground photoionized emission, consisting of contributions from both H II regions and DIG;

  • Using multiple diagnostics increased our sample size and provided a natural mechanism to robustly construct a SNR catalog. Our SNR catalog consists of all objects identified by at least two criteria, or having integrated line ratios meeting both of the criteria [S II]/Hα > 0.4 and [O I]/Hα > 0.1. The remaining objects are classified as SNR candidates;

  • We validated our method by applying it to MUSE observations of M83 (Della Bruna et al. 2022). All of the SNRs we identified using our method are found in literature catalogs (L22), along with 77% of our SNR candidates, giving us confidence that our method is robustly identifying SNRs;

  • We demonstrated that the SNRs in our sample have line ratios and line kinematics distinct from H II regions, are consistent with MAPPINGS shock models, and are similar to literature SNR samples. Variations between galaxies in our sample are consistent with the expected variations with metallicity, given the difference in galaxy stellar mass found in our sample;

  • We estimated the total number of expected SNRs based on theoretical models, given the current SFR in each galaxy, and we find this agrees within a factor of two with the total number of objects identified. This suggests that while our sample is incomplete, it is illustrative of the SNR population;

  • By using line ratios and line kinematics to identify our objects, we were able to locate SNRs that are seen in projection against bright star-forming regions. About 17% of our SNRs are blended with H II regions, while 52% of the SNR candidates are blended. This blending makes confirmation of SNR candidates as true SNRs challenging.

Building statistical samples of SNRs provides a new avenue to quantify the impact of stellar feedback; however, it remains difficult to estimate how complete these SNR catalogs are. By studying a large sample of galaxies with IFU data, we have begun exploring how automatic identification of objects allows us to systematically produce large uniform catalogs. Based on the multiple different criteria explored in this work, [O I]/Hα is the most successful at selecting objects for our parent sample (62%). BPT diagnostics select fewer objects, but they are almost all (~90%) identified as robust (isolated) SNRs.

With these techniques, it is easy to imagine that the next steps can include the injection of simulated objects with different characteristics to quantitatively determine how many can be recovered. High physical resolution simulations have also already begun to provide insights into the interpretation of these variable diagnostics, and the importance of accounting for contaminating background emission (Makarenko et al. 2023). Given the growing number of galaxies observed with multiple emission lines mapped at ~100 pc scales, these automated approaches to identifying SNRs have great potential to significantly increase our SNR samples in the coming years.

Data availability

Tables 3 and D.1 are available at the CDS via anonymous ftp to cdsarc.cds.unistra.fr (130.79.128.5) or via https://cdsarc.cds.unistra.fr/viz-bin/cat/J/A+A/690/A161.

Acknowledgements

We thank Dr. Angela Adamo for generously sharing M83 data. We extend our sincere gratitude to Dr. Eva Grebel and Dr. Andreas Sander for their valuable discussions and contributions that inspired the concept of the supernova rate (SN rate) in our research. J.L., K.K., J.E.M.-D. and O.E. gratefully acknowledge funding from the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) in the form of an Emmy Noether Research Group (grant number KR4598/2–1, PI Kreckel) and the European Research Council’s starting grant ERC StG-101077573 (“ISM-METALS”). G.A.B. acknowledges the support from the ANID Basal project FB210003. F.B. acknowledges support from the INAF Fundamental Astrophysics program 2022. S.C.O.G. and R.S.K. acknowledge financial support from the European Research Council via the ERC Synergy Grant “ECOGAL” (project ID 855130), from the German Excellence Strategy via the Heidelberg Cluster of Excellence (EXC 2181 – 390900948) “STRUCTURES”, and from the German Ministry for Economic Affairs and Climate Action in project “MAINN” (funding ID 50OO2206). Based on observations collected at the European Southern Observatory under ESO programmes 094.C-0623 (PI: Kreckel), 095.C-0473, 098.C-0484 (PI: Blanc), 1100.B-0651 (PHANGS-MUSE; PI: Schinnerer), as well as 094.B-0321 (MAGNUM; PI: Marconi), 099.B-0242, 0100.B-0116, 098.B-0551 (MAD; PI: Carollo) and 097.B-0640 (TIMER; PI: Gadotti). This research made use of astrodendro, a Python package to compute dendrograms of Astronomical data (http://www.dendrograms.org/). Other main Python Packages that have been used are ASTROPY (Robitaille et al. 2013; Price-Whelan et al. 2018, 2022), NUMPY (Harris et al. 2020), MATPLOTLIB (Hunter 2007) and PYNEB (Luridiana et al. 2015). Table 1 lists the distances that were compiled by Anand et al. (2021) from Freedman et al. (2001); Nugent et al. (2006); Jacobs et al. (2009); Kourkchi & Tully (2017); Shaya et al. (2017); Kourkchi et al. (2020) and Scheuermann et al. (2022).

Appendix A Masked regions within each galaxy

Environmental masks that applied to each criterion for 19 galaxies are described in Table A.1. Environment masks are defined in Querejeta et al. (2021), 1 is center mask, 2 is bar mask, 3 is bar ends mask, and 9 is interbar mask. They are applied to avoid over-selecting shocked regions around the center or along the stellar bar. For NGC 1566 and NGC 3351, their center is not sufficient to cover the above shocked regions, so the convolved mask of center and bar masks were applied, with corresponding Gaussian kernels as listed.

Table A.1

Environmental masks applied to each criterion for 19 galaxies.

Appendix B Atlas of supernova remnant locations in PHANGS-MUSE galaxies

A total of 2233 objects are identified in 19 PHANGS-MUSE galaxies; only three of them are shown in Fig. 7. The remaining 16 galaxies are shown in Figs. B.1–B.4. They display the distribution of identified objects across different galaxies.

thumbnail Fig. B.1

Same as in Fig. 7.
thumbnail Fig. B.2

Same as in Fig. 7.
thumbnail Fig. B.3

Same as in Fig. 7.
thumbnail Fig. B.4

Same as in Fig. 7

Appendix C Fitted emission lines and masked sky lines while running the data analysis pipeline

In Table C.1, we listed the fitted emission lines, sky lines that have been masked, ISM absorption regions that need to be masked while running the data analysis pipeline in the first column. Not all lines are fitted independently, their kinematic dependency is shown in the second column.

Table C.1

Fitted emission lines and masked sky lines while running DAP.

Appendix D Catalog of objects identified in M83

In Sec. 3.3, we validate our method by applying it to archival MUSE observations of M83. The galaxy is two to four times closer (D = 4.61 Mpc; Saha et al. 2006) than galaxies in our PHANGS-MUSE sample, and the line maps and line kinematic maps used are taken from Della Bruna et al. (2022), which used a different approach to line fitting and spatial binning. For these reasons, we do not consider this catalog completely equivalent to what was produced for the PHANGS-MUSE galaxies, but we provide part of the object locations and classifications (SNR or SNR candidate) for all 125 objects in the Table D.1.

Table D.1

Identified objects in M83 (part).

References

  1. Allen, M. G., Groves, B. A., Dopita, M. A., Sutherland, R. S., & Kewley, L. J. 2008, ApJS, 178, 20 [Google Scholar]
  2. Anand, G. S., Lee, J. C., Van Dyk, S. D., et al. 2021, MNRAS, 501, 3621 [Google Scholar]
  3. Anderson, L., Wang, Y., Bihr, S., et al. 2017, ApJ, 605, A58 [Google Scholar]
  4. Antoniucci, S., La Camera, A., Nisini, B., et al. 2014, A&A, 566, A129 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Asvarov, A. I. 2014, A&A, 561, A70 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  6. Athanassoula, E. 1992, MNRAS, 259, 345 [Google Scholar]
  7. Auchettl, K., Lopez, L. A., Badenes, C., et al. 2019, ApJ, 871, 64 [NASA ADS] [CrossRef] [Google Scholar]
  8. Bacchini, C., Fraternali, F., Iorio, G., et al. 2020, ApJ, 641, A70 [Google Scholar]
  9. Bacon, R., Accardo, M., Adjali, L., et al. 2010, in Ground-based and Airborne Instrumentation for Astronomy III, 7735, SPIE, 131 [Google Scholar]
  10. Bacon, R., Conseil, S., Mary, D., et al. 2017, ApJ, 608, A1 [Google Scholar]
  11. Badenes, C., Harris, J., Zaritsky, D., & Prieto, J. L. 2009, ApJ, 700, 727 [NASA ADS] [CrossRef] [Google Scholar]
  12. Badenes, C., Maoz, D., & Draine, B. T. 2010, MNRAS, 407, 1301 [NASA ADS] [CrossRef] [Google Scholar]
  13. Baldwin, J. A., Phillips, M. M., & Terlevich, R. 1981, PASP, 93, 5 [Google Scholar]
  14. Barnes, A., Chandar, R., Kreckel, K., et al. 2022, ApJ, 662, L6 [Google Scholar]
  15. Barnes, A. T., Watkins, E. J., Meidt, S. E., et al. 2023, ApJ, 944, L22 [NASA ADS] [CrossRef] [Google Scholar]
  16. Belfiore, F., Santoro, F., Groves, B., et al. 2022, ApJ, 659, A26 [Google Scholar]
  17. Bell, A. R., Schure, K. M., Reville, B., & Giacinti, G. 2013, MNRAS, 431, 415 [Google Scholar]
  18. Blair, W. P., & Long, K. S. 1997, ApJS, 108, 261 [NASA ADS] [CrossRef] [Google Scholar]
  19. Blair, W. P., & Long, K. S. 2004, ApJS, 155, 101 [NASA ADS] [CrossRef] [Google Scholar]
  20. Blair, W. P., Kirshner, R. P., & Chevalier, R. A. 1981, AJ, 247, 879 [NASA ADS] [Google Scholar]
  21. Blair, W. P., Winkler, P. F., & Long, K. S. 2012, ApJS, 203, 8 [NASA ADS] [CrossRef] [Google Scholar]
  22. Blair, W. P., Chandar, R., Dopita, M. A., et al. 2014, ApJ, 788, 55 [NASA ADS] [CrossRef] [Google Scholar]
  23. Blandford, R., & Eichler, D. 1987, Phys. Rep., 154, 1 [Google Scholar]
  24. Blasi, P. 2013, A&A Rev., 21, 70 [NASA ADS] [CrossRef] [Google Scholar]
  25. Bozzetto, L. M., Filipović, M. D., Sano, H., et al. 2023, MNRAS, 518, 2574 [Google Scholar]
  26. Calzetti, D. 2001, PASP, 113, 1449 [Google Scholar]
  27. Cappellari, M. 2017, MNRAS, 466, 798 [Google Scholar]
  28. Cappellari, M., & Emsellem, E. 2004, PASP, 116, 138 [Google Scholar]
  29. Chevalier, R. A. 2005, ApJ, 619, 839 [CrossRef] [Google Scholar]
  30. Chevalier, R. A., Kirshner, R. P., & Raymond, J. C. 1980, ApJ, 235, 186 [NASA ADS] [CrossRef] [Google Scholar]
  31. Chomiuk, L., & Wilcots, E. M. 2009, ApJ, 703, 370 [NASA ADS] [CrossRef] [Google Scholar]
  32. Cid Fernandes, R., Carvalho, M. S., Sánchez, S. F., de Amorim, A., & Ruschel-Dutra, D. 2021, MNRAS, 502, 1386 [NASA ADS] [CrossRef] [Google Scholar]
  33. Congiu, E., Blanc, G. A., Belfiore, F., et al. 2023, ApJ, 672, A148 [Google Scholar]
  34. Cotton, W., Kothes, R., Camilo, F., et al. 2024, ApJS, 270, 21 [NASA ADS] [CrossRef] [Google Scholar]
  35. Cronin, S. A., Utomo, D., Leroy, A. K., et al. 2021, ApJ, 923, 86 [NASA ADS] [CrossRef] [Google Scholar]
  36. Della Bruna, L., Adamo, A., Amram, P., et al. 2022, A&A, 660, A77 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  37. Díaz-Rodríguez, M., Murphy, J. W., Rubin, D. A., et al. 2018, ApJ, 861, 92 [CrossRef] [Google Scholar]
  38. Dodorico, S., Dopita, M., & Benvenuti, P. 1980, A&AS, 40, 67 [Google Scholar]
  39. Dokara, R., Brunthaler, A., Menten, K., et al. 2021, ApJ, 651, A86 [Google Scholar]
  40. Dopita, M. A., Groves, B. A., Fischera, J., et al. 2005, ApJ, 619, 755 [Google Scholar]
  41. Dopita, M. A., Blair, W. P., Long, K. S., et al. 2010, ApJ, 710, 964 [NASA ADS] [CrossRef] [Google Scholar]
  42. Dopita, M. A., Seitenzahl, I. R., Sutherland, R. S., et al. 2016, ApJ, 826, 150 [NASA ADS] [CrossRef] [Google Scholar]
  43. Drissen, L., Martin, T., Rousseau-Nepton, L., et al. 2019, MNRAS, 485, 3930 [NASA ADS] [CrossRef] [Google Scholar]
  44. Drory, N., Blanc, G. A., Kreckel, K., et al. 2024, arXiv e-prints [arXiv:2405.01637] [Google Scholar]
  45. Dubner, G. M., Holdaway, M., Goss, W. M., & Mirabel, I. F. 1998, AJ, 116, 1842 [Google Scholar]
  46. Dwek, E., & Arendt, R. G. 1992, ARA&A, 30, 11 [Google Scholar]
  47. Egorov, O. V., Kreckel, K., Glover, S. C. O., et al. 2023, A&A, 678, A153 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  48. Elwood, B. D., Murphy, J. W., & Díaz-Rodríguez, M. 2019, MNRAS, 483, 4551 [NASA ADS] [CrossRef] [Google Scholar]
  49. Emsellem, E., Schinnerer, E., Santoro, F., et al. 2022, ApJ, 659, A191 [Google Scholar]
  50. Fesen, R. A., Blair, W. P., & Kirshner, R. P. 1985, ApJ, 292, 29 [Google Scholar]
  51. Freedman, W. L., Madore, B. F., Gibson, B. K., et al. 2001, ApJ, 553, 47 [Google Scholar]
  52. Galvin, T. J., & Filipovic, M. D. 2014, Serb. Astron. J., 189, 15 [Google Scholar]
  53. Gordon, S. M., Duric, N., Kirshner, R. P., Goss, W. M., & Viallefond, F. 1999, ApJS, 120, 247 [NASA ADS] [CrossRef] [Google Scholar]
  54. Green, D. A. 2005, Mem. Soc. Astron. Italiana, 76, 534 [NASA ADS] [Google Scholar]
  55. Green, D. A. 2019, J. Astrophys. Astron., 40, 36 [Google Scholar]
  56. Green, D. A. 2014, Bull. Astron. Soc. India, 42, 47 [NASA ADS] [Google Scholar]
  57. Green, D. A., Reynolds, S. P., Borkowski, K. J., et al. 2008, MNRAS, 387, L54 [NASA ADS] [CrossRef] [Google Scholar]
  58. Groves, B., Kreckel, K., Santoro, F., et al. 2023, MNRAS, 520, 4902 [Google Scholar]
  59. Haffner, L., Dettmar, R.-J., Beckman, J., et al. 2009, Rev. Mod. Phys., 81, 969 [CrossRef] [Google Scholar]
  60. Harris, C. R., Millman, K. J., Van Der Walt, S. J., et al. 2020, Nature, 585, 357 [NASA ADS] [CrossRef] [Google Scholar]
  61. Holland-Ashford, T., Lopez, L. A., Auchettl, K., Temim, T., & Ramirez-Ruiz, E. 2017, ApJ, 844, 84 [Google Scholar]
  62. Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [NASA ADS] [CrossRef] [Google Scholar]
  63. Izotov, Y. I., Stasińska, G., Meynet, G., Guseva, N. G., & Thuan, T. X. 2006, A&A, 448, 955 [CrossRef] [EDP Sciences] [Google Scholar]
  64. Jacobs, B. A., Rizzi, L., Tully, R. B., et al. 2009, ApJ, 138, 332 [CrossRef] [Google Scholar]
  65. Jennings, Z. G., Williams, B. F., Murphy, J. W., et al. 2014, ApJ, 795, 170 [NASA ADS] [CrossRef] [Google Scholar]
  66. Katsuda, S., Morii, M., Janka, H.-T., et al. 2018, ApJ, 856, 18 [Google Scholar]
  67. Kauffmann, G., Heckman, T. M., Tremonti, C., et al. 2003, MNRAS, 346, 1055 [Google Scholar]
  68. Kavanagh, P. J., Sasaki, M., Filipović, M. D., et al. 2022, MNRAS, 515, 4099 [NASA ADS] [CrossRef] [Google Scholar]
  69. Kennicutt Jr, R. C. 1998, ARA&A, 36, 189 [NASA ADS] [CrossRef] [Google Scholar]
  70. Kewley, L. J., Dopita, M. A., Sutherland, R., Heisler, C., & Trevena, J. 2001, ApJ, 556, 121 [NASA ADS] [CrossRef] [Google Scholar]
  71. Kewley, L. J., Groves, B., Kauffmann, G., & Heckman, T. 2006, MNRAS, 372, 961 [Google Scholar]
  72. Koo, B.-C., Kim, C.-G., Park, S., & Ostriker, E. C. 2020, ApJ, 905, 35 [NASA ADS] [CrossRef] [Google Scholar]
  73. Koplitz, B., Johnson, J., Williams, B. F., et al. 2021, ApJ, 916, 58 [NASA ADS] [CrossRef] [Google Scholar]
  74. Kopsacheili, M., Zezas, A., & Leonidaki, I. 2020, MNRAS, 491, 889 [CrossRef] [Google Scholar]
  75. Kopsacheili, M., Jiménez-Palau, C., Galbany, L., Boumis, P., & González-Díaz, R. 2024, MNRAS, stae874 [Google Scholar]
  76. Kourkchi, E., & Tully, R. B. 2017, ApJ, 843, 16 [NASA ADS] [CrossRef] [Google Scholar]
  77. Kourkchi, E., Courtois, H. M., Graziani, R., et al. 2020, ApJ, 159, 67 [Google Scholar]
  78. Krause, O., Birkmann, S. M., Usuda, T., et al. 2008, Science, 320, 1195 [NASA ADS] [CrossRef] [Google Scholar]
  79. Kreckel, K., Ho, I.-T., Blanc, G., et al. 2019, ApJ, 887, 80 [NASA ADS] [CrossRef] [Google Scholar]
  80. Kroupa, P. 2002, Science, 295, 82 [Google Scholar]
  81. Lacey, C., Duric, N., & Goss, W. M. 1997, ApJS, 109, 417 [NASA ADS] [CrossRef] [Google Scholar]
  82. Lang, P., Meidt, S. E., Rosolowsky, E., et al. 2020, ApJ, 897, 122 [CrossRef] [Google Scholar]
  83. Lee, J. H., & Lee, M. G. 2014a, ApJ, 786, 130 [Google Scholar]
  84. Lee, J. H., & Lee, M. G. 2014b, ApJ, 793, 134 [NASA ADS] [CrossRef] [Google Scholar]
  85. Leroy, A. K., Schinnerer, E., Hughes, A., et al. 2021, ApJS, 257, 43 [NASA ADS] [CrossRef] [Google Scholar]
  86. Levy, R. C., Bolatto, A. D., Sánchez, S. F., et al. 2019, ApJ, 882, 84 [NASA ADS] [CrossRef] [Google Scholar]
  87. Li, W., Leaman, J., Chornock, R., et al. 2011, MNRAS, 412, 1441 [NASA ADS] [CrossRef] [Google Scholar]
  88. Long, K. S. 2017, Handbook of Supernovae (Berlin: Springer), 1 [Google Scholar]
  89. Long, K. S., Helfand, D. J., & Grabelsky, D. A. 1981, ApJ, 248, 925 [Google Scholar]
  90. Long, K. S., Blair, W. P., Winkler, P. F., et al. 2010, ApJS, 187, 495 [NASA ADS] [CrossRef] [Google Scholar]
  91. Long, K. S., Kuntz, K. D., Blair, W. P., et al. 2014, ApJS, 212, 21 [Google Scholar]
  92. Long, K. S., Blair, W. P., Milisavljevic, D., Raymond, J. C., & Winkler, P. F. 2018, ApJ, 855, 140 [NASA ADS] [CrossRef] [Google Scholar]
  93. Long, K. S., Winkler, P. F., & Blair, W. P. 2019, ApJ, 875, 85 [NASA ADS] [CrossRef] [Google Scholar]
  94. Long, K. S., Blair, W. P., Winkler, P. F., & Lacey, C. K. 2020, ApJ, 899, 14 [NASA ADS] [CrossRef] [Google Scholar]
  95. Long, K. S., Blair, W. P., Winkler, P. F., et al. 2022, ApJ, 929, 144 [NASA ADS] [CrossRef] [Google Scholar]
  96. López-Cobá, C., Sánchez, S. F., Lin, L., et al. 2022, ApJ, 939, 40 [CrossRef] [Google Scholar]
  97. Luridiana, V., Morisset, C., & Shaw, R. A. 2015, ApJ, 573, A42 [Google Scholar]
  98. Maggi, P., Haberl, F., Kavanagh, P. J., et al. 2016, ApJ, 585, A162 [Google Scholar]
  99. Maggi, P., Filipović, M. D., Vukotić, B., et al. 2019, A&A, 631, A127 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  100. Magnier, E. A., Prins, S., van Paradijs, J., et al. 1995, A&AS, 114, 215 [NASA ADS] [Google Scholar]
  101. Makarenko, E. I., Walch, S., Clarke, S. D., et al. 2023, MNRAS, 523, 1421 [NASA ADS] [CrossRef] [Google Scholar]
  102. Mannucci, F., Della Valle, M., Panagia, N., et al. 2005, ApJ, 433, 807 [Google Scholar]
  103. Maoz, D., & Badenes, C. 2010, MNRAS, 407, 1314 [NASA ADS] [CrossRef] [Google Scholar]
  104. Maoz, D., Mannucci, F., & Nelemans, G. 2014, ARA&A, 52, 107 [Google Scholar]
  105. Marcowith, A., Bret, A., Bykov, A., et al. 2016, Rep. Progr. Phys., 79, 046901 [CrossRef] [Google Scholar]
  106. Mathewson, D., & Clarke, J. 1973, ApJ, 180, 725 [NASA ADS] [CrossRef] [Google Scholar]
  107. Matonick, D. M., & Fesen, R. A. 1997, ApJS, 112, 49 [NASA ADS] [CrossRef] [Google Scholar]
  108. McKee, C. F., & Hollenbach, D. J. 1980, ARA&A, 18, 219 [NASA ADS] [CrossRef] [Google Scholar]
  109. Méndez-Delgado, J. E., Henney, W. J., Esteban, C., et al. 2021, ApJ, 918, 27 [CrossRef] [Google Scholar]
  110. Méndez-Delgado, J. E., Esteban, C., García-Rojas, J., & Henney, W. J. 2022, MNRAS, 514, 744 [CrossRef] [Google Scholar]
  111. Mendoza, C., Méndez-Delgado, J. E., Bautista, M., García-Rojas, J., & Morisset, C. 2023, Atoms, 11, 63 [NASA ADS] [CrossRef] [Google Scholar]
  112. Mesa-Delgado, A., Esteban, C., García-Rojas, J., et al. 2009, MNRAS, 395, 855 [NASA ADS] [CrossRef] [Google Scholar]
  113. Milisavljevic, D., Temim, T., De Looze, I., et al. 2024, ApJ, 965, L27 [CrossRef] [Google Scholar]
  114. Moumen, I., Robert, C., Devost, D., et al. 2019, MNRAS, 488, 803 [NASA ADS] [CrossRef] [Google Scholar]
  115. Nugent, P., Sullivan, M., Ellis, R., et al. 2006, ApJ, 645, 841 [CrossRef] [Google Scholar]
  116. Oliva, E., Moorwood, A., & Danziger, I. 1989, A&A, 214, 307 [NASA ADS] [Google Scholar]
  117. Pannuti, T. G., Duric, N., Lacey, C. K., et al. 2000, ApJ, 544, 780 [NASA ADS] [CrossRef] [Google Scholar]
  118. Pannuti, T. G., Duric, N., Lacey, C. K., et al. 2002, ApJ, 565, 966 [NASA ADS] [CrossRef] [Google Scholar]
  119. Parker, R. A., Gull, T. R., Kirschner, R. P., & Kirshner, R. P. 1979, An Emissionline Survey of the Milky Way, 434 (Scientific and Technical Information Office, National Aeronautics and Space ...) [Google Scholar]
  120. Patnaude, D., & Badenes, C. 2017, in Handbook of Supernovae, eds. A. W. Alsabti, & P. Murdin, 2233 [CrossRef] [Google Scholar]
  121. Pessa, I., Schinnerer, E., Belfiore, F., et al. 2021, ApJ, 650, A134 [Google Scholar]
  122. Points, S. D., Long, K. S., Winkler, P. F., & Blair, W. P. 2019, ApJ, 887, 66 [NASA ADS] [CrossRef] [Google Scholar]
  123. Price-Whelan, A. M., Sipőcz, B., Günther, H., et al. 2018, ApJ, 156, 123 [CrossRef] [Google Scholar]
  124. Price-Whelan, A. M., Lim, P. L., Earl, N., et al. 2022, ApJ, 935, 167 [NASA ADS] [CrossRef] [Google Scholar]
  125. Querejeta, M., Schinnerer, E., Meidt, S., et al. 2021, ApJ, 656, A133 [Google Scholar]
  126. Raymond, J. C. 1979, ApJS, 39, 1 [NASA ADS] [CrossRef] [Google Scholar]
  127. Reach, W. T., Rho, J., Tappe, A., et al. 2006, AJ, 131, 1479 [NASA ADS] [CrossRef] [Google Scholar]
  128. Reipurth, B., Yu, K. C., Heathcote, S., Bally, J., & Rodríguez, L. F. 2000, ApJ, 120, 1449 [CrossRef] [Google Scholar]
  129. Rhea, C. L., Rousseau-Nepton, L., Moumen, I., et al. 2023, RAS Tech. Instrum., 2, 345 [NASA ADS] [CrossRef] [Google Scholar]
  130. Robitaille, T. P., Tollerud, E. J., Greenfield, P., et al. 2013, ApJ, 558, A33 [Google Scholar]
  131. Rodríguez, M. 1999, A&A, 348, 222 [NASA ADS] [Google Scholar]
  132. Rodríguez, M. 2002, A&A, 389, 556 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  133. Rosolowsky, E. W., Pineda, J. E., Kauffmann, J., & Goodman, A. A. 2008, ApJ, 679, 1338 [Google Scholar]
  134. Roth, M. M., Sandin, C., Kamann, S., et al. 2018, A&A, 618, A3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  135. Russell, T. D., White, R. L., Long, K. S., et al. 2020, MNRAS, 495, 479 [NASA ADS] [CrossRef] [Google Scholar]
  136. Sabbadin, F., Minello, S., & Bianchini, A. 1977, A&A, 60, 147 [NASA ADS] [Google Scholar]
  137. Saha, A., Thim, F., Tammann, G., Reindl, B., & Sandage, A. 2006, ApJS, 165, 108 [NASA ADS] [CrossRef] [Google Scholar]
  138. Santoro, F., Kreckel, K., Belfiore, F., et al. 2022, ApJ, 658, A188 [Google Scholar]
  139. Sarangi, A., Matsuura, M., & Micelotta, E. R. 2018, Space Sci. Rev., 214, 63 [NASA ADS] [CrossRef] [Google Scholar]
  140. Sarbadhicary, S. K., Badenes, C., Chomiuk, L., Caprioli, D., & Huizenga, D. 2017, MNRAS, 464, 2326 [NASA ADS] [CrossRef] [Google Scholar]
  141. Schaefer, B. E., & Pagnotta, A. 2012, Nature, 481, 164 [NASA ADS] [CrossRef] [Google Scholar]
  142. Scheuermann, F., Kreckel, K., Anand, G. S., et al. 2022, MNRAS, 511, 6087 [NASA ADS] [CrossRef] [Google Scholar]
  143. Shaya, E. J., Tully, R. B., Hoffman, Y., & Pomarède, D. 2017, ApJ, 850, 207 [NASA ADS] [CrossRef] [Google Scholar]
  144. Slane, P., Bykov, A., Ellison, D. C., Dubner, G., & Castro, D. 2015, Space Sci. Rev., 188, 187 [NASA ADS] [CrossRef] [Google Scholar]
  145. Smartt, S. J. 2015, PASA, 32, e016 [NASA ADS] [CrossRef] [Google Scholar]
  146. Storey, P. J., & Hummer, D. G. 1995, MNRAS, 272, 41 [NASA ADS] [CrossRef] [Google Scholar]
  147. Tammann, G., Loeffler, W., & Schroeder, A. 1994, ApJS, 92, 487 [NASA ADS] [CrossRef] [Google Scholar]
  148. Thompson, T. A., Quataert, E., & Murray, N. 2009, MNRAS, 397, 1410 [NASA ADS] [CrossRef] [Google Scholar]
  149. Tüllmann, R., Gaetz, T., Plucinsky, P., et al. 2011, ApJS, 193, 31 [CrossRef] [Google Scholar]
  150. Vazdekis, A., Koleva, M., Ricciardelli, E., Röck, B., & Falcón-Barroso, J. 2016, MNRAS, 463, 3409 [Google Scholar]
  151. Veilleux, S., & Osterbrock, D. E. 1987, ApJS, 63, 295 [Google Scholar]
  152. Verner, E. M., Verner, D. A., Baldwin, J. A., Ferland, G. J., & Martin, P. G. 2000, ApJ, 543, 831 [NASA ADS] [CrossRef] [Google Scholar]
  153. Vicens-Mouret, S., Drissen, L., Robert, C., et al. 2023, MNRAS, 524, 3623 [NASA ADS] [CrossRef] [Google Scholar]
  154. Vučetić, M., Arbutina, B., & Uroršević, D. 2015, MNRAS, 446, 943 [CrossRef] [Google Scholar]
  155. White, R. L., & Long, K. S. 1991, ApJ, 373, 543 [NASA ADS] [CrossRef] [Google Scholar]
  156. White, R. L., Long, K. S., Becker, R. H., et al. 2019, ApJS, 241, 37 [NASA ADS] [CrossRef] [Google Scholar]
  157. Williams, B. J., & Temim, T. 2016, in Handbook of Supernovae (Living Reference Work Entry, 1 [Google Scholar]
  158. Williams, B. F., Peterson, S., Murphy, J., et al. 2014, ApJ, 791, 105 [NASA ADS] [CrossRef] [Google Scholar]
  159. Williams, B. F., Hillis, T. J., Blair, W. P., et al. 2019, ApJ, 881, 54 [NASA ADS] [CrossRef] [Google Scholar]
  160. Winkler, P., Smith, R., Points, S., Team, M., et al. 2015, Fifty Years of Wide Field Studies in the Southern Hemisphere: Resolved Stellar Populations of the Galactic Bulge and Magellanic Clouds, 491, 343 [NASA ADS] [Google Scholar]
  161. Winkler, P. F., Blair, W. P., & Long, K. S. 2017, ApJ, 839, 83 [Google Scholar]
  162. Winkler, P. F., Coffin, S. C., Blair, W. P., Long, K. S., & Kuntz, K. D. 2021, ApJ, 908, 80 [NASA ADS] [CrossRef] [Google Scholar]
  163. Winkler, P. F., Long, K. S., & Blair, W. P. 2023, AJ, 959, 62 [NASA ADS] [Google Scholar]
  164. Zangrandi, F., Jurk, K., Sasaki, M., et al. 2024, A&A, submitted [arXiv:2401.17307] [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  165. Zhang, K., Yan, R., Bundy, K., et al. 2017, MNRAS, 466, 3217 [NASA ADS] [CrossRef] [Google Scholar]
  166. Zurita, A., Rozas, M., & Beckman, J. E. 2000, A&A, 363, 9 [Google Scholar]
2

These regions do not have a definitive classification in Groves et al. (2023). This is principally because they are composite or shock-dominated objects in the BPT diagrams (e.g., AGN outflows, SNRs, planetary nebulae), although in some cases they simply have low S/N in one or more of the key diagnostic lines.

3

More specifically, the ratio traces the electron density in the post-shock [S II] recombination zone, and is proportional to the preshock density (Dopita et al. 2016).

All Tables

Table 1

Properties of PHANGS-MUSE galaxies used in this work.

In the text
Table 2

Number of identified SNRs and SNR candidates in 19 galaxies.

In the text
Table 3

Parent sample catalog generated in this work (see text for details).

In the text
Table A.1

Environmental masks applied to each criterion for 19 galaxies.

In the text
Table C.1

Fitted emission lines and masked sky lines while running DAP.

In the text
Table D.1

Identified objects in M83 (part).

In the text

All Figures

thumbnail Fig. 1

Zoomed-in image on four regions, showing (from left to right) log Hα flux (10−20 cm−2ergs−1), [S II]/Hα, [S II]/Hα residual, [O I]/Hα residual, and [S II] velocity dispersion maps and log [S II] flux. Objects identified by our method as SNRs and SNR candidates are marked with black solid ellipses. Objects identified in the Groves et al. (2023) nebular catalog as H II regions are enclosed by red solid lines, while regions unclassified in that catalog are outlined with blue solid lines. Region contours in the first column are from Groves et al. (2023). The second column shows region contours identified by astrodendro without any photoionization correction. The remaining columns show regions identified by astrodendro after using our final selection methods. It is apparent that many of the regions in the second column are spurious.
In the text
thumbnail Fig. 2

Log([S II]/Hα) as a function of logarithm Hα surface brightness for all pixels in each galaxy. The red solid line is the fitted correlation for each galaxy, with slope and intercept values given in the legend of each panel. The black solid line is representative of the DIG relation, as measured between 0.40 Rmax (maximum radial coverage) and 0.60 Rmax for each galaxy from Belfiore et al. (2022). The gray area indicates the 3σ range of this relation. SNRs typically lie above the fitted lines and show a looping structure toward the upper right. The horizontal black dashed line indicates the typical value of [S II]/Hα = 0.4 used to select SNRs. Many low Hα surface brightness pixels associated with the DIG have [S II]/Hα ≥ 0.4.
In the text
thumbnail Fig. 3

log([O I]/Hα) as a function of log Hα surface brightness for all pixels in all 19 galaxies. The red solid line is the fitted global correlation for all galaxies, with slope and intercept provided in the lower right corner of this plot. The black solid line is representative of the DIG relation, as measured between 0.40 Rmax (maximum radial coverage) and 0.60 Rmax for each galaxy from Belfiore et al. (2022). The gray area indicates the 3σ range of this relation. The black dashed line indicates the theoretical value of [O I]/Hα = 0.017 to select SNRs (Kopsacheili et al. 2020). The horizontal blue dotted line indicates the empirical value of [O I]/Hα = 0.1 used select SNRs. For NGC 1672, the individual pixels lying above the majority come from an SNe that just happened several days before our observation.
In the text
thumbnail Fig. 4

BPT diagrams with demarcation line from Kewley et al. (2006) with an indication of BPT distances as defined as the distance to the red solid extreme star-burst lines. The plus sign (+) indicates positive distance where shock-heating dominates, while the minus sign (-) indicates negative distance where photoionization dominates. The top left and top right plots show the distribution density of nebular regions in these 19 galaxies (Groves et al. 2023). In the bottom left and bottom right plots, SNRs and SNR candidates are indicated by red and blue dots, respectively. Nebular regions tend to occupy the photoionized region, while SNRs are concentrated in the shocked areas.
In the text
thumbnail Fig. 5

Zoomed-in map (3 kpc × 3 kpc) of the SNR population recovered in one of the galaxies in our sample, NGC 628, using our five selection criteria (described in Sec. 3). Panels from top left to bottom right show the log Hα emission, [O I]/Hα residual map, [S II]/Hα residual map, [S II]λ6716 velocity dispersion, BPT distance (OI-OIII) map, and the BPT distance (SII-OIII) map. Objects selected by each individual criterion are marked with black ellipses in the corresponding subplot and the locations of objects classified as SNR are marked with red boxes. Pixels with low S/N (≤5), and those corresponding to stellar sources, and the central region of the galaxy, have been masked (see text for details).
In the text
thumbnail Fig. 6

Example of what one SNR (in NGC 628) looks like when using the five criteria. This SNR is selected in the [O I]/Hα residual map, [S II]/Hα residual map, [S II]λ6716 velocity dispersion, BPT distance (OI-OIII), and BPT distance (SII-OIII) maps. The selection is marked with a black ellipse.
In the text
thumbnail Fig. 7

Parent sample of 2233 objects identified across the 19 PHANGS-MUSE galaxies, compared to the Hα emission. The galaxy name is in the upper left corner. SNRs are indicated by red dots while SNR candidates are in blue. The background shows the Hα map. Regions excluded from our search are outlined in solid black. Most galaxies have the center masked (see Appendix A for more details). We observe SNRs to be distributed across the full field of view, and good qualitative correspondence with Hα bright sites of star formation (see Appendix B for the rest of the galaxies).
In the text
thumbnail Fig. 8

Fraction of objects in the parent sample (left) and SNR sample (right) that are identified by specific pairwise combination of our five criteria ([O I]/Hα residual, [S II]/Hα residual, [S II]λ6716 velocity dispersion, BPT distance (OI-OIII), and BPT distance (SII-OIII)).
In the text
thumbnail Fig. 9

Spectra for four SNRs and an H II region in NGC 628. Top two: SNRs with high [S II]/Hα value > 0.8; middle two: SNRs with high-velocity dispersion > 80 km s−1; bottom one: H II region. Characteristic lines, e.g., [O I], [S II], [N II], and Hα are labeled in the plot. All SNRs have broader emission lines than the H II region, though this is most apparent in the two SNRs (NGC0628_17 and NGC0628_48) as they were selected to have particularly broad lines. The sky emission at [O I]λ6300Å, λ6363Å has been masked.
In the text
thumbnail Fig. 10

Distribution of identified objects and H II regions density (in grayscale) in the [O I]/Hα versus [S II]/Hα plane. Panel (a): SNRs (red dots); Panel (b): SNR candidates (blue dots). For objects that overlap with a H II region, the symbol is an empty circle in light red or blue. Identified objects tend to occupy different spaces than H II regions in the [O I]/Hα versus [S II]/Hα plane. However, this tendency is more obvious for a SNR than a SNR candidate.
In the text
thumbnail Fig. 11

SNRs in M83. The locations are indicated of SNRs identified by Long et al. (2022) (black boxes) and those identified by applying methods used in this work (red circles for SNRs and blue circles for SNR candidates). All the objects we identify as SNRs are also identified as SNRs in Long et al. (2022). Overall, 77% of our parent sample were previously identified as SNRs in the literature.
In the text
thumbnail Fig. 12

[S II]/Hα (top) and [O I]/Hα (bottom) versus Hα luminosity for SNRs (red dots) compared with the H II regions density distribution (in grayscale). For SNRs that overlap with a H II region, the symbol is an empty circle (in light red). Thresholds of [O I]/Hα=0.1 and [S II]/Hα=0.4 (horizontal dash-dotted lines) are shown. SNRs overlapping with H II regions are selected by other criteria so they can still be identified as SNRs. Emission lines are not corrected for the photoionization.
In the text
thumbnail Fig. 13

Radial trends of SNRs and H II regions. r/reff is the ratio of SNR distance to the effective radius of the corresponding galaxy. Panel (a): [S II]/Hα changes with the distance to their corresponding galactic centers for SNRs (red dots) and H II regions (in grayscale). Panel (b): [N II]λ6583/Hα change with the distance to their corresponding galactic centers for SNRs (red dots) and H II regions (in grayscale). For SNRs that overlap with a H II region, the symbol is an empty circle in light red. The drop in the numbers at low r/reff is likely influenced by the mask in the centers of the galaxies.
In the text
thumbnail Fig. 14

Comparison of emission-line ratios and velocity dispersion for supernova remnants (SNRs) and H II regions. Panel (a): Line ratios of [S II] versus velocity dispersion of SNRs and H II regions. The top axis shows the FWHM. The Hα velocity dispersion changes with the [S II]λ6716/[S II]λ6730 ratio for SNRs (red dots) and H II regions (in grayscale). The ratios above 1.4 are not physical and are due to noise. The horizontal dash-dotted lines indicate different electron densities with the corresponding value next to it. Panel (b): Velocity dispersion of Hα changes with [S II]λ6716,6730/Hα ratio for SNRs (red dots) and H II regions (in grayscale). For SNRs that overlap with a H II region, the symbol is an empty circle in light red. These two criteria work effectively together in distinguishing between SNRs and H II regions.
In the text
thumbnail Fig. 15

Comparison of supernova remnant properties with MAPPINGS models. Left: Observed [O III]λ5006/Hβ versus [N II]λ6583/Hα for SNRs (red dots). Right: Observed [N II]λ6583/Hα versus [S II]/Hα for SNRs (red dots). For SNRs that overlap with a H II region, the symbol is an empty light red circle. Background grids are MAPPINGS models with different metallicities. Our identified SNRs mostly lie between models for LMC and twice solar metallicities.
In the text
thumbnail Fig. 16

Relation between [S II]/Hα ratio and [N II]/Hα ratio for SNRs identified in four nearby galaxies from the literature (Long et al. 2018; Winkler et al. 2021; Long et al. 2022; Long et al. 2019) and 19 MUSE galaxies. SNRs are marked with red dots; SNRs that overlap with a H II region are marked as an empty circle in light red.
In the text
thumbnail Fig. 17

Relation between [S II]/Hα ratio and [N II]/Hα ratio for SNRs in 19 galaxies. SNRs in the same galaxy are marked in the same color and with the symbol. In the legend, from top to bottom, the stellar mass of the galaxy decreases, from 1011 M⊙ for NGC 1365 to 10941 M⊙ for NGC 5068.
In the text
thumbnail Fig. 18

Comparison between the estimated SNR number from the SFR and objects identified in this paper. The lower limit of the identified object is the number of SNRs and the upper limit is the number in the parent sample including SNR candidates. The lower limit of the estimated SNR number is given by the optically visible time of SNRs for 10 000 yr and the upper limit corresponds to 20 000 yr. The dashed black line is the one-to-one relation. The shaded area provides coverage from two times the recovery rate to 50% of the recovery rate of SNRs.
In the text
thumbnail Fig. 19

Number of SNe as a function of stellar mass (left) and SFR (right) considering the entire parent sample, colored by galaxy distance and considering galaxy integrated properties measured only within the MUSE field of view. These SN counts can also directly be converted to a SN frequency by assuming SNRs are visible in the optical for 10 000 yr. Lines of constant SN rate per unit mass (SNu) are also shown. The star markers represent the number for the parent sample while the circles are the number of SNRs. The line between them shows the range of possible true SNR numbers.
In the text
thumbnail Fig. 20

Probability distribution functions of [S II]/Hα for SNRs and SNR candidates that are isolated or overlap with H II region. More SNR candidates (bottom) overlap with H II regions (blended) and have lower [S II]/Hα values than SNRs (top).
In the text
thumbnail Fig. B.1

Same as in Fig. 7.
In the text
thumbnail Fig. B.2

Same as in Fig. 7.
In the text
thumbnail Fig. B.3

Same as in Fig. 7.
In the text
thumbnail Fig. B.4

Same as in Fig. 7
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.