// test_complex_fft.c
#include "kiss_fft.h"
#include "kiss_fftr.h"
#include "kiss_fftnd.h"
#include "kiss_fftndr.h"
// #include "kiss_fastfir.h"
#include <stdio.h>

int test_complex(); //{ printf("\n[1] 复数FFT\n"); /* 代码同示例1 */ }
int test_real(); //{ printf("\n[2] 实数FFT\n"); /* 代码同示例2 */ }
int test_2d(); //{ printf("\n[3] 2D FFT\n"); /* 代码同示例3 */ }
int test_2d_real(); //{ printf("\n[4] 2D实数FFT\n"); /* 代码同示例4 */ }
int test_conv(); //{ printf("\n[5] 快速卷积\n"); /* 代码同示例5 */ }

int main() {
    printf("=== KissFFT 功能测试集 ===\n");
    test_complex();
    test_real();
    test_2d();
    test_2d_real();
    // test_conv();
    printf("\n所有测试完成！\n");
    return 0;
}

// int test_conv() {
//     int N = 8;

//     // 创建低通滤波器（系数）
//     float coeffs[8] = {0.125, 0.125, 0.125, 0.125, 0.125, 0.125, 0.125, 0.125}; // 均值滤波

//     // 初始化快速FIR
//     kiss_fastfir_cfg cfg = kiss_fastfir_alloc(coeffs, 8, NULL, NULL);

//     // 输入信号
//     float in[8] = {1, 2, 3, 4, 5, 6, 7, 8};

//     // 输出信号
//     float out[8];

//     // 执行快速卷积
//     kiss_fastfir(cfg, in, out, 8);

//     printf("=== 快速卷积测试 ===\n");
//     printf("输入信号: ");
//     for (int i = 0; i < 8; i++) printf("%.1f ", in[i]);

//     printf("\n滤波器: ");
//     for (int i = 0; i < 8; i++) printf("%.3f ", coeffs[i]);

//     printf("\n输出信号: ");
//     for (int i = 0; i < 8; i++) printf("%.3f ", out[i]);

//     printf("\n\n说明：输出是输入与滤波器的卷积结果\n");

//     kiss_fastfir_free(cfg);
//     return 0;
// }

int test_2d_real() {
    int dims[2] = {4, 4};

    // 2D实数FFT配置
    kiss_fftndr_cfg cfg_forward = kiss_fftndr_alloc(dims, 2, 0, NULL, NULL);
    kiss_fftndr_cfg cfg_inverse = kiss_fftndr_alloc(dims, 2, 1, NULL, NULL);

    // 2D实数输入
    float in_real[16];
    for (int i = 0; i < 16; i++) {
        in_real[i] = (i % 4 == 0) ? 1.0 : 0.0;  // 第0列为1
    }

    // 输出（复数，压缩格式）
    kiss_fft_cpx out_freq[16];

    // IFFT输出
    float out_real[16];

    // 正向变换
    kiss_fftndr(cfg_forward, in_real, out_freq);

    // 逆向变换
    kiss_fftndri(cfg_inverse, out_freq, out_real);

    printf("=== 2D实数FFT测试 ===\n");
    printf("输入（实数）:\n");
    for (int i = 0; i < 4; i++) {
        for (int j = 0; j < 4; j++) {
            printf("%.1f ", in_real[i*4+j]);
        }
        printf("\n");
    }

    printf("\n重建（实数）:\n");
    for (int i = 0; i < 4; i++) {
        for (int j = 0; j < 4; j++) {
            printf("%.1f ", out_real[i*4+j]);
        }
        printf("\n");
    }

    free(cfg_forward);
    free(cfg_inverse);
    return 0;
}

