# Direct Evidence of Two-component Ejecta in Supernova 2016gkg from Nebular Spectroscopy

Hanindyo Kuncarayakti, Gaston Folatelli, Keiichi Maeda, Luc Dessart, Anders Jerkstrand, Joseph P. Anderson, Kentaro Aoki, Melina C. Bersten, Lucía Ferrari, Lluís Galbany, et al.

► To cite this version:

Hanindyo Kuncarayakti, Gaston Folatelli, Keiichi Maeda, Luc Dessart, Anders Jerkstrand, et al. Direct Evidence of Two-component Ejecta in Supernova 2016gkg from Nebular Spectroscopy. The Astrophysical Journal, 2020, 902 (2), pp.139. 10.3847/1538-4357/abb4e7. hal-02953263

HAL Id: hal-02953263

https://hal.science/hal-02953263v1

Submitted on 7 Aug 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Direct Evidence of Two-component Ejecta in Supernova 2016gkg from Nebular Spectroscopy*

Hanindyo Kuncarayakti<sup>1,2</sup>①, Gastón Folatelli<sup>3,4,5</sup>①, Keiichi Maeda<sup>5,6</sup>①, Luc Dessart<sup>7</sup>①, Anders Jerkstrand<sup>8,9</sup>①,
Joseph P. Anderson<sup>10</sup>①, Kentaro Aoki<sup>11</sup>, Melina C. Bersten<sup>3,4,5</sup>①, Lucía Ferrari<sup>3</sup>, Lluís Galbany<sup>12</sup>①, Federico García<sup>13</sup>①,
Claudia P. Gutiérrez<sup>14</sup>①, Takashi Hattori<sup>11</sup>, Koji S. Kawabata<sup>15</sup>①, Timo Kravtsov<sup>1</sup>, Joseph D. Lyman<sup>16</sup>①, Seppo Mattila<sup>1</sup>①,
Felipe Olivares E.<sup>17</sup>①, Sebastián F. Sánchez<sup>18</sup>①, and Schuyler D. Van Dyk<sup>19</sup>①

¹ Tuorla Observatory, Department of Physics and Astronomy, FI-20014 University of Turku, Finland; hanindyo.kuncarayakti@utu.fi

² Finnish Centre for Astronomy with ESO (FINCA), FI-20014 University of Turku, Finland

³ Facultad de Ciencias Astronómicas y Geofísicas, Universidad Nacional de La Plata, Paseo del Bosque S/N, B1900FWA La Plata, Argentina

⁴ Instituto de Astrofísica de La Plata (IALP), CONICET, Argentina

⁵ Kavli Institute for the Physics and Mathematics of the Universe (WPI), The University of Tokyo, 5-1-5 Kashiwanoha, Kashiwa, Chiba 277-8583, Japan

⁶ Department of Astronomy, Graduate School of Science, Kyoto University, Sakyo-ku, Kyoto 606-8502, Japan

⁷ Institut d'Astrophysique de Paris, CNRS-Sorbonne Université, 98 bis boulevard Arago, F-75014 Paris, France

⁸ Max-Planck-Institut für Astrophysik, Karl-Schwarzschild Str 1, D-85748 Garching, Germany

⁹ Department of Astronomy, Stockholm University, The Oskar Klein Centre, AlbaNova, SE-106 91 Stockholm, Sweden

¹⁰ European Southern Observatory, Alonso de Córdova 3107, Casilla 19, Santiago, Chile

¹¹ Subaru Telescope, National Astronomical Observatory of Japan, 650 North A'ohoku Place, Hilo, HI 96720, USA

¹² Departamento de Física Teórica y del Cosmo, Universidad de Granada, E-18071 Granada, Spain

¹³ Kapteyn Astronomical Institute, University of Groningen, PO BOX 800, Groningen NL-9700 AV, The Netherlands

¹⁴ Department of Physics and Astronomy, University of Southampton, Southampton, SO17 1BJ, UK

¹⁵ Hiroshima Astrophysical Science Center, Hiroshima University, 1-3-1 Kagamiyama, Higashi-Hiroshima, Hiroshima 739-8526, Japan

¹⁶ Department of Physics, University of Warwick, Coventry CV4 7AL, UK

¹⁷ Instituto de Astronomía y Ciencias Planetarias, Universidad de Atacama, Copayapu 485, Copiapó, Chile

¹⁸ Instituto de Astronomía, Universidad Nacional Autónoma de México, Circuito Exterior, Ciudad Universitaria, Ciudad de México 04510, Mexico

¹⁹ Caltech/IPAC, Pasadena, CA 91125, USA

*Received 2020 June 17; revised 2020 August 21; accepted 2020 August 26; published 2020 October 21*

## Abstract

Spectral observations of the type-Ib supernova (SN) 2016gkg at 300–800 days are reported. The spectra show nebular characteristics, revealing emission from the progenitor star’s metal-rich core and providing clues to the kinematics and physical conditions of the explosion. The nebular spectra are dominated by emission lines of [O I] λλ6300, 6364 and [Ca II] λλ7292, 7324. Other notable, albeit weaker, emission lines include Mg I] λ4571, [Fe II] λ7155, O I λ7774, Ca II triplet, and a broad, boxy feature at the location of Hα. Unlike in other stripped-envelope SNe, the [O I] doublet is clearly resolved due to the presence of strong narrow components. The doublet shows an unprecedented emission line profile consisting of at least three components for each [O I]λ6300, 6364 line: a broad component (width ~2000 km s⁻¹), and a pair of narrow blue and red components (width ~300 km s⁻¹) mirrored against the rest velocity. The narrow component appears also in other lines, and is conspicuous in [O I]. This indicates the presence of multiple distinct kinematic components of material at low and high velocities. The low-velocity components are likely to be produced by a dense, slow-moving emitting region near the center, while the broad components are emitted over a larger volume. These observations suggest an asymmetric explosion, supporting the idea of two-component ejecta that influence the resulting late-time spectra and light curves. SN 2016gkg thus presents striking evidence for significant asymmetry in a standard-energy SN explosion. The presence of material at low velocity, which is not predicted in 1D simulations, emphasizes the importance of multidimensional explosion modeling of SNe.

*Unified Astronomy Thesaurus concepts: Core-collapse supernovae (304)*

*Supporting material: data behind figure*

## 1. Introduction

Type-Ilb supernovae (SNe IIb) exhibit a change of spectral appearance from showing hydrogen lines at early times, to a typical type-Ib SN (helium-rich, hydrogen-poor), as first exemplified by

SN 1987K (Filippenko 1988) and SN 1993J (Woosley et al. 1994). The progenitor stars of this subclass are thought to retain a thin layer of hydrogen at the time of the explosion, giving rise to the SN II → SN Ib spectral transition.

SN 2016gkg is arguably one of the best-observed members of this subclass, with the detection of the progenitor star in pre-explosion archival images (Kilpatrick et al. 2017; Tartaglia et al. 2017) and the extraordinary discovery of the optical light breakout within two hours of core collapse (Bersten et al. 2018). The early light curve of SN 2016gkg following the shock breakout is exceptionally well sampled, covering the peak associated with the shock breakout cooling, as well as the subsequent main peak resulting from radioactive heating that is frequently observed in other SNe. The post shock cooling peak

<footnote>* Based on observations collected at the European Organisation for Astronomical Research in the Southern Hemisphere under ESO program 0100.D-0341; the Gemini Observatory, which is operated by the Association of Universities for Research in Astronomy, Inc., under a cooperative agreement with the NSF on behalf of the Gemini partnership: the National Science Foundation (United States), National Research Council (Canada), CONICYT (Chile), Ministerio da Ciência, Tecnologia e Inovação (Brazil), and Korea Astronomy and Space Science Institute (Republic of Korea), programs GS-2017B-Q-65, GS-2018B-FT-107 (Fast Turnaround); and the Subaru Telescope, which is operated by the National Astronomical Observatory of Japan, program S17B-162S.</footnote>
of the light curve is critical to constrain the progenitor radius and the mass of its extended envelope (see, e.g., Waxman & Katz 2017, for a general review, and references therein). With this method, the radius of the progenitor of SN 2016gkg has been constrained to be between ~200–300 $R_⊙$ and the extended envelope mass in the order of $10^{-2} M_⊙$ (Piro et al. 2017; Bersten et al. 2018). Arcavi et al. (2017) reported a similar envelope mass estimate, albeit with a smaller radius of 40–150 $R_⊙$. The detection of the pre-explosion progenitor star candidate in Hubble Space Telescope (HST) archival images points to an initial mass of 15–20 $M_⊙$ (Kilpatrick et al. 2017; Tartaglia et al. 2017). Bersten et al. (2018) further improved the identification and characterization of the progenitor star in the HST images, and concluded that the progenitor may have been a $19.5 M_⊙$ star in a binary system via binary star evolution modeling.

