Properties of the binary neutron star merger GW170817

On August 17, 2017, the Advanced LIGO and Advanced Virgo gravitational-wave detectors observed a low-mass compact binary inspiral. The initial sky localization of the source of the gravitational-wave signal, GW170817, allowed electromagnetic observatories to identify NGC 4993 as the host galaxy. In this work, we improve initial estimates of the binary's properties, including component masses, spins, and tidal parameters, using the known source location, improved modeling, and recalibrated Virgo data. We extend the range of gravitational-wave frequencies considered down to 23 Hz, compared to 30 Hz in the initial analysis. We also compare results inferred using several signal models, which are more accurate and incorporate additional physical effects as compared to the initial analysis. We improve the localization of the gravitational-wave source to a 90% credible region of $16~\mathrm{deg}^2$. We find tighter constraints on the masses, spins, and tidal parameters, and continue to find no evidence for nonzero component spins. The component masses are inferred to lie between 1.00 and 1.89 $M_\odot$ when allowing for large component spins, and to lie between 1.16 and 1.60 $M_\odot$ (with a total mass $2.73^{+0.04}_{-0.01} \, M_\odot$) when the spins are restricted to be within the range observed in Galactic binary neutron stars. Under minimal assumptions about the nature of the compact objects, our constraints for the tidal deformability parameter $\tilde \Lambda$ are $(0,630)$ when we allow for large component spins, and $300^{+420}_{-230}$ (using a 90% highest posterior density interval) when restricting the magnitude of the component spins, ruling out several equation-of-state models at the 90% credible level. Finally, with LIGO and GEO600 data, we use a Bayesian analysis to place upper limits on the amplitude and spectral energy density of a possible post-merger signal. (Abridged)


I. INTRODUCTION
On August 17, 2017 the advanced gravitational-wave (GW) detector network, consisting of the two Advanced LIGO detectors [1] and Advanced Virgo [2], observed the compact binary inspiral event GW170817 [3] with a total mass less than any previously observed binary coalescence and a matched-filter signal-to-noise ratio (SNR) of 32.4,louder than any signal to date.Followup Bayesian parameter inference allowed GW170817 to be localized to a relatively small sky area of 28 deg 2 and revealed component masses consistent with those of binary neutron star (BNS) systems.In addition, 1.7 s after the binary's coalescence time the Fermi and INTEGRAL gamma-ray telescopes observed the gamma-ray burst GRB 170817A with an inferred sky location consistent with that measured for GW170817 [4], providing initial evidence that the binary system contained neutron star (NS) matter.
According to general relativity, the gravitational waves emitted by inspiraling compact objects in a quasi-circular orbit are characterized by a chirp-like time evolution in their frequency that depends primarily on a combination of the component masses called the chirp mass [32] and secondarily on the mass ratio and spins of the components.In contrast to binary black hole (BBH) systems, the internal structure of the NS also impacts the waveform, and needs to be included for a proper description of the binary evolution.The internal structure can be probed primarily through attractive tidal interactions that lead to an accelerated inspiral.These tidal interactions are small at lower GW frequencies but increase rapidly in strength during the final tens of GW cycles before merger.Although tidal effects are small relative to other effects, their distinct behavior make them potentially measurable from the GW signal [33][34][35][36][37], providing additional evidence for a BNS system and insight into the internal structure of NSs.
In this work we present improved constraints on the binary parameters first presented in [3].These improvements are enabled by (i) re-calibrated Virgo data, (ii) a broader frequency band of 23-2048 Hz as compared to the original 30-2048 Hz band used in [3], (iii) a wider range of more sophisticated waveform models (see Table I), and (iv) knowledge of the source location from EM observations.By extending the bandwidth from 30-2048 Hz to 23-2048 Hz, we gain access to an additional ∼ 1500 waveform cycles compared to the ∼ 2700 cycles in the previous analysis.Overall, our results for the parameters of GW170817 are consistent with, but more precise than, those in the initial analysis [3].The main improvements are (i) improved 90% sky localization from 28 deg 2 to 16 deg 2 without use of EM observations, (ii) improved constraint on inclination angle enabled by independent measurements of the luminosity distance to NGC 4993, (iii) limits on precession from a new waveform model that includes both precession and tidal effects, and (iv) evidence for a nonzero tidal deformability parameter that is seen in all waveform models.Finally, we analyze the potential post-merger signal with an unmodeled Bayesian inference method [38] using data from the Advanced LIGO detectors and the GEO600 detector [39].This allows us to place improved upper bounds on the amount of post-merger GW emission from GW170817 [40].
As in the initial analysis of GW170817 [3], we infer the binary parameters from the inspiral signal while making minimal assumptions about the nature of the compact objects, in particular allowing the tidal deformability of each object to vary independently.In a companion paper [41], we present a complementary analysis assuming that both compact objects are NSs obeying a common equation of state.This results in stronger constraints on the tidal deformabilities of the NSs than we can make under our minimal assumptions, and allows us to constrain the radii of the NSs and make for novel inferences about the equation of state of cold matter at supranuclear densities.
This paper is organized as follows.Section II details the updated analysis, including improvements to the instrument calibration, improved waveform models, and additional constraints on the source location.Section III reports the improved constraints on the binary's sky location, inclination angle, masses, spins, and tidal parameters.Section IV provides upper limits on possible GW emission after the binary merger.Finally, Sec.V summarizes the results and highlights remaining work such as inference of the NS radius and equation of state (EOS).Additional results from a range of waveform models are reported in Appendix A, and an injection and recovery study investigating the systematic errors in our waveform models is given in Appendix B.

A. Bayesian method
All available information about the source parameters ϑ of GW170817 can be expressed as a posterior probability density function (PDF) p( ϑ|d(t)) given the data d(t) from the GW detectors.Through application of Bayes' theorem, this posterior PDF is proportional to the product of the likelihood L(d(t)| ϑ) of observing data d(t) given the waveform model described by ϑ and the prior PDF p( ϑ) of observing those parameters.Marginal- ized posteriors are computed through stochastic sampling using an implementation of Markov-chain Monte Carlo [42,43] available in the LALInference package [44] as part of the LSC Algorithm Library (LAL) [45].
By marginalizing over all but one or two parameters, it is then possible to generate credible intervals or credible regions for those parameters.

B. Data
For each detector we assume that the noise is additive, i.e., a data stream d(t) = h M (t; ϑ)+n(t) where h M (t; ϑ) is the measured gravitational wave strain and n(t) is Gaussian and stationary noise characterized by the one-sided power spectral densities (PSDs) shown in Fig. 1.The PSD is defined as S n ≡ (2/T ) |ñ(f )| 2 where ñ(f ) is the Fourier transform of n(t) and the angle brackets indicate a time average over the duration of the analysis T , in this case the 128 s containing the d(t) used for all results presented in Sec.III.This PSD is modeled as a sum of cubic splines, for the broad-band structure, and Lorentzians, for the line features, using the BayesWave package [38,46] which produces a posterior PDF of PSDs.
Here we approximate the full structure and variation of these posteriors as a point estimate by using a median PSD, defined separately at each frequency.
The analyses presented here use the same data and calibration model for the LIGO detectors as [3], including subtraction of the instrumental glitch present in LIGO-Livingston (cf.Fig. 2 of [3], and [47] for more details) and of other independently measured noise sources as described in [48][49][50][51].The data from Virgo has been recalibrated since the publication of [3], including the subtraction of known noise sources during post-processing of the data, following the same procedure as described in Sec.II of [52].While the assumption of Gaussian and sta-tionary noise in the detectors is not expected to hold over long timescales, the parameters reported from the study of the 128 s included in this work show no significant impact from remaining non-Gaussian or nonstationary features.
The measured strain h M (t; ϑ) may differ from the true GW strain h(t; ϑ) due to measured uncertainties in the detector calibration [53,54].We relate the measured strain to the true GW strain with the expression [55,56] where hM (f ; ϑ) and h(f ; ϑ) are the Fourier transforms of h M (t; ϑ) and h(t; ϑ) respectively.The terms δA(f ; θ cal ) and δφ(f ; θ cal ) are the frequency-dependent amplitude correction factor and phase correction factor respectively, and are each modeled as cubic splines.For each detector, the parameters are the values of δA and δφ at each of ten spline nodes spaced uniformly in log f [57] between 23 Hz and 2048 Hz.
For the LIGO detectors, the calibration parameters θ cal are informed by direct measurements of the calibration uncertainties [53], and are modeled in the same way as in [3] with 1 σ uncertainties of < 7% in amplitude and < 3 deg in phase for LIGO Hanford and < 5% in amplitude and < 2 deg in phase for LIGO Livingston, all allowing for a non-zero mean offset.The corresponding calibration parameters for Virgo follow the assumptions made in [52], with a 1 σ amplitude uncertainty of 8% and a 1 σ phase uncertainty of 3 deg.This is supplemented with an additional uncertainty in the time stamping of the data of 20 µs (to be compared to the LIGO timing uncertainty of < 1 µs [58] already included in the phase correction factor).At each of the spline nodes, a Gaussian prior is used with these 1 σ uncertainties and their corresponding means.By sampling these calibration parameters in addition to the waveform parameters, the calibration uncertainty is marginalized over.This marginalization broadens the localization posterior (Sec.III A), but does not significantly affect the recovered masses, spins, or tidal deformability parameters.

