diff options
author | Andreas Eversberg <jolly@eversberg.eu> | 2017-11-16 18:25:03 +0100 |
---|---|---|
committer | Andreas Eversberg <jolly@eversberg.eu> | 2017-11-25 19:23:57 +0100 |
commit | cb30c96bc67d195bde2586161d383cf14d1cae8c (patch) | |
tree | 8e494fc2dc38be58a3e4c2198cd14b96cdb8f08d /src/common/fft.c | |
parent | 40c24a70b2a295073428ec78c281e86a7a84f7ca (diff) |
Restructure: Move fft from common code to 'libfft'
Diffstat (limited to 'src/common/fft.c')
-rw-r--r-- | src/common/fft.c | 96 |
1 files changed, 0 insertions, 96 deletions
diff --git a/src/common/fft.c b/src/common/fft.c deleted file mode 100644 index fb819e1..0000000 --- a/src/common/fft.c +++ /dev/null @@ -1,96 +0,0 @@ -/* Fast Fourier Transformation (FFT) - * - * (C) 2016 by Andreas Eversberg <jolly@eversberg.eu> - * All Rights Reserved - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - */ - -#include <math.h> -#include <string.h> -#include "fft.h" - -/* - * Code based closely to work by Paul Bourke - * - * This computes an in-place complex-to-complex FFT - * x and y are the real and imaginary arrays of 2^m points. - * dir = 1 gives forward transform - * dir = -1 gives reverse transform - */ -void fft_process(int dir, int m, double *x, double *y) -{ - int n, i, i1, j, k, i2, l, l1, l2; - double c1, c2, tx, ty, t1, t2, u1, u2, z; - - /* Calculate the number of points */ - n = 1 << m; - - /* Do the bit reversal */ - i2 = n >> 1; - j = 0; - for (i = 0; i < n - 1; i++) { - if (i < j) { - tx = x[i]; - ty = y[i]; - x[i] = x[j]; - y[i] = y[j]; - x[j] = tx; - y[j] = ty; - } - k = i2; - while (k <= j) { - j -= k; - k >>= 1; - } - j += k; - } - - /* Compute the FFT */ - c1 = -1.0; - c2 = 0.0; - l2 = 1; - for (l = 0; l < m; l++) { - l1 = l2; - l2 <<= 1; - u1 = 1.0; - u2 = 0.0; - for (j = 0; j < l1; j++) { - for (i = j; i < n; i += l2) { - i1 = i + l1; - t1 = u1 * x[i1] - u2 * y[i1]; - t2 = u1 * y[i1] + u2 * x[i1]; - x[i1] = x[i] - t1; - y[i1] = y[i] - t2; - x[i] += t1; - y[i] += t2; - } - z = u1 * c1 - u2 * c2; - u2 = u1 * c2 + u2 * c1; - u1 = z; - } - c2 = sqrt((1.0 - c1) / 2.0); - if (dir == 1) - c2 = -c2; - c1 = sqrt((1.0 + c1) / 2.0); - } - - /* Scaling for forward transform */ - if (dir == 1) { - for (i = 0; i < n; i++) { - x[i] /= n; - y[i] /= n; - } - } -} |