another FFT amplitude question
Mostrar comentarios más antiguos
I spent a good amount of time reading/trial-error regarding FFT. Im postprocessing a measured waveform and i need to make sure my amplitude is correct from time to frequency and back to time domain but the waveform is of course quite complex, and so im not confident in the amplitude i can produce from various methods and checks.
I decided to trouble-shoot with user generated simple sin/cosines to test the methodology. In short, i have seen suggestions to multiply the FFT(x) by the sampling time (ts), and resulting IFFT(FFT(x)) by fs. However, after reading many sources on the FFT, reviewing my DFT/DTFT math and trial-error i seem to only reproduce the correct amplitude in my ideal case when normalizing the FFT(x) by the number of samples in the original waveform.
Is the method below correct for any signal? If the waveform is already sampled, and i interpolate it in order to set my own fs and # of points... do i need to multiple by fs or simply divide by nfft?
Here is my ideal test case:
N = 4096; % picked # of samples
fs = 409.6; %sampling F (determined by N and ts, since i wanted a 10second waveform and my fundamentals are adhering to nyquist condition)
ts=1/fs; % sampling time
Tstop = N*ts;
t = (0:ts:Tstop-ts); #this has 4096 points and is 10s long
x1 = 1+0.5*sin(2*pi*35*t);
x2 = 1.2*sin(2*pi*5.3*t);
x3 = 0.7*sin(2*pi*150*t);
x = x1+x2+x3;
nfft = 2^(nextpow2(length(x))); %technically not necessary since i force the signal already to a power of 2
NumUniquePts = nfft/2+1; % for single-side
f = (0:NumUniquePts-1)*(fs/nfft);
xfftnfft_dualside = fft(x,nfft)/nfft;
xfftnfft_singleside = xfftnfft_dualside(1:NumUniquePts);
xfftnfft_singleside(2:end-1) = xfftnfft_singleside(2:end-1)*2;
%attempting to use other persons idea to check energy in f and t-domain but here, despite seeing the correct amplitudes the energy isnt correct. I also dont understand the correct scaling here, or if my method above is not correct.
energy_time = sum(x.*conj(x))/nfft
energy_freq_dual = sum(xfftnfft_dualside.*conj(xfftnfft_dualside))*(fs/nfft)
energy_freq_single = sum(xfftnfft_singleside.*conj(xfftnfft_singleside))*(fs/nfft)
figure;
subplot(2,2,1)
plot(abs(xfftnfft_dualside));
legend('fft(x,nfft) dual side');
subplot(2,2,2)
plot(f,abs(xfftnfft_singleside));
legend('fft(xfftnfft_dualside(1:NumUniquePts)) single side');
subplot(2,2,3)
plot(x);
legend('signal vs index');
subplot(2,2,4)
plot(t,x);
legend('signal vs time');
Respuestas (1)
Youssef Khmou
el 12 de Ag. de 2013
0 votos
2 comentarios
Engineer Undertest
el 13 de Ag. de 2013
Youssef Khmou
el 15 de Ag. de 2013
Editada: Youssef Khmou
el 15 de Ag. de 2013
hi Engineer
Thank you for your comment, i would like to see the solution, from other contributors, that takes account of DC components, in the mean while let us try to compute the spectrum of your signal based on built in Mathworks function :
Hs=spectrum.periodogram;
psd(Hs,x,'Fs',fs)
The submission works for some kinds of signals, because there is difference between both of the solutions,
To be continued ...
Categorías
Más información sobre Transforms en Centro de ayuda y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!