课程设计任务书
学生姓名: 专业班级:
指导教师: 工作单位: 信息工程学院
题 目:FIR带阻滤波器的设计
初始条件:
具备数字信号处理的理论知识;
具备Matlab编程能力;
熟悉带阻滤波器的设计原理;
提供编程所需要的计算机一台
要求完成的主要任务:(包括课程设计工作量及其技术要求,以及说明书撰写等具体要求)
1、设计中心频率为200Hz,带宽为150Hz的FIR数字带阻滤波器,通带最大衰减 1dB 阻带最小衰减 60dB
2、独立编写程序实现
3、完成符合学校要求的设计说明书
时间安排:
一周,其中3天程序设计,2天程序调试
指导教师签名: 年 月 日
系主任(或责任教师)签名: 年 月 日
目录
摘要... 1
Abstract 2
1、 绪论... 3
2、MATLAB简介... 4
3、设计原理... 5
3.1、线性相位FIR滤波器的特点... 5
3.2、窗函数设计法原理... 6
3.3、实验中用到的的窗函数... 9
3.3.1 布莱克曼窗函数... 9
3.2.1 凯塞窗... 9
4、方案设计与分析... 11
4.1布莱克曼窗函数设计... 11
4.2凯塞窗函数设计... 14
5、心得体会... 16
6、参考文献... 17
本科生课程设计成绩评定表... 18
摘要
无限长脉冲数字滤波器的设计方法只考虑了幅度特性,没有考虑相位特性,所设的滤波器一般是某种确定的非线性相位特性。有限脉冲响应(FIR)滤波器在保证了幅度特性满足技术要求的同时,很容易做到有严格的线性相位特性。
本课题利用MATLAB软件实现。MATLAB是“矩阵实验室”(MATrix LABoratoy)的缩写,是一种科学计算软件,它使用方便,输入简捷,运算高效,内容丰富,因此利用MATLAB软件,通过一系列较为系统的函数法,根据已知的技术指标,就可以设计出满足要求的滤波器。
关键字:MATLAB 窗函数 FIR带阻数字滤波器 线性相位
Abstract
The design method of the digital filter with infinite length has only considered the amplitude characteristic, without considering the phase characteristic. The filter is usually a certain nonlinear phase characteristic. The finite impulse response (FIR) filter ensures that the amplitude characteristics meet the technical requirements, it is easy to achieve a strict linear phase characteristic.This topic uses the MATLAB software to realize. Matlab is abbreviation of "Matrix Laboratory (matrix LABoratoy), is a kind of scientific calculation software, it is convenient to use, the input is simple, computationally efficient, rich in content. Therefore, the MATLAB software, through a series of systematic function method, according to the known specifications can be designed to meet the requirements of the filter.
Keywords: MATLAB window function FIR band stop digital filter linear phase
1、 绪论
现代图像、语声、数据通信对线性相位的要求是普遍的。正是此原因,使得具有线性相位的FIR数字滤波器得到大力发展和广泛应用。
在实际进行数字信号处理时,往往需要把信号的观察时间限制在一定的时间间隔内,只需要选择一段时间信号对其进行分析。这样,取用有限个数据,即将信号数据截断的过程,就等于将信号进行加窗函数操作。而这样操作以后,常常会发生频谱分量从其正常频谱扩展开来的现象,即所谓的“频谱泄漏”。当进行离散傅立叶变换时,时域中的截断是必需的,因此泄漏效应也是离散傅立叶变换所固有的,必须进行抑制。而要对频谱泄漏进行抑制,可以通过窗函数加权抑制DFT的等效滤波器的振幅特性的副瓣,或用窗函数加权使有限长度的输入信号周期延拓后在边界上尽量减少不连续程度的方法实现。而在后面的FIR滤波器的设计中,为获得有限长单位取样响应,需要用窗函数截断无限长单位取样响应序列。另外,在功率谱估计中也要遇到窗函数加权问题。由此可见,窗函数加权技术在数字信号处理中的重要地位。
2、MATLAB简介
MATLAB是“矩阵实验室”(MATrix LABoratoy)的缩写,是一种科学计算软件,主要适用于矩阵运算及控制和信息处理领域的分析设计,它使用方便,输入简捷,运算高效,内容丰富,因此很多专家在自己擅长的领域用它编写了许多专门的MATLAB工具包,由于MATLAB功能的不断扩展,所以是科学研究中最常用必不可少的工具。
MATLAB由一系列工具组成。这些工具方便用户使用MATLAB的函数和文件,其中许多工具采用的是图形用户界面。包括MATLAB桌面和命令窗口、历史命令窗口、编辑器和调试器、路径搜索和用于用户浏览帮助、工作空间、文件的浏览器。随着MATLAB的商业化以及软件本身的不断升级,MATLAB的用户界面也越来越精致,更加接近Windows的标准界面,人机交互性更强,操作更简单。而且新版本的MATLAB提供了完整的联机查询、帮助系统,极大的方便了用户的使用。简单的编程环境提供了比较完备的调试系统,程序不必经过编译就可以直接运行,而且能够及时地报告出现的错误及进行出错原因分析。
MATLAB一个高级的矩阵/阵列语言,它包含控制语句、函数、数据结构、输入和输出和面向对象编程特点。用户可以在命令窗口中将输入语句与执行命令同步,也可以先编写好一个较大的复杂的应用程序(M文件)后再一起运行。新版本的MATLAB语言是基于最为流行的C++语言基础上的,因此语法特征与C++语言极为相似,而且更加简单,更加符合科技人员对数学表达式的书写格式。使之更利于非计算机专业的科技人员使用。而且这种语言可移植性好、可拓展性极强,这也是MATLAB能够深入到科学研究及工程计算各个领域的重要原因。
3、设计原理
3.1、线性相位FIR滤波器的特点
一个单一频率的正弦信号通过一个系统,假设它通过这个系统的时间需要t,则这个信号的输出相位落后原来信号wt的相位。从这边可以看出,一个正弦信号通过一个系统落后的相位等于它的w*t;反过来说,如果一个频率为w的正弦信号通过系统后,它的相位落后delta,则该信号被延迟了delta/w的时间。在实际系统中,一个输入信号可以分解为多个正弦信号的叠加,为了使得输出信号不会产生相位失真,必须要求它所包含的这些正弦信号通过系统的时间是一样的。因此每一个正弦信号的相位分别落后,w1*t,w2*t,w3*t。因此,落后的相位正比于频率w,如果超前,超前相位的大小也是正比于频率w。从系统的频率响应来看,就是要求它的相频特性是一条直线。在FIR滤波器的设计中,为了得到线性相位的性质,通常利用实偶对称序列的相频特性为常数0和实奇对称序列为相频特性为常数90度的特点。因此得到的是对称序列,不是因果序列,是不可实现系统,为了称为物理可实现系统,需要将它向右移动半个周期,这就造成了相移特性随时间的变化,同时也是线性变化。
单位脉冲响应h(n)(为实数)具有偶对称或奇对称性,则FIR数字滤波器具有严格的线性相位特性。
FIR滤波器的单位冲激响应h(n)是有限长的(0≤n≤N-1),其z变换为
(2.1)
显然,H(z)是z-1的(N-1)阶多项式,在有限z平面(0<|z|<∞)有(N-1)个零点,而有(N-1)阶极点全部位于z平面的原点z=0处。
h(n)的频率响应H(ejω)可表示为
(2.2)
当h(n)为实序列时,可将H(ejω)表示成
(2.3)
其中是真正的幅度响应,而实函数H(ω)称为幅度函数,θ(ω)称为相位函数。有两类准确的线性相位,要求分别满足
(2.4)
(2.5)
其中τ,β均为常数,表示相位是通过坐标原点ω=0或是通过θ(0)=β的斜直线,
二者的群延时都是常数。一般将满足式(2.4)的相位称为第一类线
性相位,满足式(2.5)的相位称为第二类线性相位。
3.2、窗函数设计法原理
数字滤波器可以理解为是一个计算程序或算法,将代表输入信号的数字时间序列转化为代表输出信号的数字时间序列,并在转化过程中,使信号按预定的形式变化。数字滤波器有多种分类,根据数字滤波器冲激响应的时域特征,可将数字滤波器分为两种,即无限长冲激响应(IIR)滤波器和有限长冲激响应(FIR)滤波器。IIR数字滤波器具有无限宽的冲激响应,与模拟滤波器相匹配。所以IIR滤波器的设计可以采取在模拟滤波器设计的基础上进一步变换的方法。FIR数字滤波器的单位脉冲响应是有限长序列。它的设计问题实质上是确定能满足所要求的转移序列或脉冲响应的常数问题,设计方法主要有窗函数法、频率采样法和等波纹最佳逼近法等。
因此设计FIR滤波器的方法之一可以从时域出发,截取有限长的一段冲击响应作为H(z)的系数,冲击响应长度N就是系统函数H(z)的阶数。只要N足够长,截取的方法合理,总能满足频域的要求。一般这种时域设计、频域检验的方法要反复几个回合才能成功。要设计一个线性相位的FIR数字滤波器,首先要求理想频率响应。是w的周期函数,周期为,可以展开成傅氏级数
= (2.6)
其中是与理想频响对应的理想单位抽样响应序列。但不能用来作为设计FIR-DF用的h(n),因为一般都是无限长、非因果的,物理上无法实现。为了设计出频响类似于理想频响的滤波器,可以考虑用h(n)来近似。
窗函数的基本思想:先选取一个理想滤波器(它的单位抽样响应是非因果、无限长的),再截取(或加窗)它的单位抽样响应得到线性相位因果FIR滤波器。这种方法的重点是选择一个合适的窗函数和理想滤波器。
设x(n)是一个长序列,w(n)是长度为N的窗函数,用w(n)截断x(n),得到N点序列an(n),即
an(n) = x(n)w(n) (2.7)
在频域上则有
(2.8)
由此可见,窗函数w(n)不仅仅会影响原信号x(n)在时域上的波形,而且也会影响到频域内的形状。
MATLAB信号工具箱主要提供了以下几种窗函数,如表2.1所示。
加矩形窗后的频谱和理想频谱可得到以下结论:
加窗使过渡带变宽,过渡带的带宽取决于窗谱的主瓣宽度。矩形窗情况下的过渡带宽是。N越大,过渡带越窄、越陡;
过渡带两旁产生肩峰,肩峰的两侧形成起伏振荡。肩峰幅度取决于窗谱主瓣和旁瓣面积之比。矩形窗情况下是8.95%,与N无关。工程上习惯用相对衰耗来描述滤波器,相对衰耗定义为
(2.9)
这样两个肩峰点的相对衰耗分别是0.74dB和-21dB。其中(-0.0895)对应的点的值定义为阻带最小衰耗。
以上的分析可见,滤波器的各种重要指标都是由窗函数决定,因此改进滤波器的关键在于改进窗函数。
表2.1 MATLAB窗函数
窗函数谱的两个最重要的指标是:主瓣宽度和旁瓣峰值衰耗。旁瓣峰值衰耗定义为
旁瓣峰值衰耗=20lg(第一旁瓣峰值/主瓣峰值) (2.10)
为了改善滤波器的性能,需使窗函数谱满足:
①主瓣尽可能窄,以使设计出来的滤波器有较陡的过渡带;
②第一副瓣面积相对主瓣面积尽可能小,即能量尽可能集中在主瓣,外泄少,使设计出来的滤波器的肩峰和余振小。
但上面两个条件是相互矛盾的,实际应用中,折衷处理,兼顾各项指标。
设计滤波器尽量要求窗函数满足以下两项要求:
(1)窗谱主瓣尽可能地窄,以获取较陡的过渡带。
(2)尽量减少窗谱的最大旁瓣的相对幅度。也就是能量尽量集中于主瓣,这样使尖峰和波纹减小,就可增大阻带的衰减。
但是这两项要求是不能同时满足的。当选用主瓣宽度较窄时,虽然得到陡峭的过渡带,但通带和阻带的波动明显增加;当选用最小的旁瓣幅度时,虽能得到平坦的幅度响应和较小的阻带波纹,但过渡带加宽,即主瓣会加宽。因此,实际所选用的窗函数往往是它们的折中。
3.3、实验中用到的的窗函数
3.3.1 布莱克曼窗函数
布莱克曼窗函数的时域形式可以表示为
(3-2-9)
它的频域特性为
(3-2-10)
其中,为矩形窗函数的幅度频率特性函数。
布莱克曼窗函数的最大旁瓣值比主瓣值低57dB,但是主瓣宽度是矩形窗函数的主瓣宽度的3倍,为12π/N。
Blackman函数:生成海明窗
调用方式
(1) w = blackman (n):输入参数n是窗函数的长度;输出参数w是由窗函数的值组成的n阶向量。
(2) w = blackman (n,sflag):参数sflag用来控制窗函数首尾的两个元素值;其取值为symmetric或periodic;默认值为symmetric。
3.2.1 凯塞窗
上面所讨论的几种窗函数,在获得旁瓣抑制的同时却增加了主瓣的宽度。而凯塞窗定义了一组可调的窗函数,它是由零阶贝塞尔函数构成的,其主瓣能量和旁瓣能量的比例是近乎最大的。而且,这种窗函数可以在主瓣宽度和旁瓣高度之间自由选择它们的比重,使用户的设计变得非常灵活。
凯塞窗函数的时域形式可表示为
(3-2-11)
其中,是第1类变形零阶贝塞尔函数,是窗函数的形状参数,由下式确定:
(3-2-12)
其中,为凯塞窗函数的主瓣值和旁瓣值之间的差值(dB)。改变β的取值,可以对主瓣宽度和旁瓣衰减进行自由选择。β的值越大,窗函数频谱的旁瓣值就越小,而其主瓣宽度就越宽。
Kaiser函数:生成凯塞窗
调用方式w = kaiser(n,beta):输入参数n是窗函数的长度;输入参数beta用于控制旁瓣的高度;输出参数w是由窗函数的值组成的n阶向量。n一定时,beta越大,其频谱的旁瓣就越小,但主瓣宽度相应的增加;当beta一定时,n发生变化,其旁瓣高度不会发生变化。
4、方案设计与分析
用窗函数法设计一个FIR带阻滤波器。指示如下:
下通带截至频率 ;
上通带截止频率 ;
阻带下限频率 ;
阻带上限频率 ;
通带最大衰减 ;
阻带最小衰减 =;
6中窗函数的基本参数如下
4.1布莱克曼窗函数设计
因为阻带最小衰减=,所以选择布莱克曼窗或凯塞窗都可以设计,先以布莱克曼窗进行设计。程序步骤如下:
flp=120;
fhp=280;
fls=130;
fhs=270;
fs=1000;
wlp=2*pi*flp/fs;
whp=2*pi*fhp/fs;
wls=2*pi*fls/fs;
whs=2*pi*fhs/fs;
wc=[(wlp+wls)/(2*pi),(whp+whs)/(2*pi)];
delta1=wls-wlp;
delta2=whp-whs;
delta_w=min(delta1,delta2);
N=ceil(11*pi/delta_w); %不同的窗要选择系数不同
N=N+rem(N,2);
n=0:N-1;
window=blackman(N+1); %选择窗函数
[h1,w]=freqz(window,1);
subplot(2,2,1)
stem(window,'.');
xlabel('n');
title(' blackman窗函数');
subplot(2,2,2)
plot(w*fs/(2*pi),20*log(abs(h1)/abs(h1(1))));
grid;
xlabel('f/Hz');
ylabel('幅度(dB)');
title(' blackman窗函数的频谱');
hn=fir1(N,wc,'stop',window);
[h2,w]=freqz(hn,1,512);
subplot(2,2,3)
stem(hn,'.');
xlabel('n');
ylabel('h(n)');
title(' blackman窗函数的单位脉冲响应');
subplot(2,2,4)
plot(w*fs/(2*pi),20*log(abs(h2)/abs(h2(1))));
grid;
xlabel('f/Hz');
ylabel('幅度(dB)');
title(' blackman带阻滤波器的幅频特性');
经过上述仿真运行结果如下:
图4.1布莱克曼窗
4.2凯塞窗函数设计
设计程序如下:
flp=120;
fhp=280;
fls=130;
fhs=270;
fs=1000;
wlp=2*pi*flp/fs;
whp=2*pi*fhp/fs;
wls=2*pi*fls/fs;
whs=2*pi*fhs/fs;
wc=[(wlp+wls)/(2*pi),(whp+whs)/(2*pi)];
delta1=wls-wlp;
delta2=whp-whs;
delta_w=min(delta1,delta2);
N=ceil(10*pi/delta_w); %不同的窗要选择系数不同
N=N+rem(N,2);
n=0:N-1;
window=kaiser(N+1); %选择窗函数
[h1,w]=freqz(window,1);
subplot(2,2,1)
stem(window,'.');
xlabel('n');
title(' kaiser窗函数');
subplot(2,2,2)
plot(w*fs/(2*pi),20*log(abs(h1)/abs(h1(1))));
grid;
xlabel('f/Hz');
ylabel('幅度(dB)');
title(' kaiser窗函数的频谱');
hn=fir1(N,wc,'stop',window);
[h2,w]=freqz(hn,1,512);
subplot(2,2,3)
stem(hn,'.');
xlabel('n');
ylabel('h(n)');
title(' kaiser窗函数的单位脉冲响应');
subplot(2,2,4)
plot(w*fs/(2*pi),20*log(abs(h2)/abs(h2(1))));
grid;
xlabel('f/Hz');
ylabel('幅度(dB)');
title(' kaiser带阻滤波器的幅频特性');
图4.2凯塞窗
5、心得体会
这次课设我们主要用布莱克曼窗和凯塞窗进行了FIR带阻滤波器的设计和制作,这次课设中主要用到matlab这种仿真软件,这个课设要用到很多新的函数,运用起来有一定的困难,不过最后通过查询一些资料,能较好地掌握这些知识。主要的困难在后面的理论部分。
设计带通滤波器时首先要计算出过渡带,然后查表得到不同窗函数所需要的阶数,不同的窗函数所设计的滤波器的形状各有差异,尤其在主瓣宽度、旁瓣的形状以及主瓣与旁瓣的高度差上有比较明显得差别,实际应用中应根据实际情况,折衷处理,兼顾各项指标。
为了这次课程设计,自己自学了数字信号处理领域中窗函数的有关知识。实际中遇到的离散时间信号总是有限长的,因此不可避免地要遇到数据截断问题。而在信号处理中,对离散序列的数据截断是通过序列与窗函数相乘来实现的。而且,有关滤波器的设计、功率谱估计等基本概念也要用到窗函数。本次课程设计对经常用到的下面6窗函数:矩形窗函数、三角窗函数、汉宁窗函数、哈明窗函函数、布莱克曼窗函数、凯塞窗函窗,先是做了基本概念上的阐释,然后对其MATLAB实现函数做出了说明,最后又结合具体的实例,对这些窗函数的频域特性等进行了介绍。
通过此次课设,我对FIR滤波器的有关知识有了更深入的认识,同时提高了读程序和写程序的能力,理论和实践都有了提高。课设做完后,也发现了自己的一些不足,平时很少自己动手写程序,以至于用的时候有很多困难,在以后的时间里,我会多看看他人写的程序,然后自己动手写一部分,提高自己的动手实践能力。
6、参考文献
[1]刘泉,厥大顺,郭志强编. 数字信号处理原理与实现. 北京:电子工业出版
社,2009.6.
[2]《MATLAB及在电子信息课程中的应用》,陈怀堔,吴大正,高西全.电子工业出版社,2006.
[3]《MATLAB 7.0从入门到精通》,求是科技.人民邮电出版社,2006.
[4]《数字信号处理(第三版)》学习指导,高西全,丁玉美.西安科技大学出版社,2001.
[5]候正信译. 数字信号处理基础. 北京:电子工业出版社,2003.8.
本科生课程设计成绩评定表
指导教师签字:
年 月 日