# Prospects for measuring off-axis spins of binary black holes with Plus-era
gravitational-wave detectors

Alan M. Knee Department of Physics & Astronomy, University of British
Columbia, Vancouver, BC V6T 1Z1, Canada  Jess McIver Department of Physics &
Astronomy, University of British Columbia, Vancouver, BC V6T 1Z1, Canada
Miriam Cabero Department of Physics & Astronomy, University of British
Columbia, Vancouver, BC V6T 1Z1, Canada

(Received 2021 September 29; Revised 2021 December 6; Accepted 2022 January 4)

###### Abstract

The mass and spin properties of binary black holes (BBHs) inferred from their
gravitational-wave signatures reveal important clues about how these systems
form. BBHs originating from isolated binary evolution are expected to have
spins preferentially aligned with their orbital angular momentum, whereas
there is no such preference in binaries formed via dynamical assembly. The
fidelity with which near-future gravitational-wave detectors can measure off-
axis spins will have implications for the study of BBH formation channels. In
this work, we examine the degree to which the Advanced LIGO Plus (A+) and
Advanced Virgo Plus (AdV+) interferometric detectors can measure both aligned
and misaligned spins. We compare spin resolution between the LIGO-Virgo
network operating at either A+/AdV+ (“Plus”) sensitivity or Advanced-era
design (“Design”) sensitivity using simulated BBH gravitational-wave signals
injected into synthetic detector noise. The signals are distributed over the
mass-spin parameter space of likely BBH systems, accounting for the effects of
precession and higher-order modes. We find that the Plus upgrades yield
significant improvements in spin estimation for systems with unequal masses
and moderate or large spins. Using simulated signals modelled after different
types of hierarchical BBH mergers, we also conclude that the Plus detector
network will yield substantially improved spin estimates for 1G+2G binaries
compared to the Design network.

##  1 Introduction

Since the first direct detection of a binary black hole (BBH), GW150914
(Abbott et al., 2016), gravitational waves from coalescing compact binaries
have revealed an emerging population of heavy stellar-mass BBHs. To date, the
LIGO Scientific and Virgo Collaborations (LVC) have reported a total of 90
gravitational-wave candidates111This tally reflects the number of events
determined to have at least a 50% probability of being astrophysical. from the
first three observing runs (Abbott et al., 2019, 2020b, 2021a, 2021b, 2021c)
of the Advanced LIGO (Aasi et al., 2015) and Advanced Virgo (Acernese et al.,
2015) interferometers. The majority of detected events are most likely
associated with BBHs, with total masses ranging from $14\,M_{\odot}$ to
$184\,M_{\odot}$. The sensitivities of gravitational-wave detectors have
steadily improved with each observing run, and the upcoming fourth (O4) and
fifth (O5) observing runs planned for this decade are expected to add hundreds
of new detections (Abbott et al., 2018; Baibhav et al., 2019).

The mass and spin properties of BBH sources inferred from their gravitational-
wave emissions test our theories of how these systems are formed. For
instance, the unusually massive component black holes (BHs) associated with
GW190521 (Abbott et al., 2020c) potentially lie in the putative upper mass
gap222GW190521 may also have an asymmetric mass ratio, which would change the
component mass estimates (Nitz & Capano, 2021; Capano et al., 2021; Estellés
et al., 2021)., approximately between $65$–$140\,M_{\odot}$. Late-stage
stellar evolution theory predicts that stars cannot produce remnant BHs in
this mass range due to the onset of pair-instability supernovae (PSN)
(Belczynski et al., 2016b; Woosley, 2017), implying that GW190521 might be of
dynamical origin (Kimball et al., 2020b; Romero-Shaw et al., 2020a; Bustillo
et al., 2021; Gayathri et al., 2020). This hypothesis is also supported by
weak evidence of orthogonal (in-plane) spin components relative to the orbit.
Although it is currently uncertain whether GW190521 points to a novel
population of BBHs (Abbott et al., 2020d), it has highlighted how the mass and
spin information encoded in gravitational waves can be used to understand how
BBHs form.

In general, BBHs can be grouped into one of two proposed categories of
formation channels: isolated binary evolution, and dynamical assembly. In the
isolated scenario, a BBH is formed from a highly evolved binary star,
independent of the gravitational influence of any outside sources. This
channel may involve a common envelope phase (Bethe & Brown, 1998; Voss &
Tauris, 2003; Belczynski et al., 2002, 2016a), or proceed via chemically
homogeneous evolution (Mandel & de Mink, 2016; de Mink & Mandel, 2016).
Alternatively, a BBH can be formed through repeated dynamical interactions.
Unlike the isolated channel, the dynamical channel is viable only in dense
stellar environments where close encounters are frequent, such as in the cores
of globular clusters and nuclear star clusters (Gerosa & Fishbach, 2021;
Sigurdsson & Hernquist, 1993; Zwart & McMillan, 2002; Rodriguez et al., 2016b,
a; Belczynski & Banerjee, 2020). Accretion disks around active galactic nuclei
are also expected to drive BBH mergers, as BBHs can have their inspirals
accelerated by interacting with the gaseous disk (McKernan et al., 2012;
Bartos et al., 2017). Hierarchical triples represent another possible channel,
in which a binary merger is driven by a perturbing tertiary mass through the
Lidov-Kozai mechanism (Lim & Rodriguez, 2020; Martinez et al., 2020).

Formation channels are distinguished by the imprints they leave on the
properties of BBHs, including the masses and spins of their components, as
well as orbital eccentricity. BBHs formed via the isolated channel are
expected to have spins preferentially aligned with the orbital angular
momentum of the system as a result of tidal interactions or mass transfer
(Kalogera, 2000; Bogdanovic et al., 2007). In dynamical binaries, the spins
should be oriented isotropically, resulting in misaligned systems (Vitale et
al., 2017c; Talbot & Thrane, 2017; Farr et al., 2017, 2018). Spin misalignment
also induces general relativistic spin-precession, causing the orbital plane
to precess about the total angular momentum of the system. Precession
introduces further periodic structure into the gravitational-wave signal by
modulating its amplitude and phase evolution (Apostolatos et al., 1994;
Kidder, 1995). Dynamical binaries may also have eccentric orbits if there is
insufficient time for the orbit to circularize via emission of gravitational
waves before the merger (Romero-Shaw et al., 2019), providing another method
for constraining formation channels (Zevin et al., 2021).

The dynamical channel can facilitate hierarchical mergers, in which the
remnant BHs of previous mergers proceed to form a new (dynamical) binary and
merge once again. This process could potentially stratify the BH population
into different “generations”, resulting in various classes of dynamical
binaries with unique mass and spin properties depending on how the generations
are paired up (Gerosa et al., 2018; Rodriguez et al., 2018, 2019; Sedda et
al., 2020; Callister et al., 2021). BHs that originated from collapsed stars
are thus termed “first generation” (1G) BHs, while the remnant of the merger
of two 1G BHs is “second generation” (2G), and so on. Hierarchical mergers are
expected to have higher average masses and spins, since these higher-
generation BHs inherit much of the mass and angular momentum of their
progenitor binaries (Fishbach et al., 2017; Gerosa & Berti, 2017), and could
theoretically populate the upper mass gap (Kimball et al., 2020a; Gerosa &
Fishbach, 2021). Precise measurements of both the masses and spins will
therefore play a critical role in the study of these formation channels with
gravitational waves.

Although our ability to measure the spins of gravitational-wave sources is
limited by current detector sensitivities, future technological upgrades to
these detectors present an opportunity to achieve higher-quality spin
measurements, potentially allowing us to make more definitive statements about
the formation histories of detected events. Work is currently underway to
bring the Advanced LIGO and Advanced Virgo detectors up to their design
sensitivity settings. This iteration of the LIGO-Virgo network is planned to
be superseded by Advanced LIGO Plus (A+) and Advanced Virgo Plus (AdV+). The
proposed A+/AdV+ upgrades feature a suite of modifications, including improved
mirror coatings and the implementation of frequency-dependent light squeezing,
which will achieve a broad-spectrum reduction in noise across the LIGO-Virgo
observing band (Tse et al., 2019; McCuller et al., 2020), greatly enhancing
the overall sensitivity of the detectors. The A+ network is projected to have
a detection range of 2.5 Gpc for a $30\,M_{\odot}+30\,M_{\odot}$ BBH; about
1.6 times the range of Advanced LIGO at design sensitivity, and twice that of
the last observing run. Similarly, AdV+ will have a BBH range of roughly a
factor of two times its design sensitivity (see Abbott et al. (2018) for an
outline of future observing prospects).

Motivated by the expected sensitivity improvements provided by these upgrades,
we analyze multiple sets of simulated gravitational-wave signals from BBH
sources and compare spin resolution between the “Design” and “Plus” networks,
where “Design” refers to the three-detector LIGO-Virgo network operating at
Advanced-era design sensitivity, and “Plus” refers to the same network with
A+/AdV+ sensitivity. By comparing spin resolution between these two networks,
as well as examining any systematic correlations with their source properties,
we evaluate the implications that these near-future detector upgrades have for
the study of BBH formation channels. Previous studies have helped us to
understand how current and future gravitational-wave detectors can resolve the
spins of BBH systems and the signatures of misalignment/precession (Biscoveanu
et al., 2021; Vitale et al., 2014, 2017b; Vitale & Evans, 2017; Calderón
Bustillo et al., 2017; Stevenson et al., 2017; Gerosa et al., 2018; Pratten et
al., 2020a; Green et al., 2020; Kalaghatgi & Hannam, 2021; Kimball et al.,
2020a). This work expands on this body of literature by explicitly comparing
the ability of these two networks to recover the spin information of a wide
range of sources, and determining the degree to which our spin measurements
improve as a result of the Plus upgrades.

The remainder of this paper is organized as follows: In Section 2, we review
how the spin properties of BBHs are parameterized in gravitational-wave
analyses; in Section 3, we outline our parameter estimation and injection
methods; in Section 4, we present and discuss the results of our injection
analyses, and finally we conclude in Section 5.

##  2 Precessing binary black holes

In general, inspiralling BBHs are fully characterized by eight intrinsic
quantities, assuming quasi-circular orbits. These are the primary (heavier)
and secondary (lighter) component masses, $m_{1}\geq m_{2}$, as measured in
the source frame333The masses measured in the frame of the detector are
cosmologically redshifted, and are related to the “true” source-frame masses
by a redshift factor, $m_{1,2}^{\rm det}=(1+z)m_{1,2}$. In this work, any time
a mass is written without the superscript, it is implied to be a source-frame
mass., plus six parameters for the two three-dimensional spin vectors of the
component BHs, $\mathbf{S}_{1,2}$ (Farr et al., 2014). It is common practice
to parameterize the binary in terms of: the chirp mass,
$\mathcal{M}=(m_{1}m_{2})^{3/5}/(m_{1}+m_{2})^{1/5}$; mass ratio,
$q=m_{1}/m_{2}\geq 1$; dimensionless spin magnitudes444Note some authors write
the spin magnitudes as $\chi_{1,2}$.,
$a_{1,2}=c|\mathbf{S}_{1,2}|/(Gm_{1,2}^{2})\leq 1$; tilt angles,
$t_{1,2}=\arccos(\hat{\mathbf{L}}\cdot\hat{\mathbf{S}}_{1,2})$, between the
spins and the orbital angular momentum, $\mathbf{L}$; and two additional
angles, $\phi_{\rm JL}$ and $\phi_{12}$, that define the in-plane orientations
of the spins. Since the spin tilts and orientation angles are not conserved as
the system evolves, they must be defined at a reference frequency. In this
work, we quote the values of these parameters at $f=20$ Hz, near the lower
limit of the LIGO-Virgo sensitivity band. Additionally, seven extrinsic
parameters are needed to describe the location and orientation of the source
relative to the detector coordinate system. These include the luminosity
distance to the source, $D_{L}$, sky location (right ascension, $\alpha$, and
declination, $\delta$), polarization angle, $\psi$, and an inclination
parameter (measured at a reference frequency), $\theta_{\rm JN}$, between the
total system angular momentum,
$\mathbf{J}=\mathbf{L}+\mathbf{S}_{1}+\mathbf{S}_{2}$, and the detector line-
of-sight vector, $\mathbf{N}$. Finally, we also need to specify the time of
coalescence (in the geocentric frame), $t_{c}$, and the orbital phase at the
chosen reference frequency, $\phi_{\rm ref}$.

