用c语言实现的FFT.docx
- 文档编号:23269641
- 上传时间:2023-05-15
- 格式:DOCX
- 页数:10
- 大小:312.92KB
用c语言实现的FFT.docx
《用c语言实现的FFT.docx》由会员分享,可在线阅读,更多相关《用c语言实现的FFT.docx(10页珍藏版)》请在冰豆网上搜索。
用c语言实现的FFT
用c语言实现的FFT
一、对FFT的介绍
1.FFT(FastFourierTransformation),即为快速傅里叶变换,是离散傅里叶变换的快速算法,它是根据离散傅里叶变换的奇、偶、虚、实等特性,对离散傅里叶变换的算法进行改进获得的。
2.FFT算法的基本原理
FFT算法是把长序列的DFT逐次分解为较短序列的DFT。
按照抽取方式的不同可分为DIT-FFT(按时间抽取)和DIF-FFT(按频率抽取)算法。
按蝶形运算的构成不同可分为基2,基4,基8,以及任意因子的类型。
3.迭代关系
由
按照上述公式的规律进行逐级分解,直到2点DFT,如下是N=8时的蝶形算法分析图:
四、FFT算法中蝶形算法的基本思想分析
(1)我们知道N点FFT运算可以分成log2(N)级,每一级都有N/2个碟形,FFT的基本思想是用3层循环完成全部运算(N点FFT)。
(2)第一层循环:
由于N=2^m需要m级计算,第一层循环对运算的级数进行控制。
(stages)
(3)第二层循环:
由于第L级有2^(L-1)个蝶形因子(乘数),第二层循环根据乘数进行控制,保证对于每一个蝶形因子第三层循环要执行一次,这样,第三层循环在第二层循环控制下,每一级要进行2^(L-1)次循环计算.(选择W)
(4)第三层循环:
由于第L级共有N/2^L即2^(n-L)个群,并且同一级内不同群的乘数分布相同,当第二层循环确定某一乘数后,第三层循环要将本级中每个群中具有这一乘数的蝶形计算一次,即第三层循环每执行完一次要进行N/2^L个碟形计算。
(执行不同group中具有相同W的蝶形运算)
(5)可以得出结论:
在每一级中,第三层循环完成N/2^L个碟形计算;第二层循环使第三层循环进行2^(L-1)次,因此,第二层循环完成时,共进行2^(L-1)*N/2^L=N/2个碟形计算。
实质是:
第二、第三层循环完成了第L级的计算。
五、用c语言实现的FFT算法如下:
18px;">#include #include #include #defineN1000 /*定义复数类型*/ typedefstruct{ doublereal; doubleimg; }complex; complexx[N],*W;/*输入序列,变换核*/ intsize_x=0;/*输入序列的大小,在本程序中仅限2的次幂*/ doublePI;/*圆周率*/ voidfft();/*快速傅里叶变换*/ voidinitW();/*初始化变换核*/ voidchange();/*变址*/ voidadd(complex,complex,complex*);/*复数加法*/ voidmul(complex,complex,complex*);/*复数乘法*/ voidsub(complex,complex,complex*);/*复数减法*/ voidoutput();/*输出快速傅里叶变换的结果*/ intmain() { inti;/*输出结果*/ system("cls"); PI=atan (1)*4; printf("输出DIT方法实现的FFT结果\n"); printf("Pleaseinputthesizeofx: \n");//输入序列的大小 scanf("%d",&size_x); printf("Pleaseinputthedatainx[N]: \n");//输入序列的实部和虚部 for(i=0;i { printf("请输入第%d个序列: ",i); scanf("%lf%lf",&x[i].real,&x[i].img); } printf("输出倒序后的序列\n"); initW();//调用变换核 fft();//调用快速傅里叶变换 printf("输出FFT后的结果\n"); output();//调用输出傅里叶变换结果函数 return0; } /*快速傅里叶变换*/ voidfft() { inti=0,j=0,k=0,l=0; complexup,down,product; change();//调用变址函数 for(i=0;i (2);i++)/*一级蝶形运算stage*/ { l=1< for(j=0;j { for(k=0;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; } } } } /*初始化变换核,定义一个变换核,相当于旋转因子WAP*/ voidinitW() { inti; W=(complex*)malloc(sizeof(complex)*size_x);//生成变换核 for(i=0;i { W[i].real=cos(2*PI/size_x*i);//用欧拉公式计算旋转因子 W[i].img=-1*sin(2*PI/size_x*i); } } /*变址计算,将x(n)码位倒置*/ voidchange() { complextemp; unsignedshorti=0,j=0,k=0; doublet; for(i=0;i { k=i;j=0; t=(log(size_x)/log (2)); while((t--)>0)//利用按位与以及循环实现码位颠倒 { j=j<<1; j|=(k&1); k=k>>1; } if(j>i)//将x(n)的码位互换 { temp=x[i]; x[i]=x[j]; x[j]=temp; } } output(); } /*输出傅里叶变换的结果*/ voidoutput() { inti; printf("Theresultareasfollows: \n"); for(i=0;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; } 六、FFT原理的理解和程序设计中遇到的相关问题及解决方法 1、遇到的相关问题: (1)首先一开始不知道什么是FFT,以及FFT原理是什么 (2)不理解FFT中迭代关系的推导以及缘由 (3)不理解变址计算的原理 (4)对蝶形运算的推导与原理的理解不透彻 (5)编程过程中对变址计算即对按位与的变换形式的不理解 2、解决的方法: (1)到图书馆借相关的书籍理解相关的原理和过程 (2)到XX收索相关的资料促进理解相应的原理过程 (3)向学长学姐请教或老师的指导 七、总结 从这次的考核中我学到了一些有关数字信号处理的相关知识,即快速傅里叶变换的原理,虽然并不是非常深入的去学习,但却深刻的领悟到了,弄懂最基本的原理是理解的关键,只有弄懂了最基本的原理并具备一定的编程语言基础,才能顺利完成一个项目。 在这次的学习过程中,我对原理的理解是通过对大量书籍和上网查阅资料以及老师的指导才渐渐的理解的。 所以我觉得做好一件事,需要有坚持不懈的努力和一点一滴的积累。
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 语言 实现 FFT
![提示](https://static.bdocx.com/images/bang_tan.gif)