Summer Research Fellowship Programme of India's Science Academies

Understanding the cosmos with gravitational modelling

Preet Baxi

Dwarkadas J. Sanghvi College of Engineering, Vile Parle, Mumbai 400056

Dr. Chandra Kant Mishra

Indian Institute of Technology Madras, Chennai 600025


The production and detection of the Gravitational Waves is a marvellous blend of science and engineering creating opportunities and expanding the horizons to the unknowns of the universe. When any mass in space undergoes acceleration it creates ripples in space-time fabric known as Gravitational Waves. These ripples are tidal forces from the sources billions of light years away traveling transparent through the universe which are recorded as strain on the test masses in the detectors. The detector data containing the signal on further analysis leads to detection of gravitational wave and understanding properties of the sources that generated it. This data when verified with the mathematical derivations of Einstein's General Relativity equations; theorized almost a century ago, potentially provides answers to questions that have puzzled astronomers for decades. The fellowship began with reading various review papers for understanding the properties, sources and solutions of Einstein’s equations. The components of detectors and the engineering involved in it were studied in detail. The LIGO open data workshop was followed to understand the step-by-step processing of the data from the detectors aided with the lectures, tutorials and exercises uploaded. From the monochromatic and chirp waveforms, the basic astrophysical quantity such as luminosity distance to the binary was calculated. Amplitude, luminosity distance and time of chirp was derived for spinning neutron stars from its basic equations. Further RMS response of antenna pattern functions for the interferometer was calculated. Also, frequency domain waveform was derived using the stationary phase approximation. Gaining interest in working on deducing and understanding the properties of a waveform with respect to its source of generation, I decided to work on parameter estimation of gravitational wave sources. The aim of the project is to compare errors in extracting source properties by considering different parameterizations of expected gravitational wave signals.

Keywords: space-time, tidal forces, relativity, LIGO, parameter estimations

Abbreviations and Symbols

GW Gravitational Waves
LIGO Laser Interferometer Gravitational-Wave Observatory
LISA Laser Interferometer Space Antenna
NS Neutron Star
BH Black Hole
BBH Binary Black Hole
SPA Stationary Phase Approximation
CBC Compact Binary Coalescence
PN Post Newtonian
SNR Signal-to-Noise Ratio
PSD Power Spectral Density
PE  Parameter Estimations
G Gravitational Constant
c Speed of light
M Mass of sun
Mc Chirp Mass
η Symmetric Mass Ratio
μ Reduced Mass Ratio
M Test Mass
ν Velocity of wave
h Strain signal
θ Polar Angle
ϕ Azimuthal angle
ψ Polarization angle
ι Inclination Angle
 A Amplitude of Wave
f Fourier Frequency
F+ , Fx Detectors' antenna pattern functions
Γij Fisher Information Matrix


The project focuses primarily on the parameter estimations for BH-BH, NS-BH and NS-NS sources detected through Advanced LIGO, VIRGO and Initial LIGO. The source parameters: Amplitude (A), coalescence time (tc), phase at coalescence (ϕc), chirp mass (Mc) and symmetric mass ratio (η) were extracted using Fisher Matrix Approach; details are mentioned further in section 7.1. Following the discussion in [5], I further tried to re-parametrize the waveform used for PE in this report, details of this parametrization appears in section 7.2.


Physics of Gravitational Radiation

Gravity is one of the fundamental forces that is produced by anything that has mass. The mass causes the curvature in space which in turn governs the motion of other masses in the neighborhood. The idea is central to Einstein’s famous equation of General Relativity. Gravity follows the principles of Einstein’s relativity, mainly; that nothing can travel faster than the speed of light. Just like acceleration of charge produces propagating electromagnetic fields, accelerating masses produce waves of gravity. The amplitude for these waves is calculated as:

h    2  ϵ  GMrc2h\;\sim\;2\;\epsilon\;\frac{GM}{rc^2} (1)​​

​​where, ϵ is a dimensionless quantity which, in some cases is also written as ν2N.S./c2 where νN.S is the non-spherical velocity inside the source.

Unlike the Electromagnetic waves, GWs are produced by bulk motions of matter and not by individual atoms, carrying a different kind of information which travels almost unabsorbed due to their weakness. The information carried by them has conceptually more data as the sources are individual masses adding the total effect compared to the ones produced by charges which cancel out each other due to their opposite polarities. The Gravitational Waves hence gives us the complete data about the sources that produce them, enabling us to make their direct observation; leading to verification of Einstein’s General Relativity Equation. To signify the importance and future scope of research in Gravitational Waves we can consider the following comparison: the observations from the Cosmic Microwave Background (Electromagnetic Waves) give us information of the universe around 105 years after the Big Band; but the GWs, if detected would give us information of around 1024 seconds after the Big Bang. Quadrupole Approximation

