课程设计（论文）基于matlab的fft算法设计

1、目目 录录1 引言.12 课程设计内容及要求.22.1 课程设计内容 .22.2 课程设计要求 .22.3 课程设计目的 .22.3 课程设计平台.23 基于 matlab 的 fft 算法设计原理.33.1 总体设计流程图 .33.2 语音信号的采集 .33.3 语音信号的时频分析 .33.4 快速傅里叶变换 .63.4.1 fft的运算规律.83.4.2 基于 matlab 的 fft 所编写程序的框图 .123.5 自编算法与机带算法仿真波形比较 .134 设计总结.16参考文献.17附录.18课程设计说明书课程设计说明书11 引言随着信息时代，数字时代的到来，数字信号处理已经成为一门极

2、其重要的学科和技术领域。以 dsp 为核心芯片的处理系统日益变成了数字信号处理系统的主流。它广泛用于电子信息、通信、图像处理、语音处理、生物医学、自动控制、地质探测等领域，受到工程设计和使用人员的青睐。matlab，它是美国 math works 公司推出的一种面向工程和科学计算的交互式计算软件。它以矩阵运算为基础，把计算、可视化、程序设计融合在一个简单易用的交互式工作环境中，是一款数据分析和处理功能都非常强大的工程适用软件。通过本次实习我们学会了分析和处理音频信号，首先要对声音信号进行采集，matlab 的数据采集工具箱提供了一整套命令和函数，通过调用这些函数和命令，可直接控制声卡进行数据采

3、集。window 自带的录音机程序也可驱动声卡来采集语音信号，并能保存为 wav 格式文件，供 matlab 相关函数直接读取、写入或播放。matlab 语言是一种数据分析和处理功能十分强大的计算机应用软件，它可以将声音文件变换位离散的数据文件，然后利用其强大的矩阵运算能力处理数据，如数据滤波、傅立叶变换、时域和频域分析、声音回放以及各种图的呈现等，它的信号处理与分析工具箱位语音信号分析提供了十分丰富的功能函数，利用这些功能函数可以快捷而又方便的完成语音信号的处理和分析以及信号的可视化，是人机交互更加便捷。信号处理是 matlab 重要应用的领域之一。语音信号处理是研究用数字信号处理技术和语音

4、学知识对语音信号进行处理的新兴的学科，是目前发展最为迅速的信息科学研究领域的核心技术之一。通过语音传递信息是人类最重要、最有效、最常用和最方便的交换信息形式。语音信号的处理与滤波的设计主要是用 matlab 作为工具平台，设计中涉及到声音的录制、播放、存储和读取，语音信号的抽样、频谱分析，滤波器的设计及语音信号的滤波，通过数字信号处理课程的理论知识的综合运用。从实践上初步实现对数字信号的处理。 课程设计说明书课程设计说明书22 课程设计内容及要求2.1 课程设计内容录制一段个人自己的语音信号，并对录制的信号进行采样；画出采样后语音信号的时域波形和频谱图；在 matlab 环境下编写基 2 di

5、t-fft 算法；利用自己编写的算法对已采集的语音信号进行频谱分析，并画出语音信号的时域与频谱图，并与 matlab 数字信号处理工具箱中的 fft 函数进行对比研究，验证自编算法的正确性。2.2 课程设计要求1.完成语音信号的采集，利用 windows 自带的录音机或其他软件，录制一段语音，时间在 1s 以内；2.在 matlab 中编写程序，实现输入信号的倒序；3.编写程序，实现蝶形运算；4.画出语音信号的频谱图，与 matlab 数字信号处理工具箱中的 fft 函数进行对比研究，并对设计结果进行独立思考和分析；2.3 课程设计目的1.学会 matlab 的使用，掌握 matlab 的程序

6、设计方法。 2.掌握在 windows 环境下语音信号采集的方法。3.掌握数字信号处理的基本概念、基本理论和基本方法。 4.掌握 matlab 设计 fir 和 iir 数字滤波器的方法。5.学会用 matlab 对信号进行分析和处理。2.4 课程设计平台matlab7.1 软件matlab 是由美国 mathworks 公司发布的主要面对科学计算、可视化以及交互式程序设计的高科技计算环境。它将数值分析、矩阵计算、科学数据可视化以及非线性动态系统的建模和仿真等诸多强大功能集成在一个易于使用的视窗环境中，为科学研究、工程设计以及必须进行有效数值计算的众多科学领域提供了一种全面的解决方案，并在很大

