Two-harmonic approximation for gravitational waveforms from precessing binaries

TL;DR

This paper introduces a two-harmonic approximation method for modeling gravitational waveforms from precessing binary black holes, simplifying analysis and aiding detection of precession effects in gravitational-wave signals.

Contribution

It proposes a novel two-harmonic approximation for precessing binary waveforms, enabling easier identification of precession effects and improving search strategies.

Findings

Waveforms can be approximated by two dominant harmonics.

Precession effects are identifiable through harmonic beating.

The method simplifies detection and characterization of precessing binaries.

Abstract

Binary-black-hole orbits precess when the black-hole spins are mis-aligned with the binary's orbital angular momentum. The apparently complicated dynamics can in most cases be described as simple precession of the orbital angular momentum about an approximately fixed total angular momentum. However, the imprint of the precession on the observed gravitational-wave signal is yet more complicated, with a non-trivial time-varying dependence on black-hole dynamics, the binary's orientation and the detector polarization. As a result, it is difficult to predict under which conditions precession effects are measurable in gravitational-wave observations, and their impact on both signal detection and source characterization. We show that the observed waveform can be simplified by decomposing it as a power series in a new precession parameter $b = \tan(\beta/2)$, where $\beta$ is the opening angle between the orbital and total angular momenta. The power series is made up of five harmonics, with frequencies that differ by the binary's precession frequency, and individually do not exhibit amplitude and phase modulations. In many cases, the waveform can be well approximated by the two leading harmonics. In this approximation we are able to obtain a simple picture of precession as caused by the beating of two waveforms of similar frequency. This enables us to identify regions of the parameter space where precession is likely to have an observable effect on the waveform, and to propose a new approach to searching for signals from precessing binaries, based upon the two-harmonic approximation.

Equations

$${\Omega }_p=\dot{\alpha }.$$

$$\tan \beta =\frac{S_{\perp }}{\mu \sqrt{Mr}+S_{\parallel }},$$

$${\Omega }_p\approx \left(2+\frac{3}{2q}\right)\frac{J}{r^3},$$

$${\chi }_{\text{eff}}=\frac{1}{M}\left(\frac{{\mathbf{S}}_1}{m_1}+\frac{{\mathbf{S}}_2}{m_2}\right)\cdot \hat{\mathbf{L}},$$

$${\chi }_{\text{p}}=\frac{1}{A_1{m}_1^2}\text{max}\left(A_1{S}_{1\perp },A_2{S}_{2\perp }\right),$$

$$h(t)$$

$${\Phi }_S(t)={\varphi }_{orb}(t)-ϵ(t)$$

$$\dot{ϵ}(t):=\dot{\alpha }(t)\cos \beta (t)$$

$$C_{+}$$

$$S_{+}$$

$$C_{\times }$$

$$S_{\times }$$

$$b=\tan (\beta /2),$$

$$\frac{1+{\cos }^2\beta }{2}$$

$$\cos \beta$$

$$\frac{\sin 2\beta }{2}$$

$$\sin \beta$$

$${\sin }^2\beta$$

$$\left(\frac{d_o}{d_L}\right)\left(C_{+}-i{S}_{+}\right)$$

$$\left(\frac{d_o}{d_L}\right)\left(C_{\times }-i{S}_{\times }\right)$$

$${\mathcal{A}}_0^{+}={\mathcal{A}}_4^{+}$$

$${\mathcal{A}}_0^{\times }=-{\mathcal{A}}_4^{\times }$$

$${\mathcal{A}}_1^{+}=-{\mathcal{A}}_3^{+}$$

$${\mathcal{A}}_1^{\times }={\mathcal{A}}_3^{\times }$$

$${\mathcal{A}}_2^{+}$$

$${\mathcal{A}}_2^{\times }$$

$$h(t)$$

$$F_{+}$$

$$F_{\times }$$

$$\Phi (t)$$

$$h=\sum _{k=0}^4{w}_{+}(h_0^k{\mathcal{A}}_k^1+h_{\frac{\pi }{2}}^k{\mathcal{A}}_k^3)+w_{\times }(h_0^k{\mathcal{A}}_k^2+h_{\frac{\pi }{2}}^k{\mathcal{A}}_k^4),$$

$$h_0^k(t)$$

$$h_{\frac{\pi }{2}}^k(t)$$

$${\mathcal{A}}_k^1$$

$${\mathcal{A}}_k^2$$

$${\mathcal{A}}_k^3$$

$${\mathcal{A}}_k^4$$

$${\varphi }_k=2{\varphi }_o+(2-k){\alpha }_o.$$

$$h_{{\varphi }_o=0;\psi =0}$$

Full text

I Introduction

When the spins of black holes in a binary system are mis-aligned with the binary’s orbital angular momentum, both the spins and orbital angular momentum will precess Apostolatos:1994mx ; Kidder:1995zr ; Apostolatos:1995pj ; Buonanno:2002fy . We therefore expect that most astrophysical binaries will undergo precession, but to date there has been no evidence of precession in gravitational-wave (GW) observations from the Advanced LIGO and Virgo detectors LIGOScientific:2018mvr ; LIGOScientific:2018jsj . This is not necessarily surprising, because precession often leaves only a weak imprint on the observable signal, particularly when the black holes are of comparable mass and the binary’s orbit is face-on to the detector, which are the most likely configurations that have been observed so far. Despite this heuristic picture, there is no simple means to estimate the measurability of precession of a given binary configuration, and as such it is difficult to predict when precession effects will be conclusively observed in GW events.

Detailed parameter estimation techniques have been developed, which enable the reconstruction of the parameters of observed signals vanderSluys:2008qx ; Veitch:2009hd ; Veitch:2014wba ; Biwer:2018osg ; Smith:2016qas , in addition to approximate Fisher-matrix methods Poisson:1995ef ; Lang:1900bz . In parallel, techniques have been developed that provide an intuitive understanding of the measurement accuracy of certain parameters (or parameter combinations)  Fairhurst:2009tc ; Baird:2012cu ; Hannam:2013uu ; Singer:2015ema ; Vallisneri:2012qq ; Usman:2018imj . These have typically involved either approximations (such as leading order, Fisher Matrix type calculations), restriction to a subset of system parameters (for example masses and spins; timing and sky location; binary orientation). Combined, these give an understanding of the accuracy of parameter estimation for non-precessing systems.

In parallel, there have been significant developments in understanding the implications of precession, starting with the early work in Refs. Apostolatos:1994mx ; Apostolatos:1995pj ; Buonanno:2002fy which provided insights into the impact of precession on the gravitational waveform emitted during the inspiral of compact binaries. Subsequently, black hole binary waveforms which incorporate precession through merger have been developed Hannam:2013oca ; Taracchini:2013rva ; Khan:2019kot ; Varma:2019csw ; Pratten:2020ceb ; large scale parameter estimation studies of precession have been performed to identify the regions of parameter space where precession will be observable Vitale:2014mka ; Trifiro:2015zda ; Vitale:2016avz ; OShaughnessy:2014shr ; Farr:2014qka ; Littenberg:2016dgj ; Farr:2015lna ; Vitale:2016icu ; and new theoretical insights into the impact of precession on both detection and parameter estimation have been obtained Brown:2012gs ; o2020semianalytic ; Lundgren:2013jla . Complementary to this, there have been several efforts to understand the impact of precession on searches Harry:2013tca ; Brown:2012gs , and to implement searches for precessing signals Apostolatos:1995pj ; Fazi:2009ifa ; Harry:2011qh ; Harry:2016ijz ; Buonanno:2002fy . This has led to an increasingly clear picture of the impact of precession: it is most significant for binaries with large mass ratios, where the in-plane spin components are large and for systems where the total angular momentum is mis-aligned with the line of sight.

At leading order, the gravitational waveform emitted by a precessing binary is composed of five harmonics, which are offset by multiples of the precession frequency Hannam:2013oca ; Lundgren:2013jla . We show that these harmonics form a natural hierarchy with the amplitude of the sub-leading harmonics suppressed by a factor that depends upon the opening angle (the angle between the orbital and total angular momenta). Using this approximation, and restricting to the two leading harmonics, we are able to obtain relatively simple expressions for the precession waveform. Each harmonic takes the form of a non-precessing-binary waveform (i.e., with monotonic amplitude and frequency evolution during the inspiral of non-eccentric systems), and the amplitude and phase modulations of the complete precessing-binary waveform arise as beating between the two harmonics.