As a consequence of their relatively subtle effects on waveform morphology,
the spin magnitudes $a_{1,2}$ are usually difficult to constrain individually,
and are sensitive to prior assumptions (Vitale et al., 2017a; Biscoveanu et
al., 2020; Abbott et al., 2019). The leading-order effect of the spins on the
gravitational waveform is characterized by the mass-weighted sum of the
aligned projections of the spins, known as the effective aligned spin (Damour,
2001; Racine, 2008; Ajith et al., 2011):

$\chi_{\rm eff}=\frac{qa_{1}\cos t_{1}+a_{2}\cos
t_{2}}{1+q}\,,\quad-1\leq\chi_{\rm eff}\leq 1\,.$ (1)

Due to its effect on the inspiral phase (Ajith et al., 2011; Baird et al.,
2013), $\chi_{\rm eff}$ is typically the most well-measured combination of the
spins (Vitale et al., 2017b). Qualitatively, the spins can either accelerate
or delay the onset of the merger (Campanelli et al., 2006; Healy & Lousto,
2018; Ng et al., 2018). Aligned-spin sources ($\chi_{\rm eff}>0$) complete a
greater number of orbits prior to merging than non-spinning sources due to the
need to dissipate more angular momentum, thereby extending the inspiral. Anti-
aligned-spin sources ($\chi_{\rm eff}<0$) instead complete fewer orbits, and
thus have a shorter observable inspiral.

Other spin effects arise when the spins are misaligned with the orbit, which
causes the system to precess. The orbital angular momentum will rotate about
the total angular momentum, which remains roughly fixed under the assumption
that the orbit precesses on a timescale significantly shorter than that of the
inspiral (Kidder, 1995). The approximate effect of precession on the waveform
is captured by the effective precession spin, defined as (Schmidt et al.,
2015):

$\chi_{\rm p}={\rm max}\bigg{(}a_{1}\sin t_{1},\frac{4+3q}{q(4q+3)}a_{2}\sin
t_{2}\bigg{)}\,,\quad 0\leq\chi_{\rm p}\leq 1\,.$ (2)

It is clear from this formula that $\chi_{\rm p}$ will be non-zero if either
of the spins are misaligned ($\sin t_{1,2}\neq 0$). In practice, this
parameter does not capture all of the precession behaviour. Alternative
methods of describing precession are given in e.g. Thomas et al. (2020);
Fairhurst et al. (2020a, b); Gerosa et al. (2021); Akcay et al. (2021).
Precession introduces periodic modulations of the phase and amplitude of the
emitted gravitational-wave signal (Apostolatos et al., 1994; Kidder, 1995), as
it causes the orientation of the orbit to change relative to the line-of-
sight. This breaks the equatorial symmetry of the system, modifying the mode
structure via mode-mixing (Schmidt et al., 2011; Khan et al., 2020) and
breaking parameter degeneracies (Lang & Hughes, 2006), including a well-known
mass-spin degeneracy that afflicts non-precessing systems (Baird et al., 2013;
Chatziioannou et al., 2015; Tiwari et al., 2018; Pratten et al., 2020a).

##  3 Method

In order to systematically compare spin resolution between the Design and Plus
detector networks, we generated a multitude of simulated gravitational-wave
signals (“injections”) from BBH sources, and attempted to recover their
parameters using these two networks. To create the simulated signals, we used
the IMRPhenomXPHM waveform approximant (Pratten et al., 2020b; García-Quirós
et al., 2020); a fully-precessing, quasi-circular inspiral-merger-ringdown
model of coalescing BBHs that supports higher-order multipoles beyond the
dominant $(\ell,|m|)=(2,2)$ quadrupole modes. The same approximant is used for
recovering the parameters later. The signals were injected into realizations
of stationary Gaussian noise555This is a reasonably close approximation to
real LIGO-Virgo noise. coloured by the noise power spectral density (PSD) of
the detector. The analyzed strain data is given by $d(t)=n(t)+h(t)$, where
$n(t)$ and $h(t)$ correspond to the noise and signal strain, respectively
(Abbott et al., 2020a). We assumed a three-detector network consisting of the
LIGO Hanford, LIGO Livingston, and Virgo interferometers, operating at their
Design and Plus settings.

###  3.1 Parameter estimation

We performed parameter estimation using the formalism of Bayesian statistical
inference (Christensen & Meyer, 2001; Thrane & Talbot, 2019). With this
approach, constraints on waveform parameters are derived from their posterior
probability distributions (“posteriors”) conditioned on the data. Given a
waveform model $h$ dependent on source parameters
$\boldsymbol{\theta}=(\mathcal{M},q,\ldots)$, and interferometric strain data
$\mathbf{d}$, the posteriors are inferred using Bayes’ theorem,

$P(\boldsymbol{\theta}|\mathbf{d},h)=\frac{L(\mathbf{d}|\boldsymbol{\theta},h)\pi(\boldsymbol{\theta}|h)}{Z(\mathbf{d}|h)}\,.$
(3)

Here, the likelihood function, $L(\mathbf{d}|\boldsymbol{\theta},h)$, is the
probability of observing the data conditioned on the parameter values, while
$\pi(\boldsymbol{\theta}|h)$ is the prior distribution, representing our
initial assumptions about the parameter distributions. The denominator is a
normalizing factor called the evidence, and is defined as

$Z(\mathbf{d}|h)=\int
L(\mathbf{d}|\boldsymbol{\theta},h)\pi(\boldsymbol{\theta}|h)\,d\boldsymbol{\theta}\,.$
(4)

For the likelihood, we employ a stationary Gaussian noise model weighted by
the detector PSD, expressed in the frequency domain as (Veitch & Vecchio,
2008; Veitch et al., 2015)

$\ln
L(\mathbf{d}|\boldsymbol{\theta},h)\propto-\sum_{k}\frac{2|\tilde{d}_{k}-\tilde{h}_{k}(\boldsymbol{\theta})|^{2}}{S_{n}(f_{k})T}\,,$
(5)

where $f_{k}$ are the frequency bins, $\tilde{d}_{k}$ and $\tilde{h}_{k}$ are
respectively the discrete Fourier transforms of the strain data and waveform
model, $S_{n}(f)$ is the noise PSD, and $T$ is the duration of the analyzed
data segment. Once the posterior distribution is known, constraints on
individual or pairs of parameters can be obtained by marginalizing the
posterior over the other remaining parameters. For a network of multiple
detectors, the full log-likelihood involves another sum over each detector.

Figure 1: Noise curves for the Advanced LIGO and Virgo detectors at their
Design and Plus sensitivity settings. The vertical axis is the PSD of the
noise in units of Hz-1. A downward shift in the noise curve implies higher
sensitivity. The solid (dashed) curves make up what we call the Design (Plus)
network.

The creation of the injections and subsequent analyses were carried out using
the publicly-available Bilby pipeline (Ashton et al., 2019; Romero-Shaw et
al., 2020b), which supports a variety of Monte Carlo and nested sampling
methods to explore the parameter space. In particular, we infer the posteriors
using a Bilby-specific implementation of the nested sampling algorithm Dynesty
(Speagle, 2020; Skilling, 2004, 2006). To improve convergence times, we use a
distance-marginalized likelihood using the techniques of Singer & Price
(2016). Each injection is analyzed assuming Design sensitivity and then again
with Plus so that we have two sets of posteriors to compare. The PSDs for the
Design and Plus networks are obtained from O’Reilly et al. (2018), and are
shown in Figure 1.

###  3.2 Injections

The injection sets are divided into two main parts. In the first part, we
create $10\times 10$ grids of injections spaced discretely in two parameters,
which we list here.

  * •

_Mass-spin_ : This grid surveys a broad region of the mass and spin parameter
space, providing a general picture of how the Design and Plus networks compare
in terms of resolving BBHs with in-plane spins, and how this depends on the
mass ratios and spins. Injections are spaced discretely in mass ratio,
$q\in[1,10]$, and spin magnitude, $a_{1,2}\in[0,0.9]$, with $a_{1}=a_{2}$. The
spin tilts, $t_{1,2}$, are both fixed to $90^{\circ}$, and there is a fixed
total mass of $M=65\,M_{\odot}$ (similar to GW150914).

  * •

_Spin tilt 1G+2G_ : This grid examines the full range of spin tilt angles,
with fixed mass and spins chosen to roughly emulate the properties of a merger
between 1G and 2G BHs. The 2G primary component is assumed to be twice the
mass of the 1G secondary. Additionally, the primary component has high spin
inherited from the angular momentum of the progenitor binary, while the
secondary component is a stellar remnant with low spin. Injections are again
spaced discretely, but in the two cosine spin tilts $\cos t_{1,2}\in[-1,1]$,
at fixed mass ratio, $q=2$, and total mass, $M=100\,M_{\odot}$. The spin
magnitudes are fixed at $a_{1}=0.7$ and $a_{2}=0.1$, the former of which is a
reasonable estimate based on the predicted spin distribution of hierarchically
formed 2G BHs (Fishbach et al., 2017; Gerosa & Berti, 2017). We do not
experiment with different total masses, as preliminary testing showed little
systematic correlation between total mass and spin recovery for
$M>80\,M_{\odot}$.

  * •

_Spin tilt 2G+2G_ : This grid is similar to the 1G+2G grid, but simulating
mergers between two 2G BHs. The two components have comparable masses, and
both have large spins by virtue of being the products of previous mergers.
Injections are spaced discretely in the two cosine spin tilts, $\cos
t_{1,2}\in[-1,1]$, at fixed mass ratio, $q=1.1$, and total mass,
$M=100\,M_{\odot}$. The spin magnitudes are fixed at $a_{1}=0.7$ and
$a_{2}=0.8$.

All extrinsic and remaining intrinsic parameters are randomly sampled from
appropriate distributions: the two intrinsic spin angles, $\phi_{\rm JL}$ and
$\phi_{12}$, are drawn from a uniform distribution between $0$ and $2\pi$; the
polarization angle, $\psi$, is uniform between $0$ and $\pi$; the orbital
reference phase, $\phi_{\rm ref}$, is uniform between $0$ and $2\pi$; the
inclination angle, $\theta_{\rm JN}$, is isotropic on the sphere; the sky
location, $(\alpha,\delta)$, is isotropic on the sky; and the luminosity
distance is drawn from a distribution such that the redshift is uniform in
comoving volume,

$P(z)\propto\frac{1}{1+z}\frac{dV_{c}}{dz}\,,$ (6)

where the extra redshift factor is to account for the transformation between
source and detector-frame time. Distances are drawn out to $\sim 2$ Gpc so
that the bulk of the simulated events have Design network SNRs in the range of
$8$–$60$. The Plus network SNRs are $1.8$–$2.1$ times larger.

The analyzed bandwidth is truncated at a minimum frequency of 15 Hz, which we
deemed appropriate for Design/Plus sensitivity. We also ensure each injection
is at least “detectable” by applying an optimal signal-to-noise ratio (SNR)
(Del Pozzo et al., 2014) cutoff to reject quiet signals. The cutoff imposes
that each injection has a network SNR666The network SNR is the SNR observed in
each detector added in quadrature. across the three detectors of at least 8,
and an SNR of at least 5 in two of the detectors at Design sensitivity.

