A hybrid geometric-random template placement algorithm for gravitational wave searches from compact binary coalescences

A hybrid geometric-random template placement algorithm for gravitational wave searches from compact binary coalescences

Soumen Roy soumen.roy@iitgn.ac.in    Anand S. Sengupta asengupta@iitgn.ac.in    Nilay Thakor thakor.nilaysinh@iitgn.ac.in Indian Institute of Technology Gandhinagar
Gujarat 382355, India.

Astrophysical compact binary systems consisting of neutron stars and blackholes are an important class of gravitational wave (GW) sources for advanced LIGO detectors. Accurate theoretical waveform models from the inspiral, merger and ringdown phases of such systems, are used to filter detector data under the template based matched filtering paradigm. An efficient grid over the parameter space at a fixed minimal match has a direct impact on the overall time taken by these searches. We present a new hybrid geometric-random template placement algorithm for signals described by two masses and one spin magnitude parameters. Such template banks could potentially be used in GW searches from binary neutron stars and neutron star-blackhole systems. The template placement is robust and is able to automatically accommodate curvature and boundary effects with no fine tuning. We also compare these banks against vanilla-stochastic template banks and show that while both are equally efficient in the fitting-factor sense, the bank sizes are larger in the stochastic method. Further, we show that the generation of the proposed hybrid banks can be sped-up by nearly an order of magnitude over the stochastic bank. Generic issues related to optimal implementation are discussed in detail. These improvements are expected to directly reduce the computational cost of gravitational wave searches.


I Introduction

The direct detection of gravitational waves (GW) was vigorously pursued by the LIGO Scientific Collaboration (LSC) consisting of several hundred scientists which culminated with the discovery of the first gravitational wave event GW150914 Abbott et al. (2016a) in the twin advanced-LIGO (advLIGO) Harry and the LIGO Scientific Collaboration (2010) detectors. This event was determined to have been caused by the inspiral and merger of a spinning binary blackhole system of component masses located nearly 1.3 billion light years away from Earth. It was also the first observational evidence for the existence of stellar mass blackholes. Subsequently, advLIGO detectors also detected a second event GW151226 Abbott et al. (2016b) from the inspiral and merger of lighter compact objects. These detections mark the transition to the era of gravitational wave astronomy.

Several other kilometer scale detectors are under upgradation or construction at present around the globe - these include the French-Italian advanced Virgo detector Acernese et al. (2015) and the Japanese cryogenic KAGRA detector Somiya (2012); Aso et al. (2013). In India, the proposal for the advanced LIGO-India detector Iyer et al. (2011); Unnikrishnan (2013) has been approved and is expected to be built over the next few years. The network of these advanced detectors is expected to improve their overall science potential and herald a new wave in astronomy with the potential to observe the very early Universe and compliment information gathered by electromagnetic observations.

Binary blackholes and neutron stars are considered one of the most promising sources for the advanced terrestrial detectors. Precise theoretical model waveforms for GW emitted from the inspiral, merger and ringdown phases of such compact binary coalescences (CBC) are now available. These models are parametrized by the system’s intrinsic properties such as the component masses, spin etc. and extrinsic properties such as sky location, distance to the source, time at coalescence etc. Accurate theoretical models allow the use of matched filtering technique to search for weak GW signals buried in detector noise. As the signal parameters are not known a priori, one filters the data using a set of expected signals spanning the deemed parameter space. Each one of these expected signal corresponds to a single point in the parameter space, and are collectively known as the template bank. Coverage of the full range of search parameters using a finite grid of discrete points leads to an inevitable loss in the signal-to-noise ratio (SNR) which can be controlled by fixing the minimal match of the bank. The latter is often decided by striking a balance between desired detection efficiency and computational cost of carrying out the search.

Advanced detectors are reaching unprecedented sensitivities at low frequencies. Over the last few years, the development of theoretical spinning waveform models have also reached a mature stage. The combined effect of these factors is that one now needs to search over a significantly larger volume compared to initial LIGO era, in a parameter space spanning three or more dimension - posing new challenges for data analysis, which include devising efficient grid placement strategies. The practical issue of optimal template placement in matched-filter searches is an important open optimization problem Manca and Vallisneri (2010); Messenger et al. (2009).

There are at least two approaches to solving the problem of efficient template placement - (a) via deterministically placed points through the tessalations of a geometrical lattice Cokelaer (2007); Babak et al. (2006) or (b) via stochastically placed points by choosing them at random Harry et al. (2009) over the parameter space. The geometric method requires the metric over the signal manifold and has been used extensively for CBC searches over 2D mass parameters by arranging the grid points in a hexagonal lattice. A variant of the geometrical placement method for aligned spin CBC systems has been explored Brown et al. (2012) for the TaylorF2 Buonanno et al. (2009); Poisson and Will (1995) signal model. Lack of availability of metric and also the intricate fine tuning required to avoid uncovered regions arising from variations in curvature across the parameter space also make it difficult to generalize the geometric placement methods to higher dimension.

The stochastic template bank is constructed from random proposals drawn from a uniform distribution over the deemed parameter space that are accepted as a new template point only if the new proposal is far (in minimal match sense) from existing templates in the bank. Such banks are easy to implement, robust and can work even without the explicit knowledge of the metric on the signal manifold by using brute force match calculations. This approach has been demonstrated to be more space efficient than a square lattice in 2D but less efficient than a hexagonal lattice. Such banks can be calculated for higher dimensions as well Ajith et al. (2014).

An earlier attempt to combine the geometric and stochastic methods Capano et al. (2016) by seeding a stochastic bank with a pre-fabricated geometric bank has also been tried and has been demonstrated to improve the efficiency marginally.

We present a new hybrid algorithm for template placement in 3D (two mass components and one reduced spin magnitude) parameter space for gravitational wave searches from compact binary coalescence by combining the efficiency of optimal geometrical placement and the robustness and ease of stochastic placement algorithm. This geometric-random bank placement method uses a local truncated octahedral lattice to place the templates and requires the metric over the signal manifold.

