Explore DFT Resolution, need help with answer interpretation
6 次查看(过去 30 天)
显示 更早的评论
Hello there, I'm relatively new to Matlab and Digital Signal Processing. I'm trying to analyze the frequency resolution of the DFT as a function of the DFT size assuming a rectangular window. Using the a sum of sinusoids sampled at 5 kHz vary the frequency spacing
of the sinusoids to show the frequency resolution for at least 4 different DFT lengths, 256, 512, 1024, and 4096.
I got two versions of code from chatgpt. But I don't know how to interpret them.
1st version
% MATLAB code to analyze the frequency resolution of the DFT using dftmtx and a rectangular window
fs = 5000; % Sampling frequency (5 kHz)
t_duration = 1; % Signal duration in seconds
t = 0:1/fs:t_duration-1/fs; % Time vector
% DFT sizes to analyze
N_dft = [256, 512, 1024, 4096];
% Frequencies of the sum of sinusoids (vary frequency spacing)
f1 = 1000; % Frequency of the first sinusoid (1 kHz)
f_spacing = [1, 5, 10, 20]; % Frequency spacing between the two sinusoids
f_end = f1 + f_spacing; % Frequency of the second sinusoid
% Prepare figure
figure;
hold on;
for N = N_dft
figure; % Create a new figure for each DFT size
hold on;
for spacing = f_spacing
f2 = f1 + spacing; % Second sinusoid frequency
% Generate the sum of two sinusoids with frequencies f1 and f2
x = sin(2*pi*f1*t) + sin(2*pi*f2*t);
% Apply rectangular window (by taking the first N samples)
x_windowed = x(1:N); % Select the first N samples for DFT calculation
% Generate DFT matrix for size N using dftmtx
DFT_matrix = dftmtx(N);
% Manually compute the DFT using the DFT matrix
X = DFT_matrix * x_windowed(:); % Compute DFT of the windowed signal
% Compute the frequency axis for the current DFT
freq_axis = (0:N-1)*(fs/N);
% Plot the magnitude of the DFT
plot(freq_axis, abs(X), 'DisplayName', ['Spacing = ', num2str(spacing), ' Hz']);
end
% Add labels and legend
xlabel('Frequency (Hz)');
ylabel('Magnitude');
title(['DFT Magnitude Spectrum (N = ', num2str(N), ')']);
legend('show');
grid on;
hold off;
end
which gives me these resulted figures.
All the peaks are so similar that I can't distinguish them at all.
2nd version of code is here:
% Parameters
fs = 5000; % Sampling frequency (Hz)
T = 1; % Signal duration (seconds)
t = 0:1/fs:T-1/fs; % Time vector
% Sinusoids: Frequencies and amplitudes
f1 = 400; % Frequency of first sinusoid (Hz)
f2 = 450; % Frequency of second sinusoid (Hz)
f3 = 500; % Frequency of third sinusoid (Hz)
% Create the signal: sum of three sinusoids
x = sin(2*pi*f1*t) + sin(2*pi*f2*t) + sin(2*pi*f3*t);
% Define DFT lengths
DFT_lengths = [256, 512, 1024, 4096];
% Plot frequency resolution for each DFT length using a rectangular window
figure;
for i = 1:length(DFT_lengths)
N = DFT_lengths(i); % DFT size
rectangular_window = rectwin(N)';
% Apply rectangular window (implicitly applied since no windowing function)
% Select only the first N samples for each DFT calculation
x_windowed = x(1:N).*rectangular_window;
% Manually compute the DFT
% Compute the N-point DFT
X = (dftmtx(N)*x_windowed(:)).'; % Compute the N-point DFT
X_mag = abs(X); % Magnitude of DFT
% Frequency vector
f = (0:N-1) * (fs / N);
% Plot the DFT magnitude (first half)
subplot(length(DFT_lengths), 1, i);
plot(f(1:N/2), X_mag(1:N/2));
title(['DFT Magnitude Spectrum with Rectangular Window, N = ', num2str(N)]);
xlabel('Frequency (Hz)');
ylabel('|X(f)|');
grid on;
end
% Adjust plot layout
sgtitle('Frequency Resolution for Different DFT Sizes Using Rectangular Window');
%Part 5
% MATLAB code to analyze zero-padding in the DFT using a rectangular window
fs = 5000; % Sampling frequency (5 kHz)
t_duration = 1; % Signal duration in seconds
t = 0:1/fs:t_duration-1/fs; % Time vector
% Window sizes to analyze
window_sizes = [256, 512];
% Zero-padded DFT sizes to analyze
N_dft = [1024, 2048, 4096];
% Frequencies of the sum of sinusoids (vary frequency spacing)
f1 = 1000; % Frequency of the first sinusoid (1 kHz)
f_spacing = [5, 10]; % Frequency spacing between the two sinusoids
f_end = f1 + f_spacing; % Frequency of the second sinusoid
% Prepare figure
for window_size = window_sizes
figure; % Create a new figure for each window size
hold on;
for N = N_dft
for spacing = f_spacing
f2 = f1 + spacing; % Second sinusoid frequency
% Generate the sum of two sinusoids with frequencies f1 and f2
x = sin(2*pi*f1*t) + sin(2*pi*f2*t);
% Apply rectangular window (by taking the first window_size samples)
x_windowed = x(1:window_size); % Select the first window_size samples
% Zero-pad the signal if DFT size is larger than window size
x_padded = [x_windowed, zeros(1, N - window_size)];
% Generate DFT matrix for size N using dftmtx
DFT_matrix = dftmtx(N);
% Manually compute the DFT using the DFT matrix
X = DFT_matrix * x_padded(:); % Compute DFT of the windowed and zero-padded signal
% Compute the frequency axis for the current DFT
freq_axis = (0:N-1)*(fs/N);
% Plot the magnitude of the DFT
plot(freq_axis, abs(X), 'DisplayName', ['Spacing = ', num2str(spacing), ' Hz, N = ', num2str(N)]);
end
end
% Add labels and legend
xlabel('Frequency (Hz)');
ylabel('Magnitude');
title(['Zero-Padded DFT Magnitude Spectrum (Window Size = ', num2str(window_size), ')']);
legend('show');
grid on;
hold off;
end
with the second version of code I can see the differences among different sample sizes, as more sampling yields more sharp signals, but why is that and what differences it will make?
0 个评论
采纳的回答
Paul
2024-10-7
Hi Zhen,
The first code gives the expected result, you just need to zoom in to see it.
If we look at how freq_axis is defined, we see that the difference between successive frequencies is:
fs/N
fs = 5000; % Sampling frequency (5 kHz)
t_duration = 1; % Signal duration in seconds
t = 0:1/fs:t_duration-1/fs; % Time vector
% DFT sizes to analyze
N_dft = [256, 512, 1024, 4096];
% Frequencies of the sum of sinusoids (vary frequency spacing)
f1 = 1000; % Frequency of the first sinusoid (1 kHz)
f_spacing = [1, 5, 10, 20]; % Frequency spacing between the two sinusoids
f_end = f1 + f_spacing; % Frequency of the second sinusoid
The ratio fs/N is
fs./N_dft
which, I think, explains the plots below (along with the fact that as N increases we'll have points in freq_axis that get closer to f1 and f2).
% Prepare figure
for N = N_dft
figure; % Create a new figure for each DFT size
hold on;
for spacing = f_spacing
f2 = f1 + spacing; % Second sinusoid frequency
% Generate the sum of two sinusoids with frequencies f1 and f2
x = sin(2*pi*f1*t) + sin(2*pi*f2*t);
% Apply rectangular window (by taking the first N samples)
x_windowed = x(1:N); % Select the first N samples for DFT calculation
% Generate DFT matrix for size N using dftmtx
DFT_matrix = dftmtx(N);
% Manually compute the DFT using the DFT matrix
X = DFT_matrix * x_windowed(:); % Compute DFT of the windowed signal
% Compute the frequency axis for the current DFT
freq_axis = (0:N-1)*(fs/N);
% Plot the magnitude of the DFT
plot(freq_axis, abs(X),'-o', 'DisplayName', ['Spacing = ', num2str(spacing), ' Hz']);
end
% Add labels and legend
xlabel('Frequency (Hz)');
ylabel('Magnitude');
title(['DFT Magnitude Spectrum (N = ', num2str(N), ')']);
legend('show');
grid on;
hold off;
xlim([f1-10, f2+10])
end
更多回答(0 个)
另请参阅
产品
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!