用MATLAB结合窗函数法设计数字带通FIR滤波器

1、武汉理工大学Matlab课程设计报告目录摘要IAbstractII1 原理说明11.1 数字滤波技术11.2 FIR滤波器11.3 窗函数21.4 MATLAB简介31.5 MATLAB结合窗函数设计法原理32 滤波器设计22.1 滤波器设计要求22.2 设计函数的选取22.3 窗函数构造32.4 设计步骤52.5 利用MATLAB自带函数设计63 滤波器测试103.1 滤波器滤波性能测试103.2 滤波器时延测量113.3 滤波器稳定性测量125 参考文献15附件一：16 摘要现代图像、语声、数据通信对线性相位的要求是普遍的。数字滤波器是一种用来过滤时间离散信号的数字系统，通过对抽样数据进行

2、数学处理来达到频域滤波的目的。根据其单位冲激响应函数的时域特性可分为两类：无限冲激响应（IIR）滤波器和有限冲激响应（FIR）滤波器。与IIR滤波器相比，FIR的实现是非递归的，总是稳定的；更重要的是，FIR滤波器在满足幅频响应要求的同时，可以获得严格的线性相位特性。因此，具有线性相位的FIR数字滤波器在高保真的信号处理，如数字音频、图像处理、数据传输、生物医学等领域得到广泛应用。在实际进行数字信号处理时，往往需要把信号的观察时间限制在一定的时间间隔内，只需要选择一段时间信号对其进行分析。这样，取用有限个数据，即将信号数据截断的过程，就等于将信号进行加窗函数操作。而这样操作以后，常常会发生频谱

3、分量从其正常频谱扩展开来的现象，即所谓的“频谱泄漏”。当进行离散傅立叶变换时，时域中的截断是必需的，因此泄漏效应也是离散傅立叶变换所固有的，可以通过窗函数加权抑制DFT的等效滤波器的振幅特性的副瓣，或用窗函数加权使有限长度的输入信号周期延拓后在边界上尽量减少不连续程度的方法实现抑制。FIR滤波器的设计中，为获得有限长单位取样响应，需要用窗函数截断无限长单位取样响应序列。另外，在功率谱估计中也要遇到窗函数加权问题。由此可见，窗函数加权技术在数字信号处理中的重要地位。此次课设基于MATLAB环境下，用窗函数法设计FIR滤波器的过程和设计实例。通过利用不同的窗函数方法设计FIR滤波器，对所设计的滤波

4、器进行分析比较，得出各种方法设计的滤波器的优缺点及其不同的使用场合，从而可以在设计滤波器时能够正确的选择FIR数字滤波器的窗函数的选取及设计方法。关键词：数字滤波器，MATLAB，有限冲激响应，窗函数，线性相位13AbstractA modern image, Voice and data communications to the linear phase requirement is universal. Digital filter is a time discrete signal which is used to filter digital system, through the

5、mathematical treatment of the sample data to achieve the frequency-domain filtering purposes. In accordance with its unit impulse response function of the time-domain characteristics can be divided into two categories: Infinite Impulse Response (IIR) filter and finite impulse response (FIR) fil

6、ters. Compared with the IIR filters, FIR realization of a non-recursive, is always stable; more importantly, FIR filters to meet the requirements of the amplitude-frequency response at the same time, get a strictly linear phase characteristics. Therefore, it is high-fidelity signal processing,

7、such as digital audio, image processing, data transmission, biomedical areas such as widely used.In practice, digital signal processing, they often need to observe the time limit for the signal in a certain time interval, only need to select a period of time to analyze the signals. In this way, acce

8、ss to a finite number of data about the process of signal data truncation is equivalent to the signal, the additional window function operation. And this after the operation, often occur from the normal component of the spectrum spread spectrum open to the phenomenon of so-called "spectral leak

9、age." When it comes to discrete Fourier transform, the time domain truncation is necessary, therefore, leakage effects are also inherent in discrete Fourier transform, and must be inhibited.Under the MATLAB environment, with the window function method for FIR filter design process and design ex

10、amples. Through the use of different window functions designed FIR filter, the filter designed to analyze and compare a variety of methods designed to filter obtained the advantages and disadvantages of the use of their different occasions, so you can design a filter to correct the choice