Plan of the Paper:

The paper is organized as follows:

In Section II we recapitulate the fundamentals of templated matched filtering technique used for gravitational wave searches from CBC and review the definition of the metric on the signal manifold with a view to set the notation used in the paper and also define key terms. Numerous papers on template placement strategies have been written over the last few years and many new ideas have emerged Brown et al. (2012); Capano et al. (2016); Harry et al. (2009); Ajith et al. (2014): we provide a concise review of these efforts for the benefit of the reader. In particular, we elucidate the stochastic algorithm by casting it in two different ways which are algorithmically equivalent. This sets the stage for Section III, where we present a new hybrid geometric-random template placement algorithm and give a detailed explanation of the issues concerning its optimal implementation. We argue that by construction this new method can not be less space efficient than the vanilla-stochastic placement method. We also construct an explicit template bank for neutron star-blackhole searches using this new algorithm and compare it against a vanilla-stochastic template bank, noting the improvement in overall bank size as well as the time taken to generate the bank. The signal model Ajith (2011) was used for this purpose. The template banks presented in this work have been calculated using the Shoemaker (2010) noise sensitivity curve for advLIGO.

In Section IV, the template banks generated by the hybrid method are tested and validated against the vanilla stochastic bank. We present the fitting factor results using signal injections using and TaylorF2 aligned spin waveform models using software implemented in the LALApps package of the LIGO Algorithm Library LIGO Data Analysis Software Working Group () and show that the two banks are nearly identical in efficiency. We also compare the hybrid bank against a lattice seeded vanilla-stochastic bank and show that the former is more efficient.

Finally, in Section V, we summarize the main results and make some comments related to several key issues related to the new method, indicating a possible way for extending it to higher dimensions.

Ii Metric on the signal manifold

In this section we shall quickly summarize the basics of template based matched filtering technique used in GW searches from CBC. A basic assumption in matched filtering based searches is that the astrophysical GW signal buried in the detector noise is faithfully represented by the signal model used in the search over the range of search parameters.

The signal manifold is the set of all possible GW signals characterized by the parameter vector . It is customary to represent the corresponding frequency domain signal as . The detector output consists of detector noise and a possible gravitational wave signal of unknown parameters. The additive noise throws the signal out from this manifold to the space of all possible functions. In order to find the point in signal manifold closest to the detector output , the latter is projected over the signal manifold by calculating the maximum likelihood over which serves as the detection statistic. For additive Gaussian noise, the likelihood is given by Finn (1992)


The inner product is the complex cross-correlation Sathyaprakash and Dhurandhar (1991) between the detector output and the gravitational wave signal weighed inversely by the noise power spectral density of the detector Babak et al. (2013):


where is the one-sided noise power spectral density defined by , asterisk denotes the complex conjugation operator, the frequency range marks the effective bandwidth of the detector and is the time delay between these two signals. The signal-to-noise ratio (SNR) after filtering is defined as:


Without any loss of generality, we assume that all template waveforms are normalized such that . From Eq. (1) it is clear that this allows us to use the log-likelihood function (or equivalently, the SNR) maximized over the parameters, as the detection statistic.

The log-likelihood function can be maximized over all time lags () by using Fast Fourier Transform (FFT) based convolution as shown in Eq. (2). It can be maximized over other extrinsic parameters analytically Sathyaprakash and Dhurandhar (1991). On the other hand, as one cannot maximize the log-likelihood function over intrinsic parameters analytically, a brute force approach is needed. In this case, a discrete grid of points is placed to cover the intrinsic parameter space. One evaluates the log-likelihood surface at each of these points and the maximum is suitably determined.

The template bank consisting of these discrete set of points on is constructed using a control parameter ; commonly known as minimal match in GW literature. This parameter is chosen such that the minimum overlap of an arbitrary vector in the signal manifold (within the deemed parameter space) and at least one template in the bank never drops below this value. The art of template placement lies in maximizing the inter-template separation without violating this constraint with the aim of achieving the smallest bank size. In this way, we can map the template placement to the sphere covering problem Conway and Sloane (1999); Prix (2007) with spherical cells of radius equal to . One wants the smallest number of overlapping templates (i.e. metric spheres) to fully cover the space (i.e. leave no holes). The equation of such a spherical cell, centered at , is given by


From the normalization of waveforms, it is clear that . For high values of this parameter, the LHS can be Taylor series expanded upto leading order in the small parameter :


where are the components of the vector and


is the metric over the signal manifold which is essentially the Fischer information matrix projected on the intrinsic parameter space and calculated using standard covariance matrix method Poisson and Will (1995); Owen (1996); Owen and Sathyaprakash (1999). The match between two nearby points in the parameter space can be calculated easily using the metric. Eq. (5) can be re-arranged as , and identified to be the equation of an ellipsoid in 3D centered at a point . For higher dimensions it represents a hyper-ellipsoid. We shall refer to this as the minimal-match ellipsoid elsewhere in the paper. In GW searches from CBC sources, template banks are usually constructed at Babak et al. (2013). This corresponds to a loss in detection rate of assuming uniform distribution of such sources.

ii.1 State of the art in template placement

As mentioned earlier, the geometric and stochastic template placement algorithms are two broad class of methods used in searches for gravitational wave signals from compact binary coalescences.

