Early Spectra of the Gravitational Wave Source GW170817: Evolution of a Neutron Star Merger

References

Acknowledgments

We thank John Mulchaey (Carnegie Observatories director), Leopoldo Infante (Las Campanas Observatory director), and the entire Las Campanas staff for their dedication, professionalism, and excitement, which were critical for obtaining the observations used in this study.

B.J.S., M.R.D., K.A.A., and A.P.J. were supported by NASA through Hubble Fellowships awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., for NASA, under contract NAS 5-26555. B.J.S. and E.B. were supported by Carnegie-Princeton Fellowships. M.R.D. was supported by a Carnegie-Dunlap Fellowship and acknowledges support from the Dunlap Institute at the University of Toronto. T.W.-S.H. was supported by a Carnegie Fellowship. Support for J.L.P. was in part provided by FONDECYT through the grant 1151445 and by the Ministry of Economy, Development, and Tourism’s Millennium Science Initiative through grant IC120009, awarded to The Millennium Institute of Astrophysics. DK is supported in part by a Department of Energy (DOE) Early Career award DE-SC0008067, a DOE Office of Nuclear Physics award DE-SC0017616, and a DOE SciDAC award DE-SC0018297, and by the Director, Office of Energy Research, Office of High Energy and Nuclear Physics, Divisions of Nuclear Physics, of the U.S. Department of Energy under Contract No.DE-AC02-05CH11231. V.M.P. acknowledges partial support for this work from grant PHY 14- 30152 from the Physics Frontier Center/JINA Center for the Evolution of the Elements (JINA-CEE), awarded by the US National Science Foundation.

The UCSC group was supported in part by NSF grant AST–1518052, the Gordon & Betty Moore Foundation, the Heising-Simons Foundation, generous donations from many individuals through a UCSC Giving Day grant, and from fellowships from the Alfred P. Sloan Foundation (R.J.F), the David and Lucile Packard Foundation (R.J.F. and E.R.) and the Niels Bohr Professorship from the DNRF (E.R.). T.B. acknowledges support from the CONACyT Research Fellowships program. G.M. acknowledges support from CONICYT, Programa de Astronomía/PCI, FONDO ALMA 2014, Proyecto No 31140024. A.M.B. acknowledges support from a UCMEXUS-CONACYT Doctoral Fellowship. C.A. was supported by Caltech through a Summer Undergraduate Research Fellowship (SURF) with funding from the Associates SURF Endowment.

This paper includes data gathered with the 6.5 meter Magellan Telescopes located at Las Campanas Observatory, Chile. Part of this work is based on a comparison to observations of GRB130603B that were obtained from the ESO Science Archive Facility under request number vmplacco308387.

We thank Antonino Cucchiara for sending us the spectrum of GRB130603B and our thoughts go out to all those in the Virgin Islands impacted by the recent series of hurricanes. We thank Jennifer van Saders for useful suggestions. We thank the University of Copenhagen, DARK Cosmology Centre, and the Niels Bohr International Academy for hosting D.A.C., R.J.F., A.M.B., E.R., and M.R.S. during a portion of this work. R.J.F., A.M.B., and E.R. were participating in the Kavli Summer Program in Astrophysics, “Astrophysics with gravitational wave detections.” This program was supported by the Kavli Foundation, Danish National Research Foundation, the Niels Bohr International Academy, and the Dark Cosmology Centre.

This research has made use of the NASA/IPAC Extragalactic Database (NED), which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration.

The data presented in this work and code used to perform the analysis are available at ftp://ftp.obs.carnegiescience.edu/pub/SSS17a. Calibrated rest-wavelength spectra are also made available at WISeREP (https://wiserep.weizmann.ac.il/).

Supplementary Materials: www.sciencemag.org Materials and Methods Figure S1 Table S1, S2 References ($48$–$62$)

The First Spectra of a Gravitational Wave Source: Dramatic Evolution of a Neutron Star Merger

B. J. Shappee,1,2∗ J. D. Simon,1 M. R. Drout,1 A. L. Piro,1

N. Morrell,3 J. L. Prieto,4,5 D. Kasen,6,7 T. W.-S. Holoien,1 J. A. Kollmeier,1

