Why is my hamming window and rectangular window values the same for sigma_omega when I make L = 21?
2 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
I have this code that produces the same value for both the hamming window and the rectangular window when I try to solve for sigma_omega for given L (L = 21, 51, and 101). For the sake of the code I will only be testing for L = 21.
Here is the equation for sigma_omega:
Here is the code I have:
% Set window length L = 21
L = 21;
% Generate rectangular and Hamming windows
rect_window = window(@rectwin, L);
hamming_window = window(@hamming, L);
% Calculate the frequency spread for each window
sigma_omega_rect = calc_freq_spread(rect_window);
sigma_omega_hamming = calc_freq_spread(hamming_window);
% Display the frequency spread for rectangular and Hamming windows for L = 21
fprintf('L = 21: sigma_omega (Rectangular) = %.4f\n', sigma_omega_rect);
fprintf('L = 21: sigma_omega (Hamming) = %.4f\n\n', sigma_omega_hamming);
% Function to calculate frequency spread
function sigma_omega = calc_freq_spread(w)
L = length(w);
% Calculate DTFT (using FFT as an approximation)
N = 1024; % Length of the FFT to improve frequency resolution
W = fft(w, N); % Compute the FFT of the window
% Frequency axis (normalized frequency ω from -π to π)
omega_hat = linspace(-pi, pi, N);
% Magnitude squared of DTFT
W_mag_sq = abs(W).^2;
% Numerator: Sum of ω^2 * |W(e^jω)|^2 (approximating the integral)
numerator = sum((omega_hat.^2) .* W_mag_sq);
% Denominator: Sum of |W(e^jω)|^2 (approximating the integral)
denominator = sum(W_mag_sq);
% Calculate the frequency spread (σ_ω)
sigma_omega = sqrt(sum(numerator) / sum(denominator)); % Summing the whole thing to a scalar
end
My question is what is wrong with my code for producing the same output? The same exact output for both rectangular and hamming window is also the same 58.0983 when I change L to either 51 or 101. Please advise.
Thank you.
0 commentaires
Réponse acceptée
Matt J
le 6 Oct 2024
Modifié(e) : Matt J
le 6 Oct 2024
% Set window length L = 21
L = 21;
% Generate rectangular and Hamming windows
rect_window = window(@rectwin, L);
hamming_window = window(@hamming, L);
% Calculate the frequency spread for each window
sigma_omega_rect = calc_freq_spread(rect_window);
sigma_omega_hamming = calc_freq_spread(hamming_window);
% Display the frequency spread for rectangular and Hamming windows for L = 21
fprintf('L = 21: sigma_omega (Rectangular) = %.4f\n', sigma_omega_rect);
fprintf('L = 21: sigma_omega (Hamming) = %.4f\n\n', sigma_omega_hamming);
% Function to calculate frequency spread
function sigma_omega = calc_freq_spread(w)
% Calculate DTFT (using FFT as an approximation)
N = 1024; % Length of the FFT to improve frequency resolution
w=w(:).'; %<---Matt J
omega_hat = ((0:N-1)-ceil((N-1)/2))*2/N*pi; %<---Matt J
W = fftshift(fft(w, N)); % <--- Matt J
% Magnitude squared of DTFT
W_mag_sq = abs(W).^2;
% Numerator: Sum of ω^2 * |W(e^jω)|^2 (approximating the integral)
numerator = sum((omega_hat.^2) .* W_mag_sq);
% Denominator: Sum of |W(e^jω)|^2 (approximating the integral)
denominator = sum(W_mag_sq);
% Calculate the frequency spread (σ_ω)
sigma_omega = sqrt(sum(numerator) / sum(denominator)); % Summing the whole thing to a scalar
end
Plus de réponses (1)
Paul
le 6 Oct 2024
Modifié(e) : Paul
le 6 Oct 2024
Another approach is to integrate the DTFT directly. The DTFT can be evaluted with freqz
% Set window length L = 21
L = 21;
% Generate rectangular and Hamming windows
rect_window = window(@rectwin, L);
hamming_window = window(@hamming, L);
dtft = @(x,w) freqz(x,1,w);
sigma_w = @(x) sqrt(integral(@(w) w.^2.*abs(dtft(x,w)).^2,-pi,pi)/integral(@(w) abs(dtft(x,w)).^2,-pi,pi));
sigma_w(rect_window)
sigma_w(hamming_window)
Taking advantage of Parseval for the denominator:
sigma_w = @(x) sqrt(integral(@(w) w.^2.*abs(dtft(x,w)).^2,-pi,pi)/(2*pi*sum(abs(x).^2)));
sigma_w(rect_window)
sigma_w(hamming_window)
0 commentaires
Voir également
Catégories
En savoir plus sur Hamming dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!