数字信号处理课程设计FFT快速傅里叶变换程序设计

1、快速傅里叶变换程序设计中文摘要数字信号处理 (digital signal processing，dsp)是一门应用十分广泛的学科。数字信号处理是指利用计算机技术，以数字形式对信号进行采样、变换、滤波、估值、增强、压缩、识别等处理，以得到人们需要的信号形式。数字信号处理器也称dsp芯片，是一种用于进行数字信号处理运算的微处理器，其突出特点是采用多组总线技术实现并行机制，有独立的家发起和乘法器，灵活的寻址方式，实时快速地实现各种数字信号处理算法及各种复杂运算。傅里叶变换是一种将时域信号变换为频域信号的变换形式。在频域分析中，信号的频率及对应的幅值、相位（统称为频谱）反映了系统的性能。快速傅里叶变

2、换(fast fourier transform)是实现离散傅里叶变换(dft)的一种快速高效的运算方法，是数字信号处理中最为重要的算法之一。fft算法的关键在于利用了蝶形因子的内在对称性和周期性，从而加快了运算的速度，使运算时间缩短1至2个数量级。此次快速傅里叶变换的dsp程序设计，根据书上提供的c编程序在相对应的dsp工作环境下进行调试、图形分析、输入比较，来加深对程序的理解以及对所学内容的融合和巩固。程序设计所选择的头文件为，完成的是采样次数n为128的fft运算。为了更好的观察到快速傅里叶运算的结果，还另外设计了一组为方波的输入，从而得到设计要求的fft运算结果图形。关键词： 数字信号

3、处理，dsp芯片，快速傅里叶变换fft，c语言程序1 设计任务描述1.1 设计题目：快速傅里叶变换程序设计1.2 设计要求1.2.1 设计目的1）理解fft的算法以及利用dsp实现的方法。2）能熟练的调试程序并能观察其结果。3）熟悉tms320c54x系列dsp芯片的软件设计方法。1.3 基本要求 1）研究fft原理以及利用dsp实现的方法。2）编写fft程序。3）调试程序，观察结果。2 设计思路根据此次课程设计的要求，采用的是c语言程序。本次快速傅里叶变换程序设计主要包括三大部分：初始化定义部分，主函数部分，子程序部分。其中调用的子程序由四个功能函数组成：fft初始化函数、计算功率谱函数、波

4、形发生函数、倒序运算函数。这四个调用函数在主函数运行到相应位置时进入操作中，实现一个完整的快速傅里叶变换。由于在程序代码中调用了pow、log、cos、sin函数，该函数所在的c文件应包含头文件math.h。初始化定义部分主要是对程序需要用到的函数或者数据进行定义，这就需要我们熟悉前面知识所学到的各种基本数据类型的格式，如数组、结构、联合等构造类型数据。主函数部分就需要对fft变换的整体过程有熟悉的了解，这样才能够知道应该调用的子程序的顺序。因为一旦程序开始运行，除了一开始的初始化阶段，直接进入的就是主函数部分。函数调用部分是整个程序的重点部分，在这里实现了波形的输入、fft变换，可以通过cc

5、s软件调试该程序，并用其中viewgraphtime/frequency菜单功能，显示变量input与data图形，观察fft的效果。在此次的程序设计中，我设计了正弦波和方波两种输入情况。对于不同的输入，经过fft变换之后就会得到不同的频谱图，在整个程序中，fft的变换过程是最大的难点，完成这部分需要对数字信号处理课程中的快速傅里叶变换的知识有很大程度的掌握，这样才能知道函数中倒序产生的原理与方法，并且对蝶形图中的运算也能相当了解，这样才能理清思路，利用c语言的字符语句达到倒序处理的目的。3设计方框图设计方框图显现的是在设计时主要实现的功能和流程，简单易懂，能清晰的体现整个设计的思路。结束程序

6、初始化开始运行main()主函数调用initforfft初始化函数调用波形发生函数跳转到fft程序在128个采样值内否输入波形是否在取点值范围内频谱值输出图 3-1 主程序流程 4 快速傅里叶变换的的算法实现4.1变换原理 若给定由个信号样本（0），（1），（-1）组成的信号序列（），dft可用式2-1给出： =0,1，-1 （2-1）式2-1中，称为旋转因子或蝶形因子，=。从中可以看出：当信号样本为复数时，计算单个需经过次复数乘法和-1次复数加法运算，相当于4次实数乘法和2（2-1）次实数加法。完成全部点dft共需次复数乘法和（-1）复数加法运算。可见，随着不断增加，整个dft运算量是相当庞