For the analyses of the two spin tilt grids, uniform priors are adopted for
the detector-frame chirp mass, $\mathcal{M}^{\rm det}\in[30,70]\,M_{\odot}$,
and for the mass-spin grid we use $\mathcal{M}^{\rm det}\in[5,50]\,M_{\odot}$.
For all three grids we use the following priors: inverted mass ratio,
$q^{-1}\in[0.05,1]$; spin magnitudes, $a_{1,2}\in[0,0.95]$; and coalescence
time, $t_{c}\in[-0.1,0.1]$ (where each injection has $t_{c}=0$). Isotropic
priors are used for the remaining angular parameters, as shown in Table I of
Ashton et al. (2019).

In the second part of our study, we analyze an injection set consisting of 200
signals with masses and spins sampled from the astrophysical population models
described in Abbott et al. (2020e). The primary masses are drawn from the
“power law + peak” model; a compound distribution consisting of a descending
power law with a smooth taper and a Gaussian peak:

$\displaystyle P(m_{1})$ $\displaystyle=[(1-\lambda)\xi(m_{1}|-\alpha,m_{\rm
max})$ $\displaystyle+\lambda G(m_{1}|\mu_{m},\sigma_{m}^{2})]S(m_{1}|m_{\rm
min},\delta_{m})\,,$ (7)

The secondary masses are obtained assuming a power law for the mass ratio,
$P(m_{2})\propto(m_{2}/m_{1})^{\beta_{q}}S(m_{2}|m_{\rm min},\delta_{m})$.
Since we also want to examine how the detector networks resolve the spins of
both aligned and misaligned systems, the spin properties are sampled from the
“default spin” model. Here, the two spin magnitudes are drawn from identical
beta distributions, and the spin tilts are drawn from a mixture of two
populations:

$P(\cos t_{1},\cos t_{2})=\frac{1-\zeta}{4}+\zeta G_{t}(\cos t_{1},\cos
t_{2}|\sigma_{t}^{2})\,.$ (8)

The first term corresponds to an isotropic tilt distribution that models
dynamical binaries, where any combination of spin angles is equally likely. In
isolated binaries, the spins are favoured to be in a nearly aligned
configuration, which is modelled by a truncated Gaussian with mean $(\cos
t_{1},\cos t_{2})=(1,1)$ and variance $\sigma_{t}^{2}$ for both spin tilts.
For the definitions and astrophysical motivations behind these models, see
e.g. Kovetz et al. (2017); Talbot & Thrane (2017); Wysocki et al. (2019);
Talbot & Thrane (2018); Abbott et al. (2020e).

We use the median posterior values for the parameters reported in Abbott et
al. (2020e). The extrinsic and remaining intrinsic parameters are sampled in
the same manner as the injection grids described earlier, and we apply the
same SNR cuts to ensure that all the signals qualify as detectable with Design
sensitivity. The minimum frequency is still $15$ Hz, but we extend the chirp
mass prior to $\mathcal{M}^{\rm det}\in[3,100]\,M_{\odot}$ to accommodate the
wider range of masses. Events are simulated out to $\sim 1$ Gpc which achieves
roughly the same SNR distribution as the injection grids777The Plus network
will be able to observe far more distant events, even beyond $z=1$, but for
purposes of determining how well we can measure off-axis spins it is not
necessary to simulate farther signals..

###  3.3 Constraint statistics

The information contained in a posterior distribution can be condensed down to
a variety of metrics, depending on what aspects of the posteriors are
considered most interesting. One commonly used metric is the Bayesian credible
interval (CI), which is defined as the interval $I$ such that some parameter
$\theta$ falls in that interval with probability $\mathcal{C}$,

$\int_{I}P(\theta|\mathbf{d},h)\,d\theta=\mathcal{C}\,.$ (9)

Note that the posterior here has been marginalized over the other parameters.
One can define the CI as an interval about the median. In what follows, we
quote the 90% ($\mathcal{C}=0.9$) CIs, so the intervals are bounded between
the 5th and 95th percentiles of the posteriors. To quantify the precision of
the CI, we can use their widths. We denote these widths (at 90% credibility)
by $\Sigma_{90}(\theta)$. For two posteriors obtained at Design and Plus
sensitivity, we can compare them using the relative change in
$\Sigma_{90}(\theta)$, given by $\Delta_{90}(\theta)=[\Sigma^{\rm
Plus}_{90}(\theta)-\Sigma^{\rm Design}_{90}(\theta)]/\Sigma^{\rm
Design}_{90}(\theta)$.

The CI only tells us about the variance of the posterior, and contains no
information about whether the analysis is recovering the correct parameter
values. Given that we work with injections, where the true values of the
parameters are known, we would like to incorporate this knowledge when
describing the constraints. Furthermore, we inject into Gaussian noise which
can cause posteriors to be offset from the simulated values, especially for
low-SNR signals near the threshold, so we have reason to expect that the
analysis will not always recover the signal parameters accurately. We thus
introduce a generalized variance quantity, referred to as a “cost function”,
which averages the squared difference from the injected value (rather than the
mean) over the posterior distribution:

$C(\theta)=\bigg{[}\int_{-\infty}^{\infty}P(\theta|\mathbf{d},h)(\theta-\theta_{0})^{2}\,d\theta\bigg{]}^{1/2}\,,$
(10)

where $\theta_{0}$ is the injected value. If the posterior mean is equal to
the injected value, the cost function is equal to the posterior standard
deviation. Should the posterior be offset from the injected value, it will be
penalized with a higher cost. Again, we can compare the cost at the two
sensitivities by looking at the relative change in the cost function, $\Delta
C(\theta)=[C^{\rm Plus}(\theta)-C^{\rm Design}(\theta)]/C^{\rm
Design}(\theta)$. Further discussion about the cost function is given in
Appendix A. To facilitate a more faithful comparison between the two
sensitivities, the noise realizations are generated with a seed that is unique
to each injected signal. Thus, the noise realizations across the two networks
retain the same likelihood, but differ through a frequency-dependent amplitude
scaling.888Specifically, this rescaling depends on the relative difference in
the noise PSDs between the Design and Plus networks.

In the following analyses, we use both the 90% CIs and $C(\theta)$ to
characterize parameter recovery. Most of our results are discussed in terms of
the latter, but we also include the posterior widths for selected results in
Appendix B.

##  4 Results

In this section, we present the results of our parameter estimation analyses.
We focus first on the mass-spin and two spin tilt injection grids, noting any
trends between the injected source properties and differences in spin recovery
between the Design and Plus networks. We then turn our attention towards the
injections sampled from population models.

The precision of parameter estimation depends largely upon the SNR of the
signal. The SNR is determined by the signal amplitude (a function of the mass,
distance, and inclination of the system) and detector sensitivity. Intrinsic
properties, namely the mass ratio and the spins, also play a significant role.
For instance, the effective aligned spin, $\chi_{\rm eff}$, will modify the
duration of the inspiral phase, potentially yielding better spin estimates.
Spin estimation also benefits from higher mass ratios, as found in e.g. Vitale
et al. (2014, 2017b). Systems with misaligned spins will also precess, de-
correlating the mass ratio and effective spin and leading to improved
constraints on both the masses and spins. Meanwhile, the relative improvement
between the Design and Plus networks will depend on differences in the PSDs,
as well as the frequency content of the signals.

###  4.1 Mass-spin grid

Figure 2: Distribution of constraints on the primary spin magnitude, $a_{1}$,
as a function of the mass ratio, $q$, for the mass-spin injection grid. Top:
Widths of the 90% CIs at Plus sensitivity, $\Sigma_{90}^{\rm Plus}(a_{1})$.
The histogram shows how the widths are distributed. The colourbars indicate
the injected value of $a_{1}=a_{2}$. The red line marks the median value of
$\Sigma_{90}^{\rm Plus}(a_{1})$. Bottom: Relative change in CI widths between
Design and Plus sensitivity, $\Delta_{90}(a_{1})$, also with a histogram. The
median $\Delta_{90}(a_{1})$ is marked by a red line.

Starting with the mass-spin grid, we show the constraints on the primary spin
magnitude, $a_{1}$, and primary tilt, $t_{1}$, in Figures 2–3 as functions of
the mass ratio. These figures include the widths of the 90% CIs,
$\Sigma_{90}$, for Plus sensitivity and the relative change $\Delta_{90}$ in
these widths compared to Design sensitivity. We can see that the uncertainties
in the primary spin decrease as the mass ratio increases. The largest
uncertainties are found in the equal-mass cases, while for $q\geq 2$ nearly
all the injections have $\Sigma_{90}^{\rm Plus}(a_{1})<0.5$ and
$\Sigma_{90}^{\rm Plus}(t_{1})<90^{\circ}$. This is consistent with how higher
mass ratios improve spin estimation as described above. This effect is
significantly enhanced when higher-order modes are present, which are
amplified by mass asymmetry (Abbott et al., 2020f), breaking a degeneracy
between the mass ratio and spins. The tilt of the primary spin is also better
measured when the spins are larger, as the high-spin injections typically have
the smallest $t_{1}$ uncertainties at fixed $q$. The largest uncertainties for
$t_{1}$ similarly occur in equal-mass systems. Note that even at high mass
ratios, we still encounter a few cases with very large $t_{1}$ uncertainties
($\Sigma_{90}^{\rm Plus}(t_{1})\sim 150^{\circ}$). These injections have
$a_{1}=a_{2}=0$, where the value of the tilt becomes arbitrary and is thus
unconstrained. Otherwise, the uncertainty in $t_{1}$ for injections with non-
zero spins decreases with higher $q$.

Figure 3: Same as Figure 2, now showing CIs for the primary tilt, $t_{1}$.
Figure 4: Recovery of the primary tilt, $t_{1}$, as a function of luminosity
distance, $D_{L}$, for the mass-spin injection grid. The ten non-spinning
injections are not shown because their tilts are not well-defined. Top:
Posterior medians (squares/diamonds) with errorbars indicating the extent of
the 90% CIs (vertical lines). The histogram shows the distributions of the
posterior medians for the Design and Plus networks. The injected value
$t_{1}=90^{\circ}$ is indicated by the solid horizontal line. Middle: Cost
functions, $C(t_{1})$, calculated for each $t_{1}$ posterior, with histograms
showing the distributions for the two networks. The solid horizontal lines are
the medians of the distributions. Bottom: Relative change in cost, $\Delta
C(t_{1})$, between Design and Plus sensitivity for each injected signal. A
horizontal line marks the median.

The medians of $\Delta_{90}(a_{1})$ and $\Delta_{90}(t_{1})$ indicate that the
constraints on the primary spin are generally around 50% narrower at Plus
sensitivity compared to Design. The bottom panels of Figures 2–3 show that the
CI widths change less for equal mass systems compared to unequal masses. We
also find that $t_{1}$ is generally better measured with the Plus network for
higher spin systems. Of the injections with spins $a_{1,2}\geq 0.3$, the
median value of $\Delta_{90}(t_{1})$ is $-0.5$, compared to injections with
spins of 0.1 or 0.2, where the median is $-0.3$. The constraints between
systems with mass ratios $2\leq q\leq 10$ and spins $0.3\leq a_{1,2}\leq 0.9$
are otherwise similar, in that $\Sigma_{90}$ and $\Delta_{90}$ are roughly
stable for injections within these ranges. Thus, we conclude that the Plus
network will yield the largest improvements in spin estimation for systems
with unequal masses and moderate-to-high spins, reducing the uncertainty in
our measurements by on average a factor of two.

