The windowed periodogram is a widely used technique for estimating the Power Spectral Density (PSD) of a signal. It enhances the classical periodogram by mitigating spectral leakage through the application of a windowing function. This technique is essential in signal processing for accurate frequency-domain analysis.
Power Spectral Density (PSD)
The PSD characterizes how the power of a signal is distributed across different frequency components. For a discrete-time signal, the PSD is defined as the Fourier Transform of the signal’s autocorrelation function:
Sx(f) = FT{Rx(Ï„)}
Here, Rx(Ï„)}is the autocorrelation function.
FT : Fourier Transform
Classical Periodogram
The periodogram is a non-parametric PSD estimation method based on the Discrete Fourier Transform (DFT):
Px(f) = \(\frac{1}{N}\) X(f)2
Here:
X(f): DFT of the signal x(n)
N: Signal length
However, the classical periodogram suffers from spectral leakage due to abrupt truncation of the signal.
Windowing to Mitigate Spectral Leakage
Spectral leakage can be minimized by applying a window function to the signal before computing the DFT. The resulting PSD estimate is called the windowed periodogram:
Pw(f) = \(\frac{1}{NW}\) Xw(f)2
Here:
w(n): Window function
W: Window normalization factor
Common Window Functions
Rectangular Window: Equivalent to the classical periodogram.
w[n]=1, 0≤n≤N−1
w[n]=0, otherwise
Where, N is the window length
Hamming Window: Reduces sidelobe amplitudes, improving frequency resolution.
w[n]=0.5(1−cos(\(\frac{\ 2\pi n}{N - 1}\ \))), 0≤n≤N−1
Where, N is the window length
Hanning Window: Similar to Hamming but with less sidelobe attenuation.
w[n]=0.54 – 0.46cos(\(\frac{\ 2\pi n}{N - 1}\ \)), 0≤n≤N−1
Where, N is the window length
Blackman Window: Offers even greater sidelobe suppression but at the cost of wider main lobes.
w[n]=0.42 – 0.5(cos(\(\frac{\ 2\pi n}{N - 1}\ \)) + 0.08(cos(\(\frac{\ 4\pi n}{N - 1}\ \)), 0≤n≤N−1
Where, N is the window length
Implementation Steps
Segment the Signal: Divide the signal into overlapping or non-overlapping segments of length N.
Apply a Window Function: Multiply each segment by a window function w(n).
Compute the DFT: Calculate the DFT of the windowed segments.
Average the Periodograms: For overlapping segments, average the periodograms to reduce variance.
Properties of the Windowed Periodogram
Bias: Windowing introduces bias in the PSD estimate as the window modifies the signal spectrum.
Variance: Averaging periodograms (Welch method) reduces variance but decreases frequency resolution.
Trade-Off: The choice of window affects the trade-off between spectral resolution and leakage suppression.
MATLAB Code
clc;
clear;
close all;
fs = 48000;
t = 0:1/fs:0.02;
f_ping = 12000;
% Base sine wave
sine_wave = sin(2*pi*f_ping*t)';
% Apply windows
w_rect = ones(size(sine_wave));
w_hann = hann(length(sine_wave));
w_hamming = hamming(length(sine_wave));
w_blackman = blackman(length(sine_wave));
% Windowed signals
s_rect = sine_wave .* w_rect;
s_hann = sine_wave .* w_hann;
s_hamming = sine_wave .* w_hamming;
s_blackman = sine_wave .* w_blackman;
% FFT
Nfft = 4096;
f = fs*(0:Nfft/2-1)/Nfft;
% Function to compute and normalize spectrum
get_norm_fft = @(sig) abs(fft(sig, Nfft))/max(abs(fft(sig, Nfft)));
S_rect = get_norm_fft(s_rect);
S_hann = get_norm_fft(s_hann);
S_hamming = get_norm_fft(s_hamming);
S_blackman = get_norm_fft(s_blackman);
% Mainlobe power (±2 bins around peak)
mainlobe_bins = 2;
% Function to compute power ratio
compute_power_ratio = @(S) ...
deal( ...
sum(S.^2), ... % Total power
max(1, find(S == max(S), 1)), ... % Peak bin
@(peak_bin) sum(S(max(1,peak_bin-mainlobe_bins):min(Nfft,peak_bin+mainlobe_bins)).^2), ...
@(total, main) 10*log10((total-main)/main) ... % dB sidelobe/mainlobe ratio
);
% Calculate ratios
[total_r, peak_r, get_main_r, get_slr_r] = compute_power_ratio(S_rect);
main_r = get_main_r(peak_r); slr_r = get_slr_r(total_r, main_r);
[total_h, peak_h, get_main_h, get_slr_h] = compute_power_ratio(S_hann);
main_h = get_main_h(peak_h); slr_h = get_slr_h(total_h, main_h);
[total_ham, peak_ham, get_main_ham, get_slr_ham] = compute_power_ratio(S_hamming);
main_ham = get_main_ham(peak_ham); slr_ham = get_slr_ham(total_ham, main_ham);
[total_b, peak_b, get_main_b, get_slr_b] = compute_power_ratio(S_blackman);
main_b = get_main_b(peak_b); slr_b = get_slr_b(total_b, main_b);
% Display Results
fprintf('Window | Mainlobe Power | Sidelobe Power | Sidelobe/Main (dB)\n');
fprintf('------------|----------------|----------------|--------------------\n');
fprintf('Rectangular | %14.4f | %14.4f | %18.2f\n', main_r, total_r - main_r, slr_r);
fprintf('Hann | %14.4f | %14.4f | %18.2f\n', main_h, total_h - main_h, slr_h);
fprintf('Hamming | %14.4f | %14.4f | %18.2f\n', main_ham, total_ham - main_ham, slr_ham);
fprintf('Blackman | %14.4f | %14.4f | %18.2f\n', main_b, total_b - main_b, slr_b);
% Plot
figure;
plot(f, 20*log10(S_rect(1:Nfft/2)), 'k'); hold on;
plot(f, 20*log10(S_hann(1:Nfft/2)), 'r');
plot(f, 20*log10(S_hamming(1:Nfft/2)), 'g');
plot(f, 20*log10(S_blackman(1:Nfft/2)), 'b');
legend('Rectangular','Hann','Hamming','Blackman');
xlim([f_ping-3000 f_ping+3000]); ylim([-100 5]);
xlabel('Frequency (Hz)'); ylabel('Magnitude (dB)');
title('Windowing Effects on Spectrum');
grid on;Output
Window | Mainlobe Power | Sidelobe Power | Sidelobe/Main (dB)
------------|----------------|----------------|--------------------
Rectangular | 3.5771 | 4.9562 | 1.42
Hann | 4.3630 | 8.4370 | 2.86
Hamming | 4.2367 | 7.3928 | 2.42
Blackman | 4.4940 | 10.2410 | 3.58
Applications
Signal Processing: Analyzing frequency content of time-varying signals.
Communications: Evaluating spectrum occupancy in wireless systems.
Bioinformatics: Investigating periodicities in biological signals (e.g., EEG, ECG).
Seismology: Characterizing seismic wave frequencies.
Further Reading