Spectral leakage due to FFT is caused by: mismatch between desired tone and chosen frequency resolution, time limiting an observation. Understand the concept using hands-on examples.
Limits of frequency domain studies
Frequency Transform is used to study a signal’s frequency domain characteristics. When using FFT to study the frequency domain characteristics of a signal, there are two limits : 1) The detect-ability of a small signal in the presence of a larger one; 2) frequency resolution – which distinguishes two different frequencies.
In practice, the measured signals are limited in time and the FFT calculates the frequency transform over a certain number of discrete frequencies called bins.
Spectral Leakage:
In reality, signals are of time-limited nature and nothing can be known about the signal beyond the measured interval. For example, if the measurement of a never ending continuous train of sinusoidal wave is of interest, at some point of time we need to terminate our observation to do further analysis. The limit on the time is also posed by limitations of the measurement system itself (example: buffer size), besides other factors.
It is often said that the FFT implicitly assumes that the signal essentially repeats itself after the measured interval and hence the FFT assumes the signal to be continuous (conceptually, juxtapose the measured signal repetitively). This lead to glitches in the assumed signal (see Figure 1). When the measurement time is purposefully made to be a non-integral multiple of the actual signal rate, these sharp discontinuities will spread out in the frequency domain leading to spectral leakage. This explanation for spectral leakage need to be carefully investigated.
Experiment 1: Effect of FFT length and frequency resolution
Consider a pure sinusoidal signal of frequency \(f_x = 10 \;Hz\) and to represent in computer memory, the signal is observed for 1 second and sampled at frequency \(F_s=100 \; Hz\). Now, there will be 100 samples in the buffer and the buffer will contain integral number of waveform cycles (10 cycles in this case). The signal samples are analyzed using N-point DFT. Two cases are considered here for investigation : 1) The FFT size \(N\) is same as the length of the signal samples, i.e, N=100 and 2) FFT size set to next power of 2 that fits the signal samples i.e, N=128. The result are plotted next.
Fx=10; %Frequency of the sinusoid
Fs=100; %Sampling Frequency
observationTime = 1; %observation time in seconds
t=0:1/Fs:observationTime-1/Fs; %time base
x=sin(2*pi*Fx*t);%sampled sine wave
N=100; %DFT length same as signal length
X1 = 1/N*fftshift(fft(x,N));%N-point complex DFT of x
f1=(-N/2:1:N/2-1)*Fs/N; %frequencies on x-axis, Fs/N is the frequency resolution
N=128; %DFT length
X2 = 1/N*fftshift(fft(x,N));%N-point complex DFT of x
f2=(-N/2:1:N/2-1)*Fs/N; %frequencies on x-axis, Fs/N is the frequency resolution
figure;
subplot(3,1,1);stem(x,'r')
title('Time domain');xlabel('Sample index (n)');ylabel('x[n]');
subplot(3,1,2);stem(f1,abs(X1)); %magnitudes vs frequencies
xlim([-16,16]);title(['FFT, N=100, \Delta f=',num2str(Fs/100)]);xlabel('f (Hz)'); ylabel('|X(k)|');
subplot(3,1,3);stem(f2,abs(X2)); %magnitudes vs frequencies
xlim([-16,16]);title(['FFT, N=128, \Delta f=',num2str(Fs/128)]);xlabel('f (Hz)'); ylabel('|X(k)|');
One might wonder, even though the buffered samples contain integral number of waveform cycles, why the frequency spectrum registered a distinct spike at \(10 \; Hz\) when \(N=100\) and not when \(N=128\). This is due to the different frequency resolution, the measure of ability to resolve two different adjacent frequencies
For case 1, the frequency resolution is \(\Delta f = f_s/N = 100/10 = 1 Hz\). This means that the frequency bins are spaced 1 Hz apart and that is why it is able to hit the bull’s eye at 10 Hz peak.
For case 2, the frequency resolution is \(\Delta f = f_s/128 = 100/128 = 0.7813 Hz\). With this frequency resolution, the x-axis of the frequency plot cannot have exact value of 10 Hz. Instead, the nearest adjacent frequency bins are 9.375 Hz and 10.1563 Hz respectively. Therefore, the frequency spectrum cannot represent 10 Hz and the energy of the signal gets leaked to adjacent bins, leading to spectral leakage.
Experiment 2: Effect of time-limited observation
In the previous experiment, the signal wave observed for 1 second duration and that fetched whole 10 cycles in the signal buffer. Now, reduce the observation time to 0.91 second and re-run the same code, results below. In this case, the signal buffer will have 9.1 cycles of the sinewave, which is not a whole number. For case 1, the frequency resolution is 1 Hz and the FFT plot has registered a distinct peak at 10 Hz. Careful investigation of the plot will reveal very low spectral leakage even in case 1 (observe the non-zero amplitude values for the rest of the bins). This is primarily due to the change in the observation interval leading to non-integral number of cycles within the observed window. The spectral leakage in case 2, when N=128, is predominantly due to mismatch in the frequency resolution.
From these two experiments, we can say that
1) The mismatch between the tone of the signal and the chosen frequency resolution (result of sampling frequency and the FFT length) leads to spectral leakage (experiment 1).
2) Time-limiting an observation (at inappropriate times), may lead to spectral leakage (experiment 2).
3) Hence, the spectral leakage from a larger signal component, if present, may significantly overshadow other smaller signals making them difficult to identify or detect.
Rate this article:
Reference
Similar articles
[1] Understanding Fourier Series
[2] Introduction to digital filter design
[3] Design FIR filter to reject unwanted frequencies
[4] FIR or IIR ? Understand the design perspective
[5] How to Interpret FFT results – complex DFT, frequency bins and FFTShift
[6] How to interpret FFT results – obtaining magnitude and phase information
[7] Analytic signal, Hilbert Transform and FFT
[8] FFT and spectral leakage
[9] Moving average filter in Python and Matlab
[10] Window function – figure of merits
[11] Equivalent noise bandwidth of window functions
Books by the author
Wireless Communication Systems in Matlab Second Edition(PDF) Note: There is a rating embedded within this post, please visit this post to rate it. | Digital Modulations using Python (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. | Digital Modulations using Matlab (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. |
Hand-picked Best books on Communication Engineering Best books on Signal Processing |
how can FFT be done with 100 points? does not it need number of points as power of 2?
The claim that FFT only functions with powers of two is untrue.
Refer the original work by Cooley & Tukey “An Algorithm for the Machine Calculation of Complex Fourier Series”. In this work, they indicate that the number of points N=r1.r2 for FFT algorithm can be a highly composite number (need not be a power of 2). However, they say that usage of power of 2 FFT has special advantage when implementing the computation on a binary computer.
How can it be explained, that there is a kind of leakage effect if we take 9 8 7 down to 1 cycle(s) of the signal? Or is it the fourier transform of a rect window, what can be seen (so in spectral dimension it will be a sinc)?
Lets do the experiment. Set observationTime=0.1 (equivalent to 1 cycle, fetches 10 samples). Now you have the choice of setting the DFT length.
Case 1: DFT length set equal to integral multiple of number of samples in the observation window. Result shows, perfect spikes at +/-10Hz frequencies. Same will be true for the cases where the DFT length is set to any multiple of 10
https://uploads.disquscdn.com/images/167cc0396eeca3547bc4bf8780402f7d5f345d28f5c95824e160f61163fa9329.png
Case 2: Set DFT length to some arbitrary number that is not a multiple of samples in the observation window. Examples shown for N=64,128,256. The frequency leaks to all other bins. Examples below
https://uploads.disquscdn.com/images/70c6061cb0aadfc8518a81a3d8193c312583a2b3ffef52af1316f60adbd62b06.png
https://uploads.disquscdn.com/images/18ee735bf0c82ac8e0ce72ccd33e2ce49188f6ce073ba97f1488f8c0ed7e6537.png
https://uploads.disquscdn.com/images/744ae6bbaad6dc9995f966a6d2cc51158546b603ba7bd6a3820723d0d003e4ff.png
If I do this experiment i get the result as seen below (setting N=100 and observationTime=0.1) . https://uploads.disquscdn.com/images/e249e89cca959d3de20e14ce5185ca98cbcfdc6f9e0121c1826de5e667d2bcdc.png
The result as can bee seen in the first picture of yours, i only get if i set N to 10 (Though the title says its N=100, because it’s not a ‘variable’in the title).
Here is the updated code, where I have replaced all the hard-coded values with variables. Now, you should see similar results as shown here
https://uploads.disquscdn.com/images/aeb4cb608c74ae51fbede9d5791b8451cd9e6182f85484049d1182285939e0b0.png
clear all; clc;
Fx=10; %Frequency of the sinusoid
Fs=100; %Sampling Frequency
observationTime = 0.1; %observation time in seconds
t=0:1/Fs:observationTime-1/Fs; %time base
x=sin(2*pi*Fx*t);%sampled sine wave
N=length(x); %DFT length same as signal length
X1 = 1/N*fftshift(fft(x,N));%N-point complex DFT of x
f1=(-N/2:1:N/2-1)*Fs/N; %frequencies on x-axis, Fs/N is the frequency resolution
figure;
subplot(3,1,1);stem(x,’r’)
title(‘Time domain’);xlabel(‘Sample index (n)’);ylabel(‘x[n]’);
subplot(3,1,2);stem(f1,abs(X1)); %magnitudes vs frequencies
xlim([-Fx-20,Fx+20]);title([‘FFT, N=’,num2str(N),’, Delta f=’,num2str(Fs/N)]);xlabel(‘f (Hz)’); ylabel(‘|X(k)|’);
N=128; %DFT length
X2 = 1/N*fftshift(fft(x,N));%N-point complex DFT of x
f2=(-N/2:1:N/2-1)*Fs/N; %frequencies on x-axis, Fs/N is the frequency resolution
subplot(3,1,3);stem(f2,abs(X2)); %magnitudes vs frequencies
xlim([-Fx-20,Fx+20]);title([‘FFT, N=’,num2str(N),’, Delta f=’,num2str(Fs/N)]);xlabel(‘f (Hz)’); ylabel(‘|X(k)|’);
Thanks for the updated code. Is it in this case correct, that only the same length of the signal and FFT provide a leakage free spectrum? Because if I try integral multiples (e.g. N=100 -> frequency res. is 1 Hz), i also get leakage kind spectrum (as mentioned in the first post of mine).
And for example if i set the observation time to 0.6 (6 cycles), and set N=60 i will get a leakage free spectrum.
I’m sorry for asking this, but i’m struggeling right now with this kind of topic.
Thanks for your time.
the two sources of spectral leakage have confused me for a while. This is exactly the article I’m look for. excellent work! thank you!
Thanks for your comment.
Solid discussion…thnk u so much for sharing…….just rock it………
Thanks for your comment
I have a question regarding the way you are calculating FFT here. It’s actually not straight forward to 100% understand every step regarding the FFT calculation, as sources seem to slightly disagree, mainly on scaling factor applied to FFT. From your code:
NFFX=2.^(ceil(log(length(x))/log(2)));
FFTX=fft(x,NFFX);%pad with zeros
NumUniquePts=ceil((NFFX+1)/2);
FFTX=FFTX(1:NumUniquePts);
MY=abs(FFTX);
MY=MY*2;
MY(1)=MY(1)/2;
MY(length(MY))=MY(length(MY))/2;
MY=MY/length(x);
f1=(0:NumUniquePts-1)*Fs/NFFX;
I can see that you multiply by 2 to account for the half FFT that you threw away. Is this correct?
What I do not get is why you then “undo” that 2x scaling factor to the first and last value of vector MY ?
Thanks
As a complement to what you just explained, you might as well have mentioned coherent sampling, which is used to improve FFT visualization a little bit more. Excellent explanation anyway. Thanks
Excellent article, thanks a lot! I don’t understand why in the fig3, the FFT has spectral leakage, given the measured signal is an integral number of cycles.
Even though it has integral number of cycles, the observation time is much shorter. The spectral leakage in this case is due to “limited observation interval#
I kindly do not agree with that. It is necessary just 1 period to entirely determine its DFT, making sure that we do not add spurious points at the edges. Cheers
He is right Seba. Observing the signal for a finite time we effectively apply a square window on it. This leads to the window response overlap with the signal spectrum. Different shape windows can be used depending on need but every one introduces some leakage into the spectrum: https://en.wikipedia.org/wiki/Window_function
Thanks for the nice post Math! Neatly explained.
I think Seba is right. If you see the code supposedly for sine wave with 1 second window (7 exact cycles of the sine wave), you can see that it is infact not the case. For it to be exactly 7 cycles, the line should be,
t=0:1/Fs:(observationTime – 1/Fs).
In the present situation, the code is adding an extra zero at the end (effectively, it is padding a zero) that makes it a non-perfect sinusoid, and that is the reason for leakage.
Good catch. Yeah I missed the -1/Fs part in the code and was oblivious to the fact that 1 sec duration should entirely contain integral number of cycles without any glitches. Thanks for finding that. !!!
Very interesting, cheers !
Thanks for your comment.