The purpose of this paper is to introduce this decomposition (Sec. III), with an alternative derivation given in the Appendix, and the two-harmonic approximation (Sec. IV), and to identify its range of validity and accuracy (Sec. V). We then discuss a proposed search for precessing binaries using the two-harmonic approximation (Sec. VI) and finally introduce the notion of a “precession SNR” that can be used to determine whether precession effects are observable in a given system (Sec. VII). We begin in the next section with a summary of precession in black-hole binaries.

II Black hole Spin Induced Precession

In the general theory of relativity a binary consisting of two objects of masses, $m_{1}$ and $m_{2}$ (where we choose $m_{1}\geq m_{2}$ and denote $q=m_{1}/m_{2}$, so that $q\geq 1$), with spin angular momenta $\mathbf{S_{1}}$ and $\mathbf{S_{2}}$, orbiting each other with angular momentum $\mathbf{L}$, will slowly inspiral due to the loss of energy and momentum through the emission of gravitational waves. If $\mathbf{S}_{1}\parallel\mathbf{S}_{2}\parallel\mathbf{L}$, then the plane of the orbit remains fixed and in non-eccentric binaries the amplitude and frequency of the emitted gravitational wave increases as the orbital separation decreases. The system eventually merges and forms a single perturbed black hole that emits gravitational radiation as a superposition of quasinormal ringdown multipoles, until the system settles down to its final state Kokkotas1999 .

For the case where the total spin is not aligned with the total orbital angular momentum, $\left(\mathbf{S}_{1}+\mathbf{S}_{2}\right)\times\mathbf{L}\neq 0$, in most cases the orbital plane of the binary will precess around the approximately constant total angular momentum $\mathbf{J}=\mathbf{S_{1}}+\mathbf{S_{2}}+\mathbf{L}$, i.e., $\mathbf{L}$ precesses around $\mathbf{J}$, and the spins precess such that $\dot{\mathbf{S}}=-\dot{\mathbf{L}}$ Apostolatos:1994mx . For configurations where $\mathbf{J}\approx 0$, the system undergoes “transitional precession” Apostolatos:1994mx ; Kidder:1995zr , but this is expected to be rare in LIGO-Virgo detections. The angle between $\mathbf{L}$ and $\mathbf{J}$ is denoted by $\beta$. In simple precession cases $\beta$ slowly increases during inspiral as $L$ decreases (recall that in the Newtonian limit $L\propto\sqrt{r}$, where $r$ is the orbital separation), but the spin magnitudes $S_{1}$ and $S_{2}$ remain fixed, and, to a good approximation, so do their orbit-averaged components parallel and perpendicular to $L$, $S_{i||}$ and $S_{i\perp}$. The opening angle $\beta$ typically varies very little over the portion of a binary’s inspiral that is visible in a GW detector, and so it is often possible to make the approximation that $\beta$ is constant. This approximation has been used to good effect in Ref. Brown:2012gs , and we will also use it in some of the discussion in this paper.

Adopting the notation that the inclination angle of the binary as seen by an observer, $\iota$, is the angle between the orbital angular momentum and the line of sight (see Fig.1), $\cos\iota=\mathbf{\hat{L}}\cdot\mathbf{\hat{N}}$, where a caret denotes a unit vector (e.g. $\hat{\mathbf{a}}=\mathbf{a}/|\mathbf{a}|$), the binary’s orbital inclination becomes times dependent. As a result the energy emitted in GWs in the $\mathbf{\hat{N}}$ direction will also be time dependent, where the maximum instantaneous energy emission is approximately in the direction of $\mathbf{\hat{L}}$. If $\mathbf{\hat{N}}$ is aligned with $\mathbf{\hat{J}}$, then $\iota\approx\beta$ and varies slowly and with minimal oscillations due only to orbital nutation. If $\mathbf{\hat{N}}$ is in some other direction, then the energy emission will be modulated on the precession timescale. In the following we will not use the inclination $\iota$, but rather combinations of $\beta$ and the angle between $\mathbf{J}$ and $\mathbf{\hat{N}}$, denoted by $\theta_{\rm JN}$. As noted previously, $\mathbf{\hat{J}}$ is approximately constant for simple precession cases, and we will treat it as a constant in the analysis in Sec. III.

The signal measured in a detector will exhibit modulations in phase and amplitude that depend on $\beta$, $\theta_{\rm JN}$, the precession angle of $\mathbf{L}$ around $\mathbf{J}$, denoted by $\alpha$, and the polarisation $\psi$ of the observed signal. These angles are illustrated in Fig. 1, and discussed further in Sec. III. For now we note several well-known features of precession waveforms Apostolatos:1994mx ; Kidder:1995zr , which will be further sharpened in the discussion later in the paper. The strength of precession in a system is characterised by the degree of tilt of the binary’s orbit, given by $\beta$, and by the precession frequency $\Omega_{P}$ of $\mathbf{L}$ around $\mathbf{J}$, which is given by

[TABLE]

The angle $\beta$ is determined primarily by the total spin in the plane, and binary’s mass ratio and separation. At leading order we can write the orbital angular momentum of the system as $L=\mu\sqrt{Mr}$, where $\mu$ is the reduced mass, $\mu=m_{1}m_{2}/M=qM/(1+q)^{2}$, and so to first approximation,

[TABLE]

which provides us with the basic dependence of $\beta$ on the binary configuration. At leading order the precession frequency can be written as,

[TABLE]

meaning that to first approximation it does not depend on the spins (or therefore the opening angle $\beta$), but only on the binary’s total mass, mass-ratio, and separation (or equivalently orbital frequency). The number of precession cycles over a certain time or frequency range (e.g., over the course of an observation), depends on the total mass and mass-ratio of the binary. In a GW observation there is a partial degeneracy between the mass ratio and the aligned spin $S_{||}$ Cutler:1994ys ; Poisson:1995ef ; Baird:2012cu , meaning that one of the chief effects of a measurement of precession will be to improve the measurement of these two physical properties OShaughnessy:2014shr .

In the remainder of this paper we choose to describe the gravitational wave signal, precessing or non-precessing, with the IMRPhenomPv2 phenomenological model presented in Ref. Hannam:2013oca . This model exploits the phenomenology of simply precessing binaries described earlier, with the additional approximation that a precessing-binary waveform can be factorised into an underlying non-precessing waveform, and the precessional dynamics Schmidt:2012rh . The underlying non-precessing-binary model is IMRPhenomD Husa:2015iqa ; Khan:2015jqa , using only the spin components aligned with $\mathbf{L}$. In IMRPhenomD both aligned spin components are used to generate an approximate post-Newtonian phasing and amplitude, with corrections provided by fits to numerical-relativity waveforms, that are parameterised by two different combinations of the two spin components. Although IMRPhenomD has been found to model well two-spin systems Kumar:2016dhh , its dominant spin dependence can be characterised well by the effective spin,

[TABLE]

which takes values between $-1$ (both maximal anti-aligned spins) and $+1$ (both maximal aligned spins) to describe the magnitude of spin aligned with the total angular momentum. For a given configuration IMRPhenomPv2 uses the corresponding IMRPhenomD waveform, but with the final spin modified to take into account the in-plane spin components. A frequency-dependent rotation is then applied to the non-precessing waveform to introduce the precession dynamics, which are modelled by frequency-domain post-Newonian expressions for the precession angles for an approximately equivalent single-spin system Schmidt:2014iyl ; Hannam:2013oca , where the large black hole has spin,

[TABLE]

where $A_{1}=2+3q/2$ and $A_{2}=2+3/(2q)$ and $\mathbf{S}_{i\perp}$ is the component of the spin perpendicular to $\mathbf{L}$. The effective precession spin parameter is obtained by averaging the relative in-plane spin orientation over a precession cycle, and so more accurate for a system that undergoes many precession cycles.

