|
FFT 代码
#include "iostream.h"
#include "math.h"
#include "memory.h"
#include "conio.h"
#define pi 3.1415926
typedef struct tagNumber
{
double re;
double im;
}Number;
Number Add(Number a,Number b)
{
Number c={a.re+b.re,a.im+b.im};
return c;
}
Number Sub(Number a,Number b)
{
Number c={a.re-b.re,a.im-b.im};
return c;
}
Number Mul(Number a,Number b)
{
Number c={a.re*b.re-a.im*b.im,a.re*b.im+a.im*b.re};
return c;
}
void fft(Number* t,Number* f,int r)
{
long count;//傅立叶变换点数
int i,j,k,p,bfsize;
Number *w,*x,*a,*b;//复数结构类型的指针变量,其中w指向加权系数
double angle;//计算加权系数所用的角度值
count=1<<r;//计算傅立叶变换点数
//分配所需运算空间
w=new Number[count/2];
a=new Number[count];
b=new Number[count];
//计算加权系数
for(i=0;i<count/2;i++)
{
angle=-i*pi*2/count;
w.re=cos(angle);
w.im=sin(angle);
}
memcpy(a,t,sizeof(Number)*count);
//采用频率分解法进行蝶形运算
for(k=0;k<r;k++)
{
for(j=0;j<1<<k;j++)
{
bfsize=1<<(r-k);
for(i=0;i<bfsize/2;i++)
{
p=j*bfsize;
b[i+p]=Add(a[i+p],a[i+p+bfsize/2]);
b[i+p+bfsize/2]=Mul(Sub(a[i+p],a[i+p+bfsize/2]),w[i*(1<<k)]);
}
}
x=a;
a=b;
b=x;
}
//将乱序的变换序列重新排序
for(j=0;j<count;j++)
{
p=0;
for(i=0;i<r;i++)
{
if(j&(1<<i))
p+=1<<(r-i-1);
}
f[j]=a[p];
}
//释放存储器空间
delete w;
delete a;
delete b;
}
int main()
{
Number t[8]={{1,1},{1,1},{1,0},{1,0},{1,1},{1,1},{1,1},{1,1}};
Number *p=t;
Number* f;
f=new Number[8];
cout<<"The data sequence in time domain:"<<endl<<endl;
for(int i=0;i<8;i++)
cout<<(*(p+i)).re<<"+("<<(*(p+i)).im<<"i) ";
cout<<endl<<endl<<endl<<endl;
cout<<"The data sequence in frequency domain after fft:"<<endl<<endl;
fft(t,f,3);
for(i=0;i<8;i++)
cout<<(*(f+i)).re<<"+("<<(*(f+i)).im<<"i) ";
cout<<endl;
delete []f;
getch();
return 0;
} |
评分
-
1
查看全部评分
-
|