Bartlett Method in MATLAB (with Simulator)

Estimate PSD by segmenting the signal into \(K\) non-overlapping segments, computing periodograms, and averaging:

\[ P_x(f) = \frac{1}{K \cdot M} \sum_{m=0}^{K-1} \left| \sum_{n=0}^{M-1} x_k[n] e^{-j 2 \pi f n} \right|^2 \]

Where \(x_k[n]\) is the m-th segment, \(K\) is number of segments, and \(M\) is the segment length.

Steps to calculate Spectral power density using Bartlett Method

1. 'M' is the length of each segment for the Bartlett method, set to 100 samples.
2. 'K' is the number of segments obtained by dividing the total number of samples N by the segment length 'M'.
3. psd_bartlett_broadband is initialized to store the accumulated periodogram.
4. For each segment k, x_k extracts the k-th segment of the broadband signal.
5. P_k computes the periodogram of the k-th segment using the FFT.
6. The periodograms are accumulated and averaged over all segments.
7. The PSD is plotted in dB/Hz by converting the power values to decibels using 10 * log10.

 

MATLAB Script

clc;
clear;
close all;

% Parameters
fs = 1000; % Sampling frequency
t = 0:1/fs:1-1/fs; % Time vector
N = length(t); % Number of samples

% Generate synthetic broadband ARMA process
arma_order = [2, 2]; % ARMA(p,q) order
a = [1, -0.75, 0.5]; % AR coefficients
b = [1, 0.4, 0.3]; % MA coefficients
%broadband_signal = filter(b, a, randn(size(t)));
% Generate sinusoids with different frequencies
frequencies = [50, 150, 300, 450]; % Frequencies in Hz
amplitudes = [1, 0.8, 0.6, 0.4]; % Amplitudes
broadband_input = zeros(size(t));

for i = 1:length(frequencies)
    broadband_input = broadband_input + amplitudes(i) * sin(2*pi*frequencies(i)*t);
end

% Add noise for realism (optional)
broadband_signal = broadband_input + 0.2*randn(size(t));
broadband_signal = filter(b, a, broadband_signal);

% Generate synthetic narrowband process
f0 = 50; % Center frequency of narrowband process
narrowband_signal = sin(2*pi*f0*t) + 0.5*randn(size(t));

% Parameters for Bartlett method
M = 100; % Length of each segment
K = N / M; % Number of segments

% Initialize the PSD estimate for broadband signal
psd_bartlett_broadband = zeros(1, M);

% Loop over each segment for broadband signal
for k = 1:K
    % Extract the k-th segment
    x_k = broadband_signal((k-1)*M + (1:M));

    % Compute the periodogram of the k-th segment
    P_k = abs(fft(x_k, M)).^2 / M;

    % Accumulate the periodogram
    psd_bartlett_broadband = psd_bartlett_broadband + P_k;
end

% Average the periodograms
psd_bartlett_broadband = psd_bartlett_broadband / K;

% Initialize the PSD estimate for narrowband signal
psd_bartlett_narrowband = zeros(1, M);

% Loop over each segment for narrowband signal
for k = 1:K
    % Extract the k-th segment
    x_k = narrowband_signal((k-1)*M + (1:M));

    % Compute the periodogram of the k-th segment
    P_k = abs(fft(x_k, M)).^2 / M;

    % Accumulate the periodogram
    psd_bartlett_narrowband = psd_bartlett_narrowband + P_k;
end

% Average the periodograms
psd_bartlett_narrowband = psd_bartlett_narrowband / K;

% Frequency axis
f = (0:M-1)*(fs/M);

% Plot the PSD for broadband signal
figure;
plot(f(1:50), 10*log10(psd_bartlett_broadband(1:50)));
title('Power Spectral Density (Bartlett Method) - Broadband Signal');
xlabel('Frequency (Hz)');
ylabel('Power/Frequency (dB/Hz)');

% Plot the PSD for narrowband signal
figure;
plot(f(1:50), 10*log10(psd_bartlett_narrowband(1:50)));
title('Power Spectral Density (Bartlett Method) - Narrowband Signal');
xlabel('Frequency (Hz)');
ylabel('Power/Frequency (dB/Hz)');
web('https://www.salimwireless.com/search?q=spectral%20estimation', '-browser');

 

Output

Bartlett Method PSD for Narrowband Signal

 

Bartlett Method PSD for Broadband Signal

 

Further Reading

1. Periodogram in MATLAB
2. Welch Method for Spectral Estimation in MATLAB
3. Correlogram in MATLAB
4. Power Spectral Density Calculation Using Simple FFT in MATLAB
5. Spectral Estimation Methods - Periodogram, Correlogram, Welch, Bartlett ... (Theory)