C. Waveform models for binary neutron stars
In this paper we use four different frequency-domain waveform models which are fast enough to be used as templates in LALInference.These waveforms incorporate point-particle, spin, and tidal effects in different ways.We briefly describe them below.Each waveform's key features are stated in detail in Table I, and further tests of the performance of the waveform models can be found in [69].In addition to these frequency domain models, we employ two state-of-the-art timedomain tidal EOB models that also include spin and tidal effects [74,75].These tidal EOB models have shown good agreement in comparison with NR simulations [74][75][76][77] in the late inspiral and improve on the post-Newtonian (PN) dynamics in the early inspiral.However, these implementations are too slow for use in LALInference.We describe these models in Sec.III D when we discuss an alternative parameter-estimation code [78,79].
The TaylorF2 model used in previous work is a purely analytic PN model.It includes point-particle and aligned-spin terms to 3.5PN order as well as leadingorder (5PN) and next-to-leading-order (6PN) tidal effects [34,[59][60][61][62][63][64][65][80][81][82][83][84][85][86][87].The other three waveform models begin with point-particle models and add a fit to the phase evolution from tidal effects, labeled NR-Tidal [68,69], that fit the high-frequency region to both an analytic EOB model [74] and NR simulations [68,88].The SEOBNRT model is based on the aligned-spin pointparticle EOB model presented in [66] using methods presented in [67] to allow fast evaluation in the frequency domain.PhenomDNRT is based on an aligned-spin point-particle model [70,71] calibrated to untuned EOB waveforms [89] and NR hybrids [70,71].Finally, Phe-nomPNRT is based on the point-particle model presented in [72] that includes an effective description of precession effects.In addition to tidal effects, PhenomPNRT also includes the spin-induced quadrupole moment that enters in the phasing at the 2PN order [90].For alignedspin systems, PhenomPNRT differs from PhenomDNRT only in the inclusion of the spin-induced quadrupole moment.We include the EOS dependence of each NS's spininduced quadrupole moment by relating it to the tidal parameter of each NS using the quasi-universal relations of [91].Although this 2PN effect can have a large phase contribution, even for small spins [37], it enters at similar PN order as many other terms.We therefore expect it to be degenerate with the mass ratio and spins.
These four waveform models have been compared to waveforms constructed by hybridizing BNS EOB inspiral waveforms [74,77] with NR waveforms [68,74,76,88] of the late-inspiral and merger.Since only the PhenomP-NRT model includes the spin-induced quadrupole moment, it was found that it has smaller mismatches than PhenomDNRT and SEOBNRT [69].In addition, because PhenomPNRT is the only model that includes precession effects, we use it as our reference model throughout this paper.
In Fig. 2 we show differences in the amplitude and phase evolution between the four models for an equalmass, non-spinning BNS system.The top panel shows the fractional difference in the amplitude ∆A/A between each model and PhenomPNRT, while the bottom panel shows the absolute phase difference |∆Φ| between each model and PhenomPNRT.Because none of the models have amplitude corrections from tidal effects, the amplitude differences between the models are entirely due to the underlying point-particle models.For the nonprecessing system shown here, PhenomPNRT and Phe-nomDNRT agree by construction, and the difference with SEOBNRT is also small.On the other hand, the purely  analytic TaylorF2 model that has not been tuned to NR simulations deviates by up to 30% from the other models.For the phase evolution of non-spinning systems, PhenomDNRT, PhenomPNRT, and SEOBNRT have the same tidal prescription, so the small < ∼ 2 rad phase differences are due to the underlying point-particle models.For non-spinning systems PhenomDNRT and Phe-nomPNRT are the same, but for spinning systems, the spin-induced quadrupole moment included in PhenomP-NRT but not in PhenomDNRT will cause an additional phase difference.For TaylorF2 the difference with respect to PhenomPNRT is due to both the underlying point-particle model and the tidal prescription, and is ∼ 5 rad for non-spinning systems.
For reference, we also show in Fig. 2 the tidal contribution to the phase for the NRTidal models (∆Φ NRTidal ) and the TaylorF2 model (∆Φ TaylorF2 Tides ).For the system here with tidal deformability Λ = 400 (Eq.( 5)), the tidal contribution is larger than the differences due to the underlying point-particle models.

D. Source parameters and choice of priors
The signal model for the quasi-circular inspiral of compact binaries is described by intrinsic parameters which describe the components of the binary, and extrinsic parameters which determine the location and orientation of the binary with respect to the observer.The intrinsic parameters include the component masses m 1 and m 2 , where we take as convention that m 1 ≥ m 2 .The best measured parameter for systems displaying a long inspiral is the chirp mass [32,80,92,93], Meanwhile, ground-based GW detectors actually measure redshifted (detector-frame) masses, and these are  I) with respect to PhenomPNRT after alignment within the frequency interval [30, 30.25]Hz.Note that in particular the alignment between SEOBNRT and PhenomPNRT is sensitive to the chosen interval due to the difference in the underlying BBH-baseline models at early frequencies.We show as an example an equalmass, non-spinning system with a total mass of 2.706M and a tidal deformability of Λ = 400.In the bottom panel, we also show the tidal contribution to the phasing for the Tay-lorF2 and the NRTidal models.This contribution can be interpreted as the phase difference between the tidal waveform models and the corresponding BBH models.The TaylorF2 waveform terminates at the frequency of the innermost stable circular orbit, which is marked by a small dot.
the quantities we state our prior assumptions on.
Detector-frame masses are related to the astrophysically relevant source-frame masses by m det = (1 + z)m, were z is the redshift of the binary [92,94].Dimensionless quantities such as the ratio of the two masses, q = m 2 /m 1 ≤ 1, are thus the same in the detector frame and the source frame.When exploring the parameter space ϑ we assume a prior PDF p( ϑ) uniform in the detector-frame masses, with the constraint that 0.
where m det 1 ≥ m det 2 , and with an additional constraint on the chirp mass, 1.184 M ≤ M det ≤ 2.168 M .These limits were chosen to mimic the settings in [3] to allow for easier comparisons, and were selected originally for technical reasons.The posterior does not have support near those limits.Despite correlations with the prior on the distance to the source, the source masses also have an effectively uniform prior in the region of parameter space relevant to this analysis.
When converting from detector-frame to source-frame quantities, we use the MUSE/VLT measurement of the heliocentric redshift of NGC 4993, z helio = 0.0098 reported in [95,96].We convert this into a geocentric redshift using the known time of the event, yielding z = 0.0099.
The spin angular momenta of the two binary components S i represent six additional intrinsic parameters, and are usually represented in their dimensionless forms For these parameters we have, following [3], implemented two separate priors for the magnitudes of the dimensionless spins, |χ| = χ, of the two objects.In both cases, we assume isotropic and uncorrelated orientations for the spins, and we use a uniform prior for the spin magnitudes, up to a maximum magnitude.In the first case we enforce χ ≤ 0.89 to be consistent with the value used in [3].This allows us to explore the possibility of exotic binary systems.The exact value of this upper limit does not significantly affect results.Meanwhile, observations of pulsars indicate that while the fastest-spinning neutron star has an observed χ < ∼ 0.4 [97], the fastest-spinning BNSs capable of merging within a Hubble time, PSR J0737-3039A [98] and PSR J1946+2052 [99], will at most have dimensionless spins of χ ∼ 0.04 or χ ∼ 0.05 when they merge.Consistent with this population of BNS systems, in the second case we restrict χ ≤ 0.05.
For the waveforms in Table I that do not support spinprecession, the components of the spins aligned with the orbital angular momentum χ 1z and χ 2z still follow the same prior distributions, which are marginalized over the unsupported spin components.We use the labels "highspin" and "low-spin" to refer to analyses that use the prior χ ≤ 0.89 and χ ≤ 0.05, respectively.
The dimensionless parameters Λ i governing the tidal deformability of each component, discussed in greater detail in Sec.III D, are given a prior distribution uniform within 0 ≤ Λ i ≤ 5000 where no correlation between Λ 1 , Λ 2 , and the mass parameters is assumed.If we assume the two components are NSs that obey the same EOS, then Λ 1 and Λ 2 must have similar values when m 1 and m 2 have similar values [100][101][102].This additional constraint is discussed in a companion paper that focuses on the NS EOS [41].
The remaining signal parameters in ϑ are extrinsic parameters which give the localization and orientation of the binary.When we infer the location of the binary from GW information alone (in the Localization section), we use an isotropic prior PDF for the location of the source on the sky.For most of the results presented here, we restrict the sky location to the known position of SSS17a/AT 2017gfo as determined by electromagnetic observations [12].In every case, we use a prior on the distance which assumes a homogeneous rate density in the nearby Universe, with no cosmological corrections applied; in other words, the distance prior grows with the square of the luminosity distance.Meanwhile we use EM observations to reweight our distance posteriors when investigating the inclination of the binary in Sec.III A, and we use the measured redshift factor to the host galaxy NGC 4993 in order to infer source-frame masses from detector frame masses in Sec.III B. For the angle cos θ JN = Ĵ • N, defined for the total angular momentum J and the line of sight N, we assume a prior distribution uniform in cos θ JN [103].To improve the convergence rate of the stochastic samplers, the analyses with the non-precessing waveform models implement a likelihood function where the phase at coalescence is analytically marginalized out [44].