7、程度上摆脱了传统非交互式程序设计语言（如c、fortran）的编辑模式，代表了当今国际科学计算软件的先进水平。课程设计说明书课程设计说明书33 基于 matlab 的 fft 算法设计原理3.1 总体设计流程图 在一个相对较安静的环境下，录下 1s 左右的 wav 声音信号，然后对声音进行采样，画出其时域波形和频谱图，其流程图如图 1 所示： 图 1 设计流程图3.2 语音信号的采集在实际工作中，我们可以利用 windows 自带的录音机录制语音文件。采集到语音信号之后，需要对语音信号进行分析，如语音信号的时域分析、频谱分析、语谱图分析。在 matlab 中，我们可以通过y,fs,bits=w

8、avread(语音信号路径,n1 n2)语句。用于读取语音，采样值放在向量 y 中，fs 表示采样频率(hz)，bits表示采样位数。n1 n2表示读取从 n1 点到 n2 点的值（若只有一个 n 的点则表示读取前 n 点的采样值）。向量 y 则就代表了一个信号（也即一个复杂的“函数表达式”）也就是说可以像处理一个信号表达式一样处理这个声音信号。3.3 语音信号的时频分析利用matlab中的“wavread”命令来读入（采集）语音信号，将它赋值给某一开始输入声音信号对声音信号采样蝶形运算原信号 fft 运算结束课程设计说明书课程设计说明书4向量。再对其进行采样，记住采样频率和采样点数。下面介绍

9、wavread 函数几种调用格式。1.y=wavread（file）功能说明：读取 file 所规定的 wav 文件，返回采样值放在向量 y 中。2.y,fs,nbits=wavread(file) 功能说明：采样值放在向量 y 中，fs 表示采样频率（hz），nbits 表示采样位数。3.y=wavread（file，n）功能说明：读取钱 n 点的采样值放在向量 y 中。4.y=wavread（file，n1,n2）功能说明：读取从 n1 到 n2 点的采样值放在向量 y 中。接下来，对语音信号 speech off.wav 进行采样。其程序如下：y,fs,nbits=wavered (sp

10、eech off.wav);功能说明：把语音信号加载入 matlab 仿真软件平台中然后，画出语音信号的时域波形，再对语音信号进行频谱分析。matlab 提供了快速傅里叶变换算法 fft 计算 dft 的函数 fft,其调用格式如下：xk=fft(xn,n)参数 xn 为被变换的时域序列向量，n 是 dft 变换区间长度，当 n 大于 xn 的长度时，fft 函数自动在 xn 后面补零。，当 n 小于 xn 的长度时，fft 函数计算 xn 的前n 个元素，忽略其后面的元素。原始信号的时域波形图如图 3 所示：课程设计说明书课程设计说明书5图 3 原始信号的时域波形图 原始信号的频域特性图如图

11、 4 所示：课程设计说明书课程设计说明书6图 4 原始信号的频域特性图3.4 快速傅里叶变换快速傅里叶变换（fft）是为提高 dft 运算速度而采用的一种算法。对一个有限长度序列 x(n)的 n 点的 dft 为：x（k）=x（n）wknn (k=0,1,，n-1；n=0,1,，n-1;w=e-j2/n)当 n=4 时，x（k）可展开为：x（0）= x（0）w0*4+ x（1）w0*4 +x（2）w0*4+ x（3）w0*4x（1）= x（0）w0*4+ x（1）w1*4 +x（2）w2*4+ x（3）w3*4x（2）= x（0）w0*4+ x（1）w2*4 +x（2）w4*4+ x（3）w6

12、*4x（3）= x（0）w0*4+ x（1）w3*4 +x（2）w6*4+ x（3）w9*4从上式可以看出,要求 4 点的 dft,需要 16 次的复数乘法运算,12 次复数乘法运算算。由此类推，要求出 n 点的 dft，需要 n2 次复数乘法运算，n*(n-1)次复数加法运算。当 n 值较大时，要完成的复数乘法运算和复数加法运算得次数都非常多，无论是用通用计算机还是用 dsp 芯片，都需要消耗大量的时间，不适合于对实时处理要求高的场合。为了能实时处理 dft,要想减少 dft 的运算量可以有两课程设计说明书课程设计说明书7个途径：第一是降 n，n 的值减小了，运算量就减少了；第二是利用旋转因

13、子的周期性和对称性，可约性。利用这两个途径实现 dft 的快速傅里叶变换（fft） ，fft 算法基本上可分为时域抽取法和频域抽取法。w=e-j2/n 的性质：（1）周期性 （2）共轭对称性 （3）可约性 本程序是用基 2 的按时间抽取的 fft 算法（dit-fft） ，设序列 x(n)的长度为n，且 n 满足 n=2m,m 为正整数。若 n 不能满足上述关系，可以将序列 x(n)补零实现，则 x(n)的 n 点 dft 为：x（k）=x（n）wknn (k=0,1,，n-1；n=0,1,，n-1;w=e-j2/n)将 n 分为奇数与偶数两部分。按时间抽取基 2-fft 算法的基本思路是将