---

Online Interactive Simulator for PSD Estimation using Bartlett Method

Parameters

Base Signal

Select an input Sine Wave Cosine Wave Rectangular Pulse Triangular Pulse

Sampling Frequency (Hz) Length of the Signal (Samples) Input Signal Frequency (Hz) Window Type: Rectangular Window Hann Window Hamming Window Window Size Pulse Width (s)

Input Signal

AWGN Noise

SNR (in dB):

Operations

Select an operation None Addition Subtraction Multiplication Convolution

Add another Signal: Select an input Sine Wave Cosine Wave Rectangular Pulse Triangular Pulse

Input Frequency (Hz):

Output

Upload CSV and See Output

Sampling Frequency (Hz):

Input Signal

Output

Access full simulation →
Other simulations →

Spectral Density Estimation Methods

Estimating Spectral Density: The Windowed Periodogram

To accurately determine the Power Spectral Density (PSD) of a discrete signal, engineers frequently utilize the windowed periodogram. This approach applies a specific weighting function to raw data to mitigate "spectral leakage"—a phenomenon where energy from a primary frequency "seeps" into adjacent bins, distorting the analysis.

Core Periodogram Principles

The most basic form of PSD estimation is the traditional periodogram, derived directly from the Discrete-Time Fourier Transform (DTFT). In its raw state, the periodogram is defined as:

\[ P_x(f) = \frac{1}{N} \left| \sum_{n=0}^{N-1} x[n] e^{-j 2 \pi f n} \right|^2 \]

Where \(x[n]\) represents the discrete samples and \(N\) is the total sample size. Without modification, this method is highly susceptible to leakage errors due to the abrupt truncation of the signal.

Improving Precision via Windowing

By multiplying the signal by a window function \(w[n]\) before processing, we taper the edges of the data, leading to a more accurate spectrum:

\[ P_x(f) = \frac{1}{N \cdot U} \left| \sum_{n=0}^{N-1} x[n] w[n] e^{-j 2 \pi f n} \right|^2 \]

To keep the total power consistent, we use a normalization constant \(U\):

\[ U = \frac{1}{N} \sum_{n=0}^{N-1} |w[n]|^2 \]

Popular Windowing Functions

• Rectangular: A simple truncation; treats all samples equally but offers minimal leakage protection.
• Hamming: Uses a cosine-based curve (\(0.54 - 0.46 \cos\)) to suppress the "sidelobes" of the frequency response.
• Hann: Similar to Hamming but provides a smoother fade-out to zero at the boundaries.
• Blackman: Adds extra terms to the sequence to further minimize sidelobes at the cost of a wider central peak.

Alternative PSD Estimation Techniques

1. The Correlogram Method

Based on the Wiener-Khinchin theorem, this technique computes the Fourier transform of the signal’s estimated autocorrelation sequence \(R_x[k]\):

\[ P_x(f) = \sum_{k=-(N-1)}^{N-1} R_x[k] e^{-j 2 \pi f k} \]

To ensure the PSD never yields negative values, a biased estimate (dividing by \(N\)) is standard practice.

2. Bartlett’s Method (Averaging)

Bartlett’s approach reduces statistical noise (variance) by splitting the signal into \(M\) non-overlapping segments and averaging their periodograms:

\[ P_x(f) = \frac{1}{M \cdot N} \sum_{m=0}^{M-1} \left| \sum_{n=0}^{N-1} x_m[n] e^{-j 2 \pi f n} \right|^2 \]

Trade-off: While it reduces variance by a factor of \(M\), it lowers frequency resolution because each segment is shorter than the original data.

3. Blackman-Tukey Method

This method applies a window to the autocorrelation sequence rather than the raw signal. By smoothing the correlation lags before the Fourier transform, it produces a much cleaner estimate.

\[ P_x(f) = \sum_{k=-K}^{K} R_x[k] w[k] e^{-j 2 \pi f k} \]

4. Welch’s Method: The Industry Standard

An evolution of Bartlett’s method, Welch’s technique allows segments to overlap (usually by 50%) and applies a window function to each segment before processing.

\[ P_x(f) = \frac{1}{K \cdot L \cdot U} \sum_{k=0}^{K-1} \left| \sum_{n=0}^{L-1} x_k[n] w[n] e^{-j 2 \pi f n} \right|^2 \]

Welch’s method is widely considered the best balance between noise reduction and spectral leakage prevention.

Return to Spectral Estimation Methods Main Page →