We do not notice similar behaviour in the secondary spin, as both its
magnitude and tilt are typically unconstrained. This is largely because the
spin of the lighter component has a smaller effect on the waveform, and is
thus harder to measure. We find that constraints on $a_{2}$ and $t_{2}$ do not
improve with higher mass ratio and spins for the Plus network, and in the few
cases where the secondary spin is better measured, it is for the sources with
the smallest distances/highest SNRs.

The 90% CIs for the primary tilt $t_{1}$ against the luminosity distance to
the source are shown in Figure 4. The injections with zero spins
($a_{1}=a_{2}=0$) are excluded from this figure. We find that the constraints
at both detector sensitivities are generally consistent with the injected
value, with $t_{1}=90^{\circ}$ contained in the 90% CI in roughly 90% of
cases. The middle panel of Figure 4 shows the cost function $C(t_{1})$ for the
primary tilt posteriors at Design and Plus sensitivity. The cost function
reveals that the medians for the Plus network are usually closer to the
injected value. We can see in the lower panel of Figure 4 that $\Delta
C(t_{1})<0$ for most (91%) of the non-spinning injections, confirming that we
generally obtain more accurate estimates with the Plus network compared to
Design.

###  4.2 Spin tilt 1G+2G grid

We find that the Plus network can recover high $a_{1}$ at approximately three
times the rate for Design for BBH systems consistent with a 1G+2G hierarchical
merger. The constraints obtained for $a_{1}$ are shown in Figure 5. In the top
panel, we show the posterior medians and 90% CIs, along with the distributions
of the medians for the two networks. Although both networks have distributions
that peak near the injected value of $a_{1}=0.7$, Design suffers from a longer
tail due to a higher occurrence of constraints that underestimate $a_{1}$. For
the Plus network, we can note qualitatively that the constraints are more
heavily clustered around the true injected value than Design.

Figure 5: Recovery of the primary spin magnitude, $a_{1}$, as a function of
luminosity distance, $D_{L}$, for the spin tilt 1G+2G injection grid. Top:
Posterior medians (squares/diamonds) with errorbars indicating the extent of
the 90% CIs (vertical lines). The histogram shows the distributions of the
posterior medians for the Design and Plus networks. The injected value
$a_{1}=0.7$ is indicated by the solid horizontal line. and the dashed
horizontal lines show the median 5th percentiles of the CIs. Middle: Cost
functions, $C(a_{1})$, calculated for each $a_{1}$ posterior, with histograms
showing the distributions for the two networks. The solid horizontal lines are
the medians of the distributions. Bottom: Relative change in cost, $\Delta
C(a_{1})$, between Design and Plus sensitivity for each injected signal. A
horizontal line marks the median.

In order to compare the amount of support for a high primary spin, we have
drawn the median position of the 5th posterior percentiles (horizontal dashed
lines) in the top panel of Figure 5. At Plus sensitivity this median sits at
0.54, i.e half of the injections have $a_{1}$ constrained above 0.54 with at
least 95% credibility, indicating strong support for a high primary spin. In
the context of formation channels, spin magnitudes $a_{1}\gtrsim 0.5$ are
astrophysically significant, as they coincide with the predicted spin
distribution of hierarchical BHs (Fishbach et al., 2017). The Design network
yields comparable constraints for only 16% of the injections.

The middle panel of Figure 5 shows the cost function of the $a_{1}$ posteriors
for the 1G+2G grid. The downward shift in the cost distribution of the Plus
network relative to Design reinforces that the constraints are generally
narrower at Plus sensitivity for these systems, and that the posteriors peaks
are indeed closer to the injected value. The lower panel in Figure 5 shows the
relative change in cost $\Delta C(a_{1})$ across the two sensitivities, and we
find the cost function is of order 40% lower at Plus sensitivity than Design.

Figure 6: Recovery of the primary spin magnitude, $a_{1}$, as a function of
the primary tilt, $t_{1}$, for the spin tilt 1G+2G injection grid. Transparent
points correspond to individual posteriors. For each set of injections with
the same value of $t_{1}$, we indicate the median of the quantity on the
vertical axis with an opaque point. Top: Cost functions, $C(a_{1})$, for
$a_{1}$. Bottom: Relative change in cost, $\Delta C(a_{1})$, between Design
and Plus sensitivity for each injected signal.

Section 2 explained how the inspiral duration is affected by having
(anti-)aligned spins. Since much spin/precession information is derived from
the inspiral, we expect the spin constraints to be sensitive to the spin
tilts. We illustrate this in Figure 6, where we show the $a_{1}$ cost function
$C(a_{1})$ against the injected values of the tilt $t_{1}$, and highlight the
median values of $C(a_{1})$ across injections with the same $t_{1}$. In the
bottom panel of Figure 6, we see that there is greater improvement in the
constraints for systems with a primary spin in the range
$90^{\circ}<t_{1}<150^{\circ}$, with $|\Delta C(a_{1})|$ being twice as large
for these cases compared to aligned spins.

1G+2G mergers are characterized by having two components that formed via
different channels. Support for large $a_{1}$ is consistent with the
interpretation of the primary component as a 2G BH, and the case for a
hierarchical merger could be strengthened if this were combined with
convincing evidence of a small secondary spin. However, the secondary spin
magnitudes $a_{2}$ are usually unconstrained at both sensitivities, owing to
its small magnitude and the asymmetric mass configuration. At Plus
sensitivity, the median 95th percentile for $a_{2}$ is 0.78. Such a large
upper limit suggests we often cannot tell whether the secondary spin is small
or not.

Figure 7: Violin posterior plots of the two spin tilts, $t_{1}$ and $t_{2}$,
for the spin tilt 1G+2G injection grid. The posteriors are arranged by their
injected $t_{1}$, given on the horizontal axis. The green/purple pluses are
the posterior medians, and the black crosses are the injected values.
Horizontal bars mark the 5th and 95th percentiles of the posteriors, so that
they enclose the 90% CI. Top: Violin plots for $t_{1}$, showing the injections
with equal component tilts, $t_{1}=t_{2}$. Bottom: Violin plots for $t_{2}$,
showing the same subset of the injection set.

Constraints on the two spin tilts are shown in Figure 7 for a subset of the
injections with equal spin tilts $t_{1}=t_{2}$, which are representative of
the wider set of results. For the aligned injections, we are able to recover
the correct alignment at both sensitivities. The posteriors necessarily tend
to zero at the prior boundaries because of the prior used for the spin tilts.
However, because this prior places more support on $t_{1,2}\sim 90^{\circ}$,
the fact that the posteriors have shifted towards $t_{1}\sim 0^{\circ}$
indicates strong support for an aligned primary spin. For the misaligned
cases, the primary tilt is typically well-resolved with both networks, and the
posteriors peak near the injected values. Plus sensitivity yields narrower
constraints on $t_{1}$ compared to Design for both aligned and misaligned
systems.

The case of anti-aligned spins appears to be the most challenging
configuration, as it is here where we find broader posteriors with
comparatively large offsets from the injected values. The weak tilt
constraints for this case are another consequence of the shorter inspirals of
these types of signals.

Similar to the secondary magnitudes, the secondary spin tilts $t_{2}$ are
often unconstrained, and their posteriors show little deviation from the prior
in most cases.

Figure 8: Recovery of the primary tilt, $t_{1}$, as a function of the
luminosity distance, $D_{L}$, for the spin tilt 1G+2G injection grid. The
solid horizontal lines indicate the median of the quantity on the vertical
axis. Top: Cost functions, $C(t_{1})$, calculated for each posterior, with a
histogram showing the distributions for the two networks. Bottom: Relative
change in cost, $\Delta C(t_{1})$, between Design and Plus sensitivity for
each injected signal.

The full set of $t_{1}$ posteriors are represented using their cost functions,
as shown in Figure 8, where we find that $C(t_{1})$ is of order 40% smaller at
Plus sensitivity compared to Design, which is similar to the $a_{1}$
constraints. We also note that there are several injections with $D_{L}<1$ Gpc
that have their primary tilts constrained at the level of $C(t_{1})\sim
3-5^{\circ}$, representing pristine resolution of the tilts. Additionally, 38%
of the injections have their primary tilt constrained at the level of
$C(t_{1})<10^{\circ}$, compared to just 11% with Design. Therefore, the Plus
network improves the rate of events with these exceptionally well-constrained
spins by about a factor of three.

Figure 9: Recovery of the primary tilt, $t_{1}$, as a function of the injected
value of $t_{1}$ for the spin tilt 1G+2G injection grid. Transparent points
correspond to individual posteriors. For each set of injections with the same
value of $t_{1}$, we indicate the median of the quantity on the vertical axis
with an opaque point. Top: Cost functions, $C(t_{1})$, calculated for each
posterior. Bottom: Relative change in cost, $\Delta C(t_{1})$, between Design
and Plus sensitivity for each injected signal.

The relation between the tilt constraints and the injected values is shown in
Figure 9, where we see the tilts are typically better measured when they are
close to orthogonal. At Plus sensitivity, the injections with a misaligned
primary tilt have median $C(t_{1})\sim 10^{\circ}$, and for the aligned and
anti-aligned injections the cost is about twice this amount. Similar to
$a_{1}$, the largest improvements in the measurement of $t_{1}$ occur for
signals with $90^{\circ}<t_{1}<150^{\circ}$, where $\Delta C(t_{1})$ reaches a
minimum.

###  4.3 Spin tilt 2G+2G grid

Figure 10: Recovery of the primary spin magnitude, $a_{1}$, as a function of
luminosity distance, $D_{L}$, for the spin tilt 2G+2G injection grid. Top:
Posterior medians (squares/diamonds) with errorbars indicating the extent of
the 90% CIs (vertical lines). The histogram shows the distributions of the
posterior medians for the Design and Plus networks. The injected value
$a_{1}=0.7$ is indicated by the solid horizontal line. and the dashed
horizontal lines show the median 5th percentiles of the CIs. Middle: Cost
functions, $C(a_{1})$, calculated for each $a_{1}$ posterior, with histograms
showing the distributions for the two networks. The solid horizontal lines are
the medians of the distributions. Bottom: Relative change in cost, $\Delta
C(a_{1})$, between Design and Plus sensitivity for each injected signal. A
horizontal line marks the median.  Figure 11: Same set of plots as in Figure
10, now for the secondary spin magnitude, $a_{2}$, where the injected value is
0.8.

We have seen that spin recovery is more effective when the masses are unequal,
and that larger spins are better constrained. Figures 10–11 show the
constraints on the two spin magnitudes $a_{1}$ and $a_{2}$ for systems
consistent with a 2G+2G hierarchical merger, where the masses are comparable
and each has a large spin. The median values of the 5th posterior percentiles
are shown to ascertain the typical level of support for a high primary or
secondary spin. Starting with the primary spin, the median 5th percentiles for
Design and Plus sensitivity do not suggest particularly strong evidence for a
large primary spin, especially when compared to the 1G+2G grid. Furthermore,
at Plus sensitivity we find that 29% of the injections recover $a_{1}>0.5$
with 95% credibility, compared to 50% in the 1G+2G grid for the same network.
The $a_{2}$ constraints similarly do not suggest a consistent preference for
large secondary spin, and only 26% of the injections have constrained
$a_{2}>0.5$ at the 95% credible level with the Plus network.

