Simple Matlab FFT problem

Discussion in 'Programmer's Corner' started by bumclouds, Mar 13, 2010.

1. bumclouds Thread Starter Active Member

May 18, 2008
82
0
Hey folks,

I've got a Matlab FFT problem. Basically all I'm trying to do is generate a sinusoidal wave which is 20Hz, and generate the FFT of it and display that in a plot.

I've done this:

Code ( (Unknown Language)):
1. clc;
2. clear;
3.
4. % Generate Sinusoidal of 20Hz
5. fo = 20;   %frequency of the sine wave
6. Fs = 100; %sampling rate
7. Ts = 1/Fs; %sampling time interval
8. t = 0:Ts:1-Ts; %sampling period  (samples 1-100, from 0sec to 0.99sec)
9. y = sin(2*pi*fo*t); %the sine curve
10.
11. N=length(y); %get the number of points
12. k=0:N-1;     %create a vector from 0 to N-1
13. T=N/Fs;      %get the frequency interval
14. freq=k/T;    %create the frequency range
15. X=fft(y,N);
16.
17.
18. figure(1)
19. plot(t,y);
20.
21. figure(2)
22. plot(freq,abs(X))
and the results I'm getting are these:

Figure 1

Figure 2

What on earth is that thingy at 80Hz? And when I change fo to 30Hz, or 40Hz, the spike on the left in Figure 1 moves up to 30Hz and 40Hz respectively, no worries, but that spike on the right moves closer in.. What on earth is going on? I thought these frequency spikes were meant to be symmetrically distributed around the 0Hz point on the frequency axis.

Thanks guys, you rock
Andrew

Last edited: Mar 13, 2010
2. t_n_k AAC Fanatic!

Mar 6, 2009
5,448
783
Looks like the image frequency equidistant from the Nyquist frequency. See ideas on sampling theory.

You have a sampling frequency fs of 100Hz. So the Nyquist frequency will be fs/2=50Hz.

The raw FFT gives two components lying equidistant from the Nyquist frequency.

It's 30 Hz from 20Hz to 50Hz. So the image lies 30Hz above 50Hz - i.e. 80Hz

Normally we are only interested in frequencies below the Nyquist frequency and take steps to prevent input signal frequencies above the Nyquist value 'aliasing' back into the lower region.

3. guitarguy12387 Active Member

Apr 10, 2008
359
12
Thats just how Matlab plots the 'negative' frequencies because the DFT is periodic with the sampling frequency mapping to pi. There's no aliasing.

Plot using fftshift(). That should take care of it.

4. Nanophotonics Active Member

Apr 2, 2009
365
3
Yes, there's no aliasing, the sampling frequency is high enough to avoid aliasing.

5. bumclouds Thread Starter Active Member

May 18, 2008
82
0
What does aliasing actually look like? It's not that extra spike at 80Hz is it?

FFTSHIFT() is not bringing my spikes back to be about the y-axis (0Hz), it's just leaving them exactly where they are >.< .. the function description says this also:

So the way I'm reading it, is that this function moves the "negative frequency" components into the +ve x-axis..

6. t_n_k AAC Fanatic!

Mar 6, 2009
5,448
783
Your 80 Hz term isn't an alias - it's an artifact of the FFT algorithm.

Aliasing occurs when the signal frequency to be analysed has frequency components which exceed the Nyquist frequency.

If you were to input a single frequency sinusoid at say 95Hz (rather than 20Hz) using your Matlab analysis you would get an alias component showing at 5Hz.

In FFT any input signal component below the Nyquist frequency will show an image frequency folded symmetrically about the Nyquist (half the sampling) frequency - that's just an artifact of the FFT algorithm. So you keep in mind that point that you are only really interested in results below the Nyquist frequency.

Image frequencies aren't a problem in the analysis since you ignore terms above the Nyquist value - BUT aliased frequencies are a problem because the observer "sees" frequency components in the lower half spectrum which aren't real values - hence the term 'alias'. So to avoid aliasing occurring you put the input signal to be analysed through an anti-aliasing filter to limit the signal bandwidth to a range less than the Nyquist frequency.

7. bumclouds Thread Starter Active Member

May 18, 2008
82
0
t_n_k, thanks mate!!

I've modified my code a little, it's listed below. The new code now produces these images (with sampling rate of 100Hz):

for 20Hz;

for 40Hz;

for 95Hz;

the code is:

Code ( (Unknown Language)):
1. clc;
2. clear;
3.
4. % Generate Sinusoidal of 20Hz
5. fo = 20;   %frequency of the sine wave
6. Fs = 100; %sampling rate
7. Ts = 1/Fs; %sampling time interval
8. t = 0:Ts:1-Ts; %sampling period  (samples 1-100, from 0sec to 0.99sec)
9. y = sin(2*pi*fo*t); %the sine curve
10.
11. N=length(y); %get the number of points
12. k=0:N-1;     %create a vector from 0 to N-1
13. T=N/Fs;      %get the frequency interval
14. freq = 0:Fs-T; %create the frequency range vector
15. X=(fft(y,N))\N; % normalize the data with / N
16.
17.
18. figure(1)
19. plot(t,y);
20.
21. figure(2)
22. plot(freq,abs(X)), axis([0 100 -2 2]);

Does this look right to you guys?

8. macellan New Member

Dec 12, 2009
2
0
Hi, i am trying to see alising effect on a sample data that i got from Labjack by using its own programme.

Your programme seems true. But if you try to find the signal components of the original signal it gives "wrong" result or the previos figure you showed, the double peak one.

Probably there must be another variable in this calculation : Scan rate!. Because i am sampling the data with 4096 Hz. But i am taking 2048,1024,512.. number of data. I chose 1024. Matlab "fft" function doesnot indicate this.
If i find how to implement it to the code , i'll write it here also. Bye...