comparison src/vamp-sdk/FFTimpl.cpp @ 460:b409560a805b

Merge from branch vampipe

 
-/* Public domain FFT implementation from Don Cross. */
+// Override C linkage for KissFFT headers. So long as we have already
+// included all of the other (system etc) headers KissFFT depends on,
+// this should work out OK
+#undef __cplusplus
 
-static void
-fft(unsigned int n, bool inverse,
-    const double *ri, const double *ii,
-    double *ro, double *io)
-{
-    if (!ri || !ro || !io) return;
+namespace Kiss {
 
-    unsigned int bits;
-    unsigned int i, j, k, m;
-    unsigned int blockSize, blockEnd;
+#undef KISS_FFT_H
+#undef KISS_FTR_H
+#undef KISS_FFT__GUTS_H
+#undef FIXED_POINT
+#undef USE_SIMD
+#undef kiss_fft_scalar
 
-    double tr, ti;
-
-    if (n < 2) return;
-    if (n & (n-1)) return;
-
-    double angle = 2.0 * M_PI;
-    if (inverse) angle = -angle;
-
-    for (i = 0; ; ++i) {
-	if (n & (1 << i)) {
-	    bits = i;
-	    break;
-	}
-    }
-
-#ifdef _MSC_VER
-    int *table = (int *)_malloca(n * sizeof(int));
+#ifdef SINGLE_PRECISION_FFT
+#pragma message("Using single-precision FFTs")
+typedef float kiss_fft_scalar;
+#define kiss_fft_scalar float
 #else
-    int table[n];
+typedef double kiss_fft_scalar;
+#define kiss_fft_scalar double
 #endif
 
-    for (i = 0; i < n; ++i) {
-        m = i;
-        for (j = k = 0; j < bits; ++j) {
-            k = (k << 1) | (m & 1);
-            m >>= 1;
-        }
-        table[i] = k;
-    }
+inline void free(void *ptr) { ::free(ptr); }
+#include "ext/kiss_fft.c"
+#include "ext/kiss_fftr.c"
 
-    if (ii) {
-	for (i = 0; i < n; ++i) {
-	    ro[table[i]] = ri[i];
-	    io[table[i]] = ii[i];
-	}
-    } else {
-	for (i = 0; i < n; ++i) {
-	    ro[table[i]] = ri[i];
-	    io[table[i]] = 0.0;
-	}
-    }
+#undef kiss_fft_scalar // leaving only the namespaced typedef
 
-    blockEnd = 1;
-
-    for (blockSize = 2; blockSize <= n; blockSize <<= 1) {
-
-	double delta = angle / (double)blockSize;
-	double sm2 = -sin(-2 * delta);
-	double sm1 = -sin(-delta);
-	double cm2 = cos(-2 * delta);
-	double cm1 = cos(-delta);
-	double w = 2 * cm1;
-	double ar[3], ai[3];
-
-	for (i = 0; i < n; i += blockSize) {
-
-	    ar[2] = cm2;
-	    ar[1] = cm1;
-
-	    ai[2] = sm2;
-	    ai[1] = sm1;
-
-	    for (j = i, m = 0; m < blockEnd; j++, m++) {
-
-		ar[0] = w * ar[1] - ar[2];
-		ar[2] = ar[1];
-		ar[1] = ar[0];
-
-		ai[0] = w * ai[1] - ai[2];
-		ai[2] = ai[1];
-		ai[1] = ai[0];
-
-		k = j + blockEnd;
-		tr = ar[0] * ro[k] - ai[0] * io[k];
-		ti = ar[0] * io[k] + ai[0] * ro[k];
-
-		ro[k] = ro[j] - tr;
-		io[k] = io[j] - ti;
-
-		ro[j] += tr;
-		io[j] += ti;
-	    }
-	}
-
-	blockEnd = blockSize;
-    }
-
-    if (inverse) {
-
-	double denom = (double)n;
-
-	for (i = 0; i < n; i++) {
-	    ro[i] /= denom;
-	    io[i] /= denom;
-	}
-    }
-
-#ifdef _MSC_VER
-    _freea(table);
-#endif
 }