Owing to its larger magnitude, the median reduction in $C(a_{2})$ is roughly
20% between the detector two networks, which is an improvement over the 1G+2G
grid. Conversely, $a_{1}$ is not measured as well. The medians of $C(a_{1})$
at Design and Plus sensitivity are roughly 1.5 times their respective values
for the 1G+2G grid. This is in part because of a spin-spin degeneracy between
$a_{1}$ and $a_{2}$ that is more significant for (near) equal-mass systems.
For the 2G+2G injections, Plus sensitivity tends to yield less support for
small secondary spins compared to Design. If this is not accompanied with
reduced posterior support for $q=1$, then this leads to a decrease in support
for a large $a_{1}$ so that $\chi_{\rm eff}$ is conserved. As a result, the
recovered $a_{1}$ is a slight underestimate of the injected value. An example
of this degeneracy playing out is shown in Figure 12, in which the Plus
network does not exclude $q=1$, but instead slightly increases support for
equal masses. In this particular case, the shifting in the $a_{1}$ posterior
results in a higher cost with the Plus network than Design ($\Delta
C(a_{1})>0$). This degeneracy was not as important in the 1G+2G grid because
it is easier to resolve the more asymmetric mass ratios of those injections,
which helps to disentangle the two spin terms in $\chi_{\rm eff}$.

Figure 12: Corner plot showing posteriors for an injection from the 2G+2G grid
with $\Delta C(a_{1})=0.17$, demonstrating the $a_{1}$–$a_{2}$ degeneracy. The
contours represent the 50% (darker shading) and 90% (lighter shading) credible
regions, and the vertical lines in the histograms show the 5th and 95th
percentiles for the one-dimensional posteriors. The injected values are marked
in red.  Figure 13: Median posterior estimates for the primary spin
magnitudes, $a_{1}$ and $a_{2}$, for the 2G+2G injection grid, with histograms
showing how they are distributed. The dashed horizontal and vertical lines are
the median values of the 5th posterior percentiles. The injected values are
indicated by the black square.

An important question is whether we can accurately recover the large spins of
both components simultaneously, as this is a characteristic signature of 2G+2G
mergers. In Figure 13, we plot the median posterior estimates of $a_{1}$ and
$a_{2}$ on the same plane. For 17% of the injections, both magnitudes are
constrained above 0.5 with 95% credibility at Plus sensitivity, compared to 5%
at Design. We show posteriors for one such injection in Figure 14. We see that
the Plus network resolves the slightly unequal mass ratio, and hence the
inference does not suffer as much from the $a_{1}$–$a_{2}$ degeneracy as
discussed earlier. The shifting in the posterior away from $q=1$ also lets us
rule out most of the negative $\chi_{\rm eff}$ region that fell into the 90%
CI at Design sensitivity.

Figure 14: Corner plot showing posteriors for an injection with nearly zero
effective aligned spin, $\chi_{\rm eff}$, from the 2G+2G grid, showing a case
where both spin magnitudes are accurately measured. Contours and colours are
the same as in Figure 12.  Figure 15: Violin posterior plots of the two spin
tilts, $t_{1}$ and $t_{2}$, for the spin tilt 2G+2G injection grid. Shown here
are injections with equal tilts, $t_{1}=t_{2}$. The posteriors are arranged by
their injected $t_{1}$, given on the horizontal axis. The green/purple pluses
are the posterior medians, and the black crosses are the injected values.
Horizontal bars mark the 5th and 95th percentiles of the posteriors, so that
they enclose the 90% CI. Top: Violin plots for $t_{1}$. Bottom: Violin plots
for $t_{2}$.  Figure 16: Recovery of the primary tilt, $t_{1}$, as a function
of effective aligned spin, $\chi_{\rm eff}$, for the spin tilt 2G+2G injection
grid. The solid horizontal lines indicate the median of the quantity on the
vertical axis. Top: Cost functions, $C(t_{1})$, calculated for each posterior,
with a histogram showing the distributions for the two networks. Bottom:
Relative change in cost, $\Delta C(t_{1})$, between Design and Plus
sensitivity for each injected signal.

Next we investigate the constraints on the two spin tilts, shown in Figure 15.
Similar to the spin magnitudes, we find that constraints on the primary tilt
are typically not as well-measured as in the 1G+2G grid. We also notice that
there are larger offsets in some of the estimates, especially when the
injected values of the tilts are greater than $90^{\circ}$. The secondary tilt
is better measured than in the 1G+2G grid, but constraints tend to broaden
again as the tilts become anti-aligned, where we also see that there is less
improvement between the two networks.

Because both of the tilts are varied in this grid, and both components of the
binary contribute greatly to the waveform due to their comparable masses and
large spins, it is difficult to identify any correlations by plotting the cost
function over the injected tilts as for the 1G+2G grid. Rather, we find it
more convenient to show how the constraints vary with the effective aligned
spin $\chi_{\rm eff}$, shown in Figure 16. We only include the primary tilt
since the overall distribution for $t_{2}$ appears similar and conveys the
same trend that we describe here. We see that $\Delta C(t_{1})$ tends towards
lower values for systems with negligible aligned spin, $\chi_{\rm eff}\sim 0$.
The spin parameters are most degenerate around zero $\chi_{\rm eff}$, since
there are a greater number of combinations of $a_{1,2},t_{1,2}$ that yield
zero $\chi_{\rm eff}$ as opposed to a large positive or negative $\chi_{\rm
eff}$, especially when $q\sim 1$ (this is clear from the shape of the
$\chi_{\rm eff}$ distribution assuming isotropic spin priors, see e.g.
Callister (2021)). The improved sensitivity of the Plus network clearly helps
to mitigate this degeneracy, in part by resolving the slight mass inequality
which leads to a more tightly constrained $\chi_{\rm eff}$ compared to Design.

In contrast to the 1G+2G grid, only 9% of the injections have $t_{1}$ measured
with $C(t_{1})<10^{\circ}$ at Plus sensitivity, roughly four times lower than
for the 1G+2G grid, and almost none (1%) have $t_{1}$ constrained at this
level with Design sensitivity. The same is true of $t_{2}$, where we find
similar percentages. The median reduction in cost between the two detector
networks is on the order of 20% for both spin tilts, which for $t_{1}$ is a
factor of two smaller than in the 1G+2G grid. We thus conclude that measuring
spins is significantly more challenging for 2G+2G mergers, as the lower mass
ratio and similar spin magnitudes suppresses the accuracy of spin recovery
both in absolute terms and in terms of the relative improvement between the
Design and Plus networks.

###  4.4 Population study

In this section, we simulate 200 events drawn from the “power law + peak” mass
and “default” spin population models introduced in Section 3.2. We examine how
well the spins can be constrained in this injection set and to what degree the
behaviours observed in the previous grids still apply when simulating a more
realistic population of events.

Figure 17: Source parameter distributions of our population injections. Each
injection is represented by a point in the $M$–$q$, where $M$ is the source-
frame total mass and $q$ is the mass ratio, and a point in the $\chi_{\rm
eff}$–$\chi_{\rm p}$ plane, where $\chi_{\rm eff}$ and $\chi_{\rm p}$ are the
effective aligned and precession spins. The points are coloured by the Design
network SNR.

The intrinsic parameters of the resulting injection set are shown in Figure
17. The set contains total (source-frame) masses ranging from $11\,M_{\odot}$
to $139\,M_{\odot}$. The vast majority (93% of the injections) have mass
ratios between $1\leq q\leq 2$, with about a dozen $q>2$ sources. The
$\chi_{\rm eff}$ distribution is skewed towards positive effective spins due
to a subpopulation with preferentially aligned spin. Systems with both spins
tilted beyond $90^{\circ}$ are rare, hence there are relatively few samples
with negative $\chi_{\rm eff}$. Nonetheless, the injections incorporate
varying degrees of precession, as shown in Figure 17. Most of the injections
(95%) have network SNR between 8 and 60 at Design sensitivity, with just a few
high-SNR injections above 100.

Figure 18: Recovery of the primary spin magnitude, $a_{1}$, plotted against
the mass ratio, $q$, for the population injection set. The solid horizontal
lines indicate the median of the quantity on the vertical axis. Top: Cost
functions $C(a_{1})$, calculated for each posterior, with a histogram showing
the distributions for the two networks. Bottom: Relative change in cost,
$\Delta C(a_{1})$, between Design and Plus sensitivity for each injected
signal.  Figure 19: Same set of plots as Figure 18, now for the primary tilt,
$t_{1}$.

The population’s tendency toward equal masses and low spins will hinder spin
estimation. Indeed, a significant fraction of the injections return
uninformative posteriors for the primary spin, as seen in Figures 18–19 where
we show the distribution of constraints on $a_{1}$ and $t_{1}$. The median
reduction in cost for $a_{1}$ is on the order of 20%, which is comparable to
the 2G+2G grid. For $t_{1}$, the median reduction is $-0.08$, so there is
typically minimal change in the primary tilt constraints between the two
networks. Nonetheless, we start to find more informative spin constraints as
the mass ratio increases, which is consistent with our previous discussions
thus far. When including only the injections with $q>2$, the Plus network
consistently performs better than Design ($\Delta C<0$) in terms of
constraining the primary spin, and the median reduction in cost rises to 40%
and 25% for $a_{1}$ and $t_{1}$, respectively.

As the two spin magnitudes are drawn from the same distribution, some
injections will have $a_{2}$ greater than or comparable to $a_{1}$, which can
allow us to place stronger constraints on the secondary spin (e.g. unlike the
1G+2G grid where $a_{2}$ is small). The median cost functions for both $a_{1}$
and $a_{2}$ are similar to the 2G+2G results ($\sim 20\%$ reduction in cost),
whereas for $t_{2}$ there is essentially no change between the two networks,
even for injections with $q>2$.

It is also interesting to look at how the constraints correlate with source
parameters that were previously fixed in the earlier grids. To make the
correlations clearer, we binned the injections by parameters of interest: the
mass ratio, source-frame total mass, primary spin magnitude, luminosity
distance, and inclination, and calculated the median of the relative change in
the cost function $\Delta C$ for the injections in each bin. These results are
shown in Figure 20, where we show $\Delta C$ for the primary spin magnitude
and tilt, and the two effective spin parameters.

Figure 20: Constraints on the primary spin magnitude, $a_{1}$, tilt, $t_{1}$,
effective aligned spin, $\chi_{\rm eff}$, and effective precession spin,
$\chi_{\rm p}$, for the population set, as functions of $a_{1}$, the injected
mass ratio, $q$, source-frame total mass, $M$, luminosity distance, $D_{L}$,
and inclination of the total angular momentum relative to the line-of-sight,
$\theta_{\rm JN}$. The cost functions for individual posteriors are shown as
the transparent points. The injections are also binned by the source parameter
on the horizontal axis. We calculate the median relative change in cost
$\Delta C$ between Design and Plus for each bin, indicated with opaque points,
which also mark the left-edges of the bins. The $q$-bins are logarithmically
spaced, and the rest are linearly spaced. Some bins were adjusted to ensure
each one contains no less than 10 injections.

The top two rows of Figure 20 show the constraints for the primary spin as a
function of the chosen injected source properties. Constraints on the primary
spin greatly improve with increasing mass ratio, as we expect. The total mass
of the system will affect spin estimation in a variety of ways, being an
important factor in determining the SNR and what fraction of the signal power
comes from the inspiral versus the merger and ringdown. Lower-mass systems are
inspiral-dominated, which is understood to enhance spin estimation. However,
it is for the higher-mass systems, particularly those with $M\gtrsim
50\,M_{\odot}$, that Plus sensitivity recovers the primary spin even more
accurately relative to Design.

Plus sensitivity also recovers $a_{1}$ more accurately for low-spin systems,
and gains in constraints gradually decrease as the primary spin increases
until there is minimal difference between the two networks. This is a
departure from the mass-spin grid, where larger spins led to greater
improvement between the two networks, and so this trend does not appear robust
against more realistic variations in the other intrinsic source properties.

