语音信号滤波去噪--使用双线性变换法设计的级联型椭圆滤波器

1、第 26 页 共 25 页Su 语音信号滤波去噪-使用双线性变换法设计的级联型椭圆滤波器 语音信号滤波去噪使用双线性变换法设计的级联型椭圆滤波器学生姓名：Su 指导老师： 摘 要 本课程设计主要内容是利用双线性变换法设计一个级联型的椭圆IIR滤波器，对一段含噪语音信号进行滤波去噪处理并根据滤波前后的波形和频谱分析滤波性能。本课程设计仿真平台为MATLAB7.0，开发工具是M语言编程。首先在windows下用录音机工具录制一段语音信号，并人为加入一单频噪声，然后对信号进行频谱分析以确定所加噪声频率，并设计滤波器进行滤波去噪处理，最后比较滤波前后的波形和频谱并进行分析。由分析结果可知，滤波器后的语

2、音信号与原始信号基本一致，即设计的IIR滤波器能够去除信号中所加单频噪声，达到了设计目的。关键词 滤波去噪；IIR椭圆滤波器；双线性变换法；级联型； MATLAB1 引 言信号处理是科学研究和工程技术许多领域都需要进行的一个重要环节，传统上对信号的处理大都采用模拟系统实现。随着人们对信号处理要求的日益提高，以及模拟信号处理中一些不可克服的缺点，对信号的许多处理而采用数字的方法进行。数字信号处理系统无论在性能、可靠性、体积、耗电量、成本等诸多方面都比模拟信号处理系统优越的多，使得许多以往采用模拟信号处理的系统越来越多地被数字处理系统所代替,数字信号处理技术在通信、语音、图像、自动控制、雷达、军事

3、、航空航天、医疗和家用电器等众多领域得到了广泛的应用。滤波器是一种用来消除干扰杂讯的器件，凡是可以使信号中特定的频率成分通过，而极大地衰减或抑制其他频率成分的装置或系统称之为滤波器。在数字信号处理中，数字滤波器设计在电子工程、应用数学和计算机科学领域都是非常重要的内容。设计滤波器的方法有多种，在各种滤波器中，椭圆滤波器相比其他类型的滤波器，在阶数相同的条件下有着最小的通带和阻带波动，它在通带和阻带的波动相同。本课程设计主要解决在含噪情况下对语音信号的滤波去噪处理，处理时采用的是利用双线性变换法设计的级联型的椭圆IIR滤波器。 1.1 课程设计目的本课程设计主要是用录音机采集一段语音信号，利用M

4、ATLAB，加入一个带外单频噪声，使用双线性变换法设计的级联型椭圆滤波器，对该含噪语音信号进行滤波去噪处理，画出级联型滤波器结构图。比较滤波前后的波形和频谱并进行分析，根据结果和学过的理论得出合理的结论，从而加深对所学知识的理解，加强我们的动手操作能力。另一方面也是对课堂所学的理论知识应用于实践中，熟悉和巩固滤波器的设计原理与方法，掌握双线性变换法级联型椭圆滤波器的设计方法和有关滤波器设计的相关经典算法，熟练地利用MATLAB语言进行设计滤波器。通过本次课程设计，可以运用所学的东西解决一些实际的问题，提高自身对所学知识的综合运用能力以及分析和解决问题的能力。1.2 课程设计的要求（1）滤波器指

5、标必须符合工程实际。（2）设计完后应检查其频率响应曲线是否满足指标。（3）处理结果和分析结论应该一致，而且应符合理论。（4）独立完成课程设计并按要求编写课程设计报告书。1.3 设计平台本次课程设计中所用到的软件平台是是由美国mathworks公司发布的MATLAB软件。MATLAB是矩阵实验室（MatrixLaboratory）的简称，它和Mathematical、Maple并称为三大数学软件。除具备卓越的数值计算能力外，它还提供了专业水平的符号计算，文字处理，可视化建模仿真和实时控制等功能。MATLAB可以进行矩阵运算、绘制函数和数据、实现算法、创建用户界面、连接其他编程语言的程序等，主要应