Previous searches for GW signals from nonspinning compact binaries in initial-LIGO and Virgo Abbott et al. (2009a, b); Abadie et al. (2010) data have used the metric based geometrical hexagonal template placement in two dimensions Cokelaer (2007); Babak et al. (2006); Owen (1996); Owen and Sathyaprakash (1999); Sathyaprakash and Dhurandhar (1991). Templates are placed in chirp time coordinates instead of component masses , since the templates are almost uniformly spaced in the former. This process starts by initializing a template point and then finding neighboring points in a hexagonal lattice, where the centre of the hexagons represent the position of individual templates. Hexagonal tiling offers the most efficient space filling in two dimensions, which optimizes the number of template points and in turn reduces the total computational cost of the search. To construct this geometric bank one require the semi-analytic metric on the signal manifold assumed to be slowly varying over the parameter space. The curvature effects leads to some loss of efficiency in this strategy. This geometric placement was initially demonstrated Cokelaer (2007) for 2PN SPA family of waveforms. At present, the metric for 3.5PN SPA waveforms Keppel et al. (2013) is available which allows template placement for such waveforms as well.

Geometric template placement in higher dimensional intrinsic parameter space (e.g. component masses and spins) has several problems. First, the metric may not be available for such signal manifolds and secondly, (unlike the hexagonal packing in 2D), optimal geometrical placement in higher dimensions with curvature are not well known. It is further complicated from curvature and boundary effects. To mitigate the curvature issues which leads to rapidly changing metric components over the parameter space, new coordinates are being explored in which the metric is slowly varying Ajith et al. (2014), but they do not solve the problem completely. Recent studies Brown et al. (2012); Harry et al. (2014) have explored geometrical placement in higher dimensions for aligned-spin BNS and NSBH systems for some specific waveform families. In this method, one constructs a metric on the parameterized coordinates - taken to be the coefficients of the 3.5PN TaylorF2 expansion of the orbital phase Pai and Arun (2013), instead of the usual chirp time coordinates. Since the metric is globally flat in these new coordinates, one can globally transform it into a Euclidean coordinate system. Finally, a principal coordinate analysis facilitates the projection to an effective lower dimensional parameter space, which can be covered by a grid placed in a hexagonal lattice. For NSBH Harry et al. (2014) systems, the parameter space is covered by stacking several two dimensional hexagonal lattices along the direction of the minor axis. This method is available in the PyCBC software package Nitz et al. (2016); Usman et al. (2016); Dal Canton et al. (2014); Allen et al. (2012).

An alternative approach of template placement is the so-called stochastic method Harry et al. (2009) where one starts with an empty template bank (or a set of seed points) to which random points, drawn from a uniform distribution over the deemed parameter space, are appended in an iterative fashion. At every iterative step, a new random point is proposed to be included in the template bank: this proposal is rejected if it happens to lie too close to the points that are already in the existing list otherwise it is accepted. For each accepted proposal, the rejection rate is determined and the process terminates when this rate exceeds a certain threshold averaged over the last few acceptances. We shall call this the bottom-up approach in building the list of templates by considering random proposals one by one, retaining only the valid ones. This has been encoded in the LSC Algorithm Library (LAL) software suite via the program lalapps_cbc_sbank.

The stochastic bank placement algorithm can also be cast in a top-down fashion. In this alternate implementation, one starts with an empty template list and a list of very large number of proposals distributed uniformly over the deemed parameter space. One picks a random point and appends it to , following which all points from that lie within the minimal match distance from are removed. The process continues in an iterative manner until all points from are exhausted. We call this the top-down approach as the template-bank is microfabricated out of a large block of random proposals by paring it down to the desired shape and order.

Both these methods are algorithmically equivalent - however, the top-down approach is more useful in projecting a geometric structure over the stochastic template bank leading to the hybrid geometric-random placement algorithm described in this paper. In Table 2 we demonstrate that the top-down approach is also faster by a factor over the traditional bottom-up approach as one is able to eliminate many proposals for a single accepted proposal using efficient computational data structures such as binary search trees (BST Hibbard (1962)). The number of templates generated by both these methods is nearly identical.

The stochastic method is relatively easier to implement and does not require the metric over the signal manifold per se. The distance between the proposed point and the elements of the current template bank can be directly calculated by evaluating the match inner product Eq. (2). This brute-force approach for match calculation allows it to be extended easily to higher dimensional parameter spaces and overcome irregular boundary effects. The disadvantages of stochastic method include the requirement of high computational time as several million proposals have to be processed to guarantee adequate coverage and the fact that it generates substantially more templates than the geometric bank. The computational time can be reduced, if we use the metric (if available) for match calculation. But the intrinsic stochastic nature of the algorithm leads to the inefficiency in grid placement.

Another instance of template placement developed for aligned spin binary black hole (BBH) searches has explored a combination of geometric and stochastic approaches Capano et al. (2016); Privitera et al. (2014). In this method, at first one generates an aligned-spin geometric hexagonal lattice template bank upto some valid range of parameters which is then used as a “seed” bank for the stochastic placement - thereby accelerating the placement. This method generates fewer template points than the stochastic method and has been used as part of the uber template bank used in CBC searches in the data from the first observational run of aLIGO Abbott et al. (2016c).

Iii A new geometric random algorithm for template placement

In this section we present the metric based hybrid geometric-random template placement algorithm in three dimensions using a truncated octahedral lattice. Such lattices are the Dirichlet-Voronoi polytope of body-centred cubic lattice Schürmann and Vallentin (2006). The latter provide optimal coverage for conformally flat spaces where the metric-coefficients are constant Prix (2007). It is interesting to note that this is in line with Lord Kelvin’s conjecture Thomson (1887) according to which, truncated octahedron based space filling is optimal in flat 3D space. A truncated octahedron is a 14-sided space filling polyhedron and has the highest volumetric quotient (Eq. 18) which makes it suitable for the lattice structure. Other geometric properties are tabulated in Appendix B.

Such a template bank would be applicable for gravitational wave searches from the compact binaries described adequately by their component masses and a single “reduced spin” parameter using signal model. This signal model is constructed using post-Newtonian (PN) template family of gravitational waveforms from inspiralling compact binaries with non-precessing spins. Here, the spin effects are captured by the parameter defined as a mass-weighted, linear combination of individual dimensionless spin magnitudes. The details of the signal model can be looked up in Appendix A.

