数字信号处理实习

1、 数字信号处理课程综合设计（十六） 数 字 信 号 处 理课程设计 题目： 数字信号处理课程综合设计（十六） 学生姓名：朱孟达学 号：201484250216班 级: 电信1402班专 业：电子信息工程所在院(系): 电气系指导教师：王创新 谢文彪 文卉起止日期 2016年12月26日30日2016年 12 月第 1 页 共 22 页目录1 本课程设计目的11.1 基本要求11.2 设计达到效果12 设计的详细原理22.1 采样定理22.2 卷积定理22.3 IIR滤波器原理32.4 FIR滤波器的原理32.5 窗函数原理43 设计的步骤和过程54 设计程序的调试和运行结果114.1 采样信号

2、的调试114.2 IIR滤波器的设计的调试114.3 FIR滤波器的调试125 课程设计的思考与体会12参考文献13附录：14附录A 采样信号MATLAB 程序代码14附录B IIR 滤波器MATLAB 程序代码16附录C FIR滤波器MATLAB 程序代码19第 1 页 共 23 页1 本课程设计目的1.1 基本要求(1)、掌握数字信号处理的基本概念，基本理论和基本方法。(2)、熟悉离散信号和系统的时域特性。(3)、掌握序列快速傅里叶变换方法。(4)、学会 MATLAB 的使用，掌握 MATLAB 的程序设计方法。(5)、掌握利用 MATLAB 对语音信号进行频谱分析。(6)、掌握滤波器的网

3、络结构。(7)、掌握 MATLAB 设计 IIR 数字滤波器的方法和对信号进行滤波的方法。1.2 设计达到效果学生通过课程设计加强实践能力，同时巩固数字信号处理的知识，学会设计设计滤波器，熟练掌握matlab软件的使用。通过matlab完成课程设计的仿真，最后提升自己的整个学习能力。2 设计的详细原理2.1 采样定理 对模拟信号进行采样可以看做一个模拟信号通过一个电子开关S。设电子开关每隔周期T合上一次，每次合上的时间为T，在电子开关输出端得到其采样信号。该电子开关的作用等效成一宽度为,周期为T的矩形脉冲T(t),采样信号就是与T(t)想乘的结果。 （2-1） （2-2）2.2 卷积定理1、

4、离散时间序列f1(k)和f2(k)的卷积和定义： （2-3） 2、 在离散信号与系统分析中有两个与卷积和相关的重要结论： a、 （2-4）即离散序列可分解为一系 列幅度由f(k)决定的单位序列(k)及其平移序列之积。 b、对线性时不变系统，设其输入序列为f(k)，单位响应为h(k)，其零状 态响应为y(k)，则有： （2-5）2.3 IIR滤波器原理 “递归滤波器”。递归滤波器，也就是IIR数字滤波器，顾名思义，具有反馈。IIR数字滤波器的设计方法有两类：间接设计法和直接设计法。间接设计法是借助模拟滤波器设计方法进行设计的，先根据数字滤波器设计指标设计相应的过渡模拟滤波器，再将过渡模拟滤波器转

5、换为数字滤波器。直接设计法师在时域或频域直接设计数字滤波器。 由于模拟滤波器设计理论非常成熟，而且有很多性能优良的典型滤波器可供选择（如，巴特沃斯滤波器、切比雪夫滤波器、椭圆滤波器等），设计公式和图表完善，而且许多实际应用需要模拟滤波器的数字仿真，所以间接设计法得到广泛的应用。N阶IIR数字滤波器的系统函数： （2-6）2.4 FIR滤波器的原理 FIR(Finite Impulse Response)滤波器：有限长单位冲激响应滤波器，又称为非递归型滤波器，是数字信号处理系统中最基本的元件，它可以在保证任意幅频特性的同时具有严格的线性相频特性，同时其单位抽样响应是有限长的，因而滤波器是稳定的系

