another FFT amplitude question
Show older comments
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');
Answers (1)
Youssef Khmou
on 12 Aug 2013
0 votes
2 Comments
Engineer Undertest
on 13 Aug 2013
Youssef Khmou
on 15 Aug 2013
Edited: Youssef Khmou
on 15 Aug 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 ...
Categories
Find more on Transforms in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!