Here we report late-time spectral observations of SN 2016gkg, resulting in spectra between +300 and 800 days after the light curve maximum. At this late stage, the SN ejecta have expanded and become optically thin, rendering the emission from the deeper layers and inner core visible. As will be demonstrated in the following sections, this condition enables the detection of very slow-moving material, in addition to the typically seen fast-moving material, in the ejecta of SN 2016gkg, which uniquely characterizes the object. The observations and data reduction are described in the next section, followed by discussions of the results and interpretations. Herein light-curve maximum refers to the peak in the light curve that resulted from radioactive heating, which occurred on MJD 57668.4 or 2016 October 7.4 (UT), in B and V bands (Bersten et al. 2018).

## 2. Observations and Data Reductions

SN 2016gkg was observed with the Gemini South telescope and the GMOS-S instrument (Hook et al. 2004; Gimeno et al. 2016), on 2017 August 20 and 26 (UT dates are used throughout). At the median date of observation, August 23, the phase corresponds to +320 days after light-curve maximum. The observations were conducted using the longslit mode with the R400 grism, with a 4 × 4 binning and a total exposure time of 5.3 hr. Another GMOS-S longslit observation was performed on 2018 November 10 (+764 d) through the Gemini Fast Turnaround program, using the R400 grism and 4 × 2 binning, with a total exposure time of 3.0 hr. Spectral dithering was performed to eliminate cosmic rays and cover the gaps between the GMOS CCDs. Data reduction was done using the Gemini package in IRAF.<fnref n="20" /> All the spectra were wavelength-calibrated, and flux-calibrated using spectrophotometric standard stars.

On 2017 November 30, the field of SN 2016gkg was observed using the MUSE instrument (Bacon et al. 2014) at the ESO Very Large Telescope. MUSE was used in Wide Field Mode, with a total exposure time of 2400 s. SN 2016gkg was positioned at the center of the field of view and the observation captured the 1' × 1' field surrounding SN 2016gkg in integral-field spectroscopy, which includes the nebular spectrum of the object at +419 days past maximum. The MUSE spaxel size is 0.2' × 0.2', corresponding to ~25 × 25 pc at the assumed

distance of 26.4 Mpc to SN 2016gkg (Bersten et al. 2018). The observations were part of the All-weather MUSe Supernova Integral-field Nearby Galaxies (AMUSING; see Galbany et al. 2016) survey. MUSE data reduction was performed using MUSE data reduction pipelines run under ESO Reflex (Freudling et al. 2013). The datacube was sky-subtracted using blank sky pointings, and corrected for atmospheric effects using the Zurich Atmospheric Package (ZAP; Soto et al. 2016). The measured image quality in the final datacube is typically 0.8 FWHM across the wavelengths, and this value is used as the radius of the circular aperture for 1D spectrum extraction of SN 2016gkg with QFitsView (Ott 2012).

SN 2016gkg was also observed using FOCAS (Kashikawa et al. 2002) at the Subaru telescope on 2018 January 1 (+451 d) in the 0.8' off-center slit mode with the B300 grism + Y47 order-cut filter setting. The total exposure time was 2400 s, taken in two exposures of 1200 s each, with a slight positional shift (AB pattern). The standard data reduction procedure, which includes wavelength and flux calibration was done using IRAF, with cosmic-ray removal using L.A.Cosmic (van Dokkum 2001).

The final reduced spectra from these observations were then corrected for the redshift of $z = 0.0049$ as measured from the underlying H II region emission lines, that agrees with the reported redshift of the host galaxy NGC 613 (Meyer et al. 2004). This redshift is assumed to be the rest reference of the SN throughout the paper.

The spectral resolutions measured from the sky emission lines in the spectra are ~3 Å FWHM for MUSE (corresponding to resolving power $R \approx 2200$ and velocity resolution $\Delta v \approx 140 \text{ km s}^{-1}$ at H$\alpha$ wavelength), and ~10 Å for GMOS and FOCAS ($R \approx 700$ and $\Delta v \approx 430 \text{ km s}^{-1}$).

## 3. Results and Discussions

### 3.1. Spectrum and Line Identifications

The four nebular spectra of SN 2016gkg at +320, +419, +451, and +764 days are presented in Figure 1. These epochs are considerably late in the SN evolution as most stripped-envelope (SE) SNe enter the nebular phase as early as a few months after maximum. The majority of nebular spectra available in the literature correspond to phases earlier than one year post maximum light, while only a few SESNe have been observed at the phase of one year or later. Typical for a hydrogen-poor core-collapse SN observed at the nebular late-phase stage, the spectra of SN 2016gkg show broad nebular emission lines, dominated by the [O I] $\lambda\lambda6300,6364$ doublet and [Ca II] $\lambda\lambda7292,7324$ (e.g., Maeda et al. 2008; Taubenberger et al. 2009; Fang et al. 2019). The [O I] $\lambda\lambda6300,6364$ doublet shows an extraordinary line profile with strong narrow emission components. These narrow components (FWHM ~ 300 km s$^{-1}$) are fully resolved with the highest resolution spectrum from MUSE ($\Delta v \approx 140 \text{ km s}^{-1}$). The [Ca II] $\lambda\lambda7292,7324$ line exhibit a sloping red wing (see Figure 1, and a zoom-in on Figure 2), possibly arising from contamination by [Ni II] $\lambda7377$ (Jerkstrand et al. 2015), or opacity and scattering effects (Jerkstrand 2017). Other weaker emission lines such as Mg I], He I, Na I, Fe II are also present, as well as a few lines from permitted transitions of O I and Ca II. The spectral lines are largely symmetric and no asymmetries that might have resulted from newly formed dust in the ejecta are detected (see, e.g., Meikle et al. 2011).

<footnote n="20">²⁰ IRAF is distributed by the National Optical Astronomy Observatory, which is operated by the Association of Universities for Research in Astronomy, Inc., under cooperative agreement with the National Science Foundation.</footnote>
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Nebular spectra of SN 2016gkg. The species responsible for the emission lines are shown, with their rest wavelengths indicated by vertical dashed lines. Narrow $H\alpha$+[N II] are likely lines arising in the interstellar medium (ISM); the telluric absorption region is indicated with the Earth symbol (⊕, ~7600 Å). Fluxes have been scaled to allow comparison between the spectra. A zoom-in on the [O I] doublet is shown on the right-hand panel inset. Line identification throughout the paper follows that of Gröningsson et al. (2008) and Jerkstrand et al. (2015). The spectra are available as the data behind the figure (DbF).
(The data used to create this figure are available.)</figcaption></figure>

The narrow emission lines are checked by constructing line maps from the MUSE datacube at their respective wavelengths. All the detected lines in the SN spectrum are confirmed to be emitted by the SN, as opposed to having originated from the diffuse background or noise, from the appearance of the SN as a point source in each linemap (Figure 3). This is one advantage of integral-field spectroscopy.

There is a relatively weak bump with a broad, boxy profile at the wavelength of $H\alpha$ in the spectra (Figure 4), typical of SNe IIb at late times (see, e.g., Fang et al. 2019). The bump is visible in all epochs. This bump has been suggested to be contributed by [N II] $\lambda\lambda$6548,6583 (Jerkstrand et al. 2015), although CSM interaction may also give rise to a boxy, flat-topped $H\alpha$ emission as seen in the case of the type-IIb SN 1993J (Matheson et al. 2000). A slight blueshift in the SN emission line peaks with respect to the rest wavelengths (~ - 500 km s⁻¹) can be discerned (also see Figure 2). This phenomenon appears to be commonplace in SESNe, and may be explained by various factors such as ejecta geometry and opacity effects (Taubenberger et al. 2009; also see radiative transfer models presented in Dessart & Hillier 2005).

