1、数字信号处理中南大学现代信号处理课程设计报告 学 院：信息科学与工程学院 专业班级： 姓 名： 学 号 指导老师：赵亚湘、郭丽梅 设计时间：2012年9月 目 录第一章 课程设计题目及设计要求1.1 课程设计题目1.2 课程设计目的及要求第二章 设计思想和系统功能分析2.1 第一题的设计思想2.2 第二题的设计思想2.3 第三题的设计思想2.4 第四题的设计思想2.5 第五题的设计思想2.6 第六题的设计思想第三章 关键部分的设计思路3.1 数字滤波器的设计思路第四章 调试及结果分析4.1 第一题的调试结果及分析4.2 第二题的调试结果及分析4.3 第三题的调试结果及分析4.4 第四题的调试结
2、果及分析4.5 第五题的调试结果及分析4.6 第六题的调试结果及分析第五章 课程设计总结5.1 课程设计中遇到的问题及解决及心得体会附 录源程序清单及参考文献第一章 课程设计题目及设计要求1.1 课程设计题目1. 给定模拟信号： 1）选择采样频率Fs = 5000Hz和合适的信号长度，采样得到序列 x1(n)。求并画出x1(n)及其序列傅里叶变换 |X1(ejw)|。 2）选择采样频率Fs = 1000Hz和合适的信号长度，采样得到序列 x2(n)。求并画出x2(n)及其序列傅里叶变换 |X2(ejw)|。 3） 说明|X1(ejw)|与|X2(ejw)|间的区别，为什么？2. 已知两系统分别
3、用下面差分方程描述：试分别写出它们的传输函数，并分别打印曲线。说明这两个系统的区别。3. 已知已调信号，其中调制信号，载波为。1） 选择合适的采样频率及信号长度，使用FFT分析该已调信号的频谱并画出其幅频特性和相频特性曲线图2） 对该已调信号进行解调，恢复原调制信号。4. 已知三角波序列和反三角波序列： 用N=8点FFT分析序列x1(n)和x2(n)的幅频特性，观察两者的序列形状和频谱特性曲线有什么异同？绘出两序列及其幅频特性曲线。在x1(n)和x2(n)的末尾补零，用N=16点FFT分析这两个信号的幅频特性，观察幅频特性发生了什么变化？两情况的FFT频谱还有相同之处吗？这些变化说明了什么？5
4、.设有一信号，设计各种IIR数字滤波器以实现：低通滤波器，滤除的成分，保留成分高通滤波器，滤除的成分，保留成分带通滤波器，滤除的成分，保留成分带阻滤波器，滤除的成分，保留成分要求：1) 求出各个滤波器的阶数，设计各滤波器。画出各滤波器的幅频和相频特性，计算滤波器的系统函数H（z） 2) 画出滤波前后信号的时域、频域波形6. (1) 用Hanning窗设计一线性相位带通数字滤波器，要求：N=15, 。观察它的实际3dB和20dB带宽。N=45，重复这一设计，观察幅频和相位特性的变化，注意长度N变化的影响；下面两题选做：7.音乐信号处理： 1）获取一段音乐或语音信号，设计单回声滤波器，实现信号的单
5、回声产生。给出加入单回声前后的信号频谱。 2）设计多重回声滤波器，实现多重回声效果。给出加入多重回声后的信号频谱。 3）设计均衡器，使得得不同频率的混合音频信号，通过一个均衡器后，增强或削减某些频率区域。8. 对混有噪声的音乐信号进行滤波： 1）在一段音乐信号中混入两个频率的正弦型干扰信号，利用FFT计算该混合信号的频谱并确定干扰信号的频谱； 2）利用二阶带阻滤波器设计一个能滤出干扰信号的梳状滤波器； 3）利用梳状滤波器滤除信号中的噪声，播放处理前后的信号，并比较处理前后的结果1.2 课程设计目的及要求一、课程设计目的：1全面复习课程所学理论知识，巩固所学知识重点和难点，将理论与实践很好地结合
6、起来。2. 掌握信号分析与处理的基本方法与实现3提高综合运用所学知识独立分析和解决问题的能力；4熟练使用一种高级语言进行编程实现。二、课程设计要求1、使用MATLAB（或其它开发工具）编程实现上述内容，写出课程设计报告。滤波器设计题目应尽量避免使用现成的工具箱函数。为便于分析与观察，设计中所有频谱显示中的频率参数均应对折叠频率归一化。2、课程设计报告的内容包括：(1)课程设计题目和题目设计要求；(2)设计思想和系统功能分析；(3)设计中关键部分的理论分析与计算，关键模块的设计思路； (4)测试数据、测试输出结果，及必要的理论分析和比较；(5)总结，包括设计过程中遇到的问题和解决方法，设计心得与
7、体会等；(6)参考文献；(7)程序源代码清单。3、课程设计考核要求： 课程考核分三部分，一部分是上机率，占20%；第二部分是检查成绩，最后两次上机为检查时间，占50%；第三部分为课程设计报告，占30%。注意：（1）演示系统使用GUI界面或混合编程实现集成打包发布，酌情加分；（2）若发现程序或课程设计报告雷同，一律不及格。4、课程设计进度安排见表1.1所示：序号阶 段 内 容合计（天）一布置设计任务，准备相关资料1二查阅资料，弄清题目要求，提出解决方案1三方案选择及初步设计2四目标项目设计实现及调试2五撰写课程设计报告2六上机检查成绩2总计（2周）10表1.1 课程设计进度安排表第二章 设计思想
8、和系统功能分析 2.1 第一题的设计思想 第一题主要的考查点是对连续时间信号的抽样，对抽样信号进行傅里叶变换以及不同采样频率对幅频特性曲线的影响等。主要的理论依据如下：对连续时间信号的抽样：设连续时间信号为x(t)，则对其采样得到的离散时间信号为x(n)=x(nT)，T为采样周期；离散时间信号的傅里叶变换：。 由于抽样后的离散时间信号题目已给出，故无需再做抽样，可直接对离散时间信号进行傅里叶变换，在MATLAB中，这一变换可以直接用FFT来实现，再绘出变换后得到的幅频特性曲线即可。2.2 第二题的设计思想 所以传输函数为： (-jw) -jw H1(e )=1+e ; (-jw) -jw H2
9、(e )=1-e ;2.3 第三题的设计思想 由题得本题第一问是对已调信号选择合适的采样频率及信号长度，使用FFT分析该已调信号的频谱并画出其幅频特性和相频特性曲线图，本程序选择的20Hz的采样频率对信号进行采样，然后进行fft分析，并分别画出其幅频特性和相频特性曲线图即可。第二问解调是先用载波与此信号相乘，再将得到的信号通过滤波器滤除其他频率的信号便可得到原调制信号。2.4 第四题的设计思想 首先先用subplot和stem等函数画出x1(n)和x2(n)的序列，得到两个序列。用N=8点FFT分析序列x1(n)和x2(n)的幅频特性，绘出两序列及其幅频特性曲线。在x1(n)和x2(n)的末尾
10、补零，用N=16点FFT分析这两个信号的幅频特性，进行FFT之前要对两个信号进行周期延拓，然后取其主值序列，再分别对其进行FFT变换画出其幅频特性和相频特性曲线图即可。 2.5 第五题的设计思想IIR滤波器的设计通常是借助于模拟滤波器的设计方法来进行的。这是因为当离散时间滤波器出现之前连续时间滤波器的设计方法已常成熟，而且实现IIR滤波器的非迭代直接设计法比较困难。设计IIR数字滤波器的任务就是寻求一个因果物理可实现的系统函数：使它的响应H(exp(j)=H(z)满足所希望得到的频域指标，即符合给定的通带衰耗、阻带衰耗、通带截频和阻带截频。设计IIR数字滤波器的方法有间接法、直接法和计算机辅助
11、法等。这里常用的是根据频响指标，设计出相应的模拟滤波器H(s)，再通过脉冲响应不变法或双线性法转换成数字滤波器H(z)。其过程如图： 双线性不变法 频域指标 设计模拟滤波器 数字滤波器 脉冲响应不变法1 设计模拟滤波器：首先把数字指标转换成相应的模拟指标，然后设计符模拟指标的滤波器。2 从模拟滤波器变换到数字滤波器：当完成模拟滤波器设计，求得归一化系统函数Ha(s)以后，可以根据给定的技术指标恰当的选取从模拟域到数字域的映射方法。将系统函数Ha(s)从s平面转换到z平面的方法有多种，常用的是脉冲响应不变法与双线性变换法。使用脉冲响应不变法设计数字滤波器时，很有可能会导致数字滤波器频率响应出现混
12、叠现象，它只适用于带限滤波器。而采用双线性变换法就能解决这个问题则可按双线性变换是变量s和z之间的代数变换，它将S平面的整个j轴映射成Z平面的整个单位圆周。因为-+映射成-，所以连续时间的频率变量和离散时间的频率变量的变换必定是非线性的。因此这种方法只能用于相应的频率轴的畸变在允许范围之内的情况。双线性变换法相当于用如下表示来代替s： 可以选取任何方便的T值，得到变量和之间的关系式： tan函数实现的是一种非线性变换，它使得S平面中的虚轴j映射成Z平面单位圆时，会产生频率畸变。但当这种非线性压缩在允许范围内加以补偿时，比如,在滤波器具有近似理想的分段恒定幅度特性的情况下，它能有效地实现稳定的模
13、拟滤波器向稳定的数字滤波器变换，克服了冲激响应不变法的多值对应而产生的混叠现象。2.6 第六题的设计思想第六题的主要考查点是几种主要的窗函数（矩形窗、汉宁窗、哈明窗等）对离散时间信号的截取，并用FFT对其频谱进行分析。其主要理论依据如下：（1）矩形窗（Rectangle Window） 其幅度函数为：（2）汉宁（Hanning）窗升余弦窗 其频谱函数为： 其幅度函数为：（3）哈明（Hamming）窗改进的升余弦窗 其频谱函数为： 当时，其幅度函数可近似表示为：设计流程如图所示：第三章 关键部分的设计思路 在这里值得我们重提一次的是数字滤波器的设计思路，因为它是这次设计的重点，也是难点。 上课的
14、时候，这一块知识学得很仓促，不扎实。一些概念也很模糊，这次设计正好可以系统温习一遍关于滤波器设计的知识，所以，在此将滤波器的设计思路再更加详尽地复述一次。3.1 数字滤波器的设计思路 实际应用中的IIR数字滤波器有低通、高通、带通和带阻等类型。设计各种数字滤波器，通常可以把一个归一化的原型模拟低通滤波器经模拟频带变换成所需类型的模拟滤波器，再通过转换成为所需类型的数字滤波器，如图3.1(a)所示。更简便的方法是从模拟低通归一化原型通过一定的频率变换关系，一步完成各类数字滤波器的设计。常用模拟低通滤波器有巴特沃兹(Butterworth，简称为BW型)，切比雪夫(Chebyshev，简称为CB型
15、)及椭圆型或考尔型(Cauer，简称为C型) 滤波器，而模拟高通、带通、带阻滤波器则可以利用变量变换方法，由低通滤波器变换得到。从模拟滤波器映射到数字滤波器有两种方法：脉冲响应不变法和双线性变换法。脉冲响应不变法又称时域取样法，它是使数字滤波器的单位脉冲响应序列h (n)逼近模拟滤波器的冲激响应ha (t)，且使得h (n) = ha (nT)。它的缺点是产生了频响的混叠。为了克服这一缺点，可采用双线性变换法。双线性变换法是先根据给定的技术指标设计一个系统函数为H(s)的模拟滤波器，再通过适当的数学变换，将无限宽的频带，变换成频带受限的系统函数H( s1)，然后再通过z变换求得数字滤波器的系统
16、函数H(z)。第四章 调试及结果分析 4.1 第一题的调试结果及分析 第一题的调试结果图结果分析：分别以所要求采样频率分别为f1=5000Hz,f2=1000Hz,对进行等间隔采样，得到，（为采样周期），可知，采样信号的频谱函数是在原模拟信号频谱函数的周期延拓；若以频率为自变量，则以采样频率为延拓周期；对频带限于的模拟信号，根据采样定理，只有当时，采样后才不会发生频谱混叠失真。但是在Matlab上是无法计算连续函数，只有在当足够大的时候，我们才将频谱混叠忽略不计，从而可对采样序列进行傅里叶变换；最后应用subplot()命令实现画图，并注意到要归一化。由以上实验得到的图形结果可以看到：，当采样
17、频率越大的时候，采样信号的频谱越陡峭，而其失真情况也越来越小，根据采样定理“采样频率越接近信号频率，其失真情况就越小”可以得到上述结论。4.2 第二题的调试结果及分析 第二题的调试结果图结果分析：运行结果可以看出，H1和H2都可以构成梳状滤波器，但是滤去的波的频率不同。而且其对应的零极点分布图也有所不同，零极点的位置有所差别。4.3 第三题的调试结果及分析 第三题的调试结果图结果分析：由运行结果可以看出已调信号的幅频特性曲线在-5和+5附近峰值较大，其他位置较小。而相频特性曲线则是在两处较松散，说明这两处的相位变化较快。解调时让已调信号乘以载波，再通过滤波器滤除其他的信号即可得到。通过图像可知
18、解调效果很好4.4 第四题的调试结果及分析第四题的调试结果图结果分析：由N=8时所得到的两幅频谱特性曲线知道，正反三角波序列的频谱图像时一样的。因为进行DFT之前要对两个信号进行周期延拓，然后取其主值序列，而正反三角函数的主值序列是相同的，只相差一个相位，因此他们的幅频特性也是相同的。而N=16时，由于后面补0，所以周期延拓后的主值序列不相同了，正三角波周期延拓后仍然是连续的，而反三角波周期延拓之后在N=8的地方出现了陡峰，因此频域中的频谱分量会增多。两图比较可知，末尾补0可以使谱线变密，从而减小栅栏效应。增加频率采样点数，原来漏掉的某些频谱分量就肯恩被检测出来。4.5 第五题的调试结果及分析
19、原序列低通滤波器 高通滤波器 带通滤波器 带阻滤波器结果分析： 实际应用中的IIR数字滤波器有低通、高通、带通和带阻等类型。设计各种数字滤波器，通常可以把一个归一化的原型模拟低通滤波器经模拟频带变换成所需类型的模拟滤波器，再通过转换成为所需类型的数字滤波器，如图3.1(a)所示。更简便的方法是从模拟低通归一化原型通过一定的频率变换关系，一步完成各类数字滤波器的设计，如图3.1(b)所示。常用模拟低通滤波器有巴特沃兹(Butterworth，简称为BW型)，切比雪夫(Chebyshev，简称为CB型)及椭圆型或考尔型(Cauer，简称为C型) 滤波器，而模拟高通、带通、带阻滤波器则可以利用变量变
20、换方法，由低通滤波器变换得到。本题的关键是选着合适的技术指标来设计滤波器，分别对应低通、高通、带通和带阻四种类型。 从模拟滤波器映射到数字滤波器有两种方法：脉冲响应不变法和双线性变换法。脉冲响应不变法又称时域取样法，它是使数字滤波器的单位脉冲响应序列h (n)逼近模拟滤波器的冲激响应ha (t)，且使得h (n) = ha (nT)。它的缺点是产生了频响的混叠。为了克服这一缺点，可采用双线性变换法。双线性变换法是先根据给定的技术指标设计一个系统函数为H(s)的模拟滤波器，再通过适当的数学变换，将无限宽的频带，变换成频带受限的系统函数H( s1)，然后再通过z变换求得数字滤波器的系统函数H(z)
21、。 4.6 第六题的调试结果及分析汉宁窗运行结果分析结论 由上面两图所得图像分析可得，随着N的增大，相频主瓣宽度变窄，同时，取样值增大，幅频宽度增加，变化较缓，相位特性急剧变化。 分别改用矩形窗和Blackman窗，设计(1)中的带通滤波器，观察并记录窗函数对滤波器幅频特性的影响，比较三种窗的特点；总结窗的不同长度和不同窗对滤波器的影响矩形窗运行结果分析结论 由上面两图得出，当加了矩形窗之后，随着N的增大，主瓣宽度变窄，同时主瓣幅度增大，旁瓣幅度也增大，但主瓣旁瓣的幅度相对值不变。幅频变化剧烈；相位变化变缓；blackman窗:分析结论： 由上面两图比较可知，随着N的增大，相频变化速度变慢，幅
22、频变化剧烈总结： 由上面三重窗在不同N之情况下的图象分析可知，当N增大时，幅频主瓣宽度变窄，同时相频变化速度变慢，幅频变化剧烈。第五章 课程设计总结5.1 课程设计中遇到的问题及解决及心得体会这次课程设计过程中遇到的主要问题是频率归一化的问题，开始对归一化的概念不是很清楚。在查阅图书馆资料及问了班上能力较强的同学后，对这个概念深入理解了，并在程序中成功地实现了。通过这两周的现代信号处理课程设计实践，我复习了MATLAB编程语言的基本概念、语法、语义和数据类型的使用特点，加深了对课堂所学理论知识的理解，掌握了运用结构化程序设计的基本思想和方法，更重要的是培养了自己的自学能力。因为这是我们第一次接
23、触MATLAB编程语言，在编写程序以及调试的过程中遇到了很多困难，但是我通过去图书馆查找资料，请教同学老师，再自己一点点改善程序，最终编写出一个比较完善的程序，实现了所有要求功能，这是最值得我欣慰的一点。 以下是我的几点切身感受： （一）编写程序需要一个清醒的头脑，明确的思路，同时也要有耐心毅力。刚拿到程序设计课题时，我感觉一片茫然，因为在之前的信号处理学习中，只是在一些例题中接触过MATLAB语言，甚至没有看过一些基础的书籍，更没有上过课，所以初次遇到一个实际问题，感觉无从下手。这是由于自己对MATLAB的模块设计不够理解，同时对MATLAB的基本语句一无所知，不过通过请教老师懂得了首先要设
24、计guide目录，再根据用户输入,执行swtich语句，在每个case后调用一个函数，来实现要求的功能，这样一下子豁然开朗，掌握了基本设计思路之后，后面的编程就顺利多了。至此，我真真体会到清晰地思路对成功编写一个程序的重要性。 当然成功编写一个程序绝非易事，之前，我总以为程序能够正常运行，就代表着编程成功，后来我才发现我大错特错了。我用了三天时间，完成了程序的编写、改错，但我立刻发现尽管程序能够正常运行，部分功能却不完善，甚至不能实现。经过一次又一次调试、修改又修改，一点一点发现问题并改正，我才真正发现编程远没有想象中的简单。它需要的不仅是清晰地编程思路、编程技巧，还需要有耐心有毅力，不要放弃
25、。（二）我们在大学需要学习的不仅是基础知识、专业知识，更重要的是一种学习能力。正如老师所说学习是终生的，因此学习能力也就能让我们受益终生。由于课堂时间有限，MATLAB课程只能靠我们自学，然而guide画图功能我们必须用到，因此只能通过自习，实践也证明我们是有自学的潜能的，通过自学课本，不理解的知识，通过查找资料，请教老师，最终掌握知识，因此，这次课程设计时间培养了我们的自学能力。（三）要学会向他人请教，与他人合作。在编程过程中，我遇到了许多问题，通过查找资料也未能找到原因，每次我想老师请教，老师总能不厌其烦地一点点讲解。有时，老师为了帮我在那几百行的程序中找到问题，花上一个多小时，这让我非常
26、感动。大学中的老师是每个大学生的好资源，他们的博学、重教，向他们请教会有很大收获。在编程过程中，我与我班同学共同讨论，为对方查找错误，提示建议，在这个过程中我收获甚多，也培养了自己的合作精神。 （四）要学会总结，提高。在编程过程中，我发现自己的程序还有很多地方可以完善。程序的界面可以做得更加人性化，这样使用户使用更加方便，画图函数也还可以改进，例如：在使用figure分图显示的时候，我尽量选择部分图，都在一个figure中显示，这样既符合了，老师的要求尽量在一个图中完成设计，同时也使整个设计的显得干净利落。附 录源程序清单第一题：1）Dt=0.00005;t=-0.005:Dt:0.005;x
27、a=exp(-1000*abs(t);Ts=0.0002;n=-25:1:25;x=exp(-1000*abs(n*Ts);K=500;k=0:1:K;w=pi*k/K;X=x*exp(-j*n*w);X=real(X);w=-fliplr(w),w(2:K+1);X=fliplr(X),X(2:K+1);subplot(1,1,1)subplot(2,1,1);plot(t*1000,xa);xlabel(t/ms);ylabel(xl(n)title(离散信号);hold onstem(n*Ts*1000,x);gtext(Ts=0.2ms);hold offsubplot(2,1,2);
28、plot(w/pi,X);xlabel(频率：pi);ylabel(Xl(w)title(离散时间傅里叶变换)2）Dt=0.00005;t=-0.005:Dt:0.005;xa=exp(-1000*abs(t);Ts=0.001;n=-5:1:5;x=exp(-1000*abs(n*Ts);K=100;k=0:1:K;w=pi*k/K;X=x*exp(-j*n*w);X=real(X);w=-fliplr(w),w(2:K+1);X=fliplr(X),X(2:K+1);subplot(1,1,1)subplot(2,1,1);plot(t*1000,xa);xlabel(t毫秒);ylabe
29、l(xl(n)title(离散信号);hold onstem(n*Ts*1000,x);gtext(Ts=1毫秒);hold offsubplot(2,1,2);plot(w/pi,X);xlabel(以pi为单位的频率);ylabel(Xl(w)title(离散时间傅里叶变换)第二题：A=1;B=1,1;H,w=freqz(B,A,whole);subplot(2,1,1);plot(w/pi,abs(H);xlabel(omega/pi);ylabel(|H(ejomega)|);title(y1(n)=x(n)+x(n-1) 幅频响应特性);axis(0,2,0,2.5);A=1;B=1
30、,-1; H,w=freqz(B,A,whole);subplot(2,1,2);plot(w/pi,abs(H);xlabel(omega/pi);ylabel(|H(ejomega)|);title(y2(n)=x(n)-x(n-1) 幅频响应特性);axis(0,2,0,2.5);第三题：fs=20;fc=9/2;k=0:199;t=k/fs;N=256;f=linspace(0,fs*(N-1)/N,N)-fs/2;x=cos(pi*t);y=x.*cos(9*pi*t);%已调信号y1=y.*cos(9*pi*t);%解调信号Y=abs(fftshift(fft(y,256);Z=angle(fftshift(fft(y,256);b,a=butter(4,9*pi,s);xb=lsim(tf(b,