11、of FIR digital filter function selection window and the design method.Key words: Digital Filters, MATLAB, Finite Impulse Response, Window FunctionMatlab课程设计报告1 原理说明随着信息时代到来，数字信号处理已成为一门极其重要的学科和技术，尤其是在通信、语音、图像、自动控制等领域得到了广泛应用。在数字信号处理中，数字滤波器占有极其重要的地位，它具有精度高、可靠性好、灵活性大等特点。现代数字滤波器可以用软件或硬件两种方式实现。软件方式实现的优点是可

12、以通过滤波器参数的改变调整滤波器的性能。1.1 数字滤波技术数字滤波器可以理解为是一个计算程序或算法，将代表输入信号的数字时间序列转化为代表输出信号的数字时间序列，并在转化过程中使信号按预定的形式变化。根据数字滤波器冲激响应的时域特征，可将数字滤波器分为两种:无限长冲激响应（IIR）滤波器和有限长冲激响应（FIR）滤波器。IIR数字滤波器具有无限宽的冲激响应，与模拟滤波器相匹配。所以IIR滤波器的设计可以采取在模拟滤波器设计的基础上进一步变换的方法。本文主要对FIR滤波器加以介绍。1.2 FIR滤波器FIR（Finite Impulse Response）滤波器，即有限长单位冲激响应滤波器，是

13、数字信号处理系统中最基本的元件。它可以在保证任意幅频特性的同时具有严格的线性相频特性，而且其单位冲激响应是有限的，没有输出到输入的反馈，是稳定的系统。因此，FIR滤波器在通信、图像处理、模式识别等领域有着广泛应用。FIR滤波器的硬件实现主要有数字集成芯片、DSP芯片FIR滤波器、可编程FIR滤波器，后两者主要通过MATLAB软件进行设计。设计方法多样，设计形式灵活，能够满足各种要求，而且不受数字集成芯片规格的限制。FIR滤波器的设计方法主要有窗函数法、多带和过渡带、约束最小二乘法、任意相应法、升余弦法，其中最常用的是窗函数法。1.3 窗函数窗函数法是设计FIR滤波器的最主要方法之一。实际中遇到

14、的离散时间信号总是有限长的，因此不可避免要遇到数据截短问题。在信号处理中，对离散序列的截短是通过序列与窗函数相乘来实现的。常用的窗函数：矩形窗、巴特利特（Bartlett）窗、三角窗、海明（Hamming）窗、汉宁（Hanning）窗、布莱克曼（Blackman）窗、切比雪夫（Chebyshev）窗、凯泽（Kaiser）窗。下面就矩形窗函数、三角窗函数和广义余弦窗函数三种窗函数展开介绍。（1）矩形窗函数矩形窗(Rectangular Window)函数的时域形式可以表示为： （1-1）频域特性为 （1-2）Boxcar函数：生成矩形窗调用方式w = boxcar (n)：输入参数n是窗函数的长

15、度；输出参数w是由窗函数的值组成的n阶向量。从功能上讲，该函数又等价于w = ones(n,1)。（2）三角窗三角窗（Triangular Window）函数时域形式可表示为： （1-3）窗谱为： （1-4）式中，当N远大于1时，此时，窗谱主瓣宽度为8/N。（3）广义余弦窗汉明窗和布莱克曼窗，都可以用一种通用的形式表示，这就是广义余弦窗。这些窗都是广义余弦窗的特例，汉明窗又被称为改进的升余弦窗，而布莱克曼窗又被称为二阶升余弦窗。采用这些窗可以有效地降低旁瓣的高度，但是同时会增加主瓣的宽度。这些窗都是频率为0、2/(N1)和4/(N1)的余弦曲线的合成，其中N为窗的长度。通常采用下面的命令来生成