7、大的，而fft算法通过对计算过程的深入分析，利用旋转因子具有的周期性与对称性，实现了降低运算复杂度的目的。当序列长度为偶数时，信号序列（）可被分解为奇、偶两个子序列，相应的点dft被分解为两个/2点的dft： =0,1, ，/2-1 （2-2） =0,1, ，/2-1 （2-3）式（2-2）和（2-3）中，和分别表示（）分解后得到的/2点偶序列点奇序列的dft。式（2-2）和式（2-3）表明，只要求出和，（）前/2点和后/2点的dft就得到了，整个序列的dft也就得到了。这样做的好处是计算点dft只需要约/2次复数乘法，总运算量约为直接dft运算量的一半。同理，当/2为偶数时，每个/2点的df

8、t又可被分解成两个/4点的dft，进一步减少了dft运算的复杂度。依次类推，直到不能继续分解为止。分解结束时，最小dft的点数称为称为基数，当=（为正整数）时，经过-1次分解，点dft最终可被分解为/2个两点的dft，即得到基数为2的fft运算，使得dft所需复数乘法次数降至。4.2蝶形运算 基于2fft的蝶形运算流图如下图4-1所示，图中以n=8为例。图4-1 蝶形运算采用的是原位运算原理和倒位序规律，其中原位运算指的是某一列的n个数据送到存储器后，经蝶形运算，其结果为下一列数据，它们以蝶形为单位仍存储在这同一组存储器中，直到最后输出，中间无需其他存储器。也就是蝶形的两个输出值仍放回蝶形的两

9、个输入所在的存储器中。 每列的n/2 个蝶形运算全部完成后， 再开始下一列的蝶形运算。而倒位序规律则是根据二进制的自然顺序和倒序得来的，其变换规律如图4-2。 自然顺序（i） 二进制数 倒位序二进制数 倒位序（j） 0123456700000101001110010111011100010001011000110101111104261537 图4-25 各部分程序设计5.1初始化程序 #includemath.h /数学函数头文件#define pi 3.1415926#define n 128 /采样次数n void initforfft(); /fft初始化函数void makewave

10、(); /波形发生函数void finv(int n1,float *xr,float *xi); /倒序运算函数，对输入序倒序int inputn,datan; float fwavern,fwavein,wn;float sin_tabn,cos_tabn; /正余弦函数表 int mum; /蝶形运算的级数这一段初始化程序是对需要用到数据和子程序段进行定义，其中包含了采样的次数、蝶形运算的级数，倒序运算正余弦函数表。5.2 fft的蝶形运算 在进行蝶形变换的时候，其中还引入了倒序运算函数，即运算完之后的fft的输出x(k)按正常顺序排列在存储单元中，即按x(0)，x(1)，x(7)的顺序

