Comparison of Spectral Estimation Methods in MATLAB: Periodogram, Welch, Bartlett, and More

MATLAB Code 

% MATLAB script to compare different Spectral Estimation Methods using Built-in Functions

% Sample signal: Sine wave with noise
fs = 1000; % Sampling frequency (Hz)
T = 1; % Signal duration (seconds)
t = 0:1/fs:T-1/fs; % Time vector
f_signal = 50; % Signal frequency (Hz)
signal = sin(2*pi*f_signal*t) + 0.5*randn(size(t)); % Sine wave + noise

% Method 1: Periodogram (Using built-in function)
[pxx_periodogram, f_periodogram] = periodogram(signal, [], [], fs);

% Method 2: Bartlett Method (Using built-in pwelch with default settings)
[pxx_bartlett, f_bartlett] = pwelch(signal, [], [], [], fs);

% Method 3: Welch Method (Using built-in pwelch with segment overlap)
segmentLength = 256; % Length of each segment
overlap = 128; % Overlap between segments
[pxx_welch, f_welch] = pwelch(signal, segmentLength, overlap, [], fs);

% Method 4: Blackman-Tukey Method (Using built-in cpsd function with autocorrelation)
% Compute the Cross Power Spectral Density (CPSD) using autocorrelation
[pxx_bartlett_tukey, f_bartlett_tukey] = cpsd(signal, signal, segmentLength, overlap, [], fs);

% Method 5: Correlogram Method (Using xcorr and FFT)
% Compute the full autocorrelation of the signal (unbiased estimator)
autocorr_full = xcorr(signal, 'unbiased'); % Unbiased autocorrelation (lags from -N to +N)
autocorr_full = autocorr_full(length(signal):end); % Keep only positive lags

% Apply a window function (e.g., Hamming) to the autocorrelation
windowed_autocorr_full = autocorr_full .* hamming(length(autocorr_full))';

% Perform FFT of the windowed autocorrelation to get the PSD
N_full = length(windowed_autocorr_full);
fft_autocorr_full = fft(windowed_autocorr_full, N_full);
pxx_correlogram = abs(fft_autocorr_full(1:N_full/2+1)).^2 / fs;

f_correlogram = (0:N_full/2) * (fs/N_full); % Frequency axis for Correlogram

% Plot all the methods in a single figure
figure;

subplot(3,2,1);
plot(f_periodogram, 10*log10(pxx_periodogram));
title('Periodogram');
xlabel('Frequency (Hz)');
ylabel('Power/Frequency (dB/Hz)');

subplot(3,2,2);
plot(f_bartlett, 10*log10(pxx_bartlett));
title('Bartlett Method');
xlabel('Frequency (Hz)');
ylabel('Power/Frequency (dB/Hz)');

subplot(3,2,3);
plot(f_welch, 10*log10(pxx_welch));
title('Welch Method');
xlabel('Frequency (Hz)');
ylabel('Power/Frequency (dB/Hz)');

subplot(3,2,4);
plot(f_bartlett_tukey, 10*log10(pxx_bartlett_tukey));
title('Blackman-Tukey Method');
xlabel('Frequency (Hz)');
ylabel('Power/Frequency (dB/Hz)');

subplot(3,2,5);
plot(f_correlogram, 10*log10(pxx_correlogram));
title('Correlogram Method');
xlabel('Frequency (Hz)');
ylabel('Power/Frequency (dB/Hz)');

% Adjust layout for the figure
sgtitle('Comparison of Spectral Estimation Methods');

 

Output

 

Further Reading

1.