D. D. Kelson,1 D. A. Coulter,8 R. J. Foley,8 C. D. Kilpatrick,8

M. R. Siebert,8 B. F. Madore,1 A. Murguia-Berthier,8 Y.-C. Pan,8

J. X. Prochaska,8 E. Ramirez-Ruiz,8,9 A. Rest,10,11 C. Adams,12

K. Alatalo,1,10 E. Bañados,1 J. Baughman,4,13 R. A. Bernstein,1 T. Bitsakis,14

K. Boutsia,3 J. R. Bravo,3 F. Di Mille,3 C. R. Higgs,15,16 A. P. Ji,1,11

G. Maravelias,17 J. L. Marshall,18 V. M. Placco,19 G. Prieto,3 Z. Wan20

This PDF file includes: Materials and Methods Figure S1 Table S1, S2 References ($48$–$62$)

S1 Data Acquisition & Reductions

We obtained spectroscopic observations of SSS17a with the Magellan/Clay and Magellan/Baade telescopes beginning 11.75 hours after the neutron star merger, and continued observing SSS17a spectroscopically for 8 days. Below we describe the data acquisition, reduction, and calibration. A log of all spectroscopic observations is given in Table S1.

We observed SSS17a with the Low Dispersion Survey Spectrograph (LDSS-3) on the Magellan/Clay telescope on 2017 Aug. 18-26 (UT). On the first night we obtained spectra in multiple spectrograph configurations, with the Volume Phase Holographic (VPH)-All, VPH-Blue, and VPH-Red grisms. The three grisms cover wavelength ranges of $3800-6200$ Å, $4250-10000$ Å, and $6000-10000$ Å with resolving powers of $R=1400$, $R=650$, and $R=1400$, respectively. The later LDSS-3 spectroscopy, when SSS17a was much fainter and redder, employed only the VPH-All or VPH-Red grisms. All observations were made with a slit width of 1 arcsec.

We reduced and calibrated the LDSS-3 spectra using IRAF following standard procedures, including bias subtraction, flat-fielding, 1-D spectral extraction, and wavelength calibration by comparison to an arc lamp. Flux calibration and telluric correction was performed using a set of custom IDL scripts based on a spectroscopic standard star observed on the same night. The statistical uncertainties on the spectra were calculated by standard error propagation.

For the first LDSS-3 spectra of SSS17a, obtained on 2017 Aug. 17-18, the spectrograph slit was oriented \sim 22\mbox{{}^{\circ}} away from parallactic angle. Given the relatively high airmass (${\rm airmass}\approx 2$) at the time, some flux was lost as a result of differential atmospheric refraction. A similar offset from parallactic angle was used for LDSS-3 observations on 2017 Aug. 18-19, with smaller flux losses because of the lower airmass. To correct for this effect, we first calculated the offset of the source from the center of the slit due to the differential refraction, assuming that the target was in the center of the slit at the effective central wavelength of the $g$ filter that was used for target acquisition. We then measured the atmospheric seeing as a function of wavelength, modeled the source with a two-dimensional Gaussian profile, and computed the fraction of light from the model source that fell in the slit as a function of wavelength. Finally, the calibrated spectrum of SSS17a was scaled to account for the light that landed outside the slit.

The magnitude of the differential refraction and seeing corrections depends on the position of the source within the slit. A source position different than the one we have assumed will result in systematic errors in the scaled spectrum. Because we do not have a way of measuring this position during the spectroscopic exposure, we quantify the resulting systematic uncertainty by calculating the change in the differential refraction and seeing corrections for source positions offset by 0.1 arcsec in either direction (over a 300 s exposure the telescope guiding is expected to be at least this accurate) from the nominal position. At each wavelength, we define the systematic uncertainty to be the average of the changes in the correction factor between the $+0.1$ arcsec and $-0.1$ arcsec offsets.

S1.2 MagE Observations

