View Question
 Question
 Subject: Fast Fourier Transform Category: Science > Math Asked by: ilya-ga List Price: \$20.00 Posted: 29 Oct 2002 22:16 PST Expires: 28 Nov 2002 22:16 PST Question ID: 92858
 ```How to apply the Fast Fourier Transform to a vector of general length N? I konw how it works if N is a power of two and I konw about the mixed-radix algorithm for cases when in N can be reprsented as a product of primes, it's not what I am looking for. I believe there is a way to pad the vector by zeros to the nearest power of two and then do some magic with FFTs and convolutions to get the DFT of the original vector. (And I'm not looking for approximations either, I need an O(N*log(N)) exact algorithm that works for general vector size N even when N is a prime.)```
 Subject: Re: Fast Fourier Transform Answered By: rbnn-ga on 31 Oct 2002 23:39 PST Rated:
 ```(This is kind of long; of course, if you would like any additional clarification please use the "request clarification" button to solicit it.) Thank you for a particularly thought-provoking question! Now, oddly, my dissertation had a chapter on Fourier transforms, particularly over non-abelian groups, and I've written my own implementation of FFT and variants. And indeed, when I first saw your question, I thought "we just pad with zeros". Indeed, I think that if you were to ask the majority of computer scientists how to derive a general FFT from a power-of-two FFT, they would just say "pad with zeroes". So I was a little surprised to see that just padding with zeroes doesn't quite work. In fact, the only way I see to do it involves a couple of extra FFTs . I will provide an outline of an algorithm, then I will describe it in more detail, and then I will show a Matlab implementation. OUTLINE OF ALGORITHM -------------------- To compute the DFT of a vector of length n, we first rewrite the DFT as a convolution. We then compute the convolution by padding both inputs and using a power-of-two FFT to compute that convolution. (I will use DFT and FFT somewhat interchangably) DESCRIPTION OF ALGORITHM ------------------------ Preliminary definitions ---------------------- In the following, n will be a positive integer. All vectors will be over C, the field of complex numbers. w will be a primitive n'th root of unity, for example, exp(-2 pi i/n) . Let x be a vector of length n whose j'th element is x_j. The DFT of x, DFT(x) is the vector of length n whose k'th element is the sum over all 0<=j> x=rand(1,1000); % Compare myfft with the built-in fft >> max(abs(myfft(x)-fft(x))) 2.2440e-011 ----------------- SEARCH STRATEGY I did not find anything useful in an internet search. The standard work on algorithms, Introduction to Algorithms, by Thomas Cormen, Charles Leiserson, Ronald Rivest, and Clifford Stein. MIT Press, 2001. ISBN 0-262-03292-7 contains the lines in reference to the DFT "Without loss of generality, we assume that n is a power of 2, since a given degree-bound can always be raised---we can always add new high-order zero coefficients as necessary". [ibid, p. 833, section 30.2] I rewrote the FFT as a chirp transform; the chirp transform algorithm is given as a starred exercise, 30.2-8, in this text. This does not seem to be the case, and I have sent in a bug report to the authors . I also referenced the text: Discrete-Time Signal Processing, by Alan Oppenheim and Ronald Schafer. Prentice Hall, 1999 . However, there are many annoying notational differences between the electrical engineering viewpoint of FFT's and the computer scientist's . Finally, I would like to mention that Fourier analysis can be carried out over an arbitrary group, even over a non-abelian group. The resulting theory is quite beautiful. An introduction is: Fourier Analysis On Finite Groups and Applications. Audry Terras. Cambridge University Press, 1999. ISBN: 0521457181``` Clarification of Answer by rbnn-ga on 01 Nov 2002 08:53 PST ```ilya-ga: I wanted to add three small comments to my answers. First, several commenters have pointed you in the direction of tailored FFTs for prime length. From the wording of your question (namely that you already knew about the prime-length FFT algorithms, and that you thought there was a way to solve this using power-of-two FFT) I am presuming that want a way to compute general FFT from power-of-2 FFT in n log n time. Obviously there are specialized and more complex FFT algorithms for arbitrary n that compute the FFT directly (see the comments). Still, the method given has applications in situations where prime-FFTs are two complicated, for example, computing an FFT on a hypercube. Second, the function definition line in the matlab program should say "myfft" and not "myfft2" at the top . (The function is in a file called "myfft.m" and so it does actually work there). Third, due to my editing error the sentence beginning "This does not seem to be the case" should occur before the sentence that precedes it in the text.```
 ilya-ga rated this answer: `Thanks for a great answer!`

 ```If you're looking for an exact answer, then I think convolution is probably not the way to go: the exact resampling convolution is probably as expensive as the FFT itself. One good approach is what FFTW does. They break n down into its factors, and create a "plan" with each stage handling one prime. For small prime factors, they have separate "twiddler" functions. If a large prime factor remains, they use Rader's algorithm, which is O(n log n), but with a worse constant factor than the usual Cooley-Tukey approach. Rader's original 1968 paper has a cite page here: http://citeseer.nj.nec.com/context/303893/0 Of course, depending on your licensing constraints, maybe you can just use fftw: http://www.fftw.org/ It certainly isn't going to be easy to do a better job than they do. Have fun!```
 ```Thanks to ilya-ga, for posting the very interesting question, and to raphlevien-ga for the perceptive answer/comment. I found a detailed example of Rader's algorithm for odd primes (independently discovered and generalized by S. Winograd to odd prime powers) in this paper (available on-line in PDF format): [Lesser Known FFT Algorithms, by R. Tolimieri and M. An] http://www.cs.umb.edu/~asi/analysis2000/papers/tolimieri.pdf The Rader-Winograd algorithm uses the fact that for prime p, Z/pZ has multiplicative group (Z/pZ)* which is cyclic of order p-1. A Fourier transform can then be converted into convolution on p-1 points, together with a permutation and reverse permutation of nodes (according to the powers of some choice of generator for (Z/pZ)*), and the convolution can be carried out by means of two Fourier transforms and a point multiplication. I'd be interested to know the target application for ilya-ga's FFT. As I read over the paper above, it reminded me of the recent Agarwal-Kayal-Saxena paper on "Primes in P". Thus I was intrigued by the references at the end of Tolimieri and An's paper, which included some by R. C. Agarwal and C. S. Burrus. But apparently this is just coincidence of names; the AKS paper Agarwal is Manindra Agarwal. Nonetheless it appears that looking at the FFT/convolution computation as a polynomial ring computation is valid, and the parallel with the AKS primality test (which requires numerous multiplications mod x^p - 1 for primes p) still seems striking. regards, mathtalk-ga```
 ```rbnn-ga, thanks again for your answer, that was exactly the information I needed. I kept looking on the Internet after I received your answer and I also found a similar algorithm: http://www.eptools.com/tn/T0001/PT11.HTM raphlevien-ga and mathtalk-ga, thanks for your comments. I know about FFTW and Winograd's algorithm. I was looking for a geral length N FFT algorithm without any assumptions on the primeness or decomposition of N and rbnn-ga's answer was right on the target. Thanks again.```