16、这些窗： （1-5） （1-6）其中，A、B、C适用于自己定义的常数。根据它们取值的不同，可以形成不同的窗函数，分别是：汉明窗 A=0.54，B=0.54，C=0；布莱克曼窗 A=0.5，B=0.5，C=0.08。1.3.1 窗函数绘图比较在MATLAB中运行以下代码：%常用窗函数绘图比较n=40;x=1:40;juxing=boxcar(n); %构造矩形窗sanjiao=triang(n); %构造三角窗hanming=hamming(n); %构造汉明窗bulaikeman=blackman(n); %构造布莱克曼窗kaize=kaiser(n); %构造凯泽窗plot(x,juxing

17、,'g.',x,sanjiao,'rx',x,hanming,'b+',x,bulaikeman,'kd',x,kaize,'c*');legend('矩形窗','三角窗','汉明窗','布莱克曼窗','凯泽窗');运行结果如图2-1所示：图2-15种窗函数绘图比较1.4 MATLAB简介MATLAB是矩阵实验室（MatrixLaboratory）之意。除具备卓越的数值计算能力外，它还提供了专业水平的符号计算，文字处理，可视化建模仿真和

18、实时控制等功能。使用MATLAB这一软件来设计滤波器,与传统的设计思路相比是非常简便的。在MATLAB环境下,使用者无需了解程序内部具体编写情况,即可输入期望的滤波器性能指标以实现系统仿真，完成对数字滤波器频率响应的处理和波形的检测。1.5 MATLAB结合窗函数设计法原理目前FIR滤波器的设计方法主要有三种：窗函数法、频率取样法和切比雪夫等波纹逼近的最优化设计方法。常用的是窗函数法和切比雪夫等波纹逼近的最优化设计方法。因此，设计FIR滤波器的方法之一可以从时域出发，截取有限长的一段冲击响应作为H(z)的系数，冲击响应长度N就是系统函数H(z)的阶数。只要N足够长，截取方法合理，总能满足频域要

19、求。通常情况下，这种时域设计、频域检验的方法要反复几个回合才能成功。要设计一个线性相位的FIR数字滤波器，首先要求理想频率响应。是w的周期函数，周期为，可以展开成傅氏级数： （1-7）使用上述的传递函数去逼近，一个理想的频率响应的傅立叶反变换： （1-8）其中是与理想频响对应的理想单位抽样响应序列，但不能用来作为设计FIR DFT用的，因为一般都是无限长、非因果的，物理上无法实现。为了设计出频响类似于理想频响的滤波器，可以考虑用来近似。窗函数的基本思想：先选取一个理想滤波器（它的单位抽样响应是非因果、无限长的），再截取（或加窗）它的单位抽样响应，从而得到线性相位因果FIR滤波器。这种方法的重点

20、是选择一个合适的窗函数和理想滤波器。设是一个长序列，是长度为N的窗函数，用截断，得到N点序列，即 （1-9）在频域上则有 （1-10）由此可见，窗函数不仅仅会影响原信号在时域上的波形，而且也会影响到频域内的形状。MATLAB主要提供了以下几种窗函数，如表1-1所示：表1-1 MATLAB常见窗函数窗类型窗函数旁瓣峰值主瓣峰值最小阻带衰减矩形窗Boxcar13dB4/M21dB三角窗Triang25dB8/M25dB汉宁窗Hamming31dB8/M44dB海明窗Blackman41dB8/M53dB凯泽窗Kaiser57dB12/M74dB加矩形窗后的频谱和理想频谱可得到以下结论：加窗使过渡带

21、变宽，过渡带的带宽取决于窗谱的主瓣宽度。矩形窗情况下的过渡带宽是。N越大，过渡带越窄、越陡。过渡带两旁产生肩峰，肩峰的两侧形成起伏振荡。肩峰幅度取决于窗谱主瓣和旁瓣面积之比。矩形窗情况下是8.95，与N无关。工程上习惯用相对衰耗来描述滤波器，相对衰耗定义为： （1-11）这样两个肩峰点的相对衰耗分别是0.74dB和-21dB。其中（-0.0895）对应的点的值定义为阻带最小衰耗。以上分析可见，滤波器各项重要指标都是由窗函数决定，因此改进滤波器的关键在于改进窗函数。窗函数谱两个最重要的指标是：主瓣宽度和旁瓣峰值衰耗。旁瓣峰值衰耗定义为： 旁瓣峰值衰耗20*lg(第一旁瓣峰值主瓣峰值) （1-12

