《课程设计（论文）-快速傅里叶变换程序设计.doc》由会员分享，可在线阅读，更多相关《课程设计（论文）-快速傅里叶变换程序设计.doc（16页珍藏版）》请在三一文库上搜索。
1、快速傅里叶变换程序设计1 设计任务描述1.1 设计题目 快速傅里叶变换程序设计1.2 设计要求1.2.1 设计目的1）理解FFT的算法以及利用DSP实现的方法。2）能熟练的调试程序并能观察其结果。3）熟悉TMS320C54x系列DSP芯片的软件设计方法。1.3 基本要求1）研究FFT原理以及利用DSP实现的方法。2）编写FFT程序。3）调试程序，观察结果。2 设计思路2.1 FFT算法原理若给定由个信号样本（0），（1），（-1）组成的信号序列（），DFT可用式2-1给出： =0,1，-1 （2-1）式2-1中，称为旋转因子或蝶形因子，=。从中可以看出：当信号样本为复数时，计算单个需经过次复数
2、乘法和-1次复数加法运算，相当于4次实数乘法和2（2-1）次实数加法。完成全部点DFT共需次复数乘法和（-1）复数加法运算。可见，随着不断增加，整个DFT运算量是相当庞大的，而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点的
3、DFT就得到了，整个序列的DFT也就得到了。这样做的好处是计算点DFT只需要约/2次复数乘法，总运算量约为直接DFT运算量的一半同理，当/2为偶数时，每个/2点的DFT又可被分解成两个/4点的DFT，进一步减少了DFT运算的复杂度。依次类推，直到不能继续分解为止。分解结束时，最小DFT的点数称为称为基数，当=（为正整数）时，经过-1次分解，点DFT最终可被分解为/2个两点的DFT，即得到基数为2的FFT运算，使得DFT所需复数乘法次数降至。2.2 基2FFT的蝶形运算流图基2FFT的蝶形运算过程可用图2-1所示，此时=8，=3。图2-1 8点基2FFT运算过程观察图2-1，根据DFT的基2FF
4、T算法，可以总结出以下几条规律：（1）点FFT运算从输入端开始，逐级进行，共需经过级运算；在第（=1，2，）级中存在个相似的蝶形运算组（除输入数据不同外）；每个组内蝶形运算的个数为，参与每个蝶形运算的两个输入数据相距个点。（2）中间数据的存储，可采用原位存储法。即每次蝶形运算的结果存储在与原数据相同的内存单元内。（3）为了保证输出数据按自然数序排列，在进行FFT之前输入数据需要按照特定的顺序存放，通过位倒序寻址可以满足这种要求。2.2.1 输入、输出序列的倒位序规律由流程图2-1可以看出，当进行原位运算时，发现当运算完成后，FFT的输出按自然顺序排列在存储单元中，即按，的顺序排列；但是这是输入
5、却不是按自然顺序存储的，而是按，的顺序存入存储单元，这种方式就称之为倒位序。当用二进制表示这个顺序时，它正好是“码位倒置”的顺序。例如，原来的自然顺序应是的地方，现在存放，用二进制码表示这一规律时，则是在处存放，出存放，即将自然顺序的二进制码位倒置过来，第一位码变成最末位码，这样倒置以后的顺序正是输入所需要的顺序。表2-1中列出的是=8时按码位倒置规律所得的顺序，其结果与按时间抽取算法FFT流程图中的输入顺序是一致的。同理，当=256时，亦可以采用同样的方法进行位倒序操作。表2-1 码位倒置顺序自然顺序二进码表示码位倒置倒位序0000000010011004201001023011110641
6、0000115101101561100113711111172.2.2 蝶距的计算设=，则整个运算流图中包含级蝶形运算，每一级则有个蝶形单元。蝶距即每个蝶形单元两个输入（出）节点的序号差。以为例，结合图2-1可知共包含3级蝶形运算，每一级蝶形单元的蝶距如下：第一级，蝶距为1，可以看作由得到：第二级，蝶距为2，可以看作由得到；第三级，蝶距为4，可以看作由得到。因此得：第级蝶形单元的蝶距为：。2.2.3 旋转因子的计算由FFT算法原理过程可知，若=，则共有级蝶形运算，各级蝶形运算中旋转因子分别如下：第级的旋转因子为（=0,1，）；第-1级的旋转因子为（=0,1，）；第一级的旋转因子为（=0,1，）
7、。由此可见， 第级蝶形运算中旋转因子为，=0，1，。2.3 FFT算法的DSP的实现方法设FFT运算的输入数据为实数，则2点实数FFT算法的实现步骤为：第一步，把2实数输入序列组合成点的复数序列。然后把该复数序列进行位倒序操作后存储在输入区。第二步，进行的FFT运算。第三步，把点FFT输出拆成2的复数序列，这2的复数序列对应于2点时实数输入序列的DFT输出。第四步，结果输出及功率谱计算。3 软件流程图结束程序初始化开始送入数据调入系数表输入数据位码倒置FFT的蝶形运算是否发生溢出出？归一化输入数据结束？各图形输出4 各部分程序设计及参数计算4.1 初始化程序#includemath.h#def
8、ine 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;初始化程序是对采样点数N、各个子程序、数据实部和虚部、输入输出以及正余弦表进行定义，在后面的程序中运用时就不用再进行定义了。4.2 子程序4.2.1 旋转因子、蝶距运算void FFT(float XrN,float XiN) int S,B;int
9、 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.5)该子程序是时间抽取法FFT程序，要求采样点数为2的整数幂次方，Xr和Xi分别为输入序列的实部和虚部。S为旋转因子的幂数，B为蝶形图运算输入数据的距离，也即各级旋转因子的个数。finv(N,Xr,Xi)语句是实现倒序运算函数，对输入序列倒序。B=(int)(pow(2,m-1)+0.5)语句计算出蝶距
10、 B=2(m-1)，S=j*(int)pow(2,Mum-m)+0.5语句计算出旋转因子的幂数S=2(Num-1)。4.2.2 蝶形图变换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; 该段子程序是实现蝶形图的变换，蝶形运算展开，结果的实部和虚部分别存储在原实部和虚部位置。4.2.3 波形发生函数void MakeWave() int i;int j;int l; for(i=0;iN;i+) j=(int)i/20; l=
11、j%2; if(l=0)INPUTi=1024; else INPUTi=0; 该程序是实现输入波形为方波，从0点开始20个点一个高电平，20个点一个低电平。4.2.4 倒序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; 该程序段是倒序运算函数finv(N1,Xr,Xi),对输入序列倒序。N1为序列长度；Xr,Xi分别为输入序列的实部和虚部。倒序原理：倒序列的加1实在最高位加
12、1，满2向次高位进1，最高位变0，依次往下，从当前倒序值可求得下一倒序值。4.3 主程序main()int 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);主程序的功能是调用各个子程序，Mum=(int)(0.5+log(N)/log(2)语句中Mum为蝶形运算的级数，N=2Num。Num=7。 小结在这次设计中，我完成了这次任务
13、，在一周的课设中，我弄懂了快速傅里叶变换程序设计，并且调试出该程序的输入输出波形以及功率谱。在这次课程设计中,我学到了很多上课时没有学到的知识。这次设计的程序是关于快速傅里叶变换的，在设计程序之前得先弄明白快速傅里叶变换的理论知识，我翻阅了上学期学的数字信号处理，弄明白了该程序要实现的功能。然后把整个程序分成好几个模块，先读懂每个模块，再把程序结合在一起，最后终于弄懂了整个程序。在程序的调试过程中我遇到了许多困难。首先，在编译时有两个错误是认为因素造成的，有一个是因为把指令的拼音打错了，还有一个是因为在程序的最开始多留出一行，由于DSP的汇编指令要求非常严格，所以在编写程序时应该注意书写的格式
14、。通过课程设计我收获很多，不仅对这个课程有了更深的理解，而且也学会了团队精神的重要性，个人的能力是有限的，团结才能有力量，我们都尽自己所能来完成这次课程设计。通过这次课程设计,我懂得了理论与实际结合是很必要的。只有理论知识是远远不够的，只有把所学的理论知识与实际结合起来，从理论中得出结论，才能真正为社会服务，从而提高自己的实际动手能力和独立思考能力。重要的是我们在设计过程中发现了自己的不足，对前面学过的知识理解的不够深刻，掌握的不够牢固。通过这次课程设计，我们相当于把前面所学过的知识又重新温故了一遍。致谢课程设计结束了，在课程设计的过程中我学到了很多课堂上所不能体会的知识。我们做的是快速傅里叶
15、变换程序设计，在这次设计当中，我学会了团体合作还有很多关于dsp的相关知识。一周的课程设计结束了，以前在上课的时候老师讲解的东西以为明白就可以，但是在这次实践中，我明白了，那是远远不够的。 在老师、同学的帮助下，终于完成了这次实验。查阅资料，请教同学，相互讨论，还好有老师和同学，一步一步的修改和讲解让我觉得非常的暖心。可能在以后步入社会以后这种感觉可能不复存在。所以我珍惜现在的每一次的设计。在这里要非常感谢我的老师。每天都来实验室为同学们解答各种问题，当我们出错误的时候，老师都会在身边所以我们感到非常的安心。这也是一次跟老师更加近距离的学习机会。在这里还要感谢我同组的同学以及帮助我们的同学。在
16、每天的课程设计中，都有你们的陪伴，在我有困难的时候，你们永远在我的身边。非常谢谢你们对我的帮助。这次课程设计对我以后的学习乃至工作都有很大的帮助，我对dsp有了更深一层次的了解。虽然对程序设计还是有很多困难，很多要学的。但我想只要我们珍惜每次设计的经过，从中我们肯定会有一些所得。再次感谢我们的指导老师，还有曾经一直帮助我的同学。参考文献1 俞一彪，孙兵. 数字信号处理理论与应用.南京：东南大学出版社，20052 李行一. 数字信号处理.重庆：重庆大学出版社，20023 郭森楙，颜允圣. 数字信号处理器体系结构、实现与应用.北京：清华大学出版社，20054姜沫岐，许涵，俞鹏，段国强 . DSP原
17、理与应用从入门到提高. 北京: 机械工业出版社， 20075 俞卞章. 数字信号处理. 西安: 西北工业大学出版社， 20026 王安国. 数字信号处理基础.北京: 电子工业出版社， 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;vo
18、id 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.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; fo
19、r(m=0;m=N/2;m+) wm=sqrt(Xrm*Xrm+Xim*Xim);main()int 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);void InitForFFT() int i; for(i=0;iN;i+) sin_tabi=sin(PI*2*i/N); cos_tabi=cos(PI*2*i/N); void MakeWave() int i;int j;int l; for(i=0;iN;i+) j=(int)i/20; l=j%2; if(l=0)INPUTi=1024; else INPUTi=0; 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; 附录A2 程序图形图一 输入信号时域波形图二 输入信号频域波形图三 输出信号功率谱16