The cost functions for the primary spin reach a minimum at edge-on
inclinations ($\theta_{\rm JN}\sim 90^{\circ}$), as we expect. Gravitational
radiation from edge-on systems can observed from above and below the orbital
plane if it is precessing, which is known to improve spin estimation (Vitale
et al., 2014). We find that, in addition to the cost functions reaching a
minimum for (near) edge-on systems, the largest improvements between Design
and Plus are also found for such systems. This behaviour is present in the
measurement of the primary spin, as well as the effective precession spin,
$\chi_{\rm p}$, showing that stronger measurements of spin/precession are in
part being driven by systems with edge-on orientation.

Figure 21: Fraction of injections in two subsets of the population set where
an aligned and anti-aligned primary spins are simultaneously excluded at 90%
credibility by some symmetric margin $t$, plotted on the horizontal axis. The
colours and shapes distinguish the network and in-plane spin subset.

We also include two more rows of plots for $\chi_{\rm eff}$ and $\chi_{\rm
p}$, as these values are often reported parameter estimation studies.
$\chi_{\rm eff}$ is consistently well-measured with a median reduction in cost
of roughly 50% between the two detector networks throughout the full range of
masses and spin. The trends in $\chi_{\rm p}$ mirror the trends we observed
for the primary spin.

Lastly, we study the rate at which we can confidently resolve misaligned spins
in this injection set. In particular, we look at how often we can exclude an
aligned ($t_{1}=0^{\circ}$) or anti-aligned ($t_{1}=180^{\circ}$) primary spin
as a function of an increasingly aggressive criterion. First, we divide the
population set into two subsets: those with “moderate” in-plane spin,
$45^{\circ}\leq t_{1}\leq 75^{\circ}$ or $105^{\circ}\leq t_{1}\leq
135^{\circ}$, and with “high” in-plane spin, $75^{\circ}\leq t_{1}\leq
105^{\circ}$. Given a margin $t$, we compute the fraction of injections in
each subset for which the 5th percentile of the $t_{1}$ posterior is greater
than $t$, and the 95th percentile is less than $180^{\circ}-t$. These
fractions are shown for each detector network as a function of $t$ in Figure
21. The horizontal axis stops before $t=45^{\circ}$ because this is the
minimum allowable tilt in the moderate subset. We can see that the points at
$t=35^{\circ}$ tell us that at Plus sensitivity, 40% of the high in-plane spin
injections have their primary tilt constrained within the range
$35^{\circ}\leq t_{1}\leq 145^{\circ}$ at the 90% credible level or higher.
Considering the same subset at Design sensitivity, only $\sim 22\%$ have
$t_{1}$ constrained in this same range and credibility.

The fractions tend to $1$ for all the subsets as the margin becomes smaller,
since we have chosen injections with $t_{1}$ well away from the prior
boundaries and it is quite easy to exclude these small regions of the
parameter space with low prior probability. As we increase the margin, we
begin to exclude the less informative posteriors that span most of the prior
range, and the different spin subsets and detector networks begin to diverge
from each other. For margins $30^{\circ}<t<40^{\circ}$, we find that the Plus
network can exclude an aligned or anti-aligned primary spin at $1.5-2$ times
higher rates than Design. Very aggressive margins of around $45^{\circ}$
result in nearly all the injections being excluded, which is why the rates for
Design and Plus converge again. We generally find better improvements between
detector networks for the high in-plane spin subset.

##  5 Summary

Spin misalignment and precession in detected BBHs are interesting potential
indicators of the formation of these systems. The ability of near-future
gravitational-wave detectors to accurately measure the signatures of off-axis
spins will have important implications for the study of these formation
channels. In order to understand what new spin information might be revealed
by the expected improved sensitivity of these detectors, we have performed an
extensive comparison of spin recovery between the Design and Plus three-
detector networks using a large-scale injection campaign, targeting a wide
range of source properties. In the first part of our study, we analyzed grids
of injections in which we systematically varied pairs of parameter in order to
determine how they affect spin estimation. Then, we analyzed an injection set
sampled from astrophysical population models to see which trends from the
injection grids applied to a more realistic set of events.

We found the Plus network yields the largest improvements in spin estimation
over Design in systems with high mass ratios and high spin magnitudes, which
we attributed to the breaking of degeneracy between the masses and spins due
to the higher-order modes and precession. In particular, for systems with mass
ratios $q\geq 2$ we found that the uncertainty in spin measurements is halved
as a result of the improved sensitivity of the Plus network. We also found
that the primary tilt is well-constrained at both sensitivities given
favourable spin configuration (in-plane spins, high precession).

Motivated by the prospects of using spins to infer the formation history of
BBH mergers, we analyzed two injection grids in which we varied the spin tilts
while fixing the masses and spin magnitudes to values reflecting two classes
of hierarchical mergers: 1G+2G and 2G+2G. In the 1G+2G grid, we found that the
Plus network can provide much more meaningful spin constraints for the 2G
primary component, especially with regards to the spin magnitude. The Plus
network constrained $a_{1}\gtrsim 0.5$ at 95% credibility for around half of
the injections. With Design, we could only obtain similar constraints for 16%
of the injections. The primary tilt was consistently well-measured at both
sensitivities for aligned and misaligned orientations. The combination of a
well-constrained magnitude and tilt gained with the Plus network will allow us
to make more confident statements about how individual binaries formed.
Furthermore, we applied a cost function to gauge the precision and accuracy of
our spin measurements, finding that the Plus network can yield constraints
with on average 40% lower cost than Design across the full set of injected
1G+2G signals. The largest improvements come from binaries with primary spin
tilted at or slightly beyond $90^{\circ}$.

The results for the 2G+2G grid show that the spins of such systems are
considerably more difficult to accurately recover, even at Plus sensitivity.
Although constraints on the secondary spin improve relative to the 1G+2G grid,
the primary spin is less well-measured, as the two component spins become
degenerate with each other for near-equal masses. We conclude that the 1G+2G
binaries are better candidates for extracting stronger spin measurements using
the Plus network relative to Design, and we expect these results to be
generalizable to other types of mergers between BHs of different generations.

Additional work is needed to understand what conditions are needed to measure
both spins accurately. Most likely, this will require more sensitive
detectors. Though we only considered near-future detector upgrades that will
come online during this decade, third-generation gravitational-wave detectors,
including LIGO Voyager (Adhikari et al., 2020), Cosmic Explorer (Reitze et
al., 2019), and Einstein Telescope (Maggiore et al., 2020), will provide even
greater signal fidelity. We leave the topic of spin resolution using these
next-generation detectors to future work, complementing Vitale & Evans (2017)
and similar studies.

Finally, we analyzed a simulated population sampled from astrophysical BBH
population models. The correlation between high mass ratio and spin recovery
with the Plus network is robust to more realistic variations in the masses and
spins, however we do not find greater improvement with the Plus network among
high-spinning binaries like we did for the mass-spin grid. We find that the
Plus network consistently provides better constraints on the primary spin for
$q>2$ sources, which is consistent our mass-spin results. Although spin
uncertainties are overall larger compared to the injection grids, we still
find that the Plus network allows us to confidently identify spin misalignment
at $1.5-2$ higher rates than Design for injections with substantial in-plane
spin.

We considered a relatively simple population in which the mass and spin
distributions were independent. However, because the masses and spins are both
influenced by formation channels, a more detailed population might include
multiple subpopulations of sources with different mass and spin properties,
and would serve as another interesting extension of our study.

We thank Scott Oser for suggesting the use of a cost function to measure the
accuracy of parameter recovery. We also thank the LIGO parameter estimation
working group for helpful discussions, including Patricia Schmidt for
providing comments on this manuscript. This material is based upon work
supported by NSF’s LIGO Laboratory (www.ligo.caltech.edu) which is a major
facility fully funded by the National Science Foundation (www.nsf.gov). We
gratefully acknowledge the support of both the LIGO Laboratory and Compute
Canada (www.computecanada.ca) for provision of computing resources. This work
was also supported by the Natural Sciences and Engineering Research Council of
Canada (www.nserc-crsng.gc.ca) through the Discovery Grants and CGS Master’s
programs.

##  Appendix A Cost function

To provide further intuition regarding the cost function introduced in Section
3.3, we show illustrative examples of synthesized posterior distributions with
different cost functions in Figure 22. The posteriors are sampled from
Gaussian distributions with certain means and standard deviations. The cost
function is computed with respect to $\theta_{0}=0$ (taken to be the injected
value in this example). In the case where the posterior mean is the same as
the injected value, the cost function simply returns the standard deviation of
the posterior. If the posterior is offset from the injected value, as is the
case for the green distribution with mean $\theta=0.05$ (so an offset of
$0.05$), the cost function is larger than the standard deviation. In fact, the
cost function is equivalent to the root sum of the variance and the offset
squared, $C(\theta)=\sqrt{\sigma^{2}+{\rm offset}^{2}}\approx 0.07$.

In terms of our ability to infer the source properties of binary coalescences,
the cost function presents different information about our inferences compared
to the widths of the posteriors/CIs. Even if the population were drawn from a
distribution that exactly matched the prior, we would naturally expect some
fraction of the events (for a given credible level) to not have their
parameters accurately recovered due to the occurrence of statistically
improbable noise realizations. For our analyses, we assigned each injected
signal a numeric seed that was used to generate the noise at both Design and
Plus sensitivity. This was done to ensure that the noise realizations had the
same probabilities across the two networks. Thus, if a posterior is offset
from the injected value at Design sensitivity, the cost function will be able
to capture whether that offset improves, worsens, or stays the same as a
result of the enhanced sensitivity provided by the Plus network, and weigh
that alongside any changes in the precision of the measurement (i.e. the width
of the posterior). This way, the cost function can be a useful tool for
flagging cases where the spin parameters are not more accurately recovered
with the Plus network compared to Design, even when the precision of the
measurement increases (see Figure 12 and the surrounding discussion).

Figure 22: Illustrative synthesized posteriors and their cost functions. Each
distribution is sampled from a Gaussian with some mean, $\mu$, and standard
deviation, $\sigma$. The red distribution has $\sigma=0.05$ and $\mu=0$, green
has $\sigma=0.05$ and $\mu=0.05$, and blue has $\sigma=0.1$ and $\mu=0$. The
cost functions are computed using $\theta_{0}=0$.  Figure 23: Distribution of
constraints on the primary spin magnitude, $a_{1}$, as a function of
luminosity distance, $D_{L}$, for the 1G+2G spin tilt grid. Top: Widths of the
90% CIs, $\Sigma_{90}(a_{1})$. The horizontal lines mark the medians of
$\Sigma_{90}(a_{1})$. Bottom: Relative change in CI widths between Design and
Plus sensitivity, $\Delta_{90}(a_{1})$, also with a histogram.  Figure 24:
Distribution of constraints on the primary spin tilt, $t_{1}$, as a function
of luminosity distance, $D_{L}$, for the 1G+2G spin tilt grid. Top: Widths of
the 90% CIs, $\Sigma_{90}(t_{1})$. The horizontal lines mark the medians of
$\Sigma_{90}(t_{1})$. Bottom: Relative change in CI widths between Design and
Plus sensitivity, $\Delta_{90}(t_{1})$, also with a histogram.  Figure 25:
Distribution of constraints on the primary spin magnitude, $a_{1}$, as a
function of luminosity distance, $D_{L}$, for the 2G+2G spin tilt grid. Top:
Widths of the 90% CIs, $\Sigma_{90}(a_{1})$. The horizontal lines mark the
medians of $\Sigma_{90}(a_{1})$. Bottom: Relative change in CI widths between
Design and Plus sensitivity, $\Delta_{90}(a_{1})$, also with a histogram.
Figure 26: Same as Figure 25, showing constraints on the secondary spin
magnitude, $a_{2}$.  Figure 27: Distribution of constraints on the primary
spin magnitude, $a_{1}$, as a function of mass ratio, $q$, for the population
injection set. Top: Widths of the 90% CIs, $\Sigma_{90}(a_{1})$. The
horizontal lines mark the medians of $\Sigma_{90}(a_{1})$. Bottom: Relative
change in CI widths between Design and Plus sensitivity, $\Delta_{90}(a_{1})$,
also with a histogram.  Figure 28: Distribution of constraints on the primary
spin tilt, $t_{1}$, as a function of mass ratio, $q$, for the population
injection set. Top: Widths of the 90% CIs, $\Sigma_{90}(t_{1})$. The
horizontal lines mark the medians of $\Sigma_{90}(t_{1})$. Bottom: Relative
change in CI widths between Design and Plus sensitivity, $\Delta_{90}(t_{1})$,
also with a histogram.

