Add inverse FFT

1 files changed, 27 insertions, 4 deletions

diff --git a/fft.c b/fft.c
index b8af1f4..f903dd4 100644
--- a/fft.c
+++ b/fft.c
@@ -4,12 +4,13 @@
 #include "fft.h"
 
 static void
-do_fft (complex_t *buffer, int n, complex_t principal)
+do_fft (complex_t *buffer, int n, complex_t principal, int is_inverse)
 {
     complex_t *temp, *even, *odd;
     complex_t alpha;
     complex_t p2;
     int i;
+    double inv;
     
     if (n == 1)
 	return;
@@ -26,17 +27,28 @@ do_fft (complex_t *buffer, int n, complex_t principal)
 
     p2 = complex_mul (principal, principal);
     
-    do_fft (even, n/2, p2);
-    do_fft (odd, n/2, p2);
+    do_fft (even, n/2, p2, FALSE);
+    do_fft (odd, n/2, p2, FALSE);
 
     alpha.re = 1.0;
     alpha.im = 0;
+
+    if (is_inverse)
+	inv = 1.0 / n;
     
     for (i = 0; i < n/2; ++i)
     {
 	buffer[i]	= complex_add (even[i], complex_mul (alpha, odd[i]));
 	buffer[i + n/2] = complex_sub (even[i], complex_mul (alpha, odd[i]));
 
+	if (is_inverse)
+	{
+	    buffer[i].re *= inv;
+	    buffer[i].im *= inv;
+	    buffer[i + n/2].re *= inv;
+	    buffer[i + n/2].im *= inv;
+	}
+	
 	alpha = complex_mul (alpha, principal);
     }
 
@@ -51,5 +63,16 @@ fft (complex_t *buffer, int n)
     p.re = cos (2 * M_PI / n);
     p.im = sin (2 * M_PI / n);
     
-    do_fft (buffer, n, p);
+    do_fft (buffer, n, p, FALSE);
+}
+
+void
+ifft (complex_t *buffer, int n)
+{
+    complex_t p;
+
+    p.re = cos (- (2 * M_PI / n));
+    p.im = sin (- (2 * M_PI / n));
+
+    do_fft (buffer, n, p, TRUE);
 }