In order to construct the three-dimensional geometric random template bank, we require the metric on the parameter space. The metric for approximant varies rapidly in coordinate system where denotes the total mass and denotes the symmetric mass ratio. In order to enhance the efficiency of the algorithm, one places the templates in dimensionless chirp time coordinates , in which the metric components are slowly varying over the parameters. These new coordinates are defined as:


However it is to be noted that the curvature effects do not vanish completely due to the above coordinate transformations. As such, Kelvin’s conjecture does not hold directly in this non-flat space. In the following sections, we show that the truncated-octahedral design can still be used for spaces with slowly varying curvature. This is achieved by merging the stochastic bank placement algorithm along with a local lattice.

The method outlined in Algorithm (1) proceeds by initializing three lists as follows:

  1. : a list of uniformly distributed random points sprayed over the deemed parameters space

  2. : an empty list for template points

  3. : an empty temporary list

At first is initialized or seeded with a random point chosen from which is immediately appended to the list of templates as its first element. The possible 14 TO neighbours of this initial point are then calculated and appended to followed by the removal of the first element. Only those TO neighbors which fall within the parameter space and are farther than from existing points in the and can be appended. This completes one iteration of the inner loop of the algorithm. We continue these steps until becomes empty. The latter happens when no new valid TO neighbors are appended in successive iterations of the inner loop of the algorithm which ultimately exhausts all the points in . At this stage, all points from that are within the minimal match distance from the elements in are deleted; following which is re-seeded with a new random point from . Termination of template placement algorithm occurs when is exhausted.

4:while  do
5:         random point
6:        while  do
7:                append
8:               find all possible TO neighbours of
9:                append valid TO neighbours
10:               delete
11:        end while
12:        delete all minimal match neighbours of from
13:end while
Algorithm 1 Geometric-Random Template Placement

The advantage of casting the stochastic template placement algorithm (see Sec II.1) in a top-down fashion is clear when we consider the extreme case where no TO lattice neighbours can be found for any proposal point in . In this case, it is clear from Algorithm 1 that the geometric-random algorithm naturally falls back to the vanilla-stochastic placement algorithm. This brings in the robustness of the latter to the proposed new algorithm. This contingency arises for a small fraction () of points in - due to boundary edge effects and small uncovered patches arising from curvature effects. By construction the geometric-random placement method presented here is more efficient than vanilla-stochastic algorithm.

Also note that the TO, whose neighbors are added to in Lines 8-9 of Algorithm (1) refers to the one inscribed inside the minimal match ellipsoid given in Eq. (5) above. The gradual change in curvature (changing metric components) leads to changing orientation and size of these ellipsoids which affects the size and orientation of the inscribed truncated-octahedron. This in turn, affects the coordinates of its 14 neighbors. In this way the placement algorithm automatically responds to the curvature effects. This is equivalent to assuming flat local patches of the signal manifold which can be optimally covered by lattice.

iii.1 Implementation

We now highlight some salient features for efficient implementation of the algorithm.

iii.1.1 Initializing with uniform random points

In order to generate a random list of points over the parameter space, at first we calculate the minimum and maximum possible values of dimensionless chirp times and for the range of parameters over which the bank is to be placed. A set of random points are generated in space along the boundaries marking the constraints for component masses, chirp masses and mass ratios. Afterwards, a coordinate transformation is taken to the dimensionless chirp time coordinate system which allows the calculation of the extremities in these coordinates. The extreme values of are evaluated using in Eq. (7). Once these extremities are known, uniform random points are proposed within these values - out of which only those which lie within the specified mass and spin ranges are retained.

The final banksize depends strongly on the initial density of points in upto a critical value, beyond which it does not change very much. Empirically, we find this number to be about of the total volume calculated in the dimensionless coordinates. The latter depends on the metric and varies across the parameter space due to which we estimate an average cell volume by sampling different points across this space. The total number of cells needed to cover the space is given by the ratio of the total volume divided by the average volume of a cell. Monte Carlo integration method can be used to estimate the total volume.

iii.1.2 Deletion of points within minimal-match ellipsoid

We discuss an efficient way to delete points from that lie within minimal match distance from the template bank points. For a fixed minimal match value, Eq. (5) represents the surface of an ellipsoid with semi axes length:


where ’s are eigenvalues of the metric evaluated at the ellipsoid centre. The orientation of the ellipsoid is defined by matrix composed of eigenvectors of the metric , where is component of the eigenvector. To visualize the deletion process, let us consider a cell centred around a template. The points in lying inside the cell can be determined by calculating the match between this template with the all other points in but such a brute force approach will not be computationally efficient. We need to find a smaller set of points against which we can identify these neighbors more efficiently. To this end, we can first use binary search tree (BST) Hibbard (1962); Douglas (1959); Windley (1960) data structure to identify a smaller list of points from that lie within a sphere of radius equal to the largest semi axes of the cell. From this smaller list of points, we then proceed to delete the ones that satisfy .

This can be further refined by binning the template points in such that each bin contains about templates. For each such subset, the removal of points from proceeds by identifying those points that lie in the same bin (within some acceptable margin in ) and then applying the above strategy to eliminate the points. This refinement is made possible due to the fact that match values decrease monotonically with the difference in parameters and is most sensitive to changes in . Additional binning in other two coordinates may possibly improve the computational efficiency of removing the points from even more.

iii.1.3 Global co-ordinate transformation

In dimensionless chirp time coordinates, most of the ellipsoidal cells have semi-axes ratio around . This implies that the number of points inside this flat and elongated ellipsoid is a small fraction of the total number of points contained in the sphere with radius equal to semi-major axis. This undermines the efficiency of deleting the points using the technique outlined in Section III.1.2. We use a conformal coordinate transformation in such a way that one of the cells transforms to a unit sphere. Note that due to curvature effects, the same transformation when applied to other cells does not guarantee them to change to unit spheres but mitigates the problem of highly asymmetric semi-axes ratio to a large extent. We use the eigenvalues and eigenvectors of the metric at a fiducial point to transform the coordinates () as,