We observed SSS17a with the Magellan/MagE spectrograph on the nights of 2017 August 17-19. We used a $0.7~{\rm arcsec}\times 10~{\rm arcsec}$ slit to provide a spectral resolving power of $R=5800$ on August 17-18 and a $1.0~{\rm arcsec}\times 10~{\rm arcsec}$ slit to provide a spectral resolving power of $R=4100$ on August 18-19. On the first night we obtained a single 322 s exposure and on the second night we obtained three 1000 s exposures. We reduced the MagE spectra using an IDL pipeline based on the techniques described in , with updates to improve the flux calibration from . Observations on both nights were flux calibrated with a standard star spectrum obtained on 2017 Aug. 18-19.

The spectrum from Aug. 17-18 was taken with the slit oriented at the parallactic angle, so the effects of differential atmospheric refraction should be negligible despite the high airmass (3.0) at the time of observation. The seeing measured from the spatial profile of the source spectrum varied as a function of wavelength, from 0.9 arcsec at the red end of the spectrum to 1.3 arcsec in the blue. We corrected for wavelength-dependent loss of light from this seeing variation by modeling the source with a two-dimensional Gaussian profile, calculating the fraction of the flux that fell in the slit as a function of wavelength, and adjusting the calibrated spectrum accordingly.

For readers who are interested in making use of the Aug. 17-18 MagE spectrum, we urge caution in interpreting the data at wavelengths redder than $\sim 7000$ Å. As a specific example, the apparent step in the spectrum at $\sim 7200$ Å is not a real feature. It occurs at the breakpoint between two spectral orders and at a wavelength where there is significant telluric absorption, which makes matching the continuum levels of the neighboring orders difficult.

The Aug. 18-19 spectra were obtained with the slit 22∘ away from parallactic angle, causing a loss of flux at short wavelengths, which we corrected as described above for LDSS-3.

The statistical uncertainties on the MagE spectra were calculated by standard error propagation. We added these in quadrature with additional uncertainties based on the seeing losses, telluric absorption corrections, and the overlap between adjacent spectral orders. To be conservative we applied generous uncertainties for each of these effects. For wavelengths within 20 Å of where orders overlap we assumed a 30% uncertainty on the measured fluxes. We also assumed that the seeing loss and telluric corrections each had an uncertainty of 30%.

S1.3 IMACS Observations

We observed SSS17a with IMACS on the Magellan/Baade telescope approximately 48 hours after its discovery on 2017 Aug. 19-20 using the f/2 camera and the 300 lines/mm grism at a blaze angle of 17.5\mbox{{}^{\circ}}. The spectrum was obtained through a 0.9 arcsec-wide slit providing a spectral resolving power of $R\sim 1000$. We reduced and extracted the IMACS spectrum using standard routines in IRAF, including a telluric correction.

S1.4 MIKE Observations

We obtained spectra of SSS17a totaling 1.03 hours of integration time with the MIKE spectrograph using a $2~{\rm arcsec}\times 5~{\rm arcsec}$ slit beginning at UT 00:18 on 2017 Aug 19. The wide slit ensured that we captured as much light as possible from the fading transient. For light that fills a 2 arcsec slit the resulting resolving power is $R=13000$ in the red ($\lambda>5000$ Å) and $R=16000$ in the blue ($\lambda<5000$ Å). However, the seeing of $\sim 1$ arcsec during the observations provides resolving power a factor of $\sim 2$ higher for the SSS17a spectrum. We reduced these data using the Carnegie Python pipeline . The spectrum has a signal-to-noise ratio of 8 per pixel at 5900 Å and 12 per pixel at 6600 Å and is shown in Figure S1.

S1.5 Calibrating and Dereddening Spectra of SSS17a

We futher calibrate the non-MIKE spectra (LDSS-3, MagE, and IMACS) against photometric measurements of SSS17a by extracting synthetic photometric magnitudes for each filter that was completely contained in the wavelength range covered by the spectrum and for which we could either interpolate the photometric light curves or extrapolate them by no more than 1 hour. The exception is the final (8.46 day) spectrum, where the $i$ band magnitude was linearly extrapolated by 1 day. Then we determined the best-fitting line to the difference between the observed and synthetic photometry as a function of central wavelength, and scaled each spectrum by this fit. The VPH-blue and VPH-red observations on the first night only covered the wavelengths of $g$ and $i$, respectively. Thus, for these two spectra we could not correct any wavelength-dependent flux calibration issues and instead only applied a zero-point correction. The bands used to calibrate each spectrum are listed in Table S1.

