How do I normalize the output of nufft?
11 ビュー (過去 30 日間)
古いコメントを表示
I am trying to figure out what normilization of the output of nufft(x,t) satisfied Parseval's theorem and am unable to figure it out. I tested uniformly sampled data and can see that nufft clearly is not normalized the same way as fft, but am unable to figure out correct normalization.
To illustrate my confusion let me first define a signal that I can calculate the energy of analytically.
%% define a signal
dt = 1/1000; % [Hz]
T = 3; % [s]
t = 0:dt:T; % [s]
f1 = 123; % [Hz]
A1 = 1/sqrt(2); % [-]
f2 = 40; % [Hz]
A2 = 1/sqrt(2); % [-]
sig = sqrt(2)*(A1.*cos(t*2*pi*f1) +A2.*cos(t*2*pi*f2));
L = length(sig);
% True energy in the time domain
E_true = A1.^2 + A2.^2;
Now, lets compute the non-uniform FFT of this uniform data useing nufft
%% Do the FFT
% Take the non-uniform FFT
NUFT = nufft(sig,t);
% Take the standard FFT
UFT = fft(sig);
Once I have computer the approximate Fourier Transform, let's calculate the power spectrum and normalize it the way FFT would require.
% calculate power spectrum
P_NU = abs(NUFT(round(1:L/2+1))).^2;
P_U = abs(UFT(round(1:L/2+1))).^2;
% create the freqs in freq domain
f = (0:L/2)/(L*dt);
% normalize both as if it were FFT
P_NU = P_NU*2/(L^2);
P_U = P_U*2/(L^2);
Now, lets look at Parsevals theroem
%% check Parseval's theorem
% Parseval's theorem is just energy conservation
% Energy in the time domain = energy in the freq domain
% The analytical Energy
disp(['The true energy = ',num2str(E_true)])
% Energy in the time signal
E_t = trapz(t,sig.^2)/(dt*(L-1));
disp(['Energy in time domain = ',num2str(E_t)])
% Energy in the freq domain
E_f_U = trapz(f,P_U)*(dt*(L-1));
disp(['Energy in freq domain (fft)= ',num2str(E_f_U)])
% Energy in the freq domain
E_f_NU = trapz(f,P_NU)*(dt*(L-1));
disp(['Energy in freq domain (nufft)= ',num2str(E_f_NU)])
disp(' ')
Running this code yeilds:
The true energy = 1
Energy in time domain = 1
Energy in freq domain (fft)= 1.0007
Energy in freq domain (nufft)= 0.0006591
How am I supposed to normalize the output of nufft?
2 件のコメント
Chris Turnes
2021 年 11 月 4 日
編集済み: Chris Turnes
2021 年 11 月 4 日
I suspect that the answer I posted on the question linked below is essentially the answer to your question:
The short version is that the default arguments for the sample/frequency points of nufft are 0:(N-1) and (0:(N-1))/N, respectively. If your time vector has a different scale than "sample units", then you have to manually scale your time units (as you've discovered).
回答 (0 件)
参考
カテゴリ
Help Center および File Exchange で Fourier Analysis and Filtering についてさらに検索
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!