There is generally little evolution in the spectrum of SN 2016gkg between +320 to +450 days and the narrow+broad

line profile persists (Figure 1). Despite the different spectral resolutions, all spectra visibly display nearly identical sets of emission lines. The last spectrum at +764 days shows that most of the weaker lines have disappeared beyond our detection limit, leaving mostly the [O I] and [Ca II] doublets and Mg I]. Narrow $H\alpha$ emission line accompanied by [N II] $\lambda$6584 and the [S II] $\lambda\lambda$6717,6731 doublet are present (Figure 4), but these are likely to be lines originating from the background interstellar medium (ISM). These lines are better subtracted in the other spectra, and their line maps do not reveal a point source at the SN position (Figure 3).

## 3.2. Line Profiles and Comparison to Other SNe

The [O I] $\lambda\lambda$6300,6364 doublet in the spectra of SN 2016gkg shows a narrow component superposed on a broad component, and while less prominent, such a composite profile is also seen in [Ca II] $\lambda\lambda$7292,7324 (Figure 2). Comparing with available spectra in the literature, the line profile containing clear broad +narrow components of [O I] and [Ca II] in SN 2016gkg is unprecedented. Most SNe show line profiles of either a single Gaussian, double peak, or asymmetric/multipeak, while those with a narrow core, like that seen in SN 2016gkg, are relatively
<figure><img src="image_4.png" /><figcaption>**Figure 2.** [O I] λ6300 (teal), Mg I] λ4571 (magenta), and [Ca II] λ7292 (orange) lines in velocity space. The Mg I] line is from the +320 day spectrum, while the others are from the +419 day one. Black vertical dashed line indicates the zero velocity for each individual line, and colored dashed lines correspond to the rest wavelengths of [O I] λ6364 and [Ca II] λ7324, and additionally [Fe II] λ7155 and possibly [Ar III] λ7135, seen near [Ca II].</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Figure 3.** Images of the SN environment in the host galaxy, generated from the MUSE datacube, spanning around 1 arcmin (300 spaxels) on each side. The SN position is indicated with crosshairs, and images are oriented north-up and east-left. The host galaxy center is off the field, due north. Clockwise from top left: integrated (white) light, and narrow emission line maps of Hα (6563 Å), [O I]λ6364 (at 6353 Å), and [O I]λ6300 (at 6289 Å). Emission line maps are continuum subtracted.</figcaption></figure>

rare (Maeda et al. 2008; Taubenberger et al. 2009) and the flux associated with the narrow core is never seen dominating over the broad component.

<figure><img src="image_6.png" /><figcaption>**Figure 4.** Spectra of SN 2016gkg around the Hα region. Dashed horizontal lines indicate the average continuum level between 6850 and 6900 Å. Likely host galaxy ISM lines of Hα, [N II] λ6584, and [S II] λ6717,6731 are indicated with vertical dashed-dotted lines. These lines are clearly visible in emission in the +764 day spectrum, and appear as background subtraction artifacts in the +451 day spectrum.</figcaption></figure>

Furthermore, SN 2016gkg displays a very prominent narrow component clearly in each line of the [O I] doublet at λλ6300 and 6364, which is not seen in other SNe thus far published. The type-Ib SN 2009jf (Sahu et al. 2011; Valenti et al. 2011) is possibly the closest analog to SN 2016gkg in terms of narrow [O I] emission as its nebular spectra show a strong narrow component superposed on the broad [O I] line, although without clear doublet identifications. The single narrow line in SN 2009jf is accompanied by weaker peaks in both the blue and red directions, and altogether the [O I] structure is interpreted as having originated from a number of clumps superposed on the bulk high-velocity ejecta along the line of sight. Figure 5 shows the comparison between the spectra of SN 2016gkg and SN 2009jf, at photospheric and nebular phases. While the narrow [O I] λ6300 is clearly evident in SN 2009jf, the λ6364 component is considerably weaker and ambiguous with the other smaller peaks. The narrow components in [Ca II] and [Fe II] are not visible. As the +361 day spectrum of SN 2009jf and the +320 day spectrum of SN 2016gkg have similar resolutions, the presence of those narrow components in SN 2009jf can therefore be ruled out.

The narrow [O I] line in SN 2016gkg possibly appeared as early as three months after the light curve peak. Figure 5 shows the ~+90 and +400 day spectra of SN 2016gkg, as compared to SN 2009jf, which also shows a narrow [O I] line core. Similarly, SN 2009jf also shows the emergence of this narrow line at around three months post-maximum. At this epoch, there is also an indication of the narrow core in O I λ7774 in SN 2016gkg; however, this is not clearly seen in SN 2009jf. The narrow cores in the other lines (e.g., [Ca II] and [Fe II]) do not appear to be present in the early phase. Note that in SESNe the ejecta become optically thin in the continuum relatively early—generally within a few weeks after the light-curve maximum. Therefore, conditions are nebular much sooner (i.e., the full ejecta is visible); forbidden lines simply require a longer time to become stronger.

Matheson et al. (2000) noted that SN 1993J showed small-scale structures at the top of the O and Mg lines, which were interpreted as multiple clumps with distinct velocities. The Ca
<figure><img src="image_7.png" /><figcaption>**Figure 5.** (Left panel) Emergence of the [O I] narrow-line core in SN 2016gkg and SN 2009jf. Early-phase spectra of the two SNe are from Bersten et al. (2018) and Valenti et al. (2011), respectively. The +419 day spectrum of SN 2016gkg is plotted with the ordinate in log scale due to the extreme intensity of the narrow component. (Right panels) Comparison of the +361 day spectrum of SN 2009jf with the spectra of SN 2016gkg at +320 days and +419 days, around the [O I] and [Ca II] doublets. All spectra are plotted in linear flux scale. The resolutions of the spectra at [O I] are indicated.</figcaption></figure>

lines of SN 1993J, however, do not show similar profiles, and therefore it was interpreted that the Ca distribution in the ejecta of SN 1993J is more uniform. For comparison, SN 2009jf has been claimed to have a more mixed ejecta due to the similar profile seen in the Mg I] and [Ca II] lines modified with an additional doublet component. This is likely to be the case for SN 2016gkg as well, since the similar broad line + narrow core profile is also clearly seen in [Ca II] (Figure 2). Ca is the strongest coolant in SE SNe, and thus [Ca II] emission may arise from all regions where Ca is present. In contrast, [O I] may only emit from the layers where O is the most abundant species.

Careful examination of the +419 day spectrum of SN 2016gkg shows that a blueshift is present in most of the emission lines (Figures 2 and 6). In addition to the two most dominant doublets [O I] λλ6300,6364 and [Ca II] λλ7292,7324, it is also present in the other species such as Mg I] λ4571 (Figure 2), [O I] λ5577, O I λ7774, and O I λ9263, and also in the Ca II lines of λ8542 and λ8662 (Figure 6). These lines do not appear to have the broad components. As the [O I] λ5577, O I λ7774, O I λ9263, and Ca II λ8542 and λ8662 lines are associated with high-density regions, this suggests that the narrow-line components originated from the inner ejecta. Accompanying the blue component, the red component is visible also in O I λ7774, Ca II λ8542, and Ca II λ8662 (Figure 6), and as previously discussed in the [O I] λλ6300,6364 doublet. It is possible that such structure is also present in [Ca II] λλ 7292,7324; however, since the narrow red components are

considerably fainter relative to the blue ones they are not readily visible.

The nebular spectra of SESNe usually show broad emission lines of [O I] and [Ca II], with FWHM of around a few 1000 km s⁻¹ (Taubenberger et al. 2009). This occurs because the energy deposition due to γ-rays prevails nonlocally and thus excites the emission lines over a large volume (i.e., velocities, given homologous expansion). SN 2016gkg shows this typical line width characteristic coexisting with narrow emission lines. To produce narrow emission lines, a certain amount of material at low velocities needs to be powered by some source.

