拉弦音乐滤波去噪――使用脉冲响应不变法设计的切比雪夫I型滤波器

拉弦音乐滤波去噪——使用脉冲响应不变法设计的切比雪夫I型滤波器指导老师：胡双红摘要本课程设计主要内容是使用脉冲响应不变法设计一个切比雪夫I型滤波器，对一段含噪拉弦音乐信号进行滤波去噪处理并根据滤波前后的波形和频谱分析滤波性能。本课程设计仿真平台为MATLAB7.0，开发工具是M语言编程。首先利用网络下载一段拉弦音乐信号，并人为加入一单频噪声，然后对信号进行频谱分析以确定所加噪声频率，并设计滤波器进行滤波去噪处理，最后比较滤波前后的波形和频谱并进行分析。由分析结果可知，滤波后的拉弦音乐信号与原始拉弦信号基本一致，即设计的IIR滤波器能够去除信号中所加单频噪声，达到了设计目的。关键词课程设计；滤波去噪；切比雪夫I型滤波器；脉冲响应不变法；MATLAB1引言本课程设计主要是在网上采集一段8000Hz，8k的单声道拉弦音乐信号，并绘制波形观察其频谱，再在MATLAB平台上利用脉冲响应不变法设计一个切比雪夫I滤波器，将该拉弦信号进行滤波去噪处理。1.1课程设计目的课程设计是教学的最后一个步骤，课程设计有利于基础知识的理解，我们掌握了基础知识和基本技能，但是要真正接触才能真正理解课程的深入部分；还有利于逻辑思维的锻炼，在许多常规学科的日常教学中，我们不难发现这样一个现象，不少学生的思维常常处于混乱的状态，写起作文来前言不搭后语，解起数学题来步骤混乱，这些都是缺乏思维训练的结果，所以我们可以通过实践来分析问题、解决问题、预测目标等目的；同时也有利于与其他学科的整合，例如我们这次的课程设计就要运用MATLAB软件的帮助才能实现；最重要的有利于治学态度的培养，在课程设计中，我们可能经常犯很多小错误，可能要通过好几次的反复修改、调试才能成功，但这种现象会随着学校的深入而慢慢改观。这当中就有一个严谨治学、一丝不苟的科学精神的培养，又有一个不怕失败、百折不饶品格的锻炼。《数字信号处理》课程设计是在学生完成数字信号处理和MATLAB的结合后的基本实验以后开设的。本课程设计的目的是为了让学生综合数字信号处理和MATLAB并实现一个较为完整的小型滤波系统。这一点与验证性的基本实验有本质性的区别。开设课程设计环节的主要目的是通过系统设计、软件仿真、程序安排与调试、写实习报告等步骤，使学生初步掌握工程设计的具体步骤和方法，提高分析问题和解决问题的能力，提高实际应用水平。1.2课程设计的要求（1）滤波器指标必须符合工程设计。（2）设计完后应检查其频率响应曲线是否满足指标。（3）处理结果和分析结论应该一致，而且应符合理论。（4）独立完成课程设计并按要求编写课程设计报告。1.3设计平台MATLAB是美国MathWorks公司生产的一个为科学和工程计算专门设计的交互式大型软件，是一个可以完成各种精确计算和数据处理的、可视化的、强大的计算工具。它集图示和精确计算于一身，在应用数学、物理、化工、机电工程、医学、金融和其他需要进行复杂数值计算的领域得到了广泛应用。它不仅是一个在各类工程设计中便于使用的计算工具，而且也是一个在数学、数值分析和工程计算等课程教学中的优秀的数学工具，在世界各地的高和大型计算机上运行，适用于Windows、UNIX等多种系统平台。MATLAB作为一种科学计算的高级语言之所以受欢迎，就是因为它有丰富的函数资源和工具箱资源，编程人员可以根据自己的需要选择函数，而无需再去编写大量繁琐的程序代码，从而减轻了编程人员的工作负担，被称为第四代编程语言的MATLAB最大的特点就是简洁开放的程序代码和直观实用的开发环境[1]。2设计原理在网上采集一段拉弦音乐信号，绘制波形并观察其频谱，给定相应技术指标，用脉冲响应不变法设计的一个满足指标的切比雪夫I型IIR滤波器，对该拉弦信号进行滤波去噪处理，比较滤波前后的波形和频谱并进行分析。2.1IIR滤波器IIR滤波器具有无限长脉冲响应，因此能够与模拟滤波器相匹敌；一般来说，所有的模拟滤波器都有无限长脉冲响应。因此，IIR滤波器设计的基本方法是利用复值映射将大家熟知的模拟滤波器变换为数字滤波器。这一方法的优势在于各种模拟滤波器设计（AFD）表格和映射在文献中普遍都能获得。这个基本方法称为A/D（模拟-数字）滤波器变换。然而，AFD表格仅对低通滤波器适用，而同时要得到设计其他频率选择性滤波器（高通、带通、带阻等等）。为此，需要对低通滤波器实行频带变换。这种IIR滤波器设计的基本方法存在两种途径[2]：途径1：期望的IIR滤波器途径2：期望的IIR滤波器在MATLAB中采用第1种途径设计IIR滤波器。这些MATLAB函数的直接使用并没有给出有关任何设计方法的细节。因此，本课程设计将研究途径2的设计基本方法，因为它涉及数字域的频带变换。这样，在这种IIR滤波器设计方法将按下列步骤进行：（1）设计模拟低通滤波器。（2）研究并实行滤波器变换以得到数字低通滤波器。（3）研究并实行频带变换以便从数字低通滤波器得到其他数字滤波器[2]。这些途径存在的主要问题是在IIR滤波器的相位特性上没有一点控制能力，因此IIR滤波器设计将仅作为幅度设计对待。2.2切比雪夫I型滤波器切比雪夫I型滤波器是电子滤波器的一种。它的特点是在通带具有等波纹响应。其设计模拟低通滤波器实行频带变换S→S实行滤波器变换S→Z设计模拟低通滤波器实行滤波器变换S→Z实行频带变换Z→Z振幅平方函数如2-1式所示：)(11|)(|222cNaTjH（2-1）式中，N为整数，称为滤波器的阶，ε是通带波纹因子（它与Rp有关）。其中N阶切比雪夫多项式如2-2式所示：xxxxNXTN1)),(cosh(cosh10)),(coscos()(11（2-2）其中cx切比雪夫滤波器存在两种可能的幅度平方响应的形状，当N分别为基数和偶数时如图2-1所示：图2-1切比雪夫I型滤波器的幅度特性从图可以看出以下性质：在x=0(或Ω=0);|Ha(j0)|2=1;N为奇数时;|Ha(j0)|2=1/(1+ε2);N为偶数时在x=1(或Ω=Ωc);|Ha(j1)|2=1/(1+ε2)对所有N而言.对0=x=1(或0=Ω=Ωc)，|Ha(jx)|2在1和1/(1+ε2)之间起伏。对x1(orΩΩc),|Ha(jx)|2单调下降到0在x=Ωr,|Ha(jx)|2=1/(A^2).2.3脉冲响应不变法如果从模拟到数字滤波器我们想要保留脉冲响应的形状，那么就得到一种方法称为脉冲不变响应法的变换方法[2]。脉冲响应不变法是从滤波器的脉冲响应出发，使数字滤波器的单位脉冲响应序h(n)模仿模拟滤波去的冲击响应ha(t)，使h(n)正好等于ha(t)的采样值，如式（2-3）所示：h(n)=ha(nT)（2-3）T为采样周期。如以Ha(s)及H(z)分别表示ha(t)的拉式变换及h(n)的z变换，如式（2-4），（2-5）所示：Ha(s)=L[ha(t)]（2-4）H(z)=Z[h(n)]（2-5）则根据采样序列z变换与模拟信号拉式变换的关系如式（2-6）所示：)2(1|)(mTjsHTzHmaezst（2-6）上式表明，采样脉冲响应不变法将模拟滤波器变换为数字滤波器时，它所完成的S平面到Z平面的变换，正是以前讨论的拉式变换到Z变换的标准变换，即首先对Ha(s)作周期严拓，然后再经过z=est的映射关系映射到Z平面上。脉冲响应不变法变换前后系统函数的关系如式（2-7），（2-8）所示：NkkkapsRsH1)(（2-7）NkTpkzeRzHk111)(（2-8）应当指出，Z=est的映射关系表明，S平面上每一条宽为2pi/T的横带部分，都将重叠地映射到Z平面的整个全部平面上。每一横带的左半部分映射到Z平面单位圆以内，每一横带的右半部分映射到Z平面单位圆以外，jΩ轴映射在单位圆上，但jΩ轴上的每一段2pi/T都应于绕单位圆一周，如下图2-2所示：图2-2脉冲响应不变法的映射关系Z=est的映射关系反映的是Ha(s)的周期严拓与H(z)的关系，而不是Ha(s)b本身与H(z)的关系，因此，使用脉冲响应不变法时，从Ha(s)到H(z)并没有一个由S平面到Z平面的简单代数映射关系，即没有一个s=f(z)的代数关系式。另外，数字滤波器的频响也不是简单的重现模拟滤波器的频响，而是模拟滤波器频响的周期严拓，周期为ΩS=2π/T=2πfs，如式（2-9）所示：)2(1)2(1)(TmwjHTTmjjHTeHmamajw（2-9）3.设计步骤3.1设计流程图拉弦音乐信号滤波去噪——使用脉冲不变响应法设计的切比雪夫I型滤波器的设计流程如图3-1所示：图3-1脉冲响应不变法切比雪夫I对吹管信号去噪流程图3.2拉弦音乐信号的采集在网上下载一段吹管音乐信号，属性设置为8000Hz，8位单声道的音频格式。如图3-2所示：开始拉弦音乐信号的采集（wavread函数），加入单频噪声，画加噪前后时域图进行快速傅里叶变换，并且画加噪前后频谱图设定滤波器性能指标，单频噪声频率fn=2000，通带截止频率fp=fn-200=1800，阻带截止频率fc=fn-50=1950，通带波纹Rp=1，阻带波纹As=16.用脉冲响应不变法法设计切比雪夫I滤波器验证滤波器是否符合指标并进行频谱分析设计好的滤波器进行滤波处理，回放滤波后的拉弦音乐比较滤波前后吹管音乐信号的波形及频谱结束回放语音信号图3-2拉弦音乐信号设置然后在MATLAB软件平台下，利用函数wavread对语音信号进行采样，源程序为：[x,fs,bits]=wavread('e:\yinyue.wav')，记住采样频率和采样点数，MATLAB实现过程如图3-3所示：图3-3拉弦音乐信号的采样频率和编码位数3.3拉弦音乐信号的频谱分析在MATLAB平台上实现频谱分析的源程序如下：N=length(x);%计算信号x的长度t=0:1/fs:(N-1)/fs;%计算时间范围，样本数除以采样频率deltaf=fs/N;%计算频谱的谱线间隔f=0:deltaf:fs/2-deltaf;%计算频谱频率范围X=abs(fft(x));%对原始信号进行fft变换,取幅度谱X=X(1:N/2);%截取前半部分fn=2000%单频噪声x=x(:,1)';y=x+0.1*sin(fn*2*pi*t);%形成加干扰的拉弦信号sound(y,fs,bits);%应该可以明显听出有尖锐的单频啸叫声Y=abs(fft(y));%对加干扰信号进行fft变换,取幅度谱Y=Y(1:N/2);%截取前半部分程序运行结果得到加噪前后信号的时域波形和频率波形，波形图如图3-4所示：00.51-0.2-0.100.10.2时间（单位：s）幅度原始拉弦信号010002000300040000204060频率（单位:HZ）幅度谱原始拉弦信号幅度谱图00.51-0.2-0.100.10.2时间（单位：s）幅度加入单频干扰后的语音信号01000200030004000050100频率（单位：HZ）幅度谱加入干扰后的语音信号幅度谱图3-4加噪前后时域和频域波形由图3-4可以看出频谱图在频率值为2000Hz处，很明显加入了一个单频噪声。加噪后的时域波形图与原时域波形图也有比较明显的差别。3.4切比雪夫I滤波器设计将数字滤波器的设计指标设为通带截止频率fn=2000,fp=fn-200=1800Hz,阻带频率fc=fn-50=1950HZ,通带波纹Rp=1，阻带波纹As=16dB,要求确定H(z)。设计步骤如下：（1）选取某一T。这是任意的，我们选T=1。并确定模拟频率如式3-1所示：Tpp和Tss（3-1）（2）利用设计参数p，s，pA和sA设计一个模拟滤波器)(sHa。此课程设计需要用切比雪夫I型原型滤波器。在MATLAB中用以下语句实现：[cs,ds]=afd_chb1(OmegaP,OmegaS,Rp,As);%计算滤波器阶数和截止频率[b,