FFT的C语言算法实现

FFT的C语言算法实现程序如下：/************FFT***********/#includestdio.h#includemath.h#includestdlib.h#defineN1000typedefstruct{doublereal;doubleimg;}complex;voidfft();/*快速傅里叶变换*/voidifft();/*快速傅里叶逆变换*/voidinitW();voidchange();voidadd(complex,complex,complex*);/*复数加法*/voidmul(complex,complex,complex*);/*复数乘法*/voidsub(complex,complex,complex*);/*复数减法*/voiddivi(complex,complex,complex*);/*复数除法*/voidoutput();/*输出结果*/complexx[N],*W;/*输出序列的值*/intsize_x=0;/*输入序列的长度，只限2的N次方*/doublePI;intmain(){inti,method;system(cls);PI=atan(1)*4;printf(Pleaseinputthesizeofx:\n);/*输入序列的长度*/scanf(%d,&size_x);printf(Pleaseinputthedatainx[N]:(suchas:56)\n);/*输入序列对应的值*/for(i=0;isize_x;i++)scanf(%lf%lf,&x[i].real,&x[i].img);initW();/*选择FFT或逆FFT运算*/printf(UseFFT(0)orIFFT(1)?\n);scanf(%d,&method);if(method==0)fft();elseifft();output();return0;}/*进行基-2FFT运算*/voidfft(){inti=0,j=0,k=0,l=0;complexup,down,product;change();for(i=0;ilog(size_x)/log(2);i++){l=1i;for(j=0;jsize_x;j+=2*l){for(k=0;kl;k++){mul(x[j+k+l],W[size_x*k/2/l],&product);add(x[j+k],product,&up);sub(x[j+k],product,&down);x[j+k]=up;x[j+k+l]=down;}}}}voidifft(){inti=0,j=0,k=0,l=size_x;complexup,down;for(i=0;i(int)(log(size_x)/log(2));i++)/*蝶形运算*/{l/=2;for(j=0;jsize_x;j+=2*l){for(k=0;kl;k++){add(x[j+k],x[j+k+l],&up);up.real/=2;up.img/=2;sub(x[j+k],x[j+k+l],&down);down.real/=2;down.img/=2;divi(down,W[size_x*k/2/l],&down);x[j+k]=up;x[j+k+l]=down;}}}change();}voidinitW(){inti;W=(complex*)malloc(sizeof(complex)*size_x);for(i=0;isize_x;i++){W[i].real=cos(2*PI/size_x*i);W[i].img=-1*sin(2*PI/size_x*i);}}voidchange(){complextemp;unsignedshorti=0,j=0,k=0;doublet;for(i=0;isize_x;i++){k=i;j=0;t=(log(size_x)/log(2));while((t--)0){j=j1;j|=(k&1);k=k1;}if(ji){temp=x[i];x[i]=x[j];x[j]=temp;}}}voidoutput()/*输出结果*/{inti;printf(Theresultareasfollows\n);for(i=0;isize_x;i++){printf(%.4f,x[i].real);if(x[i].img=0.0001)printf(+%.4fj\n,x[i].img);elseif(fabs(x[i].img)0.0001)printf(\n);elseprintf(%.4fj\n,x[i].img);}}voidadd(complexa,complexb,complex*c){c-real=a.real+b.real;c-img=a.img+b.img;}voidmul(complexa,complexb,complex*c){c-real=a.real*b.real-a.img*b.img;c-img=a.real*b.img+a.img*b.real;}voidsub(complexa,complexb,complex*c){c-real=a.real-b.real;c-img=a.img-b.img;}voiddivi(complexa,complexb,complex*c){c-real=(a.real*b.real+a.img*b.img)/(b.real*b.real+b.img*b.img);c-img=(a.img*b.real-a.real*b.img)/(b.real*b.real+b.img*b.img);}