A. Localization
For most of the analyses in this work we assume a priori that the source of GW170817 is in NGC 4993.However, the improved calibration of Virgo data enables better localization of the source of GW170817 from GW data alone.To demonstrate the improved localization we use results from the updated TaylorF2 analysis (the choice of model does not meaningfully affect localization [105]), shown in Fig. 3.We find a reduction in the 90% localization region from 28 deg 2 [3] to 16 deg 2 .This improved localization is still consistent with the associated counterpart SSS17a/AT 2017gfo (see Fig. 3).
For the remainder of this work we incorporate our knowledge of the location of the event.
While fixing the position of the event to the known location within NGC 4993, we infer the luminosity distance from the GW data alone.Combining this with the redshift associated with the Hubble flow at NGC 4993, we measure the Hubble parameter as in [106].Using the PhenomPNRT waveform model, we find that H 0 = 70 +13 −7 km s −1 Mpc −1 (we use maximum a posteriori and 68.3% credible interval for only H 0 in this work) in the highspin case and H 0 = 70 +19 −8 km s −1 Mpc −1 in the lowspin case; both measurements are within the uncertainties seen in Extended Data Table 1 and Extended Data Fig. 2 of [106].As noted in [106,107], when only measuring one polarization of GW radiation from a binary merger, in the absence of strong precession there is a degeneracy between distance and inclination of the binary.When using GW170817 to measure the Hubble constant this degeneracy is the main source of uncertainty.The slightly stronger constraints on H 0 in the high-spin case Errors quoted as x +z −y represent the median, 5% lower limit, and 95% upper limit.Errors quoted as (x, y) are one-sided 90% lower or upper limits, and are used when one side is bounded by a prior.For the masses, m1 is bounded from below and m2 is bounded from above by the equal mass line.The mass ratio is bounded by q ≤ 1.For the tidal parameter Λ, we quote results using a constant (flat) prior in Λ.In the high-spin case we quote a 90% upper limit for Λ, while in the low-spin case we report both the symmetric 90% credible interval and the 90% highest posterior density (HPD) interval, which is the smallest interval that contains 90% of the probability.
The improved localization of GW170817, with the location of the associated counterpart SSS17a/AT 2017gfo.The darker and lighter green shaded regions correspond to 50% and 90% credible regions respectively, and the gray dashed line encloses the previously-derived 90% credible region presented in [3].
arise because under that prior our weak constraint on precession (see Sec. III C) helps to rule out binary inclinations which are closer to edge-on and where preces-sion effects would be measurable, and hence increases the lower bound on the luminosity distance.Meanwhile, the upper bound on the luminosity distance is achieved with face-off binary inclinations, and is nearly the same for both high-spin and low-spin cases.This same weak constraint on precession leads to a tighter constraint on the inclination angle in the highspin case when using the precessing signal model Phe-nomPNRT, θ JN = 152 +21 −27 deg, as compared to the lowspin case.The inclination measurement in the low-spin case, θ JN = 146 +25 −27 deg, agrees with the inferred values for both the high-and low-spin cases of our three waveform models that treat only aligned-spins (see Table IV in Appendix A).This gives further evidence that it is the absence of strong precession effects in the signal, which can only occur in the high-spin case of the precessing model, that leads to tighter constraints on θ JN .This tighter constraint is absent for systems restricted to the lower spins expected from Galactic NS binaries.
Conversely, EM measurements of the distance to the host galaxy can be used to reduce the effect of this degeneracy, improving constraints on the luminosity distance of the binary and its inclination, which may be useful for constraining emission mechanisms.Figure 4 compares our posterior estimates for distance and inclination with no a priori assumptions regarding the distance to the binary (i.e., using a uniform-in-volume prior) to the improved constraints from an EM-informed prior for the distance to the binary.For the EM-informed results we have reweighted the posterior distribution to use a prior in distance following a normal distribution with mean 40.7 Mpc and standard deviation 2.36 Mpc [104].This leads to improved measurements of the inclination an-  .This measurement is consistent for both the high-spin and low-spin cases, since the EM measurements constrain the source of GW170817 to higher luminosity distances and correspondingly more face-on inclination values.They are also consistent with the limits reported in previous studies using afterglow measurements [108] and combined GW and EM constraints [104,109,110] to infer the inclination of the binary.