where, is a diagonal scaling matrix and the rotational matrix is constructed from the component of the eigenvector of the metric calculated at the fiducial point. This is carried out by first binning the templates in . The transformation in Eq. (9) is carried out using the metric of the centre point in each bin.

The metric at all other points transform as:


It is trivial to check that the metric at the fiducial points transforms to a identity matrix.

This global coordinate transformation Eq. (9) causes other nearby cells to become almost spherical with semi-axes in the ratio . In Fig. 1 we show the effect of global coordinate transformation on nearby cells. This is doubly advantageous: not only does the ball-point query volume for the BST searches decrease (leading to fewer points to deal with); the points that lie within the inscribed minimal match ellipsoid now occupy a large fraction of its volume leading to increased efficiency in removal of points.

(a) Before global transformation
(b) After global transformation
Figure 1: An example of global transformation Eq. 9) to speed up the point elimination part using BST algorithm. Here we construct the global transformer using the metric of an ellipse which is inscribed in blue dashed circle of the left figure as shown in the left panel. After the transformation the same cell becomes a unit sphere and nearby cells also become more spherical.

iii.1.4 Locally placed truncated octahedral lattice

The template placement problem in 3D can be mapped to the tiling problem such that a minimum number of similar cells is used without leaving any region uncovered. As such, truncated octahedral lattice becomes a natural choice for this problem.

As stated earlier, the TO shares 14 faces with its neighbours in the lattice. Assuming a TO inscribed in a sphere of radius , the coordinates of its 14 lattice neighbors are available in Table 4. Here is an index on the coordinates of the neighbor. When mapped to the template placement problem, we need to consider TO’s inscribed within the minimal match ellipsoids. In this case, the coordinates of the lattice neighbors can be calculated using the rotation and scaling matrices:


Fig. 2 shows the 14 neighbours of a TO inscribed in a elliptical cell.

Figure 2: Lattice neighbors of a truncated octahedron inscribed in a elliptical cell. Here only 9 neighbours are visible and the remaining 5 are on the opposite side.

Due to boundary effects, all 14 TO neighbors of a point need not necessarily be part of the template bank. The following conditions must be checked for:

  1. The point is inside the deemed parameter space and also satisfies . The latter corresponds to the condition that ’s can be inverted to yield physical masses.

  2. The point is not within the minimal match distance of existing points in and .

Check (b) above ensures that we do not double count the neighbors.

As shown in Algorithm (1), we start from a random point in the parameter space (by seeding ) and tessellate with local TO lattice. Because of curvature and boundary effects, it is not guaranteed that these tessellations cover the entire parameter space. This is marked by the exhaustion of as the placement proceeds. At this stage, we need to re-seed and continue the process iteratively until all points in have been used up.

iii.1.5 Choice of initial point and variations in banksize

The template placement algorithm starts from a randomly chosen point in the parameter space by seeding which is copied over as the first element of the template bank list .

One can start from any point in parameter space: in the present work we have started from “mid-point” corresponding to component masses with individual spins . We have checked that starting from other points (e.g. the extremities of the parameter space) also work quite well. This initial choice results in a minor fluctuation of the template bank size and it is recommended that we start from a point well inside the deemed parameter space where the local variations in curvature are less. This ensures maximum tiling before we need to re-seed .

We demonstrate these fluctuations in template bank size by constructing several template banks for compact binary systems whose bank parameters are given in Set II of Table 1 starting from different locations in the deemed parameter space. The random seed used to initialize was kept the same to eliminate bias. The five different starting points were taken to be the centre and extremities of the parameters space respectively and as expected, different choices of the initial seed point resulted in slightly different number of templates in the final bank. We generated an average of templates with fluctuation which is quite insignificant.

Iv Test of the Algorithm

In this section we construct a template bank using our algorithm and demonstrate its validity for CBC searches. We also compare its performance against a bank generated using the vanilla stochastic method available in the LSC Algorithm Library Suite (LALSuite) LIGO Data Analysis Software Working Group (). We need a metric over the dimensionless chirp coordinates and at present two such models are available for use, namely; and IMRPhenomB signal model Kalaghatgi et al. (2015). Our demonstration makes use of the model.

iv.1 Construction of the template banks for reduced-spin binary system

We generate a template bank using metric in parameter space for aligned spin compact binary system. The range of parameters is chosen such that the mass of the first object lies between with dimensionless spin magnitude in the range . The mass of the second object is taken between with dimensionless spin magnitude in the range . The NS-BH boundary mass condition is satisfied i.e any object with mass is considered to be a neutron star with corresponding limit on spin magnitude. This template bank can be used for binary neutron star (BNS) and neutron-star-blackhole (NSBH) searches.

We construct several template banks by varying the number of uniformly sprayed random points over the parameter space, . We consder sizes of varying between points. In all the cases, the placement proceeds from a point corresponding to individual masses and individual spin magnitudes . The full specification of input arguments for template generation are given in Set-I of Table 1.

The corresponding final template bank size are listed in Table 2. In Fig. 3 we show the two dimensional projections of such a bank along and directions respectively. The hat over the symbols denote unit normal vector along the specified direction. For comparison, we also generated a template bank using the vanilla stochastic method for the same parameters as listed in Set-I of Table 1. The vanilla stochastic bank placement was terminated at a point when there were 1000 rejected proposals for each accepted proposal averaged over the last 10 acceptances. This stochastic bank was found to contain templates which is larger than the geometric-random bank. The computational run time was also recorded and we observed that the geometric-random method took minutes while the vanilla stochastic method took minutes to execute on a single, unloaded processor, which is times faster.

