A large fraction of hydrogen-rich supernova progenitors experience elevated mass loss shortly prior to explosion

The Zwicky Transient Facility (ZTF) is a wide-field, high cadence, multiband survey that started operating in March 2018 (Bellm et al., 2019; Graham et al., 2019). ZTF imaging is obtained using the Samuel Oschin 48” Schmidt telescope at Palomar observatory (P48). ZTF observing time is divided among three programs: the public (MSIP) 3-day all-sky survey, partnership surveys, and Caltech programs. This paper is based on data obtained by the high-cadence partnership survey. As part of this program, during 2018, extra-galactic survey fields were observed in both the ZTF g- and r-bands 2-3 times per night per band. New images were processed through the ZTF pipeline (Masci et al., 2019) and reference images built by combining stacks of previous ZTF imaging in each band were then subtracted using the Zackay et al. (2016) image subtraction algorithm (ZOGY). A 30s integration time was used in both $g-$ and $r-$band exposures, and a 5$\sigma$ detection limit down to $\sim$20.5 mag in $r$ can be reached in a single observation.

We conducted our year-1 ZTF survey for infant SNe following the methodology of Gal-Yam et al. (2011). We selected potential targets via a custom filter running on the ZTF alert stream using the GROWTH Marshal platform (Kasliwal et al., 2019). The filter scheme was based on the criteria listed in Table 1.

