# How to Compute the PSD

March 29, 2018

Back to: Random Testing

When computing the power spectral density (PSD) of a signal, it is essential to determine the frequency spectrum. Engineers can do so using the Fourier series, which is named after the mathematician Joseph Fourier who invented the method in 1807.

### Fourier Series

The Fourier series is based on the principle that any signal *x* of length *T* can be constructed by the sum of a series of sine and cosine waves, each with an integer number of cycles in *T*. Mathematically, this is written as:

(1)

Equation 1

where *f _{1}* = 1 /

*T*and

*f*=

_{i}*i*

*f*.

_{1}*C*and

_{i}*S*are the Fourier coefficients to be determined. Figure 2.6 illustrates this graphically.

_{i}In practice, the signal *x* is digitized into a sequence of *N* numbers (*x _{n}*,

*n*= 1 to

*N*) with a time interval of Δ

*t*between samples (Figure 2.7).

### Fourier Coefficients

We can determine the value of each Fourier coefficient by multiplying both sides of Eq. (1) by the corresponding sine or cosine wave and computing the mean value over *N* = *T*/Δ*t*. Equation 2, for example, determines the value of *C*_{37}.

(2)

(3)

Equation 2

The following summation computes the mean value and *t _{n}* =

*n*Δ

*t*:

(4)

This procedure works because the different sine and cosine functions are independent of each other. That is, the mean value of the product of any two different functions is equal to zero, as illustrated in Figure 2.3. The terms in Equation 2 have a mean value of zero except for *C*_{37}. Then:

(5)

(6)

(7)

Equation 3

In Equation 3, the mean-square value of the cosine wave is 1/2 of what has been used (as illustrated in Figure 2.2). In general, the coefficients are determined by the following equation (Equation 4):

(8)

(9)

(10)

Equation 4

### Fast Fourier Transform

Historically, the computation of the summations in Equation 4 was very time-consuming until 1965, when Cooley and Tukey published an efficient algorithm for this known as the fast Fourier transform (FFT).

Their algorithm uses a sequence of *N* digital samples where *N* is a power of 2 (e.g. 2, 4, 8, 16, 32, etc.). The highest frequency that can be resolved in the FFT is *f _{q}* = 1 / (2 Δ

*t*) Hz (called the Nyquist frequency), so there are

*N*/ 2 frequency values. The signal

*x*must be conditioned by a low pass filter before it is digitized to make sure that frequencies above

*f*have been removed.

_{q}The FFT computes the values of the Fourier coefficients *C _{i}* and

*S*. The mean square amplitude is given by (

_{i}*C*

_{i}^{2}+

*S*

_{i}^{2}) / 2. Since the frequency resolution is Δ

*f*= 1 /

*T*= 1 / (

*N Δ*

*t*), the PSD is normalized to a 1Hz bandwidth by dividing the mean square amplitude by Δ

*f.*

(11)

Equation 5

The process of generating a PSD must address a number of FFT computational details. These FFT details include digital time sampling, aliasing, windowing, frequency resolution, and averaging.