请教一下大家,谁能告诉我64点fft的概念和c程序
的有关信息介绍如下:fft的概念最好看书,这个细细分析一下还是能够理解的。你要是没有相关的书的话,这个网站上的也可以看看:
c程序我以前写过一个,注释得还是挺清楚的,你看看吧
(我这里写的是1024个点,你简单修改为64个点就行。程序思路比较简单清晰,建议你自己看懂吸收了比较好,这才是学习)
/**************************************************************/
#include
#include
#ifndef pi
#define pi 3.14159265
#endif
int N=1024; //采集的数据点数,注意为2的n次方
int n=10; //2^n=N
double input; //原始输入信号
double IO; //经过重新排序的信号(实数)
int code; //0~(2^n-1)的倒置码自然数组
struct complex c_IO ; //经过重新排序的信号(复数)
struct complex out ; //输出变换的结果(复数)
struct complex temp ; //保存中间变换结果(复数)
void fft(struct complex *c_IO,int m); //fft函数,输入c_IO[]为经过了码位转换的信号数据;m=n-1表示需要分解的层数
void trans_code(int n ); //产生0~(2^n-1)的倒置码自然数组
struct complex complex_mult(struct complex a, struct complex b); //结构体复数相“乘”函数,输入为两个结构体类型的复数
struct complex complex_add(struct complex a, struct complex b); //结构体复数相“加”函数,输入为两个结构体类型的复数
struct complex complex_remove(struct complex a, struct complex b); //结构体复数相“减”函数,输入为两个结构体类型的复数
struct complex W_n_k (int n,int k); //产生权值的函数W(n,k)=exp(-j*2*k*pi/N)
struct complex doule2complex (double input); //把实数转换为结构体复数的函数
void main()
{
FILE *fp;
///////////////////////////////// 获取输入信号,这里是手动(程序)输入,也可以简单修改为 “读入文本文件”
int k=0;
for(k=0;k
{
input[k]=sin(pi*k/3)+cos(pi*k/20);
}
////////////////////////////////// 把输入数据进行码位倒置
for(k=0;k
{
code[k]=0;
}
trans_code(n);
for(k=0;k
{
IO[k]=input[ code[k] ];
}
///////////////////////////////// 把输入数据转化为复数,为下面的FFT做准备
for(k=0;k
{
c_IO[k]=doule2complex(IO[k]);
}
///////////////////////////////// 进行FFT变换,结果保存在全局变量 out[N]中
fft(c_IO, n-1);
for(k=0;k
{
printf("%f +%f j\n",out[k].x,out[k].y);
}
////////////////////////////////////////////////把结果输入到文本文件中去,其中实部输入到dx.txt文件中,虚部输入到dy.txt中
fp=fopen("dx.txt","w");
for(k=0;k
{
fprintf(fp,"%f\n ",out[k].x);
}
fclose(fp);
fp=fopen("dy.txt","w");
for(k=0;k
{
fprintf(fp,"%f\n ",out[k].y);
}
fclose(fp);
///////////////////////////////////////////////////
} //main函数截至
////////////////////////////////////////////
////////////////////////////////////////////
// 以下是子函数的实现部分 //
////////////////////////////////////////////
////////////////////////////////////////////
void fft(struct complex c_IO[],int m)
{
int group,i,p,q,k;
if(m==-1)
{
for(k=0;k
{
out[k]=c_IO[k];
}
}
else
{
fft(c_IO, m-1); //递归调用
group=(int)( N/pow(2,m+1) ); //每一级的分组数
for(i=0;i
{
for(p=(int)( i*pow(2,m+1) ) ; p<=(int)( i*pow(2,m+1)+pow(2,m)-1 ) ; p++ )
{
q=p+(int)pow(2,m);
temp[p]=complex_add(out[p],complex_mult(W_n_k( (int)pow(2,m+1),(int)(p-i*pow(2,m+1)) ),out[q] ) );
temp[q]=complex_remove(out[p],complex_mult(W_n_k( (int)pow(2,m+1),(int)( p-i*pow(2,m+1) ) ),out[q]) );
}
}
for(k=0;k
{
out[k]=temp[k];
}
}
}
void trans_code(int n ) //change code order from 1 to 2^n
{
int bit={0};
int i=0;
int j=0;
int k=0;
int mytemp=0;
int num;
num=(int)pow(2,n);
for(j=0;j
{
mytemp=j;
////////////////////////////////////////////////////////////
for(i=0;i
{
bit[i]=mytemp%2;
mytemp=mytemp/2;
}
////////////////////////////////////////////////////////////
for(k=0;k<(n/2);k++) //码位倒置
{
mytemp=bit[k];
bit[k]=bit[n-1-k];
bit[n-1-k]=mytemp;
}
///////////////////////////////////////////////////////////
for(i=0;i
{
code[j]=code[j]+(int)pow(2,i)*bit[i];
}
}
}
struct complex complex_mult(struct complex a, struct complex b)
{
struct complex c;
c.x=a.x*b.x-a.y*b.y;
c.y=a.x*b.y+a.y*b.x;
return c;
}
struct complex complex_add(struct complex a, struct complex b)
{
struct complex c;
c.x=a.x+b.x;
c.y=b.y+a.y;
return c;
}
struct complex complex_remove(struct complex a, struct complex b)
{
struct complex c;
c.x=a.x-b.x;
c.y=a.y-b.y;
return c;
}
struct complex W_n_k (int n,int k)
{
struct complex c;
c.x=cos(2*pi/n*k);
c.y=-sin(2*pi/n*k);
return c;
}
struct complex doule2complex (double input)
{
int i=0;
struct complex c;
c.x=input;
c.y=0;
return c;
}