There are several important features which are not incorporated in the IMRPhenomPv2 waveform. These include two-spin effects Taracchini:2013rva ; Babak:2016tgq ; Khan:2018fmp , gravitational wave multipoles other than the leading 22 mode Khan:2019kot , significant precession during merger OShaughnessy:2012iol , and spin alignment due to spin-orbit resonances during inspiral Gerosa:2014kta ; Gerosa:2015tea . Some of these effects will have an impact upon the distributions of black hole spin orientations when the binaries enter the LIGO or Virgo sensitivity band while others can leave imprints on the waveform which may be observable, particularly close to the merger. Nonetheless, the IMRPhenomPv2 has been used in the analysis of all LIGO-Virgo observations during the first two observing runs Abbott:2016blz ; TheLIGOScientific:2016pea ; LIGOScientific:2018jsj ; LIGOScientific:2018mvr , and it captures much of the dominant phenomenology of precessing-binary waveforms. In addition, the decomposition presented in the next section is in no way tied to the particular waveform used and could be equally well applied to other waveform models for precessing binaries which, for example, incorporate two-spin effects and precession during merger. The current formalism does not include additional gravitational wave multipoles, and we will investigate this in a future work. We expect the broad features of many of the results presented in the remainder of the paper to be relatively unaffected by the specific waveform choice, but the details for any specific signal could change.

III Harmonic decomposition of the waveform from a precessing binary

The gravitational waveform emitted by a precessing system, as observed at a gravitational wave detector, can be expressed approximately as Buonanno:2002fy ; Brown:2012gs

[TABLE]

Here, $A_{o}(t)$ denotes the amplitude of the gravitational wave signal in a (time-varying) frame aligned with the orbital angular momentum of the binary, and depends upon the masses and spins of the binary. Since the amplitude scales linearly with the luminosity distance, we have chosen to introduce a fiducial normalization $A_{o}(t)$ for a waveform at a distance $d_{o}$ and explicitly extract the distance dependence.111Of course, the observed waveform is also affected by the redshifting of frequencies. For the calculation discussed here, we work in the detector frame and consider the observed masses, which are $(1+z)$ times the source frame masses. $\Phi_{S}$ is the phase evolution in the frame aligned with the total angular momentum $\mathbf{J}$ of the binary. The phase evolution, $\Phi_{S}$, is related to the orbital phase, $\phi_{orb}$, as

[TABLE]

where Boyle:2011gg

[TABLE]

and, as before, $\beta$ is the opening angle and $\alpha$ gives the phase of the precession of $\mathbf{L}$ around $\mathbf{J}$ as shown in Fig. 1. $F_{+}$ and $F_{\times}$ give the detector response relative to the $\mathbf{J}$-aligned frame and $C_{+,\times}$, $S_{+,\times}$ encode the time-varying response to the gravitational wave due to the evolution of the binary’s orbit relative to the detector. They depend upon the three angles introduced previously: the precession opening angle $\beta$ and phase $\alpha$ and the angle between the total orbital angular momentum and the line of sight $\theta_{\rm JN}$. In terms of these angles, we can express $C_{+,\times}$ and $S_{+,\times}$ as222We have re-written the $C_{+}$ term relative to what is normally given in the literature, e.g. Buonanno:2002fy ; Brown:2012gs , to group terms with the same $\alpha$ dependence.

[TABLE]

The non-precessing expressions can be recovered in the limit of $\beta\rightarrow 0$ and $\alpha\rightarrow\mathrm{constant}$ (which is then degenerate with the polarization of the system). When $\beta$ is non-zero, the effect of precession is to modulate the detector response at frequencies $\Omega_{P}$ and $2\Omega_{P}$. To make the harmonic content of $C_{+,\times}$ and $S_{+,\times}$ more explicit, we first introduce the parameter,

[TABLE]

and write the response functions in terms of it. The terms involving $\beta$ can be expressed as

[TABLE]

Substituting the trigonometric identities from Eq. (11) into the expressions for $C_{+}$ and $S_{+}$ in Eq. (9) we obtain,

[TABLE]

where we have introduced $\mathcal{A}_{k}^{+}$ and $\mathcal{A}_{k}^{\times}$ as

[TABLE]

In the approximation where the direction of total angular momentum is constant, the $\mathcal{A}^{+,\times}_{k}$ are time independent amplitudes, and the time dependence of the amplitude functions is captured as a power series in the parameter $b=\tan(\beta/2)$.

Finally, we can use the harmonic decomposition in Eq. (12) to obtain a decomposition of the waveform, Eq. (6),

[TABLE]

This allows us to clearly identify the impact of precession on the waveform. First, precession leads to an additional phase evolution at frequency $2\Omega_{P}$ and a decrease in the amplitude by a factor $(1+b^{2})^{2}$. The precessing waveform contains five harmonics that form a power series in $b$, whose amplitude depends upon the detector response, distance and viewing angle of the binary. The frequency of each harmonic is offset from the next by the precession frequency $\Omega_{P}$. Similar results have been obtained previously, by manipulating the spin-weighted spherical harmonic decomposition of the waveform, e.g. o2020semianalytic ; Lundgren:2013jla . However, it was not previously observed that the relative amplitudes of the harmonics were related in a straightforward manner. In the Appendix, we present an alternative derivation of the result in Eq. (14) in terms of this spin-weighted spherical harmonic decomposition of the waveform, as is customary when producing waveform models for precessing binaries Hannam:2013oca .

As a final step, we would like to explicitly extract three more time-independent angles that characterize the waveform, namely the polarization angle $\psi$, the initial phase $\phi_{o}$ and the initial polarization phase $\alpha_{o}$.333The initial polarization phase $\alpha_{o}$ is sometimes denoted in the literature as $\phi_{JL}$.

The unknown polarization $\psi$ is currently folded into the detector response functions $F_{+,\times}$. It is more useful to extract $\psi$ and then consider the detector response to be a known quantity dependent upon only the details of the detector and the direction to the source. Thus, we write the detector response as,

[TABLE]

where $w_{+}$ and $w_{\times}$ are the detector response functions in a fixed frame — for a single detector it is natural to choose $w_{\times}=0$ and for a network to work in the dominant polarization for which $w_{+}$ is maximized Harry:2010fr . The unknown polarization of the source relative to this preferred frame is denoted $\psi$.

To isolate the initial orbital and precession phases, we explicitly extract them from the binary’s phase evolution by introducing,

[TABLE]

Thus $\Phi(t)$ vanishes at $t=0$ and evolves as the sum of the orbital phase and an additional, precession dependent, contribution.

We then substitute the expressions for $F_{+,\times}$, Eq. (15), and $\Phi$, Eq. (16), into the expression for $h(t)$ given in Eq. (14), and isolate the time-varying terms from the constant, orientation dependent angles. The waveform can be written as the sum of five precessing harmonics, the amplitudes of which are constants that depend upon the binary’s sky location, distance and orientation:

[TABLE]

where $h^{k}_{0,\frac{\pi}{2}}$ are the waveform harmonics and $\mathcal{A}_{k}^{\mu}$ are constants. The waveform harmonics are

[TABLE]

The harmonics form a simple power series in $be^{-i\alpha}$, so the amplitude of each successive harmonic is reduced by a factor of $b$, and the frequency is reduced by $\Omega_{P}$.

The amplitudes for the harmonics are given by

[TABLE]

where the $\mathcal{A}^{+,\times}$ were introduced in Eq. (13), $\psi$ is the polarization and the phase angle for each harmonic is,

[TABLE]

These amplitudes form a generalization of the $\mathcal{F}$-statistic decomposition of the non-precessing binary waveform (see e.g. Harry:2010fr ). In the limit that $b\rightarrow 0$, the precessing decomposition reduces to the standard expression for the non-precessing waveform as the amplitude for all harmonics other than $k=0$ vanish.

The precessing waveform can equally well be written in the frequency domain by performing a Fourier transform of the time-domain expressions given above Droz:1999qx . In this case, Eq. (17) is unchanged, as are the constant amplitude terms in Eq. (19). The frequency dependent harmonics are simply the Fourier transform of the time-domain modes given in Eq. (18), and naturally satisfy $h^{k}_{\frac{\pi}{2}}=ih^{k}_{0}$.

The expansion above is most natural when $b<1$, which corresponds to opening angles of $\beta<90^{\circ}$. In cases where the opening angle is greater than $90^{\circ}$ it is natural to re-express the waveform in terms of $c=b^{-1}=\cot(\beta/2)$ in which case the waveform can be expressed as a power series in $c$. We will not discuss the large opening angle calculation further in this paper, but note that many of the arguments presented below would extend in a straightforward manner to this case.

a Obtaining the harmonics

Here, we give an explicit prescription to obtain the five harmonics for the waveform, introduced in Eq. (17). To do so, we generate waveforms for orientations that contain only a subset of the harmonics, and combine them to isolate a single harmonic. For simplicity, we restrict attention to the $+$ polarization by fixing $w_{+}=1$, $w_{\times}=0$ and consider a binary at a distance $d_{L}=d_{o}$.

