/*(tw1+jtw2)(x_r[k]+jx_i[k]) *

* real : tw1*x_r[k] - tw2*x_i[k] * imag : tw1*x_i[k] + tw2*x_r[k] * 我想不抽象出一个 * typedef struct { * double real; // 实部 * double imag; // 虚部

* } complex; 以及针对complex的操作 * 来简化复数运算是否是因为效率上的考虑! */

/* 蝶形算法 */

x_r[k+p] = tmp_real + ( tw1*x_r[k+p+step] - tw2*x_i[k+p+step] );

x_i[k+p] = tmp_imag + ( tw2*x_r[k+p+step] + tw1*x_i[k+p+step] );

/* X[k] = A(k)+WB(k)

* X[k+N/2] = A(k)-WB(k) 的性质可以优化这里*/ // 旋转因子, 需要优化...

tw1 = cos(2*PI*(p+step)/pow(2, cur_layer)); tw2 = -sin(2*PI*(p+step)/pow(2, cur_layer)); x_r[k+p+step] = tmp_real + ( tw1*temp - tw2*x_i[k+p+step] );

x_i[k+p+step] = tmp_imag + ( tw2*temp + tw1*x_i[k+p+step] );

printf(\, k+p, x_r[k+p], x_i[k+p]);

printf(\, k+p+step, x_r[k+p+step], x_i[k+p+step]); }

/* 开跳!:) */ k += 2*step;

} } } /* * 后记:

* 究竟是颗粒在外层循环还是样本输入在外层, 好象也差不多, 复杂度完全一样. * 但以我资质愚钝花费了不少时间才弄明白这数十行代码.

* 从中我发现一个于我非常有帮助的教训, 很久以前我写过一部分算法, 其中绝大多数都是递归.

* 将数据量减少, 减少再减少, 用归纳的方式来找出数据量加大代码的规律 * 比如FFT

* 1. 先写死LayerI的代码; 然后再把LayerI的输出作为LayerII的输入, 又写死代码; ......

* 大约3层就可以统计出规律来. 这和递归也是一样, 先写死一两层, 自然就出来了! * 2. 有的功能可以写伪代码, 不急于求出结果, 降低复杂性, 把逻辑结果定出来后再添加. * 比如旋转因子就可以写死, 就写1.0. 流程出来后再写旋转因子. * 寥寥数语, 我可真是流了不少汗! Happy! */

void dft( void ) {

int i, n, k, tx1, tx2; float tw1,tw2;

float xx_r[N],xx_i[N]; /*

* clear any data in Real and Imaginary result arrays prior to DFT */

for(k=0; k<=N-1; k++)

xx_r[k] = xx_i[k] = x_i[k] = 0.0;

// caculate the DFT for(k=0; k<=(N-1); k++) {

for(n=0; n<=(N-1); n++) {

tx1 = (n*k); tx2 = tx1+(3*N)/4; tx1 = tx1%(N); tx2 = tx2%(N); if(tx1 >= (N/2))

tw1 = -twiddle[tx1-(N/2)]; else

tw1 = twiddle[tx1]; if(tx2 >= (N/2))

tw2 = -twiddle[tx2-(N/2)]; else

tw2 = twiddle[tx2];

xx_r[k] = xx_r[k]+x_r[n]*tw1; xx_i[k] = xx_i[k]+x_r[n]*tw2; }

xx_i[k] = -xx_i[k]; }

// display

for(i=0; i

printf(\, xx_r[i], xx_i[i]); } //

---------------------------------------------------------------------------

int main( void ) {

fft_init( ); bitrev( ); // bitrev2( ); //fft1( ); fft2( ); display( );

system( \ ); // dft(); return 1; }

#include #include

/********************************************************************* 快速福利叶变换C函数

函数简介：此函数是通用的快速傅里叶变换C语言函数，移植性强，以下部分不依 赖硬件。此函数采用联合体的形式表示一个复数，输入为自然顺序的复 数（输入实数是可令复数虚部为0），输出为经过FFT变换的自然顺序的 复数

使用说明：使用此函数只需更改宏定义FFT_N的值即可实现点数的改变，FFT_N的 应该为2的N次方，不满足此条件时应在后面补0 函数调用：FFT(s);