Adding eccentricity to quasicircular binary-black-hole waveform models

The detection of gravitational-wave signals from coalescing eccentric binary black holes would yield unprecedented information about the formation and evolution of compact binaries in specific scenarios, such as dynamical formation in dense stellar clusters and three-body interactions. The gravitational-wave searches by the ground-based interferometers, LIGO and Virgo, rely on analytical waveform models for binaries on quasicircular orbits. Eccentric merger waveform models are less developed, and only few numerical simulations of eccentric mergers are publicly available, but several eccentric inspiral models have been developed from the post-Newtonian expansion. Here we present a novel method to convert the dominant quadrupolar mode of any circular analytical binary-black-hole model into an eccentric model. First, using numerical simulations, we examine the additional amplitude and frequency modulations of eccentric signals that are not present in their circular counterparts. Subsequently, we identify suitable analytical descriptions of those modulations and interpolate key parameters from twelve numerical simulations designated as our training dataset. This allows us to reconstruct the modulated amplitude and phase of any waveform up to mass ratio 3 and eccentricity 0.2. We find that the minimum overlap of the new model with numerical simulations is around 0.98 over all of our test dataset that are scaled to a 50M$_\odot$ black-hole binary starting at 35 Hz with aLIGO A+ design sensitivity. A Python package \pyrex easily carries out the computation of this method.