14、n 点序列按时间下标的奇偶分为两个 n/2 点序列，计算这两个 n/2 点序列的 n/2 点 dft，计算量可减小约一半；每一个 n/2 点序列按照同样的划分原则，可以划分为两个 n/4 点序列，最后，将原序列划分为多个 2 点序列，将计算量大大降低。1. 按时间下标的奇偶将 n 点 x(n)分别抽取组成两个 n/2 点序列，分别记为x1(n)和 x2(n)，将 x(n)的 dft 转化为 x1(n)和 x2(n)的 dft 的计算。)()(nnknnnknknnwww*)(*)(nknnknknnwwwmknmnknnmknmnknnwwww/,12, 1 , 0,)() 12()()2(2

15、1nrrxrxrxrx 1n021nn0,2,4.1,3,5.1122212nn0,10,111222121n2n0,10,1221nnknnnnknknnnnrkrkrrnnrkrkrrx kx n wx n wx n wxr wxrwx r wxr w2j2j2222eerknrkrkrknnnww课程设计说明书课程设计说明书8用蝶形运算可表式为：以此类推，还可以把 x1(n)和 x2(n)按 n 值得奇偶分为两个序列，这样就达到了降 n 得目的，从而减少了运算量。fft 对 dft 的数学运算量改进：直接采用 dft 进行计算，运算量为 n2 次复数乘法和 n*(n-1)次复数乘法。当采

16、用 m 次 fft 时，由 n=2m 求得 m=logn，运算流图有 m 级蝶形，每一级都由 n/2 个蝶形运算构成，这样每一级蝶形运算都需要 n/2 次复数乘法和 n 次复数加法。m 级运算共需要复数乘法次数为 c=n/2*m,复数加法次数为 c=n*m。当 n 值较大时，fft 减少运算量的特点表现的越明显。3.4.1 fft 的运算规律（1）原位运算n=2m的fft共m级运算，每级有n/2蝶形原位计算，当数据输入到存储器以后，每一组蝶形运算后，结果仍然存放在这同一组存储器的同一位置，不需要另辟存储空间，直接最后输出。同一级的蝶形运算每个蝶形运算的输入数据对其他级输入没有影响。（2）倒序运

17、算的规律输入序列先按自然顺序存入存储单元,然后经变址运算来实现倒位序排列,用 j表示倒序的十进制数，对 n=2m，m 位的二进制数从左到右各位数权值位n/2，n/4，n/82，1。因此，最高位加 1 相当于 j+n/2。.如果最高位为 0，则直接得到下一个倒序值，j+n/2；.如果最高位为 1，则最高位为 0(j-n/2)，次高位加 1(j+n/4)。.以此类推，直到最后一位二进制数字。例如 ，n=8 时如下图 5 所示： 112212220012, 01nnrkkrknnnrrknx kx r wwxr wx kw xkn（）（k）课程设计说明书课程设计说明书9顺序倒序十进制二进制二进制十进

18、制0000000010011004201001023011110641000011510110156110011371111117图 5 码位倒序（n=8）倒序的流程图如图 3 所示：课程设计说明书课程设计说明书10图 3 倒序的流程图（3）蝶形运算的运算规律设序列x（n）倒序后存放在数组a中，如果蝶形运算两个输入相距b=2(l-1).lh=n/2n1=n-1j=0t=x(j+1)x(j+1)=x(i+1) x(i+1)=ti=0：n1i=jk=lhjkj=j+kj=j-kk=k/2ynyn课程设计说明书课程设计说明书11采用原位运算：蝶形运算可表示为： x(k)=x(k)+x(k+b)*wn

19、p; x(k+b)=x(k)-x(k+b)*wnp.其中p=j*2（m-1）,(j=0,1,2,3,2(l-1)-1) 旋转因子的确定：第l级共有2(l-1)个旋转因子，以(n=23=8为例)：l=1时，w=w j=0；pnjl2l=2时，w=w j=0，1；pnjl2l=3时，w=w j=0，1，2，3；pnjl2w=w=w p=j*2(m-l)，用来确定第l级旋转因子pnjl2)(2*lmjn从输入端开始，共进行 m 级运算，在进行第 l 级运算时，依次求出 2(l-1)个旋转因子，然后计算每个旋转因子所对应的 2(m-l)个蝶形元素。第 l 级的蝶形运算中，每个蝶形运算的两个输入相距 b

20、=2l-1,同一旋转因子对应的蝶形运算相隔2l 个，同一旋转因子对应的蝶形运算有 2m-l 个。注：1：控制第 l 级顺序运算2：控制不同种的旋转因子3：控制同种旋转因子所对应的蝶形运算蝶形运算的流程图如图 4 所示：课程设计说明书课程设计说明书12 1 2 3图 4 蝶形运算的流程图3.4.2 基于 matlab 的 fft 所编写程序的框图 基于 matlab 的 fft 所编写程序的框图如图 5 所示：开始开始开始j=0：b-1n=2mb=2l/2l=1:m输入 x(n)k=j+1:2l:np=j*2(m-l)t=x(k)+x(k+b)* wnp x(k+b)=x(k)-x(k+b)*

21、wnp x(k)=t倒序结束输出课程设计说明书课程设计说明书13图 5 程序框图3.5 自编算法与机带算法仿真波形比较 我们知道 matlab 软件自带 fft 算法，我们可以通过比较自编算法仿真结果与机带算法仿真结果来检验自编算法的正确性。自编算法与 fft 算法幅值比较图如图 6 所示：语音信号采集完成语音信号时域图完成语音信号频率特性图编写 fft 程序，画出语音信号频谱图实现输入信号的倒序实现一级中不同种蝶形算运实现一级中相同种蝶形运算与matlab自带的fft比较课程设计说明书课程设计说明书14图 6 幅值比较图自编算法 fft 算法频谱分贝计较图如图 7 所示：课程设计说明书课程设

22、计说明书15图 7 分贝比较图由以上两图可以看出，经过蝶形运算得出的频谱图和信号直接 fft 得出的频谱图一致。该程序严格按程序框图编写，思路清晰、容易理解，程序的运行过程在命令窗中一目了然。通过与 fft 函数运算的结果比对，程序编写正确。课程设计说明书课程设计说明书164 设计总结为期四天的课设很快接近尾声了，基于 matlab 的 fft 算法实现设计已按计划如期全部完成，通过这次 dsp 课程设计，我对课堂上所学到的理论知识的理解加深了许多， 自己动脑、动手设计的能力也得到了较大提高。在这次课程设计的过程中，我对 matlab 语言有了更深的认识。现在仔细想想，这次课程设计使得我对 m

23、atlab 语言的理解与应用能力得到了较大的提升，也让我认识到只要深入学习，提升的空间永远是存在的。在设计的过程中我遇到了一些问题，如：编写源程序中出现了语法错误等。通过查阅书本和以前设计的程序我发现了产生错误的原因并解决了问题完成了设计。经过反思我发现较大一部分错误是因为操作的不熟练造成的，这也让我明白了要保持设计的高效率 必须经常练习。另一方面我也发现了动手实践的重要性。动手实践是理论知识得以灵活运用的必要前提，也是今后走上工作岗位之后能够很好的完成设计工作的技术保证。只有遇到实际问题并根据自己对课堂上获得的专业知识的理解来解决才能真正的提高自己的能力。这次程序设计让我获益匪浅，对 mat

24、lab 也有了进一步的认识：matlab 功能强、使用灵活方便等。matlab 是在国内外广泛使用的一种数字信号处理系统，相信除了以上优点，还有许多我还未发现，希望能在以后的学习中有更深入的认识。在此，我也要感谢一下董翠英老师对我这次课设的大力支持和指导，如果没有老师的帮助我很难在这么短的时间内完成这次课设。课程设计说明书课程设计说明书17参考文献1范寿康 dsp 技术与 dsp 芯片.北京：电子工业出版社2程佩青.数字信号处理教程.北京：清华大学出版社出版，20013高西全 丁玉美等 数字信号处理. 北京：电子工业出版社，20094 余成波，陶红艳。数字信号处理及 matlab 实现，北京：

25、清华大学出版社，2008 5 曹弋，赵阳。matlab 实用教程，北京：电子工业出版社，2007 课程设计说明书课程设计说明书18附录源程序：clear all;x,fs,bits=wavread(gl.wav,2048); %读取声音x1=reshape(x,1,4096);sound(x1,fs,bits); %播放语音信号y1=fft(x1)n=length(x1); %计数读取信号的点数t=(1:n)/fs; %信号的时域采样点figure(1)plot(t, x1); %画出声音采样后的时域波形title(原声音信号的时域波形); %给图形加注标签说明xlabel(时间/t);yla

26、bel(振幅/a);grid ; %添加网格m=nextpow2(x1); % 求 x 的长度对应的 2 的最低幂次m n=2m;if length(x1)n x1=x1,zeros(1,n-length(x1); % 若 x 的长度不是 2 的幂，补零到2 的整数幂 end nv2=n/2;nm1=n-1;i=0;j=0;while inm1 if ij t=x1(j+1); x1(j+1)=x1(i+1); x1(i+1)=t; end k=nv2; 课程设计说明书课程设计说明书19 while k=j j=j-k; k=k/2; end j=j+k; i=i+1;end %x1; y=x1; %