《基于mltlab的各种数字模拟滤波器课程设计.doc》由会员分享,可在线阅读,更多相关《基于mltlab的各种数字模拟滤波器课程设计.doc(30页珍藏版)》请在三一办公上搜索。
1、目录1 引言1 1.1 数字滤波器的背景及意义1 1.2 设计要求1 1.3 设计目的12 数字滤波器的设计22.1 数字滤波器的基本概念22.1.1 数字滤波器的分类22.1.2 数字滤波器的技术要求22.2 利用双线性变换法设计IIR数字滤波器22.2.1 变换原理22.2.2 变换常数C的选择32.2.3 模拟滤波器的数字化42.3 利用窗函数法设计FIR数字滤波器43 基于MATLAB的数字滤波器的设计63.1 matlab软件介绍63.2 MATLAB常用函数简介73.3 IIR滤波器设计83.3.1 IIR低通滤波器设计83.3.2 IIR高通滤波器设计93.3.3 IIR带通滤波
2、器设计113.3.4 IIR滤波器设计整体流程133.4 FIR滤波器设计143.4.1 FIR低通滤波器设计143.4.2 FIR高通滤波器设计153.4.3 FIR高通滤波器设计163.4.4 FIR滤波器设计整体流程184 总结19参考文献20附录211 引言1.1 数字滤波器的背景及意义数字滤波器是指完成信号滤波处理功能的,用有限精度算法实现的离散时间线性非时变系统,其输入是一组数字量,其输出是经过变换的另一组数字量。因此,它本身即可以是用数字硬件装配成的一台完成给定运算的专用数字计算机,也可以是将所需运算编成程序,让通用计算机来执行。数字滤波器具有稳定性高、精度高、灵活性大等优点。随
3、着数字技术的发展,用数字技术实现滤波器的功能越来越受到人们的注意和广泛的应用。1.2 设计要求利用MATLAB仿真软件系统结合窗函数法和双线性变换法设计以下三种数字滤波器:(1)低通滤波器:,最大衰减 ,最小衰减;(2)高通滤波器:,最大衰减 ,最小衰减; (3)带通滤波器:, ,最大衰减 ,最小衰减;1.3 设计目的掌握双线性变换法及窗函数设计IIR FIR数字滤波器的具体设计方法及其原理,观察双线性变换及脉冲响应不变法设计的滤波器的频域特性,了解双线性变换法及窗函数法的特点。2 数字滤波器的设计2.1 数字滤波器的基本概念2.1.1 数字滤波器的分类数字滤波器从实现的网络结构或者从单位脉冲
4、响应分类,可以分成无限脉冲响应(IIR)滤波器和有限脉冲响应(FIR)滤波器。它们的系统函数分别为: 2.1.2 数字滤波器的技术要求我们通常用的数字滤波器一般属于选频滤波器。假设数字滤波器的传输函数用下式表示: 通带内和阻带内允许的衰减一般用dB数表示,通带内允许的最大衰减用Ap表示,阻带内允许的最小衰减用As表示,Ap和As分别定义为: 如将归一化为1后得:2.2 利用双线性变换法设计IIR数字滤波器2.2.1 变换原理双线性变换法是使数字滤波器的频率响应与模拟滤波器的频率响应相似的一种变换方法。为了克服脉冲响应不变法的多值映射这一缺点,我们首先把整个s平面压缩变换到某一中介的平面的一横带
5、里(宽度为,即从到),然后再通过上面讨论过的标准变换关系将此横带变换到整个z平面上去,这样就使s平面与z平面是一一对应的关系,消除了多值变换性,也就消除了频谱混叠现象,基本原理如图所示。将s平面整个平面压缩到平面的到,可采用以下的变换关系 其中C为常数;这样变为,变为,可将上式写成图1 双线性变换法的映射关系 令,则可得 再将平面通过以下标准变化关系映射到z平面后可得 2.2.2 变换常数C的选择为了使模拟滤波器与数字滤波器在低频处有较确切的对应关系,即在低频处有,当较小时有由上式可知因而得 则Z和S的表达式可重新写成: 即 2.2.3 模拟滤波器的数字化由于双线性变换法中,s与z之间有简单的
6、代数关系,故可由模拟系统函数通过代数置换直接得到数字滤波器的系统函数。即 可见数字滤波器的极点数等于模拟滤波器的极点数。频率响应也可用直接置换得到 这一公式可用于将滤波器的数字域指标,转换为模拟域指标。2.3 利用窗函数法设计FIR数字滤波器窗函数设计法是一种通过截短和计权的方法使无限长非因果序列成为有限长脉冲响应序列的设计方法。通常在设计滤波器之前,应该先根据具体的工程应用确定滤波器的技术指标。在大多数实际应用中,数字滤波器常常被用来实现选频操作,所以指标的形式一般为在频域中以分贝值给出的相对幅度响应和相位响应。用窗函数法设计FIR滤波器的步骤如下:(1)根据过渡带宽及阻带衰减要求,选择窗函
7、数的类型并估计窗口长度N(或阶数M=N-1)。窗函数类型可根据最小阻带衰减AS独立选择,因为窗口长度N对最小阻带衰减AS没有影响。在确定窗函数类型以后,可根据过渡带宽小于给定指标确定所拟用的窗函数的窗口长度N。设待求滤波器的过渡带宽为,它与窗口长度N近似成反比。窗函数类型确定后,其计算公式也确定了,不过这些公式是近似的,得出的窗口长度还要在计算中逐步修正。原则是在保证阻带衰减满足要求的情况下,尽量选择较小的N。在N和窗函数类型确定后,即可调用MATLAB中的窗函数求出窗函数。(2)根据待求滤波器的理想频率响应求出理想单位脉冲响应。如果给出待求滤波器的频率响应为,则理想的单位脉冲响应可以用傅里叶
8、反变换式得到,在一般情况下,是不能用封闭公式表示的,需要采用数值方法表示。从=0到=2采样N点,采用离散傅里叶反变换(IDFT)即可求出。(3)计算滤波器的单位脉冲响应h(n)。它是理想单位脉冲响应和窗函数的乘积,即,在MATLAB中用点乘命令表示为。(4)验算技术指标是否满足要求。为了计算数字滤波器在频域中的特性,可调用freqz子程序,如果不满足要求,可根据具体情况,调整窗函数类型或长度,直到满足要求为止。使用窗函数法设计时要满足以下两个条件:(1)窗谱主瓣尽可能地窄,以获得较陡的过渡带;(2)尽量减少窗谱的最大旁瓣的相对幅度,也就是使能量尽量集中于主瓣,减小峰肩和纹波,进而增加阻带的衰减
9、。在实际工程中常用的窗函数有五种,即矩形窗、三角窗、汉宁窗、汉明窗和凯泽窗。3 基于MATLAB的数字滤波器的设计3.1 MATLAB软件介绍Matlab是由matrix和laboratory两个英文单词的前三个字母组合而成,是“Matrix Laboratory”的缩写,意为“矩阵实验室”,它是MathWorks公司于1984年推出的一套高性能的数值计算可视化软件,集数值分析、矩阵运算、信号处理和图形显示于一体,被誉为“巨人肩上的工具”,是当今美国很流行的科学计算软件。信息技术、计算机技术发展到今天,科学计算在各个领域得到了广泛的应用在许多诸如控制论、时间序列分析、系统仿真、图像信号处理等方
10、面产生了大量的矩阵及其相应的计算问题。自己去编写大量的繁复的计算程序,不仅会消耗大量的时间和精力,减缓工作进程,而且往往质量不高美国Mathwork软件公司推出的Matlab软件就是为了给人们提供一个方便的数值计算平台而设计的。Matlab是一个交互式的系统,它的基本运算单元是不需指定维数的矩阵,按照IEEE的数值计算标准(能正确处理无穷数Inf(Infinity)、无定义数NaN(not-a-number)及其运算)进行计算。系统提供了大量的矩阵及其它运算函数,可以方便地进行一些很复杂的计算,而且运算效率极高。Matlab命令和数学中的符号、公式非常接近,可读性强,容易掌握,还可利用它所提供
11、的编程语言进行编程完成特定的工作。除基本部分外,Matlab还根据各专门领域中的特殊需要提供了许多可选的工具箱,如应用于自动控制领域的Control System工具箱和神经网络中Neural Network工具箱等。MATLAB主要特点:(1)语言简洁紧凑,使用方便灵活,库函数极其丰富。MATLAB程序书写形式自由,利用起丰富的库函数避开繁杂的子程序编程任务,压缩了一切不必要的编程工作。由于库函数都由本领域的专家编写,用户不必担心函数的可靠性。可以说,用MATLAB进行科技开发是站在专家的肩膀上。(2)运算符丰富。由于MATLAB是用C语言编写的,MATLAB提供了和C语言几乎一样多的运算符
12、,灵活使用MATLAB的运算符将使程序变得极为简短。(3)MATLAB既具有结构化的控制语句(如for循环,while循环,break语句和if语句),又有面向对象编程的特性。(4)程序限制不严格,程序设计自由度大。例如,在MATLAB里,用户无需对矩阵预定义就可使用。(5)程序的可移植性很好,基本上不做修改就可以在各种型号的计算机和操作系统上运行。(6)MATLAB的图形功能强大。在FORTRAN和C语言里,绘图都很不容易,但在MATLAB里,数据的可视化非常简单。MATLAB还具有较强的编辑图形界面的能力。(7)MATLAB的缺点是,它和其他高级程序相比,程序的执行速度较慢。由于MATLA
13、B的程序不用编译等预处理,也不生成可执行文件,程序为解释执行,所以速度较慢。(8)功能强大的工具箱是MATLAB的另一特色。MATLAB包 含两个部分:核心部分和各种可选的工具箱。核心部分中有数百个核心内部函数。其工具箱又分为两类:功能性工具箱和学科性工具箱。功能性工具箱主要用来扩充 其符号计算功能,图示建模仿真功能,文字处理功能以及与硬件实时交互功能。功能性工具箱用于多种学科。而学科性工具箱是专业性比较强的,control,toolbox,signl proceessing toolbox,commumnication toolbox等。这些工具箱都是由该领域内学术水平很高的专家编写的,所以
14、用户无需编写自己学科范围内的基础程序,而直接进行高、精、尖的研究。(9)源程序的开放性。开放性也许是MATLAB最受人们欢迎的特点。除内部函数以外,所有MATLAB的核心文件和工具箱文件都是可读可改的源文件,用户可通过对源文件的修改以及加入自己的文件构成新的工具箱。3.2 MATLAB常用函数简介(1)wavread函数:y,fs,bits=wavread(Blip)用于读取语音,采样值放在向量y中,fs表示采样频率(Hz),bits表示采样位数。sound(x,fs,bits); 用于对声音的回放。向量y则就代表了一个信号(也即一个复杂的“函数表达式”)也就是说可以像处理一个信号表达式一样处
15、理这个声音信号。(2)FFT函数:FFT用于序列快速傅立叶变换。其调用格式为:y=fft(x)。其中,x是序列,y是序列的FFT,x可以为一向量或矩阵:若x为一向量,y是x的FFT,且和x相同长度。若x为一矩阵,则y是对矩阵的每一列向量进行FFT。(3)Fir1函数:fir1函数是基于窗函数的FIR滤波器设计标准频率响应形状。其调用格式为b=fir1(N,we,ftype,window)。ftype和window可以默认。b=fir1(N,wc)可得到截止频率为wc且满足线性相位条件的N阶FIR低通滤波器,window默认选用hamming窗。其单位脉冲响应h(n)为:h(n)=b(n+1),
16、n=0,1,2,N。当wc=wc1,wc2时,得到的是通带为wc1wwc2的带通滤波器。b=fir1(N,wc,ftype),可设计高通和带阻滤波器。当ftype=high时,设计高通FIR滤波器;当ftype=stop时,设计带阻FIR滤波器。 (4)sound(x,fs,bits):该函数用于播放生硬,向量y就代表了一个信号(也即一个复杂的“函数表达式”),也就是说可以像处理一个信号表达式一样处理这个声音信号。3.3 IIR滤波器设计3.3.1 IIR低通滤波器设计根据第二节数字滤波器的设计原理设计滤波器,首先将数字域的指标转化为模拟域的指标设计模拟低通滤波器,然后应用双性变换发将模拟滤波
17、器转化为数字滤波器。将指标转换成归一化模拟低通滤波器的指标,通过归一化的模拟低通滤波器阶数N和3dB截止频率的计算,将模拟域频率变换成模拟低滤波器H(s),并用双线性变换法将H(s)转换成数字低通滤波器H(z),由此得到低通数字滤波器。其相关程序为:fb=1000;fc=1200;Ap=1;As=15;Fs=22050;ws=2*pi*fc/Fs;wp=2*pi*fb/Fs;N=ceil(12*pi/(ws-wp);Wn=(fb+fc)/Fs;b=fir1(N,Wn,blackman(N+1);h1,w1=freqz(b,1);分析结果如图2所示:图2 IIR低通滤波器的性能分析结果如图3、4
18、所示 图3滤波前后信号的波形分析图4滤波前后信号的频谱分析分析滤波前后信号的变化,波形图变窄,频谱变化表明:信号的低频率段被保留,高频率段被滤除。分析滤波前后的声音变化,滤波后声音明显变低而沉闷,这说明滤波器设计基本符合指标要求。3.3.2 IIR高通滤波器设计与低通滤波器的设计原理相似,首先将数字域的指标转化为模拟域的指标设计模拟低通滤波器,然后应用双性变换发将模拟滤波器转化为数字滤波器。将指标转换成归一化模拟低通滤波器的指标,通过归一化的模拟低通滤波器阶数N和3dB截止频率的计算,将模拟域频率变换成模拟低滤波器H(s),并用双线性变换法将H(s)转换成数字低通滤波器H(z),由此得到低通数
19、字滤波器。其相关程序为:fb=5000;fs=4800;Ap=1;As=15;Fs=22050;ws=2*pi*fb/Fs;wp=2*pi*fs/Fs;N=2*ceil(12*pi/(ws-wp);Wn=(fb+fs)/Fs;b=fir1(N,Wn,high,blackman(N+1);subplot(111);h1,w1=freqz(b,1);plot(w1/pi,20*log10(abs(h1);grid;分析滤波器显示如图5所示:。图5IIR高通滤波器的性能分析读入信号对信号傅里叶变换进行频谱分析,利用filter函数滤波,对得到的信号傅里叶变换进行频谱分析,与滤波前进行比较。并在滤波前
20、后用sound函数播放语音信号,观察滤波前后的变化。结果如图6、7所示图6滤波前后信号的波形 图7滤波前后信号的频谱分析分析滤波前后信号的变化,波形图变窄,频谱变化表明:信号的高频率段被保留,低频率段被滤除。分析滤波前后的声音变化,滤波后声音明显变高而尖锐,这说明滤波器设计基本符合指标要求。3.3.3 IIR带通滤波器设计设计程序为:fb1=1200;fb2=3000;fc1=1000;fc2=3200;fs =8000;as =15;ap = 1;wb =fb1 fb2*2*pi/fs;wc =fc1 fc2*2*pi/fs;Wb =2*fs*tan(wb/2);Wc =2*fs*tan(w
21、c/2);n,wn=buttord(Wb,Wc,ap,as,s);b,a=butter(n,wn,s);bd,ad=bilinear(b,a,fs);h,w=freqz(bd,ad);subplot(111);plot(w*fs/(2*pi),20*log10(abs(h);grid; 显示结果如图8所示:图8IIR带通滤波器的性能分析对比滤波前后的波形与频谱图,如图9、10所示:图9滤波前后信号的波形分析图10滤波前后信号的频谱分析分析滤波前后信号的变化,波形图变窄。频谱变化表明:信号的低频率段和高频率段被滤除,中间频率段被保留;分析滤波前后的声音变化,滤波后声音一定程度上变的尖锐,与高通滤
22、波器滤波后的声音相比较低。这说明滤波器设计基本符合指标。3.3.4 IIR滤波器设计整体流程 开始 读入数字滤波器技术指标将指标转换成归一化模拟滤波器的指标计算归一化的模拟滤波器阶数N和3db截止频率模拟域频率变换,将信号变换成模拟滤波器H(s)用双线性变换法将H(s)转换成数字滤波器H(z)输入信号后显示相应结果 结束图11 程序序流程图3.4 FIR滤波器设计3.4.1 FIR低通滤波器设计分析性能指标,确定所使用的窗函数。经调试,应用blackman窗比较接近所要求的滤波器指标。确定窗函数的系数b和a(FIR数字滤波器a=1),并显示滤波器的频率特性。基本程序为:fb=1000;fc=1
23、200;Ap=1;As=15;Fs=22050;ws=2*pi*fc/Fs;wp=2*pi*fb/Fs;N=ceil(12*pi/(ws-wp);Wn=(fb+fc)/Fs;b=fir1(N,Wn,blackman(N+1);h1,w1=freqz(b,1);subplot(111);plot(w1*Fs/(2*pi),20*log10(abs(h1);grid;分析显示如图12所示:图12IR低通滤波器的性能分析读入信号对信号傅里叶变换进行频谱分析,利用fftfilt函数滤波,对得到的信号傅里叶变换进行频谱分析,与滤波前进行比较。滤波前后的波形与频谱图分别为图13波前后信号的波形分析图14波
24、前后信号的频谱分析分析滤波前后信号的变化,波形图变窄,频谱变化表明:信号的低频率段被保留,高频率段被滤除。分析滤波前后的声音变化,滤波后声音明显变低而沉闷,这说明滤波器设计基本符合指标要求。 3.4.2 FIR高通滤波器设计fb=5000;fs=4800;Ap=1;As=15;Fs=22050;ws=2*pi*fb/Fs;wp=2*pi*fs/Fs;N=2*ceil(12*pi/(ws-wp);Wn=(fb+fs)/Fs;b=fir1(N,Wn,high,blackman(N+1);subplot(111);h1,w1=freqz(b,1);plot(w1/pi,20*log10(abs(h1
25、);grid;分析显示如图15所示。图14FIR高通滤波器的性能分析读入信号对信号傅里叶变换进行频谱分析,利用fftfilt函数滤波,对得到的信号傅里叶变换进行频谱分析,与滤波前进行比较。并在滤波前后用sound函数播放语音信号,观察滤波前后的变化。结果如图16、图17所示。图15波前后信号的波形分析图16波前后信号的频谱分析图14滤波前后信号的波形分析分析滤波前后信号的变化,波形图变窄,频谱变化表明:信号的高频率段被保留,低频率段被滤除。分析滤波前后的声音变化,滤波后声音明显变高而尖锐,这说明滤波器设计基本符合指标要求。3.4.3 FIR带通滤波器设计fb1=1200;fb2=3000;fc
26、1=1000;fc2=3200;Ap=1;As=15;Fs=22050;ws1=2*pi*fb1/Fs;wp1=2*pi*fc1/Fs;ws2=2*pi*fb2/Fs;wp2=2*pi*fc2/Fs;N=ceil(12*pi/(ws1-wp1);Wn1=(fb1+fc1)/Fs;Wn2=(fb2+fc2)/Fs;Wn=Wn1 Wn2;b=fir1(N,Wn,blackman(N+1);subplot(111);h1,w1=freqz(b,1);plot(w1*Fs/(2*pi),20*log10(abs(h1);grid;分析显示如图17示:图17IR带通滤波器的性能分析滤波前后的时域波形与频
27、谱如图18、19所示:图18波前后信号的波形分析图19波前后信号的频谱分析分析滤波前后信号的变化,波形图变窄。频谱变化表明:信号的低频率段和高频率段被滤除,中间频率段被保留;分析滤波前后的声音变化,滤波后声音一定程度上变的尖锐,与高通滤波器滤波后的声音相比较低。这说明滤波器设计基本符合指标要求。 开始3.4.4 FIR滤波器设计整体流程读入数字滤波器技术指标将指标转换成模拟滤波器的指标选择合适的窗函数设计滤波器将语音信号输入滤波器中对信号进行FFT变换并显示相应结果 结束图21 程序流程图4 总结通过这次课设,我认识到了MATLAB功能非常的强大,使得我们在使用的时候用户直接调用这些库函数并赋
28、予实际参数就能解决实际问题,具有极高的变成效率。我也熟悉了MATLAB的工作环境,可以很熟练的对MATLAB进行常规的操作,快速进行程序编辑和仿真。本文次课设通过一个设计实例,利用MATLAB实现 IIR数字滤波器设计与滤波的三种方法,从仿真结果可以看出它们均可以达到技术指标要求,而且方法简单、快捷,大大减轻了工作量。滤波器的设计工作完成后,可以借助于MATLAB的export操作导出所设计滤波器的系统函数H(z)。由于MATLAB具有强大的接口功能,仿真后的结果可以很方便的移植到 DSP、CPLD。在实际应用中,只需按要求修改滤波器参数,并对程序作较少的改动,即可实现不同截止频率的FIR滤波
29、器,实用性较强。完成任务后才发现这次训练的内容并不算难。首先要做的工作就是如何把设计条件写入到MATLAB中,其次就是对条件进行分析,主要是得出一些参数。最后就是对程序经行仿真,这里采用的是FDATool, Simulink仿真。以上的内容全都要通过程序的编写来完成,不过主要的操作是对函数的调用,这就要我们熟悉MATLAB的一些函数的基本功能以及它们的基本语句。当然,经过了课程设计,我也发现了自己的很多不足。但是通过自己的动手动脑,既增加了知识,又给了我专业知识以及专业技能上的提升,我也会更加努力,认真学习,争取在以后的课程中做得更好!参考文献1范寿康.DSP技术与DSP芯片. 北京:电子工业
30、出版社,20082万永革.数字信号处理的MATLAB实现. 北京:科学出版社, 20074赵健,李毅.数字信号处理学习辅导. 北京:清华大学出版社出版, 20065丁玉美.数字信号处理学习指导与题解. 北京:电子工业出版社出版, 20076丛玉良.数字信号处理原理及其MATLAB实现. 北京:电子工业出版社,20058王家文,王皓等.MATLAB7.0编程基础. 北京:机械工业出版社,2005附录%输入信号;x=wavread(d.wav);sound(x);%pause;fs=22050;ms2=floor(fs*0.002);ms10=floor(fs*0.001);ms20=floor(
31、fs*0.002);ms30=floor(fs*0.003);t=(0:length(x)-1)/fs;subplot(111);plot(t,x);legend(Waveform);xlabel(time(s);ylabel(Amplitude);sound(x);pause;%截取一段语音信号fs=22050;n1=10000;n2=30000;x,fs=wavread(d.wav,n1 n2);t=(0:length(x)-1)/fs;subplot(111);plot(t,x);legend(Waveform);xlabel(time(s);ylabel(Amplitude);soun
32、d(x);pause;%IIR 低通滤波器:clear all;fb = 1000;fc = 1200;fs = 22050;as = 15;ap = 1;wb = fb*2*pi/fs;wc = fc*2*pi/fs;Wb = 2*fs*tan(wb/2);Wc = 2*fs*tan(wc/2);n,wn=buttord(Wb,Wc,ap,as,s);b,a=butter(n,wn,s);bd,ad=bilinear(b,a,fs);h,w=freqz(bd,ad,1024); subplot(111);plot(w*fs/(2*pi),20*log10(abs(h);grid;title(
33、IIR低通滤波器的性能分析); pause;x1,fs,bits=wavread(d.wav);x=x1(10000:20000);n=length(x);f=fs*(0:(n/2-1)/n;sound(x);X=fft(x);z=filter(bd,ad,x);subplot(211);plot(x);title(原始信号的波形);subplot(212);plot(z);title(IIR低通滤波后信号的波形);pause;sound(z);待添加的隐藏文字内容1subplot(211);plot(f,abs(X(1:n/2);title(原始信号的频谱);xlabel(Hz);Z=fft
34、(z);subplot(212);plot(f,abs(Z(1:n/2); title(IIR低通滤波后的信号频谱);xlabel(Hz);pause;%IIR 高通滤波器:clear all;fb = 5000;fs = 4800;Fs = 22050;as = 15;ap = 1;wb = fb*2*pi/Fs;ws = fs*2*pi/Fs;Wb = 2*Fs*tan(wb/2);Ws = 2*Fs*tan(ws/2);n,wn=buttord(Wb,Ws,ap,as,s);b,a=butter(n,wn,high,s);bd,ad=bilinear(b,a,Fs);h,w=freqz(
35、bd,ad,1024); subplot(111);plot(w*Fs/(2*pi),20*log10(abs(h);grid;title(IIR高通滤波器的性能分析); pause;x1,Fs,bits=wavread(d.wav);x=x1(10000:20000)n=length(x);f=Fs*(0:(n/2-1)/n;sound(x);X=fft(x);z=filter(bd,ad,x);Z=fft(z);subplot(211);plot(x);title(原始信号的波形);subplot(212);plot(z);title(IIR高通滤波后的信号波形);pause;sound(
36、z);subplot(211);plot(f,abs(X(1:n/2);title(原始信号的频谱);xlabel(Hz);subplot(212);plot(f,abs(Z(1:n/2); title(IIR高通滤波后的信号频谱);xlabel(Hz);pause;%IIR 带通滤波器:clear all;fb1=1200;fb2=3000;fc1=1000;fc2=3200;fs =8000;as =15;ap = 1;wb =fb1 fb2*2*pi/fs;wc =fc1 fc2*2*pi/fs;Wb =2*fs*tan(wb/2);Wc =2*fs*tan(wc/2);n,wn=but
37、tord(Wb,Wc,ap,as,s);b,a=butter(n,wn,s);bd,ad=bilinear(b,a,fs);h,w=freqz(bd,ad);subplot(111);plot(w*fs/(2*pi),20*log10(abs(h);grid; title(IIR带通滤波器的性能分析); pause;x1,ft,bits=wavread(d.wav);x=x1(10000:20000);n=length(x);f=fs*(0:(n/2-1)/n;sound(x);X=fft(x);z=filter(bd,ad,x);Z=fft(z);subplot(211);plot(x);t
38、itle(原始信号的波形);subplot(212);plot(z);title(IIR带通滤波后的信号波形);pause;sound(z);subplot(211);plot(f,abs(X(1:n/2);title(原始信号的频谱);xlabel(Hz);subplot(212);plot(f,abs(Z(1:n/2); title(IIR带通滤波后的信号频谱);xlabel(Hz);pause;%FIR 低通滤波器:clear all;fb=1000;fc=1200;Ap=1;As=15;Fs=22050;ws=2*pi*fc/Fs;wp=2*pi*fb/Fs;N=ceil(12*pi/
39、(ws-wp);Wn=(fb+fc)/Fs;b=fir1(N,Wn,blackman(N+1);h1,w1=freqz(b,1);subplot(111);plot(w1*Fs/(2*pi),20*log10(abs(h1);grid;title(FIR低通滤波器的性能分析);pause;x1,Fs,bits=wavread(d.wav);x=x1(10000:20000)n=length(x);sound(x);X=fft(x);z=fftfilt(b,x);Z=fft(z);f=Fs*(0:(n/2-1)/n;subplot(211);plot(x);title(原始信号的波形);subp
40、lot(212);plot(z);title(FIR低通滤波后信号的波形);pause;sound(z);subplot(211);m=length(X)plot(f,abs(X(1:n/2);title(原始信号的频谱);xlabel(Hz);subplot(212);plot(f,abs(Z(1:n/2); title(FIR低通滤波后的信号频谱);xlabel(Hz);pause;%FIR 高通滤波器:clear all;fb=5000;fs=4800;Ap=1;As=15;Fs=22050;ws=2*pi*fb/Fs;wp=2*pi*fs/Fs;N=2*ceil(12*pi/(ws-w
41、p);Wn=(fb+fs)/Fs;b=fir1(N,Wn,high,blackman(N+1);subplot(111);h1,w1=freqz(b,1);plot(w1*Fs/(2*pi),20*log10(abs(h1);grid;title(FIR高通滤波器的性能分析); pause;x1,Fs,bits=wavread(d);x=x1(10000:20000)n=length(x);f=Fs*(0:(n/2-1)/n;sound(x);X=fft(x);z=fftfilt(b,x);Z=fft(z);subplot(211);plot(x);title(原始信号的波形);subplot(212);plot(z);title(FIR高通滤波后的信号的波形);pause;sound(z);subplot(211);plot(f,abs(X(1:(n-1)/2);title(原始的信号频谱);xlabel(Hz);subplot(