fftw

原文：http://bbs.csdn.net/topics/330133751

#include <stdio.h> #include <stdlib.h> #include <string.h> #include <math.h> #include <time.h> #include "fftw3.h" #define NNx 1025 #define NNy 1025 // ================= Main ================= // int main(int argc, char *argv[]) {     fftw_complex *in, *out, *in2;     fftw_plan p,q;     clock_t start, end;     int i=0, j=0, k=0;     in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) *NNx*NNy );     in2 = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) *NNx*NNy );     out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) *NNx*NNy );     for(  i=0; i < NNx*NNy; i++)     {         in[i][0] = i;         in[i][1] = 0.0;     }     start=clock();     for( k=0; k<100; k++)     {     p=fftw_plan_dft_2d( NNx, NNy, in, out, FFTW_FORWARD, FFTW_ESTIMATE);     q=fftw_plan_dft_2d( NNx, NNy, out, in2, FFTW_BACKWARD, FFTW_ESTIMATE);     fftw_execute( p ); // repeat as needed//     fftw_execute( q );     }     end=clock();     printf(" FFTW : %f second\n", (double)  (end-start)/CLOCKS_PER_SEC  );     fftw_destroy_plan(p);     fftw_destroy_plan(q);     fftw_free(in);     fftw_free(out);  system("PAUSE");  return 0;

}

建议优化：

 // 这里初始化plan可能要花费几秒钟

     p=fftw_plan_dft_2d( NNx, NNy, in, out, FFTW_FORWARD, FFTW_MEASURE);       q=fftw_plan_dft_2d( NNx, NNy, out, in2, FFTW_BACKWARD, FFTW_MEASURE);        for (  i=0; i < NNx*NNy; i++)       {           in[i][0] = i;          in[i][1] = 0.0;      }         // 以下执行部分将获得最佳性能      start= clock ();       for ( k=0; k <100; k++)       {       fftw_execute_dft( p, in, out );       fftw_execute_dft( q, out, in2 );       }      end =  clock (); 备注：

如果in和out指针相同为原位运算，否则为非原位运算。

sign可以为正变换FFTW_FORWARD(-1)，也可以为逆变换FFTW_BACKWORD(+1)，实际上就是变换公式中指数项的符号。需注意FFTW的逆变换没有除以N，即数据正变换再反变换后是原始数据的N倍。