2
$\begingroup$

I've written my own FFT but wanted to know if there was perhaps a way to combine the Cooley and Turkey methods. So far I've got N/2 log N real multiplies for transforming a real signal.

My current implementation has N log N multiplies (well a little more because I've not optimised out the case of the complex sinusoid at 0 radians), but I have an optimisation that reduces the number of multiplies by 50%.

Note: this is absolutely not homework, I've written this algorithm because I found the explanations of the Cooley/Turkey method far too convoluted, and at the same time I wanted the ability to select which bins I processed in both the FFT and IFFT as well as being able to output the IFFT to multiple channels. So I just decided to figure out my own radix-2 while at the same time understanding how and why the transform works (which I have).

ORIGINAL QUESTION (edited for clarity) ---

Example code: kfft.c

I wanted to know whether there is a similar FFT algorithm to mine.

I basically perform FFT by factorising calculations by recursively subtracting the second half of a signal from the first, meaning x[n] only needs to be multiplied by the complex sinusoid for only half of the duration, resulting in N/2 log N real multiplies (shown below). The summed signal is for multiplying with complex sinusoids that are multiples of the sampling period. For example groups {2, 6} are multiples of 2, and therefore have a periodicity of 4 samples.

// For bins { 1, 3, 5 and 7 }

x1[0] = x0[0] - x0[4]
x1[1] = x0[1] - x0[5]
x1[2] = x0[2] - x0[6]
x1[3] = x0[3] - x0[7]
// Calculate sum (required for the next level of bins)
x0[0] = x0[0] + x0[4]
x0[1] = x0[1] + x0[5]
x0[2] = x0[2] + x0[6]
x0[3] = x0[3] + x0[7]

// For bins { 2 and 6 }

x2[0] = x0[0] - x0[2]
x2[1] = x0[1] - x0[3]
// Calculate sum (required for the next level of bins)
x0[2] = x0[0] + x0[2]
x0[3] = x0[1] + x0[3]

// For bin { 4 }
x3[0] = x0[0] - x0[1]
x3[1] = x0[0] + x0[1] // <-- This is the 0hz bin (sum)
$\endgroup$
  • 1
    $\begingroup$ Your scheme is hardly any different from what the FFT does in the case when $n = 2^m$. $\endgroup$ – Dilip Sarwate Sep 15 '14 at 14:19
  • $\begingroup$ Well that's what I mean. It is FFT, but not quite radix 2, or is it? I'm really not sure as I haven't quite wrapped my head around radix-2 completely (instead I just wrote this). Also does what I said make much sense, such as considering e^ (-i * 2pi * k * n) as a resonator? $\endgroup$ – Keldon Alleyne Sep 15 '14 at 14:27
  • $\begingroup$ Why the down vote with no helpful comment or anything? Some strange people??? The only algorithm I've seen is radix-2/4, so surely it's understandable that i might not know this same process is the more general form when I've not seen it described in this way anywhere. $\endgroup$ – Keldon Alleyne Sep 15 '14 at 15:02
  • $\begingroup$ most of us, during or since our college daze, have written a simple Tukey-Cooley radix-2 FFT in some programming language. i have done so in C (can send you the code if i can find it) and in MC68000 asm and have modified the DSP56000 asm FFT to be a either unitary ($\frac{1}{\sqrt{N}}$) gain or a decent block-floating-point (mot's code sucked). the only thing i can think of is to write the code at the beginning of each group so that multiplication by $1$ or by $j$ is avoided in the butterfly. $\endgroup$ – robert bristow-johnson Sep 18 '14 at 14:32
  • $\begingroup$ Yes that is one of my planned optimisations, but with an upcoming optimisation I'm now down to N/2 log N real multiplies, so I'm wondering whether there are any further optimisation left by perhaps borrowing from optimisations used in other algorithms $\endgroup$ – Keldon Alleyne Sep 18 '14 at 19:50
1
$\begingroup$

Indeed, this is more or less the usual butterfly used in radix-2 FFT (can be extended to any prime number): http://en.wikipedia.org/wiki/Butterfly_diagram

$\endgroup$
  • $\begingroup$ Thank you Matthieu and Dilip. And I can see clearly how and why the butterfly works. $\endgroup$ – Keldon Alleyne Sep 15 '14 at 15:16

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.