6、用于工程计算、控制设计、信号处理与通讯、图像处理、信号检测、金融建模设计与分析等领域。 MATLAB 的应用范围非常广，包括信号和图像处理、通讯、控制系统设计、测试和测量、财务建模和分析以及计算生物学等众多应用领域。附加的工具箱（单独提供的专用 MATLAB 函数集）扩展了 MATLAB 环境，可以解决这些应用领域内特定类型的问题。现在MATLAB已经成为国际上最流行的科学与工程计算软件工具，它已不仅仅是一个“矩阵实验室”了，它已经成为一种最具有广泛应用前景的全新的计算机高级语言，在世界上很多国家，MATLAB已经成为线性代数、自动控制理论、数理统计、数字信号处理、时间序列分析、动态系统仿真等

7、高级课程的基本教学工具；成为攻读学位的大学生、硕士生、博士生必须掌握的基本技能。在设计研究单位和工业部门，MATLAB也被广泛用于科学研究和解决各种具体问题当中1。2 设计原理2.1 IIR滤波器IIR滤波器，又名“无限脉冲响应滤波器”，或“递归滤波器”。递归滤波器，也就是IIR数字滤波器，顾名思义，具有反馈性，一般认为具有无限的脉冲响应，因此IIR滤波器能够与模拟滤波器相匹敌2。一个IIR滤波器的系统函数给出为 （2.1） 其中和是滤波器系数。不失一般性已假定。如果。N就是这个IIR滤波器的阶。一个IIR滤波器的差分方程表示是1 （2.2）IIR滤波器具有系统的冲击响应是无穷长的、系统函数在

8、有限z平面上和结构上存在递归的特点。一个IIR滤波器可以用三种结构实现，这三种结构包括直接型、级联型和并联型，其中直接型又可分为直接型和直接型两种。级联型中，将系统函数H(z)写成具有实系数的二阶节的乘积 (2.3) 式中，和都是代表实数的二阶节系数。这些二阶节是 (2.4)称为第k个双二阶节。每个级联型结构的滤波器都可以由k个双二阶节用直接型实现，双二阶节结构图如图2-1所示。图2-1 双二阶节结构一般说来，所有的模拟滤波器都有无限长的脉冲响应，因此IIR滤波器的设计基本方法是利用复值映射将模拟滤波器变换为数字滤波器。目前，已有许多成熟的模拟滤波器应用于信号处理中，如巴特沃兹滤波器、切比雪夫

9、（型和型）滤波器、椭圆滤波器等，我们将这些各有特色的模拟滤波器成为原型滤波器。在设计IIR滤波器时，可以直接借用这些典型的模拟滤波器进行转换，这种方法称为A/D（模拟数字）滤波器变换。因为AFD表格仅对低通滤波器适用，但在实际应用中也会需要设计其他频率选择性滤波器（高通、带通、带阻滤波器等等）3。为此需要对低通滤波器实行频带变换，基本方法存在两种途径，如图2-2所示。途径1：途径2：图2-2 IIR滤波器的两种设计途径在MATLAB中采用第一种途径设计IIR滤波器。2.2 椭圆低通滤波器椭圆滤波器又称考尔滤波器，是通带和阻带等波纹的一种滤波器。在阶数相同的条件下，与巴特沃兹滤波器和切比雪夫滤波

10、器相比，椭圆滤波器有着最小的通带和阻带波动。就这点而言, 椭圆滤波器是最优的。从传递函数来看，椭圆函数滤波器在有限频率上既有零点又有极点。极零点在通带内产生等波纹, 阻带内的有限传输零点减少了过渡区, 可获得极为陡峭的衰减曲线。但它陡峭的过渡带特性是用通带和阻带的起伏为代价来换取的，并且在通带和阻带的波动相同，这一点区别于在通带和阻带都平坦的巴特沃斯滤波器，以及通带平坦、阻带等波纹或是阻带平坦、通带等波纹的切比雪夫滤波器。椭圆滤波器的传输函数是一种较复杂的逼近函数, 利用传统的设计方法进行电路网络综合要进行繁琐的计算, 还要根据计算结果进行查表, 整个设计, 调整都十分困难和繁琐。为此，常常需

