A data-analysis driven comparison of analytic and numerical coalescing binary
waveforms: nonspinning case
Yi Pan,1 Alessandra Buonanno,1 John G. Baker,2 Joan Centrella,2 Bernard
J. Kelly,2 Sean T. McWilliams,1 Frans Pretorius,3 and James R. van Meter2, 4
Department of Physics, University of Maryland, College Park, MD 20742
Gravitational Astrophysics Laboratory, NASA Goddard Space Flight Center, 8800 Greenbelt Rd., Greenbelt, MD 20771
Department of Physics, Princeton University, Princeton, NJ 08544
Center for Space Science & Technology, University of Maryland Baltimore County,
Physics Department, 1000 Hilltop Circle, Baltimore, MD 21250
(Dated: October 22, 2018)
We compare waveforms obtained by numerically evolving nonspinning binary black holes to post-
Newtonian (PN) template families currently used in the search for gravitational waves by ground-
based detectors. We find that the time-domain 3.5PN template family, which includes the inspiral
phase, has fitting factors (FFs) ≥ 0.96 for binary systems with total massM = 10–20M⊙. The time-
domain 3.5PN effective-one-body template family, which includes the inspiral, merger and ring-down
phases, gives satisfactory signal-matching performance with FFs ≥ 0.96 for binary systems with total
mass M = 10–120M⊙. If we introduce a cutoff frequency properly adjusted to the final black-hole
ring-down frequency, we find that the frequency-domain stationary-phase-approximated template
family at 3.5PN order has FFs ≥ 0.96 for binary systems with total mass M = 10–20M⊙. However,
to obtain high matching performances for larger binary masses, we need to either extend this family
to unphysical regions of the parameter space or introduce a 4PN order coefficient in the frequency-
domain GW phase. Finally, we find that the phenomenological Buonanno-Chen-Vallisneri family
has FFs ≥ 0.97 with total mass M = 10–120M⊙. The main analyses use the noise spectral-density
of LIGO, but several tests are extended to VIRGO and advanced LIGO noise-spectral densities.
PACS numbers: 04.25.Dm, 04.30.Db, 04.70.Bw, x04.25.Nx, 04.30.-w, 04.80.Nn, 95.55.Ym
I. INTRODUCTION
The search for gravitational-waves (GWs) from coa-
lescing binary systems with laser interferometer GW de-
tectors [1, 2, 3, 4, 5] is based on the matched-filtering
technique, which requires accurate knowledge of the
waveform of the incoming signal. In the last couple of
years there have been several breakthroughs in numeri-
cal relativity (NR) [6, 7, 8], and now independent groups
are able to simulate the inspiral, merger and ring-down
phases of generic black- hole (BH) merger scenarios, in-
cluding different spin orientations and mass ratios [9].
However, the high computational cost of running such
simulations makes it difficult to generate sufficiently long
inspiral waveforms that cover the parameter space of as-
trophysical interest.
References [10, 11] found good agreement between ana-
lytic (based on the post-Newtonian (PN) expansion) and
numerical waveforms emitted during the inspiral, merger
and ring-down phases of equal-mass, nonspinning binary
BHs. Notably, the best agreement is obtained with 3PN
or 3.5PN adiabatic waveforms [12] (henceforth denoted
as Taylor PN waveforms) and 3.5PN effective-one-body
(EOB) waveforms [13, 14, 15, 16, 17, 18, 19]. In addi-
tion to the inspiral phase the latter waveforms include
the merger and ring-down phases. Those comparisons
suggested that it should be possible to design hybrid nu-
merical/analytic templates, or even purely analytic tem-
plates with the full numerics used to guide the patching
together of the inspiral and ring-down waveforms. This is
an important avenue to template construction as eventu-
ally thousands of waveform templates may be needed to
extract the signal from the noise, an impossible demand
for NR alone. Once available, those templates could be
used by ground-based laser interferometer GW detectors,
such as LIGO, VIRGO, GEO and TAMA, and in the fu-
ture by the laser interferometer space antenna (LISA) for
detecting GWs emitted by solar mass and supermassive
binary BHs, respectively.
This paper presents a first attempt at investigating the
closeness of the template families currently used in GW
inspiral searches to waveforms generated by NR simula-
tions. Based on this investigation, we shall propose ad-
justments to the templates so that they include merger
and ring-down phases. In contrast, Ref. [21] examined
the use of numerical waveforms in inspiral searches, and
compared numerical waveforms to the ring-down tem-
plates currently used in burst searches. Similar to the
methodology presented here, fitting factors (FFs) [see
Eq. (2) below] are used in Ref. [21] to quantify the ac-
curacy of numerical waveforms for the purpose of de-
tection, as well as the overlap of burst templates with
the waveforms. Reference [21] found that by computing
FFs between numerical waveforms from different reso-
lution simulations of a given event, one can recast the
numerical error as a maximum FF that the numerical
waveform can resolve. In other words, any other tem-
plate or putative signal convolved with the highest reso-
lution numerical simulation that gives a FF equal to or
larger than this maximum FF is, for the purpose of de-
http://arxiv.org/abs/0704.1964v2
tection, indistinguishable from the numerical waveform.
We will explore this aspect of the problem briefly. The
primary conclusions we will draw from the analysis do
not crucially depend on the exactness of the numerical
waveforms. What counts here is that the templates can
capture the dominant spectral characteristics of the true
waveform.
For our analysis we shall focus on two nonspinning
equal-mass binary simulation waveforms which differ in
length, initial conditions, and the evolution codes used
to compute them: Cook-Pfeiffer quasi-equilibrium ini-
tial data built on Refs. [22, 23, 24, 25, 26] evolved with
Pretorius’ generalized harmonic code [6], and Brandt-
Brügmann puncture data [27] evolved using the Goddard
group’s moving-puncture code [8]. We also consider two
nonspinning unequal-mass binary simulations with mass
ratios m2/m1 = 1.5 and m2/m1 = 2 produced by the
Goddard group.
The paper is organized as follows. In Sec. II we dis-
cuss the phase differences between PN inspiraling tem-
plates. In Sec. III we build hybrid waveforms by stitch-
ing together PN and NR waveforms. We try to under-
stand how many NR cycles are needed to obtain good
agreement between NR and PN waveforms, to offer a
guide for how long PN waveforms can be used as accu-
rate templates. In Sec. IV we compute the FFs between
several PN template families and NR waveforms. We
first focus on low-mass binary systems with total mass
M = 10–30M⊙, then high-mass binary systems with to-
tal mass M = 30–120M⊙. Finally, Sec. V contains our
main conclusions. In Appendix A we comment on how
different representations of the energy-balance equations
give GW frequencies closer to or farther from the NR
ones.
II. PHASE DIFFERENCES IN
POST-NEWTONIAN INSPIRALING MODELS
Starting from Ref. [28], which pointed out the impor-
tance of predicting GW phasing with the highest possi-
ble accuracy when building GW templates, many subse-
quent studies [14, 18, 19, 20, 29, 31, 32, 33] (those ref-
erences are restricted to the nonspinning case) focused
on this issue and thoroughly tested the accuracy of those
templates, proposing improved representations of them.
These questions were motivated by the observation that
comparable-mass binary systems with total mass higher
than 30M⊙ merge in-band with the highest signal-to-
noise ratio (SNR) for LIGO detectors, It follows that the
corresponding templates demand an improved analysis.
In the absence of NR results and under the urgency
of providing templates to search for comparable-mass bi-
nary BHs, the analytic PN community pushed PN cal-
culations to higher PN orders, notably 3.5PN order [12],
and also proposed ways of resumming the PN expan-
sion, either for conservative dynamics (the EOB ap-
proach [13, 16, 17]), radiation-reaction effects (the Padé
resummation [19]), or both [14, 18]. Those results lead to
several conclusions: (i) 3PN terms improve the compar-
ison between analytic and (numerical) quasi-equilibrium
predictions [23, 26, 34, 35]; (ii) Taylor expanded and re-
summed PN predictions for equal-mass binary systems
are much closer at 3.5PN order than at previous PN
orders, indicating a convergence between the different
schemes [18, 20, 31]; (iii) the two-body motion is quasi-
circular until the end of a rather blurred plunge [14], (iv)
the transition to ring-down can be described by an ex-
tremely short merger phase [14, 18]. Today, with the NR
results we are in a position to sharpen the above con-
clusions, and to start to assess the closeness of analytic
templates to numerical waveforms.
Henceforth, we restrict the analysis to the three time-
domain physical template families which are closest to
NR results [10, 11]: the adiabatic Taylor PN model
(Tpn) [see, e.g., Eqs. (1), (10), and (11)–(13) in Ref. [30]]
computed at 3PN and 3.5PN order, and the nonadia-
batic EOB model (Epn) [see e.g., Eqs. (3.41)–(3.44) in
Ref. [14]] computed at 3.5 PN order. We shall denote
our models as Tpn(n) and Epn(n), n being the PN or-
der. The Tpn model is obtained by solving a particular
representation of the balance equation. In Appendix A
we briefly discuss how time-domain PN models based on
different representations of the energy-balance equation
would compare with NR results.
The waveforms we use are always derived in the so-
called restricted approximation which uses the amplitude
at Newtonian order and the phase at the highest PN
order available. They are computed by solving PN dy-
namical equations providing the instantaneous frequency
ω(t) and phase φ(t) = φ0 +
ω(t′)dt′, thus
h(t) = Aω(t)2/3 cos[2φ(t)] , (1)
where t0 and φ0 are the initial time and phase, respec-
tively, and A is a constant amplitude, irrelevant to our
discussion. The inclusion of higher-order PN corrections
to the amplitude can be rather important for certain
unequal-mass binary systems, and will be the subject of
a future study.
When measuring the differences between waveforms
we weight them by the power spectral-density (PSDs) of
the detector, and compute the widely used fitting factor
(FF) (i.e., the ambiguity function or normalized over-
lap), or equivalently the mismatch defined as 1-FF. Fol-
lowing the standard formalism of matched-filtering [see,
e.g., Refs. [19, 31, 36]], we define the FF as the overlap
〈h1(t), h2(t)〉 between the waveforms h1(t) and h2(t):
〈h1(t), h2(t)〉 ≡ 4Re
h̃1(f)h̃
Sh(f)
FF ≡ max
t0,φ0,λi
〈h1, h2(t0, φ0, λ
〈h1, h1〉〈h2(t0, φ0, λi), h2(t0, φ0, λi)〉
where h̃i(f) is the Fourier transform of hi(t), and Sh(f)
is the detector’s PSD. Thus, the FF is the normalized
overlap between a target waveform h1(t) and a set of tem-
plate waveforms h2(t0, φ0, λ
i) maximized over the initial
time t0, initial phase φ0, and other parameters λ
i. Some-
times we are interested in FFs that are optimized only
over t0 and φ0; we shall denote these as FF0. For data-
analysis purposes, the FF has more direct meaning than
the phase evolution of the waveforms, since it takes into
account the PSDs and is proportional to the SNR of the
filtered signal. Since the event rate is proportional to
the cube of the SNR, and thus to the cube of the FF, a
FF= 0.97 corresponds to a loss of event rates of ∼ 10%.
A template waveform is considered a satisfactory repre-
sentation of the target waveform when the FF is larger
than 0.97.
When comparing two families of waveforms, the FF is
optimized over the initial phase of the template wave-
form, and we also need to specify the initial phase of the
target waveform. Since there is no preferred initial phase
of the target, two options are usually adopted: (i) the ini-
tial phase maximizes the FF or (ii) it minimizes the FF.
The resulting FFs are referred to as the best and mini-
max FFs [29]. All FFs we present in this paper are min-
imax FFs. Although the FF of two waveform families is
generally asymmetric under interchange of the template
family [31], the best and the minimax FF0s are symmet-
ric (see Appendix B of Ref. [29] for details). Henceforth,
when comparing two waveform families using FF0, we do
not need to specify which family is the target.
We shall consider three interferometric GW detectors:
LIGO, advanced LIGO and VIRGO. The latter two have
better low-frequency sensitivity and broader bandwidth.
For LIGO, we use the analytic fit to the LIGO de-
sign PSD given in Ref. [20]; for advanced LIGO we use
the broadband configuration PSD given in Ref. [37]; for
VIRGO we use the PSD given in Ref. [20].
In Fig. 1, we show the FF0s as functions of the ac-
cumulated difference in the number of GW cycles be-
tween waveforms generated with different inspiraling PN
models and for binary systems with different component
masses. We first generate two waveforms by evolving
two PN models, say, “PN1” and “PN2” which start at
the same GW frequency fGW = 30Hz and have the same
initial phase. The two waveforms are terminated at the
same ending frequency fGW = fend up to a maximum
fend,max = min(fend,PN1 , fend,PN2), where fend,PN is the
frequency at which the PN inspiraling model ends. (For
Tpn models this is the frequency at which the PN energy
has a minimum; for Epn models it is the EOB light-ring
frequency.) Then, we compute the difference in phase
and number of GW cycles accumulated until the ending
frequency
∆NGW =
[φPN1(fend)− φPN2(fend)] . (3)
By varying fend (up to fend,max) ∆NGW changes, though
not necessarily monotonically. Although there seems to
be a loose correlation between the FF0s and ∆NGW, it
is hard to quantify it as a one-to-one correspondence.
0 0.2 0.4 0.6 0.8 1.0 1.2
Tpn(3.5) & Epn(3.5), (3+3)
Tpn(3.5) & Epn(3.5), (15+3)
Tpn(3.5) & Epn(3.5), (15+15)
Tpn(3) & Tpn(3.5), (3+3)
Tpn(3) & Tpn(3.5), (15+3)
Tpn(3) & Tpn(3.5), (15+15)
FIG. 1: We show FF0s between waveforms generated from
the three PN models Tpn(3), Tpn(3.5) and Epn(3.5) versus
∆NGW [see Eq. (3)]. The FF0s are evaluated with LIGO’s
PSD. Note that for Tpn(3.5) and Epn(3.5) and a (15+3)M⊙
binary, the lowest FF0 is 0.78 and the difference in the number
of GW cycles ∆NGW ≃ 2. In the limit ∆NGW → 0, the FF0
goes to unity.
For example, a phase difference of about half a GW cy-
cle (∆NGW ≃ 0.5) is usually thought to be significant.
However, here we find relatively high FF0s between 0.97
and > 0.99, depending on the masses of the binary and
the specific PN model used. This happens because the
FF between two waveforms is not determined by the to-
tal phase difference accumulated, but rather by how the
phase difference accumulates across the detector’s most
sensitive frequency band. The relation between FFs and
phase differences is also blurred by the maximization
over the initial time and phase: shifting the phase by
half a cycle from the most sensitive band to a less sen-
sitive band can increase the matching significantly. We
conclude that with LIGO’s PSD, after maximizing only
on initial phase and time, Epn(3.5) and Tpn(3.5) tem-
plates are close to each other for comparable-mass bi-
nary systems M = 6–30M⊙ with FF0>∼ 0.97, but they
can be different for mass ratios m2/m1 ≃ 0.3 with FF0
as low as ≃ 0.8. Tpn(3) and Tpn(3.5) templates have
FF0>∼ 0.97 for the binary masses considered. Note that
for m2/m1 = 1 [≃ 0.3] binary systems, Tpn(3.5) is closer
to Epn(3.5) [Tpn(3)] than to Tpn(3) [Epn(3.5)]. Note
also that when maximizing on binary masses the FFs can
increase significantly, for instance, for a (15 + 3)M⊙ bi-
nary, the FF between Tpn(3.5) and Epn(3.5) waveforms
becomes > 0.995, whereas FF0 ≃ 0.8.
III. BUILDING AND COMPARING HYBRID
WAVEFORMS
Recent comparisons [10, 11] between analytic and nu-
merical inspiraling waveforms of nonspinning, equal-mass
binary systems have shown that numerical waveforms are
in good agreement with Epn(3.5), Tpn(3) and Tpn(3.5)
waveforms. Those results were assessed using eight and
sixteen numerical inspiral GW cycles. Can we conclude
from these analyses that Epn(3.5), Tpn(3.5) and Tpn(3)
can safely be used to build a template bank for detecting
inspiraling GW signals? A way to address this question
is to compute the mismatch between hybrid waveforms
built by attaching either Epn or Tpn waveforms to the
same numerical waveform, and varying the time when
the attachment is made. This is equivalent to varying
the number of numerical GW cycles n in the hybrid tem-
plate. The larger n the smaller the mismatch, as we are
using the same numerical segment in both waveforms.
For a desired maximum mismatch, say 3%, we can then
find the smallest number n of numerical cycles that is
required in the hybrid waveform. This number will, of
course, depend on the binary mass and the PSD of each
detector.
A. Hybrid waveforms
We build hybrid waveforms by connecting PN wave-
forms to NR waveforms at a chosen point in the late
inspiral stage. As mentioned before, we use NR wave-
forms generated with Pretorius’ [10] code and the God-
dard group’s [38] code. Pretorius’ waveform is from an
equal-mass binary with total mass M , and equal, co-
rotating spins (a = 0.06). The simulation lasts ≃ 671M ,
and the waveform has ≃ 8 cycles before the formation of
the common apparent horizon. The Goddard waveform
refers to an equal-mass nonspinning binary. The simula-
tion lasts about ≃ 1516M , and the waveform has ≃ 16
cycles before merger.
Since we will present results from these two waveforms
it is useful to first compare them by computing the FF0.
Although the binary parameters considered by Pretorius
and Goddard are slightly different, we expect the wave-
forms, especially around the merger stage, to be fairly
close. Comparisons between (shorter) waveforms com-
puted with moving punctures and generalized-harmonic
gauge were reported in Ref. [39], where the authors dis-
cussed the different initial conditions, wave extraction
techniques, and compared the phase, amplitude and fre-
quency evolutions. Since the two simulations use dif-
ferent initial conditions and last for different amounts of
time we cut the longer Goddard waveform at roughly the
frequency where the Pretorius waveform starts. In this
way we compare waveforms that have the same length
between the initial time and the time at which the wave
amplitude reaches its maximum. In Fig. 2, we show the
FF0 as function of the total binary mass. Despite dif-
ferences in the two simulations the FF0s are rather high.
The waveforms differ more significantly at lower frequen-
cies. Indeed, as the total mass decreases the FF0s also
decrease as these early parts of the waveform contribute
more to the signal power given LIGO’s PSD.
Any waveform extracted from a numerical simulation
will inherit truncation errors, affecting both the wave-
form’s amplitude and phase [10, 21, 38]. To check
whether those differences would change the results of the
comparisons between NR and PN waveforms, we plot
in Fig. 2 the FF0s versus total binary mass between
two Goddard waveforms generated from a high and a
medium resolution run [38]. The FF0s are extremely high
(> 0.995).
Based on the comparisons between high and medium
resolution waveforms, we can estimate the FFs between
high resolution and exact waveforms. If we have several
simulations with different resolutions, specified by the
mesh-spacings xi, and xi are sufficiently small, we can
assume that the waveforms hi are given by
hi = h0 + x
i hd , (4)
where n is the convergence factor of the waveform, h0 is
the exact waveform generated from the infinite resolution
run (x0 → 0), and hd is the leading order truncation error
contribution to the waveform and is independent of the
mesh spacing xi. We find that the mismatch between the
waveforms hi and hj, 1− FFij , scales as
1− FFij ∝ (x
i − x
2 . (5)
In the Goddard simulations, the high and medium reso-
lution runs have mesh-spacing ratio xh/xm = 5/6, and
the waveform convergence rate is n = 4 [38]. The FF
between the high resolution and exact waveforms hh and
h0 is given by
FF0h = 1− 0.87(1− FFhm) , (6)
where FFhm is the FF between the high and medium
resolution waveforms hh and hm. That is to say, the
mismatch between hh and h0 is slightly smaller than that
between hh and hm, where the latter can be derived from
the FFs shown in Fig. 2. Henceforth, we shall always use
high-resolution waveforms. A similar calculation for Pre-
torius’ waveform gives FF0h = 1−0.64(1−FFhm), though
here xh/xm = 2/3 and n = 2. See Fig. 6 of Ref. [21] for a
plot of FFhm calculated from the evolution of the Cook-
Pfeiffer initial data 1; there FF0 ranges from ≈ 0.97 for
M/Ms = 30 to ≈ 0.99 for M/Ms = 100. In other words,
the mismatch between Goddard’s and Pretorius’ wave-
form shown in Fig.2 is less than the estimated mismatch
from numerical error in the latter waveform.
We build hybrid waveforms by stitching together the
PN and NR waveforms computed for binary systems with
the same parameters. At the point where we connect the
1 The plot in Ref. [21] is for “d=16” corotating Cook-Pfeiffer initial
data, whereas the results presented here are from “d=19” initial
data. However, the resolutions used for both sets were the same,
and thus the mismatches should be similar, in particular in the
higher mass range.
30 40 50 60 70 80 90 100
Total mass of the binary (      )
0.980
0.985
0.990
0.995
1.000
Goddard high and medium resolutions
Pretorius and Goddard waveforms
FIG. 2: FF0 between NR waveforms as a function of the
binary total-mass M . The solid curve are generated for wave-
forms from Pretorius and the Goddard group. The longer
Goddard waveform is shortened such that both waveforms last
≃ 671M and contain ≃ 8 cycles. The dashed curve is gen-
erated for waveforms from the high-resolution and medium
resolution simulations of the Goddard group. All FFs are
evaluated using LIGO’s PSD.
two waveforms, we tune the initial time t0 so that the
frequency of the PN waveform is almost the same as the
frequency of the NR waveform (there is a subtlety trying
to match exactly the frequencies that is discussed at the
end of this section). The initial phase φ0 is then chosen
so that the strain of the hybrid waveform is continuous
at the connecting point.
In Fig. 3, we show two examples of hybrid waveforms of
an equal-mass binary. We stitch the waveforms at points
where effects due to the initial-data transient pulse are
negligible. We find an amplitude difference on the order
of ∼ 10% between the Goddard waveform and the re-
stricted PN waveform. This difference is also present in
Pretorius’ waveform, but it is somewhat compensated for
by amplitude modulations caused by eccentricity in the
initial data. In Ref. [38] it was shown that PN waveforms
with 2.5PN amplitude corrections give better agreement
(see e.g., Fig. 12 in Ref. [38]). However, the maximum
amplitude errors in the waveforms are also on the or-
der of 10% [10, 38]. Since neither 2PN nor other lower
PN order corrections to the amplitudes are closer to the
2.5PN order, we cannot conclude that 2.5PN amplitude
corrections best approximate the numerical waves. Thus,
we decide to use two sets of hybrid waveforms: one con-
structed with restricted PN waveforms, and the other
with restricted PN waveforms rescaled by a single ampli-
tude factor, which eliminates amplitude differences with
the NR waveforms. We shall see that the difference be-
tween these two cases is small for the purpose of our
tests.
The amplitude difference between PN and NR wave-
forms is computed at the same connecting-point GW fre-
-500-1000
-1000-1500
-1500-2000
-2000-2500
-2500
Time (M)
Time (M)
FIG. 3: We show two examples of hybrid waveforms, start-
ing from 40Hz. The PN waveforms are generated with the
Tpn(3.5) model, and the NR waveforms in the upper and
lower panels are generated from Pretorius’ and Goddard’s
simulations, respectively. We mark with a dot the point where
we connect the PN and NR waveforms.
quency. There is another effect which causes a jump
in the hybrid-waveform amplitude. This is a small fre-
quency difference between PN and NR waveforms at the
connecting point. All our NR waveforms contain small
eccentricities [10, 38]. As a consequence, the frequency
evolution ω(t) oscillates. To reduce this effect we fol-
low what is done in Ref. [38] and fit the frequency to a
monotonic quartic function. When building the hybrid
waveform, we adjust the PN frequency to match the quar-
tic fitted frequency (instead of the oscillatory, numerical
frequency) at the connecting point. Since the restricted
PN amplitude is proportional to ω2/3(t) [see Eq. (1)],
this slight difference between ωs at the connecting point
creates another difference between the NR and PN am-
plitudes. Nevertheless, this difference is usually smaller
(for Goddard’s waveform) or comparable (for Pretorius’)
to the amplitude difference discussed above.
-25 -20 -15 -10 -5 0
Time (100M)
Original
Whitened
-25 -20 -15 -10 -5 0
Time (100M)
-800 -600 -400 -200 0
Time (M)
-800 -600 -400 -200 0
Time (M)
-200 -150 -100 -50 0 50
Time (M)
-200 -150 -100 -50 0 50
Time (M)
(5+5)
(10+10) (20+20)
(15+15) (30+30)
(50+50)
FIG. 4: Distribution of GW signal power. In each panel, we plot a hybrid waveform (a Tpn waveform stitched to the Goddard
waveform) in both its original form (blue curve) and its “whitened” form (red curve) [44]. We show waveforms from six
binary systems with total masses 10M⊙ 20M⊙, 30M⊙, 40M⊙, 60M⊙ and 100M⊙. The vertical lines divide the waveforms into
segments, where each segment contributes 10% of the total signal power.
B. Distribution of signal power in gravitational
waveforms
To better understand the results of the FFs between
hybrid waveforms, we want to compute how many sig-
nificant GW cycles are in the LIGO frequency band. By
significant GW cycles we mean the cycles that contribute
most to the signal power, or to the SNR of the filtered
signal. Since GW frequencies are scaled by the total bi-
nary mass, the answer to this question depends on both
the PSD and the total mass of a binary.
In Fig. 4, we show the effect of the LIGO PSD on the
distribution of signal power for several waves emitted by
coalescing binary systems with different total masses. In
each panel, we plot a hybrid waveform (a Tpn waveform
stitched to the Goddard waveform) in both its original
form and its “whitened” form [44]. The whitened wave-
form is generated by Fourier-transforming the original
waveform into the frequency domain, rescaling it by a fac-
tor of 1/
Sh(f), and then inverse-Fourier-transforming
it back to the time domain. The reference time t = 0 is
the peak in the amplitude of the unwhitened waveform.
The amplitude of a segment of the whitened waveform
indicates the relative contribution of that segment to the
signal power and takes into account LIGO’s PSD. Both
waveforms are plotted with arbitrary amplitudes, and the
unwhitened one always has the larger amplitude. The
absolute amplitude of a waveform, or equivalently the
distance of the binary, is not relevant in these figures un-
less the redshift z becomes significant. In this case the
mass of the binary is the redshifted mass (1+ z)M . Ver-
tical lines in each figure divide a waveform into segments,
where each segment contributes 10% of the total signal
power. In each plot, except for the 10M⊙-binary one, we
show all 9 vertical lines that divide the waveforms into
10 segments. In the 10M⊙-binary plot we omit the early
part of the inspiral phase that accounts for 50% of the
signal power, as it would be too long to show.
The absolute time-scale of a waveform increases lin-
early with total mass M ; equivalently the waveform is
shifted toward lower frequency bands. For a M = 10M⊙
binary, the long inspiral stage generates GWs with fre-
quencies spanning the most sensitive part of the LIGO
band, around 150Hz, while for an M = 100M⊙ binary,
only the merger signal contributes in this band. Thus,
for low-mass binary systems, most of the contribution to
the signal power comes from the long inspiral stage of
the waveform, while for high-mass binary systems most
of the contribution comes from the late inspiral, merger,
and ring-down stages. Understanding quantitatively the
distribution of signal power will let us deduce how many,
and which, GW cycles are significant for the purpose of
data analysis. We need accurate waveforms from either
PN models or NR simulations for at least those signifi-
cant cycles.
From Fig. 4 we conclude that:
• For an M = 10M⊙ binary, the last 25 inspiral cy-
cles, plus the merger and ring-down stages of the
waveform contribute only 50% of the signal power,
and we need 80 cycles (not shown in the figure) of
accurate inspiral waveforms to recover 90% of the
signal power. For an M = 20M⊙ binary, the last
23 cycles, plus the merger and ring-down stages of
the waveform contribute> 90% of the signal power,
and current NR simulations can produce waveforms
of such length;
• For an M = 30M⊙ binary, the last 11 inspiral cy-
cles, plus the merger and ring-down stages of the
waveform contribute > 90% of the signal power,
which means that, for binary systems with to-
tal masses higher than 30M⊙, current NR simula-
tions, e.g., the sixteen cycles obtained in Ref. [38],
can provide long enough waveforms for a matched-
filter search of binary coalescence, as also found in
Ref. [21];
• For an M = 100M⊙ binary, > 90% of the signal
power comes from the last inspiral cycle, merger
and ring-down stages of the waveform, with two
cycles dominating the signal power. It is thus pos-
sible to identify this waveform as a “burst” signal.
Similar analyses can be also done for advanced LIGO and
VIRGO.
C. Comparing hybrid waveforms
We shall now compute FF0s between hybrid wave-
forms. We fix the total mass of the equal-mass binary
in each comparison, i.e., we do not optimize over mass
parameters, but only on phase and time. We use the
mismatch, defined as 1− FF0, to measure the difference
between waveforms and we compute them for LIGO, ad-
vanced LIGO, and VIRGO. Note that by using FF0, we
test the closeness among hybrid waveforms that are gen-
erated from binary systems with the same physical pa-
rameters; in other words, we test whether the waveforms
are accurate enough for the purpose of parameter estima-
tion, rather than for the sole purpose of detecting GWs.
In the language of Ref. [19] we are studying the faithful-
ness of the PN templates 2.
Since at late inspiral stages PN waveforms are partly
replaced by NR waveforms, differences between hybrid
waveforms from two PN models are smaller than those
between pure PN waveforms. In general, the more NR
cycles we use to generate hybrid waveforms, the less the
2 Following Ref. [19], faithful templates are templates that have
large overlaps, say >∼ 96.5%, with the expected signal maximiz-
ing only over the initial phase and time of arrival. By contrast
when the maximization is done also on the binary masses, the
templates are called effectual.
difference is expected to be between these hybrid wave-
forms. This is evident in Figs. 5, 6 where we show mis-
matches between hybrid waveforms for binary systems
with different total masses as a function of the number
of NR cycles n. Specifically, the mismatches are taken
between two hybrid waveforms generated from the same
NR waveform (from the Goddard group, taking the last
n cycles, plus merger and ring-down) and two different
PN waveforms generated with the same masses.
The mismatches are lower for binary systems with
higher total masses, since most of their signal power
is concentrated in the late cycles close to merger (see
Fig. 4). Comparing results between LIGO, advanced
LIGO and VIRGO, we see that for the same waveforms
the mismatches are lowest when evaluated with the LIGO
PSD, and highest when evaluated with the VIRGO PSD.
This is due to the much broader bandwidth of VIRGO,
especially at low frequency: the absolute sensitivity is not
relevant; only the shape of the PSD matters. In VIRGO,
the inspiral part of a hybrid waveform has higher weight-
ing in its contribution to the signal power. As already
observed at the end of Sec. II, we can see also that the
difference between the Epn(3.5) and Tpn(3.5) models is
smaller than that between the Tpn(3) and Tpn(3.5) mod-
Figures 5, 6 show good agreement among hybrid wave-
forms. In Sec. IV, as a further confirmation of what was
found in Refs. [10, 11], we shall see that PN waveforms
from Tpn and Epn models have good agreement with the
inspiral phase of the NR waveforms. Therefore, we argue
that hybrid waveforms are likely to have high accuracy.
In fact, for the late evolution of a compact binary, where
NR waveforms are available, the PN waveforms are close
to the NR waveforms, while for the early evolution of the
binary, where we expect the PN approximations to work
better, the PN waveforms (from Tpn and Epn models)
are close to each other. Based on these observations,
we draw the following conclusions for LIGO, advanced
LIGO, and VIRGO data-analysis:
• For binary systems with total mass higher than
30M⊙, the current NR simulations of equal-mass
binary systems (16 cycles) are long enough to re-
duce mismatches between hybrid waveforms gen-
erated from the three PN models to below 0.5%.
Since these FFs are achieved without optimizing
the binary parameters, we conclude that for these
high-mass binary systems, the small difference be-
tween hybrid waveforms indicates low systematic
error in parameter estimation, i.e., hybrid wave-
forms are faithful [19].
• For binary systems with total mass around
10–20M⊙, 16 cycles of NR waveforms can reduce
the mismatch to below 3%, which is usually set as
the maximum tolerance for data-analysis purpose
(corresponding to ∼ 10% loss in event rate). By a
crude extrapolation of our results, we estimate that
with 30 NR waveform cycles, the mismatch might
4 6 8 10 12 14
Number of NR waveform cycles
Restricted PN
Rescaled PN
4 6 8 10 12 14
Number of NR waveform cycles
Quartic-fit frequency
Oscillatory frequency
FIG. 5: We show the mismatch between hybrid waveforms as a function of the number of NR waveform cycles used to generate
the hybrid waveforms. The LIGO PSD is used to evaluate the mismatches. In the left panel, we compare the Epn(3.5) and
Tpn(3.5) models. In the right panel, we compare the Tpn(3) and Tpn(3.5) models. From top to bottom, the four curves
correspond to four equal-mass binary systems, with total masses 10M⊙, 20M⊙, 30M⊙, and 40M⊙. The dots show mismatches
taken between hybrid waveforms that are generated with different methods. In the left panel, we adjust the amplitude of
restricted PN waveforms, such that they connect smoothly in amplitude to NR waveforms. In the right panel, to set the
frequency of PN waveforms at the joining point, we use the original orbital frequency, instead of the quartic fitted one. (See
Sec. IIIA for the discussion on amplitude scaling and frequency fitting).
4 6 8 10 12 14
Number of NR waveform cycles
Advanced LIGO PSD
VIRGO PSD
4 6 8 10 12 14
Number of NR waveform cycles
Advanced LIGO PSD
VIRGO PSD
FIG. 6: Mismatch between hybrid waveforms as a function of the number of NR waveform cycles used to generate the hybrid
waveforms. Following the settings of Fig. 5, we show comparisons between Epn(3.5) and Tpn(3.5), and Tpn(3) and Tpn(3.5)
models in the left and the right panels, respectively. The solid and dashed sets of curves are generated using the PSDs of
advanced LIGO and VIRGO. In each set, from top to bottom, the three curves correspond to three equal-mass binary systems,
with total masses 20M⊙, 30M⊙, and 40M⊙.
be reduced to below 1%.
• For binary systems with total mass lower than
10M⊙, the difference between the Tpn(3) and
Tpn(3.5) models is substantial for Advanced LIGO
and VIRGO. Their mismatch can be > 4% and
> 6% respectively (not shown in the figure). In this
mass range, pursuing more NR waveform cycles in
late inspiral phase does not help much, since the
signal power is accumulated slowly over hundreds
of GW cycles across the detector band. Neverthe-
less, here we give mismatches for FF0s which are
not optimized over binary masses. For the purpose
of detection only, optimization over binary param-
eters leads to low enough mismatches (see also the
end of Sec. II). In the language of Ref. [19] hybrid
waveforms for total mass lower than 10M⊙ are ef-
fectual but not faithful.
IV. MATCHING NUMERICAL WAVEFORMS
WITH POST-NEWTONIAN TEMPLATES
In this section, we compare the complete inspiral,
merger and ring-down waveforms of coalescing compact
binary systems generated from NR simulations with their
best-match PN template waveforms. We also compare
hybrid waveforms with PN template waveforms for lower
total masses, focusing on the late inspiral phase pro-
vided by the NR waveforms. We test seven families of
PN templates that either have been used in searches
for GWs in LIGO (see e.g., Refs. [40, 41]), or are
promising candidates for ongoing and future searches
with ground-based detectors. We evaluate the perfor-
mance of PN templates by computing the FFs maxi-
mized on phase, time and binary parameters. As we
shall see, for the hybrid waveforms of binary systems
with total mass M ≤ 30M⊙, both the time-domain fam-
ilies Tpn(3.5) and Epn(3.5), which includes a superposi-
tion of three ring-down modes, perform well, confirming
what found in Refs. [10, 11]. The standard stationary-
phase-approximated (SPA) template family in the fre-
quency domain has high FFs only for binary systems with
M < 20M⊙. After investigating in detail the GW phase
in frequency domain, and having understood why it hap-
pens (see Sec. IVB2), we introduce two modified SPA
template families (defined in Sec. IVB 2) for binary sys-
tems with total mass M ≥ 30M⊙. Overall, for masses
M ≥ 30M⊙, the Epn(3.5) template family in the time do-
main and the two modified SPA template families in the
frequency domain exhibit the best-match performances.
A. Numerical waveforms and post-Newtonian
templates
For binary systems with total mass M ≥ 30M⊙, the
last 8–16 cycles contribute more than 80–90% of the sig-
nal power, thus in this case we use only the NR wave-
forms. By contrast, for binary systems with total mass
10 ≤ M ≤ 30M⊙, for which the merger and ring-down
phases of the waveforms contribute only ∼ 1–10%, we
use the hybrid waveforms, generated by stitching Tpn
waveforms to the Goddard NR waveforms.
We want to emphasize that FFs computed for differ-
ent target numerical waveforms can not directly be com-
pared with each other. For instance, the Goddard wave-
form is longer than the Pretorius waveform, and the FFs
are sometime slightly lower using the Goddard waveform.
This is a completely artificial effect, due to the fact that
it is much easier to tune the template parameters and
obtain a large FF with a shorter target waveform than a
longer one.
We consider seven PN template families. The two
time-domain families introduced in Sec. II are:
• Tpn(3.5) [30, 31]: The inspiral Taylor model.
• Epn(3.5) [10, 13, 14, 16, 19]: The EOB model which
includes a superposition of three quasi-normal modes
(QNMs) of the final BH. These are labeled by three in-
tegers (l,m, n) [42]: the least damped QNM (2, 2, 0) and
two overtones (2, 2, 1) and (2, 2, 2). The ring-down wave-
form is given as:
hQNM(t) =
−(t−tend)/τ22n cos [ω22n(t− tend) + φn] ,
where ωlmn and τlmn are the frequency and decay time of
the QNM (l,m, n), determined by the mass Mf and spin
af of the final BH. The quantities An and φn in Eq. (7)
are the amplitude and phase of the QNM (2, 2, n). They
are obtained by imposing the continuity of h+ and h×,
and their first and second time derivatives, at the time of
matching tmatch. Besides the mass parameters, our Epn
model contains three other physical parameters: ǫt, ǫM
and ǫJ . The parameter ǫt takes into account possible dif-
ferences between the time tend at which the EOB models
end and the time tmatch at which the matching to ring-
down is done. More explicitly, we set tmatch = (1+ǫt)tend,
and if ǫt > 0, we extrapolate the EOB evolution, and set
an upper limit for the ǫt search where the extrapolation
fails. The parameters ǫM and ǫJ describe possible differ-
ences between the values of the mass Mend ≡ Eend and
angular momentum âend ≡ Jend/M
end at the end of the
EOB inspiral and the final BH mass and angular mo-
mentum. (The end of the EOB inspiral occurs around
the EOB light-ring.) The differences are due to the fact
that the system has yet to release energy and angular
momentum during the merger and ring-down phase be-
fore settling down to the stationary BH solution. If the
total binary mass and angular momentum at the end of
the EOB inspiral are Mend and Jend, we set the total
mass and angular momentum of the final stationary BH
to be Mf = (1 − ǫM )Mend and Jf = (1 − ǫJ)Jend, and
use af ≡ Jf/Mf to compute ωlmn and τlmn. We consider
the current Epn model with three parameters ǫt, ǫM and
ǫJ , as a first attempt to build a physical EOB model for
matching coherently the inspiral, merger and ring-down
phases. Since the ǫ-parameters are related to physical
quantities, e.g., the loss of energy during ring-down, they
are functions of the initial physical parameters of the bi-
nary, such as masses, spins, etc. In the near future we
expect to be able to fix the ǫ-values by comparing NR and
(improved) EOB waveforms for a large range of binary
parameters.
We also consider five frequency-domain models, in
which two (modified SPA models) are introduced later
in Sec. IVB 2, and three are introduced here:
• SPAc(3.5) [33]: SPAc PN model with an appropriate
cutoff frequency fcut [30, 31];
(5 + 5)M⊙ (10 + 10)M⊙ (15 + 15)M⊙
Signal Power (%) (30, 0.2) (80, 2) (85, 10)
〈hNR−hybr, hTpn(3.5)〉 0.9875 0.9527 0.8975
(M/M⊙, η) (10.18, 0.2422) (19.97, 0.2500) (29.60, 0.2499)
Mωorb 0.1262 0.1287 0.1287
〈hNR−hybr, hEpn(3.5)〉 0.9836 0.9522 0.9618
(M/M⊙, η) (10.15, 0.2435) (19.90, 0.2500) (29.49,0.2488)
(ǫt, ǫM , ǫJ )(%) (-0.02, 12.19, 30.87) (-0.02, 75.03, 95.00) (0.05, 2.38, 92.06)
Mωorb 0.1346 0.1345 0.1345
〈hNR−hybr, hSPAc(3.5)〉 0.9690 0.9290 0.8355
(M/M⊙, η) (10.16, 0.2432) (19.93, 0.2498) (29.08, 0.2500)
(fcut/Hz) 1566.8 263.9 529.6
TABLE I: FFs between hybrid waveforms [Tpn(3.5) waveform stitched to the Goddard waveform] and PN templates. In the first
row, the two numbers in parentheses are the percentages of the signal-power contribution from the 16 inspiraling NR cycles and
the NR merger/ring-down cycles. (The separation between inspiral and merger/ring-down is obtained using the EOB approach
as a guide, i.e., we match the Epn(3.5) model and use the EOB light-ring position as the beginning of the merger phase.) In
the PN-template rows, the first number in each block is the FF, and the numbers in parentheses are template parameters that
achieve this FF. The last number in each block of the Tpn(3.5) and Epn(3.5) models is the ending orbital frequency of the
best-match template. For the Epn model, the ending frequency is computed at the point of matching with the ring-down phase,
around the EOB light ring.
• BCV [31]: BCV model with an amplitude correction
term (1−αf2/3) and an appropriate cutoff frequency fcut.
• BCVimpr [31]: Improved BCV model with an am-
plitude correction term (1 − αf1/2) and an appropri-
ate cutoff frequency fcut. We include this improved
BCV model because Ref. [10] found a deviation of the
Fourier-transform amplitude from the Newtonian predic-
tion f−7/6 during the merger and ring-down phases (see
Fig. 22 of Ref. [10]). Here we shall assume n = −2/3
in the fn power law to get the (1 − αf1/2) form of the
amplitude correction. While it was found [10] that the
value of n is close to −2/3 for the l = 2,m = 2 wave-
form, this value varies slightly if other multiple moments
are included and if binary systems with different mass ra-
tios are considered. Finally, the α parameter is expected
to be negative, but in our actual search it can take both
positive and negative values.
B. Discussion of fitting-factor results
In Table I, we list the FFs for hybrid target waveforms
and three PN template families: Tpn(3.5), Epn(3.5),
and SPAc(3.5), together with the template parameters at
which the best match is obtained. As shown in the first
row, in this relatively low-mass range, i.e. 10M⊙ < M <
30M⊙, the merger/ring-down phases of the waveforms
contribute only a small fraction of the total signal power,
while the last 16 inspiraling cycles of the NR waveform
contribute a significant fraction. Therefore, confirming
recent claims by Refs. [10, 11], we can conclude that
the PN template families Tpn(3.5) and Epn(3.5) have
good agreement with the inspiraling NR waveforms. The
Tpn(3.5) model gives a low FF for M = 30M⊙ because
0.0 0.2 0.4 0.6 0.8 1.0
af /Mf
FIG. 7: Frequencies and decay times of the least damped
QNM 220, and two overtones 221 and 222. The scales of the
frequency and the decay time are listed on the left and right
sides of the plot, respectively.
for these higher masses the merger/ring-down phases,
which the Tpn model does not include, start contribut-
ing to the signal power. Note that both time-domain
templates give fairly good estimates of the mass parame-
ters. The SPAc(3.5) template family gives FFs that drop
substantially when the total binary mass increases from
10M⊙ to 30M⊙, indicating that this template family can
only match the early, less relativistic inspiral phase of
the hybrid waveforms. Nevertheless, it turns out that by
slightly modifying the SPA waveform we can match the
NR waveforms with high FFs (see Sec. IVB2).
In Table II, we list the FFs for full NR waveforms
and five PN template families: Epn(3.5), SPAextc (3.5),
SPAYc (4), BCV, and BCVimpr, together with the tem-
(15 + 15)M⊙ (20 + 20)M⊙ (30 + 30)M⊙ (50 + 50)M⊙
〈hNR−Pretorius, hEpn(3.5)〉 0.9616 0.9599 0.9602 0.9787
(M/M⊙, η) (27.93, 0.2384) (35.77, 0.2426) (52.27, 0.2370) (96.60, 0.2386)
(ǫt, ǫM , ǫJ )(%) (-0.08, 0.63, 99.70) (-0.03, 0.48, 94.38) (-0.12, 0.00, 64.14) (0.04, 0.01, 73.01)
〈hNR−Pretorius, hSPA
c (3.5)〉 0.9712 0.9802 0.9821 0.9722
(M/M⊙, η) (19.14, 08037) (24.92, 0.9097) (36.75, 0.9933) (58.06, 0.9986)
(fcut/Hz) (589.6) (476.9) (318.9) (195.9)
〈hNR−Pretorius, hSPA
c (4)〉 0.9736 0.9824 0.9874 0.9851
(M/M⊙, η) (29.08, 0.2460) (38.63, 0.2461) (57.58, 0.2441) (96.55, 0.2457)
(fcut/Hz) (666.5) (501.2) (332.5) (199.4)
〈hNR−Pretorius, hBCV〉 0.9726 0.9807 0.9788 0.9662
(ψ0/10
4, ψ1/10
2) (2.101, 1.655) (1.178, 1.744) (0.342, 2.385) (-0.092, 3.129)
(102α, fcut/Hz) (-1.081, 605.5) (-0.834, 461.7) (0.162, 320.4) (1.438, 204.3)
〈hNR−Pretorius, hBCVimpr〉 0.9727 0.9807 0.9820 0.9803
(ψ0/10
4, ψ1/10
2) (2.377, 0.930) (1.167, 1.762) (0.431, 2.077) (-0.109, 3.158)
(102α, fcut/Hz) (-3.398, 571.9) (-2.648, 458.3) (-1.196, 319.1) (-3.233, 196.0)
(15 + 15)M⊙ (20 + 20)M⊙ (30 + 30)M⊙ (50 + 50)M⊙
〈hNR−Goddard, hEpn(3.5)〉 0.9805 0.9720 0.9692 0.9671
(M/M⊙, η) (29.25, 0.2435) (38.27, 0.2422) (56.66, 0.2381) (83.52, 0.2233)
(ǫt, ǫM , ǫJ )(%) (0.05, 0.03, 99.90) (0.05, 0.27, 99.17) (0.09, 0.01, 54.56) (0.10, 1.71, 79.75)
〈hNR−Goddard, hSPA
(3.5)〉 0.9794 0.9785 0.9778 0.9693
(M/M⊙, η) (21.41, 0.5708) (27.27, 0.6695) (37.67, 0.9911) (60.90, 0.9947)
(fcut/Hz) (552.7) (444.4) (318.5) (191.7)
〈hNR−Goddard, hSPA
(4)〉 0.9898 0.9905 0.9885 0.9835
(M/M⊙, η) (30.28, 0.2456) (40.23, 0.2477) (60.54, 0.2455) (100.00, 0.2462)
(fcut/Hz) (674.6) (506.6) (330.5) (195.0)
〈hNR−Goddard, hBCV〉 0.9707 0.9710 0.9722 0.9692
(ψ0/10
4, ψ1/10
2) (3.056, -1.385) (1.650, -0.091) (0.561, 1.404) (-0.113, 3.113)
(102α, fcut/Hz) (0.805, 458.3) (0.559, 412.6) (0.218, 309.2) (1.063, 198.7)
〈hNR−Goddard, hBCVimpr〉 0.9763 0.9768 0.9782 0.9803
(ψ0/10
4, ψ1/10
2) (2.867, -0.600) (1.514, 0.448) (0.555, 1.425) (-0.165, 3.373)
(102α, fcut/Hz) (0.193, 578.0) (-1.797, 441.1) (-4.472, 308.1) (-4.467, 193.4)
TABLE II: FFs between NR waveforms and PN templates which include merger and ring-down phases. The upper table uses
Pretorius’ waveform, and the lower table uses Goddard’s high-resolution long waveform. The first number in each block is the
FF, and numbers in parentheses are template parameters that achieve this FF.
plate parameters at which the best match is obtained.
The SPAextc (3.5) and SPA
c (4) families are modified ver-
sions of the SPA family, defined in Sec. IVB 2.
We shall investigate these results in more detail in the
following sections.
1. Effective-one-body template performances
The Epn model is the only available time-domain
model that explicitly includes ring-down waveforms. It
achieves high FFs ≥ 0.96 for all target waveforms, con-
firming the necessity of including ring-down modes and
proving that the inclusion of three QNMs with three
tuning parameters ǫt, ǫM and ǫJ is sufficient for detec-
tion. As we see in Table II, the values of the tuning
parameters ǫM and ǫJ , where the FFs are achieved, are
different from their physical values. For reference, the
Goddard numerical simulation predicts Mf ≃ 0.95M
and âf ≡ Jf/M
f ≃ 0.7 [38], and Epn(3.5) predicts
Mend = 0.967 and âend ≡ Jend/M
end = 0.796, so the two
tuning parameters should be ǫM ≃ 1.75% and ǫJ ≃ 11%.
In our search, e.g., for M = 30M⊙, ǫJ tends to be tuned
to its lowest possible value and ǫt tends to take its high-
est possible value, indicating that pushing the end of the
Epn(3.5) inspiral to a later time gives higher FFs.
Since the parameters ǫM and ǫJ depend on the QNM
frequency and decay time, we show in Fig. 7 how ωlmn
and τlmn vary as functions of af [42] for the three modes
used in the Epn(3.5) model. The frequencies ωlmn of the
three modes are not really different, and grow monotoni-
cally with increasing af . The decay times τlmn, although
different for the three modes, also grow monotonically
with increasing af . Thus, the huge loss of angular mo-
mentum ǫJ , or equivalently the small final BH spin re-
quired in the Epn(3.5) model to achieve high FFs, indi-
cates that low ring-down frequencies and/or short decay
times are needed for this model to match the numerical
0.10 0.12 0.14 0.16 0.18
Time (s)
-0.02
NR frequency & waveform
Epn frequency & waveform
(15+15)
0.40 0.45 0.50 0.55 0.60
Time (s)
(50+50)
FIG. 8: Frequency evolution of waveforms from the Epn(3.5) model, and the NR simulations of the Goddard group. In the left
and right panels, we show frequency evolutions for two equal-mass binary systems with total mass 30M⊙ and 100M⊙. In each
panel, there are two nearly monotonic curves and two oscillatory curves, where the former are frequency evolutions and the
latter are binary coalescence waveforms. The solid curves (blue) are from the NR simulations, while the dashed curves (red)
are from the Epn(3.5) model. The vertical line in each plot shows the position where the three-QNM ring-down waveform is
attached to the EOB waveform.
merger and ring-down waveforms.
In Fig. 8, we show Goddard NR and Epn(3.5) wave-
forms, as well as their frequency evolutions, for two
equal-mass binary systems with total masses 30M⊙ and
100M⊙. In the low-mass case, i.e., M = 30M⊙, since the
inspiral part contributes most of the SNR, the Epn(3.5)
model fits the frequency and phase evolution of the NR
inspiral well, with the drawback that at the joining point
the EOB frequency is substantially higher than that of
the NR waveform. Then, in order to fit the early ring-
down waveform which has higher amplitude, the tuning
parameters have to take values in Table II such that the
ring-down frequency is small enough to get close to the
NR frequency during early ring-down stage, as indicated
in Fig. 8. The late ring-down waveform does not con-
tribute much to the SNR, and thus it is not too sur-
prising that waveform optimizing the FF does not ade-
quately represent this part of the NR waveform. In the
higher mass case,M = 100M⊙, the Epn(3.5) model gives
a much better, though not perfect, match to the merger
and ring-down phases of the NR waveform, at the ex-
pense of misrepresenting the early inspiral part. Again,
this is not unexpected considering that in this mass range
the merger and ring-down waveforms dominate the con-
tribution to the SNR.
Comparing the two cases discussed above, we can see
that with the current procedure of matching the inspiral
and ring-down waveforms in the EOB approach it is not
possible to obtain a perfect match with the entire NR
waveform. However, due to the limited detector sensi-
tivity bandwidth, the FFs are high enough for detection.
The large systematic error in estimating the physical pa-
rameters will be overcome by improving the EOB match-
ing procedure during the inspiral part, and also by fixing
the ǫ-parameters to physical values obtained by compar-
ison with numerical simulations.
Finally, in Figs. 9, 10 we show the frequency-domain
amplitude and phase of the NR and EOB waveforms.
Quite interestingly, we notice that the inclusion of three
ring-down modes reproduce rather well the bump in the
NR frequency-domain amplitude. The EOB frequency-
domain phase also matches the NR one very well.
2. Stationary-phase-approximated template performances
Figures 9, 10 also show the frequency-domain phases
and amplitudes for the best-match SPAc(3.5) waveforms.
We see that at high frequency the NR and SPAc(3.5)
phases rise with different slopes 3. Based on this obser-
vation we introduce two modified SPA models:
• SPAextc (3.5): SPAc PN model with unphysical values
of η and an appropriate cutoff frequency fcut. The range
of the symmetric mass-ratio η = m1m2/(m1+m2)
2 is ex-
tended from its physical range 0 ∼ 0.25 to the unphysical
3 By looking in detail at the PN terms in the SPAc(3.5) phase, we
find that the difference in slope is largely due to the logarithmic
term at 2.5PN order.
100 200 500
NR phase
SPAc(3.5)
ext(3.5)
Epn(3.5)
100 200 500
Frequency (Hz)
NR amplitude
Epn(3.5)
FIG. 9: For M = 30M⊙ equal-mass binary systems, we com-
pare the phase and amplitude of the frequency-domain wave-
forms from the SPAc models and NR simulation (Goddard
group). We also show the amplitude of the waveform from
the Epn(3.5) model.
range 0 ∼ 1.
• SPAYc (4): SPAc PN model with an ad hoc 4PN order
term in the phase, and an appropriate cutoff frequency
fcut. The phase of the SPA model is known up to the
3.5PN order (see, e.g., Eq. (3.3) of Ref. [33]):
ψ(f) = 2πft0 − φ0 −
128ηv5
k , (8)
where v = (πMf)1/3. The PN coefficients αks, k =
0, . . . , N , (with N = 7 at 3.5PN order) are given by Eqs.
(3.4a), (3.4h) of Ref. [33]. We add the following term at
4PN order:
α8 = Y log v , (9)
where Y is a parameter which we fix by imposing high
matching performances with NR waveforms. Note that a
constant term in α8 only adds a 4PN order term that is
linear in f , which can be absorbed into the 2πft0 term.
Thus, to obtain a nontrivial effect, we need to introduce
a logarithmic term. The coefficient Y could in princi-
ple depend on η. We determine Y by optimizing the
FFs of equal and unequal masses. We find that in the
equal-mass case Y does not depend significantly on the
binary total mass and is given by Y = 3923. The latter
is also close to the best match value obtained for unequal
masses. More specifically, it is within 4.5% for binary sys-
tems of mass ratiom2/m1 = 2. To further explore the de-
pendence of Y on η, we need a larger sample of waveforms
50 100 200
NR phase
SPAc(3.5)
ext(3.5)
Epn(3.5)
M=100
50 100 200
Frequency (Hz)
NR amplitude
Epn(3.5)
FIG. 10: For M = 100M⊙ equal-mass binary systems, we
compare the phase and amplitude of the frequency-domain
waveforms from the SPAc models and NR simulation (God-
dard group). We also show the amplitude of the waveform
from the Epn(3.5) model.
for unequal-mass binary systems 4. As seen in Table II,
the two modified SPAc template families have FF> 0.97
(except for one 0.9693) for all target waveforms, even
though no explicit merger or ring-down phases are in-
cluded in the waveform. The SPAYc (4) model provides
also a really good estimation of parameters.
In Fig. 11 we plot Goddard NR and SPAYc (4) wave-
forms for two equal-mass binary systems with total
massesM = 30M⊙ andM = 100M⊙. We can clearly see
tail-like ring-down waveforms at the end of the SPAYc (4)
waveforms, which result from the inverse Fourier trans-
form of frequency domain waveforms that have been cut
at f = fcut. This well-known feature is called the Gibbs
phenomenon. At first glance, it may appear surpris-
ing that the often inconvenient Gibbs phenomenon [44]
can provide reasonable ring-down waveforms in the time
domain. However, by looking at the spectra of these
waveforms in the frequency domain (see the amplitudes
in Figs. 9 & 10), we see that the SPAYc (4) cuts off at
the frequency fcut (obtained from the optimized FF)
where the NR spectra also start to drop. Thus, even
though the frequency-domain SPAc waveforms are dis-
continuous, while the frequency-domain NR waveforms
are continuous (being combinations of Lorentzians), the
SPAc time-domain waveforms contain tails with frequen-
cies and decay rates similar to the NR ring-down modes.
We expect that the values of the cutoff frequency fcut
at which the FFs are maximized are well-determined by
4 Note that the auxiliary phase introduced in Eq. (239) of Ref. [43]
also gives rise to a term in the SPA phase of the kind f log v,
except an order of magnitude smaller than Y .
the highest frequency of the NR waveforms, i.e. by the
frequency of the fundamental QNM. In the next section,
we shall show quantitative results to confirm this guess.
3. Buonanno-Chen-Vallisneri template performances
In Table II we see that the BCV and BCVimpr fam-
ilies give almost the same FFs for relatively low-mass
binary systems (M = 30, 40M⊙), while the BCVimpr
family gives slightly better FFs for higher mass binary
systems (M = 60, 100M⊙). For higher-mass binary sys-
tems, we find that the α parameter takes negative values
with reasonable magnitude. This is because the ampli-
tude of the NR waveforms in the frequency domain de-
viates from the f−7/6 power law only near the merger,
which lasts for about one GW cycle. This merger cycle
is important only when the total mass of the binary is
high enough (see Fig. 4). [See also Ref. [46] where similar
tests have been done.]
The BCV and BCVimpr template families give FFs
nearly as high as those given by the SPAYc (4) family,
but the latter has the advantage of being parametrized
directly in terms of the physical binary parameters, and
it gives fairly small systematic errors.
C. Frequency-domain templates for inspiral,
merger and ring-down
In this section, we extend our comparisons between
the SPAc families and NR waveforms to higher total-
mass binary systems (40M⊙ to 120M⊙) and to unequal-
mass binary systems with mass-ratios m2/m1 = 1.5 and
2. The numerical simulations for unequal-mass binary
systems are from the Goddard group. They last for ≃
373M and ≃ 430M , respectively, and the NR waveforms
have ≃ 4 cycles before the merger.
In Figs. 12 and 13 we show the FFs for SPAextc (3.5)
and SPAYc (4) templates, and the values of fcut that
achieved these FFs 5. For all mass combinations (ex-
cept for M = 40M⊙ for artificial reasons) the FFs of
SPAc(3.5) templates are higher than 0.96, and the FFs
of SPAYc (4) templates are higher than 0.97, confirming
that both families of templates can be used to search for
GWs from coalescing binary systems with equal-masses
as large as 120M⊙ and mass ratios m2/m1 = 2 and 1.5.
Figure 13 shows that all the fcut values from our searches
5 Note that because of the short NR waveforms for unequal-mass
binary systems, we need to search over the starting frequency of
templates with a coarse grid, and this causes some oscillations in
our results. The oscillations are artificial and will be smoothed
out in real searches. For instance, the drop of FFs at 40M⊙
for unequal-mass binary systems happens because the NR wave-
forms are too short and begin right in the most sensitive band
of LIGO.
are within 10% larger than the frequency of the funda-
mental QNM ω220 of an equal-mass binary. We have
checked that if we fix fcut = 1.07ω220/2π, the FFs drop
by less than 1%.
In Fig. 15, we show the same information as in Fig. 7,
except that here we draw ωlnm and τlnm as functions of
the mass-ratio η of a nonspinning binary. We compute
the spin of the final BH in units of the mass of the fi-
nal BH using the quadratic fit given by Eq. (3.17a) of
Ref. [47]:
≃ 3.352η − 2.461η2 . (10)
As Fig. 15 shows, ω220 does not change much, confirming
the insensitivity of the fcut on η.
However, in real searches we might request that the
template family have some deviations from the wave-
forms predicted by NR. For example, a conservative tem-
plate bank might cover a region of fcut ranging from the
Schwarzschild innermost stable circular orbit (ISCO) fre-
quency, or the innermost circular orbit (ICO) frequency
determined by the 3PN conservative dynamics, up to a
value slightly higher than the frequency of the fundamen-
tal QNM. The number of templates required to cover
the fcut dimension depends on the binary masses. We
find that to cover the fcut dimension from the 3PN ICO
frequency to the fundamental QNM frequency with an
SPAextc (3.5) template bank, imposing a mismatch < 0.03
between neighboring templates, we need only two (∼20)
templates if M = 30M⊙ (M = 100M⊙) and η = 0.25.
In the latter case, the match between templates is more
sensitive to fcut since most signal power comes from
the last two cycles, sweeping through a large frequency
range, right in LIGO’s most sensitive band. The num-
ber of templates directly affects the computational power
needed, and the false-alarm rate. Further investigations
are needed in order to determine the most efficient way
to search over the fcut dimension.
For the purpose of parameter estimation, Fig. 14 shows
that the SPAYc (4) templates are rather faithful, giving
reasonable estimates of the chirp mass: systematic errors
less than about 8% in absolute value for binary systems
with M = 40M⊙ up to M = 120M⊙. A difference of
≃ 8% may seem large, but the SPAYc (4) templates are
not exactly physical, and more importantly, for large-
mass binary systems, most of the information on the
chirp mass comes only from the last cycle of inspiral. We
notice that when the total binary mass is higher than
120M⊙, the FFs are relatively high (from 0.93 to 0.97),
and the estimates of the chirp mass are still good (within
10%). However, for binary systems with such high to-
tal masses, the ring-down waveform dominates the SNR,
and the SPAYc (4) template family becomes purely phe-
nomenological. A direct ring-down search might be more
efficient.
All results for unequal-mass binary systems are ob-
tained using the C22 component of Ψ4 [10], which is
the leading order quadrupole term contributing to the
0.06 0.08 0.10 0.12 0.14 0.16 0.18
Time (s)
NR waveform
(15+15)
0.45 0.50 0.55 0.60
Times (s)
(50+50)
FIG. 11: Binary coalescence waveforms from the SPAYc (4) model, and the NR simulations of the Goddard group. In the left
and right panels we show waveforms for two equal-mass binary systems with total mass 30M⊙ and 100M⊙. The solid lines
show the waveforms from the NR simulation, and the dashed lines give the best-matching waveforms from the SPAYc (4) model.
40 60 80 100 120
Total mass ( )
ext(3.5): equal-mass
ext(3.5): mass-ratio 2:1
ext(3.5): mass-ratio 3:2
(4): equal-mass
(4): mass-ratio 2:1
(4): mass-ratio 3:2
FIG. 12: FFs as functions of the total binary mass. The FFs
are computed between either the SPAextc (3.5) or the SPA
c (4)
templates and the NR waveforms for equal-mass and unequal-
mass binary systems.
GW radiation. For unequal-mass binary systems, higher-
order multipoles can also be important, and we need to
test the performance of the template family directly using
Ψ4. For Ψ4 extracted in the direction perpendicular to
the binary orbit, we verified that higher-order multipoles
do not appreciably change the FFs.
A natural way of improving the SPAc models would be
to replace the discontinuous frequency cut with a linear
combination of Lorentzians. We show here a first attempt
40 60 80 100 120
Total mass ( )
ext(3.5): equal-mass
ext(3.5): mass-ratio 2:1
ext(3.5): mass-ratio 3:2
(4): equal-mass
(4): mass-ratio 2:1
(4): mass-ratio 3:2
Frequency of the fundamental QNM
FIG. 13: Cutoff frequencies as functions of the total binary
mass. We show the best-match fcut for SPA
c (3.5) and
SPAYc (4) templates of Fig. 12. The solid black curve is the
fundamental QNM frequency ω220/2π. The frequencies are in
units of Hz.
at doing so. The Lorentzian L is obtained as a Fourier
transform of a damped sinusoid, e.g., for the fundamental
QNM we have
ei2πft
e±iω220t−|t|/τ220
2/τ220
1/τ2220 + (2πf ± ω220)
≡ 2L±220(f) (11)
40 60 80 100 120
Total mass ( )
-0.08
-0.06
-0.04
-0.02
(4): equal-mass
(4): mass-ratio 2:1
(4): mass-ratio 3:2
FIG. 14: Systematic errors of the chirp mass as functions
of the total binary mass when SPAYc (4) templates are used.
We show errors of the chirp masses that optimize the FFs of
Fig. 12.
0.00 0.05 0.10 0.15 0.20 0.25
FIG. 15: Frequencies and decay times of the least damped
QNM 220, and two overtones 221 and 222. The scales of the
frequency and the decay time are listed on the left and right
sides of the plot, respectively.
and the (inverse) Fourier transform of Eq. (7) reads
h̃QNM(f) = An
L+22n(f) e
iφn + L−22n(f) e
Restricting to positive frequencies we only keep the
L−22n(f) terms. In the frequency domain we attach the
fundamental mode continuously to the SPAYc (4) wave-
form at the ring-down frequency ω220 by tuning the am-
plitude and phase A0 and φ0. We denote this model
SPAL1 (note that we also need to introduce the mass-
parameter of the final BH as a scale for ω220 and τ220).
Similarly, we define the SPAL3 model where all three
QNMs are combined. With the three amplitudes and
phases as parameters, this model is similar to the spin-
BCV template family [30] and we can optimize automat-
ically over the 6 parameters. As an example, we compute
the FFs between the SPAL1(4) or SPAL3(4) and the NR
waveform of an equal-mass M = 100M⊙ binary. Using
the LIGO PSD, we obtain 0.9703 and 0.9817, respec-
tively. Those FFs are comparable to the FFs obtained
with the simpler SPAc model, shown in Fig. 12. It is
known that adding more parameters increases the FFs
but also increases the false-alarm probability. By fur-
ther investigation and comparison with NR waveforms
our goal is to express the phase and amplitude parame-
ters of the Lorentzian in terms of the physical binary pa-
rameters, relating them to the amplitudes and phases of
the QNMs and the physics of the merger. Those parame-
ters are somewhat similar to the ǫ-parameters introduced
above for the EOB model when modeling the merger and
ring-down phases.
We wish to emphasize that the results we presented in
this section are preliminary, in the sense that we consid-
ered only a few mass combinations and the NR waveforms
of unequal-mass binary systems are quite short. Never-
theless, these results are interesting enough to propose a
systematic study of the efficiency of these template fam-
ilies through Monte Carlo simulations in real data.
V. CONCLUSIONS
In this paper we compared NR and analytic waveforms
emitted by nonspinning binary systems, trying to under-
stand the performance of PN template families developed
during the last ten years and currently used for the search
for GWs with ground-based detectors, suggesting possi-
ble improvements.
We first computed FF0s (maximized only on time
and phase) between PN template families which best
match NR waveforms [10, 11], i.e., Tpn(3), Tpn(3.5) and
Epn(3.5). We showed how the drop in FF0s is not simply
determined by the accumulated phase difference between
waveforms, but also depends on the detector’s PSD and
the binary mass. Thus, waveforms which differ even by
one GW cycle can have FF0 ∼ 0.97, depending on the
binary masses (see Fig. 1).
We then showed that the NR waveforms from the
high-resolution and medium-resolution simulations of the
Goddard group are close to each other (FF0 around 0.99,
see Fig. 2). We also estimated that the FF0 between
high-resolution and exact NR waveforms is even higher,
based on the numerical convergence rates of the Goddard
simulations.
Second, by stitching PN waveforms to NR waveforms
we built hybrid waveforms, and computed FF0s (max-
imized only on time and phase) between hybrid wave-
forms constructed with different PN models, notably
Tpn(3), Tpn(3.5) and Epn(3.5) models. We found that
for LIGO’s detectors and equal-mass binary systems with
total mass M > 30M⊙, the last 11 GW cycles plus
merger and ring-down phases contribute > 90% of the
signal power. This information can be used to set the
length of NR simulations.
The FF0s between hybrid waveforms are summarized
in Figs. 5, 6. We found that for LIGO’s detectors and
binary systems with total mass higher than 10M⊙, the
current NR simulations for equal-mass binary systems
are long enough to reduce the differences between hybrid
waveforms built with the PN models Tpn(3), Tpn(3.5)
and Epn(3.5) to the level of < 3% mismatch. For GW
detectors with broader bandwidth like advanced LIGO
and VIRGO, longer NR simulations will be needed if the
total binary masses M < 10M⊙. With the current avail-
able length of numerical simulations, it is hard to esti-
mate from the FFs between hybrid waveforms how long
the simulations should be. Nevertheless, from our study
of the distribution of signal power, we estimate that for
M < 10M⊙ binary systems, at least ∼ 80 NR inspiraling
cycles before merger are needed.
Finally, we evaluated FFs (maximized on binary
masses, initial time and phase) between full NR (or hy-
brid waveforms, depending on the total binary mass)
and several time and frequency domain PN template
families. For time-domain PN templates and binary
masses 10M⊙ < M < 20M⊙, for which the merger/ring-
down phases do not contribute significantly to the to-
tal detector signal power, we confirm results obtained in
Refs. [10, 11], notably that Tpn(3.5) and Epn(3.5) mod-
els have high FFs with good parameter estimation, i.e.,
they are faithful. We found that the frequency-domain
SPA family has high FFs only for binary systems with
M < 20M⊙, for which most of the signal power comes
from the early stages of inspiral. Furthermore, we found
that it is possible to improve the SPA family by either
extending it to unphysical regions of the parameter space
(as done with BCV templates) or by introducing an ad
hoc 4PN-order constant coefficient in the phase. Both
modified SPA families achieve high FFs for high-mass
binary systems with total masses 30M⊙ < M < 120M⊙.
For time-domain PN templates and binary masses
M >∼ 30M⊙, we found that if a superposition of ring-
down modes is attached to the inspiral waveform, as nat-
urally done in the EOB model, the FFs can increase from
∼ 0.8 to > 0.9. We tested the current Epn(3.5) template
family obtained by attaching to the inspiral waveform
three QNMs [10] around the EOB light-ring. In order
to properly take into account the energy and angular-
momentum released during the merger/ring-down phases
we introduced [10] two physical parameters, ǫM and ǫJ ,
whose dependence on the binary masses and spins will
be determined by future comparisons between EOB and
NR waveforms computed for different mass ratios and
spins. We found high FFs >∼0.96. Due to small differ-
ences between EOB and NR waveforms during the fi-
Total mass ( )
50 100 150 200
Enhanced LIGO
Template: NR
Template: SPAc
ext(3.5)
Template: Epn(3.5)
FIG. 16: The sky-average SNR for LIGO and Enhanced or
mid LIGO detector versus total mass for an equal-mass binary
at 100Mpc.
nal cycles of the evolution, the best-matches are reached
at the cost of large systematic error in the merger–ring-
down binary parameters. Thus, the Epn(3.5) template
family can be used for detection, but for parameter es-
timation it needs to be improved when matching to the
ring-down, and also during the inspiral phase. The re-
finements can be achieved (i) by introducing deviations
from circular motion, (ii) adding higher-order PN terms
in the EOB dynamics, (iii) using in the EOB radiation-
reaction equations a GW energy flux closer the the NR
flux, (iv) designing a better match to ring-down modes,
etc.. The goal would be to achieve dephasing between
EOB and NR waveforms of less than a few percent in
the comparable-mass case, as obtained in Ref. [48] in the
extreme mass-ratio limit. Indeed, with more accurate nu-
merical simulations, especially those using spectral meth-
ods [49], it will be possible to improve the inspiraling
templates by introducing higher-order PN terms in the
analytic waveforms computed by direct comparison with
NR waveforms.
Frequency-domain PN templates with an appropriate
cutoff frequency fcut provide high FFs (> 0.97), even
for large masses. This is due to oscillating tails (Gibbs
phenomenon) produced when cutting the signal in the
frequency domain. We tested the SPAextc (3.5) and the
SPAYc (4) template families for total masses up to 120M⊙,
and three mass ratios m2/m1 = 1, 1.5, and 2. We always
get FFs > 0.96, even when using a fixed cutoff frequency,
fcut = 1.07ω220/2π. Because of its high efficiency, faith-
fulness, i.e., low systematic error in parameter estima-
tion, and simple implementation, the SPAYc (4) template
family (or variants of it which include Lorentzians) is, to-
gether with the EOB model, a good candidate for search-
ing coherently for GWs from binary systems with total
masses up to 120M⊙.
In Fig. 16, we show the sky averaged SNRs of a sin-
gle LIGO and Enhanced or mid LIGO [50] detector,
for an equal-mass binary at 100Mpc. The SNR peaks
at the total binary mass M ≃ 150M⊙ and shows the
importance of pushing current searches for coalescing
binary systems to M > 100M⊙. In the mass range
30M⊙ < M < 120M⊙, the SNR drops only slightly if
we filter the GW signal with SPAextc (3.5) or Epn(3.5)
instead of using NR waveforms. The difference be-
tween Epn(3.5) and SPAextc (3.5) is almost indistinguish-
able. When M > 120M⊙, although the SPA
c (3.5) and
Epn(3.5) template families give fairly good SNRs, it is
maybe not a good choice to use them as the number of
cycles reduces to a few. The key problem in detecting
such GWs is how to veto triggers from non-Gaussian,
nonstationary noise, instead of matching the effectively
short signal. This is a general problem in searches for
short signals in ground-based detectors.
Acknowledgments
A.B. and Y.P. acknowledge support from NSF grant
PHY-0603762, and A.B. also from the Alfred Sloan Foun-
dation. The work at Goddard was supported in part
by NASA grants O5-BEFS-05-0044 and 06-BEFS06-19.
B.K. was supported by the NASA Postdoctoral Program
at the Oak Ridge Associated Universities. S.T.M. was
supported in part by the Leon A. Herreid Graduate Fel-
lowship. Some of the comparisons with PN and EOB
models were obtained building on Mathematica codes de-
veloped in Refs. [14, 30, 31, 45]
APPENDIX A: COMMENT ON WAVEFORMS
OBTAINED FROM THE ENERGY-BALANCE
EQUATION
In adiabatic PN models, like the Tpn model used in
this paper, waveforms are computed under the assump-
tion that the binary evolves through an adiabatic se-
quence of quasi-circular orbits. More specifically, one
sets ṙ = 0 and computes the orbital frequency ω from the
energy-balance equation dE(ω)/dt = F(ω), where E(ω)
is the total energy of the binary system and F(ω) is the
GW energy flux. Both E(ω) and F(ω) are computed
for circular orbits and expressed as a Taylor expansion
in ω. The adiabatic evolution ends in principle at the
innermost circular orbit (ICO) [35], or minimum energy
circular orbit (MECO) [30], where (dE/dω) = 0.
By rewriting the energy-balance equation, ω(t) can be
integrated directly as
ω̇(t) =
dE(ω)/dω
. (A1)
The RHS of Eq. (A1) can be expressed as an expansion
in powers of ω. The expanded version is widely used in
generating adiabatic PN waveforms [20, 30, 31, 45], it
is used to generate the so-called Tpn template family. It
turns out that Tpn(3) and Tpn(3.5) are quite close to the
NR inspiraling waveforms [10, 11]. We wonder whether
using the energy-balance in the form of Eq. (A1), i.e.,
without expanding it, might give PN waveforms closer to
or farther from NR waveforms. In principle the adiabatic
sequence of circular orbits described by Eq. (A1) ends at
the ICO, so the adiabatic model should work better until
the ICO and start deviating (with ω going to infinity)
from the exact result beyond it.
In Fig. 17 we show the NR orbital frequency ω(t) to-
gether with the PN orbital frequency obtained by solv-
ing the unexpanded and expanded form of the energy-
balance equation. The frequency evolution in these two
cases is rather different, with the orbital-frequency com-
puted from the expanded energy-balance equation agree-
ing much better with the NR one. When many, extremely
accurate, GW cycles from NR will be available, it will be
worthwhile to check whether this result is still true.
[1] A. Abramovici et al., Science 256, 325 (1992);
http://www.ligo.caltech.edu.
[2] H. Lück et al., Class. Quant. Grav. 14, 1471 (1997);
http://www.geo600.uni-hannover.de.
[3] M. Ando et al., Phys. Rev. Lett. 86, 3950 (2001);
http://tamago.mtk.nao.ac.jp.
[4] B. Caron et al., Class. Quant. Grav. 14, 1461 (1997);
http://www.virgo.infn.it.
[5] http://www.lisa-science.org/resources/
talks-articles/science/lisa_science_case.pdf
[6] F. Pretorius, Phys. Rev. Lett. 95, 121101 (2005).
[7] M. Campanelli, C.O. Lousto, P. Marronetti, and Y. Zlo-
chower, Phys. Rev. Lett. 96, 111101 (2006).
[8] J. Baker, J. Centrella, D. Choi, M. Koppitz, and J. van
Meter, Phys. Rev. Lett. 96, 111102 (2006).
[9] M. Campanelli, C.O. Lousto, and Y. Zlochower, Phys.
Rev. D 74, 041501 (2006); ibid. D 74, 084023 (2006); U.
Sperhake, gr-qc/0606079; J. González, U. Sperhake, B.
Brügmann, M. Hannam, and S. Husa, Phys. Rev. Lett.
98, 091101 (2007); B. Szilagyi, D. Pollney, L. Rezzolla,
http://www.ligo.caltech.edu
http://www.geo600.uni-hannover.de
http://tamago.mtk.nao.ac.jp
http://www.virgo.infn.it
http://www.lisa-science.org/resources/
talks-articles/science/lisa_science_case.pdf
gr-qc/0606079
100 200 300 400 500 600 700
Time(M)
Expanded energy-balance equation
Unexpanded energy-balance equation
NR simulation result
FIG. 17: Orbital frequency evolution. The dotted and dashed curves are calculated from the unexpanded and expanded
energy-balance equations. The continuous curve refers to the really long Goddard NR simulation.
J. Thornburg and J. Winicour, gr-qc/0612150; F. Pre-
torius and D. Khurana, gr-qc/0702084.
[10] A. Buonanno, G. Cook, and F. Pretorius, Phys. Rev. D
75 (2007) 124018.
[11] J. Baker, J. van Meter, S. McWilliams, J. Centrella, and
B. Kelly (2006), gr-qc/0612024.
[12] P. Jaranowski, and G. Schäfer, Phys. Rev. D 57, 7274
(1998); Erratum-ibid D 63 029902; L. Blanchet, and G.
Faye, Phys. Rev. D 63, 062005 (2001); V. C. de An-
drade, L. Blanchet, and G Faye, Class. Quant. Grav. 18,
753 (2001); T. Damour, P. Jaranowski, and G. Schäfer,
Phys. Lett. B 513, 147 (2001); L. Blanchet, G. Faye,
B.R. Iyer, and B. Joguet, Phys. Rev. D 65, 061501(R)
(2002); L. Blanchet, and B.R. Iyer, Class. Quant. Grav.
20, 755 (2003); Erratum-ibid D 71, 129902 (2005); L.
Blanchet, T. Damour, G. Esposito-Farese, and B.R. Iyer,
Phys. Rev. Lett. 93, 091101 (2004).
[13] A. Buonanno, and T. Damour, Phys. Rev. D 59, 084006
(1999).
[14] A. Buonanno, and T. Damour, Phys. Rev. D 62, 064015
(2000).
[15] A. Buonanno, and T. Damour, Proceedings of
IXth Marcel Grossmann Meeting (Rome, July 2000),
gr-qc/0011052.
[16] T. Damour, P. Jaranowski, and G. Schäfer, Phys. Rev.
D 62, 084011 (2000).
[17] T. Damour, Phys. Rev. D 64, 124013 (2001).
[18] A. Buonanno, Y. Chen, and T. Damour, Phys. Rev. D
74, 104005 (2006).
[19] T. Damour, B.R. Iyer, and B.S. Sathyaprakash, Phys.
Rev. D 57, 885 (1998).
[20] T. Damour, B. Iyer, and B. Sathyaprakash, Phys. Rev.
D 63, 044023 (2001); ibid. D 66, 027502 (2002).
[21] T. Baumgarte, P. Brady, J.D.E. Creighton, L. Lehner, F.
Pretorius, and R. DeVoe (2006), gr-qc/0612100.
[22] J. W. York, Jr., Phys. Rev. Lett. 82, 1350 (1999).
[23] E. Gourgoulhon, P. Grandclément, and S. Bonazzola,
Phys. Rev. D 65, 044020 (2002).
[24] P. Grandclément, E. Gourgoulhon, and S. Bonazzola,
Phys. Rev. D 65, 044021 (2002).
[25] H. P. Pfeiffer L. E. Kidder, M. S. Scheel, and S. A. Teukol-
sky, Comp. Phys. Comm. 152, 253 (2003).
[26] G. B. Cook, and H. P. Pfeiffer, Phys. Rev. D 70, 104016
(2004); M. Caudill, G.B. Cook, J.D. Grigsby, and H.
Pfeiffer, Phys. Rev. D 74, 064011 (2006).
[27] S. Brandt and B. Brügmann, Phys. Rev. Lett. 78, 3606
(1997).
[28] C. Cutler et al., Phys. Rev. Lett. 70, 2984 (1993).
[29] T. Damour, B.R. Iyer, and B.S. Sathyaprakash, Phys.
Rev. D 67, 064028 (2003).
[30] A. Buonanno, Y. Chen, and M. Vallisneri, Phys. Rev. D
67, 104025 (2003); Erratum-ibid. D 74, 029904 (2006).
[31] A. Buonanno, Y. Chen, and M. Vallisneri, Phys. Rev. D
67, 024016 (2003); Erratum-ibid. D 74, 029903 (2006).
[32] T. Damour, B. Iyer, P. Jaranowski, and B.
Sathyaprakash, Phys. Rev. D 67, 064028 (2003).
[33] K. G. Arun, B.R. Iyer, B.S. Sathyaprakash, and P. Sun-
dararajan, Phys. Rev. D 71, 084008 (2005); Erratum-ibid
D 72, 069903 (2005).
[34] T. Damour, E. Gourgoulhon, and P. Grandclément,
Phys. Rev. D 66, 024007 (2002); P. Grandclément, E.
Gourgoulhon, and S. Bonazzola, Phys. Rev. D 65, 044021
(2002).
[35] L. Blanchet, Phys. Rev. D 65, 124009 (2002).
[36] L.S. Finn, Phys. Rev. D 46, 5236 (1992); L. S. Finn and
D.F. Chernoff, Phys. Rev. D 47, 2198 (1993); É.E. Flana-
gan and S.A. Hughes, Phys. Rev. D 57, 4535 (1998).
[37] http://www.ligo.caltech.edu/advLIGO/scripts/ref_
des.shtml
[38] J. Baker, S. McWilliams, J. van Meter, J. Centrella, D.
Choi, B. Kelly, and M. Koppitz (2006), gr-qc/0612117.
[39] J. Baker, M. Campanelli, F. Pretorius, and Y. Zlochower
(2007), gr-qc/0701016.
[40] B. Abbott et al. (LIGO Scientific Collaboration), Phys.
Rev.D 72, 082001 (2005).
[41] B. Abbott et al. (LIGO Scientific Collaboration), Phys.
Rev. D 73, 062001 (2006).
[42] E. Berti, V. Cardoso, and C. Will, Phys. Rev. D 73,
064030 (2006).
[43] L. Blanchet, Living Rev. Rel. 9 (2006) 4.
[44] T. Damour, B. Iyer, and B. Sathyaprakash, Phys. Rev.
D 62, 084036 (2000).
gr-qc/0612150
gr-qc/0702084
gr-qc/0612024
gr-qc/0011052
gr-qc/0612100
http://www.ligo.caltech.edu/advLIGO/scripts/ref_
des.shtml
gr-qc/0612117
gr-qc/0701016
[45] Y. Pan, A. Buonanno, Y. Chen, and M. Vallisneri, Phys.
Rev. D 69, 104017 (2004).
[46] P. Ajith et al. (2007) (in preparation).
[47] E. Berti, V. Cardoso, J. González, U. Sperhake, M. Han-
nam, S. Husa, and B. Brügmann (2007), gr-qc/0703053.
[48] A. Nagar, T. Damour, and A. Tartaglia, gr-qc/0612096
T. Damour, and A. Nagar, Proceedings of XIth Marcel
Grossmann Meeting (Berlin, July 2006), gr-qc/0612151.
[49] H. P. Pfeiffer, D.A. Brown, L.E. Kidder, L. Lindblom, G.
Lovelace, and M. A. Scheel (2007), gr-qc/0702106.
[50] http://www.ligo.caltech.edu/~rana/NoiseData/S6/
DCnoise.txt.
gr-qc/0703053
gr-qc/0612096
gr-qc/0612151
gr-qc/0702106
http://www.ligo.caltech.edu/~rana/NoiseData/S6/
DCnoise.txt

  We compare waveforms obtained by numerically evolving nonspinning binary
black holes to post-Newtonian (PN) template families currently used in the
search for gravitational waves by ground-based detectors. We find that the
time-domain 3.5PN template family, which includes the inspiral phase, has
fitting factors (FFs) >= 0.96 for binary systems with total mass M = 10 ~ 20
Msun. The time-domain 3.5PN effective-one-body template family, which includes
the inspiral, merger and ring-down phases, gives satisfactory signal-matching
performance with FFs >= 0.96 for binary systems with total mass M = 10 ~ 120
Msun. If we introduce a cutoff frequency properly adjusted to the final
black-hole ring-down frequency, we find that the frequency-domain
stationary-phase-approximated template family at 3.5PN order has FFs >= 0.96
for binary systems with total mass M = 10 ~ 20 Msun. However, to obtain high
matching performances for larger binary masses, we need to either extend this
family to unphysical regions of the parameter space or introduce a 4PN order
coefficient in the frequency-domain GW phase. Finally, we find that the
phenomenological Buonanno-Chen-Vallisneri family has FFs >= 0.97 with total
mass M=10 ~ 120Msun. The main analyses use the noise spectral-density of LIGO,
but several tests are extended to VIRGO and advanced LIGO noise-spectral
densities.

A data-analysis driven comparison of analytic and numerical coalescing binary
waveforms: nonspinning case
Yi Pan,1 Alessandra Buonanno,1 John G. Baker,2 Joan Centrella,2 Bernard
J. Kelly,2 Sean T. McWilliams,1 Frans Pretorius,3 and James R. van Meter2, 4
Department of Physics, University of Maryland, College Park, MD 20742
Gravitational Astrophysics Laboratory, NASA Goddard Space Flight Center, 8800 Greenbelt Rd., Greenbelt, MD 20771
Department of Physics, Princeton University, Princeton, NJ 08544
Center for Space Science & Technology, University of Maryland Baltimore County,
Physics Department, 1000 Hilltop Circle, Baltimore, MD 21250
(Dated: October 22, 2018)
We compare waveforms obtained by numerically evolving nonspinning binary black holes to post-
Newtonian (PN) template families currently used in the search for gravitational waves by ground-
based detectors. We find that the time-domain 3.5PN template family, which includes the inspiral
phase, has fitting factors (FFs) ≥ 0.96 for binary systems with total massM = 10–20M⊙. The time-
domain 3.5PN effective-one-body template family, which includes the inspiral, merger and ring-down
phases, gives satisfactory signal-matching performance with FFs ≥ 0.96 for binary systems with total
mass M = 10–120M⊙. If we introduce a cutoff frequency properly adjusted to the final black-hole
ring-down frequency, we find that the frequency-domain stationary-phase-approximated template
family at 3.5PN order has FFs ≥ 0.96 for binary systems with total mass M = 10–20M⊙. However,
to obtain high matching performances for larger binary masses, we need to either extend this family
to unphysical regions of the parameter space or introduce a 4PN order coefficient in the frequency-
domain GW phase. Finally, we find that the phenomenological Buonanno-Chen-Vallisneri family
has FFs ≥ 0.97 with total mass M = 10–120M⊙. The main analyses use the noise spectral-density
of LIGO, but several tests are extended to VIRGO and advanced LIGO noise-spectral densities.
PACS numbers: 04.25.Dm, 04.30.Db, 04.70.Bw, x04.25.Nx, 04.30.-w, 04.80.Nn, 95.55.Ym
I. INTRODUCTION
The search for gravitational-waves (GWs) from coa-
lescing binary systems with laser interferometer GW de-
tectors [1, 2, 3, 4, 5] is based on the matched-filtering
technique, which requires accurate knowledge of the
waveform of the incoming signal. In the last couple of
years there have been several breakthroughs in numeri-
cal relativity (NR) [6, 7, 8], and now independent groups
are able to simulate the inspiral, merger and ring-down
phases of generic black- hole (BH) merger scenarios, in-
cluding different spin orientations and mass ratios [9].
However, the high computational cost of running such
simulations makes it difficult to generate sufficiently long
inspiral waveforms that cover the parameter space of as-
trophysical interest.
References [10, 11] found good agreement between ana-
lytic (based on the post-Newtonian (PN) expansion) and
numerical waveforms emitted during the inspiral, merger
and ring-down phases of equal-mass, nonspinning binary
BHs. Notably, the best agreement is obtained with 3PN
or 3.5PN adiabatic waveforms [12] (henceforth denoted
as Taylor PN waveforms) and 3.5PN effective-one-body
(EOB) waveforms [13, 14, 15, 16, 17, 18, 19]. In addi-
tion to the inspiral phase the latter waveforms include
the merger and ring-down phases. Those comparisons
suggested that it should be possible to design hybrid nu-
merical/analytic templates, or even purely analytic tem-
plates with the full numerics used to guide the patching
together of the inspiral and ring-down waveforms. This is
an important avenue to template construction as eventu-
ally thousands of waveform templates may be needed to
extract the signal from the noise, an impossible demand
for NR alone. Once available, those templates could be
used by ground-based laser interferometer GW detectors,
such as LIGO, VIRGO, GEO and TAMA, and in the fu-
ture by the laser interferometer space antenna (LISA) for
detecting GWs emitted by solar mass and supermassive
binary BHs, respectively.
This paper presents a first attempt at investigating the
closeness of the template families currently used in GW
inspiral searches to waveforms generated by NR simula-
tions. Based on this investigation, we shall propose ad-
justments to the templates so that they include merger
and ring-down phases. In contrast, Ref. [21] examined
the use of numerical waveforms in inspiral searches, and
compared numerical waveforms to the ring-down tem-
plates currently used in burst searches. Similar to the
methodology presented here, fitting factors (FFs) [see
Eq. (2) below] are used in Ref. [21] to quantify the ac-
curacy of numerical waveforms for the purpose of de-
tection, as well as the overlap of burst templates with
the waveforms. Reference [21] found that by computing
FFs between numerical waveforms from different reso-
lution simulations of a given event, one can recast the
numerical error as a maximum FF that the numerical
waveform can resolve. In other words, any other tem-
plate or putative signal convolved with the highest reso-
lution numerical simulation that gives a FF equal to or
larger than this maximum FF is, for the purpose of de-
http://arxiv.org/abs/0704.1964v2
tection, indistinguishable from the numerical waveform.
We will explore this aspect of the problem briefly. The
primary conclusions we will draw from the analysis do
not crucially depend on the exactness of the numerical
waveforms. What counts here is that the templates can
capture the dominant spectral characteristics of the true
waveform.
For our analysis we shall focus on two nonspinning
equal-mass binary simulation waveforms which differ in
length, initial conditions, and the evolution codes used
to compute them: Cook-Pfeiffer quasi-equilibrium ini-
tial data built on Refs. [22, 23, 24, 25, 26] evolved with
Pretorius’ generalized harmonic code [6], and Brandt-
Brügmann puncture data [27] evolved using the Goddard
group’s moving-puncture code [8]. We also consider two
nonspinning unequal-mass binary simulations with mass
ratios m2/m1 = 1.5 and m2/m1 = 2 produced by the
Goddard group.
The paper is organized as follows. In Sec. II we dis-
cuss the phase differences between PN inspiraling tem-
plates. In Sec. III we build hybrid waveforms by stitch-
ing together PN and NR waveforms. We try to under-
stand how many NR cycles are needed to obtain good
agreement between NR and PN waveforms, to offer a
guide for how long PN waveforms can be used as accu-
rate templates. In Sec. IV we compute the FFs between
several PN template families and NR waveforms. We
first focus on low-mass binary systems with total mass
M = 10–30M⊙, then high-mass binary systems with to-
tal mass M = 30–120M⊙. Finally, Sec. V contains our
main conclusions. In Appendix A we comment on how
different representations of the energy-balance equations
give GW frequencies closer to or farther from the NR
ones.
II. PHASE DIFFERENCES IN
POST-NEWTONIAN INSPIRALING MODELS
Starting from Ref. [28], which pointed out the impor-
tance of predicting GW phasing with the highest possi-
ble accuracy when building GW templates, many subse-
quent studies [14, 18, 19, 20, 29, 31, 32, 33] (those ref-
erences are restricted to the nonspinning case) focused
on this issue and thoroughly tested the accuracy of those
templates, proposing improved representations of them.
These questions were motivated by the observation that
comparable-mass binary systems with total mass higher
than 30M⊙ merge in-band with the highest signal-to-
noise ratio (SNR) for LIGO detectors, It follows that the
corresponding templates demand an improved analysis.
In the absence of NR results and under the urgency
of providing templates to search for comparable-mass bi-
nary BHs, the analytic PN community pushed PN cal-
culations to higher PN orders, notably 3.5PN order [12],
and also proposed ways of resumming the PN expan-
sion, either for conservative dynamics (the EOB ap-
proach [13, 16, 17]), radiation-reaction effects (the Padé
resummation [19]), or both [14, 18]. Those results lead to
several conclusions: (i) 3PN terms improve the compar-
ison between analytic and (numerical) quasi-equilibrium
predictions [23, 26, 34, 35]; (ii) Taylor expanded and re-
summed PN predictions for equal-mass binary systems
are much closer at 3.5PN order than at previous PN
orders, indicating a convergence between the different
schemes [18, 20, 31]; (iii) the two-body motion is quasi-
circular until the end of a rather blurred plunge [14], (iv)
the transition to ring-down can be described by an ex-
tremely short merger phase [14, 18]. Today, with the NR
results we are in a position to sharpen the above con-
clusions, and to start to assess the closeness of analytic
templates to numerical waveforms.
Henceforth, we restrict the analysis to the three time-
domain physical template families which are closest to
NR results [10, 11]: the adiabatic Taylor PN model
(Tpn) [see, e.g., Eqs. (1), (10), and (11)–(13) in Ref. [30]]
computed at 3PN and 3.5PN order, and the nonadia-
batic EOB model (Epn) [see e.g., Eqs. (3.41)–(3.44) in
Ref. [14]] computed at 3.5 PN order. We shall denote
our models as Tpn(n) and Epn(n), n being the PN or-
der. The Tpn model is obtained by solving a particular
representation of the balance equation. In Appendix A
we briefly discuss how time-domain PN models based on
different representations of the energy-balance equation
would compare with NR results.
The waveforms we use are always derived in the so-
called restricted approximation which uses the amplitude
at Newtonian order and the phase at the highest PN
order available. They are computed by solving PN dy-
namical equations providing the instantaneous frequency
ω(t) and phase φ(t) = φ0 +
ω(t′)dt′, thus
h(t) = Aω(t)2/3 cos[2φ(t)] , (1)
where t0 and φ0 are the initial time and phase, respec-
tively, and A is a constant amplitude, irrelevant to our
discussion. The inclusion of higher-order PN corrections
to the amplitude can be rather important for certain
unequal-mass binary systems, and will be the subject of
a future study.
When measuring the differences between waveforms
we weight them by the power spectral-density (PSDs) of
the detector, and compute the widely used fitting factor
(FF) (i.e., the ambiguity function or normalized over-
lap), or equivalently the mismatch defined as 1-FF. Fol-
lowing the standard formalism of matched-filtering [see,
e.g., Refs. [19, 31, 36]], we define the FF as the overlap
〈h1(t), h2(t)〉 between the waveforms h1(t) and h2(t):
〈h1(t), h2(t)〉 ≡ 4Re
h̃1(f)h̃
Sh(f)
FF ≡ max
t0,φ0,λi
〈h1, h2(t0, φ0, λ
〈h1, h1〉〈h2(t0, φ0, λi), h2(t0, φ0, λi)〉
where h̃i(f) is the Fourier transform of hi(t), and Sh(f)
is the detector’s PSD. Thus, the FF is the normalized
overlap between a target waveform h1(t) and a set of tem-
plate waveforms h2(t0, φ0, λ
i) maximized over the initial
time t0, initial phase φ0, and other parameters λ
i. Some-
times we are interested in FFs that are optimized only
over t0 and φ0; we shall denote these as FF0. For data-
analysis purposes, the FF has more direct meaning than
the phase evolution of the waveforms, since it takes into
account the PSDs and is proportional to the SNR of the
filtered signal. Since the event rate is proportional to
the cube of the SNR, and thus to the cube of the FF, a
FF= 0.97 corresponds to a loss of event rates of ∼ 10%.
A template waveform is considered a satisfactory repre-
sentation of the target waveform when the FF is larger
than 0.97.
When comparing two families of waveforms, the FF is
optimized over the initial phase of the template wave-
form, and we also need to specify the initial phase of the
target waveform. Since there is no preferred initial phase
of the target, two options are usually adopted: (i) the ini-
tial phase maximizes the FF or (ii) it minimizes the FF.
The resulting FFs are referred to as the best and mini-
max FFs [29]. All FFs we present in this paper are min-
imax FFs. Although the FF of two waveform families is
generally asymmetric under interchange of the template
family [31], the best and the minimax FF0s are symmet-
ric (see Appendix B of Ref. [29] for details). Henceforth,
when comparing two waveform families using FF0, we do
not need to specify which family is the target.
We shall consider three interferometric GW detectors:
LIGO, advanced LIGO and VIRGO. The latter two have
better low-frequency sensitivity and broader bandwidth.
For LIGO, we use the analytic fit to the LIGO de-
sign PSD given in Ref. [20]; for advanced LIGO we use
the broadband configuration PSD given in Ref. [37]; for
VIRGO we use the PSD given in Ref. [20].
In Fig. 1, we show the FF0s as functions of the ac-
cumulated difference in the number of GW cycles be-
tween waveforms generated with different inspiraling PN
models and for binary systems with different component
masses. We first generate two waveforms by evolving
two PN models, say, “PN1” and “PN2” which start at
the same GW frequency fGW = 30Hz and have the same
initial phase. The two waveforms are terminated at the
same ending frequency fGW = fend up to a maximum
fend,max = min(fend,PN1 , fend,PN2), where fend,PN is the
frequency at which the PN inspiraling model ends. (For
Tpn models this is the frequency at which the PN energy
has a minimum; for Epn models it is the EOB light-ring
frequency.) Then, we compute the difference in phase
and number of GW cycles accumulated until the ending
frequency
∆NGW =
[φPN1(fend)− φPN2(fend)] . (3)
By varying fend (up to fend,max) ∆NGW changes, though
not necessarily monotonically. Although there seems to
be a loose correlation between the FF0s and ∆NGW, it
is hard to quantify it as a one-to-one correspondence.
0 0.2 0.4 0.6 0.8 1.0 1.2
Tpn(3.5) & Epn(3.5), (3+3)
Tpn(3.5) & Epn(3.5), (15+3)
Tpn(3.5) & Epn(3.5), (15+15)
Tpn(3) & Tpn(3.5), (3+3)
Tpn(3) & Tpn(3.5), (15+3)
Tpn(3) & Tpn(3.5), (15+15)
FIG. 1: We show FF0s between waveforms generated from
the three PN models Tpn(3), Tpn(3.5) and Epn(3.5) versus
∆NGW [see Eq. (3)]. The FF0s are evaluated with LIGO’s
PSD. Note that for Tpn(3.5) and Epn(3.5) and a (15+3)M⊙
binary, the lowest FF0 is 0.78 and the difference in the number
of GW cycles ∆NGW ≃ 2. In the limit ∆NGW → 0, the FF0
goes to unity.
For example, a phase difference of about half a GW cy-
cle (∆NGW ≃ 0.5) is usually thought to be significant.
However, here we find relatively high FF0s between 0.97
and > 0.99, depending on the masses of the binary and
the specific PN model used. This happens because the
FF between two waveforms is not determined by the to-
tal phase difference accumulated, but rather by how the
phase difference accumulates across the detector’s most
sensitive frequency band. The relation between FFs and
phase differences is also blurred by the maximization
over the initial time and phase: shifting the phase by
half a cycle from the most sensitive band to a less sen-
sitive band can increase the matching significantly. We
conclude that with LIGO’s PSD, after maximizing only
on initial phase and time, Epn(3.5) and Tpn(3.5) tem-
plates are close to each other for comparable-mass bi-
nary systems M = 6–30M⊙ with FF0>∼ 0.97, but they
can be different for mass ratios m2/m1 ≃ 0.3 with FF0
as low as ≃ 0.8. Tpn(3) and Tpn(3.5) templates have
FF0>∼ 0.97 for the binary masses considered. Note that
for m2/m1 = 1 [≃ 0.3] binary systems, Tpn(3.5) is closer
to Epn(3.5) [Tpn(3)] than to Tpn(3) [Epn(3.5)]. Note
also that when maximizing on binary masses the FFs can
increase significantly, for instance, for a (15 + 3)M⊙ bi-
nary, the FF between Tpn(3.5) and Epn(3.5) waveforms
becomes > 0.995, whereas FF0 ≃ 0.8.
III. BUILDING AND COMPARING HYBRID
WAVEFORMS
Recent comparisons [10, 11] between analytic and nu-
merical inspiraling waveforms of nonspinning, equal-mass
binary systems have shown that numerical waveforms are
in good agreement with Epn(3.5), Tpn(3) and Tpn(3.5)
waveforms. Those results were assessed using eight and
sixteen numerical inspiral GW cycles. Can we conclude
from these analyses that Epn(3.5), Tpn(3.5) and Tpn(3)
can safely be used to build a template bank for detecting
inspiraling GW signals? A way to address this question
is to compute the mismatch between hybrid waveforms
built by attaching either Epn or Tpn waveforms to the
same numerical waveform, and varying the time when
the attachment is made. This is equivalent to varying
the number of numerical GW cycles n in the hybrid tem-
plate. The larger n the smaller the mismatch, as we are
using the same numerical segment in both waveforms.
For a desired maximum mismatch, say 3%, we can then
find the smallest number n of numerical cycles that is
required in the hybrid waveform. This number will, of
course, depend on the binary mass and the PSD of each
detector.
A. Hybrid waveforms
We build hybrid waveforms by connecting PN wave-
forms to NR waveforms at a chosen point in the late
inspiral stage. As mentioned before, we use NR wave-
forms generated with Pretorius’ [10] code and the God-
dard group’s [38] code. Pretorius’ waveform is from an
equal-mass binary with total mass M , and equal, co-
rotating spins (a = 0.06). The simulation lasts ≃ 671M ,
and the waveform has ≃ 8 cycles before the formation of
the common apparent horizon. The Goddard waveform
refers to an equal-mass nonspinning binary. The simula-
tion lasts about ≃ 1516M , and the waveform has ≃ 16
cycles before merger.
Since we will present results from these two waveforms
it is useful to first compare them by computing the FF0.
Although the binary parameters considered by Pretorius
and Goddard are slightly different, we expect the wave-
forms, especially around the merger stage, to be fairly
close. Comparisons between (shorter) waveforms com-
puted with moving punctures and generalized-harmonic
gauge were reported in Ref. [39], where the authors dis-
cussed the different initial conditions, wave extraction
techniques, and compared the phase, amplitude and fre-
quency evolutions. Since the two simulations use dif-
ferent initial conditions and last for different amounts of
time we cut the longer Goddard waveform at roughly the
frequency where the Pretorius waveform starts. In this
way we compare waveforms that have the same length
between the initial time and the time at which the wave
amplitude reaches its maximum. In Fig. 2, we show the
FF0 as function of the total binary mass. Despite dif-
ferences in the two simulations the FF0s are rather high.
The waveforms differ more significantly at lower frequen-
cies. Indeed, as the total mass decreases the FF0s also
decrease as these early parts of the waveform contribute
more to the signal power given LIGO’s PSD.
Any waveform extracted from a numerical simulation
will inherit truncation errors, affecting both the wave-
form’s amplitude and phase [10, 21, 38]. To check
whether those differences would change the results of the
comparisons between NR and PN waveforms, we plot
in Fig. 2 the FF0s versus total binary mass between
two Goddard waveforms generated from a high and a
medium resolution run [38]. The FF0s are extremely high
(> 0.995).
Based on the comparisons between high and medium
resolution waveforms, we can estimate the FFs between
high resolution and exact waveforms. If we have several
simulations with different resolutions, specified by the
mesh-spacings xi, and xi are sufficiently small, we can
assume that the waveforms hi are given by
hi = h0 + x
i hd , (4)
where n is the convergence factor of the waveform, h0 is
the exact waveform generated from the infinite resolution
run (x0 → 0), and hd is the leading order truncation error
contribution to the waveform and is independent of the
mesh spacing xi. We find that the mismatch between the
waveforms hi and hj, 1− FFij , scales as
1− FFij ∝ (x
i − x
2 . (5)
In the Goddard simulations, the high and medium reso-
lution runs have mesh-spacing ratio xh/xm = 5/6, and
the waveform convergence rate is n = 4 [38]. The FF
between the high resolution and exact waveforms hh and
h0 is given by
FF0h = 1− 0.87(1− FFhm) , (6)
where FFhm is the FF between the high and medium
resolution waveforms hh and hm. That is to say, the
mismatch between hh and h0 is slightly smaller than that
between hh and hm, where the latter can be derived from
the FFs shown in Fig. 2. Henceforth, we shall always use
high-resolution waveforms. A similar calculation for Pre-
torius’ waveform gives FF0h = 1−0.64(1−FFhm), though
here xh/xm = 2/3 and n = 2. See Fig. 6 of Ref. [21] for a
plot of FFhm calculated from the evolution of the Cook-
Pfeiffer initial data 1; there FF0 ranges from ≈ 0.97 for
M/Ms = 30 to ≈ 0.99 for M/Ms = 100. In other words,
the mismatch between Goddard’s and Pretorius’ wave-
form shown in Fig.2 is less than the estimated mismatch
from numerical error in the latter waveform.
We build hybrid waveforms by stitching together the
PN and NR waveforms computed for binary systems with
the same parameters. At the point where we connect the
1 The plot in Ref. [21] is for “d=16” corotating Cook-Pfeiffer initial
data, whereas the results presented here are from “d=19” initial
data. However, the resolutions used for both sets were the same,
and thus the mismatches should be similar, in particular in the
higher mass range.
30 40 50 60 70 80 90 100
Total mass of the binary (      )
0.980
0.985
0.990
0.995
1.000
Goddard high and medium resolutions
Pretorius and Goddard waveforms
FIG. 2: FF0 between NR waveforms as a function of the
binary total-mass M . The solid curve are generated for wave-
forms from Pretorius and the Goddard group. The longer
Goddard waveform is shortened such that both waveforms last
≃ 671M and contain ≃ 8 cycles. The dashed curve is gen-
erated for waveforms from the high-resolution and medium
resolution simulations of the Goddard group. All FFs are
evaluated using LIGO’s PSD.
two waveforms, we tune the initial time t0 so that the
frequency of the PN waveform is almost the same as the
frequency of the NR waveform (there is a subtlety trying
to match exactly the frequencies that is discussed at the
end of this section). The initial phase φ0 is then chosen
so that the strain of the hybrid waveform is continuous
at the connecting point.
In Fig. 3, we show two examples of hybrid waveforms of
an equal-mass binary. We stitch the waveforms at points
where effects due to the initial-data transient pulse are
negligible. We find an amplitude difference on the order
of ∼ 10% between the Goddard waveform and the re-
stricted PN waveform. This difference is also present in
Pretorius’ waveform, but it is somewhat compensated for
by amplitude modulations caused by eccentricity in the
initial data. In Ref. [38] it was shown that PN waveforms
with 2.5PN amplitude corrections give better agreement
(see e.g., Fig. 12 in Ref. [38]). However, the maximum
amplitude errors in the waveforms are also on the or-
der of 10% [10, 38]. Since neither 2PN nor other lower
PN order corrections to the amplitudes are closer to the
2.5PN order, we cannot conclude that 2.5PN amplitude
corrections best approximate the numerical waves. Thus,
we decide to use two sets of hybrid waveforms: one con-
structed with restricted PN waveforms, and the other
with restricted PN waveforms rescaled by a single ampli-
tude factor, which eliminates amplitude differences with
the NR waveforms. We shall see that the difference be-
tween these two cases is small for the purpose of our
tests.
The amplitude difference between PN and NR wave-
forms is computed at the same connecting-point GW fre-
-500-1000
-1000-1500
-1500-2000
-2000-2500
-2500
Time (M)
Time (M)
FIG. 3: We show two examples of hybrid waveforms, start-
ing from 40Hz. The PN waveforms are generated with the
Tpn(3.5) model, and the NR waveforms in the upper and
lower panels are generated from Pretorius’ and Goddard’s
simulations, respectively. We mark with a dot the point where
we connect the PN and NR waveforms.
quency. There is another effect which causes a jump
in the hybrid-waveform amplitude. This is a small fre-
quency difference between PN and NR waveforms at the
connecting point. All our NR waveforms contain small
eccentricities [10, 38]. As a consequence, the frequency
evolution ω(t) oscillates. To reduce this effect we fol-
low what is done in Ref. [38] and fit the frequency to a
monotonic quartic function. When building the hybrid
waveform, we adjust the PN frequency to match the quar-
tic fitted frequency (instead of the oscillatory, numerical
frequency) at the connecting point. Since the restricted
PN amplitude is proportional to ω2/3(t) [see Eq. (1)],
this slight difference between ωs at the connecting point
creates another difference between the NR and PN am-
plitudes. Nevertheless, this difference is usually smaller
(for Goddard’s waveform) or comparable (for Pretorius’)
to the amplitude difference discussed above.
-25 -20 -15 -10 -5 0
Time (100M)
Original
Whitened
-25 -20 -15 -10 -5 0
Time (100M)
-800 -600 -400 -200 0
Time (M)
-800 -600 -400 -200 0
Time (M)
-200 -150 -100 -50 0 50
Time (M)
-200 -150 -100 -50 0 50
Time (M)
(5+5)
(10+10) (20+20)
(15+15) (30+30)
(50+50)
FIG. 4: Distribution of GW signal power. In each panel, we plot a hybrid waveform (a Tpn waveform stitched to the Goddard
waveform) in both its original form (blue curve) and its “whitened” form (red curve) [44]. We show waveforms from six
binary systems with total masses 10M⊙ 20M⊙, 30M⊙, 40M⊙, 60M⊙ and 100M⊙. The vertical lines divide the waveforms into
segments, where each segment contributes 10% of the total signal power.
B. Distribution of signal power in gravitational
waveforms
To better understand the results of the FFs between
hybrid waveforms, we want to compute how many sig-
nificant GW cycles are in the LIGO frequency band. By
significant GW cycles we mean the cycles that contribute
most to the signal power, or to the SNR of the filtered
signal. Since GW frequencies are scaled by the total bi-
nary mass, the answer to this question depends on both
the PSD and the total mass of a binary.
In Fig. 4, we show the effect of the LIGO PSD on the
distribution of signal power for several waves emitted by
coalescing binary systems with different total masses. In
each panel, we plot a hybrid waveform (a Tpn waveform
stitched to the Goddard waveform) in both its original
form and its “whitened” form [44]. The whitened wave-
form is generated by Fourier-transforming the original
waveform into the frequency domain, rescaling it by a fac-
tor of 1/
Sh(f), and then inverse-Fourier-transforming
it back to the time domain. The reference time t = 0 is
the peak in the amplitude of the unwhitened waveform.
The amplitude of a segment of the whitened waveform
indicates the relative contribution of that segment to the
signal power and takes into account LIGO’s PSD. Both
waveforms are plotted with arbitrary amplitudes, and the
unwhitened one always has the larger amplitude. The
absolute amplitude of a waveform, or equivalently the
distance of the binary, is not relevant in these figures un-
less the redshift z becomes significant. In this case the
mass of the binary is the redshifted mass (1+ z)M . Ver-
tical lines in each figure divide a waveform into segments,
where each segment contributes 10% of the total signal
power. In each plot, except for the 10M⊙-binary one, we
show all 9 vertical lines that divide the waveforms into
10 segments. In the 10M⊙-binary plot we omit the early
part of the inspiral phase that accounts for 50% of the
signal power, as it would be too long to show.
The absolute time-scale of a waveform increases lin-
early with total mass M ; equivalently the waveform is
shifted toward lower frequency bands. For a M = 10M⊙
binary, the long inspiral stage generates GWs with fre-
quencies spanning the most sensitive part of the LIGO
band, around 150Hz, while for an M = 100M⊙ binary,
only the merger signal contributes in this band. Thus,
for low-mass binary systems, most of the contribution to
the signal power comes from the long inspiral stage of
the waveform, while for high-mass binary systems most
of the contribution comes from the late inspiral, merger,
and ring-down stages. Understanding quantitatively the
distribution of signal power will let us deduce how many,
and which, GW cycles are significant for the purpose of
data analysis. We need accurate waveforms from either
PN models or NR simulations for at least those signifi-
cant cycles.
From Fig. 4 we conclude that:
• For an M = 10M⊙ binary, the last 25 inspiral cy-
cles, plus the merger and ring-down stages of the
waveform contribute only 50% of the signal power,
and we need 80 cycles (not shown in the figure) of
accurate inspiral waveforms to recover 90% of the
signal power. For an M = 20M⊙ binary, the last
23 cycles, plus the merger and ring-down stages of
the waveform contribute> 90% of the signal power,
and current NR simulations can produce waveforms
of such length;
• For an M = 30M⊙ binary, the last 11 inspiral cy-
cles, plus the merger and ring-down stages of the
waveform contribute > 90% of the signal power,
which means that, for binary systems with to-
tal masses higher than 30M⊙, current NR simula-
tions, e.g., the sixteen cycles obtained in Ref. [38],
can provide long enough waveforms for a matched-
filter search of binary coalescence, as also found in
Ref. [21];
• For an M = 100M⊙ binary, > 90% of the signal
power comes from the last inspiral cycle, merger
and ring-down stages of the waveform, with two
cycles dominating the signal power. It is thus pos-
sible to identify this waveform as a “burst” signal.
Similar analyses can be also done for advanced LIGO and
VIRGO.
C. Comparing hybrid waveforms
We shall now compute FF0s between hybrid wave-
forms. We fix the total mass of the equal-mass binary
in each comparison, i.e., we do not optimize over mass
parameters, but only on phase and time. We use the
mismatch, defined as 1− FF0, to measure the difference
between waveforms and we compute them for LIGO, ad-
vanced LIGO, and VIRGO. Note that by using FF0, we
test the closeness among hybrid waveforms that are gen-
erated from binary systems with the same physical pa-
rameters; in other words, we test whether the waveforms
are accurate enough for the purpose of parameter estima-
tion, rather than for the sole purpose of detecting GWs.
In the language of Ref. [19] we are studying the faithful-
ness of the PN templates 2.
Since at late inspiral stages PN waveforms are partly
replaced by NR waveforms, differences between hybrid
waveforms from two PN models are smaller than those
between pure PN waveforms. In general, the more NR
cycles we use to generate hybrid waveforms, the less the
2 Following Ref. [19], faithful templates are templates that have
large overlaps, say >∼ 96.5%, with the expected signal maximiz-
ing only over the initial phase and time of arrival. By contrast
when the maximization is done also on the binary masses, the
templates are called effectual.
difference is expected to be between these hybrid wave-
forms. This is evident in Figs. 5, 6 where we show mis-
matches between hybrid waveforms for binary systems
with different total masses as a function of the number
of NR cycles n. Specifically, the mismatches are taken
between two hybrid waveforms generated from the same
NR waveform (from the Goddard group, taking the last
n cycles, plus merger and ring-down) and two different
PN waveforms generated with the same masses.
The mismatches are lower for binary systems with
higher total masses, since most of their signal power
is concentrated in the late cycles close to merger (see
Fig. 4). Comparing results between LIGO, advanced
LIGO and VIRGO, we see that for the same waveforms
the mismatches are lowest when evaluated with the LIGO
PSD, and highest when evaluated with the VIRGO PSD.
This is due to the much broader bandwidth of VIRGO,
especially at low frequency: the absolute sensitivity is not
relevant; only the shape of the PSD matters. In VIRGO,
the inspiral part of a hybrid waveform has higher weight-
ing in its contribution to the signal power. As already
observed at the end of Sec. II, we can see also that the
difference between the Epn(3.5) and Tpn(3.5) models is
smaller than that between the Tpn(3) and Tpn(3.5) mod-
Figures 5, 6 show good agreement among hybrid wave-
forms. In Sec. IV, as a further confirmation of what was
found in Refs. [10, 11], we shall see that PN waveforms
from Tpn and Epn models have good agreement with the
inspiral phase of the NR waveforms. Therefore, we argue
that hybrid waveforms are likely to have high accuracy.
In fact, for the late evolution of a compact binary, where
NR waveforms are available, the PN waveforms are close
to the NR waveforms, while for the early evolution of the
binary, where we expect the PN approximations to work
better, the PN waveforms (from Tpn and Epn models)
are close to each other. Based on these observations,
we draw the following conclusions for LIGO, advanced
LIGO, and VIRGO data-analysis:
• For binary systems with total mass higher than
30M⊙, the current NR simulations of equal-mass
binary systems (16 cycles) are long enough to re-
duce mismatches between hybrid waveforms gen-
erated from the three PN models to below 0.5%.
Since these FFs are achieved without optimizing
the binary parameters, we conclude that for these
high-mass binary systems, the small difference be-
tween hybrid waveforms indicates low systematic
error in parameter estimation, i.e., hybrid wave-
forms are faithful [19].
• For binary systems with total mass around
10–20M⊙, 16 cycles of NR waveforms can reduce
the mismatch to below 3%, which is usually set as
the maximum tolerance for data-analysis purpose
(corresponding to ∼ 10% loss in event rate). By a
crude extrapolation of our results, we estimate that
with 30 NR waveform cycles, the mismatch might
4 6 8 10 12 14
Number of NR waveform cycles
Restricted PN
Rescaled PN
4 6 8 10 12 14
Number of NR waveform cycles
Quartic-fit frequency
Oscillatory frequency
FIG. 5: We show the mismatch between hybrid waveforms as a function of the number of NR waveform cycles used to generate
the hybrid waveforms. The LIGO PSD is used to evaluate the mismatches. In the left panel, we compare the Epn(3.5) and
Tpn(3.5) models. In the right panel, we compare the Tpn(3) and Tpn(3.5) models. From top to bottom, the four curves
correspond to four equal-mass binary systems, with total masses 10M⊙, 20M⊙, 30M⊙, and 40M⊙. The dots show mismatches
taken between hybrid waveforms that are generated with different methods. In the left panel, we adjust the amplitude of
restricted PN waveforms, such that they connect smoothly in amplitude to NR waveforms. In the right panel, to set the
frequency of PN waveforms at the joining point, we use the original orbital frequency, instead of the quartic fitted one. (See
Sec. IIIA for the discussion on amplitude scaling and frequency fitting).
4 6 8 10 12 14
Number of NR waveform cycles
Advanced LIGO PSD
VIRGO PSD
4 6 8 10 12 14
Number of NR waveform cycles
Advanced LIGO PSD
VIRGO PSD
FIG. 6: Mismatch between hybrid waveforms as a function of the number of NR waveform cycles used to generate the hybrid
waveforms. Following the settings of Fig. 5, we show comparisons between Epn(3.5) and Tpn(3.5), and Tpn(3) and Tpn(3.5)
models in the left and the right panels, respectively. The solid and dashed sets of curves are generated using the PSDs of
advanced LIGO and VIRGO. In each set, from top to bottom, the three curves correspond to three equal-mass binary systems,
with total masses 20M⊙, 30M⊙, and 40M⊙.
be reduced to below 1%.
• For binary systems with total mass lower than
10M⊙, the difference between the Tpn(3) and
Tpn(3.5) models is substantial for Advanced LIGO
and VIRGO. Their mismatch can be > 4% and
> 6% respectively (not shown in the figure). In this
mass range, pursuing more NR waveform cycles in
late inspiral phase does not help much, since the
signal power is accumulated slowly over hundreds
of GW cycles across the detector band. Neverthe-
less, here we give mismatches for FF0s which are
not optimized over binary masses. For the purpose
of detection only, optimization over binary param-
eters leads to low enough mismatches (see also the
end of Sec. II). In the language of Ref. [19] hybrid
waveforms for total mass lower than 10M⊙ are ef-
fectual but not faithful.
IV. MATCHING NUMERICAL WAVEFORMS
WITH POST-NEWTONIAN TEMPLATES
In this section, we compare the complete inspiral,
merger and ring-down waveforms of coalescing compact
binary systems generated from NR simulations with their
best-match PN template waveforms. We also compare
hybrid waveforms with PN template waveforms for lower
total masses, focusing on the late inspiral phase pro-
vided by the NR waveforms. We test seven families of
PN templates that either have been used in searches
for GWs in LIGO (see e.g., Refs. [40, 41]), or are
promising candidates for ongoing and future searches
with ground-based detectors. We evaluate the perfor-
mance of PN templates by computing the FFs maxi-
mized on phase, time and binary parameters. As we
shall see, for the hybrid waveforms of binary systems
with total mass M ≤ 30M⊙, both the time-domain fam-
ilies Tpn(3.5) and Epn(3.5), which includes a superposi-
tion of three ring-down modes, perform well, confirming
what found in Refs. [10, 11]. The standard stationary-
phase-approximated (SPA) template family in the fre-
quency domain has high FFs only for binary systems with
M < 20M⊙. After investigating in detail the GW phase
in frequency domain, and having understood why it hap-
pens (see Sec. IVB2), we introduce two modified SPA
template families (defined in Sec. IVB 2) for binary sys-
tems with total mass M ≥ 30M⊙. Overall, for masses
M ≥ 30M⊙, the Epn(3.5) template family in the time do-
main and the two modified SPA template families in the
frequency domain exhibit the best-match performances.
A. Numerical waveforms and post-Newtonian
templates
For binary systems with total mass M ≥ 30M⊙, the
last 8–16 cycles contribute more than 80–90% of the sig-
nal power, thus in this case we use only the NR wave-
forms. By contrast, for binary systems with total mass
10 ≤ M ≤ 30M⊙, for which the merger and ring-down
phases of the waveforms contribute only ∼ 1–10%, we
use the hybrid waveforms, generated by stitching Tpn
waveforms to the Goddard NR waveforms.
We want to emphasize that FFs computed for differ-
ent target numerical waveforms can not directly be com-
pared with each other. For instance, the Goddard wave-
form is longer than the Pretorius waveform, and the FFs
are sometime slightly lower using the Goddard waveform.
This is a completely artificial effect, due to the fact that
it is much easier to tune the template parameters and
obtain a large FF with a shorter target waveform than a
longer one.
We consider seven PN template families. The two
time-domain families introduced in Sec. II are:
• Tpn(3.5) [30, 31]: The inspiral Taylor model.
• Epn(3.5) [10, 13, 14, 16, 19]: The EOB model which
includes a superposition of three quasi-normal modes
(QNMs) of the final BH. These are labeled by three in-
tegers (l,m, n) [42]: the least damped QNM (2, 2, 0) and
two overtones (2, 2, 1) and (2, 2, 2). The ring-down wave-
form is given as:
hQNM(t) =
−(t−tend)/τ22n cos [ω22n(t− tend) + φn] ,
where ωlmn and τlmn are the frequency and decay time of
the QNM (l,m, n), determined by the mass Mf and spin
af of the final BH. The quantities An and φn in Eq. (7)
are the amplitude and phase of the QNM (2, 2, n). They
are obtained by imposing the continuity of h+ and h×,
and their first and second time derivatives, at the time of
matching tmatch. Besides the mass parameters, our Epn
model contains three other physical parameters: ǫt, ǫM
and ǫJ . The parameter ǫt takes into account possible dif-
ferences between the time tend at which the EOB models
end and the time tmatch at which the matching to ring-
down is done. More explicitly, we set tmatch = (1+ǫt)tend,
and if ǫt > 0, we extrapolate the EOB evolution, and set
an upper limit for the ǫt search where the extrapolation
fails. The parameters ǫM and ǫJ describe possible differ-
ences between the values of the mass Mend ≡ Eend and
angular momentum âend ≡ Jend/M
end at the end of the
EOB inspiral and the final BH mass and angular mo-
mentum. (The end of the EOB inspiral occurs around
the EOB light-ring.) The differences are due to the fact
that the system has yet to release energy and angular
momentum during the merger and ring-down phase be-
fore settling down to the stationary BH solution. If the
total binary mass and angular momentum at the end of
the EOB inspiral are Mend and Jend, we set the total
mass and angular momentum of the final stationary BH
to be Mf = (1 − ǫM )Mend and Jf = (1 − ǫJ)Jend, and
use af ≡ Jf/Mf to compute ωlmn and τlmn. We consider
the current Epn model with three parameters ǫt, ǫM and
ǫJ , as a first attempt to build a physical EOB model for
matching coherently the inspiral, merger and ring-down
phases. Since the ǫ-parameters are related to physical
quantities, e.g., the loss of energy during ring-down, they
are functions of the initial physical parameters of the bi-
nary, such as masses, spins, etc. In the near future we
expect to be able to fix the ǫ-values by comparing NR and
(improved) EOB waveforms for a large range of binary
parameters.
We also consider five frequency-domain models, in
which two (modified SPA models) are introduced later
in Sec. IVB 2, and three are introduced here:
• SPAc(3.5) [33]: SPAc PN model with an appropriate
cutoff frequency fcut [30, 31];
(5 + 5)M⊙ (10 + 10)M⊙ (15 + 15)M⊙
Signal Power (%) (30, 0.2) (80, 2) (85, 10)
〈hNR−hybr, hTpn(3.5)〉 0.9875 0.9527 0.8975
(M/M⊙, η) (10.18, 0.2422) (19.97, 0.2500) (29.60, 0.2499)
Mωorb 0.1262 0.1287 0.1287
〈hNR−hybr, hEpn(3.5)〉 0.9836 0.9522 0.9618
(M/M⊙, η) (10.15, 0.2435) (19.90, 0.2500) (29.49,0.2488)
(ǫt, ǫM , ǫJ )(%) (-0.02, 12.19, 30.87) (-0.02, 75.03, 95.00) (0.05, 2.38, 92.06)
Mωorb 0.1346 0.1345 0.1345
〈hNR−hybr, hSPAc(3.5)〉 0.9690 0.9290 0.8355
(M/M⊙, η) (10.16, 0.2432) (19.93, 0.2498) (29.08, 0.2500)
(fcut/Hz) 1566.8 263.9 529.6
TABLE I: FFs between hybrid waveforms [Tpn(3.5) waveform stitched to the Goddard waveform] and PN templates. In the first
row, the two numbers in parentheses are the percentages of the signal-power contribution from the 16 inspiraling NR cycles and
the NR merger/ring-down cycles. (The separation between inspiral and merger/ring-down is obtained using the EOB approach
as a guide, i.e., we match the Epn(3.5) model and use the EOB light-ring position as the beginning of the merger phase.) In
the PN-template rows, the first number in each block is the FF, and the numbers in parentheses are template parameters that
achieve this FF. The last number in each block of the Tpn(3.5) and Epn(3.5) models is the ending orbital frequency of the
best-match template. For the Epn model, the ending frequency is computed at the point of matching with the ring-down phase,
around the EOB light ring.
• BCV [31]: BCV model with an amplitude correction
term (1−αf2/3) and an appropriate cutoff frequency fcut.
• BCVimpr [31]: Improved BCV model with an am-
plitude correction term (1 − αf1/2) and an appropri-
ate cutoff frequency fcut. We include this improved
BCV model because Ref. [10] found a deviation of the
Fourier-transform amplitude from the Newtonian predic-
tion f−7/6 during the merger and ring-down phases (see
Fig. 22 of Ref. [10]). Here we shall assume n = −2/3
in the fn power law to get the (1 − αf1/2) form of the
amplitude correction. While it was found [10] that the
value of n is close to −2/3 for the l = 2,m = 2 wave-
form, this value varies slightly if other multiple moments
are included and if binary systems with different mass ra-
tios are considered. Finally, the α parameter is expected
to be negative, but in our actual search it can take both
positive and negative values.
B. Discussion of fitting-factor results
In Table I, we list the FFs for hybrid target waveforms
and three PN template families: Tpn(3.5), Epn(3.5),
and SPAc(3.5), together with the template parameters at
which the best match is obtained. As shown in the first
row, in this relatively low-mass range, i.e. 10M⊙ < M <
30M⊙, the merger/ring-down phases of the waveforms
contribute only a small fraction of the total signal power,
while the last 16 inspiraling cycles of the NR waveform
contribute a significant fraction. Therefore, confirming
recent claims by Refs. [10, 11], we can conclude that
the PN template families Tpn(3.5) and Epn(3.5) have
good agreement with the inspiraling NR waveforms. The
Tpn(3.5) model gives a low FF for M = 30M⊙ because
0.0 0.2 0.4 0.6 0.8 1.0
af /Mf
FIG. 7: Frequencies and decay times of the least damped
QNM 220, and two overtones 221 and 222. The scales of the
frequency and the decay time are listed on the left and right
sides of the plot, respectively.
for these higher masses the merger/ring-down phases,
which the Tpn model does not include, start contribut-
ing to the signal power. Note that both time-domain
templates give fairly good estimates of the mass parame-
ters. The SPAc(3.5) template family gives FFs that drop
substantially when the total binary mass increases from
10M⊙ to 30M⊙, indicating that this template family can
only match the early, less relativistic inspiral phase of
the hybrid waveforms. Nevertheless, it turns out that by
slightly modifying the SPA waveform we can match the
NR waveforms with high FFs (see Sec. IVB2).
In Table II, we list the FFs for full NR waveforms
and five PN template families: Epn(3.5), SPAextc (3.5),
SPAYc (4), BCV, and BCVimpr, together with the tem-
(15 + 15)M⊙ (20 + 20)M⊙ (30 + 30)M⊙ (50 + 50)M⊙
〈hNR−Pretorius, hEpn(3.5)〉 0.9616 0.9599 0.9602 0.9787
(M/M⊙, η) (27.93, 0.2384) (35.77, 0.2426) (52.27, 0.2370) (96.60, 0.2386)
(ǫt, ǫM , ǫJ )(%) (-0.08, 0.63, 99.70) (-0.03, 0.48, 94.38) (-0.12, 0.00, 64.14) (0.04, 0.01, 73.01)
〈hNR−Pretorius, hSPA
c (3.5)〉 0.9712 0.9802 0.9821 0.9722
(M/M⊙, η) (19.14, 08037) (24.92, 0.9097) (36.75, 0.9933) (58.06, 0.9986)
(fcut/Hz) (589.6) (476.9) (318.9) (195.9)
〈hNR−Pretorius, hSPA
c (4)〉 0.9736 0.9824 0.9874 0.9851
(M/M⊙, η) (29.08, 0.2460) (38.63, 0.2461) (57.58, 0.2441) (96.55, 0.2457)
(fcut/Hz) (666.5) (501.2) (332.5) (199.4)
〈hNR−Pretorius, hBCV〉 0.9726 0.9807 0.9788 0.9662
(ψ0/10
4, ψ1/10
2) (2.101, 1.655) (1.178, 1.744) (0.342, 2.385) (-0.092, 3.129)
(102α, fcut/Hz) (-1.081, 605.5) (-0.834, 461.7) (0.162, 320.4) (1.438, 204.3)
〈hNR−Pretorius, hBCVimpr〉 0.9727 0.9807 0.9820 0.9803
(ψ0/10
4, ψ1/10
2) (2.377, 0.930) (1.167, 1.762) (0.431, 2.077) (-0.109, 3.158)
(102α, fcut/Hz) (-3.398, 571.9) (-2.648, 458.3) (-1.196, 319.1) (-3.233, 196.0)
(15 + 15)M⊙ (20 + 20)M⊙ (30 + 30)M⊙ (50 + 50)M⊙
〈hNR−Goddard, hEpn(3.5)〉 0.9805 0.9720 0.9692 0.9671
(M/M⊙, η) (29.25, 0.2435) (38.27, 0.2422) (56.66, 0.2381) (83.52, 0.2233)
(ǫt, ǫM , ǫJ )(%) (0.05, 0.03, 99.90) (0.05, 0.27, 99.17) (0.09, 0.01, 54.56) (0.10, 1.71, 79.75)
〈hNR−Goddard, hSPA
(3.5)〉 0.9794 0.9785 0.9778 0.9693
(M/M⊙, η) (21.41, 0.5708) (27.27, 0.6695) (37.67, 0.9911) (60.90, 0.9947)
(fcut/Hz) (552.7) (444.4) (318.5) (191.7)
〈hNR−Goddard, hSPA
(4)〉 0.9898 0.9905 0.9885 0.9835
(M/M⊙, η) (30.28, 0.2456) (40.23, 0.2477) (60.54, 0.2455) (100.00, 0.2462)
(fcut/Hz) (674.6) (506.6) (330.5) (195.0)
〈hNR−Goddard, hBCV〉 0.9707 0.9710 0.9722 0.9692
(ψ0/10
4, ψ1/10
2) (3.056, -1.385) (1.650, -0.091) (0.561, 1.404) (-0.113, 3.113)
(102α, fcut/Hz) (0.805, 458.3) (0.559, 412.6) (0.218, 309.2) (1.063, 198.7)
〈hNR−Goddard, hBCVimpr〉 0.9763 0.9768 0.9782 0.9803
(ψ0/10
4, ψ1/10
2) (2.867, -0.600) (1.514, 0.448) (0.555, 1.425) (-0.165, 3.373)
(102α, fcut/Hz) (0.193, 578.0) (-1.797, 441.1) (-4.472, 308.1) (-4.467, 193.4)
TABLE II: FFs between NR waveforms and PN templates which include merger and ring-down phases. The upper table uses
Pretorius’ waveform, and the lower table uses Goddard’s high-resolution long waveform. The first number in each block is the
FF, and numbers in parentheses are template parameters that achieve this FF.
plate parameters at which the best match is obtained.
The SPAextc (3.5) and SPA
c (4) families are modified ver-
sions of the SPA family, defined in Sec. IVB 2.
We shall investigate these results in more detail in the
following sections.
1. Effective-one-body template performances
The Epn model is the only available time-domain
model that explicitly includes ring-down waveforms. It
achieves high FFs ≥ 0.96 for all target waveforms, con-
firming the necessity of including ring-down modes and
proving that the inclusion of three QNMs with three
tuning parameters ǫt, ǫM and ǫJ is sufficient for detec-
tion. As we see in Table II, the values of the tuning
parameters ǫM and ǫJ , where the FFs are achieved, are
different from their physical values. For reference, the
Goddard numerical simulation predicts Mf ≃ 0.95M
and âf ≡ Jf/M
f ≃ 0.7 [38], and Epn(3.5) predicts
Mend = 0.967 and âend ≡ Jend/M
end = 0.796, so the two
tuning parameters should be ǫM ≃ 1.75% and ǫJ ≃ 11%.
In our search, e.g., for M = 30M⊙, ǫJ tends to be tuned
to its lowest possible value and ǫt tends to take its high-
est possible value, indicating that pushing the end of the
Epn(3.5) inspiral to a later time gives higher FFs.
Since the parameters ǫM and ǫJ depend on the QNM
frequency and decay time, we show in Fig. 7 how ωlmn
and τlmn vary as functions of af [42] for the three modes
used in the Epn(3.5) model. The frequencies ωlmn of the
three modes are not really different, and grow monotoni-
cally with increasing af . The decay times τlmn, although
different for the three modes, also grow monotonically
with increasing af . Thus, the huge loss of angular mo-
mentum ǫJ , or equivalently the small final BH spin re-
quired in the Epn(3.5) model to achieve high FFs, indi-
cates that low ring-down frequencies and/or short decay
times are needed for this model to match the numerical
0.10 0.12 0.14 0.16 0.18
Time (s)
-0.02
NR frequency & waveform
Epn frequency & waveform
(15+15)
0.40 0.45 0.50 0.55 0.60
Time (s)
(50+50)
FIG. 8: Frequency evolution of waveforms from the Epn(3.5) model, and the NR simulations of the Goddard group. In the left
and right panels, we show frequency evolutions for two equal-mass binary systems with total mass 30M⊙ and 100M⊙. In each
panel, there are two nearly monotonic curves and two oscillatory curves, where the former are frequency evolutions and the
latter are binary coalescence waveforms. The solid curves (blue) are from the NR simulations, while the dashed curves (red)
are from the Epn(3.5) model. The vertical line in each plot shows the position where the three-QNM ring-down waveform is
attached to the EOB waveform.
merger and ring-down waveforms.
In Fig. 8, we show Goddard NR and Epn(3.5) wave-
forms, as well as their frequency evolutions, for two
equal-mass binary systems with total masses 30M⊙ and
100M⊙. In the low-mass case, i.e., M = 30M⊙, since the
inspiral part contributes most of the SNR, the Epn(3.5)
model fits the frequency and phase evolution of the NR
inspiral well, with the drawback that at the joining point
the EOB frequency is substantially higher than that of
the NR waveform. Then, in order to fit the early ring-
down waveform which has higher amplitude, the tuning
parameters have to take values in Table II such that the
ring-down frequency is small enough to get close to the
NR frequency during early ring-down stage, as indicated
in Fig. 8. The late ring-down waveform does not con-
tribute much to the SNR, and thus it is not too sur-
prising that waveform optimizing the FF does not ade-
quately represent this part of the NR waveform. In the
higher mass case,M = 100M⊙, the Epn(3.5) model gives
a much better, though not perfect, match to the merger
and ring-down phases of the NR waveform, at the ex-
pense of misrepresenting the early inspiral part. Again,
this is not unexpected considering that in this mass range
the merger and ring-down waveforms dominate the con-
tribution to the SNR.
Comparing the two cases discussed above, we can see
that with the current procedure of matching the inspiral
and ring-down waveforms in the EOB approach it is not
possible to obtain a perfect match with the entire NR
waveform. However, due to the limited detector sensi-
tivity bandwidth, the FFs are high enough for detection.
The large systematic error in estimating the physical pa-
rameters will be overcome by improving the EOB match-
ing procedure during the inspiral part, and also by fixing
the ǫ-parameters to physical values obtained by compar-
ison with numerical simulations.
Finally, in Figs. 9, 10 we show the frequency-domain
amplitude and phase of the NR and EOB waveforms.
Quite interestingly, we notice that the inclusion of three
ring-down modes reproduce rather well the bump in the
NR frequency-domain amplitude. The EOB frequency-
domain phase also matches the NR one very well.
2. Stationary-phase-approximated template performances
Figures 9, 10 also show the frequency-domain phases
and amplitudes for the best-match SPAc(3.5) waveforms.
We see that at high frequency the NR and SPAc(3.5)
phases rise with different slopes 3. Based on this obser-
vation we introduce two modified SPA models:
• SPAextc (3.5): SPAc PN model with unphysical values
of η and an appropriate cutoff frequency fcut. The range
of the symmetric mass-ratio η = m1m2/(m1+m2)
2 is ex-
tended from its physical range 0 ∼ 0.25 to the unphysical
3 By looking in detail at the PN terms in the SPAc(3.5) phase, we
find that the difference in slope is largely due to the logarithmic
term at 2.5PN order.
100 200 500
NR phase
SPAc(3.5)
ext(3.5)
Epn(3.5)
100 200 500
Frequency (Hz)
NR amplitude
Epn(3.5)
FIG. 9: For M = 30M⊙ equal-mass binary systems, we com-
pare the phase and amplitude of the frequency-domain wave-
forms from the SPAc models and NR simulation (Goddard
group). We also show the amplitude of the waveform from
the Epn(3.5) model.
range 0 ∼ 1.
• SPAYc (4): SPAc PN model with an ad hoc 4PN order
term in the phase, and an appropriate cutoff frequency
fcut. The phase of the SPA model is known up to the
3.5PN order (see, e.g., Eq. (3.3) of Ref. [33]):
ψ(f) = 2πft0 − φ0 −
128ηv5
k , (8)
where v = (πMf)1/3. The PN coefficients αks, k =
0, . . . , N , (with N = 7 at 3.5PN order) are given by Eqs.
(3.4a), (3.4h) of Ref. [33]. We add the following term at
4PN order:
α8 = Y log v , (9)
where Y is a parameter which we fix by imposing high
matching performances with NR waveforms. Note that a
constant term in α8 only adds a 4PN order term that is
linear in f , which can be absorbed into the 2πft0 term.
Thus, to obtain a nontrivial effect, we need to introduce
a logarithmic term. The coefficient Y could in princi-
ple depend on η. We determine Y by optimizing the
FFs of equal and unequal masses. We find that in the
equal-mass case Y does not depend significantly on the
binary total mass and is given by Y = 3923. The latter
is also close to the best match value obtained for unequal
masses. More specifically, it is within 4.5% for binary sys-
tems of mass ratiom2/m1 = 2. To further explore the de-
pendence of Y on η, we need a larger sample of waveforms
50 100 200
NR phase
SPAc(3.5)
ext(3.5)
Epn(3.5)
M=100
50 100 200
Frequency (Hz)
NR amplitude
Epn(3.5)
FIG. 10: For M = 100M⊙ equal-mass binary systems, we
compare the phase and amplitude of the frequency-domain
waveforms from the SPAc models and NR simulation (God-
dard group). We also show the amplitude of the waveform
from the Epn(3.5) model.
for unequal-mass binary systems 4. As seen in Table II,
the two modified SPAc template families have FF> 0.97
(except for one 0.9693) for all target waveforms, even
though no explicit merger or ring-down phases are in-
cluded in the waveform. The SPAYc (4) model provides
also a really good estimation of parameters.
In Fig. 11 we plot Goddard NR and SPAYc (4) wave-
forms for two equal-mass binary systems with total
massesM = 30M⊙ andM = 100M⊙. We can clearly see
tail-like ring-down waveforms at the end of the SPAYc (4)
waveforms, which result from the inverse Fourier trans-
form of frequency domain waveforms that have been cut
at f = fcut. This well-known feature is called the Gibbs
phenomenon. At first glance, it may appear surpris-
ing that the often inconvenient Gibbs phenomenon [44]
can provide reasonable ring-down waveforms in the time
domain. However, by looking at the spectra of these
waveforms in the frequency domain (see the amplitudes
in Figs. 9 & 10), we see that the SPAYc (4) cuts off at
the frequency fcut (obtained from the optimized FF)
where the NR spectra also start to drop. Thus, even
though the frequency-domain SPAc waveforms are dis-
continuous, while the frequency-domain NR waveforms
are continuous (being combinations of Lorentzians), the
SPAc time-domain waveforms contain tails with frequen-
cies and decay rates similar to the NR ring-down modes.
We expect that the values of the cutoff frequency fcut
at which the FFs are maximized are well-determined by
4 Note that the auxiliary phase introduced in Eq. (239) of Ref. [43]
also gives rise to a term in the SPA phase of the kind f log v,
except an order of magnitude smaller than Y .
the highest frequency of the NR waveforms, i.e. by the
frequency of the fundamental QNM. In the next section,
we shall show quantitative results to confirm this guess.
3. Buonanno-Chen-Vallisneri template performances
In Table II we see that the BCV and BCVimpr fam-
ilies give almost the same FFs for relatively low-mass
binary systems (M = 30, 40M⊙), while the BCVimpr
family gives slightly better FFs for higher mass binary
systems (M = 60, 100M⊙). For higher-mass binary sys-
tems, we find that the α parameter takes negative values
with reasonable magnitude. This is because the ampli-
tude of the NR waveforms in the frequency domain de-
viates from the f−7/6 power law only near the merger,
which lasts for about one GW cycle. This merger cycle
is important only when the total mass of the binary is
high enough (see Fig. 4). [See also Ref. [46] where similar
tests have been done.]
The BCV and BCVimpr template families give FFs
nearly as high as those given by the SPAYc (4) family,
but the latter has the advantage of being parametrized
directly in terms of the physical binary parameters, and
it gives fairly small systematic errors.
C. Frequency-domain templates for inspiral,
merger and ring-down
In this section, we extend our comparisons between
the SPAc families and NR waveforms to higher total-
mass binary systems (40M⊙ to 120M⊙) and to unequal-
mass binary systems with mass-ratios m2/m1 = 1.5 and
2. The numerical simulations for unequal-mass binary
systems are from the Goddard group. They last for ≃
373M and ≃ 430M , respectively, and the NR waveforms
have ≃ 4 cycles before the merger.
In Figs. 12 and 13 we show the FFs for SPAextc (3.5)
and SPAYc (4) templates, and the values of fcut that
achieved these FFs 5. For all mass combinations (ex-
cept for M = 40M⊙ for artificial reasons) the FFs of
SPAc(3.5) templates are higher than 0.96, and the FFs
of SPAYc (4) templates are higher than 0.97, confirming
that both families of templates can be used to search for
GWs from coalescing binary systems with equal-masses
as large as 120M⊙ and mass ratios m2/m1 = 2 and 1.5.
Figure 13 shows that all the fcut values from our searches
5 Note that because of the short NR waveforms for unequal-mass
binary systems, we need to search over the starting frequency of
templates with a coarse grid, and this causes some oscillations in
our results. The oscillations are artificial and will be smoothed
out in real searches. For instance, the drop of FFs at 40M⊙
for unequal-mass binary systems happens because the NR wave-
forms are too short and begin right in the most sensitive band
of LIGO.
are within 10% larger than the frequency of the funda-
mental QNM ω220 of an equal-mass binary. We have
checked that if we fix fcut = 1.07ω220/2π, the FFs drop
by less than 1%.
In Fig. 15, we show the same information as in Fig. 7,
except that here we draw ωlnm and τlnm as functions of
the mass-ratio η of a nonspinning binary. We compute
the spin of the final BH in units of the mass of the fi-
nal BH using the quadratic fit given by Eq. (3.17a) of
Ref. [47]:
≃ 3.352η − 2.461η2 . (10)
As Fig. 15 shows, ω220 does not change much, confirming
the insensitivity of the fcut on η.
However, in real searches we might request that the
template family have some deviations from the wave-
forms predicted by NR. For example, a conservative tem-
plate bank might cover a region of fcut ranging from the
Schwarzschild innermost stable circular orbit (ISCO) fre-
quency, or the innermost circular orbit (ICO) frequency
determined by the 3PN conservative dynamics, up to a
value slightly higher than the frequency of the fundamen-
tal QNM. The number of templates required to cover
the fcut dimension depends on the binary masses. We
find that to cover the fcut dimension from the 3PN ICO
frequency to the fundamental QNM frequency with an
SPAextc (3.5) template bank, imposing a mismatch < 0.03
between neighboring templates, we need only two (∼20)
templates if M = 30M⊙ (M = 100M⊙) and η = 0.25.
In the latter case, the match between templates is more
sensitive to fcut since most signal power comes from
the last two cycles, sweeping through a large frequency
range, right in LIGO’s most sensitive band. The num-
ber of templates directly affects the computational power
needed, and the false-alarm rate. Further investigations
are needed in order to determine the most efficient way
to search over the fcut dimension.
For the purpose of parameter estimation, Fig. 14 shows
that the SPAYc (4) templates are rather faithful, giving
reasonable estimates of the chirp mass: systematic errors
less than about 8% in absolute value for binary systems
with M = 40M⊙ up to M = 120M⊙. A difference of
≃ 8% may seem large, but the SPAYc (4) templates are
not exactly physical, and more importantly, for large-
mass binary systems, most of the information on the
chirp mass comes only from the last cycle of inspiral. We
notice that when the total binary mass is higher than
120M⊙, the FFs are relatively high (from 0.93 to 0.97),
and the estimates of the chirp mass are still good (within
10%). However, for binary systems with such high to-
tal masses, the ring-down waveform dominates the SNR,
and the SPAYc (4) template family becomes purely phe-
nomenological. A direct ring-down search might be more
efficient.
All results for unequal-mass binary systems are ob-
tained using the C22 component of Ψ4 [10], which is
the leading order quadrupole term contributing to the
0.06 0.08 0.10 0.12 0.14 0.16 0.18
Time (s)
NR waveform
(15+15)
0.45 0.50 0.55 0.60
Times (s)
(50+50)
FIG. 11: Binary coalescence waveforms from the SPAYc (4) model, and the NR simulations of the Goddard group. In the left
and right panels we show waveforms for two equal-mass binary systems with total mass 30M⊙ and 100M⊙. The solid lines
show the waveforms from the NR simulation, and the dashed lines give the best-matching waveforms from the SPAYc (4) model.
40 60 80 100 120
Total mass ( )
ext(3.5): equal-mass
ext(3.5): mass-ratio 2:1
ext(3.5): mass-ratio 3:2
(4): equal-mass
(4): mass-ratio 2:1
(4): mass-ratio 3:2
FIG. 12: FFs as functions of the total binary mass. The FFs
are computed between either the SPAextc (3.5) or the SPA
c (4)
templates and the NR waveforms for equal-mass and unequal-
mass binary systems.
GW radiation. For unequal-mass binary systems, higher-
order multipoles can also be important, and we need to
test the performance of the template family directly using
Ψ4. For Ψ4 extracted in the direction perpendicular to
the binary orbit, we verified that higher-order multipoles
do not appreciably change the FFs.
A natural way of improving the SPAc models would be
to replace the discontinuous frequency cut with a linear
combination of Lorentzians. We show here a first attempt
40 60 80 100 120
Total mass ( )
ext(3.5): equal-mass
ext(3.5): mass-ratio 2:1
ext(3.5): mass-ratio 3:2
(4): equal-mass
(4): mass-ratio 2:1
(4): mass-ratio 3:2
Frequency of the fundamental QNM
FIG. 13: Cutoff frequencies as functions of the total binary
mass. We show the best-match fcut for SPA
c (3.5) and
SPAYc (4) templates of Fig. 12. The solid black curve is the
fundamental QNM frequency ω220/2π. The frequencies are in
units of Hz.
at doing so. The Lorentzian L is obtained as a Fourier
transform of a damped sinusoid, e.g., for the fundamental
QNM we have
ei2πft
e±iω220t−|t|/τ220
2/τ220
1/τ2220 + (2πf ± ω220)
≡ 2L±220(f) (11)
40 60 80 100 120
Total mass ( )
-0.08
-0.06
-0.04
-0.02
(4): equal-mass
(4): mass-ratio 2:1
(4): mass-ratio 3:2
FIG. 14: Systematic errors of the chirp mass as functions
of the total binary mass when SPAYc (4) templates are used.
We show errors of the chirp masses that optimize the FFs of
Fig. 12.
0.00 0.05 0.10 0.15 0.20 0.25
FIG. 15: Frequencies and decay times of the least damped
QNM 220, and two overtones 221 and 222. The scales of the
frequency and the decay time are listed on the left and right
sides of the plot, respectively.
and the (inverse) Fourier transform of Eq. (7) reads
h̃QNM(f) = An
L+22n(f) e
iφn + L−22n(f) e
Restricting to positive frequencies we only keep the
L−22n(f) terms. In the frequency domain we attach the
fundamental mode continuously to the SPAYc (4) wave-
form at the ring-down frequency ω220 by tuning the am-
plitude and phase A0 and φ0. We denote this model
SPAL1 (note that we also need to introduce the mass-
parameter of the final BH as a scale for ω220 and τ220).
Similarly, we define the SPAL3 model where all three
QNMs are combined. With the three amplitudes and
phases as parameters, this model is similar to the spin-
BCV template family [30] and we can optimize automat-
ically over the 6 parameters. As an example, we compute
the FFs between the SPAL1(4) or SPAL3(4) and the NR
waveform of an equal-mass M = 100M⊙ binary. Using
the LIGO PSD, we obtain 0.9703 and 0.9817, respec-
tively. Those FFs are comparable to the FFs obtained
with the simpler SPAc model, shown in Fig. 12. It is
known that adding more parameters increases the FFs
but also increases the false-alarm probability. By fur-
ther investigation and comparison with NR waveforms
our goal is to express the phase and amplitude parame-
ters of the Lorentzian in terms of the physical binary pa-
rameters, relating them to the amplitudes and phases of
the QNMs and the physics of the merger. Those parame-
ters are somewhat similar to the ǫ-parameters introduced
above for the EOB model when modeling the merger and
ring-down phases.
We wish to emphasize that the results we presented in
this section are preliminary, in the sense that we consid-
ered only a few mass combinations and the NR waveforms
of unequal-mass binary systems are quite short. Never-
theless, these results are interesting enough to propose a
systematic study of the efficiency of these template fam-
ilies through Monte Carlo simulations in real data.
V. CONCLUSIONS
In this paper we compared NR and analytic waveforms
emitted by nonspinning binary systems, trying to under-
stand the performance of PN template families developed
during the last ten years and currently used for the search
for GWs with ground-based detectors, suggesting possi-
ble improvements.
We first computed FF0s (maximized only on time
and phase) between PN template families which best
match NR waveforms [10, 11], i.e., Tpn(3), Tpn(3.5) and
Epn(3.5). We showed how the drop in FF0s is not simply
determined by the accumulated phase difference between
waveforms, but also depends on the detector’s PSD and
the binary mass. Thus, waveforms which differ even by
one GW cycle can have FF0 ∼ 0.97, depending on the
binary masses (see Fig. 1).
We then showed that the NR waveforms from the
high-resolution and medium-resolution simulations of the
Goddard group are close to each other (FF0 around 0.99,
see Fig. 2). We also estimated that the FF0 between
high-resolution and exact NR waveforms is even higher,
based on the numerical convergence rates of the Goddard
simulations.
Second, by stitching PN waveforms to NR waveforms
we built hybrid waveforms, and computed FF0s (max-
imized only on time and phase) between hybrid wave-
forms constructed with different PN models, notably
Tpn(3), Tpn(3.5) and Epn(3.5) models. We found that
for LIGO’s detectors and equal-mass binary systems with
total mass M > 30M⊙, the last 11 GW cycles plus
merger and ring-down phases contribute > 90% of the
signal power. This information can be used to set the
length of NR simulations.
The FF0s between hybrid waveforms are summarized
in Figs. 5, 6. We found that for LIGO’s detectors and
binary systems with total mass higher than 10M⊙, the
current NR simulations for equal-mass binary systems
are long enough to reduce the differences between hybrid
waveforms built with the PN models Tpn(3), Tpn(3.5)
and Epn(3.5) to the level of < 3% mismatch. For GW
detectors with broader bandwidth like advanced LIGO
and VIRGO, longer NR simulations will be needed if the
total binary masses M < 10M⊙. With the current avail-
able length of numerical simulations, it is hard to esti-
mate from the FFs between hybrid waveforms how long
the simulations should be. Nevertheless, from our study
of the distribution of signal power, we estimate that for
M < 10M⊙ binary systems, at least ∼ 80 NR inspiraling
cycles before merger are needed.
Finally, we evaluated FFs (maximized on binary
masses, initial time and phase) between full NR (or hy-
brid waveforms, depending on the total binary mass)
and several time and frequency domain PN template
families. For time-domain PN templates and binary
masses 10M⊙ < M < 20M⊙, for which the merger/ring-
down phases do not contribute significantly to the to-
tal detector signal power, we confirm results obtained in
Refs. [10, 11], notably that Tpn(3.5) and Epn(3.5) mod-
els have high FFs with good parameter estimation, i.e.,
they are faithful. We found that the frequency-domain
SPA family has high FFs only for binary systems with
M < 20M⊙, for which most of the signal power comes
from the early stages of inspiral. Furthermore, we found
that it is possible to improve the SPA family by either
extending it to unphysical regions of the parameter space
(as done with BCV templates) or by introducing an ad
hoc 4PN-order constant coefficient in the phase. Both
modified SPA families achieve high FFs for high-mass
binary systems with total masses 30M⊙ < M < 120M⊙.
For time-domain PN templates and binary masses
M >∼ 30M⊙, we found that if a superposition of ring-
down modes is attached to the inspiral waveform, as nat-
urally done in the EOB model, the FFs can increase from
∼ 0.8 to > 0.9. We tested the current Epn(3.5) template
family obtained by attaching to the inspiral waveform
three QNMs [10] around the EOB light-ring. In order
to properly take into account the energy and angular-
momentum released during the merger/ring-down phases
we introduced [10] two physical parameters, ǫM and ǫJ ,
whose dependence on the binary masses and spins will
be determined by future comparisons between EOB and
NR waveforms computed for different mass ratios and
spins. We found high FFs >∼0.96. Due to small differ-
ences between EOB and NR waveforms during the fi-
Total mass ( )
50 100 150 200
Enhanced LIGO
Template: NR
Template: SPAc
ext(3.5)
Template: Epn(3.5)
FIG. 16: The sky-average SNR for LIGO and Enhanced or
mid LIGO detector versus total mass for an equal-mass binary
at 100Mpc.
nal cycles of the evolution, the best-matches are reached
at the cost of large systematic error in the merger–ring-
down binary parameters. Thus, the Epn(3.5) template
family can be used for detection, but for parameter es-
timation it needs to be improved when matching to the
ring-down, and also during the inspiral phase. The re-
finements can be achieved (i) by introducing deviations
from circular motion, (ii) adding higher-order PN terms
in the EOB dynamics, (iii) using in the EOB radiation-
reaction equations a GW energy flux closer the the NR
flux, (iv) designing a better match to ring-down modes,
etc.. The goal would be to achieve dephasing between
EOB and NR waveforms of less than a few percent in
the comparable-mass case, as obtained in Ref. [48] in the
extreme mass-ratio limit. Indeed, with more accurate nu-
merical simulations, especially those using spectral meth-
ods [49], it will be possible to improve the inspiraling
templates by introducing higher-order PN terms in the
analytic waveforms computed by direct comparison with
NR waveforms.
Frequency-domain PN templates with an appropriate
cutoff frequency fcut provide high FFs (> 0.97), even
for large masses. This is due to oscillating tails (Gibbs
phenomenon) produced when cutting the signal in the
frequency domain. We tested the SPAextc (3.5) and the
SPAYc (4) template families for total masses up to 120M⊙,
and three mass ratios m2/m1 = 1, 1.5, and 2. We always
get FFs > 0.96, even when using a fixed cutoff frequency,
fcut = 1.07ω220/2π. Because of its high efficiency, faith-
fulness, i.e., low systematic error in parameter estima-
tion, and simple implementation, the SPAYc (4) template
family (or variants of it which include Lorentzians) is, to-
gether with the EOB model, a good candidate for search-
ing coherently for GWs from binary systems with total
masses up to 120M⊙.
In Fig. 16, we show the sky averaged SNRs of a sin-
gle LIGO and Enhanced or mid LIGO [50] detector,
for an equal-mass binary at 100Mpc. The SNR peaks
at the total binary mass M ≃ 150M⊙ and shows the
importance of pushing current searches for coalescing
binary systems to M > 100M⊙. In the mass range
30M⊙ < M < 120M⊙, the SNR drops only slightly if
we filter the GW signal with SPAextc (3.5) or Epn(3.5)
instead of using NR waveforms. The difference be-
tween Epn(3.5) and SPAextc (3.5) is almost indistinguish-
able. When M > 120M⊙, although the SPA
c (3.5) and
Epn(3.5) template families give fairly good SNRs, it is
maybe not a good choice to use them as the number of
cycles reduces to a few. The key problem in detecting
such GWs is how to veto triggers from non-Gaussian,
nonstationary noise, instead of matching the effectively
short signal. This is a general problem in searches for
short signals in ground-based detectors.
Acknowledgments
A.B. and Y.P. acknowledge support from NSF grant
PHY-0603762, and A.B. also from the Alfred Sloan Foun-
dation. The work at Goddard was supported in part
by NASA grants O5-BEFS-05-0044 and 06-BEFS06-19.
B.K. was supported by the NASA Postdoctoral Program
at the Oak Ridge Associated Universities. S.T.M. was
supported in part by the Leon A. Herreid Graduate Fel-
lowship. Some of the comparisons with PN and EOB
models were obtained building on Mathematica codes de-
veloped in Refs. [14, 30, 31, 45]
APPENDIX A: COMMENT ON WAVEFORMS
OBTAINED FROM THE ENERGY-BALANCE
EQUATION
In adiabatic PN models, like the Tpn model used in
this paper, waveforms are computed under the assump-
tion that the binary evolves through an adiabatic se-
quence of quasi-circular orbits. More specifically, one
sets ṙ = 0 and computes the orbital frequency ω from the
energy-balance equation dE(ω)/dt = F(ω), where E(ω)
is the total energy of the binary system and F(ω) is the
GW energy flux. Both E(ω) and F(ω) are computed
for circular orbits and expressed as a Taylor expansion
in ω. The adiabatic evolution ends in principle at the
innermost circular orbit (ICO) [35], or minimum energy
circular orbit (MECO) [30], where (dE/dω) = 0.
By rewriting the energy-balance equation, ω(t) can be
integrated directly as
ω̇(t) =
dE(ω)/dω
. (A1)
The RHS of Eq. (A1) can be expressed as an expansion
in powers of ω. The expanded version is widely used in
generating adiabatic PN waveforms [20, 30, 31, 45], it
is used to generate the so-called Tpn template family. It
turns out that Tpn(3) and Tpn(3.5) are quite close to the
NR inspiraling waveforms [10, 11]. We wonder whether
using the energy-balance in the form of Eq. (A1), i.e.,
without expanding it, might give PN waveforms closer to
or farther from NR waveforms. In principle the adiabatic
sequence of circular orbits described by Eq. (A1) ends at
the ICO, so the adiabatic model should work better until
the ICO and start deviating (with ω going to infinity)
from the exact result beyond it.
In Fig. 17 we show the NR orbital frequency ω(t) to-
gether with the PN orbital frequency obtained by solv-
ing the unexpanded and expanded form of the energy-
balance equation. The frequency evolution in these two
cases is rather different, with the orbital-frequency com-
puted from the expanded energy-balance equation agree-
ing much better with the NR one. When many, extremely
accurate, GW cycles from NR will be available, it will be
worthwhile to check whether this result is still true.
[1] A. Abramovici et al., Science 256, 325 (1992);
http://www.ligo.caltech.edu.
[2] H. Lück et al., Class. Quant. Grav. 14, 1471 (1997);
http://www.geo600.uni-hannover.de.
[3] M. Ando et al., Phys. Rev. Lett. 86, 3950 (2001);
http://tamago.mtk.nao.ac.jp.
[4] B. Caron et al., Class. Quant. Grav. 14, 1461 (1997);
http://www.virgo.infn.it.
[5] http://www.lisa-science.org/resources/
talks-articles/science/lisa_science_case.pdf
[6] F. Pretorius, Phys. Rev. Lett. 95, 121101 (2005).
[7] M. Campanelli, C.O. Lousto, P. Marronetti, and Y. Zlo-
chower, Phys. Rev. Lett. 96, 111101 (2006).
[8] J. Baker, J. Centrella, D. Choi, M. Koppitz, and J. van
Meter, Phys. Rev. Lett. 96, 111102 (2006).
[9] M. Campanelli, C.O. Lousto, and Y. Zlochower, Phys.
Rev. D 74, 041501 (2006); ibid. D 74, 084023 (2006); U.
Sperhake, gr-qc/0606079; J. González, U. Sperhake, B.
Brügmann, M. Hannam, and S. Husa, Phys. Rev. Lett.
98, 091101 (2007); B. Szilagyi, D. Pollney, L. Rezzolla,
http://www.ligo.caltech.edu
http://www.geo600.uni-hannover.de
http://tamago.mtk.nao.ac.jp
http://www.virgo.infn.it
http://www.lisa-science.org/resources/
talks-articles/science/lisa_science_case.pdf
gr-qc/0606079
100 200 300 400 500 600 700
Time(M)
Expanded energy-balance equation
Unexpanded energy-balance equation
NR simulation result
FIG. 17: Orbital frequency evolution. The dotted and dashed curves are calculated from the unexpanded and expanded
energy-balance equations. The continuous curve refers to the really long Goddard NR simulation.
J. Thornburg and J. Winicour, gr-qc/0612150; F. Pre-
torius and D. Khurana, gr-qc/0702084.
[10] A. Buonanno, G. Cook, and F. Pretorius, Phys. Rev. D
75 (2007) 124018.
[11] J. Baker, J. van Meter, S. McWilliams, J. Centrella, and
B. Kelly (2006), gr-qc/0612024.
[12] P. Jaranowski, and G. Schäfer, Phys. Rev. D 57, 7274
(1998); Erratum-ibid D 63 029902; L. Blanchet, and G.
Faye, Phys. Rev. D 63, 062005 (2001); V. C. de An-
drade, L. Blanchet, and G Faye, Class. Quant. Grav. 18,
753 (2001); T. Damour, P. Jaranowski, and G. Schäfer,
Phys. Lett. B 513, 147 (2001); L. Blanchet, G. Faye,
B.R. Iyer, and B. Joguet, Phys. Rev. D 65, 061501(R)
(2002); L. Blanchet, and B.R. Iyer, Class. Quant. Grav.
20, 755 (2003); Erratum-ibid D 71, 129902 (2005); L.
Blanchet, T. Damour, G. Esposito-Farese, and B.R. Iyer,
Phys. Rev. Lett. 93, 091101 (2004).
[13] A. Buonanno, and T. Damour, Phys. Rev. D 59, 084006
(1999).
[14] A. Buonanno, and T. Damour, Phys. Rev. D 62, 064015
(2000).
[15] A. Buonanno, and T. Damour, Proceedings of
IXth Marcel Grossmann Meeting (Rome, July 2000),
gr-qc/0011052.
[16] T. Damour, P. Jaranowski, and G. Schäfer, Phys. Rev.
D 62, 084011 (2000).
[17] T. Damour, Phys. Rev. D 64, 124013 (2001).
[18] A. Buonanno, Y. Chen, and T. Damour, Phys. Rev. D
74, 104005 (2006).
[19] T. Damour, B.R. Iyer, and B.S. Sathyaprakash, Phys.
Rev. D 57, 885 (1998).
[20] T. Damour, B. Iyer, and B. Sathyaprakash, Phys. Rev.
D 63, 044023 (2001); ibid. D 66, 027502 (2002).
[21] T. Baumgarte, P. Brady, J.D.E. Creighton, L. Lehner, F.
Pretorius, and R. DeVoe (2006), gr-qc/0612100.
[22] J. W. York, Jr., Phys. Rev. Lett. 82, 1350 (1999).
[23] E. Gourgoulhon, P. Grandclément, and S. Bonazzola,
Phys. Rev. D 65, 044020 (2002).
[24] P. Grandclément, E. Gourgoulhon, and S. Bonazzola,
Phys. Rev. D 65, 044021 (2002).
[25] H. P. Pfeiffer L. E. Kidder, M. S. Scheel, and S. A. Teukol-
sky, Comp. Phys. Comm. 152, 253 (2003).
[26] G. B. Cook, and H. P. Pfeiffer, Phys. Rev. D 70, 104016
(2004); M. Caudill, G.B. Cook, J.D. Grigsby, and H.
Pfeiffer, Phys. Rev. D 74, 064011 (2006).
[27] S. Brandt and B. Brügmann, Phys. Rev. Lett. 78, 3606
(1997).
[28] C. Cutler et al., Phys. Rev. Lett. 70, 2984 (1993).
[29] T. Damour, B.R. Iyer, and B.S. Sathyaprakash, Phys.
Rev. D 67, 064028 (2003).
[30] A. Buonanno, Y. Chen, and M. Vallisneri, Phys. Rev. D
67, 104025 (2003); Erratum-ibid. D 74, 029904 (2006).
[31] A. Buonanno, Y. Chen, and M. Vallisneri, Phys. Rev. D
67, 024016 (2003); Erratum-ibid. D 74, 029903 (2006).
[32] T. Damour, B. Iyer, P. Jaranowski, and B.
Sathyaprakash, Phys. Rev. D 67, 064028 (2003).
[33] K. G. Arun, B.R. Iyer, B.S. Sathyaprakash, and P. Sun-
dararajan, Phys. Rev. D 71, 084008 (2005); Erratum-ibid
D 72, 069903 (2005).
[34] T. Damour, E. Gourgoulhon, and P. Grandclément,
Phys. Rev. D 66, 024007 (2002); P. Grandclément, E.
Gourgoulhon, and S. Bonazzola, Phys. Rev. D 65, 044021
(2002).
[35] L. Blanchet, Phys. Rev. D 65, 124009 (2002).
[36] L.S. Finn, Phys. Rev. D 46, 5236 (1992); L. S. Finn and
D.F. Chernoff, Phys. Rev. D 47, 2198 (1993); É.E. Flana-
gan and S.A. Hughes, Phys. Rev. D 57, 4535 (1998).
[37] http://www.ligo.caltech.edu/advLIGO/scripts/ref_
des.shtml
[38] J. Baker, S. McWilliams, J. van Meter, J. Centrella, D.
Choi, B. Kelly, and M. Koppitz (2006), gr-qc/0612117.
[39] J. Baker, M. Campanelli, F. Pretorius, and Y. Zlochower
(2007), gr-qc/0701016.
[40] B. Abbott et al. (LIGO Scientific Collaboration), Phys.
Rev.D 72, 082001 (2005).
[41] B. Abbott et al. (LIGO Scientific Collaboration), Phys.
Rev. D 73, 062001 (2006).
[42] E. Berti, V. Cardoso, and C. Will, Phys. Rev. D 73,
064030 (2006).
[43] L. Blanchet, Living Rev. Rel. 9 (2006) 4.
[44] T. Damour, B. Iyer, and B. Sathyaprakash, Phys. Rev.
D 62, 084036 (2000).
gr-qc/0612150
gr-qc/0702084
gr-qc/0612024
gr-qc/0011052
gr-qc/0612100
http://www.ligo.caltech.edu/advLIGO/scripts/ref_
des.shtml
gr-qc/0612117
gr-qc/0701016
[45] Y. Pan, A. Buonanno, Y. Chen, and M. Vallisneri, Phys.
Rev. D 69, 104017 (2004).
[46] P. Ajith et al. (2007) (in preparation).
[47] E. Berti, V. Cardoso, J. González, U. Sperhake, M. Han-
nam, S. Husa, and B. Brügmann (2007), gr-qc/0703053.
[48] A. Nagar, T. Damour, and A. Tartaglia, gr-qc/0612096
T. Damour, and A. Nagar, Proceedings of XIth Marcel
Grossmann Meeting (Berlin, July 2006), gr-qc/0612151.
[49] H. P. Pfeiffer, D.A. Brown, L.E. Kidder, L. Lindblom, G.
Lovelace, and M. A. Scheel (2007), gr-qc/0702106.
[50] http://www.ligo.caltech.edu/~rana/NoiseData/S6/
DCnoise.txt.
gr-qc/0703053
gr-qc/0612096
gr-qc/0612151
gr-qc/0702106
http://www.ligo.caltech.edu/~rana/NoiseData/S6/
DCnoise.txt