int test_2d() {
    int dims[2] = {4, 4};  // 4x4矩阵

    // 分配配置
    kiss_fftnd_cfg cfg_forward = kiss_fftnd_alloc(dims, 2, 0, NULL, NULL);
    kiss_fftnd_cfg cfg_inverse = kiss_fftnd_alloc(dims, 2, 1, NULL, NULL);

    // 2D输入数据（4x4 = 16个点）
    kiss_fft_cpx in[16];
    for (int i = 0; i < 16; i++) {
        in[i].r = (i % 4 == 0) ? 1.0 : 0.0;  // 第0列为1，其他为0
        in[i].i = 0;
    }

    kiss_fft_cpx out[16];
    kiss_fft_cpx reconstructed[16];

    // 2D FFT
    kiss_fftnd(cfg_forward, in, out);

    // 2D IFFT
    kiss_fftnd(cfg_inverse, out, reconstructed);

    printf("=== 2D FFT测试 (4x4) ===\n");
    printf("输入矩阵（只第0列为1）:\n");
    for (int i = 0; i < 4; i++) {
        for (int j = 0; j < 4; j++) {
            printf("%.1f\t", in[i*4+j].r);
        }
        printf("\n");
    }

    printf("\n2D FFT频谱（幅度）:\n");
    for (int i = 0; i < 4; i++) {
        for (int j = 0; j < 4; j++) {
            float mag = sqrt(out[i*4+j].r * out[i*4+j].r + out[i*4+j].i * out[i*4+j].i);
            printf("%.2f\t", mag);
        }
        printf("\n");
    }

    printf("\n重建矩阵:\n");
    for (int i = 0; i < 4; i++) {
        for (int j = 0; j < 4; j++) {
            printf("%.1f\t", reconstructed[i*4+j].r);
        }
        printf("\n");
    }

    free(cfg_forward);
    free(cfg_inverse);
    return 0;
}

int test_real() {
    int N = 8;

    // 实数FFT配置
    kiss_fftr_cfg cfg_forward = kiss_fftr_alloc(N, 0, NULL, NULL);
    kiss_fftr_cfg cfg_inverse = kiss_fftr_alloc(N, 1, NULL, NULL);

    // 实数输入信号
    float in_real[N];
    for (int i = 0; i < N; i++) {
        in_real[i] = cos(2 * M_PI * i / N) + 0.5 * sin(4 * M_PI * i / N);
    }

    // FFT输出（复数，N/2+1个点）
    kiss_fft_cpx out_freq[N/2+1];

    // IFFT输出（实数）
    float out_real[N];

    // 正向实数FFT
    kiss_fftr(cfg_forward, in_real, out_freq);

    // 逆向实数FFT
    kiss_fftri(cfg_inverse, out_freq, out_real);

    printf("=== 实数FFT测试 (N=%d) ===\n", N);
    printf("原始实数信号 -> FFT -> IFFT\n\n");
    printf("索引\t原始\t\t频谱(复数)\t\t重建\n");
    for (int i = 0; i < N; i++) {
        if (i <= N/2) {
            printf("%d\t%.3f\t\t%.3f%+.3fj\t%.3f\n",
                   i, in_real[i], out_freq[i].r, out_freq[i].i, out_real[i]);
        } else {
            printf("%d\t%.3f\t\t-\t\t%.3f\n", i, in_real[i], out_real[i]);
        }
    }

    // 验证
    float error = 0;
    for (int i = 0; i < N; i++) {
        error += fabs(in_real[i] - out_real[i]);
    }
    printf("\n重建误差: %f\n", error / N);

    kiss_fftr_free(cfg_forward);
    kiss_fftr_free(cfg_inverse);
    return 0;
}

int test_complex() {
    int N = 8;

    // 分配配置
    kiss_fft_cfg cfg_forward = kiss_fft_alloc(N, 0, NULL, NULL);
    kiss_fft_cfg cfg_inverse = kiss_fft_alloc(N, 1, NULL, NULL);

    // 输入信号：复数指数 e^(j*2π*k/N)
    kiss_fft_cpx in[N];
    for (int i = 0; i < N; i++) {
        in[i].r = cos(2 * M_PI * i / N);  // 实部
        in[i].i = sin(2 * M_PI * i / N);  // 虚部
    }

    // 输出缓冲区
    kiss_fft_cpx out[N];
    kiss_fft_cpx reconstructed[N];

    // 正向FFT
    kiss_fft(cfg_forward, in, out);

    // 逆向FFT
    kiss_fft(cfg_inverse, out, reconstructed);

    // 打印结果
    printf("=== 复数FFT测试 (N=%d) ===\n", N);
    printf("原始信号 -> FFT -> IFFT\n\n");
    printf("原始\t\t\tFFT输出\t\t\t重建信号\n");
    for (int i = 0; i < N; i++) {
        printf("%.2f%+.2fj\t%.2f%+.2fj\t%.2f%+.2fj\n",
               in[i].r, in[i].i,
               out[i].r, out[i].i,
               reconstructed[i].r, reconstructed[i].i);
    }

    // 验证精度
    float error = 0;
    for (int i = 0; i < N; i++) {
        error += fabs(in[i].r - reconstructed[i].r) + fabs(in[i].i - reconstructed[i].i);
    }
    printf("\n重建误差: %f\n", error / N);

    free(cfg_forward);
    free(cfg_inverse);
    return 0;
}