Alerts that passed our filter (typically $50-100$ alerts per day) were then visually scanned by a duty astronomer, in order to reject various artefacts (such as unmasked bad pixels or ghosts) and false positive signals, such as flaring M stars, CVs and AGN. Most spurious sources could be identified by cross matching with additional catalogues (e.g., WISE IR photometry (Wright et al., 2010) to detect red M stars, the Gaia DR2 catalog (Gaia Collaboration et al., 2018), and catalogs from time domain surveys such as the Palomar Transient Factory (PTF; Law et al. 2009) and the Catalina Real-Time Survey (CRTS; Drake et al. 2014) for previous variability of CVs and AGN.

Due to time-zone differences, our scanning team (located mostly at the Weizmann Institute in Israel and at the Oskar Klein Center (OKC) in Sweden) could routinely scan the incoming alert stream during the California night time, with the goal of triggering spectroscopic follow-up of promising infant SN candidates within hours of discovery (and thus typically within $<2$ days from explosion), as well as Swift (Gehrels et al., 2004) Target-of-Opportunity (ToO) UV photometry.

Figure  1 shows the SN Type distribution amongst the $\sim 2500$ spectroscopically-confirmed SNe gathered by ZTF between March and December 2018. About 37% are core-collapse events, and $\sim 63\%$ of those are of Type II. Since the large majority of flash events are SNe II, we can only place statistically meaningful constraints on the frequency of this phenomenon among Type II SNe. We therefore analyze here this population only.

Our infant SN program allowed us to obtain early photometric and spectroscopic follow-up of young SNe. However, it is possible that we have missed some relevant candidates. In order to ensure the completeness of our sample, we therefore inspected all spectroscopically classified SNe II (including subtypes IIn and IIb) from ZTF ¹¹1between March 2018 and December 2018 using the ZTFquery package (Rigault, 2018). We pulled from this sample all events (the large majority) lacking a ZTF non-detection limit within 2.5 days prior to the first detection recorded on the ZTF Marshal. To include events in our final sample, we required that they show significant and rapid increase in flux, as previously observed for very young SNe (e.g., Gal-Yam et al. 2014, Yaron et al. 2017), with respect to the last non-detection. This excludes older events that are just slightly below our detection limit and are picked up by the filter when they slowly rise, or when conditions improve. We implemented a cut on the observed rise of $\Delta$r or $\Delta$g$>0.5$ mag with respect to the recent limit in the same band, and note all events that satisfy this cut as “real infant” (RI; Fig. 1, panel C).

All in all, we gathered 43 candidates which fulfilled the RI criteria. Additional inspection led us to determine that 15 candidates were spurious (see Appendix A for details). Our final sample (Table 2) thus includes a total of 28 RI Type II SNe, or about 5% of all the SNe II found by the ZTF survey for 2018. During its first year of operation (starting March 2018), ZTF obtained useful observations for our program during approximately 32 weeks, excluding periods of reference image building (initially), periods dedicated to Galactic observations, and periods of technical/weather closure. We therefore find that the survey provided about one real infant SN II per week.

Our goal was to obtain rapid spectroscopy of RI SN candidates following the methods of Gal-Yam et al. (2011). This was made possible using rapid ToO follow-up programs as well as on-request access to scheduled nights on various telescopes. During the scanning campaign, we applied the following criteria for rapid spectroscopic triggers. The robotic SEDm (see below) was triggered for all candidates brighter than a threshold limiting magnitude ($19$ mag during 2018). Higher resolution spectra (using WHT, Gemini or other available instruments) were triggered for events showing recent non-detection limits (within $2.5$ d prior to first detection) as well as a significant rise in magnitude compared to a recent limit or within the observing night.

The ZTF alert system (Patterson et al., 2018) provides on the fly photometry (Masci et al., 2019) and astrometry based on a single image for each alert. In order to improve our photometric measurements (and in particular, to test the validity of non-detections just prior to discovery) we performed forced PSF photometry at the location of each event. As shown by Yaron et al. (2019), the $95\%$ astrometric scatter among ZTF alerts is $\sim 0.44$”; for our events we have multiple detections, with typically higher signal-to-noise ratio data around the SN peak compared to the initial first detections. We therefore compute the median coordinates of all the alert packages and perform forced photometry using this improved astrometric location.

We use the pipeline developed by F. Masci and R. Laher⁷⁷7http://web.ipac.caltech.edu/staff/fmasci/ztf/forcedphot.pdf to perform forced PSF photometry at the median SN centroid on the ZTF difference images available from the IRSA database . For each light curve, we filter out measurements returned by the pipeline with non-valid flux values.

We perform a further quality cut on each light curve by rejecting observations with a data quality parameter scisigpix⁸⁸8A parameter calculated by the pipeline that measures the pixel noise in each science image that is more than 5 times the median absolute deviation (MAD) away from the median of this parameter for each light curve. We also remove faulty measurements where the infobitssci parameter is not zero. According to the Masci & Laher prescription we rescale the flux errors by the square root of the $\chi^{2}$ of the PSF fit estimate in each image. We then correct each measured forced photometry flux value by the photometric zero point of each image, as provided by the pipeline:

We determine our zero-flux baseline using forced photometry observations obtained prior to the SN explosion. We calculate the median of these observations, reject outliers that are $>3$ MAD away from the median, re-calculate the median and subtract it from our measured post-explosion flux values; these corrections are typically very small, of the order of $<0.1\%$ of the supernova flux values.

If the ratio between the measured flux and the uncertainty $\sigma$ is below $3$, we consider this measurement as a non-detection, and report a 5$\sigma$ upper limit. Otherwise (if the flux to error ratio is above $3\sigma$) we report the flux, magnitude and respective errors. Finally, we correct for Galactic extinction using the python package extinction⁹⁹9https://github.com/kbarbary/extinction, using local extinction values from Schlafly & Finkbeiner (2011) and assuming a selective extinction of $R_{V}=3.1$ and the Cardelli et al. (1989) extinction law. We also correct the light curves to restframe time according to the spectroscopic redshifts.

We recovered detections prior to the first detection by the real-time pipeline using the forced photometry pipeline in 11 cases ¹⁰¹⁰10ZTF18aarqxbw, ZTF18aavpady, ZTF18aawyjjq, ZTF18abcezmh, ZTF18abckutn, ZTF18abcptmt, ZTF18abdbysy, ZTF18abddjpt, ZTF18abokyfk, ZTF18abrlljc, ZTF18abvvmdf. We redefined the first detection and last non-detection according to the forced photometry pipeline measurements in these cases.

We present our photometry for all RI objects in Table 3.

In order to estimate the explosion time, which we define here as the time of zero-flux, we fit a general power law of the form to our flux measurements:

using the routine curvefit within the astropy python package (Astropy Collaboration et al. 2013). We fit the first 2 days of data following the first detection as well as the first 5 days (see Fig. 2 for example) in both the $g$ and $r-$bands. The estimated explosion time is taken as the weighted mean of the four fits, and we adopt the standard deviation as the error on this value. In ten cases, however, there were not enough data in either band to perform the fit. In those cases, we set the explosion date as the mean between the time of the last non detection and the first detection. In all but four of the cases the estimated explosion date is within less than a day from the first detection (Fig. 4; Table 2).

Following Khazov et al. (2016), we also test if events showing flash features are on average more luminous. As can be seen from Table 2, the relevant events to consider are only those with relatively early spectra. We therefore compute the peak magnitude of all seventeen events with a first spectrum obtained within $7$ days from explosion. We use the forced photometry lightcurves to evaluate the peak magnitude. We fit a polynomial of order 3 to the flux measurements, over several intervals of time whose lower bound is within the first few days from explosion time and upper bound between 10 to 40 days after the estimated explosion time (Fig. 3). We adopt the mean and range of peak times obtained from these fits as the peak date and its error (vertical grey band in Fig. 3) and take the mean and standard deviation of the flux value within this range to be the peak flux and error (horizontal grey band in Fig. 3). We report these values for each event in each band in Table 4.

We sort the 28 RI SNe in our sample according to the difference between the estimated explosion time and the time of first spectrum (Table 2, “First spectrum” column; Fig. 4). From previous work (Gal-Yam et al. 2014, Yaron et al. 2017, Khazov et al. 2016), we know that flash features are typically present from the time of explosion up to several days later. We therefore define a sub-sample including events with spectra obtained within 2 d from explosion (top of Table 2). For about one third of the total sample (ten objects) we have been able to secure a first spectrum within less than 2 days from the estimated explosion time.

Throughout the 2018 campaign, we find that seven infant supernovae of Type II show flash features (Table 2; Fig 5). Two additional infant objects were marked as potential flash events (Fig. 8; see below). Four of the seven confirmed flashers had their first spectrum obtained with SEDm.

The two-day sub-sample we are considering includes 6 events showing flash features (one object, SN 2018leh, shows flash features but its first spectrum was obtained only $>3$ days after explosion, Table 2), the two potential flashers, and two events have high signal to noise early spectra that show no flash features (Fig. 7).

Based on our systematic survey of infant SNe II with spectra obtained within two days of discovery, we have found that at least $60\%$, and perhaps as many as $80\%$ of the sample of ten events show evidence for flash-ionized emission. Taking into account our limited sample size and assuming binomial statistics ($\mathcal{B}(k,n,p)$), we infer the true fraction of SNe with CSM that manifests as flash features from the true probability p to observe a flash event given we the observed fraction D using a Bayesian model:

Where $p$ is the probability of observing a flash ionised event (here $p\in[0;1]\,$), D is the observation presented in this paper (i.e.: 6 out of 10 candidates are showing flash features). The probability of our observation, $P(D)$ can be calculated with the formula of total probability, i.e. $P(D)=\int_{0}^{1}\mathcal{B}(6,10,p)\times\pi(p)\,\mathrm{d}p$. We assume a uniform distribution for the prior $\pi(p)$ which allows us to write the posterior function as:

This results in a Beta distribution (see Figure 10). We can put a strict lower limit on the fraction of infant SNe II showing flash features of $>30.8\%$ ($>23.5\%$) at the $95\%$ ($99\%$) confidence level. This fraction rapidly drops when events with spectra obtained within $7$ d from explosion are considered; presumably the fraction could be even higher for events with even earlier spectra.

These results are broadly consistent with previous work by Khazov et al. (2016), which estimate that $7-36\%$ show flash features in spectra obtained within $<2$ d from explosion ($68\%$ confidence level). It is also consistent with the low observed frequency of flash features among the general population of Type II SNe reported in the literature, as these events very rarely have a spectrum obtained $<2$ d after explosion, and Table 2 shows that the fraction of flash events falls rapidly at ages $>2$ d. The unique nightly cadence of the ZTF partnership survey enabled us to routinely discover infant SNe and rapidly obtain spectra, while the systematic design of our survey allowed for a robust measurement of the frequency of this phenomenon.

Khazov et al. (2016) (see their Fig. 8) show that Type II SNe showing flash-ionized features tend to be brighter at peak than other events. We cannot confirm this is also true for our sample. We consider here the subsample of Infant Supernovae whose first spectrum was obtained within less than 7 days from the estimated explosion time, the peak magnitudes were obtained following the method described in § 3.3.2. Figure  11, top panel, shows the peak magnitudes in both g and r bands for flashers and non flashers. Flashers appear to be brighter in both bands. However, when one considers SN 2018cyg as a flasher, the average peak magnitude of both groups is inverted and non-flashers appear brighter than flashers (see Table 6, top section). Since SN 2018cyg is strongly reddened, we repeated this same analysis but with SN 2018cyg being host extinction corrected. To apply the extinction correction, we consider the spectrum from 2018 August, 4 ¹²¹²12see on WISeREP : https://wiserep.weizmann.ac.il/object/698 and apply the method described in Poznanski et al. 2012, using the line doublet of sodium. We consider the doublet not to be resolved and apply the following formula:

We estimate the EW of D₁+D₂ using the built-in tool from WISeREP by measuring it several times. The mean EW is 1.64 Å with an error of 0.17 Å. Following Eq. (5), the final peak magnitudes for SN 2018cyg are : $M_{peak,r}=-18.45\pm 0.50$ and $M_{peak,g}=-18.77\pm 0.80$. Table 6 summarises the different cases : whether SN 2018cyg is a flasher and whether SN 2018cyg was corrected for estimated host extinction. We find that flash events are not inherently brighter than non-flash events.

We also inspect in Fig. 11 (lower panel) the distribution of apparent magnitudes at discovery for our $<7$ d sample. As can be seen there, we find that the flash events were not significantly brighter at discovery than other events, and thus neither more likely to be discovered, nor to be followed, as both of these depend on the apparent magnitude of the object at discovery.

We have shown here that a significant fraction, and possibly most, Type II SN progenitors, show transient emission lines in their early spectra, that provide evidence that these stars are embedded in a compact distribution of CSM (Yaron et al., 2017). The narrow width of these emission lines indicates a slow expansion speed for the CSM ($100-800$ km s^-1, Boian & Groh 2020 ), and combined with its compact radial dimension ($<10^{15}$ cm) we have evidence that the CSM was deposited by the stars within months to a few years prior to its terminal explosion. Assuming these progenitors are mostly red supergiants (RSGs; Smartt 2015), this would suggest that most exploding RSGs experience an enhanced mass loss shortly prior to explosion.

While RSGs certainly lose mass during their final stages of evolution (Smith, 2014), such a period of enhanced mass loss shortly (months to a year) prior to explosion is not explained by standard stellar evolution models. Our work thus may indicate that additional physical processes leading to such pre-explosion instabilities (e.g., Arnett & Meakin 2011, Shiode & Quataert 2014) not only exist, but are ubiquitous among massive stars.

As we have shown that most SN II progenitors likely undergo a remarkable evolution shortly prior to explosion, it may be needed to re-examine the stellar models used as initial conditions to explosion simulations. In particular, at least some of the effects proposed to explain such pre-explosion mass loss, may render the spherical pre-explosion stellar models used in explosion simulations less realistic (Arnett & Meakin, 2016). Perhaps our work thus provides a clue how to tackle some of the problems encountered in trying to reproduce the observed distribution of SN explosions using numerical explosion models.