Harmonics $k=0$ and $k=4$. When the viewing angle of the signal is aligned with the total angular momentum, $\theta_{\rm JN}=0$, the observed waveform contains only the zeroth and fourth harmonics as $\mathcal{A}^{+,\times}_{1,2,3}$ vanish for $\theta_{\rm JN}=0$. We also fix $\alpha_{o}=0$, to obtain,

[TABLE]

From these, we can extract the $k=0$ and 4 harmonics,

[TABLE]

The $\tfrac{\pi}{2}$ phases of the harmonics can be obtained in an identical way.

Harmonics $k=1$ and $k=3$. When the signal is edge on, the $\times$ polarization contains only the first and third harmonics. Then, fixing $\theta_{\rm JN}=\tfrac{\pi}{2}$ and $\psi=\tfrac{\pi}{4}$, we have,

[TABLE]

so that,

[TABLE]

Harmonic $k=2$. Finally, from the $+$ polarization of the edge-on waveform, we can extract the second harmonic — in principle we could also get $k=0$ and $k=4$, but we have already described a method to obtain them. Fixing $\theta_{\rm JN}=\tfrac{\pi}{2}$ and $\psi=0$ we have,

[TABLE]

so that,

[TABLE]

b Precession with varying orientation

The observable effect of precession will vary significantly with the binary orientation, as has been discussed in many previous works, for example Apostolatos:1994mx ; Brown:2012gs . Interestingly, both the amplitude and frequency of the observed precession depends upon the viewing angle. The harmonic decomposition derived above provides a straightforward way to understand this effect. The observed amplitude and phase modulations can be understood as the beating of the different harmonics against each other, with the amplitude of the composite waveform being maximum when the harmonics are in phase and minimum when they are out of phase.

In Fig. 2, we show the waveform for four different orientations: a) along $\mathbf{J}$, b) $\times$ polarization at $45^{\circ}$ to $\mathbf{J}$, c/d) $+/\times$ polarization orthogonal to $\mathbf{J}$. In all cases, we show the last two seconds of the waveform (from around 25 Hz) for a $40M_{\odot}$ binary, with $q=6$, and in-plane spin on the larger black hole of $\chi_{P}=0.6$. This configuration gives an opening angle of $\beta\approx 45^{\circ}$ (and $b\approx 0.4$) which leads to significant precession effects in the waveform.

When viewed along $\mathbf{J}$, there is minimal precession as only the $k=0$ and $4$ harmonics are present in the system and the $k=4$ harmonic is down-weighted by a factor of $b^{4}\approx 0.03$ relative to the leading harmonic. Furthermore, the modulation comes from the beating of the $k=0$ and $k=4$ harmonics and occurs at four times the precession frequency. When the line of sight is orthogonal to the total angular momentum, the $k=0,2,4$ harmonics are present in the $+$ polarized waveform and $k=1,3$ in the $\times$ polarization. The $k=0$ and 2 harmonics have close to equal amplitude (although $k=2$ is down-weighted by $b^{2}\approx 0.17$, the amplitude as given in Eq. (13) is maximal). Consequently the observed waveform has maximal amplitude and phase modulation due to precession. For the $\times$ polarized signal, it is the $k=1,3$ harmonics that contribute, with $k=3$ a factor of $b^{2}\approx 0.17$ smaller than $k=1$. Consequently, precession effects are less significant. In both cases, precession occurs at twice the precession frequency as it is from the beating of $k=0$ and $k=2$ ($+$ polarization) or $k=1$ and $k=3$ ($\times$ polarization). For the $\times$ polarized signal with $\theta_{\rm JN}=45^{\circ}$, the $k=0,1,3,4$ harmonics are present, with $k=0,1$ dominating and having approximately equal amplitude. For this signal, the binary precesses from a face-on orientation, $\iota=0$ to edge-on, $\iota=90^{\circ}$, and the waveform amplitude oscillates from the maximum to zero. Here, modulations occur at the precession frequency.

c Importance of precession over parameter space

From the intuitive discussion of precession presented in Apostolatos:1994mx ; Buonanno:2002fy ; Brown:2012gs and summarized in Section II, it is straightforward to identify regions of parameter space where precession is most likely to have a significant impact upon the binary dynamics and, consequently, the observed waveform. Specifically, we expect that higher mass ratios, larger in-plane spins and negative aligned spin components will all lead to a larger opening angle and more significant precession Brown:2012gs . Here we briefly revisit this discussion, framing our results in terms of the parameter $b$ introduced earlier. Explicitly, we introduce the waveform-averaged value of $b$ as,

[TABLE]

where $h^{0,1}$ are the harmonics of the waveform introduced in Eq. (18) and $S_{n}(f)$ is the noise power spectrum of the detector. For this work, we choose $S_{n}(f)$ to be the design-sensitivity Advanced LIGO noise curve LIGOScientific:2018mvr and evaluate the integral over the frequency range $f\in[20,1024]$ Hz 444Using a realistic noise curve similar to the observed curves during 01 and O2 would change the reported values slightly, as these noise curves are less sensitive than design, particularly at low frequencies. The qualitative patterns seen in the figure would remain the same however. For binaries where the opening angle $\beta$ is approximately constant, $\overline{b}\approx\tan(\overline{\beta}/2)$.

Fig. 3 shows the value of $\overline{b}$ on several two-dimensional slices through the four dimensional parameter space of total mass $M$, mass ratio $q$, effective spin $\chi_{\mathrm{eff}}$ and precessing spin $\chi_{p}$. Keeping other quantities fixed, the value of $\overline{b}$ increases with total mass. For higher masses, the late inspiral and merger occur in the sensitive band of the detectors and, close to merger, the opening angle increases as orbital angular momentum is radiated. For masses above $40M_{\odot}$ the mass dependence of $\overline{b}$ is small, with only a $10\%$ decrease from $40M_{\odot}$ to $100M_{\odot}$. Thus, for the other figures, we fix $M=40M_{\odot}$ and investigate the dependence of $\overline{b}$ on $q,\chi_{\mathrm{eff}}$ and $\chi_{p}$. The dependence of $\overline{b}$ follows directly from Eq. (2). The opening angle will increase with mass ratio, as the orbital angular momentum decreases. The opening angle, and also $\overline{b}$, increase with $\chi_{p}$. It follows directly from the definition that $\tan{\beta}$ scales linearly with $\chi_{p}$, and hence approximately linearly for $b=\tan(\beta/2)$. Finally, the opening angle decreases as the effective spin $\chi_{\mathrm{eff}}$ increases, so that the largest value of $\overline{b}$ is obtained with significant spin anti-aligned with $\mathbf{J}$.

Over much of the parameter space we have explored, $\overline{b}\lesssim 0.3$. This includes binaries with mass ratio up to 4:1, with precessing spin $\chi_{p}\lesssim 0.6$, and zero or positive aligned spin, $\chi_{\mathrm{eff}}\geq 0$. Only a small part of parameter space has $\overline{b}>0.4$, the value used in generating the waveforms in Figure 2, and $b>0.5$ is only achieved with at least two of: a) close to maximal $\chi_{p}$, b) high mass ratio, $q\gtrsim 5$ or c) significant spin anti-aligned with the orbital angular momentum $\chi_{\mathrm{eff}}\lesssim-0.4$.

Next, we consider the importance of precession for an astrophysically motivated population. In Fig. 4, we show the distribution of $\overline{b}$ for three distributions of black hole masses and spins. For each population, we generate 100,000 binaries uniformly in co-moving distance, with masses drawn from a power law distribution — $p(m_{1})\propto m_{1}^{-\alpha}$, with $\alpha=2.35$ — and different spin distributions, which are the same as those used in Refs. Farr:2017uvj ; Tiwari:2018qch ; fairhurst2019will . We consider populations where the spins are preferentially low and aligned with the binary orbit; low and isotropically aligned or drawn from a flat distribution and preferentially leading to precession. A low spin distribution is a triangular distribution peaked at zero spin and dropping to zero at maximal spin while a flat distribution is a uniform between zero and maximal spin. The aligned distribution is strongly peaked towards aligned spins, while the isotropic distribution assumes a uniform distribution of spin orientations over the sphere. The precessing distribution is strongly peaked towards spins orthogonal to the orbital angular momentum, i.e., with significant orbital precession Antonini:2017tgo ; Rodriguez:2018jqu . To account for observational biases, we keep only those signals that would be observable above a fixed threshold in a gravitational wave detector. We find that even for the most extreme precessing population considered, the mean value of $\overline{b}$ is 0.15 with over 90% of binaries having $\overline{b}<0.3$. This result is obviously sensitive to the assumptions on the mass and spin distribution. In Ref. fairhurst2019will we investigate a larger number of spin distributions, including ones which allow for large spin magnitudes, and we find that the peak of the $\overline{b}$ distribution is below $0.2$ and that over 90% of binaries have $\overline{b}<0.4$ in all cases.