Bank Parameter Set-I Set-II
Waveform model
Noise model
Lower cut-off frequency
Higher cut-off frequency
Mass of first object
Mass of second object
Spin of first object
Spin of second object
Size of
Minimal Match 0.97 0.97
Table 1: Parameters used to generate the geometric-random and vanilla stochastic banks. The results for different sizes of are summarized in in Fig. 4. In Set-I, the parameter space is chosen by satisfying the NS-BH boundary mass i.e. components with individual mass are identified as NS with dimensionless spin magnitude in the range .
(a) Projection along
(b) Projection along
(c) Projection along
Figure 3: Area normalized histograms of the template density in various planes. The template bank was constructed using hybrid geometric random algorithm presented in this paper. Each bin of the histogram was normalized by the square-root of the determinant of the metric to ensure equal area. The metric was calculated at the bin centre. The boundary effects are clearly seen. We can also see that bank is highly elongated along direction as compared to both and directions.

iv.2 Validation of the template banks

We investigate the performance of both the geometric-random and the vanilla stochastic template banks against a set of signal injections from the reduced spin signal model. In this section we summarize the results of this comparison and demonstrate that the two banks are nearly identical in performance.

Following Apostolatos Apostolatos (1995), the “fitting factor” is defined as a measure of the maximum match over the template bank for a putative injected signal :


The parameters of the injected signal (chosen from within the deemed parameter space over which the bank is placed) is chosen at random and may not coincide with that of a template point. The mismatch indicates the fractional loss of optimal SNR. Because differences between true waveforms and waveform models will always reduce the fitting factor, banks are usually tested to achieve fitting factors slightly larger than the minimal match used in their construction. The minimal match at which the template banks are constructed ( in this case) is somewhat arbitrary - and is arrived by carefully balancing the computational cost of the search against the desired detection efficiency.

For the banks generated as per the parameters in Set-I of Table 1, we injected signals from reduced-spin waveform model and calculated the fitting factor using Eq. (12). The lalapps_cbc_sbank_sim program as implemented in LALApps package of the LSC Algorithm Library Suite (LALSuite) LIGO Data Analysis Software Working Group () was used for this calculation. The range of mass and spin parameters of the injected signals were chosen to be the same as that of the bank and were drawn from a uniform distribution in this range. Other extrinsic parameter were drawn from (i) a uniform random distribution over all possible sky locations (ii) a fixed inclination angle corresponding to edge-on orientation of the plane of the binary, and (iii) a fixed luminosity distance of 1 Mpc. For the banks (corresponding to different initial sizes of ), we found that signals were recovered with fitting factor for both geometric-random bank as well as for the stochastic bank. The results are shown in Fig. 4. The fluctuation at higher significant digit was ignored due to the limited sample size of injections. These results are graphically depicted in Fig. 4 and demonstrate the equivalence of the two banks along with the computational speed up and corresponding efficiency of template placement using the new method.

Although the injection signal model is identical to the one used to construct the template banks, we still find that the fitting factors fall below the minimal match of the bank for a small fraction of the injections. This alludes to the fact that both template banks have holes in certain regions of the deemed parameter space where neighbouring templates do not provide adequate coverage. This arises from curvature effects.

(a) Comparison of hybrid bank performance
(b) Fitting factor distribution for hybrid bank
(c) Fitting factor distribution for hybrid bank
Figure 4: Panel (a)a shows a comparison of the geometric-random bank and the vanilla-stochastic bank constructed over identical parameter ranges is made by plotting fitting factors for signals. The horizontal line depicts the percentage of such injections for which the fitting factor is above the bank minimal match in the case of vanilla stochastic bank. The solid dots correspond to the percentage for geometric-random banks constructed with different sizes of as tabulated in Table 2. The bottom horizontal axis measures the computational speed-up of the hybrid geomtric-random placement while the one on top shows the corresponding efficiency in bank-size. Panel (b)b and (c)c show the cumulative distribution of the hybrid bank fitting factors for the different initial size of . They have been split into two panel for clarity.

The fitting factor depends on the parameters of the injected signals - e.g. masses, spin components, sky position and inclination angle of the binary’s plane to the line of sight etc. For non-precessing signals including only the dominant radiation mode, the fitting factor depends only on the intrinsic parameters. To understand the systematics, it is convenient to represent it as a function of two parameters while averaging out over the remaining ones.

(a) over component masses
(b) over total mass and reduced spin parameters
Figure 5: The left figure shows the minimum fitting factor for the hybrid template bank for a set of injected aligned spin NSBH, BNS signals as a function of component masses. The figure on the right shows the minimum fitting factor as a function of the reduced spin parameter and the total mass. Both the signal and template waveforms are modeled using the approximant. We see a region showing poor coverage corresponding to high total masses and high reduced spin values. This is also seen for the stochastic bank as well. One possible reason could be that the metric could be doing worse for systems with high masses and high reduced spin.

In Fig. 5LABEL:sub@fig_ff_aLABEL:sub@fig_ff_b, we show the histogram of “minimum” fitting factors Harry et al. (2014) (for geometric-stochastic template bank) over various combinations of intrinsic parameters of the compact binary system where both the signals and template waveforms are generated from waveform model. The bank performs well to match the injected signals throughout the parameter space for NSBH and BNS systems except for those regions where both the total mass and reduced spins are high.

In order to carry out a high-precision test comparing the efficiency of the hybrid bank with that of the vanilla-stochastic bank, we need a large number of injections to calculate the fitting factors upto high significant figures. To this end, we construct both the template banks using metric and parameters given in Set-II of Table 1. The geometric random bank was constructed by initializing with uniform random points in dimensionless chirp-time coordinates whereas the stochastic bank code was set to terminate when the rejection rate reached a value averaged over the last acceptances. The geometric random bank was found to contain templates whereas the stochastic bank contained ( more) template points. In this case, the geometric-stochastic bank took times less time than the vanilla-stochastic bank on a single unloaded processor.