The weak gravitational fields are studied in post-Newtonian expansions which include higher-order corrections to Newtonian gravity originating from General Relativity. This method has yielded the most powerful insights of different emission sources, resulting in the fundamental Quadrupole Formula which gives an approximation for weakly relativistic systems. The quadrupole formula is similar to the dipole formula of electromagnetism, approximating emission of gravitational radiation from a generic system undergoing acceleration. It is in the form of a time-dependent integral. But, in case of general relativity, the time derivative of the dipole moment is calculated just as the integral of velocity as its mass-energy conservation. This dipole moment is constant and hence there is no energy radiated indicating a constant field in the dipole. Thus, to find the exact radiation through General Relativity we have to go beyond the dipole approximation known as quadrupole (second moment of mass distribution). It is defined to be a spatial tensor matrix given as:

Qjk=ρ  xj  xk3xQ_{jk}=\int\rho\;x_{j\;}x_kⅆ^3x (2)​

Here, ρ is the Mass-Energy density and xj, xk are the cartesian co-ordinates. This tensor matrix describes the source of the GWs and the waves are represented as:

hjk=2Gc4r  d2Qjkdt2h_{jk}=\frac{2G}{c^4r}\;\frac{d^2Q_{jk}}{dt^2} (3)

It is observed that the h is not a scalar quantity, instead, a tensor. Since the geometry of the gravitation waves represents stretching and squeezing of space-time fabric, this tensor contains the whole distortion information.

Luminosity in Gravitational Waves

The luminosity in GWs is the energy carried by the waves is measured as the square of time-derivative of the wave amplitude, hence, depends on the squares of the components d3Qjkdt3\frac{d^3Q_{jk}}{dt^3}(See Eq.(3)). The luminosity contains a factor of G/c5 and 1/5 concerning dimensional grounds and further General Relativity calculations respectively leading to a final equation in the quadrupole approximation is given as:

L=Gc5(jk  Q...jkQ...jk13Q...2)L=\frac G{c^5}\left({\textstyle\sum_{jk}}\;{\overset{...}Q}_{jk}{\overset{...}Q}_{jk}-\frac13\overset{...}Q^2\right) (4)

Where Q is the trace of the matrix Qjk . Its squared third time-derivative is subtracted to ensure that spherical motions do not radiate.

Radiation from sSpinning Neutron Star

According to Newtonian Gravity, an absolutely spherical body has the same gravitational field outside it over every point in space at a constant distance and the body is treated as a point mass. Thus, even when these bodies spin, they do not cause any disturbance in the gravitational field hence generating no gravitational radiation. So, to calculate radiation from a neutron star of radius R, spinning with frequency f, it is assumed to have a small bulge of mass m present on its surface, having the non-spherical velocity as νN.S =2πRf .

Thus, the gravitational field generated by the bulge. Hence on substituting the properties concerning it in Equation (1), we get the radiation amplitude as:

hbump  2  (2πRfc2)Gmrc2h_{bump}\sim\;2\;\left(\frac{2\pi Rf}{c^2}\right)\frac{Gm}{rc^2} (5)

For luminosity, we assumed roughly four comparable components of Qjk contributing to the sum in Equation (4):

Lbump  (G5c5)(2πf)6m2R4L_{bump}\sim\;\left(\frac G{5c^5}\right)(2\pi f)^6m^2R^4 (6)

The radiated energy passed on as gravitational waves hence come from the rotation of the neutron stars leading their spin-down on a timescale of:

tspindown=12mν2  (Lbump  )    54πf1(GmRc2)1(νc)3t_{spindown}=\frac12m\nu^2⁄\;(L_{bump\;})\;\sim\;5⁄4\pi f^{-1}(Gm⁄Rc^2)^{-1}(\nu⁄c)^{-3} (7)

Spinning neutron stars have the weakest of gravitational waves which are difficult to detect but they are continuous sources. So, if the gravitational waves are detected from neutron stars over a longer period, promising data will be produced.

Radiation from a Binary Star System

Consider a Binary Star System with the stars of same masses M rotating as a centrifuge in a circular orbit of radius R and having ν2N.S.=GM/4R. For this system, the gravitational wave amplitude, from Equation (1) is:

hbinary12GMRc2Gmrc2  h_{binary}\sim\frac12\frac{GM}{Rc^2}\frac{Gm}{rc^2}\; (8)

The gravitational wave luminosity for such a system, comparing to neutron star bulges, given as:

Lbinary180c5G(GMRc2)5L_{binary}\sim\frac1{80}\frac{c^5}G\left(\frac{GM}{Rc^2}\right)^5 (9)​

The energy released through the rotation orbital motion results in shrinkage of orbit. The shrinking makes the observed gravitational wave increase in time. This is called a chirp. The time-scale for this is:

tchirp=(Mν2)/Lbinary  20GM/c3(GM/Rc2)4t_{chirp}=(M\nu^2)/L_{binary}\sim\;20GM/c^3(GM/Rc^2)^{-4} (10)