6、统。设单位脉冲响应h(n)长度为N，器系统函数H(z)和差分方程分别为： （2-7） （2-8）2.5 窗函数原理在实际进行数字信号处理时，往往需要把信号的观察时间限制在一定的时间间隔内，只需要选择一段时间信号对其进行分析。这样，取用有限个数据，即将信号数据截断的过程，就等于将信号进行加窗函数操作。而这样操作以后，常常会发生频谱分量从其正常频谱扩展开来的现象，即所谓的“频谱泄漏”。当进行离散傅立叶变换时，时域中的截断是必需的，因此泄漏效应也是离散傅立叶变换所固有的， 必须进行抑制。而要对频谱泄漏进行抑制，可以通过窗函数加权抑制 DFT 的等效滤波器的振幅特性的副瓣，或用窗函数加权使有限长度的输

7、入信号周期延拓后在边界上尽量减少不连续程度的方法实现。而在后面的FIR 滤波器的设计中，为获得有限长单位取样响应，需要用窗函数截断无限长单位取样响应序列。另外，在功率谱估计中也要遇到窗函数加权问题。由此可见，窗函数加权技术在数字信号处理中的重要地位。下面介绍窗函数的基本概念。设x(n)是一个长序列， w(n)是长度为N的窗函数，用w(n)截断x(n)，得到N点序列xn(n)，即 （2-9）在频域上则有 （2-10）由此可见，窗函数 w(n)不仅仅会影响原信号 x(n)在时域上的波形，而且也会影响到频域内的形状。数字信号处理领域中所用到的基本窗函数主要有：表1-1窗函数旁瓣峰值幅度/db过渡带宽

8、阻带最小衰减/db矩形窗-134/N12三角形窗-258/N-25汉宁窗-318/N-44哈明窗-418/N-53不萊克曼窗-5712/N-74凯塞窗-5710/N-803 设计的步骤和过程3.1采样定理的验证对下面连续信号采样: A=444.128, a= 10,=222.144， A 为幅度因子， a为衰减因子，为模拟角频率。画出原信号时域图，x轴为时间，y轴为振幅，如图1：图1原信号时域图 对原信号采样，x轴为离散点，y轴为振幅，如图2：图2原信号采样时域图原信号采样频域图，x轴为频率，y轴为振幅，如图3：图3原信号采样频域图改变采样频率，频谱混叠时域图，x轴为采样点，y轴为振幅，如图4

9、：图4频谱混叠时域图频谱混叠频域图，x轴为频率，y轴为振幅，如图5：图5频谱混叠频域图3.2 IIR滤波器的设计 录自己的一段语音，画出原信号时域图频域图，时域图x轴为时间，y轴为振幅，频域图x轴为频率，y轴为振幅，如图6：图6原信号时域图频域图 加入噪音信号，得出加噪信号时域频域图，时域图x轴为时间，y轴为振幅，频域图x轴为频率，y轴为振幅，如图7：图7加噪时域图频域图同时画出椭圆IIR滤波器幅频图，x轴为频率，y轴为DB，如图8：图8幅频图然后对加噪信号进行滤波，得出滤波后信号时域频域图，时域图x轴为时间，y轴为振幅，频域图x轴为频率，y轴为振幅，如图9：图9时域频域图3.3FIR滤波器的

10、设计 根据上次技术指标，仿真出汉明窗FIR滤波器，FIR窗函数滤波器幅频图，x轴为频率，y轴为振幅，如图10：图10FIR窗函数滤波器幅频图 FIR滤波器时域频域图，时域图x轴为时间，y轴为振幅，频域图x轴为频率，y轴为振幅，如图11：图11时域频域图4 设计程序的调试和运行结果4.1 采样信号的调试 根据实习的要求计算出信号的截止频率fc,当采样频率fs大于二倍fc时，由上图可知采样频谱并未发生频谱混叠，当fs小于二倍fc时，会产生频谱混叠。4.2 IIR滤波器的设计的调试 通过自己在电脑录的一段自己的一段话，在matlab中，利用wavread函数对语音信号进行采样，同时利用快速傅里叶变换