We quantify the performance of these two template banks by computing fitting factors for two different injection families of aligned spin waveforms: and TaylorF2. injections were made in both cases where the intrinsic parameters of the injected waveforms were chosen from Set-II of Table 1 and other parameters chosen as earlier. For the case of geometric-random template bank created using local TO lattice, we found signals were recovered below a fitting factor of 0.97 for injections, and signals are recovered below this level for TaylorF2 injections. In the latter case, the injection model is different from the one used to construct the template bank - hence it is expected that the coverage for TaylorF2 will be less. The corresponding numbers for the vanilla-stochastic bank are found to be ( injections) and (TaylorF2 injections) respectively. From these numbers it is evident that the geometric-stochastic bank is equally efficient as the vanilla-stochastic bank. These results are summarized in Fig. 6.

Figure 6: Fitting factors computed for various sets of aligned-spin signal families against geometric-random (GR) / vanilla stochastic template (ST) banks (see Set-II of Table 1 for parameters). The performance of geometric random bank when both template and injected signals are generated from (TF2RS) approximant (black solid line), both generated from TaylorF2 (TF2) approximant (black dashed line) are shown. The performance of vanilla stochastic bank when both the templates injected signals are generated from approximant (black dashed dot line) and when both are generated from TaylorF2 approximant (black dot line) are also shown. In this case, the vanilla stochastic bank has 25% more templates than the geometric random bank and can be placed about times faster in time.

As mentioned earlier, the lattice provides optimal coverage for conformally flat spaces in 3D. One may be tempted therefore to use such a lattice as a seed for stochastic template placement. We have investigated this approach and compared it with the geometric-random and vanilla-stochastic algorithms. At first, the deemed parameter space (corresponding to parameters in Set-II of Table 1) was covered by a lattice using the metric at a putative point to determine the dimensions of the unit cell. Using a point in this volume for which the unit cell had the smallest dimension, we placed lattice points to entirely cover the deemed parameter space. Using these as seed points for stochastic placement, the final bank size was found to have templates. This is marginally () smaller than the vanilla-stochastic bank which has templates as reported above. The geometric-random bank outperforms this seeded stochastic bank by a good margin of more than .

As remarked earlier, the vanilla-stochastic algorithm can be cast in two different ways. The traditional bottom-up approach has been implementation in the LSC Algorithm Library and has been used in this work for comparison with the geometric-random bank. In order to compare it head-to-head with the top-down approach we implemented it in software and ran it for the exact same parameters as given in Set II of Table 1. As expected, the top-down implementations gives nearly identical bank sizes ( difference in size) but takes less than half the time as the bottom-up LAL implementation. The computational speed comes from the fact that in this implementation, one can reject many proposals that lie within the minimal match ellipsoid centred around a single accepted proposal. Efficient computational data structures like binary search trees are readily available for such queries.

The summary of various templates banks referred to in the above discussion is available in Table 2.

Bank Parameters Placement Algorithm Size of Bank Size Execution Time (min) Comments
Set-I of Table 1
more templates
Hybrid construction
Vanilla-Stochastic (lalapps_cbc_sbank)
Set-II of Table 1 Hybrid construction Geometric-Random
more templates
Vanilla-Stochastic (lalapps_cbc_sbank)
seeded Stochastic
Vanilla-Stochastic                   (top_down)  134, 426
Table 2: Summary of various template banks mentioned in this paper. The semi-analytic metric for signal model was used in all cases. The Vanilla Stochastic algorithm can also work by directly calculating matches (instead of using ) but in this head-to-head comparison, we have used the semi-analytic metric which speeds up the Vanilla Stochastic template placement significantly. The usage of the metric is compulsory for Geometric-Random placement.

V Discussion and outlook

Templated matched filtering is the mainstay of gravitational wave detection pipelines. With unprecedented improvement in low frequency sensitivity of advanced detectors, and the availability of theoretical spinning waveform models, it has become imperative to conduct these searches over increasingly larger volumes in higher dimensional parameter spaces. For such cases, the stochastic algorithm is used for template placement as it is easily scalable to higher dimensions. But it is computationally expensive and by design not the most space efficient.

This paper introduces a new template placement algorithm in 3D with an attempt to combine the space efficiency of lattice along with the robustness of stochastic placement algorithm. Such a template bank can be used in gravitational wave searches from binary neutron stars and neutron-star-blackhole compact binary systems where the waveform is described by two mass parameters and a mass-weighted spin magnitude parameter providing coverage for aligned-spin systems.

The truncated octahedron which is a Dirichlet-Voronoi polytope of the lattice, inscribed within the minimal match ellipsoid is used as a unit cell for the geometric placement. Such lattices are known to provide optimal coverage for conformally flat 3D spaces. While the signal manifold is not globally flat, one can assume local flat patches and use such cells to cover them optimally. The interface to the stochastic placement algorithm is made by spraying random points over the parameters space - which are removed if found within the minimal match ellipsoid of any template. We discuss how this merger of methods is able to handle the issues arising out of varying curvature and irregular boundaries. The nuances of its optimal implementation are discussed in detail. We show in a direct comparison with stochastic algorithm that the new method generates significantly fewer templates and is computationally more efficient in Table2. We now make a couple of comments related to the new template placement algorithm presented in this paper.

One of the key issues with the geometric template placement algorithms (e.g. geometric hexagonal placement Cokelaer (2007)) is related to the amount of fine tuning needed in the method to account for curvature effects. For example, in 2D the hexagonal cells used to place templates need to be whittled or pared to a slightly smaller size so that no holes are left due to the changes in relative orientation of neighboring hexagonal cells. This reduces the overall efficiency of placement leading to increase in template bank sizes. But more importantly, one loses the ability to make generic codes that are usable for different waveform models. In this regard, the geometric-random algorithm proposed here is robust against such fine-tuning by design. The template placement proceeds by first spraying a large number of random points over the parameters space which are later removed if found within the minimal match ellipsoid of the templates in the bank. Suppose a small portion of the deemed space is left uncovered due to curvature effects, it would lead to some residual points that are not removed from . As evident from Algorithm 1, these points are revisited in subsequent iterations where a random point (out of the residuals) is added to the template bank leading to complete coverage.