##  Appendix B Posterior widths

The cost function combines together information about the position of the
posterior peak relative to the true injected value and its width. In cases
where the posterior is significantly offset from the injected value of a
parameter, the cost function will be much larger than the standard deviation.
Regardless, the posterior width alone, representing the precision of the
measurement, may still be of interest. Here, we compile a set a of figures
showing the widths of the 90% CIs and their relative change between Design and
Plus sensitivity, $\Sigma_{90}(\theta)$ and $\Delta_{90}(\theta)$, on the
vertical axis instead of $C(\theta)$ and $\Delta C(\theta)$. When discussing
results in terms of the CI widths, our overall assessment of the injection
studies does not significantly change, though we maintain the cost function
still provides useful information for the reasons described above.

Beginning with the spin tilt 1G+2G grid, in Figures 23–24 we show the widths
of the 90% CIs for the primary spin. At Plus sensitivity, 50% of the
injections have $a_{1}$ constrained at $\Sigma_{90}<0.3$. 50% of the
injections also have $t_{1}$ constrained at $\Sigma_{90}<30^{\circ}$. Similar
to the cost function, the widths of the primary spin CIs decreases by nearly
50% on average between Design and Plus sensitivity.

Figures 25–26 show the widths of the 90% CIs for the primary spin for the
2G+2G spin tilt grid. Looking at the medians again, 50% of the injections have
$a_{1}$ constrained at $\Sigma_{90}<0.35$, whereas for $a_{2}$, 50% have it
constrained at $\Sigma_{90}<0.65$. Our previous conclusion that the spins of
1G+2G binaries are better constrained than in 2G+2G binaries, which was based
on the cost function, still holds when looking at the widths of the
constraints themselves. Despite the slightly better $a_{1}$ constraints than
for $a_{2}$, we do observe cases where Plus sensitivity produces a worse
constraint (in the sense of a wider posterior) than Design, whereas all the
Plus network $a_{2}$ constraints are narrower than or equal to the Design
constraint. For instance, this occurs for the nearest signal located at around
200 Mpc, but from the cost function in Figure 10 we can see that despite the
slightly worse precision, the peak is closer to the injected value at Plus
sensitivity than Design, resulting in a lower cost. We also note that, like
the cost function, the spin constraints for 2G+2G binaries has a weaker
correlation with distance compared to 1G+2G binaries.

Primary spin constraints for the population injections are shown in Figures
27-28. When analyzing the constraints using the cost function, we found modest
improvements for the $a_{1}$ constraints, and minimal improvements for the
$t_{1}$ constraints, which remains true when looking at the posteriors widths.
Similarly, systems with unequal masses have their constraints narrow by a
greater amount between Design and Plus sensitivity.

