Problem Set 6 for Fourier transforms/matched filters

1) Write a function that will shift an array by an arbitrary amount using

a convolution (yes, I know there are easier ways to do this). The function

should take 2 arguments – an array, and an amount by which to shift the

array. Plot a gaussian that started in the centre of the array shifted by half

the array length.

2) a)The correlation function f ?g is R

f(x)g(x+y)dx. Through a similar

proof, one can show f ? g = if t(df t(f) ∗ conj(df t(g))). Write a routine to

take the correlation function of two arrays. Plot the correlation function of

a Gaussian with itself.

b) using these results, write a routine to take the correlation function

of a Gaussian (shifted by an arbitrary amount) with itself. How does the

correlation function depend on the shift? Does this surprise you?

3) The circulant (wrap-around) nature of the dft can sometimes be problematic. Write a routine to take an FFT-based convolution of two arrays

without any danger of wrapping around. You may wish to add zeros to the

end of the input arrays.

4) DFTs work very nicely out of the box when there are an integer

number of periods of a wave in the region analyzed. Sadly, when we are

dealing with real data, we usually are forced to analyze a finite chunk of data,

and there will in general be no particular relation between the frequencies

in the data and the interval we’re analyzing. We’ll look at the effects of this

a bit now.

a) Show that:

N

X−1

x=0

exp(−2πikx/N) = 1 − exp(−2πik)

1 − exp(−2πik/N)

It may help to recognize that the sum can be re-written P α

x where α =

exp(−2πik/N) so we can treat it as the sum of a geometric series.

b) Show that this approaches N as k approaches zero, and is zero for

any integer k that is not a multiple of N.

c) We can use this to analytically write down the DFT of a non-integer

sine wave. Pick a non-integer value of k and plot your analytic estimate

of the DFT. Show that the FFT agrees (to within machine precision) with

your analytic estimate. Normally, we think of the Fourier transform of a pure

1

sine wave to be a delta function. Are we close to that? This phenomenon is

usually known as spectral leakage.

d) A common tool to get around this is the use of window functions.

The leakage essentially comes from the fact that we have a sharp jump at

the edge of the interval. If we multiply our input data by a function that

goes to zero at the edges, this cancels out the jump, and so prevents the

leakage from the jumps at the edges. Of course, since we have multiplied

by the window in real space, we have convolved by it in Fourier space. One

simple window we could use is 0.5 − 0.5 cos(2πx/N) (there are many, many

choices). Show that when we multiply by this window, the spectral leakage

for a non-integer period sine wave drops dramatically.

e) Show that the Fourier transform of the window is [N/2 N/4 0 … 0 N/4]

(either numerically or analytically). Use this to show that you can get the

windowed Fourier transform by appropriate combinations of each point in

the unwindowed Fourier transform and its immediate neighbors (you may

need to be careful with signs here, since if you work through the math, some

of the transforms need to be inverse FFTs). The choice of suitable windows

is as much art as science (and depends on the details of what you’re most

concerned about), but I hope this gives at flavor of what’s going on and will

be useful for matched filter results.

5) Matched Filter of LIGO data

We are going to find gravitational waves! Key will be getting LIGO data

from github:

https://github.com/losc-tutorial/LOSC_Event_tutorial

While they include code to do much of this, please don’t use it (although

you may look at it for inspiration) and instead write your own. You can

look at/use simple read ligo.py that I have posted for concise code to read

the hdf5 files. Feel free to have your code loop over the events and print the

answer to each part for that event. In order to make our life easy, in case we

have to re-run your code (which we should not have to do), please also have

a variable at the top of your code that sets the directory where you have

unzipped the data. LIGO has two detectors (in Livingston, Louisiana, and

Hanford, Washington) and GW events need to be seen by both detectors

to be considered real. Note that my read template functon returns the

templates for both the plus and cross polarizations, but for simplicity you

can just pick one of them to work with.

a) Come up with a noise model for the Livingston and Hanford detectors

separately. Describe in comments how you go about doing this. Please

2

mention something about how you smooth the power spectrum and how

you deal with lines (if at all). Please also explain how you window the data

(you may want to use a window that has an extended flat period near the

center to avoid tapering the data/template where the signal is not small).

b) Use that noise model to search the four sets of events using a matched

filter. The mapping between data and templates can be found in the file

BBH events v3.json, included in the zipfile.

c) Estimate a noise for each event, and from the output of the matched

filter, give a signal-to-noise ratio for each event, both from the individual

detectors, and from the combined Livingston + Hanford events.

d) Compare the signal-to-noise you get from the scatter in the matched

filter to the analytic signal-to-noise you expect from your noise model. How

close are they? If they disagree, can you explain why?

e) From the template and noise model, find the frequency from each event

where half the weight comes from above that frequency and half below.

f) How well can you localize the time of arrival (the horizontal shift of

your matched filter). The positions of gravitational wave events are inferred

by comparing their arrival times at different detectors. What is the typical

positional uncertainy you might expect given that the detectors area a few

thousand km apart?

3

Phys 512

# Problem Set 6 for Fourier transforms/matched filters

Original price was: $35.00.$30.00Current price is: $30.00.

## Reviews

There are no reviews yet.