We note here that the identification of the narrow components in SN 2016gkg is not because of the relatively high spectral resolution of MUSE (R ~ 2200) compared to the typical resolutions (R ~ 500–1000) at which nebular spectra of SNe are usually acquired. In fact, they are clearly evident in the low-resolution GMOS and FOCAS spectra (R ~ 700) as well (Figures 1 and 5). While the MUSE resolution is not needed to identify the narrow lines, it is clearly required in order to resolve them. Higher resolution does help in detecting narrow emission lines; regardless, such lines as in SN 2016gkg would have been detected if they were present in the other SNe observed thus far. Therefore the absence of such a narrow component in other SNe cannot be attributed to the low spectral resolution, and the fact that it is detectable and prominent in SN 2016gkg suggests that there must be a unique circumstance in this SN that leads to the emergence of the line.
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Close-up of narrow emission lines, in velocity space, in the MUSE spectrum at +419 d. Dashed vertical line indicates zero velocity for each of the indicated lines. The lines are scaled in flux if indicated (×2), and in the case of [O I] λ6364 (plotted in red) it is scaled to match the λ6300 component (plotted in blue). The profiles of these mirrored narrow components at ±500 km s⁻¹ match each other well.</figcaption></figure>

To measure the individual components comprising the [O I] doublet complex, spectral line decomposition was performed. This was done using MPFIT (Markwardt 2009), via PAN: Peak Analysis application.<fnref n="21" /> The MUSE spectrum at +419 days, corrected for the host redshift, was used due to the superior spectral resolution compared to the other spectra. The fit assumed a linear background and multiple Gaussian profiles in order to reproduce the observed line profile. Examining the region around the [O I] line, it appears that there are three components in the [O I] λ6300 line: a broad component, a narrow component blueward (∼10 Å offset) of the rest wavelength of 6300 Å, and a narrow component redward of 6300 Å (also with ∼10 Å offset). A similar profile is seen in the weaker [O I] λ6364 line. Using the apparent line centers, peak intensities, and FWHM as initial guesses, we fitted six Gaussians and a linear background function to the [O I] doublet line profile, without fixing any parameters. As such, the line centers, widths, and intensities, are left as free parameters. The

<figure><img src="image_9.png" /><figcaption>**Figure 7.** Fitting of the [O I] λλ6300,6364 doublet. The observed spectrum is plotted by a black continuous line, while the fitted functions in long-dashed lines (red for the free-parameter fit and blue for the FWHM-tied fit). Short-dashed vertical lines indicate zero velocities for 6300 Å and 6364 Å rest wavelengths.</figcaption></figure>

result of the fit is shown in Figure 7 and the parameters are listed in Table 1.

Comparing with other stripped-envelope SNe at similar epochs, the broad component in SN 2016gkg appears to be narrower, while looking more similar to SNe II in both profile and width (Figure 8). Taubenberger et al. (2009) reported typical FWHM of the [O I] doublet from a single-component Gaussian fit for SNe Ib/IIb to be ~5200 km s⁻¹, while Silverman et al. (2017) noted that SNe II typically show [O I] half width at half maximum (HWHM) of 1000–1200 km s⁻¹ (FWHM ~ 2000–2400 km s⁻¹). The broad components of SN 2016gkg are observed to have FWHM of ~2000 km s⁻¹, while the narrow components’ FWHM are ~300 km s⁻¹. These are much greater than the instrument resolution of 140 km s⁻¹, and thus correcting for it would not significantly reduce the derived FWHM velocities. The doublet nature of the [O I] emission line is more easily discerned in SNe II, whereas in SESNe the doublet is usually blended due to the higher expansion velocity. SN 2016gkg, while belonging to the SESN group, clearly displays the doublet nature of [O I]. Difference in the broad component probably means that SN 2016gkg differs from other SESNe in terms of different inner ejecta structure, which may be characterized by higher density overall, or more ⁵⁶Ni at low velocity.

### 3.3. Properties of the Oxygen-rich Region

Here we analyze the core oxygen properties based on the [O I] doublet complex, using the spectrum at +419 days. The measured line ratio between the 6300 Å and 6364 Å lines yields 3.6 ± 0.1:1 for the broad components, while the ratio is 2.8 ± 0.1:1 for the narrow components (blue+red combined), suggesting a relatively denser narrow-line emitting region as the ratio should be close to 3:1 in the optically thin regime (Li & McCray 1992; Jerkstrand 2017), and 1:1 in the optically thick regime. The [O I] doublet profile in late-time SESNe

<footnote n="21">²¹ https://www.ncnr.nist.gov/staff/dimeo/panweb/pan.html</footnote>
<table><caption>**Table 1**

Resulting Parameters of the [O I] Doublet Fit, +419 Day Spectrum</caption><thead><tr><th>Component</th><th>λ<sub>center</sub><br>(Å)</th><th>v<sub>center</sub><br>(km s<sup>-1</sup>)</th><th>FWHM<br>(Å)</th><th>v<sub>FWHM</sub><br>(km s<sup>-1</sup>)</th><th>Normalized Flux</th></tr></thead><tbody><tr><td>Broad-6300</td><td>6291</td><td>-428</td><td>50</td><td>2381</td><td>1.00<sup>a</sup></td></tr><tr><td>Broad-6364</td><td>6351</td><td>-613</td><td>34</td><td>1603</td><td>0.28</td></tr><tr><td>Narrow-6300, blue</td><td>6289</td><td>-524</td><td>8</td><td>381</td><td>0.36</td></tr><tr><td>Narrow-6300, red</td><td>6311</td><td>+524</td><td>6</td><td>286</td><td>0.03</td></tr><tr><td>Narrow-6364, blue</td><td>6353</td><td>-471</td><td>7</td><td>330</td><td>0.12</td></tr><tr><td>Narrow-6364, red</td><td>6375</td><td>+519</td><td>8</td><td>377</td><td>0.03</td></tr></tbody></table>

Notes. The spectral resolution is ~3 Å ≈140 km s⁻¹; nominal fitting uncertainties are typically 0.1 Å ≈5 km s⁻¹.

a. Equals $2.36 \times 10^{-15}$ erg s⁻¹ cm⁻² Å⁻¹.

has been discussed in the literature (Maeda et al. 2008; Taubenberger et al. 2009; Modjaz et al. 2008), in particular whether the double-peak features result from the λλ6300,6364 doublet components or distinct velocity components in the ejecta (Milisavljevic et al. 2010). Here, in SN 2016gkg, the [O I] line profile with the 64 Å peak separation and close to 3:1 flux ratio unambiguously point to the doublet origin.

The line ratio of λ6300 and λ6364 (combining the broad and narrow components) can be compared with theoretical calculations to derive the density of the [O I] emitting region. Using Figure 6 of Leibundgut et al. (1991), these line ratios yield an estimated [O I] density of <3 × 10⁸ cm⁻³. The [O I] λ5577 line is detected in the spectrum, and thus can be used to estimate the temperature of the emitting region (e.g., Elmhamdi 2011). In SN 2016gkg, the observed [O I]λ5577/(λλ6300 + 6364) line ratio is 0.07. Using Equation (2) of Jerkstrand et al. (2014), and assuming β<sub>5577</sub>/β<sub>6300,6364</sub> = 1.5, the estimated temperature is then 3800 K. We adopt this temperature estimate and further derive the mass of the emitting oxygen from the total luminosity of the [O I] line (accounting for both the narrow and broad components), using Equation (3) of Jerkstrand et al. (2014) with β<sub>6300,6364</sub> = 0.5. Note that the λ5577 line falls out of local thermodynamic equilibrium (LTE) conditions after around 250 days in SNe IIP models, and after 150 days in SN IIb models (Jerkstrand et al. 2015), although in SN 2016gkg this appears to be delayed due to the lower expansion velocity. Even if λ5577 is in LTE but λ6300,6364 still in LTE the method will give an upper limit to the oxygen mass, which is still useful.

Assuming a distance of 26.4 Mpc and reddening A<sub>V</sub> = 0.053 mag (Bersten et al. 2018), the total oxygen mass is therefore estimated to be in the order of 0.3 M<sub>☉</sub>, which accounts for the ~0.2 and 0.1 M<sub>☉</sub> associated with the broad and narrow components, respectively, assuming that line flux corresponds to mass in the same way for both components.

