C++实现FFT算法
#include "iostream.h"
#include "stdio.h"
#include "math.h"
#include "stdlib.h"
#include "malloc.h"
#define PI (float)3.1415926//复数结构体
typedef struct{ float re; float im;}complex;//定义旋转因子
complex W(int N,int n){
complex out; out.re=(float)cos(2*PI*((float)n/(float)N));
out.im=-(float)sin(2*PI*((float)n/(float)N)); return out;
}
//复数加法,out=a+b
complex add(complex a,complex b){
complex out;
out.re=a.re+b.re;
out.im=a.im+b.im;
return out;}
//复数减法,out=a-b
complex sub(complex a,complex b){
complex out;
out.re=a.re-b.re;
out.im=a.im-b.im;
return out;}
//复数乘法,out=a*b
complex mul(complex a,complex b){
complex out;
out.re=a.re*b.re-a.im*b.im;
out.im=a.re*b.im+a.im*b.re;
return out;
}
//复数赋值
complex comcpy(complex a){
complex out;
out.re=a.re;
out.im=a.im;
return out;
}
bool DOFFT(complex *x,complex *F,int N){
int K=(int)(log(N+1)/log(2));// cout<<K;
int i,j,m; for(m=1;m<=K;m++) {// cout<<"m="<<m<<endl;
for(i=0;i<(1<<K);i+=1<<m) {// cout<<"i="<<i<<"\t";
for(j=0;j<(1<<(m-1));j++){// cout<<i+j<<","<<i+j+(1<<(m-1))<<",";
F[i+j]=add(x[i+j],mul(x[i+j+(1<<(m-1))],W(1<<m,j)));
F[i+j+(1<<(m-1))]=sub(x[i+j],mul(x[i+j+(1<<(m-1))],W(1<<m,j)));
} // cout<<"\t";
}
// cout<<endl;
for(i=0;i<N;i++)
x[i]=comcpy(F[i]);
}
return true;
}
int reverse(int n,int j){
int i,k,power;
power=(int)(log(n+1)/log(2));
k=0;
for(i=0;i<power;i++)//检查每个为0~power-1
if(j&(1<<i))
k+=1<<(power-i-1);//如果正序数中位i的值为1,则倒序数中power-i-1的相应位置1。
return k;
}
void main(){
int m,i,j;
int N=4;
complex *x,*F;
x=new complex[N];
F=new complex[N];
for(i=0;i<N/2;i++){
x[i].re=1;
x[i].im=0;
cout<<"x["<<i<<"]="<<x[i].re<<", ";
}
for(i=N/2;i<N;i++){
x[i].re=0;
x[i].im=0;
cout<<"x["<<i<<"]="<<x[i].re<<", ";
}
cout<<endl;
for(i=0;i<N;i++){
F[i]=comcpy(x[reverse(N,i)]);
}
for(i=0;i<N;i++){
x[i]=comcpy(F[i]);
}
for(i=0;i<N;i++){
cout<<"x["<<i<<"]="<<x[i].re<<", ";
}
cout<<endl;
DOFFT(x,F,N);
for(i=0;i<N;i++)
{
cout<<"("<<F[i].re<<","<<F[i].im<<")";
}
cout<<endl<<endl;
}