We adopt reddening and extinction estimates of $E(B-V)$ = 0.106 and $A_{\rm V}=0.34$ mag for our analysis based on the far-infrared dust maps of . These measurements are consistent with the extinction value of $A_{\rm V}=0.37\pm 0.06$ mag determined by using Pan-STARRS1 stellar colors. As a consistency check, we also estimated the extinction along the line of sight to SSS17a with the Na I D absorption lines in the MIKE data. We modeled the Na lines from the Milky Way with a single Gaussian component for each of the D2 and D1 lines, which provides an accurate fit to the spectrum (Figure S1). The spectrum can also be fit with additional weaker components, but because of the modest signal-to-noise ratio of the data and the presence of residuals from the subtraction of the telluric Na D emission lines at nearly the same velocity, those features are not statistically significant. We measured a best-fitting heliocentric velocity of $4.7\pm 1.2$ km s-1 for this absorbing gas, with a full width at half maximum (FWHM) of $26\pm 3$ km s-1 and equivalent widths (EWs) of $328\pm 52$ mÅ for the D2 line and $256\pm 62$ mÅ for the D1 line. We calculated a Na I column density of $5.2\times 10^{12}$ cm-2 from the EW measurements using a standard curve-of-growth analysis. This column density corresponds to a V-band extinction $A_{\rm V}\approx 0.4$ mag . Although this value is consistent with the results of and , because of the uncertainties involved in translating Na D absorption strength into extinction we do not make further use of it in this paper.

We correct the observed spectra for this foreground reddening using a Milky Way extinction curve . The extinction curve is parameterized by the value $R_{V}\equiv A_{V}/E(B-V)$, where $E(B-V)\equiv A_{B}-A_{V}$ is the selective extinction between the $B$ and $V$ photometric bands ($A_{B}$ is the total extinction in the $B$-band). We assume $R_{V}=3.1$.

S1.6 Synthetic Photometry of SSS17a

We measured synthetic photometry in any Sloan ($griz$) or Johnson/Cousins (BVRI) photometric bandpass whose transmission window falls within the wavelength range of the observed spectrum. In order to facilitate direct comparison to broadband observations of SSS17a, this photometry was performed after correction for slit losses, but before the correction for Milky Way reddening described above. To estimate the uncertainties on the synthetic magnitude measurements we run a Monte Carlo simulation. We randomly redraw the photometric measurements 50,000 times based on the observed photometry , assuming that the photometric uncertainties are normally distributed around the measured magnitudes with a width given by the photometric uncertainties. For each draw, we recalibrate the spectrum to the drawn photometry and then perform synthetic photometry. We adopt the 16th and 84th percentiles of the distribution of synthetic magnitude measurements in the Monte Carlo simulation as the $1\sigma$ uncertainties on each measurement. For the day 8.46 spectrum we do not make any synthetic measurements because the photometry used to calibrate that spectrum was extrapolated from one day earlier. We report the synthetic photometry and associated uncertainties in Table S2.

S1.7 GRB130603B Spectra

Both GRB130603B spectra plotted in Figure 2 were presented in . The GTC/OSIRIS spectrum was taken from that paper while we retrieved and reduced the raw X-Shooter optical and near-infrared spectra using the ESO Reflex environment and the X-Shooter standard pipeline recipes. These spectra are mostly featureless within the noise, except for telluric lines and narrow Mg and Ca absorption.

S2 Constraints on Host Galaxy Absorption Lines

