Example: optimal program

0 This program reads the 40-meter data stream, and then filters it though a chirp template corresponding to a pair of inspiraling $1.4M_\odot$ neutron stars.

The correspondence between different arrays in this program, and the quantities discussed previously in this section, is given below. In these equations, $\Delta t = 1/{\tt srate}$ is the sample time in seconds, and $\Delta f = (n \Delta t)^{-1} = {\tt srate/npoint}$ is the size of a frequency bin, in Hz. Here $n={\tt npoint}$ is the number of points in the data stream which are being optimally filtered in one pass.

Chirp templates (in frequency space) for the two polarizations are related to the arrays chirp0[ ] and chirp1[ ] by

$$\tilde T_{0}(f) = {\Delta t \over {\tt HSCALE}} \;{\tt chirp0[\;]}$$ (6.30.145)

$$\tilde T_{90}(f) = {\Delta t \over {\tt HSCALE}} \; {\tt chirp1[\;]}$$ (6.30.146)

where the elements chirp0[2j] and chirp0[2j+1] are the real and imaginary parts at frequency $f= j \Delta f$ (with the exception of the Nyquist frequency, stored in chirp0[1]). Note that to ensure that quantities within the code remain within the dynamic range of floating point numbers, we have scaled up the template strain by a constant factor HSCALE; we also scale up the interferometer output by the same factor, so that all program output (such as signal-to-noise ratios) is independent of the value of HSCALE. If you're not comfortable with this, go ahead and change HSCALE to 1. It won't change anything, provided that you don't overflow the dynamic range of the floating point variables! The scaled interferometer response function is

$${\tt response[\;]} = {\tt HSCALE}/{\tt ARMLENGTH} \times R(f),$$ (6.30.147)

where the function $R(f)$ is defined by equation (). The Fourier transform $\tilde h$ of the dimensionless strain is obtained by multiplying $\Delta t$ and the FFT of channel.0 by response[ ], yielding

$$\tilde h(f) = {\Delta t \over {\tt HSCALE}} {\tt htilde[\;]}.$$ (6.30.148)

The one-sided noise power spectrum $S_h(f)$ is the average of

$$S_h(f) = {2 \over n \Delta t} \vert \tilde h(f) \vert^2 = {2... ...lta t \over n \; {\tt HSCALE}^2} \vert{\tt htilde[\;]}\vert^2.$$ (6.30.149)

The power spectrum $S_h(f)$ is averaged using the same exponential averaging technique described for the routine avg_spec(). This average is stored as

$${ S_h(f)} = {2 \Delta t \over n \; {\tt HSCALE}^2 } \langle ... ...\Delta t \over n \; {\tt HSCALE}^2 } {\tt mean\_pow\_spec[\;]}$$ (6.30.150)

Twice the inverse of this average is stored in the array twice_inv_noise[ ], so that

$${2 \over S_h(f)} = { n \; {\tt HSCALE}^2 \over \Delta t} \; {\tt twice\_inv\_noise[\;]}.$$ (6.30.151)

The expected noise-squared for the plus polarization is given by equation ():

$$\begin{eqnarray*} \langle N^2 \rangle &=& {1 \over 2} (Q,Q) = {1 \over 2} \int_{... ...ate(\cdots,chirp0[\;],chirp0[\;],twice\_inv\_noise[\;],npoint)}. \end{eqnarray*}$$

where the subscript on the inverse FFT means ``at zero lag", and ``$\rightarrow f$" means ``returned by the call to the function f". We have chosen a distance for the system producing the ``chirp" $\tilde T_(f)$ so that the expected value of $\langle N^2 \rangle=1$.

In similar fashion, the signal $S$ at lag $t_0$ is given by

$$S = ({\tilde h \over S_h}, \tilde Q)$$

$$= \int_{-\infty}^\infty df\; {\tilde h(f) \tilde T^*_{0}(f) \over S_h(f)} \; {\rm e}^{-2 \pi i f t_0}$$

$$= {1 \over n \Delta t} FFT^{-1}_{i}\left[ {\Delta t \over {\tt HSCA... ...t HSCALE}^2 \over \Delta t } \; {1 \over 2} {\tt twice\_inv\_noise[\;]} \right]$$

$$= FFT^{-1}_{i} \left[ {\tt htilde[\;]} \; {\tt chirp0[\;]} \; {1 \over 2} {\tt twice\_inv\_noise[\;]} \right]$$

$$\rightarrow {\tt correlate(\cdots,htilde[\;],chirp0[\;],twice\_inv\_noise[\;],npoint)},$$ (6.30.152)

where now the subscript on the FFT means ``at lag $t = i \; \Delta t$".

You might wonder why we have been so careful - after all, both the signal and the noise, as we've defined them, are dimensionless, so it's not surprising that all of the factors of $\Delta t$ drop out of the final formulae for the signal and the expected noise-squared. The main reason we've been so long winded is to show exactly how the units cancel out, and to demonstrate that there aren't any missing dimensionless constants, like npoint, left out of the program. Some sample output from this program is shown in the next section. Includes/optimal.tex