11、与plot画出频谱图，通过分析得到技术指标，通带截止频率，阻带截止频率。然后加上一段噪音，画出加噪信号的时域图和频谱图。然后根据技术指标仿真出滤波器，同时画出滤波器的幅频特性图，以下是滤波器的程序：Fs=10000;%采样频率wp=600*2/Fs;%通带截止频率ws=800*2/Fs;%阻带截止频率Rp=1;%通带最大衰减Rs=40;%阻带最大衰减n,Wn=ellipord(wp,ws,Rp,Rs,s); %计算阶数 边界频率B,A=ellip(n,Rp,Rs,Wn,low);%计算滤波器系数h,w=freqz(B,A,512);%求数字滤波器的频率响应plot(w*Fs/(2*pi),20

12、*log10(abs(h)/max(abs(h);axis(0,1200,-65,3);title(数字低通滤波器的幅度响应); xlabel(Hz/频率);%x轴ylabel(DB/分贝数);%y轴然后对加噪信号滤波，画出滤波后的信号频谱。最后完整滤出噪音信号。4.3 FIR滤波器的调试同样在技术指标下，仿真出滤波器，下面是滤波器设计程序：wp=600/(0.5*fs);%通带截止频率ws=800/(0.5*fs);%阻带截止频率rp=1;%通带最大衰减rs=60;%阻带最大衰减N1,wc=buttord(wp,ws,rp,rs);%求阶数B,A=butter(N1,wc);%Filter(

13、);yo=filter(B,A,x2);%滤波hn=Fir1(95,0.5*(wp+ws),low,hanning(96);%加窗口函数求系统函数h1,w=freqz(hn,1);figure(2);plot(w*fs/(2*pi),20*log10(abs(h1);xlabel(hz/频率);ylabel(DB);title(窗函数幅频特性);yof=conv(hn,x2); 最后画出滤波后的信号频谱图，看滤波效果。5 课程设计的思考与体会 在课程设计的这段时间，我获益匪浅。不但进一步掌握了数字信号处理的基础知识及MATLAB的基本操作，还使我了解了信号的产生、采样及频谱分析的方法。我进一步

14、了解到凡事都需要耐心，细心仔细是成功的保证。虽然在做的过程中遇到了一些问题，但是我都通过自己的努力解决了它们，证明了自己的能力。这次课程设计对我各方面的综合能力有了很大的提高，对我以后的工作，实践都有很大的帮助。 在此次课程设计当中，我经常把C语言的语法知识照搬到MATALAB设计中，从而导致调试失败，所以下次用此类语言做课程设计时，应事先学习下这类语言的基本语法，以免与其他语言相混淆。还有就是有些不定参数存在时，可先取定值，用于调试，这样可以节约调试时间，从而提高效率。 本次课程设计不但让我又学到了一些知识，而且也提高了我的综合能力。使我在各个方面都得到了锻炼，以后有这样的机会一定会更加的很

15、好利用，它不仅可以提高学习的针对性而且可以很好的锻炼动手能力以及自己的逻辑设计能力和处理问题的能力，希望在以后这方面的能力会很好的加强。参考文献1数字信号处理（第二版），丁玉美等，西安电子科技大学出版社；2数字信号处理试验指导书王创新、文卉编 长沙理工大学印刷（内部使用）3数字信号处理及其 MATLAB 实现，陈怀琛等译，电子工业出版社；4 MATLAB 及在电子信息课程中的应用，陈怀琛等，电子工业出版社。5数字信号处理 A.V.奥本海姆， R.W.谢弗著，北京：科学出版社6数字信号处理理论、算法与实现（第二版）胡广书编著，北京：电子工业出版社附录：附录A 采样信号MATLAB 程序代码fig

16、ure(1);t=0:0.005:1;%时域A=444.128;p=222.144;x=A*exp(-10)*t).*sin(p*t);plot(t,x);%连续画图title(信号波形201484250216zmd);xlabel(t/s);ylabel(振幅);figure(2);n=0:200;A=444.128;p=222.144;fs=200;%离散信号采样频率x=A*exp(-10)*n/fs).*sin(p*n/fs);stem(n,x);title(采样后的信号);xlabel(hz);ylabel(振幅);figure(3);y=fft(x,128);y1=abs(y);f=

17、100*(0:63)/128;plot(abs(y1(1:128);title(信号频谱);xlabel(hz);ylabel(振幅);figure(4);t1=0:0.01:1;%时域x1=A*exp(-10)*t1).*sin(p*t1);stem(x1);xlabel(hz);ylabel(振幅);title(混频信号频谱);figure(5);y=fft(x1,128);y1=abs(y);f=100*(0:63)/128;plot(y1);title(混频频谱);xlabel(hz);ylabel(振幅);附录B IIR 滤波器MATLAB 程序代码fs=8000;x,fs,bits

18、=wavread(不说再见1.wav);%读取录音文件%sound(x,fs,bits);figure(1);%创建一个窗口subplot(2,1,1);plot(x);%信号图xlabel(s/时间);%x轴ylabel(fudu);%y轴title(原始信号图);y1=fft(x,4096);%快速傅里叶变换f=fs*(0:2047)/4096;%设x=0：512代表一个一维数列0、1、2、3。512，相当于x=0,1,2,3,。,512f=fs*x/1024就是每个x值进行fs*x(i)/1024运算subplot(2,1,2);plot(abs(y1);%频谱图xlabel(hz/频率

19、);%x轴ylabel(fudu);%y轴title(原始信号频谱);axis(0 4500 0 300);%限定画图的范围figure(2);subplot(2,1,1);%加噪t1=(0:length(x)-1)/1024;x1=0.5*sin(2*pi*7000*t1);%产生噪声信号x0=x(:,1);x2=x0+x1;%加噪音%sound(x2,fs,bits);plot(x2);title(加噪信号);wavwrite(x2,fs,E:数字信号处理实习加噪.wav);%将滤波的信号生成wav文件subplot(2,1,2);y2=fft(x2,4096);plot(abs(y2);

20、%频谱图xlabel(hz/频率);%x轴ylabel(fudu);%y轴title(加噪信号频谱);axis(0 4500 0 300);figure(3);Fs=10000;%采样频率wp=600*2/Fs;%通带截止频率ws=800*2/Fs;%阻带截止频率Rp=1;%通带最大衰减Rs=40;%阻带最大衰减n,Wn=ellipord(wp,ws,Rp,Rs,s); %计算阶数 边界频率B,A=ellip(n,Rp,Rs,Wn,low);%计算滤波器系数h,w=freqz(B,A,512);%求数字滤波器的频率响应plot(w*Fs/(2*pi),20*log10(abs(h)/max(a

21、bs(h);axis(0,1200,-65,3);title(数字低通滤波器的幅度响应);xlabel(Hz/频率);%x轴ylabel(DB/分贝数);%y轴grid figure(4);subplot(2,1,1);y3=filter(B,A,x2);%滤波plot(y3);xlabel(t/时间);%x轴ylabel(fudu);%y轴title(滤波后的时域图);subplot(2,1,2);y4=fft(y3,4096);plot(abs(y4);xlabel(hz/频率);%x轴ylabel(fudu);%y轴title(滤波后的频域图);wavwrite(y3,fs,E:数字信号处理实习去噪.wav);axis(0,1200,0,300);%限定画图的范围附录C FIR滤波器MATLAB 程序代码fs=8000;%采样频率%t=0:(size(x)-1);x,fs,bits=wavread(不说再见1.wav);%读取录音文件t1=(0:length(x)-1)/1024;x1=0.5