22、）为了改善滤波器的性能，需使窗函数谱满足：(1)主瓣尽可能窄，以使设计出来的滤波器有较陡的过渡带。(2)尽量减少最大旁瓣的相对幅度，也就是能量集中于主瓣，以减小带内、带外波动的最大幅度，增大阻带衰减。一般来说，以上两点很难同时满足。当选取主瓣宽度很窄时，旁瓣的分量势必增加，从而带内、带外的波动也增加了；当选取最小的旁瓣幅度时，降低了带内、带外的波动，但是过渡带的陡度减小了。因此实际采用的窗函数的特性往往是二者的折中，即在保证主瓣宽度达到一定要求的前提下，适当牺牲主瓣宽度来换取旁瓣波动的减小。2 滤波器设计2.1 滤波器设计要求利用MATLAB仿真软件系统结合窗函数法设计一个数字带通FIR滤波器

23、。要求：分别使用矩形窗、三角形窗、汉明窗、布莱克曼窗、凯泽窗进行设计，并输出滤波器的频率特性。技术指标：1 采样频率为40kHz；2 通带边缘频率：fp1=4kHz,fp2=6kHz;3 通带峰值起伏：p<1dB;4 阻带边缘频率：fs1=2.5k,fs2=7.5k;5 最小阻带衰减：As>40dB2.2 设计函数的选取MATLAB信号处理工具箱提供了基于加窗的线性相位FIR滤波器设计函数fir1和fir2，fir1函数的调用格式为：b=fir1(n,Wn)；b=fir1(n,Wn,'ftype')；b=fir1(n,Wn,window)b=fir1(n,Wn,&#

24、39;ftype',window)；b=fir1(.,'normalization')函数参数说明如下表：表2-1 Matlab函数参数说明函数参数参数说明n滤波器的阶数'ftype'设计滤波器的类型'high'高通滤波器'stop'带阻滤波器'DC-1'多通带滤波器，第一频带为通带'DC-0'多通带滤波器，第一频带为阻带；默认时为低通或带通滤波器'window'窗函数，是长度为n+1的列向量，默认时函数自动取Hamming窗该函数实现加窗的线性相位FIR滤波器设计，可以设计

25、标准低通、带通、高通和带阻滤波器。（具有任意频率响应的加窗滤波器可以采用fir2进行设计）2.3 窗函数构造MATLAB工具箱已经提供了各种窗函数的构造函数，因而窗函数的构造十分方便。下面给出几种常用窗函数的构造方法：1矩形窗：利用w=boxcar(n)的形式得到窗函数，其中n为窗函数的长度，而返回值w为一个n阶的向量，它的元素由窗函数的值组成。w=boxcar(n)等价于w=ones(1,n).2三角窗：利用w=triang(n)的形式得到窗函数，其中n为窗函数的长度，而返回值w为一个n阶的向量，它的元素由窗函数的值组成。w=triang(N-2)等价于bartlett(N)。3海明窗：利用

26、w=hamming(n)得到窗函数，其中n为窗函数的长度，而返回值w为一个n阶的向量，包含了窗函数的n个系数。它和汉宁窗的主瓣宽度相同，但是它的旁瓣进一步被压低。4布莱克曼窗：利用w=blackman(n)得到窗函数，其中n为窗函数的长度，而返回值w为一个n阶的向量，包含了窗函数的n个系数。它的主瓣宽度是矩形窗主瓣宽度的3倍，为12*pi/N，但是它的最大旁瓣值比主瓣值低57dB。5凯泽窗：利用w=kaiser(n,beta)的形式得到窗函数。2.3.1 窗函数设计条件在使用窗函数设计滤波器时要满足以下两个条件：1 窗谱主瓣尽可能地窄，从而可以获得较陡峭的过渡带；2 尽量减少窗谱的最大旁瓣的相

27、对幅度，及尽可能是能量集中于主瓣，减少峰肩和波纹，进而增加阻带的衰减。2.4 设计步骤鉴于MATLAB内置有丰富函数，还集成了通信设计的工具箱，两者均有专用于滤波器设计的工具。从效率和方便性出发，本文采用MATLAB内置函数进行滤波器相关设计。图2-1 程序设计流程图2.5 利用MATLAB自带函数设计利用原理说明介绍的fir1函数进行设计，该设计方法只需要给出滤波器的阶数，截止频率，窗函数等参数，MATLAB即可自行完成设计，并可通过freqz函数查看滤波器的幅频响应和相频响应，以验证滤波器是否满足设计要求。下面给出利用fir1函数设计的不同窗函数的数字带通滤波器。1. 利用矩形窗进行设计%

28、矩形窗函数的数字带通滤波器fs=40000; %设定采样频率fp1=4000;fp2=6000; %第一截止频率fs1=2500;fs2=8500; %第二截止频率As=40; %最小阻带衰减Ws1=(fp1+fs1)/fs;Ws2=(fp2+fs2)/fs; %截止频率归一化处理w=(fp1-fs1)/fs; %求归一化过渡带M=ceil(As-7.95)/(14.36*w) %计算所需滤波器的阶数juxing=boxcar(M+1); %生成长度为M+1的矩形窗boxb=fir1(M,Ws1,Ws2,juxing); %生成矩形窗设计的fir滤波器freqz(boxb,1,fs,fs);

29、%绘制幅频和相频响应曲线图2-2 矩形窗设计的FIR滤波器仿真结果图从幅频响应上看，通带基本无波纹，阻带中波纹较大，因此阻带较不理想，相频响应曲线在通带内为直线，效果较好，信号失真小。 2. 利用三角窗进行设计利用三角窗进行设计时，原理与矩形窗基本相同，生成窗函数时选用triang（）函数生成三角窗，程序运行结果如下：图2-4三角窗设计的FIR滤波器幅频和相频响应曲线3. 利用汉明窗进行设计利用汉明窗进行设计时，原理与矩形窗基本相同，生成窗函数时选用hamming（）函数生成三角窗，程序运行结果如下：图2-5汉明窗设计的FIR滤波器幅频和相频响应曲线4. 利用布莱克曼窗进行设计利用布莱克曼窗进

30、行设计时，原理与矩形窗基本相同，生成窗函数时选用blackman（）函数生成三角窗，程序运行结果如下：图2-6布莱克曼窗设计的FIR滤波器幅频和相频响应曲线5. 利用凯泽窗进行设计利用凯泽窗进行设计时，滤波器的参数可以用以下函数求得：调用函数n,wn,bta,ftype=kaiserord（f，a，dev，fs）利用该函数修改代码得到凯泽窗设计FIR滤波器的代码如下：fs=40000; %设定采样频率fp1=4000;fp2=6000; %第一截止频率fs1=2500;fs2=8500; %第二截止频率n,wn,bta,ftype=kaiserord(fs1,fp1,fp2,fs2,. 0,1

31、,0,0.01 0.1087 0.01,fs) %求滤波器参数b=fir1(n,wn,ftype,kaiser(n+1,bta); %生成fir滤波器freqz(b,1,fs,fs); %绘制幅频和相频响应曲线程序运行结果如下：图2-7凯泽窗设计的FIR滤波器幅频和相频响应曲线3 滤波器测试利用MATLAB中提供的filter函数可以选择不同的滤波器对数字信号进行滤波，这里主要对滤波器的滤波性能进行简单测试，这里主要针对滤波器的滤波性能。3.1 滤波器滤波性能测试由MATLAB模拟生成含有不同频率的数字信号，然后选取设计的滤波器（汉明窗函数FIR滤波器）对数字信号进行滤波，为方便观察，模拟生成

32、的信号只含有包含在阻带的两个频率（1KHz，8KHz）一个包含于通带的频率（5KHz），测试代码如下：%汉明窗函数FIR滤波器性能测试fs=40000; %设定采样频率fp1=4000;fp2=6000; %第一截止频率fs1=2500;fs2=8500; %第二截止频率As=40; %最小阻带衰减Ws1=(fp1+fs1)/fs;Ws2=(fp2+fs2)/fs; %截止频率归一化处理w=(fp1-fs1)/fs; %求归一化过渡带M=ceil(As-7.95)/(14.36*w) %计算所需滤波器的阶数hanming=hamming(M+1); %生成长度为M+1的汉明窗窗hanb=fir1(M,Ws1,Ws2,hanming); %生成汉明窗设计的fir滤波器t=0:1/fs:2; x=sin(2*pi*2000*t)+sin(2*pi*5000*t). +sin(2*pi*8000*t); %生成混合信号xo=filter(hanb,2,x); %用滤波器对信号进行滤波figure;nn=5000:5100; %取一段信号subplot(211);tt=nn/fs;plot