Fig. 5 shows $\overline{b}$ for a range of neutron star–neutron star and neutron star–black hole binaries. For neutron star–black hole binaries, the picture is similar to that for black hole binaries, with large values of $\overline{b}$ observed for high mass ratios and large $\chi_{P}$. However, as an earlier part of the waveform is in the detector’s sensitive band, the impact of precession is less observable at fixed mass ratio than for higher mass black hole binaries. For neutron star binaries, the value of $\overline{b}$ remains below $0.15$ across the parameter space, and is less than $0.05$ for reasonable neutron star spins, $\chi\lesssim 0.4$.

IV The two-harmonic approximation

The precessing waveform can be expressed as the sum of five harmonics whose amplitudes form a power series in $b=\tan(\beta/2)$. Furthermore, over the majority of the space of binary mergers, the value of $b$ is less than $0.3$. In addition, for $\overline{b}\leq 0.4$ the dominant harmonic — the one containing the most power — must be either $k=0$ or $1$. Thus, for the vast majority of binary mergers, we expect that these two harmonics will be the most significant.

This motivates us to introduce the two-harmonic approximation, in which we generate a waveform containing only the $k=0$ and $k=1$ harmonics, i.e.,

[TABLE]

The expression for the two-harmonic waveform can be simplified by restricting to the single detector case (i.e., setting $w_{+}=1$ and $w_{\times}=0$), explicitly working with the waveform in the frequency domain, for which $h^{k}_{\tfrac{\pi}{2}}(f)=ih^{k}_{0}(f)$, and dropping the subscript [math] on the zero-phase waveform, so that $h^{k}(f):=h^{k}_{0}(f)$. The two harmonics of interest are,

[TABLE]

and the two-harmonic waveform then becomes,

[TABLE]

where,

[TABLE]

Thus, the two-harmonic waveform is composed of two components that have frequencies offset by $\Omega_{P}$, and any observed amplitude and phase modulation of the waveform is caused by the beating of one waveform against the other. The relative amplitude and phase of the two harmonics is encoded by

[TABLE]

The value of $\zeta$ depends upon the viewing angle, encoded in $\theta_{\rm JN}$ and $\psi$, and the initial precession phase $\alpha_{o}$. It is not difficult to show that $\zeta$ can take any value as the parameters $\theta_{\rm JN}$, $\psi$, $\alpha_{o}$ are varied. For example, at $\theta_{\rm JN}=0$, $\mathcal{A}_{1}$ vanishes and so does $\zeta$, while at $\theta_{\rm JN}=\pi/2$ and $\psi=\pi/4$, $\mathcal{A}_{0}$ vanishes and $\zeta\rightarrow\infty$. Since the initial precession phase $\alpha_{o}$ is a free parameter, the phase of $\zeta$ also can take any value. The overall amplitude and phase of the signal also depends upon the distance and coalescence phase so that any values of the amplitude and phase of the signal in the two harmonics are consistent with a signal.

V Validity of the two-harmonic waveform

To investigate the validity of the two-harmonic approximation, we compare the approximate waveform with the full, five-harmonic, precessing waveform across the parameter space. The error will be of order $b^{2}$, which is small over much of the parameter space, and for the majority of orientations.

Fig. 6 shows the overlap between the full waveform and a subset of the harmonics for a binary with $M=40M_{\odot}$, $q=4$ and $\chi_{\mathrm{eff}}=0$, while varying the orientation and value of $\chi_{P}$. In each case, we calculate,

[TABLE]

where,

[TABLE]

and $S(f)$ is the power spectral density of the detector data. Thus the overlap is maximized over the phase, but not over time or any of the mass and spin parameters. An overlap of close to unity shows that the two waveforms are very similar, while a lower value of overlap implies significant deviations between the waveforms. As a rule of thumb, an overlap $O(h,h^{\prime})\lesssim 1-3/\rho^{2}$ will be observable at a signal to noise ratio $\rho$ Flanagan:1997kp ; Miller:2005qu ; Baird:2012cu .

We calculate the overlap of the full waveform, $h$, against

the leading order waveform in the precession expansion, $h^{0}$; 2. 2.

the dominant harmonic, i.e. the harmonic of $h^{0}$ and $h^{1}$ which contains the largest fraction of the power in the full waveform; 3. 3.

the two-harmonic waveform with the appropriate values of $\mathcal{A}_{0}$ and $\mathcal{A}_{1}$.

For the $+$ polarized waveform (left column), the $k=0$ harmonic is dominant for all values of $\theta_{\rm JN}$ and $\chi_{P}$, so that the observed overlap with the full waveform is above $0.8$ across the parameter space. For $\theta_{\rm JN}\approx 0$ or small values of $\chi_{P}$, the other harmonics make a minimal contribution and the overlap is close to unity. For larger values of $\theta_{\rm JN}$ and $\chi_{P}$ the other harmonics are more significant and the overlap drops to $0.9$ or less. The two-harmonic waveform is a significantly better match to the full waveform, with an overlap greater than $0.99$ for much of the parameter space, and only below $0.9$ for edge-on systems with high $\chi_{P}$ where the $k=2$ harmonic contributes most strongly (and the $k=1$ contribution vanishes).

For the $\times$ polarized waveform (center column), the effect of incorporating the $k=1$ harmonic is dramatic. For $\theta_{\rm JN}=90^{\circ}$ the $k=0$ contribution vanishes and only the $k=1,3$ harmonics are present. Thus, the overlap with harmonic $k=0$ is essentially zero. Using the best of $k=0,1$ provides a good overlap with the edge-on waveform, but there is still a poor overlap at $\theta_{\rm JN}\approx 60^{\circ}$ where both the $k=0$ and $1$ harmonics contribute significantly to the waveform. This effect has been observed previously, for example in Brown:2012gs ; o2020semianalytic and a geometric understanding of its origin provided. The two-harmonic waveform matches remarkably well to the full waveform, with the largest differences for $\theta_{\rm JN}=90^{\circ}$ and $\chi_{P}\approx 1$ where the overlap drops to $0.99$ due to the contribution from the $k=3$ harmonic.

The right column shows the overlap as the orientation of the binary changes. As expected, at points where the $k=0$ harmonic vanishes ($\theta_{\rm JN}=90^{\circ}$ and $\psi=45^{\circ}$), the overlap with this harmonic drops to zero. The dominant harmonic is a good match to the waveform, except for orientations where two harmonics contribute significantly. As discussed in detail in Ref. Brown:2012gs , this corresponds to configurations where the binary orientation passes through the null of the detector response (i.e. the signal goes to zero) once per precession cycle. Thus, the radius of the circle with poor overlaps is approximately equal to the opening angle of the binary. The two-harmonic approximation provides an excellent fit to the full waveform over the majority of orientations, only dropping below $0.95$ for orientations where $\theta_{\rm JN}\rightarrow 90^{\circ}$ and $\psi\approx 0,90^{\circ}$, where the $k=2$ harmonic is most significant.

Next, we investigate the validity of the two-harmonic approximation for a population of binaries. To begin with, let us fix the masses and spins and just consider the effect of binary orientation. As before, we choose $M=40M_{\odot}$, $q=4$, $\chi_{\mathrm{eff}}=0$ and $\chi_{P}=0.6$, corresponding to $\overline{b}\approx 0.3$, with the binary orientation distributed uniformly over $\cos(\theta_{\rm JN}),\phi_{o},\alpha_{o},\psi$. Fig. 7 shows the distribution of the overlap between the full waveform and