11、要用一些程序和表格来设计椭圆滤波器，用MATLAB设计它们可以大大简化设计过程。椭圆滤波器的幅度平方响应给出为 (2.5) 式中是通带波纹(它与有关) ,是N阶雅可比(Jacobian)椭圆函数，N是阶次。对于奇数和偶数N，可得到不同的椭圆滤波器的幅度响应特性图，具体如下图2-3所示。图2-3 椭圆滤波器的幅度特性其中椭圆滤波器阶次N的计算公式给出为 (2.6)式中，和。当确定了阶次N，便可设计椭圆滤波器了。MATLAB中提供了丰富的函数，我们可将MATLAB中一个afd_elip函数用于设计一个已知设计指标的模拟低通滤波器3。2.3 双线性变换法双线性变换法保留的是从模拟到数字域的系统函数的

12、表示，它是一种稳定的变换。双线性变换的主要优点是S平面与Z平面是一种一一对应的单值关系，这种影射是最好的变换方法，具体变换公式为 (2.7)这里T为采样间隔。式(2.7)亦可转换为 (2.8)从式(2.8)可以看出，只要其中一个变量固定，这在每个变量上都是线性的，也就是说在s和z上是双线性的。在式(2.7)的制约下的复平面影射如图2-4所示。从图中可以看出，图2-4 双线性变换法中的复平面影射第一步先将整个S平面压缩影射到S1平面的之间的一条横带里；第二步再通过标准变换关系将此横带转换到整个Z平面上去3。由图2-4可以看出，S1平面的整个左半平面影射到单位圆内，S1平面的虚轴(整个j)对应于Z

13、平面单位圆的一周，S1平面的=0处对应于Z平面的=0处，所以双线性变换不存在混迭效应。为了将S平面的j轴压缩到S1平面j1轴上的到段上，可以通过如下正切变换实现 (2.9)由此可求得1。当1由经过0 变化到时，由-经过0变化到，也即影射了整个j轴。将式(2.9)写成 (2.10)令s=j,s1= j1,则可得到 (2.11)最后通过标准变换关系将S1平面映射到Z平面上去，便可得到式(2.7)所示的S平面和Z平面的单值影射关系。这便是双线性变换4。3 设计步骤3.1 设计流程图语音信号滤波去噪使用双线性变换法设计的椭圆滤波器的设计流程如图3-1所示。利用录音机录制一段语音信号对语音信号进行频谱分

14、析，画出时域和频域波形图使用双线性变换法设计级联型椭圆滤波器，用设计好的椭圆滤波器对语音信号进行滤波，并对比分析滤波前后语音信号波形 开始结束在MATLAB软件平台下，用wavread函数对该语音信号进行采样，并加入单频噪声在MATLB平台下用sound函数回放语音信号绘制级联型椭圆滤波器结构图图3-1 设计流程图3.2 录制语音信号及加噪利用Windows平台下的录音机，录制一段语音信号“大家好，我是*”，时间在1-2s左右。录制时，单击如图3-2所示的录音机界面的圆圈，便可录制语音信号。录音完毕，点击如图3-2所示的方形，录制结束。图3-2 录音机界面通过点击录音界界面上的文件，选择下拉菜

15、单中的属性，在弹出的对话框中点击立即转换，设置录音的格式为PCM编码，属性8kHz,8位，单声道7kB/s,这样的设置在相同的时间下数据最少，处理时间也就最短，如图3-3所示。图3-3 格式设置对话框格式设置好后，点击文件，在下拉菜单中选择另存为su.wav，并将其放入到MATLAB安装盘下的work文件夹中。然后在MATLAB软件平台下，利用函数wavread将语音信号读入workspace，并对其进行频谱分析，画出时域和频域波形图。函数wavreadd的调用格式为x,fs,nbits=wavread(file)其中x为采样值；fs为采样频率，单位为Hz；bits表示采样位数；file表示待

16、读取的文件名和全路径。采集完成后在信号中加入一频率为2200Hz的单频噪声，然后设计滤波器将信号中加入的单频噪声滤除，还原出原始信号。下面是语音信号采集及加入单频噪声的部分程序。x,fs,bits=wavread('e:su.wav'); % 输入参数为文件的全路径和文件名N=length(x); % 计算信号x的长度fn=2200; % 单频噪声频率t=0:1/fs:(N-1)/fs; % 计算时间范围，样本数除以采样频率x=x' y=x+0.1*sin(fn*2*pi*t);wavwrite(y,fs,bits,'su1'); %将加噪后的语音信号写