I. INTRODUCTION
Coalescing stellar-mass black-hole binaries are one of the primary sources of gravitational-wave (GW) signals detected by the ground-based interferometers, the advanced Laser Interferometer Gravitational-wave Observatory (aLIGO) [1], Virgo [2], and KAGRA [3]. In the first three observing runs (O1-O3), detection pipelines assumed binary-black-hole (BBH) mergers to have negligible eccentricity when entering the orbital frequencies to which aLIGO, Virgo, and KAGRA are sensitive [4][5][6]. BBHs formed in an isolated environment through a massive stellar evolution are expected to circularize and therefore have undetectable eccentricity by the time they enter the LIGO band [7]. However, BBHs with a detectable eccentricity can form in a dense stellar cluster through dynamical capture [8,9].
A possible scenario is that the binary gains eccentricity due to gravitational torques exchanged with a circumbinary disk [10]. Eccentric BBHs can also form from three-body interactions [9], where the BBH behaves as the inner binary. In this system, the Kozai-Lidov [11,12] mechanism triggers the oscillation that boosts the eccentricity.
Interactions of BBHs in a typical globular cluster suggest a significant eccentric BBH merger rate. As many as ∼5% of binaries may enter the LIGO detector band ( f ≥ 10 Hz) with eccentricities e > 0.1 [13][14][15]. A confident measurement of significant eccentricity in a BBH system would be strong evidence for the dynamical formation scenarios in dense stellar clusters and would boost our understanding of the dynamical evolution of compact objects.
The impact of eccentricity is more substantial during the early inspiral and therefore plays a vital role in the spacebased detector era [16]. In the LIGO band, the detection of GWs from an eccentric orbit would suggest that the binary was formed with a small initial separation and did not have time to circularize, or the binary evolved through an unknown dynamical process. Incorporating eccentric BBH simulations may also lead to an increase in the LIGO/Virgo/KAGRA detection rate [14]. Besides, the detection of eccentric BBH mergers could capture effects from the extreme-gravity regime and therefore can be used for testing the general theory of relativity [17,18].
We highlight the significance of detecting GWs from eccentric BBHs. Constructing template models for eccentric waveforms is challenging, and we aim to make progress towards this goal especially for the late inspiral and merger regimes that are most accessible with today's observations. One of the main difficulties in developing an eccentric waveform model is that only a few numerical relativity (NR) simulations with higher eccentricity are available. Thus, many studies focus on developing eccentric models from the post-Newtonian (PN) expansion. The development of full inspiral-merger-ringdown (IMR) eccentric waveform models is currently an actively researched area [19][20][21].
Huerta et al. [19] construct a time-domain eccentric nonspinning waveform model (e 0 < 0.2) up to mass ratio 5.5, where e 0 is the eccentricity 10 cycles before the merger. Their model is called ENIGMA, a hybrid waveform that has been calibrated using a set of numerical simulations and trained using Gaussian process regression (GPR). Reference [20] presents a low-eccentricity model (e 0 < 0.2) called SEOBNRE using the expansion of the effective one-body (EOB) waveform family. A more up-to-date EOB formalism is demonstrated in Refs. [22,23]. Hinder et al. [21] present a timedomain, nonspinning eccentric waveform model up to mass ratio q = m 1 /m 2 = 3 from 23 NR simulations that are publicly available in the SXS catalog. The referenced eccentricity is e ref ≤ 0.08 starting at seven cycles before the merger. Like Ref. [19], the early inspiral of this model is hybridized with a PN expansion to produce a full IMR model in a Mathematica package [24]. In addition, Ref. [25] recently developed an eccentric model NRSur2dq1Ecc for nonspinning waveforms and eccentricities up to 0.2 from 47 NR simulations. Although the model was trained for q = 1, it can be extended to mass ratio q ≈ 3. Apart from the studies above, nonspinning, loweccentricity frequency-domain models from the PN expansion are publicly available in the LIGO algorithm library (LAL) [26][27][28].
The excitement to search for an eccentric BBH motivated the following analysis. References [29][30][31] recently developed an analysis to find the signature of an eccentric BBH in the O1, O2 and several events in the O3 data using the SEOBNRE model. Additionally, Ref. [32] analyzed the heaviest BBH system during O1-O3, GW190521 [33] with 325 NR simulations. They found that this event is consistent with highly precessing, eccentric model with e ≈ 0.7.
We present a promising method to add eccentricity to quasicircular systems independent of the PN expansion. We apply this method to nonspinning, time-domain waveforms, although in principle it can be used in more general settings. Our technique focuses on a fast reconstruction of the nearmerger eccentric BBH waveform and can be applied to any analytical circular nonspinning model. We build our model from 12 NR simulations and test against further 8 NR simulations from the open SXS catalog [34]. Our method is very simple and can be applied to any circular time-domain model obtained from, e.g., the phenomenological [35][36][37] or EOB [38,39] families.
We model the deviation from circularity visible in the amplitude and phase of eccentric GW signals. This deviation is modeled across the parameter space and can be simply added to any quasicircular model, which elevates that model to include eccentric effects. This approach is inspired by the "twisting" technique that is applied for reconstructing precessing spins from an aligned-spin model to build, e.g., the IMRPhenomP family [36,[40][41][42][43]. The dynamic calibration of the waveform model is motivated by our previous study [44] and the regression techniques tested in detail in Ref. [45]. We calibrate our model for mass ratios q ≤ 3 and eccentricity e ≤ 0.2, and provide it as a Python package called pyrex [46]. Our model has been constructed for a fiducial 50 M BBH and can then be rescaled for other total masses M. We find that the overlap of all our test data against NR is above 98%. Moreover, we expand the construction to earlier regimes than the calibrated time span. Although we do not calibrate for higher mass ratios, the early inspiral, or higher orbital eccentricity, we allow the building of waveforms beyond the parameter boundaries used for calibration.
The organization of this manuscript is as follows: In Sec. II, we present the methodology to construct this model. Section III discusses the primary outcome and the faithfulness of our model. Finally, Sec. IV summarizes and concludes the prospect of our studies. Throughout this article, we use geometric units in which G = c = 1.

II. METHOD
Using NR simulations, we investigate the frequency and amplitude modulations in eccentric BBH signals and implement them in analytical waveforms to develop our model. As described by Peters [7], the orbital eccentricity in binary systems decreases over time due to energy loss through GW radiation. Pfeiffer et al. [47] investigated this in numerical simulations of the SXS catalog. The authors point out that one of the main differences in the evolution of low-eccentricity initial data compared to quasicircular binaries is an overall time and phase shift, where the quasicircular data represent the binary at a point close to merger. Following these studies, Hinder et al. [21] showed that the GW emissions from low-eccentric binaries and circular binaries are indistinguishable near the merger stage. Specifically, Hinder et al. suggest that one only loses 4% of the signal when substituting the GW emission from low-eccentricity binaries with circular orbits 30M before the peak of the amplitude (t = 0). They use this fact to build an eccentric IMR model by replacing the late inspiral eccentric model with a circular waveform. Combining the finding above, we model the decaying eccentricity as amplitude and phase modulation up to t = −29M. We then substitute the GW strain at t > −29M with the circular model for the same binary masses.

A. Data preparation
We use 20 nonspinning NR simulations from the SXS catalog up to mass ratio 3 and eccentricity 0.2 to build our model (see Table I). We follow the definition of eccentricity e comm in Ref. [21] as the eccentricity measured at the referenced frequency, x = (Mω) 2/3 = 0.075. These simulations are divided into a training data set of 12 simulations and the test datasets of 8 simulations, as shown in Fig. 1. Binaries of the test dataset fall within the training data's parameter boundaries. Hence, we do not perform extrapolation with the test data.
We combine the "+" and "×" polarization using the spinweighted spherical harmonics with the following expression [48]: where M and r are the total mass of the system and the distance from the observer, respectively; −2 Y m are the spinweighted spherical harmonics that depend on the inclination angle ι and the phase angle φ; and h m (t) can be extracted from the NR data in the corresponding catalog. We construct TABLE I. NR simulations from the SXS catalog used in this study with mass ratio q = m 1 /m 2 , eccentricity at the reference frequency e comm , and the number of orbits before the maximum amplitude of h 22 . e comm is the eccentricity at the reference frequency (Mω) 2/3 = 0.075 as described in Ref. [21]. The quasicircular waveforms (e comm = 0.000) have eccentricities lower than 10 our model for h 2±2 , the leading contribution of spherical harmonic modes with = 2, m = ±2. Reference [21] suggests that other, subdominant modes are less significant for nearly equal-mass systems with low eccentricity. Here we consider only moderately small eccentricities; therefore we only model the dominant mode. For future studies, subdominant harmonics will be important to model high-eccentricity signals accurately.
We prepare the data as follows. First, we align all the waveforms in the time domain such that the peak amplitude is at t = 0. Subsequently, we remove the first 250M from the start of the waveforms due to the junk radiation, and the last 29M before t = 0 due to circularization (see Fig. 2). Later, we use a circular waveform for t > −29M. We then decompose h 2±2 into amplitude (A), phase (Ψ), and the phase derivative, ω = dΨ dt , where the referenced frequency follows Ref [21]. We model amplitude A 22 and frequency (ω 22 ) as a simple quasicircular piece plus an oscillatory function. The final model then yields the phase (Ψ 22 ) by integrating the frequency.

B. Eccentricity estimator
In numerical simulations, eccentricity is often discussed as a consequence of imperfections in the initial data [49]. It manifests itself as small oscillations on top of the gradual binary evolution, where the oscillation's amplitude is proportional to the eccentricity (see A 22 and ω 22 plots in Figs. 2 and 3). We use this residual oscillation as a key to estimating the eccentricity evolution.
Mroué et al. [50] compare various methods to estimate eccentricity using e X (t). The orbital eccentricity is proportional  Table I. Although higher-eccentricity waveforms produce more oscillations than the lower-eccentricity waveforms, all data appear identical at t > −30M due to circularization as shown in the top panels. We employ the residual amplitude A 22 and frequency ω 22 to develop our model in the late inspiral regime.
to the amplitude of a sinusoidal function, e X (t), expressed by where X is either ω 22 or A 22 , and X c (t) is the X quantity in circular binaries instead of low-order polynomial fitting functions that are often used in the literature. We reverse this relation to convert a circular model [with given X c (t)] to an eccentric model using an analytical description of the oscillatory function e X (X c ). We apply the Savitzky-Golay filter [51] to smooth the e X (t) curves from noises caused by numerical artifacts. Savitzky-Golay is a digital filter applied to smooth the selected data points without altering the signal direction by fitting the adjacent data with a low-degree polynomial fit. We stress that the definition of the orbital eccentricity is not unique. Thus, one could use different definitions of eccentricity. In principle, any definition can be accepted if consistently applied to the study in question. The NR data we use are labeled with a value for the initial eccentricity that is based on PN initial data [21]. As we shall discuss below, these labels are similar to what we estimate for the eccentricity using Eq. (2), but not identical. However, we refrain from redefinition of the initial eccentricity of the NR data and instead identify each NR simulation with the value of eccentricity at the reference frequency (Mω) 2/3 = 0.075 determined by the original Ref. [21]. We do this because (i) we want to avoid any confusion as to what NR data we are using and what their properties are, and (ii) by making the amplitude of e X a function of the eccentricity label imposed by Ref. [21], we introduce an extra uncertainty that may be seen as representing the ambiguity in determining the initial eccentricity of the respected NR simulations. Thus, we present a conservative estimate of the approach's accuracy.
As a check, we compute the orbital eccentricity using the eccentricity estimator (e X ) and find that the results agree with a maximum relative error of roughly 10% against e comm quoted in the SXS catalog and given in Table. I. In Fig. 4, we present the eccentricity estimator e X (X c ) as a function of its circular amplitude and frequency, A c and ω c , respectively.

C. Fitting e X
Our main goal is to model an eccentric waveform by modulating the amplitude and phase of a circular model. To construct the model, we interpolate the additional oscillation of an eccentric waveform depending on its eccentricity and mass ratio, where the relationship between the circular and the eccentric model is expressed in Eq. (2). Accordingly, we look for a fitting function to model e X (X c ) that relies on the desired parameters (q, e) and reverse Eq. (2) to obtain the eccentric amplitude and frequency. We then integrate the frequency to obtain the eccentric phase and construct the eccentric h 22 .
We note that alternatives to fitting the amplitude and frequency modulations have been studied in Ref. [25]. In particular, they investigated using the phase residual instead of the frequency, or fitting the eccentric amplitude and phase (or frequency) directly instead of recasting the problem in terms of differences to noneccentric signals. Here we find that the most suitable strategy for our approach is to fit the residual amplitude and frequency oscillation defined as the eccentricity estimator (e X ) that comes from {A 22 , ω 22 } and integrate ω 22 to obtain the phase (Ψ 22 ).
In a suitable parametrization, the eccentricity estimator e X is a decaying sinusoidal function (see Fig. 4) with its ampli-tude defined by the orbital eccentricity e [50]. To model e X for various eccentricities and mass ratios, we fit e X with a set of free parameters modifying a damped sinusoidal function. These parameters are two amplitude quantities (A and B), a frequency ( f ), and phase (ϕ) with the following relation: A, B, f , and ϕ are standard damped sinusoidal parameters obtained from the optimized curve fitting. We use a X κ c instead of X c to describe the evolution of the residual oscillations of the amplitude and frequency mainly for the following reasons: X c is a rapidly evolving function. Therefore, it is more difficult to model e X with a standard sinusoidal function with a constant frequency. Although it is in principle possible to use X c directly in the model, we would have to slice the data into multiple small time windows that overlap. Thus, the results will be less smooth; one would have to blend all those individual functions defined on small intervals into one big function. Besides, we cannot guarantee our result beyond our calibration range, especially for the early inspiral. Using a power law allows us to fit the entire region with one set of free parameters. However, we note that the power law of X c induces a twist resulting in infinitely large eccentricities for the very early inspiral stage. That is a problem with assuming exponential decay, and the fact that the power law we use has a negative exponent. We fit our model e X from the starting frequency f low = 25 Hz for a circular BBH with a total mass M = 50 M . The power law for ω c is κ = −59/24, and for A c it is κ = −83/24. We emphasize that these values are customized i.e., we expect that one might need different values to calibrate with higher eccentricity, a higher mass ratio, or a different starting frequency.
By optimizing the curve fit between e X and Eq. 3, we obtain the four quantities (A, B, f, ϕ) for all training data. The relation between the mass ratio (q), eccentricity (e), and the three parameters A, B, f is shown in Fig. 5. The amplitude components A and B are strongly correlated to eccentricity, whereas the mass ratio determines the frequency squared. Hence, we perform one-dimensional linear interpolation across eccentricity to obtain the values of A and B. Similar to that, we linearly interpolate f 2 across mass ratios. We choose f 2 instead of f because the data is smoother for interpolation. The square root of f 2 gives either positive or negative values. However, this ambiguity can be absorbed by the phase parameter ϕ.
The phase parameter ϕ is an additional degree of freedom that we cannot explore sufficiently with the available NR data. For small sets of NR simulations with nearly constant values of q and e, but varying , we find that the best-fit ϕ mirrors changes in . Thus, we expect that it may correlate strongly with the mean anomaly. Because the orientation of the ellipse is astrophysically less interesting than the value of the eccentricity, we do not attempt to model the effect of varying the mean anomaly other than introducing the phenomenological nuisance parameter ϕ. We interpolate the other parameters when generating a new waveform model with different mass ratios and referenced eccentricities.
We apply a one-dimensional interpolation for each key quantity shown in Fig. 5. A and B are interpolated over different eccentricities e, f 2 is interpolated over the mass ratio q, and the phase of the oscillation ϕ can be chosen arbitrarily.
Once we obtain the eccentricity estimators e X using the interpolated quantities, we substitute the results to reconstruct A 22 and ω 22 using Eq. 2. To construct Ψ 22 , we integrate ω 22 numerically using the trapezoidal rule. We truncate the waveform at t = −50M and join it with the nonspinning circular model. We then smooth the transition with the Savitzky-Golay filter at −46M < t < −25M.
We then build h 2±2 as the combination of the amplitude and phase as follows: To reconstruct the gravitational-wave strain h = h + − h × , we compute the spin-weighted spherical harmonics Y m (ι, φ) and employ Eq. 1.

III. RESULTS
We built a new nonspinning eccentric model by modulating the residual amplitude and phase oscillations of the circular analytical models, IMRPhenomD [35] and SEOBNRv4 [52]. IMRPhenomD is an aligned-spin IMR model that was originally built in frequency domain and calibrated to numerical simulations for mass ratios q ≤ 18. SEOBNRv4 is an alignedspin time-domain IMR model [52,53] that has been calibrated to 140 NR waveforms produced with the SpEC code up to mass ratio 8 and extreme-mass-ratio signals.
As described in Sec. II, we interpolate the residual amplitude and phase oscillations of the training dataset for the given mass ratio and eccentricity. To construct a new, eccentric waveform for the intermediate to near-merger regime, we then use one of the nonspinning circular models with the desired mass ratio, compute the eccentricity estimators (e X ) from the analytical description given in Eq. (3), and reconstruct the desired eccentric waveform model for each test data. We develop a map from circular nonspinning waveforms to eccentric waveforms that can be applied to any analytical model with a relatively simple and fast function using only 20 NR simulations.
We evaluate the results by computing the overlap between the new model and the NR test data. The overlap is maximized over a time and phase shift, as well as the free phase offsets of the residual oscillations. Mathematically, we define the overlap O based on an inner product between two waveforms: where S n is the sensitivity curve of the corresponding GW interferometer,h( f ) is the unfaithfulness is defined by We investigate three sensitivity curves for the future GW detectors, aLIGO A+, the Einstein Telescope (ET), and Cosmic Explorer (CE). LIGO A+ is the future GW interferometer with 1.7 times better sensitivity than the current detector, expected to start observing in mid-2024 at the earliest [54]. The ET is a 10 km GW observatory planned to be built on the border between Belgium, Germany, and Netherlands which could be operating in the mid-2030s [55]. The ET is expected to have higher sensitivity towards the low-frequency range. CE is a 40 km third-generation GW detector which has higher sensitivity towards low redshift (z > 10) that is planned to start observing in the 2030s [56]. Since our model focuses on the late inspiral case, and because the unfaithfulness is insensitive to a change in overall signal-to-noise ratio, the values obtained for the future third-generation detectors show similar behavior [57]. Hence, we only show the overlap results for the LIGO A+ design sensitivity. A possible caveat is that our model might not fill the LIGO A+ band down to 10 Hz. Thus, there is a chunk of inspiral power missing in the signal. Figure 6 visually compares the strain h 2±2 of each NR test dataset with the new eccentric nonspinning signal built from analytical models, IMRPhenomD and SEOBNRv4 for a 50 M BBH with inclination angle ι=0 (face-on) and phase of coalescence, φ c =0. Using our method, we find that the minimum overlap between the new model and NR is ≈ 0.98 (log 10 M = −1.8) over all of our test datasets. The minimum overlap occurs at the highest eccentricity in the test dataset.
Although we calibrated the new model for limited ranges in mass ratio, eccentricity, and time, we let the production of the new model go beyond our calibration range. In Fig. 7, we show the unfaithfulness of the new model against the NR test data for various total masses with the aLIGO A+ design sensitivity curve. The left panel shows the unfaithfulness within the calibrated frequency range, between 25 Hz and the ISCO frequency scaled over the total mass. Similarly, the right panel presents the unfaithfulness beyond the calibrated frequency range, between 20 Hz and the ringdown frequency. We use the definitions of the ISCO and ringdown frequencies as follows: and Figure 7 shows that the mismatches decrease toward higher-total-mass systems. As the total mass increases, the overlap computation covers a smaller waveform regime towards merger in the frequency space. Since the eccentricity decreases over time, the near-merger regime has lower eccentricities. Thus, the overlap between the model and the corresponding NR simulation is better for the higher-mass systems compared to the lower-mass ones. For comparison, we find that mismatches between circular analytical models and the eccentric NR test data are at least 1 order of magnitude worse than the results we find for our eccentric model.
The unfaithfulness between eccentric waveforms is better for {25, ISCO} than for {20, Ringdown}. We investigate the contribution weight between the early inspiral and the ringdown in the unfaithfulness results by comparing with the {25, Ringdown} and {20, ISCO} ranges. We argue that the mismatches for the low masses are dominated by the inspiral, whereas for high masses, the mismatches are dominated by the merger or ringdown. In the mismatch computation, we FIG. 6. IMRPhenomD (orange) and SEOBNRv4 (green) circular waveforms twisted into eccentric models. log 10 M 1 is the log mismatch of IMRPhenomD against the NRs waveform (shown in blue), and log 10 M 2 gives the log mismatch of SEOBNRv4 against NR with the same mass ratio and eccentricity, respectively. The total mass of the system is M = 50M , and the mass ratio (q) and eccentricity (e) are shown in the title of each plot. We employ the A+ design sensitivity curve starting at f = 35 Hz (see text) to compute the match. The black vertical lines mark the range in which we perform the interpolation and compute the match. add padding in the ringdown area, but the early inspiral should come purely from the fitting data.
Furthermore, we test how well one can extract the parameters of an eccentric signal h(q, e) by comparing with various waveforms with different eccentricities e and mass ratios q. We generate a pyrex waveform (q = 1, e = 0.144) and compare it with various other signal parameters (q, e) using the same analytical waveform model. The results are shown in Fig. 8. We emphasize that in this study, we did not run a standard parameter estimation (PE) pipeline that stochastically explores a much greater parameter space. In particular, we do not consider varying the total mass or spin. Hence, our results are only a first indication of potential parameter ambiguities. Our results in Fig. 8 show that the mismatch between the generated waveform and other waveforms having similar mass ratios but different eccentricities is relatively low, suggesting that an accurate measurement of the eccentricity is challenging for high-mass BBH systems where only the late inspiral and merger are accessible through the GW detection.

IV. CONCLUSION AND FUTURE PERSPECTIVES
The detection of GWs from an eccentric BBH merger would be a crucial step towards understanding the physical evolution of compact binary coalescences and the nature of BBHs in globular clusters. Due to limitations in waveform modeling, the current search and parameter estimation pipelines in the LIGO/Virgo data analysis rely on analytical waveform models for circular binaries. One of the limitations to developing eccentric BBH models is the small number of eccentric NR simulations. NR simulations that are publicly available have low eccentricities (e ≤ 0.2) at Mω 2/3 = 0.075. We use 20 NR simulations from the open SXS catalog and split them into 12 training datasets and 8 test datasets to develop our method.
We presented a novel method to convert any circular nonspinning waveform model into a low-eccentricity nonspinning waveform. To develop our method, we analyzed the residual modulations in the amplitude and frequency of eccentric waveforms compared to the circular signals with the same mass ratio in the 12 NR simulations of the training dataset. We modeled the decrease of eccentricity over time, known as the eccentricity estimators, e X , using a damped sinusoidal fit, where the fitting function is built upon four key parameters. We then performed a one-dimensional interpolation for each key parameter (A, B, and f ) to build the eccentric waveform with the desired mass ratio and eccentricity. One of our model parameters, ϕ, shows no clear correlation with the physical parameters we explore. However, the small number of NR simulations used here did not allow us to model the effect of varying the mean anomaly in detail, and we expect ϕ to represent this degree of freedom. When quantifying the agreement between our model and the test data, we maximize over this nuisance parameter.
We then build a new model using the fitting values of e X and the amplitude and frequency of the circular model which here we take from IMRPhenomD and SEOBNRv4. Our new model has an overlap 0.98 O 0.999 over all NR simulations in our test dataset with the LIGO A+ design sensitivity curve. We generate an eccentric waveform model derived from a nonspinning circular model, IMRPhenomD or SEOBNRv4, and compare the signal with models for different mass ratios and eccentricities. Waveforms with higher parameter distance have lower overlap. The color bar shows the log 10 mismatch.
We hint that we need more training and test datasets for further development of this model beyond the current parameter boundaries. The computation of our method can be performed easily and quickly in the Python package pyrex [46].
Although we calibrate our model to a 50 M BBH (q ≤ 3, e ≤ 0.2) starting at frequency f low = 25 Hz, we let the computation go slightly beyond the calibrated range. The calibrated time range of the waveform is from the late inspiral up to the near-merger phase, but we can extend the model through merger and ringdown by using the circular data. For the early inspiral, an analytical PN model could be used to complete the description of the entire coalescence. This way, our approach can be adapted to develop a complete IMR eccentric model. This would be especially important for future generations of GW interferometers as they have higher sensitivity especially in the low-frequency range. Careful studies of eccentric search and parameter estimation are needed to detect eccentric compact binary coalescences and their origin.