11、排列，但是这时输入x(n)却不是按自然顺序存储的，而是按x(0)，x(4)， ， x(7)的顺序存入存储单元，看起来好像是“混乱无序”的，实际上是根据倒序原理：倒序数的加1是在最高位加1，满2向次高位进1，最高位变0，依次往下，从当前的倒序值求出下一个倒序值。其蝶形变换的主要程序如下： for(k=j;k=n-1;k+=(int)(pow(2,m)+0.5) x=xrk+b*cos_tabs+xik+b*sin_tabs; y=xik+b*cos_tabs-xrk+b*sin_tabs; xrk+b=xrk -x; /运算结果的实部和虚部分别存储在原 xik+b=xik-y; /实部和虚部的位

12、置 xrk=xrk+x; xik=xik+y; 蝶形运算展开采用的是原位运算，即某一列的n个数据送到存储器后，经蝶形运算，其结果为下一列数据，它们以蝶形为单位仍存储在这同一组存储器中，直到最后输出，中间无需其他存储器。也就是蝶形的两个输出值仍放回蝶形的两个输入所在的存储器中。 每列的n/2 个蝶形运算全部完成后， 再开始下一列的蝶形运算。 这样存储器数据只需n个存储单元。 下一级的运算仍采用这种原位方式，只不过进入蝶形结的组合关系有所不同。 5.3 旋转因子 旋转因子是复数，可表示为： （4-1） 由式（4-1）可以看出旋转因子的实部为余弦函数，虚部为正弦函数。为了获得fft运算中需要的全部旋

13、转因子，需要分别存储正弦表和余弦表，且每个表长度为，对应于0180，同时，采用循环寻址方式对正弦表和余弦表进行寻址。5.4 实现n点的fft变换实现n点的fft变换需要知道采样的次数n值，且n为2的整数幂次方；同时还需要旋转因子的幂数，蝶形运算输入数据的距离，即各级旋转因子的个数；然后再进行蝶形运算，展开结果的实部和虚部分别存储在原实部和虚部位置。fft变换的程序如下： void fft(float xrn,float xin)int s,b;int m,j,k;float x,y;finv(n,xr,xi);for(m=1;m=mum;m+) b=(int)(pow(2,m-1)+0.5);

14、 for(j=0;jb;j+) s=j*(int)(pow(2,mum-m)+0.5); for(k=j;k=n-1;k+=(int)(pow(2,m)+0.5) x=xrk+b*cos_tabs+xik+b*sin_tabs; y=xik+b*cos_tabs-xrk+b*sin_tabs; xrk+b=xrk-x; xik+b=xik-y; xrk=xrk+x; xik=xik+y; for(m=0;m=n/2;m+) wm=sqrt(xrm*xrm+xim*xim); /功率谱计算在fft变换的后，经过功率谱计算公式得到频率域的频谱，根据输入信号的不同，得到不同的频谱图，当输入信号为正弦

15、波时，输出的data数据如图5-1，当输入的信号为方波时，输出的data数据序列见图5-2。 图5-1 图5-25.5 倒序运算部分的设计 倒序运算是在fft变换的蝶形运算开始之前调用的，其作用是不断的从当前的倒序值算出下一个倒序值，同时把输入序列的实部和虚部按顺序调换，是程序中需要仔细去理解的部分，其函数代码如下： void finv(int n1,float *xr,float *xi) int m,n,n2,k; float t; n2=n1/2; n=n2; for(m=1;m=n1-2;m+) if(m=k) n=n-k; k=(int)(k/2+0.5); n=n+k; 5.6正余

16、弦函数表与波形发生部分的设计 正余弦函数表的建立的目的是为蝶形变换公式中的正余弦序列提供原始数据依据，而波形发生器的作用则是为整个程序提供一个输入序列，使系统根据输入的波形来得到与之相对应的频谱。正余弦函数表的程序代码为： void initforfft() int i; for(i=0;in;i+) sin_tabi=sin(pi*2*i/n); cos_tabi=cos(pi*2*i/n); 至于波形发生的输入，我用了正弦波和方波两种波形作为互相对比的输入，正弦波为f=3hz的序列： void makewave() int i; for(i=0;in;i+) inputi=sin(pi*2

17、*i/n*3)*1024; 与之相对应的正弦波序列的地址为0x3f90c0，波形图如图5-3所示。 图5-3若以方波作为输入，除了多定义两个变量之外则是把for循环里的语句变为j=(int)i/20; l=j%2; if(l=0)inputi=1024; else inputi=0;即可。其方波的波形图见图5-4 图5-45.7主函数部分的设计 主函数main()是每一个c语言程序中有且仅有的一个重要部分，是一个程序的入口点，可以根据自己的爱好把它放在程序的某个地方。有些程序员把它放在最前面，而另一些程序员把它放在最后面, 无论放在哪个地方。本次程序的主函数main()为： main()int

18、 i;initforfft();makewave();for(i=0;in;i+) fwaveri=inputi; fwaveii=0.0;wi=0.0; mum=(int)(0.5+log(n)/log(2);fft(fwaver,fwavei); for(i=0;in;i+)datai=wi;while(1);5 工作过程分析 根据此次课程设计的要求，为了实现快速傅里叶变换的目的，首先要做的就是让程序有相应的输入信号，为此，在程序中，我用了正弦波和方波两种波形作为输入以形成对比。再然后就是调用程序各部分分支以实现相应功能，最终达到fft变换的目的。首先，需要对程序用到的函数进行定义，这部分

19、就是程序是初始化阶段。由于傅里叶变换是一种数学运算，所以对文件类型进行定义为数学函数头文件，再来就是程序用了不少的函数，就定义了波形发生函数、倒序运算函数、fft初始化函数和正余弦函数表等变量。然后，进行main()主函数部分，这是程序的主干、入口，主函数内调用子程序的顺序决定了各部分子程序运行的顺序，所以这部分是很重要的，要根据实际fft变换的顺序来编写。再来，就是对各个子程序的编写，这些程序分支就包含fft初始化函数，顺利建立正余弦函数表，倒序运算函数，波形发生及蝶形运算展开并计算功率谱。最后，就是根据输入输出程序代码段的地址，通过ccs软件调试该程序，并用其中viewgraphtime/

20、frequency菜单功能，改变相应的起始地址和显示方式显示输入变量input与输出结果data图形，观察fft的效果。 其中，输入变量input的起始地址为0x3f90c0，输出结果data的地址为0x3f9040，这样就能显示输入输出图形，顺利观察fft的效果。小 结通过本次的dsp课程设计,让我对这个工作原理以及它编程的调试方法都类似于单片机的学科有了更好的掌握与理解。与次同时，通过此次的课程设计也让我充分地认识到了我在学习这门课程中存在的不足之处。当刚知道有dsp这门学科的程序设计时，首先闯入我脑中的信息便是这次的课设肯定不好做，但是，在真正的知道了我要设计的题目在树上能找到类似的程序

21、时，我就知道我把问题给想复杂了。在整个的设计过程中，我需要做的就是把书上的程序看明白，然后再对其中的某些部分进行小小的修改，算起来，我做的设计算是三组里最简单的一个了，程序调试很快，波形除了一些瑕疵之外算得上是不错了。这其中出现的主要问题便是对主函数调用的子程序内容，特别是蝶形运算那部分，是这次设计的瓶颈 ，经过几天的研究与讨论，让我对以前学的懵懵懂懂的fft变换有了更深层次的了解，特别是蝶形图展开及倒序运算部分。 通过这次的课程设计让我学到了很多平时课上学不到的东西，让我知道了哪些方面是自己的薄弱环节，自己需要在哪些方面下一些功夫。更重要的是能过这次的实验我熟练地掌握了程序的调试，知道如果一

22、个程序运进的时候出现了问题应该怎样找到错误并改正，也知道了该如何去利用图书馆和网络上的资源，这对我们以后的课程设计和学习是有很大帮助的。总之课程设计对我们的帮助是非常大的，在提高我们的动手能力和开发逻辑思维之于也能在一定的程度上激发我们对dsp这门课学习的兴趣。让我们的学习从单一的理论升华到实践的高度。 很庆幸能有这样一次课设的机会，让我学尝试了理论与实际的结合，让我学到了很多知识，让我认识到了自己的不足，知道光是掌握了理论知识并不代表什么，必须还要与实际结合才能更好的理解这些知道，使之掌握得更牢靠。在以后的学习生涯中我会不断的完善自己，用更多知识来丰富和充实自己，为自己以后的人生道路打下坚实

23、的基础！致谢短短一周的dsp课程设计就这样结束了，很感谢吕老师和雷老师在这一周里给我带来的帮助和指导。还要感谢在这一周时间里与我共同讨论和交流的同学，感谢他们的细心指导和他们提出的疑点。吕老师在软件方面很厉害，只要遇到跟程序里遇到的问题就可以找他，老师肯定会尽责为我们讲解，而雷老师因为本来就是数字信号处理的任课老师，在傅里叶变换这个设计中，必然会有很多请教雷老师的时候，他都会细心地为我们一一解惑，知道我们听懂为止。还要感谢为我提供课设资料的学校图书馆，图书馆的大量书籍为我提供了广阔的资源，让我熟悉掌握了定时器、外部中断的工作原理，由其是液晶显示器的原理开启方法，对我的帮助非常大。仅仅运用课堂所

24、学的知识，靠自己的力量要在两周之内完成设计的题目是难以想象的。所以这两周我能够顺利地完成课设任务有很大一部分原因是大家的帮助和图书的资源。再次，我要感谢指导老师和与我同一个课设题目的同学，因为一个人的思路是有限的，感谢大家的集思成就了我们的妙想。最后感谢学校给我们提供这种独立自主研究性学习的机会，充分开发了我们的创新能力。每个学期的课程设计都让我学会了很多，让我对学习充满浓厚的兴趣! 参考文献1 姜沫岐，许涵，俞鹏，段国强 . dsp原理与应用从入门到提高. 北京: 机械工业出版社， 20072 李行一. 数字信号处理.重庆大学出版社，20023 张东亮. dsp控制器原理与应用.机械工业出版

25、社，20114 李哲英. dsp基本理论与应用系统设计.北京航空航天大学出版社，20025 俞卞章. 数字信号处理. 西安: 西北工业大学出版社， 20026 王安国. 数字信号处理基础.北京: 电子工业出版社， 20037 何苏秦，王忠勇. tms320c2000dsp原理及实用技术.电子工业出版社，2003附录a1 程序清单#includemath.h#define pi 3.1415926#define n 128void initforfft();void makewave();void finv(int n1,float *xr,float *xi);int inputn,datan;float fwavern,fwavein,wn;float sin_tabn,cos_tabn;int mum;void fft(float xrn,float xin)int s,b;int m,j,k;float x,y;finv(n,xr,xi);for(m=1;m=mum;m+)b=(int)(pow(2,m-1)+0.5); for(j=0;jb;j+) s=j*(int)(pow(2,mum-m)+0.5); for(k=j;k=n-1;k+=(int)(pow(2,m)+0.