17、成wav文件，并命名为su1在MATLAB平台上执行，可以通过sound (y,fs,bits) 回放该语音信号，加入单频干扰后，在明显听出“大家好，我是*”的声音中伴随较尖锐的干扰啸叫声。加入单频干扰信号后，画出干扰前后语音信号的时域波形，然后进行快速傅里叶变换，得到信号的频谱特性图。快速傅里叶变换算法FFT计算DFT的函数为fft，其调用格式为X=fft(x,n)其中x为被变换的时域序列向量，n是DFT的变换区间长度，当n大于x的长度时，fft函数自动在x后面补零。当n小于xn的长度时，fft函数计算x的前n个元素，忽略其后面的元素。X=abs(fft(x); Y=abs(fft(y);

18、% 对原始信号和加噪信号进行fft变换,取幅度谱X=X(1:N/2); Y=Y(1:N/2); % 截取前半部分deltaf=fs/N; % 计算频谱的谱线间隔f=0:deltaf:fs/2-deltaf; % 计算频谱频率范围在MATLAB平台中执行后，得到原始语音信号和加入干扰后的语音信号的时域和频域波形图如图3-4所示。由图可以看出，加入单频干扰信号后，在时域波形图中，加入干扰后的语音信号的幅度明显比原始语音信号的大；在频域波形图中，在2200Hz处明显出现一尖锐单峰，与干扰频率刚好一致。图3-4 加噪前后语音信号时域和频域波形比较 3.3 滤波器设计滤波器的设计指标是通带截至频率fp=

19、2000Hz，阻带截至频率fc=2150Hz，通带波纹为1dB，阻带波纹为60dB。这里用双线性变换法设计一个满足上述指标的椭圆低通滤波器。MATLAB中提供了一个称为z,p,k=ellipap(N,Rp,AS)的函数用于设计一个阶次为N，通带波纹为Rp和阻带波纹为As的归一化椭圆模拟原型滤波器。在这里用一个不是MATLAB系统自带的称为afd_elip的函数用于设计一个已知设计指标时的模拟椭圆低通滤波器，然后再利用双线性进行模数转换得到所要求的IIR滤波器。afd_elip的调用格式为b,a=afd_elip(Wp,Ws,Rp,As)其中b为滤波器系统函数的分子，a为滤波器系统函数的分母，W

20、p为进行了预畸变的通带截至频率，Ws为进行了预畸变的阻带截至频率。设计滤波器的部分主要程序如下fn=2200; fs=8000; Rp=1;As=60;%低通滤波器设计指标fp=fn-200;fc=fn-50; %定义通带和阻带截止频率wp=fp/fs*2*pi;ws=fc/fs*2*pi; %计算对应的数字频率T=1; %定义采样间隔Omegap=2/T*tan(wp/2);Omegas=2/T*tan(ws/2); %截止频率预畸变c,d=afd_elip(Omegap,Omegas,Rp,As) % 计算滤波器系统函数分子分母系数b,a=bilinear(c,d,1/T) % 双线性变换

21、得到数字滤波器系统函数分子分母系数db,mag,pha,grd,w=freqz_m(b,a); %计算椭圆低通滤波器的幅度，幅度响应，脉冲响应和相位响应等数据图，验证滤波器是否达到指定性能delta=1,zeros(1,99);ha=filter(b,a,delta);%计算脉冲响应得到的椭圆滤波器的幅度、幅度响应、相位响应和脉冲响应图如图3-5所示。图3-5 滤波器的幅度、幅度响应、相位响应和脉冲响应图图3-5中，第一幅图中的两根横线上面的横线是Rp，下面的一根横线表示As；两根竖线中，左边的一根是Wp，右边的一根表示Ws，把这四根线通过函数line标出，可以更直观得判断设计是否达标。由图中

22、可看出，设计的滤波器为理想的滤波器。3.4 信号滤波处理用设计好了的滤波器对采样信号进行滤波处理，将加入的干扰信号滤除，还原出原始语音信号，并将还原出的语音信号用wavwrite函数写成wav文件。在MATLAB，IIR滤波器利用函数filter对信号进行滤波，其调用格式为y_fil=filter(b,a,y)其中b为进行双线性变换后得到的数字滤波器系统函数的分子系数，a为进行双线性变换后得到的数字滤波器系统函数的分母系数。滤波器对信号进行滤波的主要程序如下y_fil=filter(b,a,y); % 用设计好的IIR滤波器对信号进行滤波处理Y_fil=fft(y_fil);Y_fil=Y_f

23、il(1:N/2); % 计算频谱取前一半wavwrite(y_fil,fs,bits,'su2');%将滤波后的语音信号写成wav文件得到原始信号、加噪信号和滤波去噪信号的时域和频谱图如图3-6所示。由图3-6可以图3-6 滤波前后信号的时域和频谱图明显地看出，无论时域波形图还是频谱图，滤波后，加入的噪声被有效地抑制，达到了滤波去噪的效果，能够较好地还原出原始信号。3.5 结果分析在MATLAB中，经sound函数对经过椭圆滤波器滤波之后的语音信号进行回放，可以听出滤波之后的信号比原始信号还清晰一些，因为周围环境的声音被滤除了一些，加入的干扰信号也被滤除了。Sound函数的调

24、用格式为Sound(x,fs,bits)其中x表示语音信号文件，fs表示采样频率，bits表示文件编码时用到的样本编码位数。sound(x,fs,bits); %按指定的采样率和每样本编码位数回放sound(y,fs,bits); %可以明显听出有尖锐的单频啸叫声sound(y_fil,fs,bits); % 按指定的采样率和每样本编码位数回放在MATLAB平台上，利用前面的程序，分别插入上面的三个语句并运行程序，我们可以依次听到清晰的原始录音、伴有尖锐的啸叫声的语音和清晰的与原始语音差别不大的滤波后的语音，所得的结果证明了椭圆滤波器设计成功。下面与两个同学的使用双线性变换设计的级联型的切比雪

25、夫型滤波器和巴特沃兹滤波器的滤波效果进行比较。加入干扰后的语音信号经这三位同学的滤波器滤波后，得到的原始信号、加噪信号和滤波去噪信号的时域和频谱图如图3-7所示。(a) 切比雪夫型滤波器滤 (b)巴特沃兹滤波器图3-7 双线性变换设计的不同滤波器滤波效果图将图3-6和图3-7对比就可以看出，相同的设计指标fp=2000Hz， fc=2150Hz，Rp=1dB，As=60dB，巴特沃兹的滤波效果相对于其他三种滤波器来说滤波效果是最差的，而且需要的阶数也最高，为65阶。但从图中，很难看出椭圆滤波器和切比雪夫型滤波器滤这两种滤波器的优劣，我们可以通过MATLAB的运行结果可以知道，契比雪夫、型滤波器

26、所需要的阶数都是17阶，说明它们具有相同的性能；而椭圆滤波器具有最小的阶数8，说明它的性能在这个意义上是最优的。3.6 滤波器结构设计利用前面的3.3节滤波器的设计中经过双线性变换得到的数字滤波器系统函数的分子分母，调用函数C,B,A=dir2cas(b,a)，将IIR滤波器的直接型转换成级联型，运行的结果如下C,B,A=dir2cas(b,a)C = 0.0249B = 1.0000 1.7552 1.0000 1.0000 0.8262 1.0000 1.0000 0.3594 1.0000 1.0000 0.2106 1.0000A = 1.0000 -0.0013 0.9731 1.0

27、000 -0.1267 0.8885 1.0000 -0.4776 0.6957 1.0000 -1.0324 0.3982根据以上参数，将滤波器结构用visio绘图软件画好，并填上系数。画滤波器结构图的具体操作步骤首先是打开visio绘图软件，点击文件新建新建绘图，进入绘图的界面。点击绘图对话框中的线条工具，可以跟据自己的喜好选择绘制结构图线条的大小和线段，在绘图界面中间白色的网格中便可以绘制滤波器的结构图了。画好的级联型椭圆滤波器的结构图如图3-8所示。图3-8 级联型椭圆滤波器的结构图4 出现的问题及解决方法从一开始接触MATLAB已经有一年多的时间了，而且已经运用MATLB软件做过了一

28、次课程设计，但是到真正再次操作起来还是感觉很生疏，对MATLAB的提供的函数的调用还是不熟悉。所以，这次又是在MATLAB平台上操作的课程设计，从一开始就遇到了许多问题，但在老师和同学的帮助下，当然也在自己的努力下，问题都一一得到了解决。刚开始，因为没有电脑的录音机的麦克风坏掉了，一连换了好多台都录不进去，最后在一个麦克风是好的的同学那里好不容易录了一段音，结果因为背景有太多杂音，运行后得到的语音信号用自己设计的滤波器滤波效果不是很好，最后在下课后机房没有什么人的时候又重录了一段，这次的录音质量终于好了一些。还有就是在画语音信号的时域和频域波形图的时候，没有对波形的坐标进行修正，观察比较的效果

29、图不是很好。在老师的指导下，用axis函数进行修正，便可以得到比较好看的信号时域和频谱波形图，下面图4-1是修正前和修正后的波形图。图4-1 用函数axis修正前后的波形图在调用fld_elip等函数时，因为没有把老师提供的function文件夹放进到MATLAB工作下的work文件夹中，在运行程序的时候老是出错，提示调用的函数不存在。最后问同学，才知道系统中没有这些函数，调用时要把这些函数的脚本放到MATLAB工作的路径下。设计滤波器时，因为参数没有设置好，得到的滤波器滤波效果很不好，最后通过修改阻带波纹，改善了滤波器的滤波效果。在设计滤波器的时候还遇到了一个问题，就是在画滤波器的幅度谱图的

30、时候，需要在Wp，Ws，Rp和As处画线以便更直观判断设计得滤波器是否达标。因为当时还不太理解标注的原理，结果因为第一次用的是下面的程序xx=0 1;yy=0 0;line(xx, yy,'Color','r','LineWidth',2,'LineStyle','-');xx=0 1;yy=-4 -4;line(xx, yy,'Color','r','LineWidth',2,'LineStyle','-');xx=wp/pi wp/

31、pi;yy=-90 90;line(xx, yy,'Color','r','LineWidth',2,'LineStyle','-');xx=ws/pi ws/pi;yy=-90 90;line(xx, yy,'Color','r','LineWidth',2,'LineStyle','-');运行得到的结果得到如下图4-2所示的滤波器波形图，这个图是不正确的。图4-2 -As标错的滤波器波形图后来在老师批改了作业，并指出了自己的问题，

32、把上面的程序改为X_l=0,0,wp/pi,ws/pi;1,1,wp/pi,ws/pi;Y_l=-As,-Rp,-90,-90;-As,-Rp,3,3;% 在wp，ws，Rp，As处画线以更直观判断设计是否达标，每列参数是每个线条的端点坐标line(X_l,Y_l,'Color','r','LineWidth',2,'LineStyle','-') % 添加线条，红色，线宽为2运行得到正确的滤波器波形图如下图4-3所示。图4-3 As标错的滤波器波形图5 结束语在这短短的两个星期的数字信号处理课程设计中，因为一直在

33、不停的做老师布置的任务，感觉时间过得飞快而且充实。这次课程设计，因为在第一个星期内就得把所有设计任务完成，同时还要把论文初稿完成提交给老师批改，所以其实也不能算是两个星期。但是，时间虽然短，但感觉比自己一个学期学到的关于数字信号处理的知识还要多，而且动手能力也得到了大大的提高。首先，因为之前上数字信号处理的课程，都是只听课，而没有动手操作，即使后来有实验课也因为某种原因没有将自己所学的东西在MATLAB平台上进行实践操作，致使到做这个课程设计之前对数字信号处理的意义都是不太明朗的。如今，因为这个课程设计，我懂得了数字信号处理的许多用途，对模拟滤波器和数字滤波及它们之间的转换也有了更深入的了解，

34、而且可以自己动手设计滤波器了。在这几次课程设计中，现在再次体会到了平时学习理论知识与实践之间的差距。学习是一个过程，如果你不会把自己所学的东西应用于实践，那么这些知识永远都不会是你的。而我平时学习，总没有深入地去研究，所以到实际去操作时就会碰到各种各样的问题，然后被卡住，实践就进行不下去了。这次课程设计我的课题是语音信号滤波去噪使用双线性变换法设计级联型的椭圆滤波器，开始就因为自己对知识理解得不够透彻，遇到许多问题，但最后在老师和同学的帮助下都解决了，成功地完成了滤波器的设计任务。所以，学习就应扎扎扎实实地去学，不仅要往广度学，还要深入研究，而且要加强自己的动手操作能力。其次，在这两周的课程设

35、计中，我学到了很多东西，了解到了MATLAB强大的功能。同时，也懂得了要熟练的操作MARTLAB，就必须掌握编程所需的各种函数和语句。MATLAB的函数库与编程语句量是非常大的，要完全掌握是很困难的，所以需要我们平时自己多练习，慢慢的理解MATLAB的编程规律，一点一点地积累。最后，在这两周的课程设计中，我每天的日子都过得很忙碌，每天写程序，编译，纠错，上网查资料，拿自己设计的课题所学要设计的滤波器反复地研究。虽然有时候在自己编译的程序老是报错的时候会很苦恼，上网查资料查得自己头晕眼花的，写报告的那两天常常写凌晨困得不行的时候才爬到床上去睡觉，等老师批改完之后还要修改，心里会很郁闷，可是等到自

36、己的程序编译成功，自己设计的滤波器成功的那一刹那，心情又是极度喜悦的，这样的日子让自己感到过得很充实。在这段时间中，同学和同学之间，老师和同学之间的交流也是最多的。一群同学之间在一起共同讨论自己所遇到的问题，老师在一边耐心地指导的场面，感觉很和谐。这次课程设计让我充分认识到同学间合作的重要性，你自己的一个小小的问题怎么也发现不了，找同学帮忙，他可能一下子就可以帮你指出来了。在这个过程中我认识到坚持、耐心和细心等品质的重要性，这对我今后的无论是学习、生活，抑或是工作都有很大的帮助。总之，这次课程设计使我获益匪浅。最后，在此我由衷地感谢各位帮助过我的老师和同学们！ 6 参考文献1 张威. MATL

37、ANB基础与编程入门.2.西安.西安电子科技大学出版社,20082 王靖. 李永全.数字椭圆滤波器的Matlab设计与实现.现代电子技术，2007,6 : 2453 维纳K恩格尔，约翰G普罗克斯.刘树棠.陈志刚.数字信号处理.3.西安.西安交通大学出版社，20134 刘凤举.吴简彤.刘力.基于双线性变换法的IIR数字滤波器设计与matlab仿真.自动化技术与应用，2008，27(8)：4445 附录1：语音信号滤波去噪设计源程序清单%程序名称：su.m% 程序功能：使用双线性变换法设计的级联型椭圆滤波器对语音信号滤波去噪。% 程序作者：*% 最后修改日期：2013-3-7x,fs,bits=w

38、avread('e:su.wav'); % 输入参数为文件的全路径和文件名，输出的第一个参数是每个样本的值，fs是生成该波形文件时的采样率，bits是波形文件每样本的编码位数。sound(x,fs,bits); % 按指定的采样率和每样本编码位数回放N=length(x); % 计算信号x的长度fn=2200; % 单频噪声频率，此参数可改为2200t=0:1/fs:(N-1)/fs; % 计算时间范围，样本数除以采样频率x=x' y=x+0.1*sin(fn*2*pi*t);%sound(y,fs,bits); 应该可以明显听出有尖锐的单频啸叫声wavwrite(y,

39、fs,bits,'su1');%将加噪后的语音信号写成wav文件，并命名为su1X=abs(fft(x); Y=abs(fft(y); % 对原始信号和加噪信号进行fft变换,取幅度谱X=X(1:N/2); Y=Y(1:N/2); % 截取前半部分deltaf=fs/N; % 计算频谱的谱线间隔f=0:deltaf:fs/2-deltaf; % 计算频谱频率范围figure(1);subplot(2,2,1),plot(t,x);grid;axis(0,1.7,-0.2,0.2);xlabel('时间（单位：s）'),ylabel('幅度'),t

40、itle('原始语音信号');%原始语音信号时间幅度图subplot(2,2,2),plot(f,X);grid;axis(0,4000,0,50);xlabel('频率（单位：Hz）'),ylabel('幅度谱'),title('语音信号幅度谱图');%语音信号幅度谱图subplot(2,2,3),plot(t,y);grid;axis(0,1.7,-0.2,0.2);xlabel('时间（单位：s）'),ylabel('幅度'),title('加入单频干扰后的语音信号');%加入

41、单频干扰后的语音信号subplot(2,2,4),plot(f,Y);grid;axis(0,4000,0,50);xlabel('频率（单位：Hz）'),ylabel('幅度谱'),title('加入单频干扰后的语音信号幅度谱图');%加入单频干扰后的语音信号幅度谱图fp=fn-200;fc=fn-50; %定义通带和阻带截止频率Rp=1;As=60; % 定义通带波纹和阻带衰减wp=fp/fs*2*pi;ws=fc/fs*2*pi; %计算对应的数字频率T=1; %定义采样间隔% 第一步 预畸变Omegap=2/T*tan(wp/2);Ome

42、gas=2/T*tan(ws/2); %截止频率预畸变% 第二步 根据指标设计模拟滤波器c,d=afd_elip(Omegap,Omegas,Rp,As) % 计算滤波器系统函数分子分母系数% 第三步 模数滤波器变换b,a=bilinear(c,d,1/T) % 双线性变换得到数字滤波器系统函数分子分母系数%验证db,mag,pha,grd,w=freqz_m(b,a); %计算椭圆低通滤波器的幅度，幅度响应，群延时和相位响应等数据图，验证滤波器是否达到指定性能delta=1,zeros(1,99);ha=filter(b,a,delta);%计算脉冲响应figure(2);subplot(2

43、,2,1);plot(w/pi,db);grid;axis(0,1,-90,3);xlabel('w/pi');ylabel('dB');title('滤波器幅度响应图');X_l=0,0,wp/pi,ws/pi;1,1,wp/pi,ws/pi;Y_l=-As,-Rp,-90,-90;-As,-Rp,3,3; % 在wp，ws，Rp，As处画线以更直观判断设计是否达标，每列参数是每个线条的端点坐标line(X_l,Y_l,'Color','r','LineWidth',2,'LineStyl

44、e','-') % 添加线条，红色，线宽为2subplot(2,2,2);plot(w/pi,mag);grid;axis(0,1,0,1.1);xlabel('w/pi');ylabel('幅度mag');title('滤波器幅度响应图');subplot(2,2,3);plot(w/pi,pha);grid;axis(0,1,-4,4);xlabel('w/pi');ylabel('相位pha');title('滤波器相位响应图');subplot(2,2,4);plot

45、(ha);grid;axis(0,100,-0.12,0.12);xlabel('n');ylabel('h(n)');title('滤波器脉冲响应图');y_fil=filter(b,a,y); % IIR滤波器对信号signal进行滤波处理，其中b,a为上面设计好的滤波器参数。Y_fil=fft(y_fil);Y_fil=Y_fil(1:N/2); % 计算频谱取前一半sound(y_fil,fs,bits); %按指定的采样率和每样本编码位数回放wavwrite(y_fil,fs,bits,'su2');%将加噪后的语音信号

46、写成wav文件，并命名为su2figure(3)subplot(3,2,1);plot(t,x);grid;axis(0,1.7,-0.2,0.2);ylabel('幅度');xlabel('时间t'),title('原始语音信号时间x');subplot(3,2,2);plot(f,X);grid;axis(0,4000,0,50);ylabel('幅度');xlabel('频率f');title('原始语音信号幅度频谱X');subplot(3,2,3);plot(t,y);grid;axis(

47、0,1.7,-0.2,0.2);ylabel('幅度');xlabel('时间t');title('加入干扰语音信号时间x1');subplot(3,2,4);plot(f,Y);grid;axis(0,4000,0,50);ylabel('幅度');xlabel('频率f');title('加入干扰语音信号幅度谱X1');subplot(3,2,5);plot(t,y_fil);grid;axis(0,1.7,-0.2,0.2);ylabel('幅度');xlabel('时间t'),title('滤波后语音信号时间y');subplot(3,2