1. the $k=0$ harmonic, 2) the dominant harmonic and 3) the two-harmonic approximation. The results are shown for both a uniformly distributed population, and a population of signals observable above a fixed threshold in the detector — thereby favoring orientations that produce the largest amplitude gravitational wave. The median overlap with either the $k=0$ or dominant harmonic is $\lesssim 0.9$, while the two-harmonic approximation improves the median overlap to $0.99$. Using the dominant harmonic, there are a small fraction of signals with overlaps of $0.7$ or lower (and for the $k=0$ harmonic, this tail extends to overlaps of 0.2), while for the two harmonic approximation, the worst overlap is $0.88$.

We can use these results to obtain a rough sense of the benefits of performing a search using the two-harmonic approximation. Previous, more detailed, investigations of this question have been carried out in, e.g. Harry:2016ijz ; Brown:2012gs ; Ajith:2012mn . Current gravitational wave searches make use of spin-aligned waveforms Messick:2016aqy ; Usman:2015kfa , and a precessing waveform will naturally be identified by a spin-aligned waveform which matches well the dominant harmonic. Thus, we can use the overlaps between the precessing waveforms and dominant harmonics as a proxy for the performance of an aligned spin search. Since the median overlap is $0.9$ we would expect to recover approximately 70% as many signals ($\approx 0.9^{3}$ for a population uniform in volume) as with a full precessing search, above a fixed threshold. A search based upon the two-harmonic approximation would recover around 97% of these signals, indicating an improvement of over 30% in sensitivity to such systems.

We also show how the distribution of overlaps varies across the mass and spin parameter space, as encoded by the parameter $\overline{b}$ and plotted for three choices of spin distribution in Figure 4.555While these plots were made with fixed masses and $\chi_{\mathrm{eff}}$, they should give a reasonable indication of the accuracy of the two-harmonic waveform across the mass and spin parameter space, as a function of $\overline{b}$. For different masses and spins, the evolution of the precession angle during the coalescence can have a slight impact upon the relative importance of the modes but, as $b$ typically does not change significantly over the observable waveform, this effect is likely to be small. Furthermore, as different modes are not perfectly orthogonal, the degree to which they are not will also have a small effect upon the results. As shown in Section VII, the harmonics are close to orthogonal for $M\lesssim 40M_{\odot}$ so that the results shown here will be representative, at least at lower masses. For $\overline{b}\lesssim 0.13$ — accounting for three quarters of signals in the low-isotropic population — the median overlap between the dominant harmonic and the full waveform is above $0.97$. Thus, for the majority of expected signals, the spin-aligned search will have good sensitivity. However, even for low values of $\bar{b}$ there will be some orientations of signals where two dominant harmonic will not match the waveform well, while the two-harmonic waveform still provides an essentially perfect representation of the waveform for all orientations. At $\overline{b}\approx 0.25$ the median overlap with the dominant harmonic waveform drops to $0.9$, and it is here that a search with the two-harmonic approximation could provide a 30% improvement. We note, however, that for the low-isotropic distribution this accounts for only 5% of systems. While systems with such significant precession may be rare they would come from interesting areas of parameter space, with high mass ratios and spins. It is only at $\overline{b}=0.4$ that the median overlap for the two harmonic waveform drops to 0.97, indicating a 10% loss relative to an ideal search, but also 70% improvement over a spin-aligned search.

VI Searching for precessing binaries

The two-harmonic approximation provides an ideal basis to develop a search for binaries with precession. The typical approach to searching for binary coalescences has been to generate a template-bank of waveforms that covers the parameter space Owen:1995tm ; Owen:1998dk ; Babak:2006ty . These templates comprise discrete points in the mass and spin space chosen so that the waveform produced by a binary anywhere in the parameter space of interest has a match of at least 97% with one of the templates. The waveform for each template is then match-filtered against the data to identify peaks of high SNR, and various signal consistency and coincidence tests are used to differentiate signals from non-stationary noise transients Allen:2005fk ; Allen:2004gu ; Babak:2012zx ; Messick:2016aqy ; Usman:2015kfa . Current searches make use of a template bank covering the four dimensional mass and aligned-spin space DalCanton:2017ala ; Mukherjee:2018yra .666As we have discussed, the most significant effect on the observed waveform arises due to the effective spin $\chi_{\mathrm{eff}}$, which is a combination of the aligned spin components of the two waveforms. Thus, although the template space is four dimensional, one of the spin directions provides limited variation to the waveforms, and thus is relatively straightforward to cover. The search takes advantage of the fact that changing the sky location, distance and orientation of the binary only changes the overall amplitude and phase of the signal, and these quantities can be maximized over in a simple manner.

When developing a search for precessing binaries, the search becomes more challenging due to the increasing number of parameters. In principle, it is necessary to search over two masses and six spin components, although, in practice it will probably be sufficient to restrict to the masses, $\chi_{\rm eff}$ and $\chi_{P}$. The second complication is that the observed morphology of the waveform varies as the orientation of the binary changes, and it becomes necessary to search over binary orientation $\theta_{\rm JN}$, polarization $\psi$ and precession phase $\alpha_{o}$, although methods have been developed to straightforwardly handle a subset of these parameters Pan:2003qt ; Harry:2016ijz .

The two-harmonic waveform can be used to maximize the SNR over the binary orientation in a simple way. The two complex amplitudes $\mathcal{A}_{0}$ and $\mathcal{A}_{1}$, defined in Eq. (IV), are dependent upon five variables: the distance, $d_{L}$, binary orientation, $\theta_{\rm JN}$, $\psi$, and the initial orbital and precession phases, $\phi_{o}$, $\alpha_{o}$. Since $\mathcal{A}_{0}$ and $\mathcal{A}_{1}$ can take any value in the complex plane, it is possible to construct the two-harmonic SNR by filtering the two harmonics $h_{0}$ and $h_{1}$ against the data and then freely maximizing the amplitudes so that,

[TABLE]

If the harmonics are not orthogonal, the two-harmonic SNR should be calculated using $h^{0}$ and $h_{\perp}^{1}$ — the $k=1$ harmonic with any component proportional to $h^{0}$ removed. The extrinsic parameters of the binary (distance, sky location, orientation, orbital and precession phase) can be searched over through maximization over the amplitudes of the two harmonics, leaving only the masses and spins as dimensions to search using a bank of waveforms.

We must still construct a bank of waveforms to cover the four-dimensional parameter space of masses, the effective aligned $\chi_{\mathrm{eff}}$ and precessing $\chi_{P}$ components of the spins. The amplitude and phase evolution of a single harmonic does not carry the tell-tale amplitude and phase modulation caused by precession, but does have a different phase evolution due to precession Lundgren:2013jla ; Hannam:2013oca . Since the phase evolution of each precessing harmonic is degenerate with a non-precessing waveform with different mass-ratio or effective spin, the bank of templates will essentially be a bank of non-precessing waveforms. This may allow us to reduce the size of the template bank.

The k=0 harmonic of the precessing waveform has an additional phase (see Eq. (16)) of,

[TABLE]

For systems in which orbital angular momentum dominates over spin angular momentum, the precession frequency is inversely proportional to orbital frequency, $\Omega_{P}=\dot{\alpha}\propto f^{-1}$ Apostolatos:1994mx ; Buonanno:2002fy ; Brown:2012gs . This is the same frequency dependence as the 1PN contribution to the waveform, whose amplitude depends upon the mass ratio. Consequently, it is reasonable to expect that the precession-induced phase will be indistinguishable from a systematic offset in the binary mass ratio, or the effective spin Hannam:2013uu . Similarly, the $k=1$ harmonic has essentially the same amplitude evolution as the non-precessing waveform, but with a phase difference of,

[TABLE]

which will also, in many cases, be degenerate with a change in the mass ratio or aligned spin.

In Figure 8, we investigate the degeneracy in the spin ($\chi_{\mathrm{eff}}$–$\chi_{P}$) space of the two leading precession harmonics. We consider a system with masses, $M=40M_{\odot}$ and $q=4$, and spins $\chi_{\mathrm{eff}}=0$, $\chi_{p}=0.6$ and investigate how the two waveform harmonics vary as the spins are changed. The figure shows the match — the overlap maximized over time-offsets — between our fiducial waveform and one with the same masses but different spins. For both harmonics, there is a band in the $\chi_{\mathrm{eff}}$–$\chi_{P}$ plane where the is mismatch is small — the different phase evolution of each harmonic caused by varying $\chi_{P}$ can be offset by a suitable change in $\chi_{\mathrm{eff}}$. The relation is approximately quadratic, $\Delta\chi_{\mathrm{eff}}\propto(\Delta\chi_{P})^{2}$, which is to be expected. Recall, from Eq. (37), that change in phase due to precession is quadratic in $b$, and therefore also in $\chi_{P}$ at least for small values of $b$. Meanwhile the phasing of the waveform varies, at leading order, linearly with $\chi_{\mathrm{eff}}$.

