利用c语言写fft的算法

下面是一个简单的基于C语言的FFT算法实现：

#include <stdio.h>
#include <math.h>

#define PI 3.14159265358979323846

// 定义复数结构体
typedef struct {
    double real, imag;
} Complex;

// 复数加法
Complex add(Complex a, Complex b) {
    Complex c;
    c.real = a.real + b.real;
    c.imag = a.imag + b.imag;
    return c;
}

// 复数减法
Complex sub(Complex a, Complex b) {
    Complex c;
    c.real = a.real - b.real;
    c.imag = a.imag - b.imag;
    return c;
}

// 复数乘法
Complex mul(Complex a, Complex b) {
    Complex c;
    c.real = a.real * b.real - a.imag * b.imag;
    c.imag = a.real * b.imag + a.imag * b.real;
    return c;
}

// FFT算法
void FFT(Complex *x, int n, int inv) {
    if (n == 1) return;
    Complex odd[n / 2], even[n / 2];
    for (int i = 0; i < n / 2; i++) {
        even[i] = x[i * 2];
        odd[i] = x[i * 2 + 1];
    }
    FFT(even, n / 2, inv);
    FFT(odd, n / 2, inv);
    Complex wn, w;
    for (int i = 0; i < n / 2; i++) {
        wn.real = cos(2 * PI / n * i);
        wn.imag = inv * sin(2 * PI / n * i);
        w = mul(wn, odd[i]);
        x[i] = add(even[i], w);
        x[i + n / 2] = sub(even[i], w);
    }
    if (inv == -1) {
        for (int i = 0; i < n; i++) {
            x[i].real /= n;
            x[i].imag /= n;
        }
    }
}

// 测试
int main() {
    int n = 8;
    Complex x[n];
    for (int i = 0; i < n; i++) {
        x[i].real = i + 1;
        x[i].imag = 0;
    }
    FFT(x, n, 1);
    for (int i = 0; i < n; i++) {
        printf("%lf + %lfi\n", x[i].real, x[i].imag);
    }
    FFT(x, n, -1);
    for (int i = 0; i < n; i++) {
        printf("%lf + %lfi\n", x[i].real, x[i].imag);
    }
    return 0;
}

这里实现的是基于递归的Cooley-Tukey FFT算法，可以进行任意长度的FFT变换。在测试部分，我们将一个长度为8的实数序列进行FFT变换，然后再进行IFFT逆变换，输出变换前后的结果，可以看到变换前后的序列完全一致。