0
$\begingroup$

Consider a sine wave having an amplitude of 9 and a frequency of 0.3 Hz. I wish to perform an fft and plot the frequency on the x-axis and the real amplitude (i.e. 9) on the y-axis. I am using the following MATLAB code:

Fs = 10;
data = xlsread('Signal.xlsx');
Y(1)=0;
Y = fft(data);
N = length(Y);
freq = 0:Fs/N:Fs/2-Fs/N;
freq = freq';
amplitude = abs(Y(1:floor(N/2)))/floor(N/2);
plot(freq,amplitude);
xlabel('Frequency [Hz]');
ylabel('Amplitude');

When running the code, an amplitude close to 9 is obtained at a frequency of 0.3 Hz. In reality, I have a sinusoidal signal of a rotating motor shaft. I am measuring the angle using an encoder. The data can be found at:

http://www.dropbox.com/s/pe3g44b48ugauqp/Signal.xlsx?dl=0&m=

I know that the amplitude of the sine wave should be about 0.6 m but when conducting an fft, it is indicating that the amplitude is around 0.4 m. I am sampling at 10Hz and have about 185 000 data points. When trying to reduce the data, the amplitude of the sine wave is approximating the desired one i.e. 0.6 m. As far as I know, the more data points you have, the better will the resolution be when conducting an fft. So why am I getting the opposite? Can someone please help me? Thanks.

$\endgroup$
4
$\begingroup$

Judging from your post, you have over 5 hours of data, your are using to calculate the spectrum. This is a lot of data points to process with FFT. Let's focus on your code and modify it a little bit to generate the $0.3\; \mathtt{Hz}$ sinusoid with an amplitude of $9$:

Fs = 10;
t = 0:0.1:1000;
A = 9;
f0 = 0.3;
data = A*sin(2*pi*f0*t);
N = length(data);    

Y(1)=0;
Y = fft(data);

freq = 0:Fs/N:Fs/2-Fs/N;
freq = freq';
amplitude = abs(Y(1:floor(N/2)))/floor(N/2);

subplot(211)
plot(t, data)
xlabel('Time [s]')
ylabel('Amplitude')
title('Time domain signal')
grid on

subplot(212)
plot(freq, amplitude);
xlabel('Frequency [Hz]');
ylabel('Amplitude');
title('Amplitude spectrum')
grid on

This gives us what you expect to obtain.

enter image description here

Now let's assume that your $0.3\; \mathtt{Hz}$ signal is not always present (sometimes there is only noise present). I suspect this is how your data looks like, especially over so many hours. I will approximate this situation by adding the random noise after half of a sinusoid duration: data(5000:10000)=randn(1,5001); (consider it as someone turned off the machine halfway through the experiment). Obviously you can make it even more intelligent by having signal for 1 second, followed by 1 second of noise, and so on.

How it will affect your amplitude spectrum? Obviously you won't get $9$ as an amplitude because there is only half of the energy present before at $0.3\; \mathtt{Hz}$. Thus you will observe something like:

enter image description here

If you want to detect a true amplitude, then I suggest you to use the Time-Frequency analysis tools. Easiest one is the Short-Time Fourier Transform. You will be able to observe the amplitude of your component at a given time, and extract it later on.


After analysing your signal indeed you can expect the peak with an amplitude of $0.62$ at $0.2\; \mathtt{Hz}$ (not as you claiming $0.3\; \mathtt{Hz}$):

enter image description here

You problem relates to the leakage as you get lot's of frequency bins and likelihood of your frequency falling exactly between get's lower. Indeed, to many is also not a good idea. Below you can observe two spectra for two synthesized sinusoids with same amount of samples as your dataset and amplitude $0.62$. Frequencies are $f=0.2\; \mathtt{Hz}$, and $f=0.2003\; \mathtt{Hz}$. You can observe that for a perfect match, amplitude is estimated correctly, although for $0.2003\; \mathtt{Hz}$ you get $0.41$.

enter image description here

One workaround to deal with that problem in your case, would be to apply windowing to your signal. Although it will decrease amount of leakage, but still it won't remove it absolutely. Below you can observe example for the Hamming window.

enter image description here

In my opinion better solution would be to use less samples. For example if you expect your frequency to be around $0.2\; \mathtt{Hz}$, then you can choose $100$ samples (I know - very crude), which will give you resolution of $0.1\; \mathtt{Hz}$ and component will fall exactly inside. Here is the example:

enter image description here

One last remark, regarding your code. I do not understand what your declaration Y(1)=0 is suppose to do? If it's initialization, then it is unnecessary. If you are trying to remove the DC component, then it's not a good place nor way.


In case of using other window than rectangular one, i.e. Hamming you must normalize it properly. In case of DFT you must divide your output samples by sum of all window samples. In case of 'no windowing' which is in fact rectangular window you get sum of all samples being equal to N. Putting it into code:

win = hamming(N); % rectwin(N), hann(N), ...
data = data .* win;
Y = fft(data)/sum(win);

For more theory please refer to other posts on DSP SE or take a look into this great paper.

$\endgroup$
  • $\begingroup$ Thanks for your reply. Please have a look at my signal available at: dropbox.com/s/pe3g44b48ugauqp/Signal.xlsx?dl=0&m= If plotted, it can be seen that it closely resembles a cosine wave with an amplitude of about 0.62. However, when performing an FFT using the code that I indicated in the previous post, the amplitude of the FFT is about 0.42 and I believe that the signal does not contain a lot of noise. This is why I cannot understand why such thing is happening. I am not sure about the use of STFT since I only have 1 frequency in my signal. What do you think? thanks. $\endgroup$ – user10881 Aug 26 '14 at 12:17
  • $\begingroup$ @user10881: Answer updated. $\endgroup$ – jojek Aug 26 '14 at 20:28
  • $\begingroup$ Thanks for your answer. Very helpful. Yes, I was defining Y(1)=0 for initialization. I have some questions. When applying no window, which if I understood correctly, MATLAB uses the rectangular window by default, I am dividing the magnitude by N/2 to obtain the true amplitude of the sine wave. Can you please tell me how you're obtaining the true amplitude when applying a hamming window? Is it possible to paste the code you're using here to better understand? I was considering using another window the flat top window, since I'm mostly interested in the amplitude. What do you think? $\endgroup$ – user10881 Aug 28 '14 at 11:07
  • $\begingroup$ @user10881: Initialization of Y(1)=0 is unnecessary in that case. Matlab doesn't care about that. Answer updated with info regarding windows. $\endgroup$ – jojek Aug 28 '14 at 11:34
1
$\begingroup$

When doing an FFT of a length on the order of 100k samples, the frequency resolution will be on the order of 0.001%, which means even tiny variations in the signal frequency (from being perfectly integer periodic in the FFT length) or the sampling frequency can spread the energy of a sinusoid out into more than one FFT result bin, thus reducing the magnitude of the peak bin.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.