Fast Fourier transform

`Y = fft(x)`

Y = fft(X,n)

Y = fft(X,[],dim)

Y
= fft(X,n,dim)

The functions `Y = fft(x)`

and `y = ifft(X)`

implement the transform
and inverse transform pair given for vectors of length *N* by:

$$\begin{array}{l}X(k)={\displaystyle \sum _{j=1}^{N}x}(j){\omega}_{N}^{(j-1)(k-1)}\\ x(j)=(1/N){\displaystyle \sum _{k=1}^{N}X}(k){\omega}_{N}{}^{-(j-1)(k-1)}\end{array}$$

where

$${\omega}_{N}={e}^{(-2\pi i)/N}$$

is an *N*th root of unity.

`Y = fft(x)`

returns the
discrete Fourier transform (DFT) of vector `x`

,
computed with a fast Fourier transform (FFT) algorithm.

If the input `X`

is a matrix, `Y = fft(X)`

returns the Fourier
transform of each column of the matrix.

If the input `X`

is a multidimensional array, `fft`

operates
on the first nonsingleton dimension.

`Y = fft(X,n)`

returns the `n`

-point
DFT. `fft(X)`

is equivalent to `fft(X, n)`

where `n`

is
the size of `X`

in the first nonsingleton dimension.
If the length of `X`

is less than `n`

, `X`

is
padded with trailing zeros to length `n`

. If the
length of `X`

is greater than `n`

,
the sequence `X`

is truncated. When `X`

is
a matrix, the length of the columns are adjusted in the same manner.

`Y = fft(X,[],dim)`

and ```
Y
= fft(X,n,dim)
```

applies the FFT operation across the dimension `dim`

.

A common use of Fourier transforms is to find the frequency components of a signal buried in a noisy time domain signal. Consider data sampled at 1000 Hz. Form a signal containing a 50 Hz sinusoid of amplitude 0.7 and 120 Hz sinusoid of amplitude 1 and corrupt it with some zero-mean random noise:

Fs = 1000; % Sampling frequency T = 1/Fs; % Sample time L = 1000; % Length of signal t = (0:L-1)*T; % Time vector % Sum of a 50 Hz sinusoid and a 120 Hz sinusoid x = 0.7*sin(2*pi*50*t) + sin(2*pi*120*t); y = x + 2*randn(size(t)); % Sinusoids plus noise plot(Fs*t(1:50),y(1:50)) title('Signal Corrupted with Zero-Mean Random Noise') xlabel('time (milliseconds)')

It is difficult to identify the frequency components by looking
at the original signal. Converting to the frequency domain, the discrete
Fourier transform of the noisy signal `y`

is found
by taking the fast Fourier transform (FFT):

NFFT = 2^nextpow2(L); % Next power of 2 from length of y Y = fft(y,NFFT)/L; f = Fs/2*linspace(0,1,NFFT/2+1); % Plot single-sided amplitude spectrum. plot(f,2*abs(Y(1:NFFT/2+1))) title('Single-Sided Amplitude Spectrum of y(t)') xlabel('Frequency (Hz)') ylabel('|Y(f)|')

The main reason the amplitudes are not exactly at 0.7 and 1
is because of the noise. Several executions of this code (including
recomputation of `y`

) will produce different approximations
to 0.7 and 1. The other reason is that you have a finite length signal.
Increasing `L`

from 1000 to 10000 in the example
above will produce much better approximations on average.

`fft`

supports inputs of data types `double`

and `single`

.
If you call `fft`

with the syntax ```
y =
fft(X, ...)
```

, the output `y`

has the same
data type as the input `X`

.

[1] Cooley, J. W. and J. W. Tukey, "An
Algorithm for the Machine Computation of the Complex Fourier Series,"*Mathematics
of Computation*, Vol. 19, April 1965, pp. 297-301.

[2] Duhamel, P. and M. Vetterli, "Fast
Fourier Transforms: A Tutorial Review and a State of the Art," *Signal
Processing*, Vol. 19, April 1990, pp. 259-299.

[3] FFTW (`http://www.fftw.org`

)

[4] Frigo, M. and S. G. Johnson, "FFTW:
An Adaptive Software Architecture for the FFT,"*Proceedings
of the International Conference on Acoustics, Speech, and Signal
Processing*, Vol. 3, 1998, pp. 1381-1384.

[5] Oppenheim, A. V. and R. W. Schafer, *Discrete-Time
Signal Processing*, Prentice-Hall, 1989, p. 611.

[6] Oppenheim, A. V. and R. W. Schafer, *Discrete-Time
Signal Processing*, Prentice-Hall, 1989, p. 619.

[7] Rader, C. M., "Discrete Fourier
Transforms when the Number of Data Samples Is Prime," *Proceedings
of the IEEE*, Vol. 56, June 1968, pp. 1107-1108.

Was this topic helpful?