This degeneracy in the $\chi_{\mathrm{eff}}$–$\chi_{P}$ plane suggests that a single template waveform could be used to search over an extended region corresponding, for example, to the region of mismatch $<0.03$ in Figure 8. However, this will only work if the degenerate region for the $k=0$ and $k=1$ harmonics is the same. It is clear from Equations (37) and (38) and Figure 8 that they are not identical. Nonetheless,777Strictly, when doing this comparison, we must use the same time offset for the two harmonics, whereas the figure allows for an independent maximization of the time delay for each harmonic. Fixing a single time delay does slightly decrease the matches, but not significantly enough to change the conclusions. for the example we have considered, the two degenerate regions are similar, and along the line that traces the mid-point between best fit values of $\chi_{\mathrm{eff}}$ for the two harmonics, both harmonics have a match above 0.97 with the initial point. Thus, to an accuracy appropriate for generating a template bank, we can use the two harmonics from a single waveform to cover a band in the $\chi_{\mathrm{eff}}$–$\chi_{P}$ plane which spans all values of $\chi_{P}$. This effectively reduces the dimensionality of the parameter space to three dimensions: mass, mass ratio and one spin parameter.

Our proposal to develop a precessing search is as follows. First, generate a bank of templates to cover the space of non-precessing binaries. At each $M$, $q$, $\chi_{\mathrm{eff}}$ point in the template bank, construct the two-harmonic waveform for a fixed value of $\chi_{P}$. Then, filter the data against the two harmonics and calculate the two-harmonic SNR, as defined in Eq. (36) to identify candidate events in a single detector. It will be necessary to extend the existing $\chi^{2}$ signal consistency test Allen:2004gu to each harmonic, taking into account the presence of the other harmonics, to reduce the impact of non-stationarity in the data. Next, perform coincidence between detectors by requiring a signal in the same template at the same time, up to the allowed time delays based upon speed of propagation. For a non-precessing signal observed in two detectors, the relative amplitude and phase of the SNR in each detector can take any value, even though some are astrophysically more likely Nitz:2017svb (and this can be used to increase search sensitivity). However, for the two-harmonic waveform not every signal observed in two detectors will be compatible with an astrophysical source. This can be seen through simple parameter counting: there are ten measured quantities (two complex amplitudes and a time of arrival in each detector), which depend upon eight parameters, the five orientation parameters ($d_{L}$, $\theta_{\rm JN}$, $\psi$, $\phi_{o}$, $\alpha_{o}$), sky location and merger time. An additional coincidence test to check for consistency between parameters will likely be necessary to reduce the search background. A similar problem arises already in extending the amplitude and phase consistency of Nitz:2017svb to three or more detectors and methods developed for that purpose may be helpful for the precessing search.

We can estimate the likely sensitivity improvement from a precessing search, as we have briefly discussed in Section V. A non-precessing search will typically find the dominant harmonic of the waveform. Thus, for signals where two harmonics provide a significant contribution, a search based on the two-harmonic waveform has the potential to out-perform the non-precessing search. The two-harmonic waveform has four degrees of freedom, encoded in $\mathcal{A}_{0}$ and $\mathcal{A}_{1}$, compared to two for the non-precessing search. Thus, the noise background is higher for the two-harmonic search and, based upon a comparison of the tails of the $\chi^{2}$ distribution with 2 and 4 degrees of freedom, an increase of around 5% in SNR is required to obtain the same false alarm rate (see e.g., Ref. Harry:2016ijz for a discussion of this issue). Thus, a signal will be observed as more significant in the two-harmonic search than a non-precessing search if the SNR can be increased by 5% or more. Fig. 7 shows that this occurs for $\overline{b}\gtrsim 0.15$, and for binaries with $\overline{b}$ above this value the two-harmonic search has the potential to outperform a non-precessing search. We note, however, that a given template will cover a range of spin values and consequently a range of $\overline{b}$, so it may be more appropriate to deploy the two-harmonic search for templates with an average of $\overline{b}$ which is greater than $0.15$.

Another challenge of searches for precessing systems is the associated computational cost Harry:2016ijz , which can be prohibitive. The maximum computational cost for the two-harmonic search would be double that of a comparable non-precessing search: it becomes necessary to filter both the $k=0$ and 1 harmonics, and computational time is dominated by this matched filtering. However, since both the $k=0$ and $k=1$ harmonics are essentially non-precessing waveforms, there may be waveforms associated with the $k=1$ harmonics are already in the set of $k=0$ waveforms, but associated with different parameters. If so, this could further reduce the computational cost.

VII Observability of precession

The two-harmonic approximation allows us to easily identify regions of the binary merger parameter space for which precession will leave an observable imprint on the waveform. Since the amplitude and phase evolution of a single harmonic is generally consistent with that of a non-precessing waveform (see above and Lundgren:2013jla ; o2020semianalytic ), it is only when two harmonics can be observed that we are able to clearly identify precession in the system. We are therefore interested in deriving an expression for the precession SNR, $\rho_{p}$, defined as the SNR in the second most significant harmonic, and determining when it will be observable. If the two harmonics $h^{0}$ and $h^{1}$ in Eq. (31) are orthogonal, then the precession SNR is simply,

[TABLE]

where $\zeta$, defined in Eq. (IV), gives the ratio of the SNR in the $k=1$ and $k=0$ harmonics and $\rho_{\mathrm{2harm}}$ is the total SNR in the two-harmonic waveform.

Let us briefly examine where in parameter space the two harmonics are close to orthogonal. Where there are sufficient precession cycles we expect the two harmonics, $h^{0}$ and $h^{1}$, will be close to orthogonal, and the overlap to be close to zero Lundgren:2013jla . The overlap between the two harmonics for various two-dimensional slices through the parameter space is shown in Fig. 9. At higher masses, where the binary completes one, or fewer, precession cycles in the detector’s sensitive band, there is a larger overlap between the harmonics. At negative $\chi_{\mathrm{eff}}$ and minimal $\chi_{p}$, the overlap is also significant. However, providing the mass of the system is below $50M_{\odot}$, for the much of the parameter space the overlap is less than $0.1$ and simple expression in Eq. (39) will be applicable.

Taking into account the overlap between harmonics, the total power in the two-harmonic waveform is,

[TABLE]

where $o_{1,0}$ is complex overlap between the two harmonics:

[TABLE]

We can project the SNR onto directions parallel and perpendicular to the $h^{0}$ waveform to obtain the SNR in these two directions as,

[TABLE]

Similarly, the power parallel to and perpendicular to the $k=1$ harmonic is,

[TABLE]

The precession SNR is defined as the power orthogonal to the dominant harmonic,888In exceptional circumstances, where the overlap is large and $\zeta o_{1,0}$ is close to $-1$, there can be more power in $\rho_{\perp,i}$ than $\rho_{i}$. In such cases, it is natural to use $\rho_{i}$ to determine if precession is present, although this is not ideal as $\rho_{\perp,i}$ need not resemble a non-precessing waveform.

[TABLE]

As expected, the precession SNR scales with the total SNR of the signal, so that precession will be more easily observed for louder events. If there is significant degeneracy between the harmonics, the numerator will be reduced, making the observation of precession more difficult. Finally, in the limit that $o_{1,0}\rightarrow 0$, the expression simplifies to the one given earlier for orthogonal harmonics in (39), as expected.

What value of $\rho_{p}$ will be required to observe precession? This will happen if the evidence for a signal with $\chi_{p}\neq 0$ in the data is greater than that for a non-precessing source. This can be evaluated through Bayesian model selection, by considering the Bayes factor between the hypotheses. However, such a calculation requires a full exploration of the parameter space. We can, instead, obtain an approximate answer by considering the maximum likelihood. Since the two-harmonic waveform is more general than the non-precessing waveform, it will always give a larger maximum likelihood even in the absence of precession do to its ability to fit the detector noise. Thus, we are interested in examining the expected increase in SNR due to the inclusion of the second harmonic, in the absence of any power in it.

The two-harmonic SNR can be written as

[TABLE]