We searched the MIKE and MagE spectra for host galaxy absorption or emission lines. The host galaxy, NGC 4993 , has a redshift of $z=0.00988$ . We are unable to detect any absorption lines associated with the host galaxy.This result is in agreement with the results of VLT/X-Shooter observations . Specifically, we do not detect host galaxy Na D absorption, from which we conclude that all of the extinction along the line of sight to SSS17a is located in the Milky Way. Assuming the same linewidth as for the Milky Way Na D absorption (S1.5) and using the formula given by , we place a 2$\sigma$ upper limit on host galaxy Na D absorption of $92$ mÅ (for either the D2 or D1 lines). In comparison, the detected Na D EWs for GRB130603B were $530\pm 90$ mÅ and $590\pm 80$ mÅ for the D2 and D1 lines, respectively . The signal-to-noise ratio of the MIKE spectrum in the blue is too low to place meaningful constraints on Ca H and K absorption. The lack of any neutral gas in NGC 4993 near the site of the merger might indicate that SSS17a was on the near side of the host galaxy, or that NGC 4993 is a largely gas-free system. We also do not detect H$\alpha$ in either absorption or emission (with a 2$\sigma$ upper limit of 63 mÅ for a linewidth of 25 km s-1) or the O III $\lambda 5007$ Å emission line, although the signal-to-noise ratio of the spectrum near 5000 Å is low because it is close to the wavelength of the dichroic of the spectrograph.

S3 BlackBody Spectral Fitting

We fit blackbody models to the observed dereddened rest-frame spectra from the first night using Markov Chain Monte Carlo (MCMC) methods. For a single MCMC run, the statistical uncertainties described above for the LDSS-3 (S1.1) and MagE (S1.2) spectra define the width of the posterior probability distributions for the blackbody temperatures and radii. We carry out additional Monte Carlo simulations to incorporate the effects of systematic uncertainties on the spectra and the photometric flux measurements to which they are tied (S1.5). For the LDSS-3 spectrum we draw 500 random positions of the source within the slit from a Gaussian distribution with a FWHM of 0.1 arcsec. We scale the observed spectrum by the slit losses from differential atmospheric refraction and seeing for each of the randomly drawn source positions to create 500 Monte Carlo spectra. We also randomly draw 5 $g$ and $i$ magnitudes for each source position based on the measured photometric uncertainties and scale the Monte Carlo spectra accordingly. We then re-run the MCMC with the resulting spectra. The reported uncertainties on the BB temperatures and radii are the 90% confidence intervals from these 2500 Monte Carlo iterations. Because the MagE spectrum was obtained at parallactic angle there is no systematic component to the spectroscopic uncertainties. The uncertainties related to the photometry still apply, so we randomly draw 2000 $g$ and $i$ magnitudes, scale the spectrum, re-run the MCMC on the Monte Carlo spectra, and define the uncertainties as above.

Figure 1 shows the observed spectra as well as shaded regions with lower and upper bounds corresponding to the blackbody spectra expected for the 5th and 95th percentile temperature and radii obtained from the fits, respectively.

S4 Existing Kilonova Models

There are many theoretical models of kilonovae in the literature, but they have been developed with few observational constraints. To investigate whether any aspects of the SSS17a spectra agree with theoretical predictions, we searched the three main classes of published models: i) classic lanthanide-rich red kilonova , ii) lanthanide-poor accretion disk winds following the merger , and iii) a model including both hydrodynamical ejecta and wind from a black hole (BH)-neutron star (NS) merger . We find that no existing model we considered simultaneously produced satisfactory fits to the spectroscopic features, color, evolution, and luminosity of our spectroscopic time series. However, studying the ways in which models do agree often leads to additional physical insights. Therefore, we re-examined the models for qualitative similarities with the data by scaling the model luminosities at each epoch to match the observed spectroscopic times series. Even with that additional freedom almost all models were unsatisfactory. However, for the two models shown in Figures 4A and 4B some features resemble the spectra of SSS17a. The NS-NS merger model of has 0.1 solar masses of ejecta composed of Ca, Fe, and Nd distributed in a broken power law density profile with $v=0.2c$. The resulting spectra exhibit smooth red continua with a similar shape to the observed spectra from 4.51 days onward, with the exception of the predicted emission feature at $\sim 6000$ Å (Figure 4A). Conversely, the disc wind outflow model of with a NS with lifetime of 0 ms reproduces the spectral slopes during the first 3.5 days. However, it strongly over-predicts absorption features and under-predicts the photospheric velocity after the merger.

Figure S1. High-resolution MIKE spectrum of SSS17a. The plotted wavelength range is centered on the Na I D absorption lines from the interstellar medium of the Milky Way. The red curve is our Gaussian fit to the spectrum, and the gray spectrum below is the residuals from the fit.