拉弦音乐滤波去噪使用脉冲响应不变法设计的切比雪夫I型滤波器文档格式.docx
《拉弦音乐滤波去噪使用脉冲响应不变法设计的切比雪夫I型滤波器文档格式.docx》由会员分享，可在线阅读，更多相关《拉弦音乐滤波去噪使用脉冲响应不变法设计的切比雪夫I型滤波器文档格式.docx（18页珍藏版）》请在冰豆网上搜索。
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型滤波器是电子滤波器的一种。
它的特点是在通带具有等波纹响应。
其振幅平方函数如2-1式所示：
（2-1）
式中，N为整数，称为滤波器的阶，ε是通带波纹因子（它与Rp有关）。
其中N阶切比雪夫多项式如2-2式所示：
（2-2）
其中
切比雪夫滤波器存在两种可能的幅度平方响应的形状，当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）之间起伏。
对x>
1（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-6）
上式表明，采样脉冲响应不变法将模拟滤波器变换为数字滤波器时，它所完成的S平面到Z平面的变换，正是以前讨论的拉式变换到Z变换的标准变换，即首先对Ha（s）作周期严拓，然后再经过z=est的映射关系映射到Z平面上。
脉冲响应不变法变换前后系统函数的关系如式（2-7），（2-8）所示：
（2-7）
（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-9）
3.设计步骤
3.1设计流程图
拉弦音乐信号滤波去噪——使用脉冲不变响应法设计的切比雪夫I型滤波器的设计流程如图3-1所示：
图3-1脉冲响应不变法切比雪夫I对吹管信号去噪流程图
3.2拉弦音乐信号的采集
在网上下载一段吹管音乐信号，属性设置为8000Hz，8位单声道的音频格式。
如图3-2所示：
图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:
%截取前半部分
程序运行结果得到加噪前后信号的时域波形和频率波形，波形图如图3-4所示：
图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所示：
和
（3-1）
（2）利用设计参数
，
设计一个模拟滤波器
。
此课程设计需要用切比雪夫I型原型滤波器。
在MATLAB中用以下语句实现：
[cs,ds]=afd_chb1（OmegaP,OmegaS,Rp,As）;
%计算滤波器阶数和截止频率
[b,a]=imp_invr（cs,ds,T）;
%脉冲响应不变法得到数字滤波器系统函数分子分母系数
[db,mag,pha,grd,w]=freqz_m（b,a）;
%验证滤波器是否达到指定性能
（3）利用部分分式将
展开。
（4）现在将模拟极点{pk}变换为数字极点{epkT}，得到数字滤波器。
（5）化简得出作为z-1有理函数的H（z）。
根据上述步骤，按要求用脉冲响应不变法设计好了切比雪夫I滤波器。
其MATLAB实现源程序如图3-5所示：
图3-5得到满足要求的切比雪夫I的源程序
要验证滤波器是否满足要求就要知道wp/pi,ws/pi的值。
wp/pi
ans=
0.4500
ws/pi
0.4875
由上述程序运行的结果如图3-6所示：
图3-6设计得到的切比雪夫I滤波器的幅度、相位和脉冲响应
可以知道wp/pi=0.45,ws/pi=0.4875.由滤波器幅度响应图可以看出当w/pi大于0.4785时dB都大于16，所以设计的切比雪夫滤波器达到要求。
3.5信号滤波处理
前面已经用脉冲响应不变法设计好了我们要的切比雪夫滤波器，接着就对拉弦信号进行滤波处理，看自己设计的切比雪夫滤波器有没有对我们的拉弦信号进行了滤波处理。
所以就用filter函数进行滤波，源程序如下：
y_fil=filter（b*T,a*T,y）;
%得到滤波后的时间图
Y_fil=abs（fft（y_fil））;
%用设计好的滤波器对y进行滤波
Y_fil=Y_fil（1:
%计算频谱取前一半
我们将滤波前后的时域波形进行比较，并对其进行快速傅里叶变换，Y_fil=fft（y_fil）;
N/2），目的是对比前后的频域频谱，具体分析设计的滤波器是否达到设计要求。
滤波前后的拉弦音乐信号时域和频域对比如图3-7所示：
图3-7滤波前后时域和频域波形对比图
3.6结果分析
我们先采集拉弦信号，再按照步骤用脉冲响应不变法设计切比雪夫I型滤波器。
由图3-6可知，设计的阻带截止频率wc=（ws-wp）/pi，纵坐标差不多刚好在As处，所以设计的滤波器达到要求。
我们观察到图3-7滤波前后拉弦音乐信号的波形对比图，发现时域波形中加干扰噪声后有明显的变化，不过经过滤波后几乎没有变化，说明设计的滤波器达到要求。
但是可以看到在频域波形中，很明显地反应出设计的滤波器滤去了我们采集的拉弦信号中的噪声。
对比是一种很好的学习方法，不仅可以发现自己的问题，还可以通过对比不同的滤波方法更深入的理解所学的知识。
一个同学的课题是《弹拨音乐滤波去噪—使用脉冲响应不变法设计的切比雪夫I滤波器》，与我的对比就是信号源不一样，得到的图如图3-8所示：
图3-8弹拨音乐—脉冲响应不变法切比雪夫I型滤波器得到的对比图
与图3-7对比可知：
弹拨音乐比拉弦音乐幅度更大，加噪声后幅度变化幅度没我的变化明显。
最终的滤波效果差不多。
用sound函数听时，滤波后的效果和原始音乐信号相同。
4.出现的问题及解决方法
在这次的课程设计中我们设计的滤波器对拉弦信号的滤噪处理或多或少会因为我们的疏忽或者实践能力不够还是会遇到问题的，当遇到问题时我们不但运用自己学到的知识还得到同学和老师的帮助将其解决，不断完善和修改。
设计中出现的问题及解决方法如下：
（1）在看原始音乐信号的频谱时出现两种颜色的图形，原来是因为需要的拉弦信号，录音机的音频格式应该是PCM8000kHz，8位，单声道，而从网上下载下来的吹管音乐是立体音。
所以要改下属性就可以解决了。
（2）刚开始设计出的切比雪夫滤波器无法滤除加噪信号频谱中的噪声，而且滤除噪音后的拉弦音乐感觉噪声更大了。
向老师询问后才知道原来是因为设置T的时候，把T设置为1/fs所以声音扩大了，经更改后，设置T=1,就能滤波了，滤波后基本上与原来的拉弦音乐一样。
（3）在利用脉冲响应不变法设计切比雪夫滤波器时得不到理想的滤波器，把通带截止频率与阻带截止频率之间的差值设置的太小或太大，通过适当的调整参数，绘制出来的图像效果比较明显，基本符合设计指标。
通过不断调整As，最后取As=16dB,这样通过MATLAB运算出来的滤波器的阻带波纹达到要求。
（4）还要在滤波器性能曲线的wc处画一根竖线，这样更方便看出结果，其中wc处线的确定还需计算出wp/pi的值。
（5）所有的时间波形横坐标都要化为时间，滤波前后频谱的横坐标应是频率，这样在观察通带截止频率和阻带截止频率时更加精确，误差较小。
（6）开始设计滤波器时不知道怎样看出是否满足要求，询问老师后知道了当大与ws/pi时都大于As=16dB就说明满足要求了。
4.结束语
为期两周的DSP课程设计即将结束，这是第二次做有关专业方面的课程设计了。
因为有上一次的经验，所以这次再做课程设计就显得从容多了。
尽管如此，还是犯了不少错误，发现自己还不够认真不够细心，以前的知识学习的不是很牢固。
通过课程设计，感觉又复习了一遍书本的知识，遇到不懂得就翻书，或问老师。
不过还是要培养自己解决问题的良好习惯，所以有时候还得去图书馆借阅相关的书籍，来解决自己的疑问。
这次课程设计的平台还是MATLAB，由此可见MATLAB对本专业的重要性，所以通过又一次的学习和使用该软件，让我对MATLAB又更深的了解，同时也更加巩固了对MATLAB软件的使用。
经过这次的DSP课程设计，让我有机会将自己学到的理论知识运用到实际中，提高了自己的动手能力和思维能力。
在课程设计中发现自己的不足，所以在今后的学习和生活中我们要更加努力，学习好我们的专业知识并要能运用到实际。
这次的设计中老师指出的几个问题，在老师的指导和同学的帮助下解决了，当然也少不了自己的努力和思考。
能将自己平时学到的东西能运用到实际中，让理论和实际得以结合还是很不错的。
也让我在课程设计中找到了动手的乐趣和思考的快乐，很有成就感。
我从这次的课程设计中体会到，什么事都要我们真真正正用心去做了，才会使自己更加成长，没有学习就不可能有实践的机会，没有自己的实践就不会有所突破，希望这次的经历能让我们在以后的学习生活中不断成长。
最后，在此衷心地感谢老师和同学对我的帮助，也感谢学校给我们的机会，让我们能够将自己学到的知识运用到实际中！
参考文献
[1]张圣勤.MATLAB7.0实用教程[M].北京：
机械工程出版社，2006
[2]（美）维纳·
K·
英格尔,（美）约翰·
G·
普罗克斯（著）.刘树棠（译）.数字信号处理（MATLAB版）[M].第2版.西安：
西安交通大学出版社，2008
[3]程佩青．数字信号处理教程[M]．北京：
清华大学出版社，2002
[4]吴镇扬，数字信号处理[M].高等教育出版社，2004
[5]罗军辉，MATLAB7.0在数字信号处理中的应用[M].机械工业出版社，2005
[6]陈怀琛，数字信号处理教程：
MATLAB释义与实现[M].电子工业出版社，2004
附录1：
拉弦音乐信号滤波去噪设计源程序清单
%程序名称：
laxian.m
%程序功能：
采用脉冲响应不变法设计的切比雪夫I滤波器，对拉弦音乐滤波去噪。
%程序作者：
巢楚颉
%最后修改日期：
2012-3-2
\yinyue.wav'
）;
%输入参数为文件的全路径和文件名
%输出的第一个参数是每个样本的值
%fs是生成该波形文件时的采样率
%bits是波形文件每样本的编码位数
sound（x,fs,bits）;
%按指定的采样率和每样本编码位数回放
N=length（x）;
t=0:
deltaf=fs/N;
f=0:
X=abs（fft（x））;
X=X（1:
subplot（2,2,1）；
plot（t,x）；
xlabel（'
时间（单位：
s）'
ylabel（'
幅度'
title（'
原始拉弦信号'
）
subplot（2,2,2）；
plot（f,X）;
频率（单位:
HZ）'
幅度谱'
原始拉弦信号幅度谱图'
fn=2000%单频噪声
x=x（:
sound（y,fs,bits）;
Y=abs（fft（y））;
Y=Y（1:
subplot（2,2,3）；
plot（t,y）
加入单频干扰后的语音信号'
subplot（2,2,4）；
plot（f,Y）
频率（单位：
加入干扰后的语音信号幅度谱'
fs=8000;
fn=2000;
fp=fn-200;
fc=fn-50;
%定义通带和阻带截止频率
Rp=1;
As=16;
%定义通带波纹和阻带衰减
wp=fp/fs*2*pi;
ws=fc/fs*2*pi;
%计算对应的数字频率
T=1;
%定义采样间隔
OmegaP=wp/T;
OmegaS=ws/T;
%截止频率线性变换
Ripple=10^（-Rp/20）;
Attn=10^（-As/20）;
%转化为绝对形式
delta=[1,zeros（1,99）];
ha=filter（b,a,delta）;
%计算脉冲响应
Subplot（221）;
plot（w/pi,db）%画出滤波器的幅度响应图
w/pi'
dB'
滤波器幅度响应图'
Subplot（222）;
plot（w/pi,mag）%画出滤波器幅度响应图
幅度mag'
滤波器幅度响应'
Subplot（223）;
plot（w/pi,pha）%画出滤波器的相位响应图
相位pha'
滤波器相位响应'
Subplot（224）;
plot（ha）%画出滤波器的脉冲响应图
n'
h（n）'
滤波器脉冲响应'
）
subplot（3,2,1）;
plot（t,x）;
%画原始拉弦信号时间图
时间t'
原始拉弦信号x'
subplot（3,2,2）;
plot（f,X）;
%画原始拉弦信号幅度谱图
频率f'
原始拉弦信号幅度谱X'
subplot（3,2,3）;
plot（t,y）;
%画出加干扰拉弦信号时间
加干扰拉弦信号时间y'
subplot（3,2,4）;
plot（f,Y）;
%画出加干扰拉弦信号幅度谱图
加干扰拉弦信号幅度谱Y'
subplot（3,2,5）;
plot（t,y_fil）;
滤波后拉弦信号时间y_fil'
subplot
升级会员 