The core oxygen mass is related to the initial mass of the progenitor star. Following the method presented in Kuncarayakti et al. (2015), the initial mass of the progenitor of the SN 2016gkg is estimated by comparing the core oxygen mass to the theoretical yields from Nomoto et al. (1997) and Limongi & Chieffi (2003). The resulting estimate for the progenitor initial mass with 0.3 M<sub>☉</sub> of core oxygen mass is 13–15 M<sub>☉</sub>. This estimate is effectively a lower limit since some of the oxygen may not be emissive. In such a case, for example assuming that half of the oxygen is not emissive, doubling the oxygen mass estimate of SN 2016gkg would bring the initial mass estimate to around 16–17 M<sub>☉</sub>. In comparison, Bersten et al. (2018) suggested that the progenitor of SN 2016gkg might have been a 19.5 M<sub>☉</sub> primary star in a close binary system, based on their HST progenitor detection and binary

star evolution modeling. Whereas, Kilpatrick et al. (2017) suggested the mass to be around 15 M<sub>☉</sub>, and Tartaglia et al. (2017) estimated 15–20 M<sub>☉</sub>, from photometric analyses of the progenitor star detected in pre-explosion archival HST images.

We further estimate the oxygen filling factor in the ejecta following the method of Leibundgut et al. (1991). At 400 d, assuming an isotropic expansion with constant velocity of 2000 km s⁻¹, the ejecta would fill a sphere with a volume of 1.4 × 10⁴⁸ cm³. If this sphere was uniformly filled with oxygen with an upper limit of mean density of 3 × 10⁸ atoms cm⁻³ as previously estimated, a total of around 5.5 M<sub>☉</sub> of oxygen would be contained. As the estimated oxygen mass of 0.3 M<sub>☉</sub> is far below this total oxygen mass (~5%; however, this is a lower limit), the interpretation would be that the oxygen filling factor in the ejecta is probably low, or that a significant fraction of the oxygen is located at a density lower than 3 × 10⁸ atoms cm⁻³.

### 3.4. Two-component Ejecta and Explosion Geometry

In SN 2016gkg, the narrow [O I] lines are found to be ~300 km s⁻¹ wide, indicating that the emitting region is expanding at least five to six times slower compared to the broad component (~2000 km s⁻¹), and have a combined [O I] flux (λλ6300 + 6364) about 40% of the broad component (see Table 1). If the emitting region is assumed to be of constant density, the density is scaled by the mass and the inverse cube of velocity, where the mass of the material is proportional to the emission line flux. Considering the line flux of the narrow [O I] component is about 40% of the broad component, and the expansion velocity 20% of the broad component, the density would be 0.4/(0.2)³ = 50 times higher compared to the broad component. In general (i.e., in central explosions such as core-collapse SNe), density drops with increasing velocity and thus the broad component of [O I] should form at lower density than the narrow component. If the absorbing material is not dense, the power absorbed by the low-velocity material would raise the ionization significantly. Furthermore, a low-density inner region would not efficiently absorb power.

The above illustration is consistent with an interpretation of line profiles seen in other SNe, where a narrow core on top of a broader component indicates enhanced central density. Further, the central dense part of the ejecta may take the form of a torus-like structure, which when viewed from the side gives rise to the double-peak profile symmetric around the rest wavelength (Taubenberger et al. 2009). In principle, a bipolar structure pointing toward and away from the observer may also give rise to the blue and red pair of narrow emission lines. The line displacement from rest wavelength immediately suggests asymmetric ejecta, as in the spherical case no offset is expected.
<figure><img src="image_10.png" /><figcaption>**Figure 8.** (Left panel) Spectra SN 2016gkg compared to other SNe representative of the SN Ib/IIb and II classes, at similar phases: SN 2009jf (Ib; Valenti et al. 2011; Shivvers et al. 2019), SN 1993J (Iib; Matheson et al. 2000), SN 2011dh (Iib; Ergon et al. 2015), SN 2004et (II; Sahu et al. 2006), SN 2004dj (II; Silverman et al. 2017), and SN 1999em (II; Faran et al. 2014). The most prominent lines are indicated. Spectra are scaled to match the broad [O I] components and shifted for clarity. (Right panel) [O I] profile comparison (zoom-in of left panel), with the broad component of SN 2016gkg shown by dashed lines. SN 2016gkg appears more similar to SNe II in terms of the broad component profile and width.</figcaption></figure>

Note that as mentioned earlier in Section 2, the SN rest frame is assumed to coincide with the underlying H II region. If the rest wavelength is instead the central wavelength of the broad [O I]λ6300 emission line, then the interpretation would become more general: the broad and narrow components of [O I] lines are mostly symmetric with a minor red clump, suggesting a dense core in spherical ejecta (although an inner disk-like structure viewed face-on may in principle also be accommo- dated in this scenario).

Comparing the two GMOS spectra obtained with the identical instrument and resolution (+320 d and +764 d), the broad component is found to be diminishing in time with respect to the narrow component (Figure 9). One explanation might be that γ-rays increasingly leak from the ejecta, and thus the inner denser regions trap more efficiently relative to the faster-moving material above. The flux in the broad component then decreases faster than that of the narrow component.

Assuming that the broad component is emitted from a region defined by 0.2 $M_\odot$ and 2000 km s$^{-1}$, then the γ-ray optical depth is described by the following: $\tau \approx 2 \times 10^4 \times (M/0.5 M_\odot) \times (v/2000 \text{ km s}^{-1})^{-2} \times (\text{day})^{-2}$. This component starts to become transparent to γ-rays (below $\tau \approx 1$) at around 90 days. After that, the luminosity would follow the decrease in the deposition rate, i.e., $\tau \propto t^{-2}$. On the other hand, if the narrow component is defined by 0.1 $M_\odot$ and 300 km s$^{-1}$, then $\tau$ is larger by a factor of 20 at around 90 days. For this component, the transition time to γ-ray transparency would be around 420 days. As a consequence, the

<figure><img src="image_11.png" /><figcaption>**Figure 9.** Comparison of the two GMOS spectra with identical spectral resolution, at phases +320 days and +764 days at the [O I] doublet region. The spectra are normalized in flux to match the peak intensity of the λ6300 line.</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>**Figure 10.** Evolution of $[O\,I]/[Ca\,II]$ flux ratio in SN 2016gkg, including both the broad and narrow components, integrated over the velocity range of ±2500 km s⁻¹.</figcaption></figure>

luminosity ratio of the broad/narrow component would decrease over time, and after the narrow component transitions to be γ-ray thin, this ratio would saturate. This estimate, while rough, supports qualitatively that the narrow component is relatively denser and thus its transition to the optically thin phase is delayed. Such a phenomenon was predicted by the two-component SN explosion model of Maeda et al. (2003). Alternatively, as density increases inwards in all core-collapse SN explosions, the trapping of γ-rays is more efficient at low velocity than high velocity even in 1D spherical models.

The flux ratio of $[O\,I]/[Ca\,II]$ has been frequently used as an indicator of core mass (see applications in, e.g., Kuncarayakti et al. 2015; initially in Fransson & Chevalier 1989). Modeling of SESN nebular spectra suggests that this ratio increases slowly with time (Jerkstrand 2017). In SN 2016gkg, the line ratio is found to be increasing (Figure 10). The ratio increase might indicate that some part of the oxygen-rich region is not becoming transparent to γ-rays after around two years. At late time, a high-density oxygen-rich core selectively absorbs γ-ray energy input. If the core is described by low-density and high-density regions (i.e., not a single component; Maeda et al. 2003), the latter would stay opaque to gamma-rays and could emit for longer time. This behavior is seen in the diminishing broad component relative to the narrow component of SN 2016gkg (Figure 9).

The two-component ejecta model has been proposed in the past to describe the late-time light curves of hypernovae (Maeda et al. 2003), as a single-component model is insufficient to explain the observed late-time light curve and spectra of these objects (Mazzali et al. 2000, 2001; Sollerman et al. 2000; Izzo et al. 2019). This model describes the ejecta as consisting of two zones: the outer component that is dominating the early-phase light curve and the broad spectral features, and the dense inner component that is dominating in the later phase. In the case of hypernovae, the presence of the dense core material has been predicted by jet-driven explosion models (e.g., Maeda et al. 2002). Such a two-component shell model of the ejecta is, however, impossible to produce in 1D hydrodynamic simulations, and the only way to produce this is to break spherical symmetry (Dessart et al. 2017).

