![]() | This article is rated C-class on Wikipedia's content assessment scale. It is of interest to the following WikiProjects: | |||||||||||||||||||||||||||||||||||||||||||
|
Weighting is not trivial: beautiest way to do this is by dividing each f by in each step. Inverse transform is done by simply replacing the negative exponents with a positive one. --152.66.224.133 14:52, 20 May 2004 (UTC)lorro, 2004.05.20.
There's not.. um.. a more simple way of describing this, is there? For "laymen"? I understand Fourier transforms fine, have used FFTs plenty of times and used them to speed up computations of things like the Hilbert transform, but I really don't understand this article. Basically it sounds to me like they have created this formula by which you take, for instance, a 64-point DFT and break it down into 2 32-point DFTs, then you break those down into 4 16-point DFTs, and so on, until you are left with just 32 two-point DFTs to calculate? And then you calculate them and work your way back up, reassembling the results? Are there any graphical analogies? What is this "butterfly" they mention? - Omegatron 17:49, 10 July 2005 (UTC)
It seems that it might actually be *much* clearer if the sigma notation were turned into (1st, ellipses, last) format. I'll do this in a few days myself if there are no objections. 68.163.184.132 03:39, 2 October 2006 (UTC)
I think an 8 point picture of the butterfly network would help the reader. Does anyone have one that is fair use? Jdufresne 01:01, 5 March 2006 (UTC)
Here's a simple implementation in Python based on the description in the article (which I thought was quite good). I spent a lot of time looking for a quick little implementation like this for python but all I could find was heavy mathematical libraries with lots of baggage.
from cmath import *
def fft(x):
N=len(x)
if N==1: return x
even=fft([x[k] for k in range(0,N,2)])
odd= fft([x[k] for k in range(1,N,2)])
M=N/2
l=[ even[k] + exp(-2j*pi*k/N)*odd[k] for k in range(M) ]
r=[ even[k] - exp(-2j*pi*k/N)*odd[k] for k in range(M) ]
return l+r
—The preceding unsigned comment was added by 24.77.42.65 (talk) 03:13, 15 March 2007 (UTC).
Another interesting implementation is in Matlab. Although Matlab has it own fft function, which can perform the Discrete-time Fourier transform of arrays of any size, a recursive implementation in Matlab for array of size 2^n, n as integer (Cooley–Tukey FFT algorithm), follows.
function y = myfft(x)
N = length(x);
if N == 1
y = x;
else
par = myfft(x(1:2:N));
impar = myfft(x(2:2:N));
M = N / 2;
W = exp(-2i*pi*(0:M-1)/N);
l = par(1:M) + W.*impar(1:M);
r = par(1:M) - W.*impar(1:M);
y = [l,r];
end
end
An iterative implementation in Matlab, that also returns the intermediate steps in the matrix A, is shown below.
function [y, A] = myfft_ite(x)
n = length(x);
x = x(bitrevorder(1:n));
q = round(log(n)/log(2));
A = zeros(q,length(x));
for j = 1:q
m = 2^(j-1);
d = exp(-1i*pi*(0:m-1)/m);
for k = 1:2^(q-j)
s = (k-1)*2*m+1; % start-index
e = k*2*m; % end-index
r = s + (e-s+1)/2; % middle-index
y_top = x(s:(r-1));
y_bottom = x(r:e);
z = d .* y_bottom;
y = [y_top + z, y_top - z];
x(s:e) = y;
end
A(j,:) = x;
end
end
— Preceding unsigned comment added by FHZ-1988 (talk • contribs) 05:53, 8 June 2019 (UTC)
It seems silly to use e-sub-m to denote the outputs of the sub-transform performed on the even elements of the original sequence, given that e stands for the base of the natural logarithm, in which role it appears in the same term where e-sub-m is used!
If your programming language allowed you to have a scalar e and an unrelated array e in the same lexical scope, would you do it? :) — Preceding unsigned comment added by 192.139.122.66 (talk) 20:05, 7 April 2008 (UTC)
I was a radio-astronomy student for a time in the Cavendish, in the late 60s/early 70s. I remember one of the great duo giving us a talk and telling us how the Cambridge University Maths Lab guys had actually implemented the algorithm independently, for use by the radio-astronomers, well before its publication by CT. Linuxlad (talk) 20:48, 10 March 2009 (UTC)
Well Gauss discover a formula who predicts how many prime numbers are smaller than a given number N, as you may guess the formula is N/log(N), is quite similar!!!!, later Gauss was proven right. — Preceding unsigned comment added by 190.207.17.229 (talk) 02:28, 26 October 2012 (UTC)
I have just tried adding a link to a blog post of my own, but was automatically reverted, describing a python implementation of the Cooley-Tukey algorithm for general factorizations. Most sources on the web seem to concentrate on the radix 2 FFT, be the decimation in time or frequency, so it is hard to get a feeling of how anything else really works or is implemented. I felt (and still feel) this to be a worthwhile addition. Wouldn't mind rewriting it for wikipedia, if the problem is linking a personal site, but would like some pointers on how to proceed. Jfrio (talk) 11:38, 29 May 2009 (UTC)
"0.02 minutes"? Why not 1.2 seconds? The original paper may have quoted the time in this unusual fashion, but I think we'd be justified in using more conventional units in the article. Tevildo (talk) 23:48, 8 December 2009 (UTC)
The definition presented here is consistant with the definition in my undergraduate college textbook (Lathi), but both definitions seem to be inconsistant with how the mathamatics tools I use implement them. For example, MATHCAD uses the 1/N factor on the FFT, but not on the IFFT. It also does not produce the first value since the summation starts from 1 instead of from 0. On the other hand, MATLAB produces the first value and uses a factor of 1 on the FFT and the 1/N factor on the IFFT like here, but it uses a positive sign in the exponent in the exponential of the FFT rather than a negative sign, and the negative sign is used in the IFFT. So it is reversed. My questions are these: Are these variations simply wrong, or are they as correct as what is presented here and just dependant upon how they were derived? It seems that the Cooley-Tukey Method should only be one way. Am I wrong? Also, what is the impact of using these different tools and/or methods? I would suspect that cycling through the MATLAB version (i.e., time domain sampled function --> transform --> time domain sampled function) should give the same result as the version presented here. But what about the transforms while in the frequency domain? Is it comparing apples-to-apples? And then for the MATHCAD version, I suspect that this one would have to have some additional operations on it to yeild the same results as the other two, like adding back in the first term and changing the 1/N weighting by multiplying thorugh by N. Am I wrong in assuming the results would be very different? Why do they define it so differently across the board? I would think that a company selling a toolbox like this would want to stay consistant with what the math/engineering community expects. Does anyone have any insight into this? Thanks in advance. J.G. Lindley (talk) 21:48, 21 March 2011 (UTC)
In the section on input ordering, it is stated that the defining difference between DIT and DIF algorithms is the order of the inputs and outputs. This is not true. The DIT algorithm can be implemented with either the inputs or outputs in normal order. Same for the DIF algorithm. DIT applies the twiddle factor to both arms in the butterfly (twiddle factor is applied to odd component in the flowgraph before it is summed with the even), whereas with DIF, the twiddle factor is applied after the sum takes place. (the 2pt butterfly elements are different). See Richard Lyons' "Understanding Digital Signal Processing" 1997 edition pg 151. 128.59.159.194 (talk) 08:00, 4 December 2011 (UTC)
Gauss write the article in Neo Latin, because he didn't know anything of Chinesse, if in his hands were he would have writen the article in this lenguage, is a coincidence at his times, not because he could not learn Chinesse, we all know that he can. — Preceding unsigned comment added by 190.207.17.229 (talk) 01:28, 26 October 2012 (UTC)
I know its not easy to write articles on complex mathematical subjects in the tone of a general encyclopedia, but I can't make heads or tails of what is going on in here. A few additional paragraphs putting this in layman's terms would go a long way. ThemFromSpace 06:46, 14 December 2012 (UTC)
I have to agree with Themfromspace, this is one of the typical examples of maths/science articles in Wikipedia which fails across the board at sharing knowledge to non-maths/science experts. Although it may be accurate, the maths short-hand and the obscure references to other concepts that the author(s) assume are known to the reader mean that it is like struggling through a barrel of treacle to understand what is going on. I have just put together an FFT package in C++ for transforming WAV data into frequency data and it works, but my what a struggle to accumulate the information from various places online and piece together very patiently the actual process. I belive that an encyclopedia is meant primarily to share knowledge, not obfusticate it and condem it to the realms of references for people who are already experts in their field. Along with phrases such as "This is hard" or "very complicated process", maths shorthand in large blocks should be avoided absolutely unless it is provided with a key and can also be talked through in plain english. I hope that maths and science articles in Wikipedia improve dramatically in the future. I feel that this article is rich in information, but arrogant in it's delivery. And I do not wish to be rude. But I am (and I am sure others are) frustrated. - Sam, UK — Preceding unsigned comment added by 46.233.70.199 (talk) 10:15, 17 January 2014 (UTC)
Shouldn't the time complexity be O(N log log N) instead of O(N log N) when the radix, r, is always chosen to be sqrt(N)? With that choice of r, there are clearly log_2(log_2(N)) levels, (since it takes exactly log_2(log_2(N)) times taking the square root of N to reach the value 2), and each level clearly does O(N) work, (where N is the total N, not the N at a given level), by multiplying the twiddles to every number twice and adding in each number twice, (twice because first all of the rows are FFT'd and then all of the columns are FFT'd). Of course, sqrt(N) is only exact each time if N is of the form 2^(2^k), such as 4, 16, 256, 65,536, or 4,294,967,296, but the article specifically mentions that it takes O(N log N) time in the case when r=sqrt(N). This seems to conflict. Ndickson (talk) 04:44, 27 September 2013 (UTC) (Edited: Added missing 256 = 2^(2^3) to list.) Ndickson (talk) 05:00, 27 September 2013 (UTC)
I am expecting to add a code block to the Pseudocode section of Cooley-Tukey FFT algorithm based on my naive c++ implementation. I have not written pseudocode before and would like some other users to vet it before I submit it in the article. If a user thinks it is not a good idea to add such to this section please explain.
This is the current development of my pseudocode. Users are invited to make edits using the same guidelines as apply to articles.
x0,...,N−1 ← diffft2(x, N): DFT of {x0, ...,xN-1}: φT ← eπi/N for n ∈ {2-k × N : k ∈[1, log2(N))} φT ← φT × φT T ← conj(φT) for a ∈ [0, n/2) T ← φT × T for b ∈ {a + kn : k ∈ [0, (N/n))} t ← xb - xb+n/2 xb ← xb + xb+n/2 xb+n/2 ← T × t
Sorry, this is all the effort I have today. I think this approach has a lot of advantages.
This is a copy of my c++ implementation. It requires annotation so I apologize.
unsigned int reverse(unsigned int x, unsigned int maxbits) { x = (((x & 0xaaaaaaaa) >> 1) | ((x & 0x55555555) << 1)); x = (((x & 0xcccccccc) >> 2) | ((x & 0x33333333) << 2)); x = (((x & 0xf0f0f0f0) >> 4) | ((x & 0x0f0f0f0f) << 4)); x = (((x & 0xff00ff00) >> 8) | ((x & 0x00ff00ff) << 8)); return ((x >> 16) | (x << 16)) >> (32 - maxbits); } void decimate(vector<complex<double>> &input) { unsigned int N = input.size(), m = (unsigned int)log2(N); for (unsigned int x = 0; x < input.size(); x++) { unsigned int y = reverse(x, m); if (y > x) { complex<double> t = input[x]; input[x] = input[y]; input[y] = t; } } } void fft(vector<complex<double>> &x, bool normalize) { unsigned int N = x.size(), k = N, n; complex<double> phiT = 3.14159265358979 / N, T; phiT = complex<double>(cos(phiT.real()), sin(phiT.real())); while (k > 1) { n = k; k >>= 1; phiT = phiT * phiT; T.real(phiT.real()); T.imag(-phiT.imag()); for (unsigned int l = 0; l < k; l++) { T *= phiT; for (unsigned int a = l; a < N; a += n) { unsigned int b = a + k; complex<double> t = x[a] - x[b]; x[a] += x[b]; x[b] = t * T; } } } if (normalize) { complex<double> Factor = 1.0 / sqrt(N); for (int i = 0; i < N; i++) x[i] *= Factor; } }
184.17.219.58 (talk) 16:03, 30 June 2015 (UTC)
input : where is a power of , output : its discrete Fourier transform (DFT)
the naive algorithm takes multiplications/additions but we will do it in operations by a divide and conquer strategy
transform into two downsampled vectors
recursively apply the algorithm to calculate the DFT of those two vectors so that
-periodize the two obtained DFT so that and , then mix them to obtain the whole DFT
78.227.78.135 (talk) 18:27, 15 January 2016 (UTC)
I was having trouble understanding the pseudocode in the Radix 2 DIT case. I made this pseudocode, which I think is clearer.
complex X0...N fft(complex x0...N, integer N) if( N==1 ) complex X0 = x0 return X complex x_evens0...N/2 complex x_odds0...N/2 for( i = 0 to N ) if( is_even(i) ) x_evens.append( xi ) else x_odds.append( xi ) complex fft_odds0...N/2 = fft(x_odds, N/2) complex fft_evens0...N/2 = fft(x_evens, N/2) complex X0...N for( k = 0 to N/2 ) complex nth_rou = e( -2 * i * π * k ) / N Xk = fft_evensk + nth_rou * fft_oddsk Xk + N/2 = fft_evensk - nth_rou * fft_oddsk return X
I would also mention that the inverse of the fft can be taken by changing the sign of the power of the root of unity (i.e., e( 2 * i * π * k ) / N) and multiplying each term of the final output of the inverse fft by 1/N.
Thoughts?
Cormac596 (talk) 15:28, 6 May 2017 (UTC)
Given at the beginning of the article, it says N=N1*N2, so the example better to have different N1 & N2, rather than N1=N2=2. Jackzhp (talk) 15:33, 10 October 2018 (UTC)
algorithm iterative-fft is input: Array a of n complex values where n is a power of 2 output: Array A the DFT of a bit-reverse-copy(a,A) n ← a.length for s = 1 to log(n) m ← 2s ωm ← exp(−2πi/m) for k = 0 to n-1 by m ω ← 1 for j = 0 to m/2 – 1 t ← ω A[k + j + m/2] u ← A[k + j] A[k + j] ← u + t A[k + j + m/2] ← u – t ω ← ω ωm return A
At one point it says "ω ← ω ωm". What does it mean? There's no operator in the middle. — Preceding unsigned comment added by Braden1127 (talk • contribs) 00:39, 18 October 2018 (UTC)
Another question on this pseudocode - what is the 'i' in exp(−2πi/m)? — Preceding unsigned comment added by Jrogerz (talk • contribs) 20:22, 3 November 2019 (UTC)
At https://en.wikipedia.orghttps://demo.azizisearch.com/lite/wikipedia/page/Cooley%E2%80%93Tukey_FFT_algorithm#Data_reordering,_bit_reversal,_and_in-place_algorithms a piece of code is shown to implement an in-place FFT. I have used this code in a program I'm writing (though with at least one variable having a different name). Unfortunately the code here in the above linked Wikipedia page section, doesn't tell me anything about the license for that code. Is it considered public domain? Can I use it in any way I want in software I write? What about if I want to sell my software? Is the FFT code protected with GPL so that I will need to publicly reveal the source code for my software (defeating the entire purpose for me considering selling my software)? Does the fact that the code in question on Wikipedia is pseudo code (rather than implemented in a specific language like C or VB or Java) have an impact on its copyright regarding my actual implementation in a specific language (C++ in my case)? Benhut1 (talk) 22:45, 11 July 2019 (UTC)