《基于MATLAB的FIR数字低通滤波器设计.doc》由会员分享,可在线阅读,更多相关《基于MATLAB的FIR数字低通滤波器设计.doc(22页珍藏版)》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、摘 要FIR数字滤波器是数字信号处理一个重要组成部分,由于FIR数字滤波器具有严格线性相位,因此在信息采集与处理过程中得到了广泛应用。本文介绍了FIR数字滤波器概念与线性相位条件,分析了窗函数法、频率采样法与等波纹逼近法设计FIR滤波器思路与流程。在分析三种设计方法原理基础上,借助Matlab仿真软件工具箱中fir1实现窗函数法中哈明窗设计FIR低通滤波器。关键词:FIR数字滤波器;线性相位;窗函数法;哈明窗;Matlab目 录前言1一 基本原理21.1 FIR数字滤波器概述21.2 FIR数字滤波器线性相位定义31.3 FIR数字滤波器线性相位时域条件3二 系统设计52.1 FIR数字滤波器
2、窗函数设计方法52.1.1 窗函数法设计思路52.1.2 常见窗函数介绍62.1.3 吉布斯效应82.2 FIR数字滤波器频率采样设计法92.3 FIR数字滤波器等波纹逼近设计法10三 详细设计123.1 程序设计流程123.2 Matlab简介12窗函数法Matlab实现133.3.1 fir1函数介绍13基于fir1函数窗函数法FIR滤波器设计14总结17参考文献18附录19致谢21第 20 页前 言随着信息科学与计算机技术不断发展,数字信号处理(DSP,Digital Signal Processing)理论与技术也得到了飞速发展,并逐渐成为一门重要学科,它重要性在日常通信、图像处理、遥
3、感、声纳、生物医学、地震、消费电子、国防军事、医疗方面等显得尤为突出。在我们面临信息革命中,数字信号处理几乎涉及了所有工程技术领域1。数字信号处理是一种将信号以数字形式进行处理一种理论与技术,它目是将真实世界中一些信号进行分析并滤波,最后得出其中有用信号。数字滤波器是数字信号处理一种,一般根据单位脉冲响应h(n)分为无限脉冲响应(IIR)与有限脉冲响应(FIR)系统。IIR数字滤波器设计方法简单,特别是采用双线性变换法来设计数字滤波器不存在频域混叠现象,但是IIR滤波器存在一个较为明显缺憾,就是它相位响应一般都是非线性,而在传输频带内相位响应如果不是线性,就会造成有用信号传输失真,而FIR数字
4、滤波器不仅可以设计成任意幅度响应,而且可以设计成在通频带内具有良好线性相位响应2。FIR数字滤波器单位脉冲响应h(n)有限长,所以FIR数字滤波器是稳定,不存在稳定性问题,且可以通过快速傅里叶变换(FFT)算法来实现信号滤波,大大提高运算效率。因此,FIR数字滤波器日益引起了人们关注。FIR数字滤波器设计方法有很多,比较常用有窗函数设计法、频率采样设计法、等波纹逼近法等。本课题通过运用窗函数设计FIR数字低通滤波器,并实现对给定信号进行滤波。窗函数设计法是最基本数字滤波方法,是利用傅里叶反变换(IDTFT)计算给定频响理想单位脉冲响应,再加以窗函数进行截断与平滑 3。Matlab软件信号处理工
5、具箱提供了FIR数字滤波器设计子函数,运用Matlab软件设计可以避免繁杂数学运算,而且具有丰富绘图功能,可以方便地查看所设计数字滤波器幅度响应与相位响应是否满足设计要求。因此,本课题在理论分析各种FIR数字滤波器设计方法基础上,运用Matlab软件进行仿真分析。一 基本原理1.1 FIR数字滤波器概述一般来说一个经典数字滤波器是一个线性时不变系统,其数学模型可以用Z域系统函数来表示: (1-1)其中均为滤波器参数。在式(1-1)中,当值不全为零值时,Z域系统函数必定含有一个或一个以上极值点,此时单位脉冲响应为无限长,对于一个稳定数字滤波器来说,Z域系统函数必须在单位圆内,因而把含有极值点Z域
6、系统函数数字滤波器称为无限脉冲响应数字滤波器(Infinite Impulse Response),即IIR数字滤波器4。而当值全为零时,Z域系统函数只有一个零点,式(1-1)表示系统函数可以写成: (1-2)公式(1-2)表明,FIR滤波器系统函数是阶多项式,在有限平面上有个零点,而在平面原点处有阶极点。公式(1-2)表示系统,其单位脉冲响应可以表示为: (1-3)在式(1-3)中,只有当,才有非零值,所以数字滤波器脉冲响应是有限长,因此在数字信号处理中把这种数字滤波器称为有限脉冲响应数字滤波器(Finite Impulese Response),即FIR数字滤波器。FIR数字滤波器最突出两
7、个优点是:(1)只要对附加一定条件,就很容易获得严格线性相位。(2)由于极点位于原点处,所以FIR数字滤波器不存在稳定性问题。1.2 FIR数字滤波器线性相位定义设FIR数字滤波器脉冲响应长度为N,则其频率响应可以表示为: (1-4)上式通过欧拉恒等式展开可得到相位特性,有两种线性相位特性,通常称为第一类线性相位与第二类线性相位。第一类线性相位特性: 是一个及无关常数第二类线性相位特性: 是起始相位严格地说第二种情况时是不具有线性相位特性,但上述两种情况都满足群延迟是一个常数,仍可以视为具有线性相位,在第二类线性相位中是常用一种情况5。1.3 FIR数字滤波器线性相位时域条件对于第一类线性相位
8、,即,通过一系列运算整理之后可得到一个三角函数求与公式: (1-5)式中正弦函数为奇对称,当时,对称中心为,需要满足关于偶对称,即要求: , (1-6)对于第二类线性相位,即时,通过运算得到公式: (1-7)函数为偶对称,当时,对称中心也为。若要使上式成立,则要使关于奇对称,即要求:, (1-8)从上述分析看来,线性相位FIR数字滤波器时域约束条件是指满足线性相位时对约束条件,对于第一类线性相位,冲激响应满足(1-6)式;对于第二类线性相位,冲激响应h(n)满足(1-8)式。二 系统设计FIR数字滤波器设计方法主要有窗函数设计法、频率采样设计法以及等波纹逼近设计法三种,其中窗函数设计法是最常用
9、,其次是频率采样法,但这两种方法在设计中还会存在一些不足之处,所以需要优化设计方法,而等波纹逼近法很好弥补了窗函数法与频率采样法不足6。2.1 FIR数字滤波器窗函数设计方法 窗函数法设计思路窗函数设计法是FIR数字滤波器里最简单一种设计法,又叫傅里叶级数法,为了设计简单方便,通常选择所希望逼近滤波器频率响应函数为具有片段常数特性理想滤波器,寻找一组,确定其频率响应,然后用来逼近1。窗函数法设计FIR滤波器是在时域中进行,那么可以通过傅里叶反变换得到得到频率响应,即: (2-1)在实际中,一般是处于逐段恒定,在边界频率处有不连续点,因而单位脉冲响应是无限长非因果序列,不能直接作为FIR数字滤波
10、器单位脉冲响应,因此需要对进行截断,转换为有限长一段因果序列,也就是用一个有限长度窗函数序列来截取,即,并将非因果序列转变为一个因果序列。截取长度与加权窗函数形状都直接影响到逼近精度。以截止频率为,相位为零理想低通滤波器为例,其频率特性为: (2-2)通过傅里叶反变换得到对应为: (2-3)此时是一个无限长非因果序列,我们需要对其进行截断,变成一个有限长因果序列。可以先把向右平移个点,得到为: (2-4)相应传输函数为: (2-5)然后对截取从0到N个点,N为窗函数长度,所得结果表示为: (2-6)表示窗函数,一般用下标来表示窗函数类型,矩形窗记为。 常见窗函数介绍常见窗函数有矩形窗(Rect
11、angle Window)、三角形窗(Bartlerr Window)、汉宁(Hanning)窗升余弦窗、哈明(Hamming)窗改进升余弦窗、布莱克曼(Blackman)窗、凯塞贝塞尔窗(Kaiser-Basel Window)7。矩形窗窗函数为: (2-7)其频谱幅度函数为 (2-8)矩形窗主瓣宽度为,用矩形窗设计FIR数字滤波器过渡带宽度近似为。三角形窗窗函数为: (2-9)其频谱幅度函数为 (2-10)三角窗主瓣宽度为。汉宁窗窗函数为: (2-11)汉宁窗频谱幅度函数为 (2-12)汉宁窗主瓣宽度为,汉宁窗在其两个端点都为零,实际中这两个端点数据是不可用。哈明窗窗函数为 (2-13)其
12、幅度函数为 (2-14)哈明窗是一种改进余弦窗,能量更加集中在主瓣,是一种高效窗函数,主瓣宽度及汉宁窗相同。布莱克曼窗窗函数为 (2-15)其频谱幅度函数为 (2-16)该窗函数位移不同,幅度函数也不同,会使旁瓣进一步抵消,主瓣宽度为。凯塞窗是一种最优窗函数,不同于前面五种窗函数,凯塞窗是一种参数可调窗函数,其函数形式如下: (2-17)其中 (2-18) (2-19)一般取15-25项可以满足精度要求。参数可以控制窗形状。一般越大,主瓣越宽,而旁瓣幅度会随之减小,典型数据在4到9之间。2.1.3 吉布斯效应 用窗函数对进行直接截断,得到有限长序列,并以代替,肯定会引起误差,表现在频域就是通常
13、所说吉布斯(Gibbs)效应。对于一个在有限区间分布信号,其连续频谱在频域上分布往往是无限,而在实际信号处理时,我们通常只能在有限区间内做傅里叶分析,也就是说,我们只能用有限区间来代替理论分析中无限区间,多数情况下,我们总是选择信号低频部分,而舍弃高频部分。而信号高频部分往往是反应信号快速变化特征,如果信号本身是连续,这样做一般不会引起信号显著变化,但实际中信号往往是比较丰富,特别是信号本身存在剧烈变化,这样做必定会引起一些误差。该误差引起过渡带加宽以及通带与阻带内波动。为了减小吉布斯效应带来影响,需要调整窗口长度来控制过渡带宽度,但要减小带内波动以及增大阻带衰减,还需要从窗函数形状上寻找解决
14、方法8。为了减少序列因截断而产生Gibbs效应,窗函数在设计时需要注意:(1)频率特性主瓣要尽可能窄,并且尽量把能量都集在主瓣内。(2)窗函数频率特性旁瓣趋于过程中,其能量迅速减小为零。虽然窗函数设计法设计思路简单,但是它边界频率不容易控制,而且窗函数还有吉布斯效应,需要选择不同窗函数来减小吉布斯效应对结果影响,但无论哪种窗函数,都无法很好解决这一问题,所以我们需要通过其他设计方法来进行滤波,便于满足实际工程中不同要求。2.2 FIR数字滤波器频率采样设计法窗函数设计法是从时域出发来设计FIR数字滤波器,而频率采样法是从频域出发设计FIR数字滤波器。与窗函数设计法相同,频率采样法也需要预先构造
15、一个希望逼近滤波器频率响应函数,对其加以等间隔采样后,作为FIR数字滤波器频率响应。对在到之间等间隔采样点,得到频率采样值: (2-20)再对进行点IDFT,得到: (2-21)将作为所涉及FIR数字滤波器单位脉冲响应,其系统函数为为 (2-22)由于滤波器频率响应是理想,即有间断点,那么其单位冲激响应是无限长。这样,由于时域混叠,引起所设计h(n)与有偏差。因此,采样点处及相等,逼近误差为0,而在采样点之间,由有限项之与形成。其误差与特性平滑程度有关,特性愈平滑误差愈小;特性曲线间断点处,误差越大。误差表现形式为间断点用倾斜线取代,且间断点附近形成振荡特性,使阻带衰减减小,往往不能满足实际工
16、程中技术要求。当然,增大N值,可以减小逼近误差,但间断点附近误差仍然最大,且N太大会增加滤波器级数及成本。提高阻带衰减最有效方法是在频响间断点附近区间内插一个或几个过渡采样点,使不连续点变成缓慢过渡。过渡带采样点个数及阻带最小衰减关系以及使阻带最小衰减最大化每个过渡带采样值求解都要用优化算法解决。其基本思路是将过渡带采样值设为一个自由量,用一种优化算法改变它们,最终使阻带最小衰减最大。2.3 FIR数字滤波器等波纹逼近设计法窗函数设计法与频率采样设计法虽然设计方法简单,但都存在滤波器边缘频率不易精确控制缺点,且这两种设计方法设计出来滤波器通带与阻带波动幅度都是相等,两种设计方法都不能分别控制通
17、带与阻带波动幅度,而现实工程中往往对二者都有不同要求,需要分别进行控制。等波纹逼近法是一种优化设计方法,它克服了窗函数设计法与频率采样法缺陷,是最大误差最小化设计方法,并在整个逼近频段上均匀分布。设为希望逼近幅度特性函数,且要求设计线性相位FIR数字滤波器时,必须满足线性相位约束条件。用表示实际设计幅度特性函数,定义加权误差函数为 (2-23)式中,被称为误差加权函数,是由设计者定义,用来控制不同频段逼近精度。经过推导可把统一标示为: (2-24)式中,是系数不同余弦组合式,记;是不同常数,在设计FIR滤波器时存在四种线性相位,当且奇对称时,N为奇数,为1;N为偶数时,为;而当偶对称时,不管N
18、为奇数还是偶数,都取。利用等波纹逼近法设计FIR滤波器,其误差均匀分布在频带中,可以得优良滤波特性,它在同样过渡带较窄情况下,通带最稳定,阻带有最大化最小衰减。三 详细设计3.1 程序设计流程采用窗函数法设计FIR数字低通滤波器程序设计流程如图3-1所示。图3-1 窗函数法设计FIR滤波器流程3.2 Matlab简介Matlab是MATrixLABoratory(“矩阵实验室”)缩写,是由美国MathWorks公司开发集数值计算、符号计算与图形可视化三大基本功能于一体,功能强大、操作简单语言。Matlab是国际公认优秀数学应用软件之一。它集中了日常数学处理中各种功能,包括高校数值计算、矩阵运算
19、、信号处理与图像生成等功能。在Matlab环境下,用户可以进行程序设计、数值计算、图形绘制、输入输出、文件管理等各项操作。除此之外,Matlab易于扩充。除内部函数外,所有Matlab核心文件与工具箱文件都是可读可改源文件,用户可修改源文件与加入自己文件,它们可以及库函数一样被调用。Matlab是一种矩形运算为基础交互式程序语言,着重针对科学计算、工程计算与绘画需求,及其他机器语言相比,其特点是简单与智能化,适应科技专业人员思维方式与书写习惯,使得编程与调试效率大大提高。Matlab由一系列工具组成。这些工具方便用户使用Matlab函数与文件,其中许多工具采用是图形用户界面。包括Matlab桌
20、面、命令窗口、历史命令窗口、编辑器与调试器、路径搜索与用于用户浏览帮助、工作空间、文件浏览器。随着Matlab商业化以及软件本身不断升级,Matlab用户界面也越来越精致,更加接近Windows标准界面,人机交互性更强,操作更简单。Matlab自产生之日起就具有方便数据可视化功能,将向量与矩阵用图形表现出来,并且可以对图形进行标注与打印。Matlab具有功能强大工具箱,工具箱可分为两类:功能性工具箱与学科性工具箱。功能性工具箱主要用来扩充其符号计算功能、图示建模仿真功能、文字处理功能以及及硬件实时交互功能。而学科性工具箱是专业性比较强,如优化工具箱、统计工具箱、控制工具箱、小波工具箱、图象处理
21、工具箱、通信工具箱等。 fir1函数介绍Matlab信号处理工具箱提供了基于加窗线性相位FIR滤波器设计函数fir1。fir1调用格式为:b=fir1(n, Wc,ftype,window)函数参数说明如下:2. Wc为滤波器归一化截止频率,它是一个大于0小于1一个数。3.ftype表示所设计滤波器类型,如果ftype=high,则表示高通滤波器,如果ftype=stop,则表示带阻滤波器,如果此时没有参数,就表示低通滤波器。4.window表示是指定窗函数,如矩形窗为rectwin(n),三角窗为bartlett(n),如果缺省window参数,则fir1默认为是哈明窗hamming(n)。
22、基于fir1函数窗函数法FIR滤波器设计下面利用fir1函数哈明窗法设计数字低通滤波器。利用fir1函数进行设计,这种设计方法只需要给出滤波器阶数、截止频率、窗函数等参数,Matlab即可自行完成设计,并可通过freqz函数查看滤波器幅频响应与相频响应。,通带衰减小于等于0.1dB,阻带衰减大于等于50dB。程序参见附录。图3-1与图3-2为哈明窗法设计数字低通滤波器幅频、相频曲线图以及滤波器增益响应图。图3-1哈明窗FIR滤波器幅频与相频特性曲线图3-2哈明窗FIR低通滤波器增益响应曲线,通带衰减为0.019dB,阻带衰减为53dB。图3-3与图3-4分别为信号滤波前时域与频域图以及信号滤波
23、后时域与频域图。 图3-3 信号滤波前时域与频域图图3-4 信号滤波后时域与频域图从图3-3与图3-4图像中可以看到:输入信号是由两个不同频率正弦信号叠加而成信号频域图中位于滤波器通带内频率分量保留了下来位于滤波器阻带内频率分量被滤除,符合设计要求。总 结本次设计主要分析了FIR数字滤波器基本理论,讨论了FIR数字滤波器线性相位种类及其约束条件,分析了窗函数设计法、频率采样设计法、等波纹逼近法三种不同设计方法,并借助Matlab软件对窗函数法进行了分析及仿真。窗函数设计法对信号加窗之后会使不连续点处边沿加宽形成过渡带,其宽度(两肩峰之间宽度)等于窗函数频率响应主瓣宽度。在处出现肩峰值,两侧形成
24、起伏振荡,振荡幅度与多少取决于旁瓣幅度与多少,改变N只能改变窗谱主瓣宽度,但不能改变主瓣及旁瓣相对比例,其相对比例由窗函数形状决定。窗函数设计法设计FIR数字滤波器是傅里叶变换典型运用,而频率采样法设计指导思想是通过频域采样点实现,而等波纹逼近法可以使误差均匀分布,对于相同技术指标,这种逼近法所需要阶数最少,而对于相同滤波器阶数,等波纹逼近法误差可以达到最小。设计一个FIR滤波器不管是哪种方法,需要完成大量计算与图形绘制工作,而且从上面设计过程中可以看到,设计中只用到两个技术指标,也就是通过截止频率与阻带截止频率,其他指标:带内允许最大衰减,带外允许最小衰减指标,无法表达在设计过程中来,所以,
25、设计结果通常不可能一次计算得到,往往需要反复多次处理,对比才能最终得到符合各项技术指标设计结果,利用Matlab软件,可以减少设计中工作量,从而使FIR滤波器设计变得简单快捷。参考文献1丁玉美, 高西全.数字信号处理 M.第三版. 西安:西安电子科技大学出版社, 2008, 6. 2余成波,陶红艳,杨菁等.数字信号处理及Matlab实现M.第二版.清华大学出版社,2008,1.3郑国强,傅江涛,彭勃等.数字信号处理理论及实践M.第一版.西安电子科技大学出版社,2009,8.4 闫胜利. FIR滤波器及设计原理J. 长春工程学院学报(自然科学版), 2003, 6, 4.5杨永昌,李晨辉,王凯.
26、 FIR数字滤波器设计方法J. 桂林航天工业高等专科院校学报, 2006, 11.6李寿柏, 胡业林. MATLAB在FIR滤波器设计中应用J. 机电工程技术, 2007, 12.7 朱敏. MATLAB数字信号处理工具箱开发与应用J. 信息及电脑, 2010, 2.8 张小虹, 黄忠虎, 邱正伦等. 数字信号处理M.北京:机械工业出版社, 2008, 9.附 录f1=100;f2=200;%待滤波正弦信号频率fs=2000;%采样频率m=(0.3*f1)/(fs/2);%定义过度带宽M=round(8/m);%定义窗函数长度N=M-1;%定义滤波器阶数b=fir1(N,0.5*f2/(fs/
27、2);%使用fir1函数设计滤波器%输入参数分别为滤波器阶数与截止频率figure(1)freqz(b);%输出滤波器幅频与相频曲线figure(2)h,f=freqz(b,1,512);%滤波器幅频特性图%H,W=freqz(B,A,N)当N是一个整数时函数返回N点频率向量与幅频响应向量plot(f*fs/(2*pi),20*log10(abs(h)%参数分别是频率与幅值xlabel(频率/赫兹);ylabel(增益/分贝);title(滤波器响应增益);figure(3)subplot(211)t=0:1/fs:0.5;%定义时间范围与步长s=sin(2*pi*f1*t)+sin(2*pi
28、*f2*t);%滤波前信号plot(t,s);%滤波前信号图像xlabel(时间/秒);ylabel(幅度);title(信号滤波前时域图);subplot(212)Fs=fft(s,512);%将信号变换到频域AFs=abs(Fs);%信号频域图幅值f=(0:255)*fs/512;%频率采样plot(f,AFs(1:256);%滤波前信号频域图xlabel(频率/赫兹);ylabel(幅度);title(信号滤波前频域图);figure(4)sf=filter(b,1,s);%使用filter函数对信号进行滤波subplot(211)plot(t,sf);%滤波后信号图像xlabel(时间
29、/秒);ylabel(幅度);title(信号滤波后时域图);axis(0.2 0.5 -2 2);%限定图像坐标范围subplot(212)Fsf=fft(sf,512); AFsf=abs(Fsf); %滤波后信号频域图及信号频域图幅值f=(0:255)*fs/512;%频率采样plot(f,AFsf(1:256);%滤波后信号频域图xlabel(频率/赫兹);ylabel(幅度);title(信号滤波后频域图);致 谢这次课程设计终于顺利完成了,在设计中遇到了很多专业知识问题,最后在教师辛勤指导下,终于迎刃而解,我们学也到很多实用知识,学到了很多课内学不到东西,比如独立思考解决问题,出现差错随机应变,与及人合作共同提高,都受益非浅。经过一个星期努力与完善,这次课程设计终于顺利完成了。 在这里首先感谢我课程设计指导教师何继爱及朱宁宁教师,在这段时间一直给我支持及鼓励。认真负责监督我们课程设计进度,耐心指导我们使我们能够按时完成任务。还有感谢学校为我们提供良好实验环境以及充足实验设备,为我们设计与调试提供了很大方便。在这段时间学到了很多,虽然由于自身不足没有能够为系统提出更好解决方案。但这对我来说绝对是一个非常宝贵历练。从中我切身体会到了理论与现实差距,只有真正动手去做才能发现问题。最后,感谢所有对我提供过帮助人。谢谢。