B. Masses
Owing to its low mass, most of the SNR for GW170817 comes from the inspiral phase, while the merger and post-merger phases happen at frequencies above 1 kHz, where LIGO and Virgo are less sensitive (Fig. 1).This is different than the BBH systems detected so far, e.g.GW150914 [111][112][113][114] or GW170814 [52].The inspiral phase evolution of a compact binary coalescence can be written as a PN expansion, a power series in v/c, where v is the characteristic velocity within the system [87].The intrinsic parameters on which the system depends enter the expansion at different PN orders.Generally speaking, parameters which enter at lower orders have a large impact on the phase evolution, and are thus easier to measure using the inspiral portion of the signal.
The chirp mass M enters the phase evolution at the lowest order, thus we expect it to be the best-constrained among the source parameters [32,80,92,93].The mass ratio q, and consequently the component masses, are instead harder to measure due to two main factors: 1) they are higher-order corrections in the phase evolution, and 2) the mass ratio is partially degenerate with the component of the spins aligned with the orbital angular momentum [92,93,115], as discussed further below.
In Fig. 5 we show one-sided 90% credible intervals of the joint posterior distribution of the two component masses in the source-frame.We obtain m 1 ∈ (1.36, 1.89) M and m 2 ∈ (1.00, 1.36) M in the high-spin case, and tighter constraints of m 1 ∈ (1.36, 1.60) M and m 2 ∈ (1.16, 1.36) M in the lowspin case.These estimates are consistent with, and generally more precise than, those presented in [3].The inferred masses for the components are also broadly consistent with the known masses of Galactic neutron stars observed in BNS systems (see e.g.[116]).
As expected, the detector-frame chirp mass is measured with much higher precision, with These uncertainties are decreased by nearly a factor of two as compared to the value reported in [3] for the detectorframe chirp mass, while the median remains consistent with the 90% credible intervals previously reported.The main source of uncertainty in the source-frame chirp mass comes from the unknown velocity of the source: the line-of-sight velocity dispersion σ v = 170 km s −1 of NGC 4993 reported in [95] translates into an uncertainty on the geocentric redshift of the source z = 0.0099 ± 0.0009, and thereby onto the chirp mass.This dominates over the statistical uncertainty in M and over the sub-percent level uncertainty in the redshift measurement of NGC 4993 reported in [96].The use of the velocity dispersion to estimate the uncertainty in the radial velocity of the source is consistent with the impact of the second supernova on the center-of-mass velocity of the progenitor of GW170817 being relatively small [95,117], especially given that the probable delay time of GW170817 is much longer than the dynamical time of its host galaxy.Both the sources of uncertainty are incorporated into the values reported in Table II, which still correspond to a sub-percent level of precision on the measurement of M. This method of determining M from the detector-frame chirp mass differs from the original method used in [3], and the resulting median value of M lies at the edge of the 90% credible interval reported there, with uncertainties reduced by a factor of two or more.The fact that chirp mass is estimated much better than the individual masses is the reason why in Fig. 5 the two-dimensional posteriors are so narrow in one direction. 1Meanwhile, the unknown velocity of the progenitor of GW170817 impacts the component masses at a sub-percent level, and is neglected in the bounds reported above and in Table II.
The spins of compact objects directly impact the phasing and amplitude of the GW signal through gravitomagnetic interactions (e.g.[118][119][120]), and through additional contributions to the mass-and current-multipole moments which are the sources of GWs (e.g.[87]).This allows for the measurement of the spins of the compact objects from their GW emission.The spins produce two qualitatively different effects on the waveform.
First, the components of spins along the orbital angular momentum L have the effect of slowing down or speeding up the overall rate of inspiral, for aligned-spin components and anti-aligned spin components, respectively [121].The most important combination of spin components along L is a mass-weighted combination called the effective spin, χ eff [122][123][124], defined as This combination contributes to the gravitational wave phase evolution at the 1.5PN order, together with M, q, and an additional spin degree of freedom [92,125].This leads to a degeneracy among these quantities, especially between q and χ eff , which complicates the measurement of both of these parameters from the GW phase.The remaining aligned-spin degree of freedom at 1.5PN order is more important for systems with lower mass ratios [126], while the perpendicular components of the spins first contribute to the phasing at the 2PN order [80,125,127].
Second, the components of the spins perpendicular to the instantaneous direction of L precess due to spin-orbit and spin-spin interactions.This leads to the precession of the orbital plane itself in order to approximately conserve the direction of the total angular momentum, which modulates the GW phasing and signal amplitude measured by a fixed observer [103].One benefit of considering the effective spin χ eff is that it is approximately conserved throughout inspiral, even as the other components of spins undergo complicated precessional dynamics [128].
The precession-induced modulations of the GW amplitude and phase occur on time scales which span many orbital periods.They are most measurable for systems with large spin components perpendicular to L, for systems with smaller mass ratios q, and for systems viewed close to edge-on [129,130], where precession of the orbital plane strongly modulates the observed signal [103].Precession effects are commonly quantified by an effective spin-precession parameter χ p , which is defined as [131] where χ i⊥ are the magnitudes of the components of the dimensionless spins which are perpendicular to L. When considering precessing binaries, we must specify a reference frequency at which spin-related quantities such as As discussed previously we use two choices for priors on component spins, a prior which allows for high spins (χ ≤ 0.89) and one which restricts to lower spin magnitudes (χ ≤ 0.05).The choice of prior has a strong impact on our spin inferences, which in turn influences the inferred component masses through the q-χ eff degeneracy.
Figure 6 shows the marginalized posterior probability distributions for χ eff from the four waveform models, along with the high-spin and low-spin priors.For the high-spin case we find that negative values of χ eff are mostly excluded for all of the models, although small negative χ eff and negligible values are still allowed.Large values of χ eff are also excluded, and the 90% credible interval for PhenomPNRT is χ eff ∈ (−0.00, 0.10).The uncertainty in χ eff is reduced by nearly a factor of two as compared with the more conservative constraint χ eff ∈ (−0.01, 0.17) reported in [3] for this prior, and remains consistent with negligibly small spins.For the low-spin prior, the constraints on negative values of χ eff are nearly identical, but in this case the upper end of the χ eff marginal posterior is shaped by the prior distribu-  7. Marginalized two-dimensional posteriors for the effective spin χ eff and mass ratio q using the PhenomPNRT model for the high-spin prior (blue) and low-spin prior (orange).The 50% (dashed) and 90% (solid) credible regions are shown for the joint posterior.The 90% credible interval for χ eff is shown by vertical lines and the 90% lower limit for q is shown by horizontal lines.1-D marginal distributions have been renormalized to have equal maxima.
tion.The 90% credible interval in the low-spin case for PhenomPNRT is χ eff ∈ (−0.01, 0.02), which is the same range as reported in [3] for the low-spin case.
Figure 7 shows two-dimensional marginalized posteriors for q and χ eff for PhenomPNRT, illustrating the degeneracy between these parameters.The twodimensional posterior distributions are truncated at the boundary q = 1, and when combined with the degeneracy this causes a positive skew in the marginalized χ eff posteriors, as seen in Fig. 6 [132].Compared to the high-spin priors, the low-spin prior on χ eff cuts off smaller values of q, favoring nearly equal-mass systems.
While all of the models provide constraints on the effective spin, only the PhenomPNRT model provides constraints on the spin-precession of the binary.The top panel of Fig. 8 shows the inferred component spin magnitudes and orientations for the high-spin case.In the highspin case, Fig. 8 shows that we rule out large spin components aligned or anti-aligned with L, but the constraints on in-plane spin components are weaker.As such, we can only rule out large values for the effective precession parameter χ p , as seen in the bottom panel of Fig 8, with the upper 90th percentile at 0.53.Nevertheless, in this case we can place bounds on the magnitudes of the component spins; we find that the 90% upper bounds are χ 1 ≤ 0.50 and χ 2 ≤ 0.61, still well above the range of spins inferred for Galactic binary neutron stars.
Figure 9 shows the same quantities as Fig. 8 using the low-spin prior.In this case we primarily constrain the spins to lie in or above the orbital plane at the reference frequency.This is consistent with the inferences on χ eff , which rule out large negative values of χ eff but whose upper bounds are controlled by the prior distribution.Meanwhile, for χ p the upper 90th percentile is at 0.04, which is nearly unchanged between the prior and posterior distributions.The inability to place strong constraints on precession is consistent with an analysis reported in [3] using a precessing model which neglects tidal effects [72].