We would like to point out that the template banks shown in this paper make use of the waveform family which model only the inspiral part of the evolution of the compact objects. They are not ideal for BBH searches directly, where a significant part of the SNR is expected to be contributed from the merger and ringdown phases. However the method developed in this paper is quite general and can be used for placing templates for IMR waveform families as well. If analytical metrics for IMR approximants (IMRPhenomDKhan et al. (2016); Husa et al. (2016), SEOBNRv2Taracchini et al. (2014); Pürrer (2014), etc.) were available, the proposed method would be an elegant and efficient solution for covering the entire combined BNS+NSBH+BBH space. This will thereby mesh well with modern LIGO searches which are carried out over a large range of parameters with combined BNS+NSBH+BBH banks using a mix of inspiral- only and IMR templates.

The assumption that the metric on the signal manifold is slowly varying and is locally flat is crucial for space efficiency. One can imagine a hypothetical case where this assumption does not hold true (e.g. metric coefficients are random at every point) in which case, the geometric-random algorithm will effectively fall back to the top-down version of the vanilla-stochastic template placement by design. In other words, the number of templates in the bank from the new method will not exceed the stochastic template bank in the limiting case. We have also shown in a direct comparison that the top-down stochastic bank implementation is computationally more efficient and should be used where the metric is available. Incorporating an intelligent way of spraying the random proposals (instead of drawing them from a uniform distribution) over the parameter space may lead to further optimization of this method.

On the other hand if the metric was perfectly flat (and given in coordinates where it was constant), the hybrid construction would fall back on a perfect lattice.

Finally, the geometric-random placement method presented here for 3D is generically scalable to signal manifolds in higher dimensions by using the appropriate lattices.

Appendix A The signal model for gravitational waves from inspiraling compact binary coalescence

The TaylorF2 “reduced spin” waveform model in frequency domain is given byAjith (2011)


where, the amplitude depends on the component masses, distance to the source, sky position, orientation of the binary’s plane, and is the instantaneous phase which can be explicitly written as:


where is the time of arrival of the signal at the detector marking the epoch at which the instantaneous frequency takes a fiducial value, is the corresponding phase, is the instantaneous velocity, is the total mass and is the symmetric mass ratio of the binary and is the Euler gamma constant.

The spin effects are encoded through , and which appear at 1.5PN, 2PN and 2.5PN phase terms respectively and are given by:


where the reduced spin parameter is defined as the weighted sum of individual spins and of the component masses as:


The individual spins of the components are the projections of their spin vectors along the Newtonian orbital angular momentum vector and defined as:


Appendix B Space filling truncated octahedron

A polyhedron is a three dimensional solid that has a finite number of polygon faces. One can fill a 3D space completely without any overlap or gap through the tessellations of space-filling polyhedra. Examples of such space-filling polyhedron include cube, hexagonal prism etc. Solution to such space-filling problems find many practical applications like optimal placement of a network of communication towers (Alam and Haas, 2006). The template placement problem addressed in this paper can also be mapped to an optimal space-filling problem in curved space. The geometric properties of optimal space-filling polyhedra can be understood from the volumetric quotient defined as


where is the volume of polyhedron, is maximum distance from its center to any vertex and is the surface area of the polyhedron. Note that is the inverse of the thickness, a common measure for the efficiency of a given covering. The polyhedra with the highest is expected to provide the optimum coverage.

The TO is a 14-faced Archimedean solid, with 8 hexagonal faces, 6 square faces and has 24 vertices. It is generated by joining two regular pyramids upside down and cutting a pyramid from all six vertices in such a way that the length of all the sides generated are equal. Thus a truncated octahedron of side can be created by removing six pyramids of side from an octahedron of side . Fig. 7 shows the TO obtained from two pyramids. The geometric properties of a TO and pyramid are given in Table 3.

Figure 7: Truncated Octahedron created by truncating six pyramids from the six vertices of two square base pyramids.
Name Value
Side length of Pyramid
Side length of TO
Height of Pyramid
Height of TO along square face
 Height of TO along hexagonal face
Height of TO along vertices
Volume of TO
Surface Area of TO
Volumetric Quotient of TO
Table 3: Geometrical properties of truncated octahedron and pyramid where the truncated octahedron constructed using two truncated pyramid.

Suppose a three dimensional volume is covered by the tessellations of TO cells. Each face of such a cell can be shared by another neighboring TO cell. There are two kind of neighboring cells: ones that share the square faces which we shall call S-neighbors and the others that share the hexagonal faces which we shall refer to as H-neighbors. Each TO has a maximum of six S-neighbors and eight H-neighbors. The distance between S-neighbors is twice the height of TO along square face and for H-neighbors, twice of height of TO along the hexagonal face. When a TO is inscribed in a sphere of radius such that the -axis goes through the centre of one of the square faces, then, the sides of the squares and hexagons is given by , the distance from the centre to each of the S-neighbors is equal to , and the distance to the H-neighbors is . The coordinates of all the 14 neighbours are listed in Table 4.

Neighborhood Type Position
square faces ,
hexagonal faces
Table 4: These coordinates are the neighborhood positions of a TO where a TO is inscribed in a sphere with radius , where the center of a TO is placed at origin and the -axis goes through the centre of one of the square faces.
A.S. would like to thank LS for motivation. SR thanks IIT Gandhinagar for research fellowship and fellow graduate students (Amit Reza, Chakresh Singh, Md. Yousuf) for useful discussions and help with the manuscript. We thank Samarth Vaijanapurkar (B.Tech. 2016 student of IIT Gandhinagar) for useful discussions. The authors thank Duncan Brown for valuable feedback. Last but not least, the authors thank the anonymous referees whose careful reading of the manuscript, detailed queries and many suggestions have helped to improve the content and presentation of this work.


Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description