Hence, if the timescale for an equal mass chirping binary is known then one can determine both M and R and further if the amplitude is detected then the luminosity distance to the binary can also be known. This also holds true for binaries with different masses for which the measurable mass is known as the Chirp Mass denoted as:


Where, M1 and M2 are the individual masses, μ is the reduced mass of the binary system and MT is the total mass.

Also, for a circular binary system, where the components are treated as point masses, the gravitational waves are represented in simple form as:

h(t)=A(t)    Cos  ϕ(t)h(t)=A(t)\;\;Cos\;\phi(t) (11)

where, h(t) is the gravitational waveform (gravitational strain), A(t) is the time dependent amplitude and ϕ(t) is the gravitational wave phase. The amplitude can be expressed as:

A(t)=2(GMc)5/3c4rπPgw(t)2/3A(t)=\frac{2(GM_c)^{5/3}}{c^4r}\frac\pi{P_{gw}(t)^{2/3}} (12)

where, Pgw is the gravitational wave period and r is the luminosity distance to the binary. Further, for circularized binary systems, orbital wave period and gravitational wave period are related as: Porb(t)=2Pgw(t). The factor of two is because the lowest possible order of gravitational radiation is the quadrupole order and quadrupole moments are invariant under 180° rotation.


Following the Equivalence principle, to calculate the effect of gravitational waves we always consider more than two particles that come in the path of their propagation. The travelling tidal forces or the gravitational waves are transverse waves which create stretching and squeezing in these particles due to changes in the gravitational field, recorded as strain. This disturbance in the ensemble of particles occurs in two directions as shown; known as polarization of gravitational waves.

    Plus and Cross Polarizations

    Consider an ensemble of particles in the plane of the paper assumed in the x-y plane as shown. The gravitational waves travel perpendicular to the plane of the paper i.e. moving in z-direction hence producing degrees of freedom along x-x, x-y and y-y. In General Relativity, gravitational radiation is represented by a second rank, symmetric trace free tensor hjk having only components {hxx, hyy, hxy} and the remaining set to zero. Therefore, when there is stretching and squeezing along the x and y axes it is called as Plus polarization (h+) and when it occurs along the x=y and x=-y lines it is called as Cross polarization (hx). The polarization of GWs from a source depends on the orientation of the dynamics inside the source relative to the observer and is also a time-dependent quantity. Therefore, measuring the polarization provides information about the orientation of a particular source in space.


    Engineering of LIGO

    The basic idea that came up for the detection of GWS was to create a machine that could detect the change in distance of two freely falling bodies. As these two bodies will be in a gravitational field, so the measure of non-uniformities or changes in this field known as the tidal forces detects the presence of gravitational wave. Hence, a GW is a measure of strain (h) in the distance between two test masses. LIGO (Laser Interferometer Gravitational-Wave Observatory) was so built.


      The large ground base detectors like LIGO, Virgo, etc. are principal beam interferometers. The basic principle of interferometers is that when two lasers ‘Interfere’ with each other creating a pattern that can be analyzed. It is a complete vacuum L-shaped beam splitter with a laser input and a photodetector. A typical bar detector consists of a cylinder of aluminum with a length l~ 4m, a very narrow resonant frequency between f~500Hz and 1.5KHz, and a mass M~1000Kg. A short gravitational wave burst with h ~ 10-21 will make the bar vibrate with an amplitude δlgw  hl  1021m\delta l_{gw}\sim\;hl\sim\;10^{-21}m. The arms are made long because of the effect that GWs have on matter is almost negligible. The amplitude or the distance of the shift h of an object due to gravitational waves is fixed by nature at 10-21 . The shift in the distance is given by δl  hl\delta l\sim\;hl. So, h is a constant, we see that the resulting shift is completely dependent upon the initial length between the objects. Therefore, the longer the better.

      From the source, a laser beam is shot with hits a mirror known as a beam splitter. This mirror splits the beam into two beams. One travelling straight through and the other deflected at an angle of 9090^\circ. At the two ends, which are placed at exactly the same distance from the beam splitter, more mirrors reflect the light back and they meet or interfere at the middle mirror merging into a single beam which is then shot on a screen where the light is studied. These mirrors are known as test masses. Since, both the beams of light travel the exact same distance, the interference that will occur is fixed- the two beams with either merge to reform the same pre-split beam or there will be no light at all. Any change in this interference results in a change in the behavior of the beam. This alerts the observer to the presence of any change that may be caused by gravitational waves. In addition to the fully reflective test masses mirrors, two more partially reflective mirrors are placed to create additional reflection paths of the laser beams as shown. Effectively, the beam is reflected about 280 times before being merged. This makes the arms longer than their actual 4Km length to around 1120Km. The laser begins at 4 Watts, enters the arms at 200Watts, finally reaching the optimum power level of 750KWatts. This high power is achieved by a power recycling mirror placed in between the light source and the beam splitter. The laser is recycled in the sense that when the light reflects back from the end mirrors on the beam splitter, it is directed towards the additional mirror which throws it back into the interferometer thereby boosting the power.

      Thus, finally based on the interference pattern on the screen the gravitational waves are detected and after working on the final productive output is evolved.

      Detection through Interferometers

      The wave amplitude h(t) from detectors is the one that has exRand  eyRe_x^Rand\;e_y^Ras the axes of the ellipse formed during polarization. Thus the full wave amplitude is described by the wave tensor as:

      h(t)=h+(t)e++hX(t)eXh(t)=h_+(t)e_++h_X(t)e_X (13)

      where  e+  and  eX  arewhere\;e_+\;and\;e_{X\;}are the polarization tensors. These on the basis of the sky-plane reference with respect to the detectors is also written as:

      ϵ+=e+cos2ψ+eXsin2ψ  and  ϵX=e+sin2ψ+eXcos2ψ\epsilon_+=e_+cos2\psi+e_Xsin2\psi\;and\;\epsilon_X=-e_+sin2\psi+e_Xcos2\psi (14)

        relative orientation of sky and detector

        Following further the detectors orientation and the angles in the sky-plane as shown, the strain can be further written as:

        h=δll=F+h++FXhXh=\frac{\delta l}l=F_+h_++F_Xh_X (15)

        where F+  and  FXF_+\;and\;F_Xare the antenna pattern functions for the two polarizations defined on the sky-plane basis as:

        F+=1/2(1+cos2θ)cos2ϕcos2ψcosθsin2ϕsin2ψF_+=1/2(1+cos^2\theta)cos2\phi cos2\psi-cos\theta sin2\phi sin2\psi

        FX=1/2(1+cos2θ)cos2ϕsin2ψ+cosθsin2ϕcos2ψF_X=1/2(1+cos^2\theta)cos2\phi sin2\psi+cos\theta sin2\phi cos2\psi (16)

        Finally, the polarization amplitudes from an inspiraling binary, a rotating neutron star, or a ringing blackhole, takes a simple form as follows:

        h+=ho/2(1+cos2ι)  cosΦ(t),  hx=ho  cosι  sinΦ(t)h_+=h_o/2(1+cos^2\iota)\;cos\Phi(t),\;h_x=h_o\;cos\iota\;sin\Phi(t) (17)

        Where hoh_ois an overall time dependent amplitude, Φ(t)\Phi(t)is the signal’s phase and ι\iotais the angle made by the characteristic direction in source with the line of sight. Therefore, Equation (15) can be rewritten as:

        h(t)=Ahocos(Φ(t)Φo)h(t)=Ah_ocos(\Phi(t)-\Phi_o) (18)

        where A=(A+2+AX2)1/2,  tanΦo=A+/AX,  A+=1/2F+(1+cos2ι)  and  AX=FXcosι.A=(A_+^2+A_X^2)^{1/2},\;tan\Phi_o=A_+/A_X,\;A_+=1/2F_+(1+cos^2\iota)\;and\;A_X=F_Xcos\iota.

        Global and Spatial Approach

        To extend the horizons of research further; the global and spatial approach is adopted. This includes detection of gravitational waves using different combinations of various detectors set-up at different locations all over the world to be used as network. This increases the efficiency of the data and also the probability of viewing maximum sources in the sky. These multiple detectors include LIGO in Hanford and Livingston (USA), VIRGO (Italy), KAGRA (Japan), INDIGO (India), GEO (Germany) and TAMA (Japan).

        Further the original set of Initial LIGO and VIRGO have been made more efficient and are named as Advanced LIGO and Advanced VIRGO.

        Furthermore, there a plan-in-action to set up a space based interferometer called as LISA. The benefits of which are removal of terrestrial noises, natural vacuum and extremely large lengths of beams to be used for detection. Hence, increasing the efficiency of the data received.



          The Stationary Phase Method(SPA) is a procedure for evaluation of integrals of the form

          I=F(x)ej(x)dxI=\int_{-\infty}^\infty F(x)e^{-j\varnothing(x)}dx (19)

          Where ∅(𝑥) is a rapidly-varying function of x over most of the range of integration, and F(x) is slowly varying (by comparison). Such integrals frequently arise in radiation or scattering problems. Rapid oscillations of the exponential term is approximately zero over those regions of the integrand; the only significant non-zero contributions to the integral occur in regions of the integration range where d∅/dt=0 i.e. at points of stationary phase.

          Hence on following this and considering Fourier transform for the time-domain generic waveform:

          h(t)=A(t)  cos((t))h(t)=A(t)\;\cos(\varnothing(t)) (20)

          the Fourier transform will be:

          h~(f)=h(t)exp(2πift)dt          =1/2A(t)(exp(2πift+i(t))+exp(2πifti(t)))dt\begin{array}{l}\widetilde h(f)=\int_{-\infty}^\infty h(t)exp(2\pi ift)dt\\\;\;\;\;\;=1/2\int_{-\infty}^\infty A(t)(exp(2\pi ift+i\varnothing(t))+exp(2\pi ift-i\varnothing(t)))dt\end{array} (21)

          We consider the positive frequency, f > 0.

          When the phase ∅(𝑡) is monotonically increasing function of t, 𝑑∅/𝑑𝑡>0. We then have,

          ddt  (2πft+)=2πf+(d)dt0\frac d{dt}\;(2\pi ft+\varnothing)=2\pi f+\frac{(d\varnothing)}{dt}\neq0 (22)

          Further, we assume that at t=𝑡𝑜 following relation is satisfied.

          d/dt  (2πft)=2πf(d)/dt=0d/dt\;(2\pi ft-\varnothing)=2\pi f-(d\varnothing)/dt=0 (23)

          Furthermore, let 𝑀=𝑚1+𝑚2, 𝜇=𝑚1𝑚2𝑀⁄, 𝜈=𝜇𝑀⁄, Ω: orbital angular frequency, F:orbital frequency and Ω=2πf. Then,

          x(GM/c3)2/3=(2πGMF/c3)2/3x\equiv(GMΩ/c^3)^{2/3}=(2\pi GMF/c^3)^{2/3} (24)

          where we have:

          h~12A(t)  exp(2πiftiϕ(t)dt)\widetilde h\cong\frac12\int_{-\infty}^\infty A(t)\;exp(2\mathrm{πift}-\mathrm{iϕ}(\mathrm t)\mathrm{dt})

          12A(t0)e2πift0i(t0)(d2dt2t=t0)1/2(2π)1/2eπ4i\cong\frac12A(t_0)e^{2{\mathrm{πift}}_0-\mathrm i\varnothing({\mathrm t}_0)}\left(-{\left.\frac{d^2\varnothing}{dt^2}\right|}_{t=t_0}\right)^{-1/2}(2\mathrm\pi)^{1/2}e^{\frac{\mathrm\pi}4i}

          =12A(t0)ei(2πft0(t0)π4)(dFdtt=t0)1/2=\frac12A(t_0)e^{i(2{\mathrm{πft}}_0-\varnothing({\mathrm t}_0)-\frac{\mathrm\pi}4)}\left(-{\left.\frac{dF}{dt^{}}\right|}_{t=t_0}\right)^{-1/2} (25)​

          ​and hence the Energy and Energy emission rate equation via the gravitational wave is given as:

          2019-08-14 (2).png
            Energy and Emission rate equation
            2019-08-14 (4).png
              phase equation 


              The modelling of gravitational waves is the analysis of the data from the detectors. The waves are modelled in the General Relativistic framework to generate a data set which gives information of the source and the wave formed. An efficient data analysis scheme involves two independent aspects: first, the theoretical computational of very high accuracy templates and second, the design strategy adapted to the particular signal one is looking for. These strategies vary according to the type of signal. Gravitational waves from inspiraling binaries are transients lasting for a short duration in the sensitivity bandwidth of a ground-based detector. As the binary evolves, the waveform sweeps up in frequency and amplitude, leading to a characteristic chirp signal. As the phasing of the waves is known very accurately, it is possible to enhance their detectability by using matched filtering.

                Chirp waveform

                Among the different methods suggested for the detection of chirps for the inspiraling and merging binaries, matched filtering is the most effective technique. Matched filtering consists of passing the detector data through a linear filter, or a template, constructed from the expected signal h(t; θ). Here θ is a ‘vector’ whose components are the parameters of the template. The templates h(t; θ) generally use the restricted waveform where for binaries in quasi-circular orbits the phase is computed at the highest PN order available, but the amplitude is taken to be Newtonian, involving only the dominant signal harmonic at twice the orbital frequency. This is different from the complete waveform, which incorporates the PN corrections to the amplitude, arising from the ‘plus’ and ‘cross’ GW polarizations, and hence includes the contribution from other harmonics (both higher and lower) besides the dominant one. Till date, for non-spinning binaries, the restricted waveform is computed to 3.5PN accuracy and the complete waveform up to 2.5PN order. The best template is probably the one which consists of the phasing at 3.5PN and the amplitude at 2.5PN.

                In matched filtering, the unknown set of parameters characterizing the signal are measured by maximising the correlation of the data with a whole family of templates which correspond to different values of the parameters. The parameters of the template which maximises the output of a matched filter is an estimate of the true parameters. The parameters of a signal measured in a single experiment will be different from the actual values due to the presence of noise. Parameter estimation basically aims at computing the probability distribution for the measured values of a signal. Given a measured value from a single experiment one then uses the probability distribution function to compute the interval in which the true parameters of the signal lie at a specified confidence level.


                The main objective of the project was to understand the theory involving modelling of a waveform studying in particular the process of Matched Filtering. A particular waveform was studied and the error values for the 1PN order were derived. This calculation was performed for each Advanced LIGO, initial LIGO and VIRGO interferometers individually. All the three cases: BH-BH, NS-BH and NS-NS were studied and the values were tabulated. The main idea was to compare the efficiency of the data from each detectors and hence find improvement methods related to them. The scope in the first case was restricted to only the 1st Newtonian order and the amplitude considered while solving the frequency domain waveform was the equivalent amplitude for each case with the primary term. The SNR in this case was fixed at 10 hence the values were calculated for the fixed SNR system. The masses of the Neutron star and the Black hole considered were 1.4 Mand 10 M respectively. The table from[4] was reproduced with five point accuracy values and compared. The phasing formula used increases in steps of 0.5PN. The error bounds were derived through Cramer-Rao methods as it gives efficient results at higher SNR.

                Further, it was studied that for improved results a network configuration should be used. Another method called Monte Carlo was studied in further detail and the given frequency domain waveform was then linearized using this method. This method is mainly applied for a network of detectors. Here in this project, this waveform with further parameters was used to derive error bounds for Advanced LIGO; BH-BH case hence proving that this method is not reasonable for single interferometers but for networks as a whole. When in first case of Cramer-Rao bound, the preferred part of the waveform is only the inspiral phase whereas in the case of Monte Carlo the complete waveform can be studied together. But here we restrict our scope to the inspiral part only. Initially only five parameters are involved and later we increase it to a nine parameter space.


                The output of a gravitational wave detector contains both the signal and noise and is schematically represented as:

                x(t) = h(t) + n(t) ( 27)

                where x(t) is the signal registered and n(t) is the noise, which is assumed to be a stationary Gaussian random variable, with zero mean, i.e.,

                n(t)\overline{n(t)} = 0 (28)

                Here an overbar denotes the ensemble average (over many realisations of the noise or, equivalently, over an ensemble of detectors). Let q(t) define a linear filter and c(t) its correlation with the detector output x(t)

                c(t)=-dt'x(t')q(t+t') (29)

                Define a new quantity [q](t), such that c(t) is normalized by the square root of its variance,

                where 𝑥̃(f) and 𝑞̃(f) are the Fourier transforms of x(t) and q(t), respectively, 𝑆ℎ(f) is the real, one-sided power spectral density defined only for positive frequencies

                And 𝑛̃(𝑓) is the Fourier transform of n(t) defined as 𝑛̃(𝑓)=∫ⅆ𝑡 𝑛(𝑡)𝑒−2𝜋𝑖𝑓𝑡∞−∞. the filtered SNR is defined by the ensemble average:

                ​An optimal filter is the one which maximizes the SNR at a particular instant, say t = 0, and is given by the matched filtering theorem as:

                where 𝛾 is an arbitrary real constant. Thus, the SNR corresponding to the optimal filter is given by:

                ρ2=40h(f)~2Sh(f)df\rho^2=4\int_0^\infty\frac{\left|\widetilde{h(f)}\right|^2}{S_h(f)}\operatorname df

                For a given incident gravitational wave, different realizations of the noise will give rise to somewhat different best-fit parameters. However, if the SNR is high enough, the best-fit parameters will have a Gaussian distribution centred around the actual values.

                Let 𝜃̃ 𝑎denote the ‘true values’ of the parameters and let 𝜃̃ 𝑎+Δ𝜃𝑎 be the best-fit parameters in the presence of some realization of the noise. Then for large SNR, errors in the estimation of parameters Δ𝜃𝑎 obey a Gaussian probability distribution of the form:


                Where 𝑝(0)is a normalization constant. In the above expression Γ𝑏𝑐=(ℎ𝑎|ℎ𝑏), where ℎ𝑎≡𝜕ℎ/𝜕𝜃𝑎, is the Fisher information matrix evaluated at the measured value 𝜃̂ of the parameters 𝜃 . Here, ( | ) denotes the noise weighted inner product. Given any two functions g and h their inner product is defined as :

                (gh)20dfg~(f)h~(f)+g~(f)h~(f)Sh(f)\begin{array}{l}(g\left|h\right.)\equiv2\int_0^\infty df\frac{\widetilde g^\ast(f)\widetilde h(f)+\widetilde g(f)\widetilde h^\ast(f)}{S_h(f)}\\\end{array}

                Using the definition of the inner product one can re-express Γ𝑎𝑏 more explicitly as:

                Γab=20ha~(f)  hb~(f)+ha~(f)hb~(f)Sh(f)df\Gamma_{ab}=2\int_0^\infty\frac{\widetilde{h_a}^\ast(f)\;\widetilde{h_b}(f)+\widetilde{h_a}(f)\widetilde{h_b}^\ast(f)}{S_h(f)}df

                The variance-covariance matrix, or simply the covariance matrix, defined as the inverse of the Fisher information matrix, is given by

                    abθa  θb=(Γ1)ab\textstyle\sum_{\;\;}^{ab}\equiv\left\langle\triangle\theta^a\;\triangle\theta^b\right\rangle=(\Gamma^{-1})^{ab}

                where〈 .〉 denotes an average over the probability distribution function in Eq. (36). The root-mean-square error 𝜎𝑎 in the estimation of the parameters 𝜃𝑎 is

                σa=(θa)21/2=    aa\sigma_a=\left\langle\left(\triangle\theta^a\right)^2\right\rangle^{1/2}=\sqrt{\textstyle\sum_{\;\;}^{aa}}

                and the correlation coefficient 𝑐𝑎𝑏 between parameters 𝜃𝑎 and 𝜃𝑏 is defined as

                cab=θaθbσaσb= ab aa bb

                Since the start, we choose the set of independent parameters 𝜃 describing the GW signal to be:


                where 𝑡𝑐 refers to the coalescence time, ∅𝑐 refers to the phase at coalescence instant, 𝑓0 is a scaling frequency related to the power spectral density (PSD) of the detectors.

                The upper cut-off frequency in computing the integrals is taken to be the GW frequency at the last stable circular orbit (LSO) given, for a test mass in a Schwarzschild spacetime of mass M, to be:


                The noise PSD of advanced LIGO is given as:

                Sh(f)=S0[x4.145x2+111(1x2+x4/2)(1+x2/2)],ffsS_h(f)=S_0\left[x^{-4.14}-5x^{-2}+\frac{111(1-x^2+x^4/2)}{(1+x^2/2)}\right],f\geq f_s

                Where 𝑥=𝑓/𝑓0,𝑓0=215𝐻𝑧 ,𝑓𝑠=20𝐻𝑧 and 𝑆0=10−49. Here, 𝑓0 is the scaling frequency and 𝑓𝑠 is the lower cut-off frequency.

                The noise PSD of initial LIGO is given as:

                𝑆(𝑓)=𝑆0[(4.49𝑥)−56+0.16𝑥−4.52+0.52+0.32𝑥2] ,𝑓≥𝑓𝑠

                Where 𝑥=𝑓/𝑓0,𝑓0=150𝐻𝑧 ,𝑓𝑠=40𝐻𝑧 and 𝑆0=9 x 10−46𝐻𝑧−1.

                The noise PSD of advanced LIGO is given as:

                𝑆(𝑓)=𝑆0[(6.23𝑥)−5+2𝑥−1+1+𝑥2] ,𝑓≥𝑓𝑠

                Where 𝑥=𝑓/𝑓0,𝑓0=500𝐻𝑧 ,𝑓𝑠=20𝐻𝑧 and 𝑆0=3.24 x 10−46𝐻𝑧−1.


                The parameters considered in this method are given by 𝜗 as:


                We map four of the six extrinsic signal parameters: {𝐴, 𝜓,𝜄,𝜑0} into new parameters 𝑎𝑘 with k=1,…,4. Such that a given signal can be written as:

                h(t)=k=14akhk(t)h(t)={\textstyle\sum_{k=1}^4}a^kh_k(t) (47)​

                where the ℎ𝑘(𝑡)𝑠 are completely independent of those four extrinsic parameters. (The two remaining extrinsic parameters are the sky-position angles.) To deduce their dependencies as well as the forms of the 𝑎𝑘𝑠 we begin by noting that the antenna-pattern functions can be treated as the components of a vector that are related to two sky position dependent functions, 𝑢(𝜃,∅) and 𝑣(𝜃,∅) through a two-dimensional rotation by 2𝜓:

                (F+Fx)=(cos2ψsin2ψsin2ψcos2ψ)(uv)\begin{pmatrix}F_+\\F_x\end{pmatrix}=\begin{pmatrix}\cos2\psi&\sin2\psi\\-\sin2\psi&\cos2\psi\end{pmatrix}\begin{pmatrix}u\\v\end{pmatrix} (48)​

                With this well-known observation, one finds

                ℎ1(𝑡)∝ 𝑢(𝜃,∅)cos[𝜑(𝑡)], ℎ2(𝑡)∝ 𝑣(𝜃,∅)cos[𝜑(𝑡)],

                ℎ3(𝑡)∝ 𝑢(𝜃,∅)sin[𝜑(𝑡)], ℎ4(𝑡)∝ 𝑣(𝜃,∅)sin[𝜑(𝑡)] (49)

                where the proportionality factor is a dimensionless (mass-dependent) function of time.

                The new parameters are themselves defined as:

                Ma(a1a3a2a4)=y(ι)dLOφ0.  I  .  O2ψM_a\equiv\begin{pmatrix}a^1&a^3\\a^2&a^4\end{pmatrix}=\frac{y(\iota)}{d_L}O_{\varphi_0}.\;I\;.\;O_{2\psi} (50)​

                Where 𝑦(𝜄)≡[(1+𝑐𝑜𝑠2𝜄)2+4𝑐𝑜𝑠2𝜄]1/2, 𝑂𝛼 is the two dimensional orthonormal rotation matrix for angle 𝛼 and

                I=((1+cos2ι)/y(ι)00(2cosι)/y(ι))I=\begin{pmatrix}(1+\cos^2\iota)/y(\iota)&0\\0&(2\cos\iota)/y(\iota)\end{pmatrix} (51)

                The Fisher matrix is then computed on the space (M,𝜂,𝜃,∅,𝜓,𝜄,𝑡0,𝑎1,𝑎2,𝑎3,𝑎4) The errors in the are obtained by inverting that matrix. By using error propagation equations obtained from Eq. (49), we are able to deduce error estimates for all four extrinsic parameters.


                The calculations for the parameter estimations were performed on the WOLFRAM MATHEMATICA 11.3 platform. The code was written for the equations above and were solved keeping in mind the scope assumed. The results then derived and compared.


                The error bounds of the parameters were generated as expressed in[4] and was tabulated as:

                2019-08-14 (13).png

                  Here the values of 𝑡𝑐 and ∅𝑐 are in milliseconds and radians respectively.

                  SUMMARY AND CONCLUSION

                  The tabulated data was compared and studied and hence it was proved that the sensitivity of the VIRGO interferometers is better than the other interferometers in single detector system. The degeneracy in the parameters estimated is minimum therefore asserting its efficiency and providing newer regions of improvisations where required. But since this case is restricted to only the inspiral part of a waveform, for better efficiency, calculations should be performed using Monte-Carlo Estimations. This method on the other hand does not gives reasonable results when parameters are estimated using only one detector but gives efficient outcomes for network of interferometers. Therefore, it was inferred that to achieve extremely efficient results one should consider a network of detectors using the Monte Carlo approach and calculate for a complete chirp waveform. Further the higher modes of the waveform can also be included with a complete 3.5PN order of Equations.


                  [1] Bernard F. Schutz (2000) “Gravitational Radiation”.

                  [2] B.S. Satyaprakash and Bernard F. Schutz (2009) “Astrophysics and cosmology with Gravitational waves”, Living Reviews in Relativity.

                  [3] Hideyuki Tagoshi (2013) “Inspiral waveforms with the Stationary Phase Approximation”.

                  [4] K.G. Arun et al. (2008) “Parameter Estimation of inspiraling compact binaries using 3.5 post-Newtonian gravitational wave phasing: the non-spinning case”.

                  [5] P. Ajith and Sukanta Bose (2009) “Estimating the parameters of non-spinning binary black holes using ground-base gravitational wave detectors: statistical errors”.

                  [6] Hideyuki Tagoshi et al. (2018) “Parameter Estimation of neutron star-black hole binaries using an advanced gravitational wave detector network: Effects of the full post-Newtonian waveform”.


                  ·      Books

                  1.       A first course in General Relativity - Bernard F. Schutz

                  2.      Gravitational Waves Volume 1: Theory and Experiments – Michele Maggiore

                  3.      Gravitational Waves: How Einstein’s space-time ripples reveal the secrets of then universe – Brian Clegg

                  ·      Magazine “Gravitational Waves” (November 2017) – Digit Dmystify

                  ·      Open Data Workshop

                  LIGO Scientific Collaboration – Open data workshop #1 (2018)



                  Foremost, I would like to express my sincere gratitude to my guide Dr. Chandra Kant Mishra, for the amazing journey through space-time. The knowledge and methodologies that he has imparted in me has benefited me throughout my fellowship and will be a boon for my future. His guidance, motivation, patience and belief in me resulted in the satisfactory completion of research work and this report.

                  Besides, I would also like to put forward my hearty thanks to the Indian Academy of Sciences for organizing such a platform where research enthusiasts like me can join and set the beginning of our research career. It is indeed a great privilege to be selected for the summer fellowship(2019).

                  I would like to thank my college professor Prof. Shashikant Auti for posting my reference to the academy and for always helping me out in all possible ways to not let any opportunities that come in my way go vain. Further, I would like to thank my HOD Dr. K.N. VijayaKumar, Vice-Principal Dr. Ashish C. Daptardar and Principal Dr. Hari Vasudevan for permitting me for the fellowship and their constant support.

                  I would like to thank all of my research colleagues helping me throughout the research.

                  Above all, I would like to thank my parents, Mr. Pradip Baxi and Mrs. Vasanti Baxi as without their blessings nothing could have been achievable. Their trust and confidence in me is my key to success.

                  Finally, I would like to thank my Family and Friends for the love and enthusiasm they showed, which became my strength and made this research more exciting.

                  Written, reviewed, revised, proofed and published with