D. Tidal parameters
In the post-Newtonian formalism, matter effects for non-spinning objects first enter the waveform phase at 5PN order through the tidally induced quadrupolar ( = 2) deformation [133].The amount of deformation is described by the dimensionless tidal deformability of each NS, defined by Λ = (2/3)k 2 [(c2 /G)(R/m)] 5 , where k 2 is the dimensionless = 2 Love number and R is the NS radius.These quantities depend on the NS mass m and EOS.For spinning NSs, matter effects also enter at 2PN due to the spin-induced quadrupole moment as discussed in Sec.II C, and of the models considered here only Phe-nomPNRT implements this effect.
We show marginalized posteriors for the tidal parameters Λ 1 and Λ 2 in Fig. 10 for the four waveform models.For TaylorF2, the results in this work are in general agreement with the values reported in the detection paper that also used the TaylorF2 model [3].However, here we used a lower starting frequency of 23 Hz instead of 30 Hz, resulting in upper bounds on Λ 1 and Λ 2 that are ∼ 10% (for the high-spin prior) and ∼ 20% (for the lowspin prior) smaller than in [3].This improvement occurs because, although most of the tidal effects occur above several hundred Hz as shown in Fig. 2, the tidal parameters still have a weak correlation with the other parameters.Using more low-frequency information improves the measurement of the other parameters, and thus decreases correlated uncertainties in the tidal parameters.
The three waveform models that use the same NR-Tidal prescription produce nearly identical 90% upper limits that are ∼ 10% smaller than those of TaylorF2.This results because the tidal effect for these models is larger than for TaylorF2 as shown in Fig. 2, so the tidal parameters that best fit the data will be smaller to compensate.Including precession and the spin-induced quadrupole moment in the PhenomPNRT model does not noticeably change the results for the tidal parameters compared to the other two models with the NRTidal prescription.Overall, as already found in [3] the NRTidal models have 90% upper limits that are ∼ 20%-30% lower than the TaylorF2 results presented.
For reference, we also show Λ 1 -Λ 2 contours for a representative subset of theoretical EOS models that span 0  the range of plausible tidal parameters using piecewisepolytrope fits from [134]. 2 The values of Λ 1 and Λ 2 are calculated using the samples for the source-frame masses m 1 and m 2 contained in the 90% credible region for Phe-nomPNRT.The widths of these bands are determined  by the small uncertainty in chirp mass.The lengths of these bands are determined by the uncertainty in mass ratio.They have most of their support near the Λ 1 = Λ 2 line corresponding to the equal mass case, and end at the 90% lower limit for the mass ratio.The predicted values of the tidal parameters for the EOSs MS1, MS1b, and H4 lie well outside of the 90% credible region for both the low-spin and high-spin priors, and for all waveform models.This can be compared to Fig. 5 of [3] where H4 was still marginally consistent with the 90% credible region.
The leading tidal contribution to the GW phase evolution is a mass-weighted linear combination of the two tidal parameters Λ [135].It first appears at 5PN order and is defined such that Λ = Λ 1 = Λ 2 when m 1 = m 2 : In Fig. 11 we show marginalized posteriors of Λ for the two spin priors and four waveform models.Because we used flat priors for Λ 1 and Λ 2 , the prior for Λ, and thus the posterior for Λ, goes to zero in the limit Λ → 0. To avoid the misinterpretation that there is no evidence for Λ = 0, we reweight the posterior for Λ by dividing by the prior used, effectively imposing a flat prior in Λ.In practice, this is done by dividing a histogram of the posterior by a histogram of the prior.The resulting histogram is then resampled and smoothed with kernel density estimation.We have verified the validity of the reweighting procedure by comparing the results to runs where we fix Λ 2 = 0 and use a flat prior in Λ.This differs from the reweighting procedure only in the small, next-to-leadingorder tidal effect.
After reweighting there is still some support at Λ = 0.For the high-spin prior, we can only place a 90% upper limit on the tidal parameter, shown in Fig. 11 and listed in Tables II and IV.For the TaylorF2 model, this 90% upper limit can be directly compared to the value reported in [3].We note, however, that due to a bookkeeping error the value reported in [3] should have been 800 instead of 700.Our improved value of 730 is ∼ 10% less than this corrected value.As with the Λ 1 -Λ 2 posterior (Fig. 10), the three models with the NRTidal prescription predict 90% upper limits that are consistent with each other and less than the TaylorF2 results by ∼ 10%.For the lowspin prior, we can now place a two-sided 90% highest posterior density (HPD) credible interval on Λ that does not contain Λ = 0.This 90% HPD interval is the smallest interval that contains 90% of the probability.
The PDFs for the NRTidal waveform models are bimodal.The secondary peak's origin is the subject of further investigation, but it may result from a specific noise realization, as similar results have been seen with injected waveforms with simulated Gaussian noise (see Fig. 4 of [135]).
In Fig. 11 we also show posteriors of Λ (gray PDFs) predicted by the same EOSs as in Fig. 10, evaluated using the masses m 1 and m 2 sampled from the posterior.The sharp cutoff to the right of each EOS posterior corresponds to the equal mass ratio boundary.Again, as in Fig. 10, the EOSs MS1, MS1b, and H4 lie outside the 90% credible upper limit, and are therefore disfavored.
The differences between the high-spin prior and lowspin prior can be better understood from the joint posterior for Λ and the mass ratio q. Figure 12 shows these posteriors for the PhenomPNRT model without reweighting by the prior.For mass ratios near q = 1, the two posteriors are similar.However, the high-spin prior allows for a larger range of mass ratios, and for smaller values of q there is more support for small values of Λ.If we restrict the mass ratio to q > ∼ 0.5, or equivalently m 2 > ∼ 1 M , we find that there is less support for small values of Λ, and the two posteriors for Λ are nearly identical.
To verify that we have reliably measured the tidal parameters, we supplement the four waveforms used in this paper with two time-domain EOB waveform models: SEOBNRv4T [75,136] and TEOBResumS [74].SEOB-NRv4T includes dynamical tides and the effects of the spin-induced quadrupole moment.TEOBResumS incorporates a gravitational-self-force re-summed tidal potential and the spin-induced quadrupole moment.Both models are compatible with state-of-the-art BNS numerical simulations up to merger [77,137].
Unfortunately, these waveform models are too expensive to be used for parameter estimation with LALInference.We therefore use the parallelized, but less validated parameter estimation code RapidPE [78,79].This code uses a different procedure from the standard LALInference code for generating posterior samples and allows for parameter estimation with significantly more expensive waveform models.For each point in the intrinsic parameter space, RapidPE marginalizes over the extrinsic parameters with Monte Carlo integration.12. PDFs for the tidal parameter Λ and mass ratio q using the PhenomPNRT model for the high-spin (blue) and low-spin (orange) priors.Unlike Fig. 11, the posterior is not reweighted by the prior, so the support that is seen at Λ = 0 is due to smoothing from the KDE.The 50% (dashed) and 90% (solid) credible regions are shown for the joint posterior.The 90% credible interval for Λ is shown by vertical lines and the 90% lower limit for q is shown by horizontal lines.
For aligned-spin models, the resulting 6-dimensional intrinsic marginalized posterior is then adaptively sampled and fit with Gaussian process regression.Samples from this fitted posterior are then drawn using a Markov-chain Monte Carlo algorithm.We performed runs with RapidPE using the low-spin prior for three waveform models.The first used the PhenomDNRT waveform for a direct comparison with the LALInference result.The 90% highest posterior density credible interval for Λ is shifted downward from (70, 730) using LALInference to (20, 690) using RapidPE.Although these differences are not negligible, they are still smaller than the differences between different waveform models.The main difference, however, is that Λ has a bimodal structure using LALInference that is not seen with RapidPE.There are several possible reasons for this difference.One possibility is over-smoothing from the Gaussian process regression fit used in RapidPE.Another possibility is differences in data processing when evaluating the likelihood functions for the two codes.In addition, RapidPE does not marginalize over detector calibration uncertainties.However, comparisons using LALInference with and without calibration error marginalization show that this cannot account for the differences between LALInference and RapidPE.Unfortunately, we have not been able to resolve the differences in the shape of the pos-terior.Given its extensive previous use and testing we use LALInference for our main results, and only use RapidPE for exploratory studies, leaving detailed comparisons to future work.For the two EOB waveforms, the 90% highest posterior density credible interval for Λ is (0, 560) for SEOBNRv4T and (10, 690) for TEOBRe-sumS.We note that the lower bound for SEOBNRv4T of Λ = 0 is not the result of railing against the Λ = 0 bound.The posterior value is the same for the lower and the upper limits of the 90% credible interval.
Recently, De et al. performed an independent analysis of the GW data to measure the tidal parameters [138].Their results are broadly consistent with those presented here, but are made under the assumption that the two merging NSs have the same EOS.They assume that the two NSs have identical radii and that the tidal deformability of the individual stars are related by the approximate relation Λ 1 = q 6 Λ 2 , whereas we allow the tidal parameters to vary independently.Furthermore [138] used the TaylorF2 waveform model and restricted the spin components χ 1z and χ 2z to be uniform in [-0.05, 0.05] (similar to our low-spin prior on the spin magnitudes).Their results for their least restrictive mass prior, where they allow m 1 and m 2 to vary uniformly in [1,2]M , are most directly comparable to our results.They obtain a 90% symmetric credible interval of Λ = 310 +679 −234 (Table I of [138]).This can be compared to our equivalent 90% symmetric credible interval for the TaylorF2 waveform and low-spin prior of 340 +580 −240 (Table IV).Note, however, that a symmetric 90% interval will always (by construction) yield a nonzero lower bound on Λ, even if values close to zero are not disfavored by the data.We therefore prefer to use an HPD interval, which may have either a nonzero or zero lower bound, determined by the posterior distribution.A more direct comparison of the results can be done with our companion paper where we assume a common EOS using approximate universal relations as well as directly sampling a parameterized EOS [41, 100-102, 139, 140].