It has been established that most, if not all, stripped-envelope SN are aspherical explosions (Maeda et al. 2008; Taubenberger et al. 2009). The degree of asphericity generally range from extreme (explosions with relativistic jets, such as those of gamma-ray bursts and broad-lined type-Ic SNe) to mildly deviating from spherical geometry. In between, SN explosions may exhibit bipolar outflows (Maeda & Nomoto 2003). In a highly asymmetric explosion, the line-of-sight orientation may affect the observed nebular spectra. The broad-lined type-Ic SNe 1998bw and 2003jd are both inferred to be energetic explosions that launched relativistic jets. The jet was pointed toward Earth in the case of SN 1998bw, resulting in the gamma-ray burst (GRB) 980425, and also evident in the nebular spectra where the $[O\,I]$ line is sharp-peaked, and $[Fe\,II]$ is found at a higher velocity compared to $[O\,I]$ (Maeda et al. 2002).

In a spherical explosion, oxygen is expected to be located in the outer layers relative to iron, and thus to show higher velocities. Whereas, in an aspherical event, explosive nucleosynthesis near the jet axis produces iron that moves with higher velocity compared to oxygen, which is distributed closer to the equatorial region, perpendicular to the jet axis. Therefore, in the case of SN 2003jd, the line of sight goes through the jet equator as evidenced by the double-peaked nebular $[O\,I]$ line, and since the jet is pointing away from Earth, no GRB was observed (Mazzali et al. 2005). In SN 2016gkg, the double-peaked $[O\,I]$ line indicates that the line of sight is probing the oxygen-rich region that produces the strong $[O\,I]$ narrow emission. However, there is no low-velocity analog to a jet such as those seen in hypernovae. The asymmetry in SN 2016gkg has to arise from something completely different.

Evidence of asymmetry or bipolarity in core-collapse SNe, in addition to nebular spectroscopy, include a mildly relativistic radio jet in the energetically normal type-Ic SN 2007gr (Paragi et al. 2010). SN 2007gr is also found to match the two-component explosion model, along with a number of type-Ib SNe (Cano et al. 2014). In addition, observations of SN remnants also support this, such as the Si-rich bilateral jet of ejecta in the Vela SN remnant (García et al. 2017) and the distribution of $^{44}$Ti in the SN remnant Cassiopeia A (Grefenstette et al. 2014). SN light echo observations that probe an SN from multiple points of view, also add more evidence for explosion asymmetry (Rest et al. 2011). In the case of hydrogen-rich (type-II) SNe, signs of bipolarity and asphericity have been observed (e.g., SN 2010jp, Smith et al. 2012; SN 2016X, Bose et al. 2019; SN 2017gmr, Andrews et al. 2019; Nagao et al. 2019), and it is likely that most of them are asymmetric explosions (Leonard et al. 2006, 2016). Even in low-energy SN II explosions such as in SN 2008bk, asymmetry is clearly detected (Leonard et al. 2012), and this further suggests that asymmetry does not require high energy explosions. Since SESNe lack the hydrogen envelope, they should be even more asymmetric compared to the type-II SNe. Three-dimensional simulations of neutrino-driven core-collapse SN explosions also suggest that asymmetry is widespread (e.g., Kifonidis et al. 2000, 2003; Wongwathanarat et al. 2015).

