弹拨音乐滤波去噪—用凯塞窗设计的FIR滤波器教材

1、弹拨音乐滤波去噪一一使用凯塞窗设计的滤波器第14页 共15页弹拨音乐滤波去噪使用凯塞窗设计的FIR滤波器学生姓名：唐柯指导老师：胡双红摘 要 本课程设计主要是通过使用凯塞窗设计一个FIR滤波器以对弹拨进行滤波去噪处理。本设计首先从网下下载一段弹拨音乐，依据对该信号的频谱分析， 给定相关指标。以MATLAB件为平台，采用凯塞窗设计满足指标的FIR滤波器, 以该弹拨音乐进行滤波去噪处理。 通过对比滤波前后的波形图，深入了解滤波器 的基本方法。通过程序调试及完善，该设计基本满足设计要求。关键词 滤波去噪；FIR滤波器；凯塞窗函数；MATLAB1引言数字滤波器是一种用来过滤时间离散信号的数字系统，它是

2、通过对抽样数据 进行数学处理来达到频域滤波的目的。随着现代通信的数字化，数字滤波器变得更加重要。数字滤波器的种类很多，但总的来说可以分成两大类，一类是经典滤 波器，另一类可称为现代滤波器。从滤波特性方面考虑，数字滤波器可分成数字 高通、数字低通、数字带通和数字带阻等滤波器。从实现方法上考虑，将滤波器 分成两种，一种称为无限脉冲响应滤波器，简称IIR ( Infinite Impulse Response滤波器，另一种称为FIR( Finite Impulse Response滤波器。设计FIR数字滤 波器的方法有窗函数法、频率采样法和等波纹最佳逼近法等。1.1课程设计目的数字信号处理(Digi

3、tal Signal Processing,简称DSP)是一门涉及许多学科 而又广泛应用于许多领域的新兴学科。20世纪60年代以来，随着计算机和信息 技术的飞速发展，数字信号处理技术应运而生并得到迅速的发展。在过去的二十多年时间里，数字信号处理已经在通信等领域得到极为广泛的应用。数字信号处理是利用计算机或专用处理设备，以数字形式对信号进行采集、变换、滤波、估 值、增强、压缩、识别等处理，以得到符合人们需要的信号形式。在本次课程设计中，最主要的设计是设计FIR滤波器，FIR滤波器的设计方法主 要分为两类：第一类是基于逼近理想滤波器器特性的方法包括窗函数法、 频率采 样法、和等波纹最佳逼近法；第二

4、类是最优设计法。本次的课程设计主要采用的是第一类设计方法，是利用凯塞窗函数法设计 FIR滤波器对一段弹拨音乐进行滤波去噪，通过这一过程，对滤波前后波形进行 对比分析得到结论。此课程设计比较简单，主要是将书本中的知识运用到现实中， 并且根据自己对设计题目的理解， 运用软件编写出程序实现这一设计，也是我们 对数字信号处理的原理进行验证的一个过程。对此，也可以加深我们对所学知识 的理解，培养我们的动手能力。1.2课程设计的要求（1）通过利用各种不同的开发工具实现模拟信号数字化、信道编解码、基带数字信号编解码、数字信号的调制解调和弹拨音乐的滤波去噪等课题，掌握数字信 号的分析方法和处理方法。（2）按要

5、求编写课程设计报告书，能正确阐述设计和实验结果等等。（3） 通过课程设计培养学生严谨的科学态度，认真的工作作风和团队协作精神。（4）在老师的指导下，要求每个学生独立完成课程设计的全部内容。1.3设计平台课程设计的主要设计平台式 MATLAB 7.0。如下图1.1所示：MATLAB是矩 阵实验室（Matrix Laboratory）的简称，是美国MathWorks公司出品的商业数学 软件，用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言 和交互式环境，主要包括 MATLAB和Simulink两大部分。MATLAB和 Mathematics Mape并称为三大数学软件。它在数学类科

6、技应用软件中在数值计 算方面首屈一指。MATLAB可以进行矩阵运算、绘制函数和数据、实现算法、 创建用户界面、连接其他编程语言的程序等，主要应用于工程计算、控制设计、 信号处理与通讯、图像处理、信号检测、金融建模设计与分析等领域。Fie Edit View D电ug Desktop Vndow HelpK 冒裂 Cdrreni DlrectOty： tWATLAB7V0rk 直ISfiortcuts 3 Hem to Add 回 Whsfs blewH, v =fzeqz (b, a, Wi * rB iiLagH=abs(H), angH=argle (H) ； subp plot (wft

7、agH):.To gmt started select NALL AB Help or DgmogThe eleiuent typemust be termnat&d byCould not parse the file: f:：matlsb7toolbox0处收敛，极点全部在z=0处（因果系统）；（3）结构上主要是非递归结构，没有输出到输入的反馈，但有些结构中（例如频 率抽样结构）也包含有反馈的递归部分。FIR滤波器的系统函数用下式表示：。 H（n）就是FIR滤波器的单位脉冲响应。FIR滤波器最重要的优点就是由于不存在系统极点，FIR滤波器是绝对稳定的系统。相较于IIR滤波器,FIR滤波器有

8、以下的优点：（1）可以很容易地设计线性相位的滤波器。线性相位滤波器延时输入信号，却 并不扭曲其相位。（2）实现简单。在大多数DSP处理器，只需要对一个指令积习循环就可以完成 FIR计算。（3）适合于多采样率转换，它包括抽取（降低采样率），插值（增加采样率）操作。无 论是抽取或者插值，运用FIR滤波器可以省去一些计算，提高计算效率。相反，如果使用IIR滤波器，每个输出都要逐一计算，不能省略，即使输出要丢弃。（4） 具有理想的数字特性。在实际中，所有的 DSP滤波器必须用有限精度（有 限bit数目）实现，而在IIR滤波器中使用有限精度会产生很大的问题，由于采 用的是反馈电路，因此IIR通常用非常少

9、的bit实现，设计者就能解决更少的与 非理想算术有关的问题。（5） 可以用小数实现。不像IIR滤波器，FIR滤波器通常可能用小于1的系数 来实现。（如果需要，FIR滤波器的总的增益可以在输出调整）。当使用定点DSP 的时候，这也是一个考虑因素，它能使得实现更加地简单。FIR滤波器的缺点在于它的性能不如同样阶数的IIR滤波器，不过由于数字计算硬件的飞速发展，这一点已经不成为问题。再加上引入计算机辅助设计，FIR滤波器的设计也得到极大的简化。基于上述原因，FIR滤波器比IIR滤波器的应用更广。2.2窗口设计法窗函数设计法的基本思想是用 FIRDF逼近洗完的滤波特性。设希望逼近的 滤波器的频率响应函

10、数为m鈔 其单位脉冲响应为h/n)表示。为了设计简单 方便，通常选择叫辭巧为具有片段常数特性的理想滤波器。因此牟回是无限长 非因果序列，不能直接作为FIRDF的单位脉冲响应。窗函数设计法就是截取呃 为有限长的一段因果序列，并用合适的窗函数进行加权作为FIRDF的单位脉冲响应h(n)。用窗函数法设计FIRDF的具体设计步骤如下：(1) 构造希望逼近的频率响应函数 气妒个(2) 求出吨。(3) 加窗得到FIRDF的单位脉冲响应h(n)。h(n)二i:w(n)。式中，w(n)称为 窗函数，其长度为N。如果要求设计第一类线性相位 FIRDF，则要求h(n)关于 (N-1)/2点偶对称。而加吧关于n=t

11、点偶对称，所以要求t=(N-1)/2。同时要求w(n) 关于(N-1)/2点偶对称。用窗函数法设计第一类线性相位 FIRDF的步骤：(1) 选择窗函数类型和长度，写出才窗函数褪M；的表达式。根据阻带最小衰减选择窗函数喩时的的类型，再根据过渡带宽度确定所选窗函数 的长度N。用窗函数法设计的FIRDF通带波纹幅度近似等于阻带波纹幅度。一 般阻带最小衰减达到40dB以上，则通带最大衰减就小于 O.ldB。所以用窗函数 法设计FIRDF时，通常只考虑阻带最小衰减就可以了。(2) 构造希望逼近的频率响应函数比根据设计需要，一般选择线性相位理想滤波器(理想低通、理想高通、理想带通、 理想带阻)。理想滤波器

12、的截止频率 c近似为最终设计的FIRDF的过渡带中心频率，幅度函数衰减一半(约-6dB)。所以一般取和分别为通带边界频率和阻带边界频率(3) 计算 】期抑如的辭呱(4) 加窗得到设计结果：h( n) =：1一 w (n)。3设计步骤3.1设计流程图3.1设计流程3.2采集弹拨音乐在网络中找到自己想要的歌曲文件， 并下载下来，即可对音乐文件进行文件 格式转换及声音属性的变化。在转换过程中，将音乐格式设置为wav格式。然后再点击文件/另存为，即可将该音乐文件存在 MATLAB文件下的work文件夹 里。然后在Matlab软件平台下，利用函数 wavread对音乐信号进行采样，记住 采样频率和采样点

13、数。通过 wavread函数的使用，让我们很快理解了采样频率、采样位数等概念。采集完成后在信号中加入一个单频噪声。对应程序如下：x,fs,bits=wavread(Music.wav);soun d(x,fs,bits);N=length(x);% 计算信号x的长度fn=1900;% 单频噪声频率t=O:1/fs:(N-1)/fs;%计算时间范围，样本数除以采样频率x=x;y=x+0.1*si n(fn*2*pi*t);%加单频噪声sound(y,fs,bits);%应该可以明显听出有尖锐的单频啸叫声X=abs(fft(x);Y=abs(fft(y);%对原始信号和加噪信号进行fft变换X=X

14、(1:le ngth(X)/2);Y=Y(1:le ngth(Y)/2);%截取前半部分deltaf=fs/2/le ngth(X);灿算频谱的谱线间隔f=0:deltaf:fs/2-deltaf;%计算频谱频率范围用plot函数画出原始语音信号及加噪语音信号的时域图和频谱图，如图3.2所示:图3.2频谱分析图3.3滤波器设计在该滤波器的设计中，我们给出该滤波器的性能指标如下：fpd=1800;fsd=1850;fsu=1950;fpu=2000;Rp=1;As=100;截止频率也可以任意自选，在单频噪声干扰附近即可。我设计滤波器过程对应的程序如下：fpd=1800;fsd=1850;fsu=

15、1950;fpu=2000;Rp=1;As=100;% 带阻滤波器设计指标fcd=(fpd+fsd)/2;fcu=(fpu+fsu)/2;df=mi n( (fsd-fpd),(fpu-fsu);%计算上下边带中心频率，和频率间隔wcd=fcd/fs*2*pi;wcu=fcu/fs*2*pi;dw=df/fs*2*pi;%将Hz为单位的模拟频率换算为rad为单位的数字频率wsd=fsd/fs*2*pi;wsu=fsu/fs*2*pi;M=ceil(As-7.95)/(2.285* dw)+1)+1;%计算凯塞窗设计该滤波器时需要的阶数n=0:M-1; beta=0.1102*(As-8.7);

16、% 定义时间范围w_kaiser= kaiser(M,beta); % 产生 M 阶的凯塞窗 hd_bs=ideal_lp(wcd,M)+ideal_lp(pi,M)-ideal_lp(wcu,M);%调用自编函数计算理想带阻滤波器的脉冲响应h_bs=w_kaiser*hd_bs;%用窗口法计算实际滤波器脉冲响应db,mag,pha,grd,w=freqz_m(h_bs,1); %调用自编函数计算滤波器的频率特性 在将滤波器设计好，频率特性求出来之后，我编写如下程序，来画出该带阻滤波 器的特性图。通过该程序得到的图如下图3.3.1所示。Zile Edit jLi ew Lnsert Zools

17、 desktop Window Helpd&q ei|V|st 宴丨囤| 也qp百率O频JOOO404率频0020000O6图3.3滤波器的特性图3.4信号滤波处理在将滤波器设计好后，我们用自己设计的带阻滤波器对采集的弹拨音乐进行 滤波。在Matlab中，FIR滤波器利用函数fftfilt对信号进行滤波。我对信号进行滤波处理对应的程序如下：y_fil=fftfilt(h_bs,1,y); %用设计好的滤波器对加噪信号y进行滤波 Y_fil=fft(y_fil);Y_fil=Y_fil(1:length(Y_fil)/2); % 计算频谱取前一半3.5结果分析在将加噪信号滤波之后，我们将滤波前后

18、弹拨音乐的波形及频谱图相互比较。在同一张大图里分别绘制原始信号x，加噪信号y,滤波去噪信号y_fil的时域波形和频谱，以便比较和分析。经过这段程序画出来的三个信号的时域波形和频谱图如下图3.4所示。FigurE 3匚|石|叹Zile Edit Yi ew Insert lools Uesktop iLindow Kelp Q S 別釵霁|硬| 口目| 叵o5 0 5OL0o o o3一日fB,2时o o o3o Hu o o oO O O3OOO2O-OOO o- OOO图3.4三个信号的时域波形和频谱图比较从图3.5.1中我们可以看出，原信号与滤波去噪信号的时域图基本相似，只 有边缘部分有点

19、差异；原信号与滤波去噪信号的频谱图波形也大致相似。通过观 察可以看到，加噪信号的时域图中大部分都被加入的噪声给遮盖了，加噪信号的频谱图中，我们可以很明显地看到与原信号频谱图相比，它在1900Hz左右处有一个尖脉冲，而滤波去噪信号的频谱图中该尖脉冲已经消失，波形大致与原图相似，可见滤波去噪效果基本不错。在将三个信号的时域波形和频谱图比较之后，我们还要通过回放去滤波去噪弹拨音乐，来跟原信号相比，以检验滤波器的效果。在Matlab中，函数sou nd可以对声音进行回放。其调用格式为：sou nd (x, fs, bits)。我用sou nd(y_fil,fs,bits) 语句回放该滤波去噪信号，便可

20、以感觉到滤波后的弹拨音乐与原信号差不多，但 仍有一点点变化。4出现的问题和解决方法在这前面一周的时间里，我们每天的上午都是待在学校的机房上机，同时，老师跟我们进行讲解，然后，当我们有不懂的问题时可以随时向老师请教。 在老 师没有进行讲解时，我总以为这次的课程设计会有不小的难度， 可是当老师跟我 们讲解了之后，我却发现，我们的课程设计其实也不是很难，只要我们细心一点， 并有耐心，就能解决我们所遇到的问题。首先，我按照课程设计的要求在网络上找到一段音乐， 并对其进行格式转换 及属性变换，但是当我真正制作完成之后并存盘后，运行 MATLAB后，才感觉 音乐文件过于长，造成出现的频谱图很模糊，根本看不

21、清楚，后面看了一下任务 指导书及问了同学之后，截取其中的一小段，效果就好了许多。在此次课程设计中，虽然没有预计的那么难，可还是遇到了一些小小的问题。 庆幸的是，能够在老师和同学的帮助下解决了。我觉得最重要的就是要有耐心， 要细致，要不懂就问，这样，只要自己坚持，一切问题都会解决的。5结束语在本次课程设计中，让我知道了以前所学信号与系统与数字信号处理 的结合，并在实际运用中设计滤波器的过程。 课程设计不仅要求对滤波器理论的 研究，更重要的是实际设计中遇到的问题。因为有了这次课程设计，我不得不对其设计原理进行更深一层次的理解， 对 书中原来学到的理论，仅知道了其表面，而不知其原因。在设计中也使我对

22、一些 概念有了更深刻的认识。例如：在指标方面，我混淆了模拟指标和数字指标的概 念，经过老师的点拨，自己更加明确，而且记忆深刻。还有在课程设计中每一次 的数据输入都有其重要意义，用 MATLAB编译程序时，可以根据滤波器指标的 要求实时知道对滤波器的影响。例如，编好程序后，调试成功，任意改变输入阻 带或者通带衰减，可以看到输出波形的变化，改变截止频率wc，同样可以看到输出波形的变化。由此，对理论的理解就更加简单方便，而且记忆力深刻。除此之外，对程序的编译不是一蹴而就的，而是经过多次的编译与调试。以 前用MATLAB就是简单的输入，可以说都不是自己的劳动成果，但这次不一样， 课程设计没有别人给你编

23、好，而是自己写，出错率就大大提高了，但这是过程， 学习就是在过程中进行的，经过自己几天的脑力劳动，再加上同学们的帮助，不 仅对读程序有了很大提高，而且自己的编译水平也上了一个新台阶， 更加熟系了 MATLAB的应用，也对其中的函数有了大概的了解，对其中一些函数也相当熟 练，滤波器设计中用到了 ceil ()、freqz ()等很重要的函数，虽然可以用 firl ()、fir2 ()可以直接调用，但那样就达不到对真正理论设计过程的理解和运用。完成整个设计过程后，学到的东西已经不仅仅上面的那些东西， 还有就是同 学之间的共同努力和探讨和设计过程中的每一个细节， 也许每一个细节的错误就 可能导致结局

24、的失败，所以我认为这次收获最大的莫过于静心， 学习不能急，一 定要冷静，心无旁骛，不放过任何一个细节，就能带来凯旋的消息。在此，感谢老师一学年以来对我们不仅学习上的关心，还有生活中的照顾，我也不会辜负老师的期望，继续努力深造。3 弹拨音乐滤波去噪一一使用凯塞窗设计的滤波器第13页 共15页参考文献1 程佩青数字信号处理教程M北京：清华大学出版社，2002.2 刘敏，魏玲.Matlab通信仿真与应用M北京：国防工业出版社，2001.3 Sanjit K.Mitra著，孙洪，余翔宇 译，数字信号处理实验指导书，电子 工业出版社，2005年。4 Edward W. Kame n, Bonnie S.

25、Heck编，信号与系统基础应用 Web和MATLAB （第二版），科学出版社，2002年。 张贤达 编，现代信号处理（第二版），清华大学出版社， 2002年。弹拨音乐滤波去噪一一使用凯塞窗设计的滤波器第19页 共15页附录1 ：弹拨音乐滤波去噪一一用凯塞窗设计的滤波器源程序清单x,fs,bits=wavread(Music.wav);soun d(x,fs,bits);N=le ngth(x);fn=1900;t=O:1/fs:(N-1)/fs;x=x;y=x+ 0.1*si n(fn*2*pi*t);soun d(y,fs,bits);X=abs(fft(x); Y=abs(fft(y);X

26、=X(1:length(X)/2); Y=Y(1:length(Y)/2);deltaf=fs/N;f=0:deltaf:fs/2-deltaf;subplot(2,2,1);plot(t,x);grid on;xlabel(时间(单位：s);ylabel(幅值);title(原始音乐信号);axis(0,4,-1,1);subplot(2,2,2);plot(f,X);grid on;xlabel(频率(单位：Hz);ylabel(幅度谱);title(音乐信号幅度谱图);axis(0,4000,0,5000);subplot(2,2,3);plot(t,y);grid on;xlabel(

27、时间(单位：s);ylabel(幅值);title(加入单频干扰后音乐信号);axis(0,4,-1,1);subplot(2,2,4);plot(f,Y);grid on;xlabel(频率(单位：Hz);ylabel(幅度谱);title(加入干扰后的音乐信号幅度谱图);axis(0,4000,0,5000);fpd=1800;fsd=1850;fsu=1950;fpu=2000;Rp=1;As=100;fcd=(fpd+fsd)/2;fcu=(fpu+fsu)/2;df=mi n( fsd-fpd),(fpu-fsu);wcd=fcd/fs*2*pi;wcu=fcu/fs*2*pi;dw=df/fs*2*pi;wsd=fsd/fs*2*pi;wsu=fsu/