As I've demonstrated in earlier posts in this series, the fundamental idea behind the Fourier series & transforms is that ALL functions can be decomposed in waves of various frequencies - therefore EVERY SIGNAL will exhibit frequencies in a PSD. The question becomes which frequencies are physically relevant and not an artifact of the data selection and distribution.

An important issue in working with the PSD is to understand just how it behaves with different types of known data as input. Generally, if your data contains a clean, sinusoidal signal, the PSD will generate nice, sharp peaks at the frequency of the signal. But as we will see below, even a little bit of deviation from the perfect, infinite signal changes the PSD.

A note on normalization: For comparisons, I want to maintain a matching scale between analytically computed power spectra using Fourier transforms with the power spectra using the FFT. I have adopted the convention of the previous post (Quantized Redshifts. IV. The Fourier Transform, FFT & PSD) for analytic computations. I have also included a multiplicative factor for normalizing the PSD generated by the FFT. In Python using the numpy library, this can be expressed as

>>> import numpy as np

>>> from numpy import fft

FFT computation. nbins=number of data bins

>>> ff=fft.fft(dataf,nbins)

>>> psd=np.abs(ff)**2/nbins

Analytic computation

>>> psdref=(nbins/(2.0*np.pi))*np.abs(fftref)**2

**PSD of Pure Sinusoidal Signal:**

First we'll consider the basic case of the sine wave. The ideal sine wave has infinite in duration in time, t, with some well-defined frequency, \omega_0. It is expressed mathematically by

But in the real world, we cannot have a sine wave of infinite duration. All data runs are of finite duration in time and space. Therefore, we will define a sine wave of finite duration, starting at a time \tau and lasting for time, T. Mathematically, this can be expressed as:

We can compute the Fourier transform, and compute the FFT of this system, and we plot them below.

Click to enlarge |

In this example, the infinite sine wave produces a PSD that is a delta function (wikipedia) at the frequency of the wave, which I've designated with red in the graphic above. Note that in the above plot, we are plotting the frequency, f, not the angular frequency, \omega. But the PSD generated by the FFT (blue) and the analytic function (green) are slightly different. We see this more clearly in the plot below where we plot the amplitude of the PSD on a logarithmic scale.

Click to enlarge |

**PSD of A Couple of (Independent) Sine Waves**

Next we'll consider the superposition of two overlapping sine waves of finite duration. Because the Fourier Transform is a linear transformation, the Fourier transform of the sum of the Fourier transform of the two input waves. However, when we compute the PSD, we must still compute the product with the

*total*Fourier transform. This means that the PSD of two signals is not necessarily equal to the sum of the PSDs of the individual signals.

Click to Enlarge |

Note that while in this example, the PSD is equal to the sum of the PSDs of the two input signals, this will not always the case for reasons we will see next week.

**PSD of a Square Pulse**

Finally, we consider a rectangular pulse, which has a constant amplitude, different from zero, for some finite duration, T. To make the inquiry more general, I'll also consider offset from time t=0 by a time \tau.

For the analytic part of the analysis, we compute the Fourier Transform of this function:

And now we compute the power spectrum using the function above, and the FFT. The results are plotted below.

Click to Enlarge |

Here again are frequencies in the power spectrum from an input that is clearly perfectly uniform within its boundaries, it has no 'internal structure'.

This is a very counter-intuitive result. What generates all these peaks?

Whatever the origin of the peaks, it means we must be very careful when claiming peaks in a PSD represent actual internal structure in a dataset. To make sure we understand what our analysis techniques are doing, we must test the null hypothesis (Wikipedia) - what kind of PSD is produced if we have a uniform region subjected to the same sampling as our data?

Here's another question. Suppose we change the offset value, \tau. How will the Fourier Transform change? How will the PSD change? Suppose the signal is not a square pulse, will the answers change? How is this question relevant to cosmological power spectral analysis?

**Some General Observations**

- We see that the power spectrum computed by the Fourier Transform and the FFT are indeed different, but still exhibit many of the same qualitative features, such as the locations of local peak values. Some of this difference is due to aliasing, in the Discrete Fourier Transform (DFT).
- Narrow features in time or space are broad in frequency and wave number. Well-defined frequencies correspond to infinite waves. Broad structures in spatial dimensions generate low-frequencies in the power spectra.

Enough of the introductions. I've now covered enough material to explore some specific errors. Time to get to the meat of our story...

Next Weekend:

*Common Power Spectral Density "Bloopers"*

## No comments:

Post a Comment