## References

  * Aasi et al. (2015) Aasi, J., et al. 2015, Class. Quant. Grav., 32, 074001, doi: 10.1088/0264-9381/32/7/074001
  * Abbott et al. (2016) Abbott, B., et al. 2016, Phys. Rev. Lett., 116, 061102, doi: 10.1103/PhysRevLett.116.061102
  * Abbott et al. (2019) —. 2019, Phys. Rev. X, 9, 031040, doi: 10.1103/PhysRevX.9.031040
  * Abbott et al. (2018) Abbott, B. P., et al. 2018, Living Rev. Rel., 21, 3, doi: 10.1007/s41114-018-0012-9
  * Abbott et al. (2020a) —. 2020a, Class. Quant. Grav., 37, 055002, doi: 10.1088/1361-6382/ab685e
  * Abbott et al. (2020b) Abbott, R., et al. 2020b.  https://arxiv.org/abs/2010.14527
  * Abbott et al. (2020c) —. 2020c, Phys. Rev. Lett., 125, 101102, doi: 10.1103/PhysRevLett.125.101102
  * Abbott et al. (2020d) —. 2020d, Astrophys. J. Lett., 900, L13, doi: 10.3847/2041-8213/aba493
  * Abbott et al. (2020e) —. 2020e.  https://arxiv.org/abs/2010.14533
  * Abbott et al. (2020f) —. 2020f, Phys. Rev. D, 102, 043015, doi: 10.1103/PhysRevD.102.043015
  * Abbott et al. (2021a) —. 2021a.  https://arxiv.org/abs/2108.01045
  * Abbott et al. (2021b) —. 2021b, Astrophys. J. Lett., 915, L5, doi: 10.3847/2041-8213/ac082e
  * Abbott et al. (2021c) —. 2021c.  https://arxiv.org/abs/2111.03606
  * Acernese et al. (2015) Acernese, F., et al. 2015, Class. Quant. Grav., 32, 024001, doi: 10.1088/0264-9381/32/2/024001
  * Adhikari et al. (2020) Adhikari, R. X., et al. 2020, Class. Quant. Grav., 37, 165003, doi: 10.1088/1361-6382/ab9143
  * Ajith et al. (2011) Ajith, P., et al. 2011, Phys. Rev. Lett., 106, 241101, doi: 10.1103/PhysRevLett.106.241101
  * Akcay et al. (2021) Akcay, S., Gamba, R., & Bernuzzi, S. 2021, Phys. Rev. D, 103, 024014, doi: 10.1103/PhysRevD.103.024014
  * Apostolatos et al. (1994) Apostolatos, T. A., Cutler, C., Sussman, G. J., & Thorne, K. S. 1994, Phys. Rev. D, 49, 6274, doi: 10.1103/PhysRevD.49.6274
  * Ashton et al. (2019) Ashton, G., et al. 2019, Astrophys. J. Suppl., 241, 27, doi: 10.3847/1538-4365/ab06fc
  * Baibhav et al. (2019) Baibhav, V., Berti, E., Gerosa, D., et al. 2019, Phys. Rev. D, 100, 064060, doi: 10.1103/PhysRevD.100.064060
  * Baird et al. (2013) Baird, E., Fairhurst, S., Hannam, M., & Murphy, P. 2013, Phys. Rev. D, 87, 024035, doi: 10.1103/PhysRevD.87.024035
  * Bartos et al. (2017) Bartos, I., Kocsis, B., Haiman, Z., & Márka, S. 2017, Astrophys. J., 835, 165, doi: 10.3847/1538-4357/835/2/165
  * Belczynski & Banerjee (2020) Belczynski, K., & Banerjee, S. 2020, Astron. Astrophys., 640, L20, doi: 10.1051/0004-6361/202038427
  * Belczynski et al. (2016a) Belczynski, K., Holz, D. E., Bulik, T., & O’Shaughnessy, R. 2016a, Nature, 534, 512, doi: 10.1038/nature18322
  * Belczynski et al. (2002) Belczynski, K., Kalogera, V., & Bulik, T. 2002, Astrophys. J., 572, 407, doi: 10.1086/340304
  * Belczynski et al. (2016b) Belczynski, K., et al. 2016b, Astron. Astrophys., 594, A97, doi: 10.1051/0004-6361/201628980
  * Bethe & Brown (1998) Bethe, H. A., & Brown, G. E. 1998, Astrophys. J., 506, 780, doi: 10.1086/306265
  * Biscoveanu et al. (2021) Biscoveanu, S., Isi, M., Varma, V., & Vitale, S. 2021.  https://arxiv.org/abs/2106.06492
  * Biscoveanu et al. (2020) Biscoveanu, S., Isi, M., Vitale, S., & Varma, V. 2020.  https://arxiv.org/abs/2007.09156
  * Bogdanovic et al. (2007) Bogdanovic, T., Reynolds, C. S., & Miller, M. C. 2007, Astrophys. J. Lett., 661, L147, doi: 10.1086/518769
  * Bustillo et al. (2021) Bustillo, J. C., Sanchis-Gual, N., Torres-Forné, A., et al. 2021, Phys. Rev. Lett., 126, 081101, doi: 10.1103/PhysRevLett.126.081101
  * Calderón Bustillo et al. (2017) Calderón Bustillo, J., Laguna, P., & Shoemaker, D. 2017, Phys. Rev. D, 95, 104038, doi: 10.1103/PhysRevD.95.104038
  * Callister (2021) Callister, T. 2021.  https://arxiv.org/abs/2104.09508
  * Callister et al. (2021) Callister, T. A., Haster, C.-J., Ng, K. K. Y., Vitale, S., & Farr, W. M. 2021.  https://arxiv.org/abs/2106.00521
  * Campanelli et al. (2006) Campanelli, M., Lousto, C. O., & Zlochower, Y. 2006, Phys. Rev. D, 74, 041501, doi: 10.1103/PhysRevD.74.041501
  * Capano et al. (2021) Capano, C. D., Cabero, M., Westerweck, J., et al. 2021.  https://arxiv.org/abs/2105.05238
  * Chatziioannou et al. (2015) Chatziioannou, K., Cornish, N., Klein, A., & Yunes, N. 2015, Astrophys. J. Lett., 798, L17, doi: 10.1088/2041-8205/798/1/L17
  * Christensen & Meyer (2001) Christensen, N., & Meyer, R. 2001, Phys. Rev. D, 64, 022001, doi: 10.1103/PhysRevD.64.022001
  * Damour (2001) Damour, T. 2001, Phys. Rev. D, 64, 124013, doi: 10.1103/PhysRevD.64.124013
  * de Mink & Mandel (2016) de Mink, S. E., & Mandel, I. 2016, Mon. Not. Roy. Astron. Soc., 460, 3545, doi: 10.1093/mnras/stw1219
  * Del Pozzo et al. (2014) Del Pozzo, W., Grover, K., Mandel, I., & Vecchio, A. 2014, Class. Quant. Grav., 31, 205006, doi: 10.1088/0264-9381/31/20/205006
  * Estellés et al. (2021) Estellés, H., et al. 2021.  https://arxiv.org/abs/2105.06360
  * Fairhurst et al. (2020a) Fairhurst, S., Green, R., Hannam, M., & Hoy, C. 2020a, Phys. Rev. D, 102, 041302, doi: 10.1103/PhysRevD.102.041302
  * Fairhurst et al. (2020b) Fairhurst, S., Green, R., Hoy, C., Hannam, M., & Muir, A. 2020b, Phys. Rev. D, 102, 024055, doi: 10.1103/PhysRevD.102.024055
  * Farr et al. (2018) Farr, B., Holz, D. E., & Farr, W. M. 2018, Astrophys. J. Lett., 854, L9, doi: 10.3847/2041-8213/aaaa64
  * Farr et al. (2014) Farr, B., Ochsner, E., Farr, W. M., & O’Shaughnessy, R. 2014, Phys. Rev. D, 90, 024018, doi: 10.1103/PhysRevD.90.024018
  * Farr et al. (2017) Farr, W. M., Stevenson, S., Coleman Miller, M., et al. 2017, Nature, 548, 426, doi: 10.1038/nature23453
  * Fishbach et al. (2017) Fishbach, M., Holz, D. E., & Farr, B. 2017, Astrophys. J. Lett., 840, L24, doi: 10.3847/2041-8213/aa7045
  * García-Quirós et al. (2020) García-Quirós, C., Colleoni, M., Husa, S., et al. 2020, Phys. Rev. D, 102, 064002, doi: 10.1103/PhysRevD.102.064002
  * Gayathri et al. (2020) Gayathri, V., Healy, J., Lange, J., et al. 2020.  https://arxiv.org/abs/2009.05461
  * Gerosa & Berti (2017) Gerosa, D., & Berti, E. 2017, Phys. Rev. D, 95, 124046, doi: 10.1103/PhysRevD.95.124046
  * Gerosa et al. (2018) Gerosa, D., Berti, E., O’Shaughnessy, R., et al. 2018, Phys. Rev. D, 98, 084036, doi: 10.1103/PhysRevD.98.084036
  * Gerosa & Fishbach (2021) Gerosa, D., & Fishbach, M. 2021.  https://arxiv.org/abs/2105.03439
  * Gerosa et al. (2021) Gerosa, D., Mould, M., Gangardt, D., et al. 2021, Phys. Rev. D, 103, 064067, doi: 10.1103/PhysRevD.103.064067
  * Green et al. (2020) Green, R., Hoy, C., Fairhurst, S., et al. 2020.  https://arxiv.org/abs/2010.04131
  * Healy & Lousto (2018) Healy, J., & Lousto, C. O. 2018, Phys. Rev. D, 97, 084002, doi: 10.1103/PhysRevD.97.084002
  * Kalaghatgi & Hannam (2021) Kalaghatgi, C., & Hannam, M. 2021, Phys. Rev. D, 103, 024024, doi: 10.1103/PhysRevD.103.024024
  * Kalogera (2000) Kalogera, V. 2000, Astrophys. J., 541, 319, doi: 10.1086/309400
  * Khan et al. (2020) Khan, S., Ohme, F., Chatziioannou, K., & Hannam, M. 2020, Phys. Rev. D, 101, 024056, doi: 10.1103/PhysRevD.101.024056
  * Kidder (1995) Kidder, L. E. 1995, Phys. Rev. D, 52, 821, doi: 10.1103/PhysRevD.52.821
  * Kimball et al. (2020a) Kimball, C., Talbot, C., Berry, C. P., et al. 2020a.  https://arxiv.org/abs/2005.00023
  * Kimball et al. (2020b) Kimball, C., et al. 2020b.  https://arxiv.org/abs/2011.05332
  * Kovetz et al. (2017) Kovetz, E. D., Cholis, I., Breysse, P. C., & Kamionkowski, M. 2017, Phys. Rev. D, 95, 103010, doi: 10.1103/PhysRevD.95.103010
  * Lang & Hughes (2006) Lang, R. N., & Hughes, S. A. 2006, Phys. Rev. D, 74, 122001, doi: 10.1103/PhysRevD.75.089902
  * Lim & Rodriguez (2020) Lim, H., & Rodriguez, C. L. 2020, Phys. Rev. D, 102, 064033, doi: 10.1103/PhysRevD.102.064033
  * Maggiore et al. (2020) Maggiore, M., et al. 2020, JCAP, 03, 050, doi: 10.1088/1475-7516/2020/03/050
  * Mandel & de Mink (2016) Mandel, I., & de Mink, S. E. 2016, Mon. Not. Roy. Astron. Soc., 458, 2634, doi: 10.1093/mnras/stw379
  * Martinez et al. (2020) Martinez, M. A. S., Fragione, G., Kremer, K., et al. 2020, Astrophys. J., 903, 67, doi: 10.3847/1538-4357/abba25
  * McCuller et al. (2020) McCuller, L., et al. 2020, Phys. Rev. Lett., 124, 171102, doi: 10.1103/PhysRevLett.124.171102
  * McKernan et al. (2012) McKernan, B., Ford, K. E. S., Lyra, W., & Perets, H. B. 2012, MNRAS, 425, 460, doi: 10.1111/j.1365-2966.2012.21486.x
  * Ng et al. (2018) Ng, K. K. Y., Vitale, S., Zimmerman, A., et al. 2018, Phys. Rev. D, 98, 083007, doi: 10.1103/PhysRevD.98.083007
  * Nitz & Capano (2021) Nitz, A. H., & Capano, C. D. 2021, Astrophys. J. Lett., 907, L9, doi: 10.3847/2041-8213/abccc5
  * O’Reilly et al. (2018) O’Reilly, B., Branchesi, M., Haino, S., & Gemme, G. 2018.  https://dcc.ligo.org/LIGO-T2000012/public
  * Pratten et al. (2020a) Pratten, G., Schmidt, P., Buscicchio, R., & Thomas, L. M. 2020a, Phys. Rev. Res., 2, 043096, doi: 10.1103/PhysRevResearch.2.043096
  * Pratten et al. (2020b) Pratten, G., et al. 2020b.  https://arxiv.org/abs/2004.06503
  * Racine (2008) Racine, E. 2008, Phys. Rev. D, 78, 044021, doi: 10.1103/PhysRevD.78.044021
  * Reitze et al. (2019) Reitze, D., et al. 2019, Bull. Am. Astron. Soc., 51, 035.  https://arxiv.org/abs/1907.04833
  * Rodriguez et al. (2018) Rodriguez, C. L., Amaro-Seoane, P., Chatterjee, S., & Rasio, F. A. 2018, Phys. Rev. Lett., 120, 151101, doi: 10.1103/PhysRevLett.120.151101
  * Rodriguez et al. (2016a) Rodriguez, C. L., Chatterjee, S., & Rasio, F. A. 2016a, Phys. Rev. D, 93, 084029, doi: 10.1103/PhysRevD.93.084029
  * Rodriguez et al. (2016b) Rodriguez, C. L., Haster, C.-J., Chatterjee, S., Kalogera, V., & Rasio, F. A. 2016b, Astrophys. J. Lett., 824, L8, doi: 10.3847/2041-8205/824/1/L8
  * Rodriguez et al. (2019) Rodriguez, C. L., Zevin, M., Amaro-Seoane, P., et al. 2019, Phys. Rev. D, 100, 043027, doi: 10.1103/PhysRevD.100.043027
  * Romero-Shaw et al. (2019) Romero-Shaw, I. M., Lasky, P. D., & Thrane, E. 2019, Mon. Not. Roy. Astron. Soc., 490, 5210, doi: 10.1093/mnras/stz2996
  * Romero-Shaw et al. (2020a) Romero-Shaw, I. M., Lasky, P. D., Thrane, E., & Bustillo, J. C. 2020a.  https://arxiv.org/abs/2009.04771
  * Romero-Shaw et al. (2020b) Romero-Shaw, I. M., et al. 2020b, Mon. Not. Roy. Astron. Soc., 499, 3295, doi: 10.1093/mnras/staa2850
  * Schmidt et al. (2011) Schmidt, P., Hannam, M., Husa, S., & Ajith, P. 2011, Phys. Rev. D, 84, 024046, doi: 10.1103/PhysRevD.84.024046
  * Schmidt et al. (2015) Schmidt, P., Ohme, F., & Hannam, M. 2015, Phys. Rev. D, 91, 024043, doi: 10.1103/PhysRevD.91.024043
  * Sedda et al. (2020) Sedda, M. A., Mapelli, M., Spera, M., Benacquista, M., & Giacobbo, N. 2020, Astrophys. J., 894, 133, doi: 10.3847/1538-4357/ab88b2
  * Sigurdsson & Hernquist (1993) Sigurdsson, S., & Hernquist, L. 1993, Nature, 364, 423, doi: 10.1038/364423a0
  * Singer & Price (2016) Singer, L. P., & Price, L. R. 2016, Phys. Rev. D, 93, 024013, doi: 10.1103/PhysRevD.93.024013
  * Skilling (2004) Skilling, J. 2004, American Institute of Physics Conference Series, Vol. 735, Nested Sampling, ed. R. Fischer, R. Preuss, & U. V. Toussaint, 395–405, doi: 10.1063/1.1835238
  * Skilling (2006) Skilling, J. 2006, Bayesian Analysis, 1, 833 , doi: 10.1214/06-BA127
  * Speagle (2020) Speagle, J. S. 2020, MNRAS, 493, 3132, doi: 10.1093/mnras/staa278
  * Stevenson et al. (2017) Stevenson, S., Berry, C. P. L., & Mandel, I. 2017, Mon. Not. Roy. Astron. Soc., 471, 2801, doi: 10.1093/mnras/stx1764
  * Talbot & Thrane (2017) Talbot, C., & Thrane, E. 2017, Phys. Rev. D, 96, 023012, doi: 10.1103/PhysRevD.96.023012
  * Talbot & Thrane (2018) —. 2018, Astrophys. J., 856, 173, doi: 10.3847/1538-4357/aab34c
  * Thomas et al. (2020) Thomas, L. M., Schmidt, P., & Pratten, G. 2020.  https://arxiv.org/abs/2012.02209
  * Thrane & Talbot (2019) Thrane, E., & Talbot, C. 2019, PASA, 36, e010, doi: 10.1017/pasa.2019.2
  * Tiwari et al. (2018) Tiwari, V., Fairhurst, S., & Hannam, M. 2018, Astrophys. J., 868, 140, doi: 10.3847/1538-4357/aae8df
  * Tse et al. (2019) Tse, M., et al. 2019, Phys. Rev. Lett., 123, 231107, doi: 10.1103/PhysRevLett.123.231107
  * Veitch & Vecchio (2008) Veitch, J., & Vecchio, A. 2008, Phys. Rev. D, 78, 022001, doi: 10.1103/PhysRevD.78.022001
  * Veitch et al. (2015) Veitch, J., et al. 2015, Phys. Rev. D, 91, 042003, doi: 10.1103/PhysRevD.91.042003
  * Vitale & Evans (2017) Vitale, S., & Evans, M. 2017, Phys. Rev. D, 95, 064052, doi: 10.1103/PhysRevD.95.064052
  * Vitale et al. (2017a) Vitale, S., Gerosa, D., Haster, C.-J., Chatziioannou, K., & Zimmerman, A. 2017a, Phys. Rev. Lett., 119, 251103, doi: 10.1103/PhysRevLett.119.251103
  * Vitale et al. (2017b) Vitale, S., Lynch, R., Raymond, V., et al. 2017b, Phys. Rev. D, 95, 064053, doi: 10.1103/PhysRevD.95.064053
  * Vitale et al. (2017c) Vitale, S., Lynch, R., Sturani, R., & Graff, P. 2017c, Class. Quant. Grav., 34, 03LT01, doi: 10.1088/1361-6382/aa552e
  * Vitale et al. (2014) Vitale, S., Lynch, R., Veitch, J., Raymond, V., & Sturani, R. 2014, Phys. Rev. Lett., 112, 251101, doi: 10.1103/PhysRevLett.112.251101
  * Voss & Tauris (2003) Voss, R., & Tauris, T. M. 2003, Monthly Notices of the Royal Astronomical Society, 342, 1169, doi: 10.1046/j.1365-8711.2003.06616.x
  * Woosley (2017) Woosley, S. E. 2017, The Astrophysical Journal, 836, 244, doi: 10.3847/1538-4357/836/2/244
  * Wysocki et al. (2019) Wysocki, D., Lange, J., & O’Shaughnessy, R. 2019, Phys. Rev. D, 100, 043012, doi: 10.1103/PhysRevD.100.043012
  * Zevin et al. (2021) Zevin, M., Romero-Shaw, I. M., Kremer, K., Thrane, E., & Lasky, P. D. 2021.  https://arxiv.org/abs/2106.09042
  * Zwart & McMillan (2002) Zwart, S. F. P., & McMillan, S. L. W. 2002, The Astrophysical Journal, 576, 899, doi: 10.1086/341798