The velocity contrast between the narrow and broad components in SN 2016gkg is very large, suggesting that two different regions are responsible for the emitted profiles. If the energy sources are the same (i.e., $^{56}$Ni), since the narrow component region encompass a much smaller volume (20% of the broad component's velocity, therefore 0.8% of the volume) then the factor of 40% in flux is an order of magnitude larger compared to
the case where the oxygen density is homogeneous. SN 2016gkg has been modeled as a standard-energy explosion with kinetic energy of $1.2 \times 10^{51}$ erg and ejecta mass of $3.4 M_{\odot}$ (Bersten et al. 2018)—in SESNe, where typical explosions produce $\sim 0.25–3.0 \times 10^{51}$ erg of kinetic energy and $1–6 M_{\odot}$ of ejecta (Lyman et al. 2016; Taddia et al. 2018), the ejecta material would be expelled with velocity in the order of thousands of km s$^{-1}$, and thus one would not expect material at low velocities. This material with low velocity is not predicted in 1D explosion models and immediately implies an asymmetric explosion.

SN 2016gkg thus presents compelling evidence that significant asymmetry may be present even in a standard-energy explosion. The asymmetry is most striking at low velocity, possibly the lowest velocity material ever identified in an SESN, and not at large velocity. Despite the presence of asymmetry, the SN explosion parameters, such as the explosion energy, ejecta mass, and $^{56}$Ni mass, inferred from modeling of the early phase (around light-curve peak) with 1D explosion models would still hold as low-velocity material mostly affects the late-time light curves (Maeda et al. 2003, 2006) and spectra, as shown here. In 1D explosion models, the energy needed to explain the early light curve and spectra would yield no material at low velocity (see, e.g., Dessart et al. 2017). If multidimensional models are invoked, there could be mass at low velocity that causes the ejecta mass estimate from hydrodynamic simulations of the light curves to be slightly underestimated. This material would still have little impact to the derived explosion energy, due to its low velocity, and to the total ejecta mass, due to its small contribution ($\sim 0.1 M_{\odot}$ in the case of SN 2016gkg). In asymmetric neutrino-driven explosions, one typically observes that the explosion occurs in some directions, while continued accretion (and fallback) occurs in other directions. Low-velocity material would correspond to the material that just barely avoids fallback and stays at low velocity in the final 3D ejecta (e.g., Wongwathanarat et al. 2013; Chan et al. 2020). Such material resides in the deeper layers of the exploding star and thus would not significantly affect the early light curves and spectra.

### 3.5. Alternative Interpretations of the Narrow Emission Lines

The narrow lines of SN 2016gkg cannot be explained by interaction with circumstellar material (CSM). To produce a spectrum with strong, narrow [O I] lines, the CSM needs to be oxygen-rich, and it has to be depleted in H and He as the narrow lines of these elements are not present. SN 2016gkg clearly showed signatures of H and He in the early phase (Tartaglia et al. 2017; Bersten et al. 2018), and thus an oxygen-rich CSM poor in H/He would require a fine tuning and is therefore highly unlikely.

Thus far, there is only one instance of an SESN interacting with H/He-free CSM in the H/He-poor type-Ic SN 2010mb (Ben-Ami et al. 2014), while on the other hand a type-Ic SN interacting with H/He-rich CSM displays strong H and He emission lines (SN 2017dio, Kuncarayakti et al. 2018). The emergence of a blue continuum in SN 2010mb is also attributed to the interaction with the H/He-free CSM, and this feature is not apparent in SN 2016gkg. In addition, the type-Ib SNe characteristically display strong He emission lines resulting from interaction with He-rich CSM (e.g., Pastorello et al. 2007). None of these examples show the characteristics seen in the spectra of SN 2016gkg.

Weak interaction with H-rich CSM in SN 2016gkg is not ruled out, nevertheless, as the bump seen in the spectra around

<figure><img src="image_13.png" /><figcaption>**Figure 11.** Zoom-in of the $H\alpha$ region in the spectra of SN 2016gkg and other well-studied type-Ib SNe with CSM interaction.</figcaption></figure>

the $H\alpha$ region may be contributed by such interaction in addition to the emission from [N II]. Such interaction with H-rich CSM has been seen in a number of SNe Ib, such as SNe 1993J and 2013df (Matheson et al. 2000; Maeda et al. 2015), where shock-induced $H\alpha$ emission may reach velocities around $10,000$ km s$^{-1}$ (see Figure 11). Still, this cannot explain the observed strong narrow lines in SN 2016gkg.

Another interesting possibility is the interaction with a stripped stellar companion. However, this scenario cannot explain the symmetric blue and red narrow-line peaks. The narrow lines of SN 2016gkg are also not likely to be the result from a pulsar wind nebula in the inner ejecta, as they are of low ionization. As exemplified by the type-Ib SN 2012au (Milisavljevic et al. 2018), the [O III] $\lambda\lambda4959,5007$ and [O II] $\lambda\lambda7320,7330$ appear strong, comparable to [O I], as a result of a pulsar wind nebula inside the SN ejecta. We note the absence of strong [O II] and [O III] emission lines, which could also arise from CSM interaction or from magnetar power (Jerkstrand et al. 2017; Dessart 2019), in the spectra of SN 2016gkg. [O II] and [O III] are more easily produced in low-density conditions, cf. [O I] that is associated with high density. There are weak features in the positions of [O II] and [O III] in the +764 day spectrum, therefore similarly a pulsar wind nebula may be present at the core, but its power contribution is limited even at the latest epoch and therefore it cannot be responsible for powering the set of conspicuous emission lines seen in the earlier spectra.

### 3.6. SN Environment

NGC 613, the host galaxy of SN 2016gkg, is a large spiral galaxy with a measured metallicity $8.6 \pm 6$ north from the SN explosion site of $12+\log(O/H) = 8.7$ dex, close to the solar value (Bersten et al. 2018, using a combination of indices). The galaxy exhibits gas outflows near the center (López-Cobá et al. 2020); however, this does not seem to affect the outer parts where SN 2016gkg lies.
We extracted the background stellar population of SN 2016gkg from our MUSE datacube, integrating light within an annulus extending between 0.''9–1.''7 (~100–200 pc) from the SN. We do not detect any anomalous flows around the SN in the Hα velocity field extracted from the MUSE datacube. The redshift of the background narrow Hα line is consistent with the surrounding H II regions in the field and the host redshift. Utilizing the [N II] λ6584 and Hα emission lines, an oxygen abundance of 12+log (O/H) = 8.6 dex in the N2 scale (Marino et al. 2013) was calculated from the extracted stellar population spectrum, consistent with the near-solar metallicity estimate. The equivalent width of the Hα emission line in this background stellar population was measured to be around 7 Å. While other factors such as shocks and emissions from post-asymptotic giant branch stars may contribute to the observed emission line equivalent width (see, e.g., Cid Fernandes et al. 2011), here we assume that the Hα emission is contributed by the star-forming parent stellar population of the progenitor of SN 2016gkg and thus can be used as an age indicator. Following Kuncarayakti et al. (2018), this Hα equivalent width value is consistent with the stellar population age of around 12 Myr by comparing with Starburst99 simple stellar population models (Leitherer et al. 1999), which subsequently points to the lifetime of a single 16 M⊙ star. This estimate is consistent with the progenitor initial mass range of 15–20 M⊙ derived with other methods discussed above.

## 4. Summary

We present late-time spectral observations of the type-Ib SN 2016gkg. The nebular spectra show the characteristic [O I] λλ6300,6364 and [Ca II] λλ7292,7324 emission doublets as the strongest lines, accompanied by other lines, which are fainter.

The most striking features of the spectra are the strong narrow lines superposed on the broad component base, reaching low velocities down to ~300 km s⁻¹. The [O I] doublet is the strongest line in the spectra, and both the λλ6300 and 6364 lines feature a composite line profile consisting of at least three single components: a broad emission (FWHM ~ 2000 km s⁻¹), and blue and red emission components mirrored against the zero velocity. Such a low velocity seen in the narrow lines is not expected from spherical 1D SN explosion models and points to the presence of slow-moving material in the ejecta.

SN 2016gkg suggests that such low-velocity material might also be present in other SNe, although not observed thus far. A more detailed modeling of the power source at the origin of the narrow-line component is needed. While we made no attempt to model the explosion and spectra of SN 2016gkg, future SN modeling efforts may take advantage of these observations and incorporate more detailed ejecta structures to improve the constructed models of SNe and their progenitors.

We thank the anonymous referee and Stefan Taubenberger for useful comments and discussions. Archival data used in this article were obtained via the UC Berkeley Supernova Database<fnref n="22" /> (Shivvers et al. 2019), the Open Supernova Catalog<fnref n="23" /> (Guillochon et al. 2017), and the Weizmann Interactive Supernova Data Repository<fnref n="24" /> (WISEREП, Yaron & Gal-Yam 2012). H.K. was funded by the Academy of Finland projects

324504 and 328898. K.M. acknowledges support by JSPS KAKENHI Grant (20H00174, 20H04737, 18H04585, 18H05223, 17H02864). L.G. was funded by the European Union's Horizon 2020 research and innovation program under the Marie Skłodowska-Curie grant agreement No. 839090. F.G. acknowledges support from PIP 0102 (CONICET) and PICT-2017-2865 (ANPCyT). This work has been partially supported by the Spanish grant PGC2018-095317-B-C21 within the European Funds for Regional Development (FEDER). C.P.G. acknowledges support from EU/FP7-ERC grant No. [615929]. F.O.E. acknowledges support from the FONDECYT grant No. 11170953 and 1201223.

**Facilities:** Gemini-S (GMOS-S), Subaru (FOCAS), VLT (MUSE).

**Software:** Astropy (Astropy Collaboration et al. 2018), ESO Reflex (Freudling et al. 2013), GDL (Coulais et al. 2010), IRAF (Tody 1986, 1993), L.A.Cosmic (van Dokkum 2001), mpfit (Markwardt 2009), QFitsView (Ott 2012), ZAP (Soto et al. 2016).

## ORCID iDs

Haindyo Kuncarayakti https://orcid.org/0000-0002-1132-1366

Gastón Folatelli https://orcid.org/0000-0001-5247-1486

Keiichi Maeda https://orcid.org/0000-0003-2611-7269

Luc Dessart https://orcid.org/0000-0003-0599-8407

Anders Jerkstrand https://orcid.org/0000-0001-8005-4030

Joseph P. Anderson https://orcid.org/0000-0003-0227-3451

Melina C. Bersten https://orcid.org/0000-0002-6991-0550

Lluís Galbany https://orcid.org/0000-0002-1296-6887

Federico García https://orcid.org/0000-0001-9072-4069

Claudia P. Gutiérrez https://orcid.org/0000-0003-2375-2064

Koji S. Kawabata https://orcid.org/0000-0001-6099-9539

Joseph D. Lyman https://orcid.org/0000-0002-3464-0642

Seppo Mattila https://orcid.org/0000-0001-7497-2994

Felipe Olivares E. https://orcid.org/0000-0002-5115-6377

Sebastián F. Sánchez https://orcid.org/0000-0001-6444-9307

Schuyler D. Van Dyk https://orcid.org/0000-0001-9038-9950

## References

* Andrews, J. E., Sand, D. J., Valenti, S., et al. 2019, ApJ, 885, 43

* Arcavi, I., Hosseinzadeh, G., Brown, P. J., et al. 2017, ApJL, 837, L2

* Astropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M., et al. 2018, AJ, 156, 123

* Bacon, R., Vernet, J., Borisova, E., et al. 2014, Msngr, 157, 13

* Ben-Ami, S., Gal-Yam, A., Mazzali, P. A., et al. 2014, ApJ, 785, 37

* Bersten, M. C., Folatelli, G., García, F., et al. 2018, Natur, 554, 497

* Bose, S., Dong, S., Elias-Rosa, N., et al. 2019, ApJL, 873, L3

* Cano, Z., Maeda, K., & Schulze, S. 2014, MNRAS, 438, 2924

* Chan, C., Müller, B., & Heger, A. 2020, MNRAS, 495, 3751

* Cid Fernandes, R., Stasińska, G., Mateus, A., et al. 2011, MNRAS, 413, 1687

* Coulais, A., Schellens, M., Gales, J., et al. 2010, in ASP Conf. Ser. 434, adass XIX, ed. Y. Mizumoto, K.-I. Morita, & M. Ohishi (San Francisc, CA: ASP), 187

* Dessart, L. 2019, A&A, 621, A141

* Dessart, L., & Hillier, D. J. 2005, A&A, 437, 667

* Dessart, L., John Hillier, D., Yoon, S.-C., et al. 2017, A&A, 603, A51

* Elmhamdi, A. 2011, AcA, 61, 179

* Ergon, M., Jerkstrand, A., Sollerman, J., et al. 2015, A&A, 580, A142

* Fang, Q., Maeda, K., Kuncarayakti, H., et al. 2019, NatAs, 3, 434

* Faran, T., Poznanski, D., Filippenko, A. V., et al. 2014, MNRAS, 442, 844

* Filippenko, A. V. 1988, AJ, 96, 1941

<footnote n="22">22 http://heracles.astro.berkeley.edu/sndb/</footnote>

<footnote n="23">23 https://sne.space/</footnote>

<footnote n="24">24 https://wiserep.weizmann.ac.il/</footnote>
Fransson, C., & Chevalier, R. A. 1989, ApJ, 343, 323

Freudling, W., Romaniello, M., Bramich, D. M., et al. 2013, A&A, 559, A96

Galbany, L., Anderson, J. P., Rosales-Ortega, F. F., et al. 2016, MNRAS, 455, 4087

García, F., Suárez, A. E., Miceli, M., et al. 2017, A&A, 604, L5

Gimeno, G., Roth, K., Chiboucas, K., et al. 2016, Proc. SPIE, 9908, 99082S

Grefenstette, B. W., Harrison, F. A., Boggs, S. E., et al. 2014, Natur, 506, 339

Gröningsson, P., Fransson, C., Lundqvist, P., et al. 2008, A&A, 479, 761

Guillochon, J., Parrent, J., Kelley, L. Z., & Margutti, R. 2017, ApJ, 835, 64

Hook, I. M., Jørgensen, I., Allington-Smith, J. R., et al. 2004, PASP, 116, 425

Izzo, L., de Ugarte Postigo, A., Maeda, K., et al. 2019, Natur, 565, 324

Jerkstrand, A., Ergon, M., Smartt, S. J., et al. 2015, A&A, 573, A12

Jerkstrand, A., Smartt, S. J., Fraser, M., et al. 2014, MNRAS, 439, 3694

Jerkstrand, A., Smartt, S. J., Inserra, C., et al. 2017, ApJ, 835, 13

Jerkstrand, A. 2017, in Handbook of Supernovae, ed. A. W. Alberti & P. Murdin (Cham: Springer), 795

Kashikawa, N., Aoki, K., Asai, R., et al. 2002, PASJ, 54, 819

Kifonidis, K., Plewa, T., Janka, H.-T., et al. 2000, ApJL, 531, L123

Kifonidis, K., Plewa, T., Janka, H.-T., et al. 2003, A&A, 408, 621

Kilpatrick, C. D., Foley, R. J., Abramson, L. E., et al. 2017, MNRAS, 465, 4650

Kuncarayakti, H., Anderson, J. P., Galbany, L., et al. 2018, A&A, 613, A35

Kuncarayakti, H., Maeda, K., Ashall, C. J., et al. 2018, ApJL, 854, L14

Kuncarayakti, H., Maeda, K., Bersten, M. C., et al. 2015, A&A, 579, A95

Leibundgut, B., Kirshner, R. P., Pinto, P. A., et al. 1991, ApJ, 372, 531

Leitherer, C., Schaerer, D., Goldader, J. D., et al. 1999, ApJS, 123, 3

Leonard, D. C., Dessart, L., Hillier, D. J., et al. 2012, in AIP Conf. Ser. 1429, Stellar Polarimetry: From Birth to Death (Melville, NY: AIP), 204

Leonard, D. C., Dessart, L., Pignata, G., et al. 2016, IAUFM, 29B, 458

Leonard, D. C., Filippenko, A. V., Ganeshalingam, M., et al. 2006, Natur, 440, 505

Li, H., & McCray, R. 1992, ApJ, 387, 309

Limongi, M., & Chieffi, A. 2003, ApJ, 592, 404

López-Cobá, C., Sánchez, S. F., Anderson, J. P., et al. 2020, AJ, 159, 167

Lyman, J. D., Bersier, D., James, P. A., et al. 2016, MNRAS, 457, 328

Maeda, K., Hattori, T., Milisavljevic, D., et al. 2015, ApJ, 807, 35

Maeda, K., Kawabata, K., Mazzali, P. A., et al. 2008, Sci, 319, 1220

Maeda, K., Mazzali, P. A., Deng, J., et al. 2003, ApJ, 593, 931

Maeda, K., Mazzali, P. A., & Nomoto, K. 2006, ApJ, 645, 1331

Maeda, K., Nakamura, T., Nomoto, K., et al. 2002, ApJ, 565, 405

Maeda, K., & Nomoto, K. 2003, ApJ, 598, 1163

Marino, R. A., Rosales-Ortega, F. F., Sánchez, S. F., et al. 2013, A&A, 559, A114

Markwardt, C. B. 2009, in ASP Conf. Ser. 411, adass XVIII, 411, ed. D. A. Bohlender, D. Durand, & P. Dowler (San Francisco, CA: ASP), 251

Matheson, T., Filippenko, A. V., Ho, L. C., Barth, A. J., & Leonard, D. C. 2000, AJ, 120, 1499

Mazzali, P. A., Iwamoto, K., & Nomoto, K. 2000, ApJ, 545, 407

Mazzali, P. A., Kawabata, K. S., Maeda, K., et al. 2005, Sci, 308, 1284

Mazzali, P. A., Nomoto, K., Patat, F., & Maeda, K. 2001, ApJ, 559, 1047

Meikle, W. P. S., Kotak, R., Farrah, D., et al. 2011, ApJ, 732, 109

Meyer, M. J., Zwaan, M. A., Webster, R. L., et al. 2004, MNRAS, 350, 1195

Milisavljevic, D., Fesen, R. A., Gerardy, C. L., et al. 2010, ApJ, 709, 1343

Milisavljevic, D., Patnaude, D. J., Chevalier, R. A., et al. 2018, ApJL, 864, L36

Modjaz, M., Kirshner, R. P., Blondin, S., Challis, P., & Matheson, T. 2008, ApJL, 687, L9

Nagao, T., Cikota, A., Patat, F., et al. 2019, MNRAS, 489, L69

Nomoto, K., Hashimoto, M., Tsujimoto, T., et al. 1997, NuPhA, 616, 79

Ott, T. 2012, QFitsView: FITS file viewer, Astrophysics Source Code Library, ascl:1210.019

Paragi, Z., Taylor, G. B., Kouveliotou, C., et al. 2010, Natur, 463, 516

Pastorello, A., Smartt, S. J., Mattila, S., et al. 2007, Natur, 447, 829

Piro, A. L., Muhleisen, M., Arcavi, I., et al. 2017, ApJ, 846, 94

Rest, A., Foley, R. J., Sinnott, B., et al. 2011, ApJ, 732, 3

Sahu, D. K., Anupama, G. C., Srividya, S., et al. 2006, MNRAS, 372, 1315

Sahu, D. K., Gurugubelli, U. K., Anupama, G. C., et al. 2011, MNRAS, 413, 2583

Shivvers, I., Filippenko, A. V., Silverman, J. M., et al. 2019, MNRAS, 482, 1545

Silverman, J. M., Pickett, S., Wheeler, J. C., et al. 2017, MNRAS, 467, 369

Smith, N., Cenko, S. B., Butler, N., et al. 2012, MNRAS, 420, 1135

Sollerman, J., Kozma, C., Fransson, C., et al. 2000, ApJL, 537, L127

Soto, K. T., Lilly, S. J., Bacon, R., Richard, J., & Conseil, S. 2016, MNRAS, 458, 3210

Taddia, F., Stritzinger, M. D., Bersten, M., et al. 2018, A&A, 609, A136

Tartaglia, L., Fraser, M., Sand, D. J., et al. 2017, ApJL, 836, L12

Taubenberger, S., Valenti, S., Benetti, S., et al. 2009, MNRAS, 397, 677

Tody, D. 1986, Proc. SPIE, 627, 733

Tody, D. 1993, in ASP Conf. Ser. 52, adass II, 52, ed. R.J. Hanisch, R. J. V. Brissenden, & J. Barnes (San Francisco, CA: ASP), 173

Valenti, S., Fraser, M., Benetti, S., et al. 2011, MNRAS, 416, 3138

van Dokkum, P. G. 2001, PASP, 113, 1420

Waxman, E., & Katz, B. 2017, in Handbook of Supernovae, ed. A. W. Alberti & P. Murdin (Cham: Springer), 967

Wongwathanarat, A., Janka, H.-T., & Müller, E. 2013, A&A, 552, A126

Wongwathanarat, A., Müller, E., & Janka, H.-T. 2015, A&A, 577, A48

Woolesley, S. E., Eastman, R. G., Weaver, T. A., & Pinto, P. A. 1994, ApJ, 429, 300

Yaron, O., & Gal-Yam, A. 2012, PASP, 124, 668