fftpack库的学习

fftpack库,使用的是pffft库中下载的，不需要额外配置，只需要fftpack.h和fftpack.c两个文件即可。使用也比较简单，有实数fft和复数fft两种。

1、复数例子

/复数fft

void cfft_test()

{

    int N = 16;

    int Nfloat = 2 * N;

    int Nbytes = Nfloat * sizeof(fftpack_real);

    fftpack_real *wrk   = (fftpack_real*)malloc(2 * Nbytes + 15 * sizeof(fftpack_real));

    fftpack_real *refin = (fftpack_real*)malloc(Nbytes);





    printf("==================\ncomplex orgin data:\n");

    for (int k = 0; k < Nfloat; k += 2) {

         refin[k] = k*0.5f;

         refin[k + 1] = 0 * k*0.5f;

         printf("%d %f %f\n", k / 2, refin[k], refin[k + 1]);

    }

    cffti(N, wrk);

    cfftf(N, refin, wrk);



    printf("==================\nfft:\n");

    for (int i = 0; i < 2 * N; i += 2)

         printf("%d  %f  %f\n", i / 2, refin[i], refin[i + 1]);

    printf("======================\nifft:\n");

    cfftb(N, refin, wrk);

    for (int i = 0; i < 2 * N; i += 2)

         printf("%d  %f  %f\n", i / 2, refin[i] / N, refin[i + 1] / N);

    free(wrk);

    free(refin);

}

2、实数例子

//实数fft

void rfft_test()

{

    int N = 16;

    int Nfloat = N;

    int Nbytes = Nfloat * sizeof(fftpack_real);

    fftpack_real *wrk    = (fftpack_real*)malloc(2 * Nbytes + 15 * sizeof(fftpack_real));

    fftpack_real *refin  = (fftpack_real*)malloc(Nbytes);

    fftpack_real *refinr = (fftpack_real*)malloc(2 * Nbytes);



    printf("==================\nreal orgin data:\n");

    for (int k = 0; k < Nfloat; k++) {

       refin[k] = k*1.0f;

       printf("%d %f \n", k, refin[k]);

    }

    rffti(N, wrk);

    rfftf(N, refin, wrk);



    printf("==================\nfft:\n");

    {

       refinr[0] = refin[0];

       refinr[1] = 0;

       refinr[Nfloat + 1] = 0;

       for (int i = 2; i <= Nfloat; i++)

       {

           refinr[i] = refin[i - 1];

       }

       for (int i = Nfloat + 2; i <= 2 * Nfloat - 1; i++)

       {

           if (i % 2)

              refinr[i] = -refinr[2 * Nfloat - i + 2];

           else

              refinr[i] = refinr[2 * Nfloat - i];

       }

    }





    for (int i = 0; i < 2 * Nfloat; i += 2)

       printf("%d  %f   %f\n", i / 2, refinr[i], refinr[i + 1]);

    printf("======================\nifft:\n");

    rfftb(N, refin, wrk);

    for (int i = 0; i < Nfloat; i++)

       printf("%d  %f \n", i / 2, refin[i] / N);



    free(wrk);

    free(refin);

    free(refinr);

}

3、分析

1. 与matlab结果比对一致,运算速度还可以。
2. 支持double和float。
3. 使用简单。