IV. LIMITS ON POST-MERGER SIGNAL
Having used the inspiral phase of the GW signal to constrain the properties of the component bodies, we now place limits on the signal content after the two stars merged to make inferences about the remnant object.The outcome of a BNS coalescence depends on the progenitor masses and the NS EOS.Soft EOSs and large masses result in the prompt formation of a black hole immediately after the merger [141].Stiffer EOS and lower masses result in the formation of a stable or quasistable NS remnant [142,143].A hypermassive NS, whose mass exceeds the maximum mass of a uniformly rotating star but is supported by differential rotation and possibly thermal gradients [142], will survive for < ∼ 1 s, after which time the NS collapses into a black hole [144,145].A supramassive star, whose mass is lower but still ex-ceeds the threshold for non-rotating NSs, will spin down on longer time scales before forming a black hole [146].Finally, extremely stiff EOSs and low masses will result in a stable NS.
We use the BayesWave algorithm [38] to form frequency-dependent upper limits on the strain amplitude and radiated energy by following the approach described in [147].BayesWave models GWs as a superposition of an arbitrary number of elliptically polarized Morlet-Gabor wavelets.This signal model has been found to be capable of accurate waveform reconstruction for a variety of signal morphologies, including short duration post-merger signals [147].The priors of this analysis are expressed in terms of the individual wavelet parameters and on the SNR of each wavelet.Consequently, the priors on the signal amplitude and waveform morphology are derived from the individual wavelet priors, rather than being directly specified.The priors on the wavelet quality factor and phase are flat in (0, 200) and (0, 2π) respectively.The priors on the central frequency and time are determined by the analysis duration and bandwidth described below, while the amplitude prior is determined through the SNR of each wavelet and discussed in more detail in [38].
We use the analysis described in [147] to estimate an upper bound on the amplitude of a putative GW signal assumed to be present but at insufficiently high SNR to generate a statistically significant detection candidate.We use coincident data from the two LIGO detectors and from GEO600 [39], which has comparable sensitivity to Virgo at high frequency.Indeed, the sky-location of GW170817 is particularly favorable for the GEO600 antenna response so that any high-frequency signal component observed by GEO600 will have an SNR greater than or equal to that expected in Virgo.During this period, the Virgo data above 2 kHz suffer from an abundance of spectral lines and transient noise and, therefore, are not included in this analysis.It should also be noted that GEO600 was not in science mode due to investigations into a degraded squeezer phase error point signal leading to a reduced level of squeezing.At the time of the event, the investigations were passive observations.Otherwise, GEO600 was in nominal running condition.The calibration of the LIGO detectors is more uncertain above 2 kHz than at lower frequencies, but is still within 8% in amplitude and 4 deg in phase [53].The GEO600 calibration uncertainty is estimated to be within 15% in amplitude and 15 deg in phase in the 1-4 kHz band.GEO600 was not used in a previous search for high-frequency GW emission due to an insufficient characterization of data quality and analysis tuning, which would have been required for accurate background estimation [40].The analysis reported in this work, by contrast, is a Bayesian characterization of an underlying signal and involves only the 1 s of data around the coalescence time of the merger, which relaxes the data quality requirements somewhat.Furthermore, the analysis configuration has been chosen based on studies of the expected signal (i.e.[147]), and is  not optimized to eliminate statistical outliers in a background distribution.We use a 1 s segment of data centered around the time of coalescence and we restrict the analysis to waveforms whose peak amplitude lies within a 250 ms window at the center of the segment.This window is sufficient to account for statistical or systematic uncertainties in the time-of-coalescence measurement inferred from the inspiral signal, and the total length of segment used encompasses the duration of post-merger signals predicted by numerical simulations for hypermassive NSs that eventually collapse to black holes.The analysis is performed over the 1024-4096 Hz band, which is sufficient to contain the full post-merger spectrum.
We determine the relative evidence for the two models that the on-source data is described by Gaussian noise only, or by Gaussian noise plus a GW signal as described in [148,149].We find that the Gaussian noise model is strongly preferred, with a Bayes factor (evidence ratio) of 256.79 over the signal model.This result is consistent with both prompt collapse to a BH and with a post-merger signal which is too weak to be measurable with our current sensitivity.We further characterize the absence of a detectable signal by forming 90% credible upper limits on three measures of signal strength: (i) the network SNR, evaluated over 1-4 kHz, (ii) the strain amplitude spectral density (ASD), and (iii) the spectral energy density (SED).
We compute the 90% credible upper limit on the network SNR directly using the reconstructed waveform posterior.We exclude signal power in our analysis band with ρ net > 6.7 at the 90% level.The top panel of Fig. 13 reports the upper limits and expectations for the strain ASD induced in the LIGO-Hanford instrument.These limits are formed directly from the posterior probability distribution for the reconstructed waveform in the 1 s of data around the merger time measured from the premerger observations using a coherent analysis of data from all 3 detectors.The noise ASDs for each instrument are shown for comparison.As one would expect in the absence of a signal the upper limits on strain amplitude approximately follow the shape of the noise spectrum.We also overlay a small set of spectra obtained from simulations of BNS mergers using different EOSs with extrinsic parameters (i.e.sky-location, inclination, and distance) determined from the pre-merger analysis.Information about the simulations used is presented in Table III.Depending on the EOS the analysis frequency band might contain significant contributions from the inspiral and merger phases of the coalescence.If the simulated waveforms are truncated at peak amplitude such that we only include the postmerger phase, the network SNR of each waveform is ∼ 0.5.
Finally, the peak-like structures evident in the posterior upper limit are due to low-significance instrumental artefacts and, particularly around 2.4 kHz, a nonstationary spectral line in LIGO-Livingston.The low end of the strain ASD posterior extends to include zero, consistent with the non-detection of a post-merger signal.
We apply a similar procedure to form a frequencydependent 90% credible upper limit on GW energy (see [147] for details).The bottom panel of Fig. 13 shows the 90% credible upper limits on the SED.As with the reconstructed amplitude, the prior on the SED is imposed by the priors on individual wavelet parameters rather than any specific astrophysical argument.SEDs derived from BNS simulations in which the source is held at the distance, sky-location, and orientation of GW170817 are shown for comparison with our upper limits.Our 90% credible upper limit is still too large to make any inference about the EOSs from this part of the signal.Instead, we characterize the sensitivity improvement required to begin to probe astrophysically interesting energy regimes by comparing the peaks of the simulated SEDs to the 90% credible upper limit on the energy radiated at that frequency.
We find that our upper limits on energy are 12-215 times larger than expectations based on our choice of EOS and simulations, shown in Table III.We therefore require amplitude sensitivity to improve by a fac-tor ∼ 3.5-15 compared to our current results in order to probe realistic energy scales for an equivalent event.This should be regarded as a rather conservative estimate of the upgrade required before we can start probing the astrophysically interesting energy regime, as a number of improvements can increase the sensitivity of our analysis.The current methodology described in [147] is agnostic when it comes to the morphology of the postmerger signal.Additional information about the signal, such as its broadband structure or the finite extent of the post-merger peak, could increase the sensitivity of our analysis, making it easier to detect and characterize the post-merger signal.
As stated earlier, the analysis described here complements the previous, more generic high-frequency search in [40].The upper limits here are given by the 90% credible interval of the posterior probability distribution on the signal amplitude spectrum and its power spectral density.The analysis in [40], by contrast, reports the root-sum-squared amplitude that a number of numerical simulations would require in order that 50% of a population of those signals would produce a ranking statistic with false alarm probability of 10 −4 .Nonetheless, one can compare the amplitude sensitivity improvement required such that each analysis begins to probe astrophysically interesting energy scales.The two analyses share a subset of simulated signals: those with the H4, SLy and SFHx EOSs reported in Table III.In [40], sensitivities are quoted in terms of the root-sum-squared amplitude which scales with the square root of the gravitational wave energy.The energy scales probed by [40] for the H4, SLy and SFHx waveforms are respectively 169, 144 and 121 times higher than the values expected from merger simulations with extrinsic parameters of GW170817.In this analysis the best energy upper limits for the same waveforms are 70, 64 and 31 times greater than the peak energies of those waveforms.With the caveat that we are free to compare our limits with the dominant postmerger frequency, the analysis reported here effectively probes a factor of ∼ 2-4 smaller energies.
Sensitivity improvements may come from more stringent and accurate waveform models, serendipitously located sources, as well as improved instrumental highfrequency sensitivity.The Advanced LIGO design sensitivity, for example, is expected to be three times better at high frequencies than has been achieved to date [1,150], while squeezing is expected to improve the sensitivity by another factor of 2 [151].Similarly, the high-frequency sensitivity of Virgo may see as much as a factor ∼ 40 improvement when design sensitivity is achieved [2,150].The post-merger SNR of the simulated waveforms is about 6-8 times smaller than the SNR required for marginal reconstruction of the post-merger signal [147] depending on the EOS and its energy content [152].A similar event observed with the full LIGO-Virgo network operating at design sensitivity would, therefore, offer an opportunity to probe an astrophysically interesting energy regime and may even provide an estimate of the dominant post-merger oscillation frequency and corresponding constraints on the NS EOS.