where $\rho_{np}$ is the non-precessing SNR or, equivalently, the SNR in the dominant harmonic. In the absence of precession, $\rho_{p}$ will be $\chi^{2}$ distributed with 2 degrees of freedom, as we are able to freely maximize over the amplitude and phase of the two harmonics independently Allen:2005fk ; Babak:2012zx . Consequently, in 90% of cases, noise alone will give a value of $\rho_{p}<2.1$. Therefore, as a simple criterion, we require that,

[TABLE]

for precession to be observable. In Ref. Green:2019_prec_pe we use this definition to investigate in detail the observability of precession over the parameter space.

VIII Discussion

We have presented a new, intuitive way to understand the observability of precession in GW observations. By keeping only the leading precession term, we have derived a precession SNR and argued that this can be used to determine when precession will be observable. Before discussing applications we point out the main limitations of this analysis. As is clear from the formulation, this analysis works best for binaries where $b=\tan(\beta/2)$ is small. This typically corresponds to situations where the masses are comparable, the precessing spin is small and any aligned component of the spin is aligned (rather than anti-aligned) with the orbital angular momentum. We have shown above that this assumption is valid for a reasonable population.

We now point to several advantages and applications of this formulation:

First, it gives new understanding of the observability of precession, and also of the origin of precession as the beating of two waveform components with slightly differing frequencies (also discussed in Lundgren:2013jla ). It is difficult to identify the presence of precession in a GW observation directly from $\chi_{P}$, since the prior astrophysical expectation disfavours $\chi_{P}=0$. While the deviation from the prior can be determined through the Bayes factor, the results in this paper suggest that the precession SNR $\rho_{P}$ could provide a direct measure of whether precession has been measured in a signal. The potential applications of $\rho_{P}$ are discussed in the companion paper fairhurst2019will , and will be investigated in more detail in Ref. Green:2019_prec_pe , where we probe the measurability of precession across the gravitational wave parameter space.

There exist a number of detailed population analyses which extract the features of the underlying population of gravitational waves from the set of observed gravitational wave events, for example Talbot:2017yur ; Talbot:2019okv ; Wysocki:2018mpo ; LIGOScientific:2018jsj . These typically use the full posterior distributions recovered from the gravitational wave signal Veitch:2014wba ; Ashton:2018jfp to infer the population and, as such, naturally account for precession effects in the observed signals when inferring the black hole mass and spin populations. Nonetheless, there have been a number of studies performed which investigate the population properties using a subset of the recovered parameters, see e.g. Farr:2017gtv ; Farr:2017uvj ; Fishbach:2017zga ; Fishbach:2018edt ; Tiwari:2018qch ; LIGOScientific:2018jsj , and have been successfully used to infer interesting properties of the mass and spin distributions. The majority of these studies have restricted attention to the aligned components of the spins. The precession SNR provides a straightforward method to determine the significance of precession, and provides away to probe observability of precession in populations of binaries. In using this method we have been able to derive constraints on the preferred spin distribution including precession effects fairhurst2019will .

Both of the applications highlighted above are currently possible using other more sophisticated but computationally expensive methods such as Bayesian model comparison. This is, of course, a more general method that makes fewer assumptions than we do in computing $\rho_{p}$, however the computational costs associated with calculating the marginal likelihood over multiple, e.g. precessing and non-precessing, models per binary are not feasible for a large number of binaries. For example the analysis in fairhurst2019will involved calculating $\rho_{p}$ for 1 million binaries, and computing the Bayes factor for 1 million binaries would certainly not be practical. Similar, lightweight analyses, could also be developed using the formalism introduced in, e.g. Brown:2012gs , and if this is done,it would be interesting to compare them with the results from the two harmonic analysis.

Finally, we have outlined a method by which the two-harmonic approximation could be used to develop a search for precessing binaries. We have shown that in principle that this approach could result in a significant increase in sensitivity without the computational overheads associated with other precessing search methods. In addition, the formalism should provide a way to identify the parts of parameter space where a precessing search is likely to increase sensitivity. We plan a detailed investigation into the feasibility of a precessing search based upon the two-harmonic approximation in future work.

Acknowledgements

We would like to thank Eleanor Hamilton, Ian Harry, Andrew Lundgren, Frank Ohme, Francesco Pannarale, Vivien Raymond, Patrick Sutton and Vaibhav Tiwari for discussions. This work was supported by Science and Technology Facilities Council (STFC) grant ST/L000962/1, European Research Council Consolidator Grant 647839 and the National Science Foundation under Grant No. NSF PHY-1748958. We are grateful for computational resources provided by Cardiff University, and funded by an STFC grant supporting UK Involvement in the Operation of Advanced LIGO.

Appendix: Derivation using spin-weighted spherical harmonics

In this appendix, we provide an alternative derivation of the power series decomposition of the precessing waveform, given in Section III, based upon the spin-weighted spherical harmonic decomposition of the waveform Thorne:1980ru and its application to precession as described in Hannam:2013oca ; Boyle:2011gg . Specifically, we wish to obtain the result in Eq. (14). Throughout, we follow the notation used in Khan:2019kot .

The gravitational waveform emitted during a binary merger,

[TABLE]

can naturally be decomposed into a set of spin-weighted spherical harmonics as

[TABLE]

where $\theta$ and $\phi$ give the orientation of the observer relative to a co-ordinate system used to identify the spherical harmonics, $\vec{\lambda}$ encodes the physical parameters of the system (masses, spins, etc) and $t$ is the time.

The multipoles for a precessing system are approximated by “twisting up” Hannam:2013oca ; Boyle:2011gg the multipoles of the non-precessing counterpart based upon the orientation of the orbital angular momentum given by the opening angle $\beta$, precession angle $\alpha$ and the third Euler angle $\epsilon$ defined via

[TABLE]

Then, the precessing multipoles are given by

[TABLE]

where the Wigner D-matrix is

[TABLE]

and the Wigner d-matrix given, for example, in Brown:2007jx .

Combining these decompositions gives the waveform for a precessing binary as

[TABLE]

In performing the twisting, it’s natural that the precessing waveform is described in a coordinate system aligned with the orbital angular momentum, so that $\theta=\theta_{\mathrm{JN}}$. Furthermore, the orientation of the $x$-axis will be specified relative to the (initial) precession phase so that $\phi=-\alpha_{o}$.

In this work, we restrict attention to the case where the non-precessing model contains only the $\ell=2$ and $n=\pm 2$ modes, and require symmetry in gravitational wave emission above and below the plane of the binary so that $h_{\ell,n}=(-1)^{\ell}h^{\star}_{\ell,-n}$. This eliminates the sum over $\ell$ and $m$ from Eq. (52). Furthermore, we can expand the spherical harmonics using

[TABLE]

to obtain

[TABLE]

We now wish to re-write the above to show that the waveform can be decomposed in modes whose amplitudes form a power series in $b=\tan(\beta/2)$. To do so, we note that the Wigner d-matrices can be evaluated as powers of $\sin(\beta/2)$ and $\cos(\beta/2)$, so that if we are able to group terms with the same indices we will arrive at the desired expression. We do this by using the d-matrix identities:

[TABLE]

and relabelling the dummy index $m\rightarrow-m$ in the second term of Eq. (Appendix: Derivation using spin-weighted spherical harmonics) to obtain:

[TABLE]

Finally, we can evaluate the Wigner d-matrices as

[TABLE]

where $C_{\pm 2}=1$, $C_{\pm 1}=2$, $C_{0}=\sqrt{6}$ and, as before, $b=\tan(\beta/2)$. Similarly, we introduce $\tau=\tan\theta_{JN}/2$, and evaluate the d-matrices for the angle $\theta_{\mathrm{JN}}$. This gives

[TABLE]

This is close to the desired form and, in particular, we have obtained an decomposition where the relative strength of each mode is decreased by a factor of $b$. To obtain an expression comparable to Eq. (14) we must evaluate the waveform observed at a detector with response $F_{+}$ and $F_{\times}$ to the two gravitational polarizations.

[TABLE]

Then, to finally equate this with the desired expression, we must make the identification

[TABLE]

where $\Phi_{S}$ is defined in Eq. (7). Thus the amplitude of the waveform, $A_{o}(t)$ is the same as the scaled $22$ mode while the phase of the 22 mode is the (negative) of the orbital phase. Furthermore, it is straightforward to show that the $\mathcal{A}^{+,\times}_{k}$ coefficients are given by

[TABLE]

Substituting these identifications, we obtain the desired expression for the waveform observed at a detector,

[TABLE]