Go back to previous page
Forum URL: http://www.cut-the-knot.org/cgi-bin/dcforum/forumctk.cgi
Forum Name: College math
Topic ID: 562
Message ID: 4
#4, RE: spectral analysis question
Posted by mr_homm on Mar-11-06 at 11:17 PM
In response to message #2
What you are seeing is not an artifact of the sampling size, nor is it a peculiarity of MatLab's implementation of the discrete Fourier transform. Rather, it is an inherent feature of the DFT itself, regardless of the particular algorithmic implementation used.

Here are the relevant properties of Fourier (and Fourier-like) transforms in general:

Transforms of periodic functions result in discrete functions (this is fairly obvious, since if a periodic function can be decomposed into a sum of sinusoids at all, it must be decomposed into a sum containing only frequencies that are multiples of the original function's frequency).

Transforms of aperiodic functions result in continuous functions (or actually functions that are nonzero on a dense set, of which continuous functions are a subclass).

Inverse transforms exist, and they have the same mathematical form as the forward transforms, except for some sign or scaling constant differences.

Based on these facts, for any pair of functions that are transforms of each other, if one is periodic, the other is discrete, and if one is aperiodic, the other is continuous, and this works both ways. Therefore, there are 4 distinct versions of the Fourier transform, based on whether the functions are discrete (D), continuous (C), periodic (P) or aperiodic (A):


CA ---Fourier Transform------------------>CA
CP ---Fourier Series--------------------> DA
DA ---Discrete Time Fourier Transform---> CP
DP ---Discrete Fourier Transform--------> DP

Note that the DFT is inherently periodic, so the complete information content of the DFT is displayed within any given single period. It is your choice which period to display. However, since the DFT is a complex Fourier transform, trig functions such as sine and cosine will of course have amplitude peaks at both positive and negative frequencies in their DFT, and because they are real functions, these complex amplitudes will be conjugates of each other.

Therefore, the most natural interval on which to display the DFT is one that is symmetrical around zero. For an N point DFT, this would range from -N/2 up to N/2 - 1, assuming N is even. However, because of the way MatLab is set up, it likes to use vectors that are indexed by nonnegative integers, so the period that is chosen is from 0 to N-1. Hence the information that is located at negative frequencies will be displayed in the upper half of the spectrum. Just imagine that you cut off the upper half of the spectrum, move it over, and reattach it at the left side of the lower half of the spectrum. I think that answers your question about the graph.

There are a couple of other points as well. You probably noticed that the spikes are not perfectly sharp, but frequencies near 50, 75, and 125 also have nonzero amplitudes. This is NOT an inherent property of the DFT. As in a true Fourier transform, these peaks should be delta functions (all value concentrated at a single point), but the DFT in MatLab shows some widening of the peaks. This is an artifact of the calculation, as you can verify by computing the DFT of a pure sine wave by hand (an easy computation if you visualize the terms in the sum as vectors in the complex plane).

Also, as to the selection of the number of sample points, there is no need to make N a power of 2, because there are versions of the Fast Fourier Transform algorithm that can handle other values for N. Basically, the FFT algorithm is so fast because it divides the points into subsets and then recursively calls itself on these smaller sets. It is easy to write an FFT algorithm that divides the points into two subsets, but it is perfectly possible (and has been done) for arbitrary partitions into equal size subsets. There is an nice implementation of this version of the FFT in C language source code, free for noncommercial use. I forget its name, but if you do a Google search for the phrase "fastest FFT in the west" you will find it, as that's the slogan on the website.

It is still true that the fastest values of N are the powers of 2, but basically any N that factors into several small factors will be nearly as fast. What you want to avoid at all costs is a large prime number value for N. For example, 499 is prime, and will take several times as long to compute as 500, which factors as 2*2*5*5*5.

(My job is to coach students through this kind of stuff at the university, so if you have any further questions, please ask, it's no trouble.)

Hope this helps,

--Stuart Anderson