《双音频信号滤波去噪——使用双线性变换法设计的切比雪夫I型滤波器》第1页共19页双音频信号滤波去噪——使用双线性变换法设计的切比雪夫I型滤波器摘要本课程设计是采用双线性变换法设计的切比雪夫I型滤波器对双音频信号滤波去噪。在网上下载一段双音频信号,在MATLAB集成环境下,首先用wavread函数求出双音频信号的相关参数,对双音频信号进行读取和加噪;然后再给定相应技术指标,设计一个满足指标的切比雪夫I型滤波器,对该双音频信号进行滤波去噪处理,并绘制对比图,比较滤波前后的波形和频谱并进行分析;最后通过回放双音频信号,对比滤波前后的信号变换。本课程设计成功的对双音频信号进行滤波去噪,初步完成了设计指标。关键词双音频信号;滤波设计;MATLAB;切比雪夫I型滤波器1引言用麦克风采集一段8000Hz,8k的双音频信号,绘制波形并观察其频谱,给定通带截止频率为2000Hz,阻带截止频率为2150Hz,通带波纹为1dB,阻带波纹为35dB,用双线性变换法设计的一个满足上述指标的切比雪夫I型IIR滤波器,对该双音频信号进行滤波去噪处理。1.1课程设计目的《数字信号处理》课程设计是在学生完成数字信号处理和MATLAB的结合后的基本实验以后开设的。本课程设计的目的是为了让学生综合数字信号处理和MATLAB并实现一个较为完整的小型滤波系统。这一点与验证性的基本实验有本质性的区别。开设课程设计环节的主要目的是通过系统设计、软件仿真、程序安排与调试、写实习报告等步骤,使学生初步掌握工程设计的具体步骤和方法,提高分析问题和解决问题的能力,提高实际应用水平。《双音频信号滤波去噪——使用双线性变换法设计的切比雪夫I型滤波器》第2页共19页1.2课程设计的要求(1)学会MATLAB的使用,掌握MATLAB的程序设计方法;(2)滤波器指标必须符合工程实际,根据模拟滤波器的性能指标,确定数字滤波器指标;(3)采用双线性变换法,设计满足上述性能指标要求的ChebyshevI型数字低通滤波器;(4)设计完后应检查其频率响应曲线是否满足指标;(5)处理结果和分析结论应该一致,而且应符合理论;(6)独立完成课程设计并按要求编写课程设计报告书;1.3设计平台本次课程设计是在MATLAB软件平台上进行的。MATLAB是矩阵实验室(MATRIXLABORATORY)的简称,是美国MATHWORKS公司推出的具有强大数值分析、矩阵运算、图形绘制和数据处理等功能的软件,现已广泛应用到教学、科研、工程设计等领域[2]。随着MATLAB软件信号处理工具箱的推出,MATLAB已成为信息处理,特别是数字信号处理DSP应用中分析和设计的主要工具。就MATLAB信号处理中的滤波器设计而言,在很大程度上能快速有效地实现滤波器的分析、设计及仿真,大大节约了设计时间,相对传统设计而言,简化了滤波器设计的难度。2设计原理用麦克风采集一段双音频信号,绘制波形并观察其频谱,给定相应技术指标,用双线性变换法设计的切比雪夫I型IIR滤波器,对该双音频信号进行滤波去噪处理,比较滤波前后的波形和频谱并进行分析。2.1IIR滤波器从离散时间来看,若系统的单位抽样(冲激)响应延伸到无穷长,称之为“无限长单位冲激响应系统”,简称为IIR系统。无限长单位冲激响应(IIR)滤波器有以下几个特点:(1)系统的单位冲激响应h(n)是无限长;《双音频信号滤波去噪——使用双线性变换法设计的切比雪夫I型滤波器》第3页共19页(2)系统函数H(z)在有限z平面(0z∞);(3)结构上存在着输出到输入的反馈,也就是结构上是递归型的。IIR滤波器采用递归型结构,即结构上带有反馈环路。同一种系统函数H(z)可以有多种不同的结构,基本网络结构有直接Ⅰ型、直接Ⅱ型、级联型、并联型四种,都具有反馈回路。同时,IIR数字滤波器在设计上可以借助成熟的模拟滤波器的成果,巴特沃斯(Butterworth)滤波器、切比雪夫(Chebyshev)滤波器、椭圆(Cauer)滤波器、贝塞尔(Bessel)滤波器等,这些典型的滤波器各有特点。有现成的设计数据或图表可查,在设计一个IIR数字滤波器时,我们根据指标先写出模拟滤波器的公式,然后通过一定的变换,将模拟滤波器的公式转换成数字滤波器的公式。2.2切比雪夫I型滤器切比雪夫滤波器(又译车比雪夫滤波器)是在通带或阻带上频率响应幅度等波纹波动的滤波器。在通带波动的为“I型切比雪夫滤波器”,在阻带波动的为“II型切比雪夫滤波器”。切比雪夫滤波器在过渡带比巴特沃斯滤波器的衰减快,但频率响应的幅频特性不如后者平坦。切比雪夫滤波器和理想滤波器的频率响应曲线之间的误差最小,但是在通频带内存在幅度波动。这种滤波器来自切比雪夫多项式,因此得名,用以纪念俄罗斯数学家巴夫尼提·列波维其·切比雪夫巴特沃兹滤波器在通带内幅度特性是单调下降的,如果阶次一定,则在靠近截止处,幅度下降很多,或者说,为了使通带内的衰减足够小,需要的阶次N很高,为了克服这一缺点,采用切比雪夫多项式来逼近所希望的。切比雪夫滤波器的在通带范围内是等幅起伏的,所以在同样的通常内衰减要求下,其阶数较巴特沃兹滤波器要小。切比雪夫滤波器的振幅平方函数为(2-1)式中Ωc—有效通带截止频率《双音频信号滤波去噪——使用双线性变换法设计的切比雪夫I型滤波器》第4页共19页—与通带波纹有关的参量,大,波纹大01VN(x)—N阶切比雪夫多项式(2-2)|x|≤1时,|VN(x)|≤1|x|1时,|x|↗,VN(x)↗切比雪夫滤波器的振幅平方特性如图所示,通带内的变化范围为1(max)→211(min)时,|x|1,随↗,→0(迅速趋于零)当=0时,(2-3)N为偶数,cos2()=1,得到min,,(2-4)N为奇数,cos2(,得到max,(2-5)切比雪夫滤波器的振幅平方特性如图2.1所示。《双音频信号滤波去噪——使用双线性变换法设计的切比雪夫I型滤波器》第5页共19页图2.1切比雪夫滤波器的振幅平方特性2.3双线性变换法双线性变换法是使数字信号滤波器的频率响应与模拟滤波器的频率响应相似的一种变换方法。为了客服多值映射这一缺点,我们首先把整个s平面压缩变换到某一中介的s1平面的一条横带里(宽度为T2,即从T到T),其次再通过上面讨论过的标准变换关系Tsez1将此横带变换到这个z平面上去,这样就使s平面与z平面式一一对应的关系,消除了多值变换性,也就消除了频谱混叠现象。将s平面整个j轴压缩变换到s1平面j轴上的T到T一段,可以采用以下变换关系:)2tan(1T(2-6)这样,变换到T1,0变到01可将(6)式写成22221111TjTjTjTjeeeej(2-7)解析延拓到整个s平面和s1平面,令sj,11sj,则得TsTsTjTjTjTjeeTstheeeej11111111212222(2-8)再将s1平面通过以下标准变换关系映射到z平面:《双音频信号滤波去噪——使用双线性变换法设计的切比雪夫I型滤波器》第6页共19页Tsez1(2-9)从而得到s平面和z平面的单值映射的关系为1111zzs(2-10)ssz11(2-11)一般来说,为了使模拟滤波器的某一频率与数字滤波器的任一频率有对应的关系,可以引入待定常数c,使(6)式和(7)式变换成)2tan(1Tc(2-12)TsTseecTscths111121(2-13)仍将Tsez1代入(13)式,可得1111zzcs(2-14)scscz(2-15)(14)式和(15)式是s平面与z平面之间的单值映射关系,这种变换就称为双线性变+换。3设计步骤3.1设计流程图双音频信号滤波去噪——使用双线性变换法设计的切比雪夫I型滤波器的设计流程如图3.1所示:《双音频信号滤波去噪——使用双线性变换法设计的切比雪夫I型滤波器》第7页共19页图3.1双线性变换法切比雪夫I型滤波器对双音频信号去噪流程图3.2双音频信号的采集点击windows系统桌面的“开始”按钮,点击开始菜单栏里的“附件”,选择“录音机”选项,点击录音机“文件”选项,进入“声音选定”设置,把属性一栏设置成“8000Hz,8位,单声道,7KB/秒”(见图3.2)。点击确定,然后开始双音频信号的采集,采集时间为2—3秒左右为最佳。采集的声音文件以“.wav”格式存储(见图3.3)。双音频信号的采集(wavread函数),画时域图快速傅里叶变换,并且画频谱图设定滤波器性能指标,通带截止频率fb=2000,阻带截止频率fc=2150,通带波纹Rp=1,阻带波纹As=35双线性变换法设计切比雪夫I型滤波器验证并进行频谱分析设计好的滤波器进行滤波处理比较滤波前后双音频信号的波形及频谱回放音频信号结束开始《双音频信号滤波去噪——使用双线性变换法设计的切比雪夫I型滤波器》第8页共19页图3.2采集声音的参数设置图3.3采集声音3.3双音频信号的频谱分析在MATLAB中编辑m函数,使用wavread函数读取采集的声音文件(.wav)将它赋值给某一向量,再对其进行采样,然后使用plot语句画出相关的频谱图形。(1)Wavread函数调用格式:[y,fs,nbits]=wavread(file)功能说明:采样值放在向量x中,fs表示采样频率(Hz),bits表示采样位数。(2)快速傅里叶变换算法FFT计算DFT的函数fft,其调用格式如下:Xk=fft(x,n)参数x为被变换的时域序列向量,N是DFT变换区间长度,当n大于x的长度时,fft函数自动在x后面补零。,当n小于xn的长度时,fft函数计算x的前n个元素,忽略其后《双音频信号滤波去噪——使用双线性变换法设计的切比雪夫I型滤波器》第9页共19页面的元素。在本次课程设计中,我们利用fft函数对双音频信号进行快速傅里叶变换,就可以得到信号的频谱特性。(3)声音采样文件读取的程序(文件名:xiaomiao.wav);双音频信号的提取:[x,fs,bits]=wavread('D:\MATLAB7\work\xiaomiao.wav');%输入参数为文件的全路径和文件名,输出的第一个参数是每个样本的值,fs是生成该波形文件时的采样率,bits是波形文件每样本的编码位数。sound(x,fs,bits);%按指定的采样率和每样本编码位数回放N=length(x);%计算信号x的长度fn=2200;%单频噪声频率,此参数可改t=0:1/fs:(N-1)/fs;%计算时间范围,样本数除以采样频率x=x(:,1)';y=x+sin(fn*2*pi*t);sound(y,fs,bits);%应该可以明显听出有尖锐的单频啸叫声X=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;%计算频谱频率范围得到的原始双音频信号和加上噪音后的双音频信号的时域波形和频谱图如图3.4所示。《双音频信号滤波去噪——使用双线性变换法设计的切比雪夫I型滤波器》第10页共19页图3.4原始双音频信号和加噪后信号的时域和频谱图3.4滤波器设计设计指标:通带截止频率为2000Hz,阻带截止频率为2150Hz,通带波纹为1dB,阻带波纹为35dB,用双线性变换法设计的一个满足上述指标的切比雪夫I型滤波器双线性变换法设计切比雪夫I型滤波器fp=fn-200;fc=fn-50;%定义通带和阻带截止频率Rp=1;As=35;%定义通带波纹和阻带衰减wp=fp/fs*2*pi;ws=fc/fs*2*pi;%计算对应的数字频率T=1;Fs=1/T;%定义采样间隔Omegap=(2/T)*tan(wp/2);Omegas=(2/T)*tan(ws/2);%截止频率预畸变[cs,ds]=afd_c