V. CONCLUSIONS
This work provides the most constraining measurements of the source of GW170817 to date.Without imposing strong astrophysical priors on the masses or spins, we show that the GW data constrains the masses to the range expected for BNS systems and constrains spin components parallel to the orbital angular momentum to be small.The GW data, however, does not significantly constrain the spin components perpendicular to the orbital angular momentum.If there is significant spin, it must lie near the orbital plane of the binary.Imposing a prior on the distance to GW170817 from the known distance to the host NGC 4993 allows us to constrain the inclination angle of the binary, providing insight into the nature of gamma-ray bursts.
Our improved constraints on the tidal deformation of the binary components reduce the upper bounds on this deformation, further ruling out some of the stiffest equation-of-state models.In addition, we find evidence for finite size effects by establishing a lower bound for the tidal deformation parameter Λ when we restrict the spins to be within the ranges observed in Galactic binaries.However, when we allow for large component spins we are still unable to rule out the possibility of no tidal deformation of the component stars, as would occur for example in a surprisingly low-mass binary black hole merger.While the measured properties are consistent with what we expect for binary neutron star systems, we cannot definitively say from GW measurements alone that both components of the binary were indeed neutron stars.
Comparing results from four different waveform models provides assurance that systematic uncertainties are small compared to statistical uncertainties.Improved waveform models, as well as optimizations to the models and parameter estimation codes that allow them to be used, will further reduce systematic uncertainties.We have shown initial results with the RapidPE code and SEOBNRv4T and TEOBResumS waveform models, and found that the measured tidal parameters are consistent with the main results of the paper.Furthermore, updated instrumental calibration could improve constraints further.However, we do not expect these improvements to change the conclusions obtained here.
We have also placed new, morphology-agnostic bounds on the post-merger signal and argue that the Advanced LIGO-Virgo network at design sensitivity could have potentially reconstructed the post-merger signal.
Where there is still significant potential for improved constraints on GW170817 is in the use of additional information in the priors for tidal deformation.In this work we allowed the component tidal parameters to vary independently, implicitly allowing each neutron star to have a different equation of state.One can require the two neutron stars obey the same EOS through the use of binary universal relations [100][101][102] or a parameterized EOS [139,140].This assumption also allows one to place bounds on the radii of the two neutron stars, and results are discussed in a companion paper [41].
[154] 3342 2.1 450 H4 [155] [154] 2541 4.5 320 GNH3 [156] [157] 2522 3.1 380 SLy [158] [157] 3299 5.0 320 SFHx [159] [  13.We report the value of the spectral energy density (SED) from each simulation at the peak frequency f peak and our 90% credible upper limit on the SED at that frequency.The distance DL = 44.74Mpc and inclination θJN = 166.05deg are determined from the point of maximum posterior probability sampled from the PhenomPNRT model.Note that this is the maximum posterior probability sample drawn from the full posterior probability distribution and does not necessarily correspond to the maxima of the 1-D and 2-D marginal distributions shown in figure 4.
Research Corporation, the Ministry of Science and Technology (MOST), Taiwan and the Kavli Foundation.The authors gratefully acknowledge the support of the NSF, STFC, MPS, INFN, CNRS and the State of Niedersachsen/Germany for provision of computational resources.Data associated with the figures in this article, including posterior samples generated using the Phe-nomPNRT model, can be found at dcc.ligo.org/LIGO-P1800061/public. The GW strain data for this event are available at the LIGO Open Science Center [163].This article has been assigned the document number ligo-p1800061.
Appendix A: Source properties from additional waveform models In this Appendix we present additional results for the source properties of GW170817.Table IV presents the same inferred parameters quoted in Table II for the three additional waveform models TaylorF2, PhenomDNRT, and SEOBNRT.As expected from Figs. 5, 6, 10, and 11, the results among the four waveform models are largely consistent with each other.
One exception is the binary inclination angle θ JN : in the high-spin case, the precessing waveform PhenomP-NRT achieves tighter bounds centered around a more face-off orientation than for the low-spin case.As discussed in Sec.III A, we attribute the tighter constraints on θ JN with the fact we disfavor configurations where strong precession effects would be observable, hence prefer values of θ JN closer to face-off (180 deg).Meanwhile, the other three waveform models only treat aligned spins, and so the absence of strong precession does not help improve their inclination measurements.For all four waveforms, in the small-spin case the spins are constrained to sufficiently small values that there can be no strong precession effects, and so again the inclination measurements for the small-spin case are consistent with the alignedspin measurements in the high-spin case.Finally, when we incorporate EM information about the distance to the source of GW170817, we eliminate the portion of the posteriors at closer distances and lower θ JN , achieving consistent inclination constraints across all cases.
The upper bounds on the spin magnitudes χ 1 and χ 2 are also lower for the waveforms which treat aligned spins only.This is as expected, given that only the components χ i,z contribute to the spin magnitudes for the alignedspin runs, and these spin components are constrained by χ eff in the high-spin case and by our prior in the low-spin case.The differences between the remaining inferred parameters among the four waveforms give a sense of the possible size of systematic errors from our signal modeling, although PhenomPNRT includes the greatest number of relevant physical effects, as seen in Table I.
Table V presents the inferred intrinsic parameters of the binary as produced by RapidPE.

Appendix B: Injection and recovery study
The reliability of the parameter estimation techniques used here was studied in detail for the first BBH detection by injecting state-of-the-art numerical waveform models into the data and verifying that the waveform templates correctly recover the injected parameters [114].We perform a similar analysis for GW170817 by injecting a state-of-the-art BNS waveform model using parameters consistent with the data, then verifying that our waveform templates reliably recover the injected values.For BNS systems, we will focus on the additional tidal parameters which are particularly sensitive to errors in the waveform models [135,[164][165][166].We use as our injected waveform model the time-domain alignedspin SEOBNRv4T model [75,136] discussed in Sec.III D. The version used in this study did not include the spininduced quadrupole moment, but later implementations of SEOBNRv4T such as the one used for the results in Table V include this effect.
We inject SEOBNRv4T with the following parameters in Table VI that are consistent with the measured posterior for GW170817: (i) an approximately equal mass, nonspinning case, (ii) an unequal mass ratio (q = 0.67), nonspinning case, and (iii) an approximately equal mass case with a small spin for the primary star.For these systems we choose a reference EOS that is near the peak of the tidal parameter Λ in Fig. 11, APR4, from which to calculate the tidal parameters Λ 1 and Λ 2 .Finally, we also choose (iv) a stiffer parameterized EOS sometimes used in NR simulations, H [167], that is near the maximum allowed value of Λ.In all four cases, we use the same PSD as for GW170817 and inject the waveform with a network SNR of 32, consistent with GW170817.
We inject these waveforms into zero noise data, resulting in posteriors that are independent of the specific noise realization.We use the high-spin prior (χ i ≤ 0.89) and the three aligned-spin waveform models (TaylorF2, Phe-nomDNRT, and SEOBNRT) as templates.
We show in Fig. 14 the recovered tidal parameter Λ when using the soft APR4 EOS.The posteriors for the three templates are peaked near the injected value of Λ.As with Fig. 11, the 90% upper limits are nearly the same for the two waveforms that use the NRTidal description, while the 90% upper limit for TaylorF2 is ∼ 100 higher.
Because the TaylorF2 tidal effect is smaller than that for the NRTidal models, the TaylorF2 model will estimate a larger tidal parameter to compensate for the smaller tidal effect, cf.Fig. 2.
In Fig. 15 we show the recovered tidal parameter using the stiffer H EOS. The width of the posteriors, distance of the peaks from the injected value, and spread in the 90% credible intervals between the waveform models are larger than in Fig. 14, indicating that the statistical error and systematic waveform errors scale with the true tidal parameter.As with the APR4 injections, the credible interval for the NRTidal waveforms agree fairly well with each other, while the credible interval for the TaylorF2 waveform is ∼ 400 larger.
For GW170817 with a network SNR of 32, waveform systematic errors are important but do not dominate over statistical errors.However, as the detectors improve and results from multiple BNS observations are combined, the statistical errors will decrease.In this case systematic waveform errors may become the dominant source of error, and improved waveform modeling will be needed.FIG. 15.Same as Fig. 14, but with injection (iv) using the H EOS. Solid vertical lines represent the 90% HPD credible interval for each waveform.

FIG. 1 .
FIG.1.Power spectral densities (PSDs) of the Advanced LIGO-Advanced Virgo network.Shown, for each detector, is the median PSD computed from a posterior distribution of PSDs as estimated by BayesWave[38,46] using 128 s of data containing the signal GW170817.

FIG. 2 .
FIG.2.Relative amplitude and phase of the employed waveform models starting at 23 Hz (see TableI) with respect to PhenomPNRT after alignment within the frequency interval[30, 30.25]Hz.Note that in particular the alignment between SEOBNRT and PhenomPNRT is sensitive to the chosen interval due to the difference in the underlying BBH-baseline models at early frequencies.We show as an example an equalmass, non-spinning system with a total mass of 2.706M and a tidal deformability of Λ = 400.In the bottom panel, we also show the tidal contribution to the phasing for the Tay-lorF2 and the NRTidal models.This contribution can be interpreted as the phase difference between the tidal waveform models and the corresponding BBH models.The TaylorF2 waveform terminates at the frequency of the innermost stable circular orbit, which is marked by a small dot.
FIG.7.Marginalized two-dimensional posteriors for the effective spin χ eff and mass ratio q using the PhenomPNRT model for the high-spin prior (blue) and low-spin prior (orange).The 50% (dashed) and 90% (solid) credible regions are shown for the joint posterior.The 90% credible interval for χ eff is shown by vertical lines and the 90% lower limit for q is shown by horizontal lines.1-D marginal distributions have been renormalized to have equal maxima.

FIG. 8 .
FIG.8.Top: Inferred spin parameters using the PhenomP-NRT model, in the high-spin case where the dimensionless component spin magnitudes χ < 0.89.Plotted are the probability densities for the dimensionless spin components χ 1 and χ 2 relative to the orbital angular momentum L, plotted at the reference gravitational wave frequency of f = 100 Hz.A tilt angle of 0 • indicates alignment with L. Each pixel has equal prior probability.Bottom: The posterior for the precession parameter χp, plotted together with its prior distribution, also plotted at the reference frequency of f = 100 Hz.The vertical lines represent the 90th percentile for each distribution.

FIG. 9 .
FIG.9.Inferred spin parameters using the PhenomPNRT model as in Fig.8, but in the low-spin case where the dimensionless component spin magnitudes χ < 0.05.The posterior probability densities for the dimensionless spin components and for χp are plotted at the reference gravitational wave frequency of f = 100 Hz.

FIG. 10 .
FIG.10.PDFs for the tidal deformability parameters Λ1 and Λ2 using the high-spin (top) and low-spin (bottom) priors.The blue shading is the PDF for the precessing waveform PhenomPNRT.The 50% (dashed) and 90% (solid) credible regions are shown for the four waveform models.The seven black curves are the tidal parameters for the seven representative EOS models using the masses estimated with the Phe-nomPNRT model, ending at the Λ1 = Λ2 boundary.

FIG. 11 .
FIG. 11.PDFs of the combined tidal parameter Λ for the high-spin (top) and low-spin (bottom) priors.Unlike in Fig.6, the PDFs have been reweighted by dividing by the original prior for Λ (also shown).The 90% HPD credible intervals are represented by vertical lines for each of the four waveform models: TaylorF2, PhenomDNRT, SEOBNRT, and PhenomPNRT.For the high-spin prior, the lower limit on the credible interval is Λ = 0.The seven gray PDFs are those for the seven representative EOSs using the masses estimated with the PhenomPNRT model.Their normalization constants have been rescaled to fit in the figure.For these EOSs, a 1.36M NS has a radius of 10.4 km (WFF1), 11.3 km (APR4), 11.7 km (SLy), 12.4 km (MPA1), 14.0 km (H4), 14.5 km (MS1b), and 14.9 km (MS1).
FIG.12.PDFs for the tidal parameter Λ and mass ratio q using the PhenomPNRT model for the high-spin (blue) and low-spin (orange) priors.Unlike Fig.11, the posterior is not reweighted by the prior, so the support that is seen at Λ = 0 is due to smoothing from the KDE.The 50% (dashed) and 90% (solid) credible regions are shown for the joint posterior.The 90% credible interval for Λ is shown by vertical lines and the 90% lower limit for q is shown by horizontal lines.

FIG. 14 .
FIG.14.Marginalized PDF of Λ for the three aligned-spin waveform models using the high-spin prior of χi < 0.89.As in Fig.11the PDF is reweighted by the prior.The SEOBNRv4T model was injected into zero-noise data with a network SNR of 32.The injected tidal parameter shown by the dotted vertical line was calculated with the APR4 EOS.Top panel: Approximately equal mass and nonspinning.Middle panel: Unequal mass and nonspinning.Bottom panel: Approximately equal mass and primary component spinning.Solid vertical lines represent the 90% upper limit for each waveform.

TABLE II .
Properties for GW170817 inferred using the PhenomPNRT waveform model.All properties are source properties except for the detector frame chirp mass

TABLE III .
Numerical simulations of the 1.35 M -1.35 M binary neutron star mergers with different EOSs shown in Fig.

TABLE VI .
Parameters used for the injected SEOBNRv4T waveform.The chosen masses and spins are consistent with the measured posteriors for GW170817.The tidal parameters are calculated from the mass and chosen EOS.