Move the ALcomplex and FFT functions to a separate file

2 files changed, 124 insertions, 0 deletions

diff --git a/common/alcomplex.c b/common/alcomplex.c
new file mode 100644
index 00000000..c1a312e3
--- /dev/null
+++ b/common/alcomplex.c
@@ -0,0 +1,62 @@
+
+#include "config.h"
+
+#include "alcomplex.h"
+
+#include <math.h>
+
+
+extern inline ALcomplex complex_add(ALcomplex a, ALcomplex b);
+extern inline ALcomplex complex_sub(ALcomplex a, ALcomplex b);
+extern inline ALcomplex complex_mult(ALcomplex a, ALcomplex b);
+
+void complex_fft(ALcomplex *FFTBuffer, ALsizei FFTSize, ALdouble Sign)
+{
+    ALsizei i, j, k, mask, step, step2;
+    ALcomplex temp, u, w;
+    ALdouble arg;
+
+    /* Bit-reversal permutation applied to a sequence of FFTSize items */
+    for(i = 1;i < FFTSize-1;i++)
+    {
+        for(mask = 0x1, j = 0;mask < FFTSize;mask <<= 1)
+        {
+            if((i&mask) != 0)
+                j++;
+            j <<= 1;
+        }
+        j >>= 1;
+
+        if(i < j)
+        {
+            temp         = FFTBuffer[i];
+            FFTBuffer[i] = FFTBuffer[j];
+            FFTBuffer[j] = temp;
+        }
+    }
+
+    /* Iterative form of Danielson–Lanczos lemma */
+    for(i = 1, step = 2;i < FFTSize;i<<=1, step<<=1)
+    {
+        step2 = step >> 1;
+        arg   = M_PI / step2;
+
+        w.Real = cos(arg);
+        w.Imag = sin(arg) * Sign;
+
+        u.Real = 1.0;
+        u.Imag = 0.0;
+
+        for(j = 0;j < step2;j++)
+        {
+            for(k = j;k < FFTSize;k+=step)
+            {
+                temp               = complex_mult(FFTBuffer[k+step2], u);
+                FFTBuffer[k+step2] = complex_sub(FFTBuffer[k], temp);
+                FFTBuffer[k]       = complex_add(FFTBuffer[k], temp);
+            }
+
+            u = complex_mult(u, w);
+        }
+    }
+}
diff --git a/common/alcomplex.h b/common/alcomplex.h
new file mode 100644
index 00000000..cf4683fa
--- /dev/null
+++ b/common/alcomplex.h
@@ -0,0 +1,62 @@
+#ifndef ALCOMPLEX_H
+#define ALCOMPLEX_H
+
+#include "AL/al.h"
+
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+typedef struct ALcomplex {
+    ALdouble Real;
+    ALdouble Imag;
+} ALcomplex;
+
+/** Addition of two complex numbers. */
+inline ALcomplex complex_add(ALcomplex a, ALcomplex b)
+{
+    ALcomplex result;
+
+    result.Real = a.Real + b.Real;
+    result.Imag = a.Imag + b.Imag;
+
+    return result;
+}
+
+/** Subtraction of two complex numbers. */
+inline ALcomplex complex_sub(ALcomplex a, ALcomplex b)
+{
+    ALcomplex result;
+
+    result.Real = a.Real - b.Real;
+    result.Imag = a.Imag - b.Imag;
+
+    return result;
+}
+
+/** Multiplication of two complex numbers. */
+inline ALcomplex complex_mult(ALcomplex a, ALcomplex b)
+{
+    ALcomplex result;
+
+    result.Real = a.Real*b.Real - a.Imag*b.Imag;
+    result.Imag = a.Imag*b.Real + a.Real*b.Imag;
+
+    return result;
+}
+
+/**
+ * Iterative implementation of 2-radix FFT (In-place algorithm). Sign = -1 is
+ * FFT and 1 is iFFT (inverse). Fills FFTBuffer[0...FFTSize-1] with the
+ * Discrete Fourier Transform (DFT) of the time domain data stored in
+ * FFTBuffer[0...FFTSize-1]. FFTBuffer is an array of complex numbers, FFTSize
+ * MUST BE power of two.
+ */
+void complex_fft(ALcomplex *FFTBuffer, ALsizei FFTSize, ALdouble Sign);
+
+#ifdef __cplusplus
+} // extern "C"
+#endif
+
+#endif /* ALCOMPLEX_H */