基于MATLAB的数字滤波器设计

更新时间:2024-01-11 02:39:01 阅读量: 教育文库 文档下载

说明:文章内容仅供预览,部分内容可能不全。下载后的文档,内容与下面显示的完全一致。下载之前请确认下面内容是否您想要的,是否完整无缺。

淮北煤炭师范学院

2009届学士学位论文

基于MATLAB的数字滤波器设计

学院、 专业 物理与电子信息学院

电子信息科学与技术 研 究 方 向 基于MATLAB的数字滤波器设计 学 生 姓 名 耿 博 学 号 200513432024 指导教师姓名 邹 锋 指导教师职称 讲 师

2009 年 4 月 18

淮北煤炭师范学院2009届学士学位论文 基于MATLAB的数字滤波器设计

基于MATLAB的数字滤波器设计

耿 博

(淮北煤炭师范学院物理与电子信息学院 235000)

摘 要 随着信息时代和数字世界的到来,数字信号处理已成为今一门极其重要的学科和技术领域。数字信号处理在通信、语音、图像、自动控制、雷达、军事、航空航天、医疗和家用电器等众多领域得到了广泛的应用。在数字信号处理应用中,数字滤波器十分重要并已获得广泛应用。

数字滤波是数字信号处理的重要内容,数字滤波器可分为IIR和FIR两大类。对于IIR数字滤波器的设计,需要借助模拟原型滤波器,再将模拟滤波器转化为数字滤波器,文中采用的设计方法是脉冲响应不变法、双向性变换法和完全函数设计法;对于FIR数字滤波器的设计,可以根据所给定的频率特性直接设计,文中采用的设计方法是窗函数法。本文根据IIR滤波器和FIR滤波器的特点,在MATLAB坏境下分别用双线性变换法设计IIR和用窗函数设计FIR数字滤波器,并对采集的语音信号进行分析,最后给出了IIR和FIR对语音滤波的效果。 关键词 数字滤波器;IIR ;FIR ;MATLAB

II

淮北煤炭师范学院2009届学士学位论文 基于MATLAB的数字滤波器设计

The Design of Digital Filter based on MATLAB

Geng Bo

School of Physics and Electronics Information, Huaibei Coal Industry Teachers? College, 235000

ABSTRACT Along with the information age and the digital world arrival, the digital signal processing has become a now extremely important discipline and the area of technology.The digital signal processing in the correspondence, the multitudinous domains the pronunciation such as the image, the automatic control, the radar, the military, the aerospace, the medical service and the domestic electric appliances and so on have obtained the widespread application.In the digital signal processing application, the digital filter are extremely important and have obtained the widespread application.

The digital filter are the digital signal processing important content, the digital filter may divide into IIR and the FIR two main kinds. As for the IIR digital filter design, we need the help of analog prototype filter, and then transform analog filter into digital filter. In the paper we use the design of the pulse response invariable method, the bilinear method and full function design; as for the FIR filter, we can design it directly based on the giving frequency, in the paper it uses the design of the window function.This article according to the IIR filter and the FIR filter characteristic, uses the bilinearity method of transformation under the MATLAB bad boundary to design IIR and to design the FIR numeral filter separately with the window box number, and carries on the analysis to the gathering pronunciation signal, and finally gives IIR and FIR to the pronunciation filter effect. Keywords Digtial Filter;IIR;FIR;MATLAB

III

淮北煤炭师范学院2009届学士学位论文 基于MATLAB的数字滤波器设计

目 次

1 引言 ......................................................................................... 1 2 数字滤波器及MATLAB语言概述 .......................................... 2

2.1 数字滤波器的定义和分类 ........................................ 2 2.2 常用滤波器的性能指标 ............................................ 3 2.3 MATLAB概述 ................................................................ 6

3 IIR滤波器设计 ...................................................................... 8

3.1 双线性变换法设计IIR数字滤波器 ........................ 8 3.2 脉冲响应不变法 ....................................................... 12 3.3 完全设计函数法 ....................................................... 15 3.4 语音滤波实例 ........................................................... 16

4 FIR滤波器设计 .................................................................... 21

4.1 窗函数法 ..................................................................... 21 4.2 FIR滤波器滤波实例 ............................................... 25

5 总结 ....................................................................................... 29 参考文献 ..................................................................................... 30 致谢 ............................................................................................. 31

IV

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

1 引言

数字滤波在通信、图像编码、语音编码、雷达等许多领域中有着十分广泛的应用。目前,数字信号滤波器的设计在图像处理、数据压缩等方面的应用取得了令人瞩目的进展和成就。近年来迅速发展起来。

MATLAB是美国MathWorks公司推出的一套用于工程计算的可视化高性能语言与软件环境。MATLAB为数字滤波的研究和应用提供了一个直观、高效、便捷的利器。它以矩阵运算为基础,把计算、可视化、程序设计融合到了一个交互式的工作环境中。MATLAB推出的工具箱使各个领域的研究人员可以直观方便地进行科学研究、工程应用,其中的信号处理(signalproeessing)、图像处理(imageproeessing)、小波(wavelet)等工具箱为数字滤波研究的蓬勃发展提供了有力的工具。

数字滤波器与模拟滤波器相比,具有精度高、稳定、体积小、重量轻、灵活、不要求阻抗匹配以及能实现模拟滤波器无法进行的特殊滤波等优点[1]。本文主要介绍有限冲激响应数字滤波器(FIR)和无限冲激响应数字滤波器(IIR)的设计原理、方法、步骤以及在MATLAB中的实现,并以实例形式列出设计程序和仿真结果。

1

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

2 数字滤波器及MATLAB语言概述

2.1 数字滤波器的定义和分类

数字滤波器是指完成信号滤波处理功能的,用有限精度算法实现的离散时间线性非时变系统,其输入是一组数字量,其输出是经过变换的另一组数字量。因此,数字滤波器本身既可以是用数字硬件装配成的一台完成给定运算的专用的数字计算机,也可以将所需要的运算编成程序,让通用计算机来执行。数字滤波器,输入输出均为数字信号,通过一定的运算关系,改变输入信号中所含频率成分的相对比例,或则滤除某些频率成分的器件[2]。对于数字滤波器而言,若系统函数为H(z),其脉冲响应为h(n),输入时间序列为x(n),则它们在时域内的关系式如下:

y(n)=h(n)﹡x(n) (2-1)

在Z域内,输入和输出存在如下关系:

Y(z)= H(z)X(z) (2-2)

式中,X(z)、Y(z)分别为x(n)和y(n)的Z变换。

在频域内,输入和输出则存在如下关系:

Y(j?)=H(j?)X(j?) (2-3)

式中,H(j?)是数字滤波器的频率特性;X(j?)Y(j?)分别为x(n)和y(n)的频谱,而

?为数字角频率。

数字滤波器可以有很多种分类方法,但总体上可分为两大类。一类称为经典滤波器,即一般的滤波器,其特点是输入信号中的有用成分和希望滤除的成分占用不同的频带,通过合适的选频滤波器可以实现滤波[3]。例如,若输入信号中有干扰,信号和干扰的频带互不重叠,则可滤出信号中的干扰得到纯信号。但是,如果输入信号中信号和干扰的频带相重叠,则干扰就不能被有效的滤出。另一类称为现代滤波器,如维纳滤波器、卡尔曼滤波器等,其输入信号中有用信号和希望滤除的成分频带重叠。对于经典滤波器,从频域上也可以分为低通、高通、带通和带阻滤波器。从时域特性上看,数字滤波器还可以分为有限冲激响应数字滤波器(FIR)和无限冲激响应数字滤波器(IIR)。

对于有限冲激响应数字滤波器(FIR),其输出y(n)只取决于有限个过去和现在的输入,x(n),x(n-1),…,x(n-m),滤波器的输入输出关系可表示为

y(n)=?brx(n?r) (2-4)

r?0M2

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

对于无限冲激响应数字滤波器(IIR),它的输出不仅取决于过去和现在的输入,而且还取决于过去的输出,其差分方程为

y(n)+?aky(n?k)=?brx(n?r) (2-5)

k?1r?0NM该差分方程的单位冲激响应是无限延续的。

2.2 常用滤波器的性能指标

j?H(e)来说明,常用的性能指标主要有以下滤波器性能一般用系统频率特性

三个参数:

1. 幅度平方函数

H(ej?)?H(ej?)?H*(ej?)?H(ej?)H(e?j?) ?H(z)H*(z)z?ej?2 (2-6)

该性能指标主要用来说明系统的幅频特性。

2.相位函数

H(ej?)?Re[H(ej?)]?jIm[H(ej?)]?H(ej?)ej?(ej?) (2-7)

其中:

?Im[H(ej?)]??(e)?arctg??j?Re[H(e)]?? (2-8)

j?该指标主要用来说明系统的相位特性。

3.群延时

d[?(ej?)]?(?)??d? (2-9)

定义为相位对角频率导数的负值,说明了滤波器对不同的频率成分的平均延时。当要求在通带内的群延迟是常数时,滤波器相位响应特性应该是线性的。

实际设计中所能得到的滤波器的频率特性与理想滤波器的频率特性之间存在着一些显著的差别,现以低通滤波器的频率特性为例进行说明。

1.理想滤波器的特性:

设滤波器输入信号为x(t),信号中混入噪音u(t),它们有不同的频率成分。滤

3

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

波器的单位脉冲响应为h(t)。则理想滤波器输出为:

y(t)?[x(t)?u(t)]?h(t)?K?x(t??) (2-10)

即噪音信号被滤除u(t)?h(t)?0,而信号无失真只有延时和线性放大。对上式作傅里叶变换得:

Y(j?)?X(j?)?H(j?)?U(j?)?H(j?)?Ke?j??X(j?) (2-11)

假定噪音信号被滤除,即:

U(j?)?H(j?)?0 (2-12)

整理得:

H(j?)?Y(j?)?Ke?j??X(j?) (2-13)

h(t)

H(j?)K0

?c?

t

图1 理想低通滤波器频率特性 图2 理性滤波器的单位脉冲响应(??0)

假定信号频率成分为:低通滤波器特性是:

???c,噪音频率成分为

???c。则完成滤波的理想

?j?? |?|??cY(j?)?K?eH(j?)???X(j?)?0 |?|??c (2-14)

即:

?K |?|??c|H(j?)|???0 |?|??c (2-15)

arg(H(j?))???? (2-16)

系统的单位脉冲响应为:

1h(t)?2???c??cKe?j???ej?td??Ksin[(t??)?c]?(t??) (2-17)

理性低通滤波器的频率特性如图1所示,单位脉冲响应的波形如图2所示。 理想滤波器具有非因果、无限长的单位脉冲响应和不连续的频率特性,要用

4

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

稳定的线性时不变(LTI)系统来实现这样的特性是不可能的[4]。工程上是用脉冲响应为有限长的、因果的、稳定的线性时不变系统或具有连续频率特性的线性时不变系统来逼近理想特性。在满足一定的误差要求的情况下来实现理想滤波特性。因此实际的滤波器的频率特性如图3所示。

其中:

1??11??1H(ej?)过渡带 ?2通 带 阻带 0?c?s?图3实际滤波器的频率特性性

?c——截止频率 ?s——阻带起始频率

?s??c——过渡带宽 在通带内幅度响应以

??1的误差接近于1,即:

1??1?H(ej?)?1??1???c (2-18)

?s为阻带起始频率,在阻带内幅度响应以小于?2的误差接近于零,即:

H(ej?)??2?s???? (2-19)

为了使逼近理想低通滤波器的方法成为可能,还必须提供一带宽为不为零的过渡带。在这个频带内,幅度响应从通带平滑的下落到阻带。

?s??c的

5

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

2.3 MATLAB概述

MATLAB是美国MathWorks公司开发的一种功能极其强大的高技术计算语言和内容极其丰富的软件库,集数值计算、矩阵运算和信号处理与显示于一身。该软件最初是由美国教授Cleve Moler创立的。1980年前后,他在教线性代数课程时,发现用其他高级语言编程时极不方便,便构思开发了MATLAB,即矩阵实验室(Matrix Laboratory)。该软件利用了当时代表数值线性代数领域最高水平的EISPACK和LINPACK两大软件包,并且利用Fortran语言编写了最初的一套交互式软件系统,MATLAB的最初版本便由此产生了[5]。

最初的MATLAB由于语言单一,只能进行矩阵的运算,绘图也只能用原始的描点法,内部函数只有几十个,因此功能十分简单。1984年该公司推出了第一个MATLAB的商业版,并用C语言作出了全部改写。现在的MATLAB程序是MathWorks公司用C语言开发的,第一版由steve Bangert主持开发编译解释程序,Steve Kleiman完成图形功能的设计,John Little和Cleve Moler主持开发了各类数学分分析的子模块,撰写用户指南和大部分的M文件。接着又添加了丰富的图形图像处理、多媒体功能、符号运算和与其它流行软件的接口功能,使MATLAB的功能越来越强大。

MTALAB系统主要由以下五个部分组成[6]:

(1)MATALB语言体系。 MATLAB是高层次的矩阵/数组语言,具有条件控制、函数调用、数据结构、输入输出、面向对象等程序语言特性。利用它既可以进行小规模端程,完成算法设计和算法实验的基本任务,也可以进行大规模编程,开发复杂的应用程序。

(2)MATLAB工作环境 。这是对MATLAB提供给用户使用的管理功能的总称。包括管理工作空间中的变量据输入输出的方式和方法,以及开发、调试、管理M文件的各种工具。

(3)图形句相系统 。这是MATLAB图形系统的基础,包括完成2D和3D数据图示、图像处理、动画生成、图形显示等功能的高层MATLAB命令,也包括用户对图形图像等对象进行特性控制的低层MATLAB命令,以及开发GUI应用程序的各种工具。

(4)MATLAB数学函数库。这是对MATLAB使用的各种数学算法的总称。包括各种初等函数的算法,也包括矩阵运算、矩阵分析等高层次数学算法。

(5)MATLAB应用程序接口(API)。这是MATLAB为用户提供的一个函数库,使得用户能够在MATLAB环境中使用C程序或FORTRAN程序,包括从MATLAB中调

6

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

用于程序(动态链接),读写MAT文件的功能[7]。

除此之外,MATLAB系统还具有如下特点: (1)具有易学易用的语言体系; (2)具有交互式的工作环境; (3)具有多层面的图像处理系统; (4)具有丰富高效的MATLAB工具箱; (5)具有便利的程序接口(API); (6)应用领域广泛;

(7)嵌入了面向对象编程语言。

7

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

3 IIR滤波器设计

3.1 双线性变换法设计IIR数字滤波器

从S平面到Z平面是多值的映射关系会造成频率响应的混叠失真。为了克服这一缺点,可以采用非线性频率压缩方法,将整个频率轴上的频率范围压缩到-π/T~π/T之间,再用z=esT转换到Z平面上。也就是说,第一步先将整个S平面压缩映射到S1平面的-π/T~π/T一条横带里;第二步再通过标准变换关系z=es1T将此横带变换到整个Z平面上去。这样就使S平面与Z平面建立了一一对应的单值关系,消除了多值变换性,也就消除了频谱混叠现象,映射关系如图4所示。

图4 双线性变换的映射关系

为了将S平面的整个虚轴jΩ压缩到S1平面jΩ1轴上的-π/T到π/T段上,可以通过以下的正切变换实现

(3-1)

式中,T仍是采样间隔。

当Ω1由-π/T经过0变化到π/T时,Ω由-∞经过0变化到+∞,也即映射了整个jΩ轴。将式(3-1)写成

(3-2)

将此关系解析延拓到整个S平面和S1平面,令jΩ=s,jΩ1=s1,则得

(3-3)

再将S1平面通过以下标准变换关系映射到Z平面z=es1T 从而得到S平面和Z平面的单值映射关系为:

8

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

(3-4)

(3-5)

式(3-4)与式(3-5)是S平面与Z平面之间的单值映射关系,这种变换都是两个线性函数之比,因此称为双线性变换

式(3-1)与式(34)的双线性变换符合映射变换应满足的两点要求。 首先,把z=ejω,可得

(3-6)

即S平面的虚轴映射到Z平面的单位圆。

其次,将s=σ+jΩ代入式(3-6),得

(3-7)

因此

(3-8)

由此看出,当σ<0时,|z|<1;当σ>0时,|z|>1。也就是说,S平面的左半平面映射到Z平面的单位圆内,S平面的右半平面映射到Z平面的单位圆外,S平面的虚轴映射到Z平面的单位圆上。因此,稳定的模拟滤波器经双线性变换后所得的数字滤波器也一定是稳定的。

双线性变换法优缺点:

双线性变换法,其主要的优点是避免了频率响应的混叠现象。这是因为这里S平面与Z平面是单值的一一对应关系。S平面整个jΩ轴单值地对应于Z平面单位圆一周,即频率轴是单值变换关系。这个关系如式(3-6)所示,重写如下:

(3-9)

上式表明,S平面上Ω与Z平面的ω成非线性的正切关系,如图5所示。

9

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

由图5看出,在零频率附近,模拟角频率Ω与数字频率ω之间的变换关系接近于线性关系;但当Ω进一步增加时,ω增长得越来越慢,最后当Ω→∞时,ω终止在折叠频率ω=π处,因而双线性变换就不会出现由于高频部分超过折叠频率而混淆到低频部分去的现象,从而消除了频率混叠现象。

图5 双线性变换法的频率变换关系

但是双线性变换的这个特点是靠频率的严重非线性关系而得到的,如式(3-6)及图5所示。由于这种频率之间的非线性变换关系,就产生了新的问题。首先,一个线性相位的模拟滤波器经双线性变换后得到非线性相位的数字滤波器,不再保持原有的线性相位了;其次,这种非线性关系要求模拟滤波器的幅频响应必须是分段常数型的,即某一频率段的幅频响应近似等于某一常数(这正是一般典型的低通、高通、带通、带阻型滤波器的响应特性),不然变换所产生的数字滤波器幅频响应相对于原模拟滤波器的幅频响应会有畸变,如图8所示。

图6 双线性变换法幅度和相位特性的非线性映射

对于分段常数的滤波器,双线性变换后,仍得到幅频特性为分段常数的滤波器,但是各个分段边缘的临界频率点产生了畸变,这种频率的畸变,可以通过频率的预畸来加以校正。也就是将临界模拟频率事先加以畸变,然后经变换后正好映射到所需要的数字频率上。

在MATLAB中,双线性变换法的调用函数是bilinear。其调用格式为: a.[zd,pd,kd]= bilinear(z,p,k,fs)

10

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

b.[zd,pd,kd]= bilinear(z,p,k,fs,fp) c. [numd,dend]=bilinear(num,den,fs) d. [numd,dend]=bilinear(num,den,fs,fp) e. [Aa,Bb,Cc,Dd]=bilinear(A,B,C,D,fs) f. [Aa,Bb,Cc,Dd]=bilinear(A,B,C,D,fs,fp)

[zd,pd,kd]=bilinear(z,p,k,fs)是把模拟滤波器的零极点模型转换为数字滤波器的零极点模型,fs为采样频率,z,p,k分别为滤波器的零点、极点和增益;

[numd,dend]=bilinear(num,den,fs)是把模拟滤波器的传递函数模型转换为数字滤波器的传递模型;

[Aa,Bb,Cc,Dd]=bilinear(A,B,C,D,fs)是把模拟滤波器的状态方程模型转换为数字滤波器状态方程模型。

例如,用双线性变换法设计一个巴特沃思数字低通滤波器,技术指标如下:通带截止频率?p=2??4k rad/s,阻带截止频率?s=2??8k rad/s,通带波纹系数

Rp=0.3dB, 阻带波纹系数Rs=50dB,采样频率fs=20000Hz。

程序如下:

[N,Wn]=buttord(wp,ws,Rp,Rs,?s?) %估计滤波器最小阶数 [z,p,k]=buttap(N); [Bap,Aap]=zp2tf(z,p,k); [b,a]=lp2lp(Bap,Aap,Wn); [bz,az]=bilinear(b,a,fs) freqz(bz,az,Nn,fs)

程序在MATLAB环境下的运行及结果如图7所示: 结果如下: N =11

Wn =1.4892e+004 bz =

Columns 1 through 6

0.0110 0.1211 0.6055 1.8166 3.6333 5.0866 Columns 7 through 12

5.0866 3.6333 1.8166 0.6055 0.1211 0.0110 az =

Columns 1 through 6

11

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

1.0000 2.7098 4.6379 5.2252 4.3685 2.7207 Columns 7 through 12

1.2885 0.4561 0.1181 0.0211 0.0023 0.0001

0Magnitude (dB)-100-200-300050010001500Frequency (Hz)200025000Phase (degrees)-200-400-600-800-1000050010001500Frequency (Hz)20002500

图7巴特沃思数字低通滤波器幅频-相频特性

3.2 脉冲响应不变法

顾名思义,脉冲响应不变法就是要求数字滤波器的脉冲响应序列h(n)与模拟滤波器的脉冲响应ha(t)的采样值相等,即

h(n)=ha(t)的关系,我们知道

H(z)

z?eSTt?nT=ha(nT) (3-10)

式中,T为采样周期。根据模拟信号的拉普拉斯变换与离散序列的Z变换之间

=

1T?Ha(S?jk?ks) (3-11)

此式表明,ha(t)的拉普拉斯变换在s平面上沿虚轴,按照周期?s=2?/T延拓后,按式z=eST,进行Z变换,就可以将Ha(s)映射为H(z)。事实上,用脉冲响应不变法设计IIR滤波器,只适合于Ha(s)有单阶极点,且分母多项式的阶次高于分子多项式阶次的情况[5]。将Ha(s)用部分分式表示:

Ha(s)=LT[ha(t)]=?i?1NAis?si (3-12)

式中,LT[·]代表拉普拉斯变换,si为的单阶极点。将Ha(s)进行拉普拉斯反变换,即可得到

ha(t)=?AieSitu(t) (3-13)

i?1N12

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

N式中,u(t)是单位阶跃函数。则ha(t)的离散序列h(n)=ha(nT)=?AieSinTu(nT)

i?1对h(n)进行z变换之后,可以得到数字滤波器的系统函数H(z)

H(z)=?h(n)zn?0??n=

?1?ei?1NAiSiTz?1 (3-14)

对比Ha(s)与H(z),我们会发现:s域中Ha(s)的极点是si,映射到z平面之后,其极点变成了eSiT,而系数没有发生变化,仍为Ai。因此,在设计IIR滤波器时,我们只要找出模拟滤波器系统函数Ha(s)的极点和系数Ai,通过脉冲响应不变法,代入H(z)的表达式中,即可求出H(z),实现连续系统的离散化[8]。

但是脉冲响应不变法只适合于设计低通和带通滤波器,而不适合于设计高通和带阻滤波器。因为,如果模拟信号ha(t)的频带不是介于??/T之间,则会在

??/T的奇数倍附近产生频率混叠现象,映射到z平面后,则会在???附近产生频率混叠现象。从而使所设计的数字滤波器不同程度的偏离模拟滤波器在???附近的频率特性,严重时使数字滤波器不满足给定的技术指标。为此,希望设计的滤波器是带限滤波器,如果不是带限的,例如,高通滤波器、带阻滤波器,需要在高通滤波器、带阻滤波器之前加保护滤波器,滤出高于折叠频率?/T以上的频带,以免产生频率混叠现象[9]。但这样会增加系统的成本和复杂性。因此,高通与带阻滤波器不适合用这种方法。

将模拟滤波器转化为数字滤波器,牵涉到一个关键的问题,即寻找一种转换关系,将s平面上的Ha(s)转换成z平面上的H(z)。这里Ha(s)是模拟滤波器的传输函数,H(z)是数字滤波器的系统函数[10]。为了确保转换后的H(z)稳定且满足技术要求,转换关系要满足以下要求:

1)将因果稳定的模拟滤波器转换为数字滤波器后,仍然是因果稳定的。我们知道,当模拟滤波器的传输函数Ha(s)的极点全部位于s平面的左平面时,模拟滤波器才是因果稳定的;对于数字滤波器而言,因果稳定的条件是其传输函数H(z)的极点要全部位于单位圆内。因此,转换关系应是s平面的左半平面映射到z平面的单位圆内。

2)数字滤波器的频率响应与模拟滤波器的频率响应相对应,S平面的虚轴映射为Z平面的单位圆,而响应的频率之间是线性变换关系。

在MATLAB中,脉冲响应不变法的调用函数是impinvar,其调用格式为: a. [bz,az]=impinvar(b,a,fs) b. [bz,az]=impinvar(b,a)

13

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

c. [bz,az]=impinvar(b,a,fs,tol)

该函数的功能是将分子向量为b、分母向量为a的模拟滤波器,转换为分子向量为bz、分母向量为az的数字滤波器。fs为采样频率,单位为Hz,默认值为1Hz。tol指误差容限,表示转换后的离散系统函数是否有重复的极点。

例如,用脉冲响应不变法设计一个契比雪夫型数字低通滤波器,指标要求:通带截止频率?p=1000Hz,阻带截止频率?s=1200Hz,采样频率fs=5000Hz,通带衰减系数Rp=0.3dB,阻带衰减系数Rs=40dB。

程序如下:

wp=1000*2*pi;ws=1200*2*pi;fs=2500;Rp=0.3;Rs=40; [N,Wn]=cheb1ord(wp,ws,Rp,Rs,'s'); %估计滤波器最小阶数 [z,p,k]=cheb1ap(N,Rp); %模拟滤波器函数引用 [A,B,C,D]=zp2ss(z,p,k); %返回状态转移矩阵形式 [AT,BT,CT,DT]=lp2lp(A,B,C,D,Wn); %频率转换

[b,a]=ss2tf(AT,BT,CT,DT); %返回传递函数形式 [bz,az]=impinvar(b,a,fs); %调用脉冲相应不变法 [H,W]=freqz(bz,az); %返回频率响应 plot(W*fs/(2*pi),abs(H)); %画图 grid;

xlabel('frequency/Hz'); ylabel('magnitude'); N,Wn

程序在MATLAB环境下的运行及结果如图8所示:

1.41.21magnitude0.80.60.40.200200400600800frequency/Hz100012001400

图8 契比雪夫型数字低通滤波器幅频特性曲线

N =11 Wn =6.2832e+003

14

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

3.3 完全设计函数法

此法是根据设计要求,直接调用函数来设计数字滤波器。所用到的函数有butter、cheby1、cheb2ord、ellipd以及besself等。Butter用来直接设计巴特沃思数字滤波器,cheby1用来直接设计切比雪夫?型滤波器,cheb2ord用来设计切比雪夫

?型滤波器,ellipd用来设计椭圆滤波器,besself用来设计贝塞尔滤波器。

例如,用完全设计函数法设计一个巴特沃思数字低通滤波器,技术指标要求为:wp=1000;ws=1200;Rp=0.3;Rs=40;fs=8000;

程序如下:

wp=1000;ws=1200;Rp=0.3;Rs=40;fs=8000;

[N,Wn]=buttord(wp/(fs/2),ws/(fs/2),Rp,Rs) %估计滤波器最小阶数 [b,a]=butter(N,Wn);

[H,W]=freqz(b,a); %返回频率响应 plot(W*fs/(2*pi),abs(H)); %画图 grid;

xlabel('Frequency/Hz'); ylabel('magnitude');

程序在MATLAB环境下的运行及结果如图9所示:

1.41.21magnitude0.80.60.40.2005001000150020002500Frequency/Hz300035004000

图9巴特沃思数字低通滤波器幅频特性曲线

N =29 Wn =0.2611

15

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

3.4 语音滤波实例

根据语音信号的特点给出有关滤波器的性能指标:1)低通滤波器性能指标,fp=1000Hz,fc=1200 Hz, As=100dB,Ap=1dB;2)高通滤波器性能指标,fc=2800 Hz,fp=3000 Hz As=100dB,Ap=1dB;3)带通滤波器性能指标,fp1=1200 Hz,fp2=3000 Hz,fc1=1000 Hz,fc2=3200 Hz,As=100dB,Ap=1dB。用IIR对其分析,最后,利用MATLAB中的函数freqz画出各滤波器的频率响应。

1.语音信号的采集:

1)利用windows下的录音机(开始—程序—附件—娱乐—录音机,文件—属性—立即转换—8000KHz,8位,单声道),录制一段自己的话音“信号”, 时间控制在1秒左右,然后将音频文件保存“E:\\耿博.wav”

2)在MATLAB软件平台下,利用函数wavread对语音信号进行采样,记住采样频率和采样点数。

[z1,fs,bits]=wavread('E:\\耿博.wav') 运行结果:fs =44100 bits =16

wavread函数调用格式:

y=wavread(file),读取file所规定的wav文件,返回采样值放在向量y中。 [y,fs,nbits]=wavread(file),采样值放在向量y中,fs表示采样频率(Hz),nbits表示采样位数。

y=wavread(file,N),读取前N点的采样值放在向量y中。

y=wavread(file,[N1,N2]),读取从N1点到N2点的采样值放在向量y中。 2.语音信号的频谱分析 ①首先画出语音信号的时域波形 z1=wavread('E:\\耿博.wav'); plot(z1);图像输出如图10

②对语音信号进行频谱分析,在MATLAB中,可以利用函数fft对信号进行快速付立叶变换,得到信号的频谱特性

z1=wavread('E:\\耿博.wav'); y1=z1(1:8192); Y1=fft(y1); n=0:8191;

plot(n,Y1);图像输出如图11

16

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

0.80.66004000.40.20-0.2-2002000-0.4-0.6-0.8-40000.511.522.533.544.5x 104-6000100020003000400050006000700080009000

图10 信号时域波形 图 11 信号频谱分析图

根据语音信号的特点,设计出不同性能的数字滤波器,并用MATLAB进行仿真。设计程序及仿真结果如下: 1.设计IIR低通滤波器滤波 程序如下: clear;close all

[z1,fs,bits]=wavread('E:\\耿博.wav') y1=z1(1:8192); Y1=fft(y1);

fp=1000;fc=1200;As=100;Ap=1; ;Fs=8000; wc=2*fc/Fs;wb=2*fp/Fs; [n,wp]=cheb1ord(wc,wb,Ap,As); [b,a]=cheby1(n,Ap,wp); figure(1); freqz(b,a); x=filter(b,a,z1); X=fft(x,8192); figure(2);

subplot(2,2,1);plot(abs(Y1));axis([0,1000,0,1.0]); title('滤波前信号频谱');

subplot(2,2,2);plot(abs(X));axis([0,4000,0,0.03]); title('滤波后信号频谱'); subplot(2,2,3);plot(z1); title('滤波前信号波形');

17

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

subplot(2,2,4);plot(x); title('滤波后信号波形'); sound(x,fs,bits); 仿真结果如下:

滤波前信号频谱0滤波后信号频谱0.0310.8Magnitude (dB)-2000.60.40.20.02-4000.01-60000.10.20.30.40.50.60.70.8Normalized Frequency (?? rad/sample)0.9100500滤波前信号波形1000001000200030004000滤波后信号波形10.50-0.5010.50-0.500.10.20.30.40.50.60.70.8Normalized Frequency (?? rad/sample)0.91Phase (degrees)-500-1000-1500-2000-1024x 1064-1024x 1064

图 12 IIR数字低通滤波器幅频-相频特性 图 13 滤波前后信号频谱和波形对比

2.设计IIR高通滤波器滤波

程序设计如下: clear;close all

[z1,fs,bits]=wavread('E:\\耿博.wav') y1=z1(1:8192); Y1=fft(y1);

fc=2800 ;fp=3000 ;As=100;Ap=1; Fs=8000; wc=2*fc/Fs;wb=2*fp/Fs; [n,wp]=cheb1ord(wc,wb,Ap,As); [b,a]=cheby1(n,Ap,wp,'high'); figure(1); freqz(b,a); x=filter(b,a,z1); X=fft(x,8192); figure(2);

subplot(2,2,1);plot(abs(Y1));axis([0,1000,0,1.0]); title('滤波前信号频谱');

subplot(2,2,2);plot(abs(X));axis([0,4000,0,0.03]);

18

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

title('滤波后信号频谱'); subplot(2,2,3);plot(z1); title('滤波前信号波形'); subplot(2,2,4);plot(x); title('滤波后信号波形'); sound(x,fs,bits);

图形分析如图14、图15:

滤波前信号频谱滤波后信号频谱0.030Magnitude (dB)10.8-2000.60.40.20.02-4000.01-60000.10.20.30.40.50.60.70.8Normalized Frequency (?? rad/sample)0.9100500滤波前信号波形10000010002000300040000Phase (degrees)滤波后信号波形0.040.020-0.0210.50-0.5-500-1000-1500-200000.10.20.30.40.50.60.70.8Normalized Frequency (?? rad/sample)0.91-1

024x 1064-0.04024x 1064

图14 IIR数字低通滤波器幅频-相频特性 图15 滤波前后信号频谱和波形对比

3.设计IIR带通滤波器滤波

程序设计如下: clear;close all

[z1,fs,bits]=wavread('E:\\耿博.wav') y1=z1(1:8192); Y1=fft(y1);

fp1=1200 ;fp2=3000; fc1=1000 ;fc2=3200 ;As=100;Ap=1; Fs=8000; wc=[2*fc1/Fs,2* fc2/Fs];wb=[2*fp1/Fs,2*fp2/Fs]; [n,wp]=cheb1ord(wc,wb,Ap,As); [b,a]=cheby1(n,Ap,wp,'stop'); figure(1); freqz(b,a); x=filter(b,a,z1); X=fft(x,8192); figure(2);

subplot(2,2,1);plot(abs(Y1));axis([0,1000,0,1.0]);

19

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

title('滤波前信号频谱');

subplot(2,2,2);plot(abs(X));axis([0,4000,0,0.03]); title('滤波后信号频谱'); subplot(2,2,3);plot(z1); title('滤波前信号波形'); subplot(2,2,4);plot(x); title('滤波后信号波形'); sound(x,fs,bits);

图形分析如图16、图17:

0)Bd( -200edutinga-400M-60000.10.20.30.40.50.60.70.80.91Normalized Frequency (?? rad/sample)0)seerg-1000ed( esa-2000hP-300000.10.20.30.40.50.60.70.80.91Normalized Frequency (?? rad/sample)图 16 IIR数字低通滤波器幅频-相频特性 滤波前信号频谱滤波后信号频谱10.030.80.60.020.40.010.2005001000001000200030004000滤波前信号波形滤波后信号波形110.50.500-0.5-0.5-10246-10246x 104x 104

图17 滤波前后信号频谱和波形对比

20

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

4 FIR滤波器设计

FIR滤波器的设计问题在于寻求一系统函数H(z),使其频率响应H(ej?)逼近滤波器要求的理想频率响应Hd(ej?),其对应的单位脉冲响应hd(n)。

4.1 窗函数法

4.1.1用窗函数设计FIR数字滤波器的基本方法

设计思想:从时域从发,设计h(n)逼近理想hd(n)。设理想滤波器Hd(ej?)的单位脉冲响应为hd(n)。以低通线性相位FIR数字滤波器为例。

Hd(e)?hd(n)?j?n????h??d(n)e?jn?1j?jn?H(e)ed?d???2?hd(n)一般是无限长的,且是非因果的,不能直接作为FIR滤波器的单位脉冲

响应。要想得到一个因果的有限长的滤波器h(n),最直接的方法是截断h(n)?hd(n)w(n),即截取为有限长因果序列,并用合适的窗函数进行加权作为FIR滤波器的单位脉冲响应。按照线性相位滤波器的要求,h(n)必须是偶对称的。对称中心必须等于滤波器的延时常数,即

?h(n)?hd(n??)w(n) ????(N?1)/2用矩形窗设计的FIR低通滤波器,所设计滤波器的幅度函数在通带和阻带都呈现出振荡现象,且最大波纹大约为幅度的9%,这个现象称为吉布斯效应。

根据过渡带宽及阻带衰减要求,选择窗函数的类型并估计窗口长度N(或阶数M=N-1),窗函数类型可根据最小阻带衰减As独立选择,因为窗口长度N对最小阻带衰减As没有影响,在确定窗函数类型以后,可根据过渡带宽小于给定指标确定所拟用的窗函数的窗口长度N,设待求滤波器的过渡带宽为Δw,它与窗口长度N近似成反比,窗函数类型确定后,其计算公式也确定了,不过这些公式是近似的,得出的窗口长度还要在计算中逐步修正,原则是在保证阻带衰减满足要求的情况下,尽量选择较小的N,在N和窗函数类型确定后,即可调用MATLAB中的窗函数求出窗函数wd(n)。

根据待求滤波器的理想频率响应求出理想单位脉冲响应hd(n),如果给出待求滤波器频率应为Hd,则理想的单位脉冲响应可以用下面的傅里叶反变换式求出:

hd(n)?12?????Hd(ej?)ej?nd?21

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

在一般情况下,hd(n)是不能用封闭公式表示的,需要采用数值方法表示;从w=0到w=2π采样N点,采用离散傅里叶反变换(IDFT)即可求出。

用窗函数wd(n)将hd(n)截断,并进行加权处理,得到

根据上式中的正、 负号和长度N的奇偶性又将线性相位FIR滤波器分成四类。要根据所设计的滤波特性正确选择其中一类。 例如, 要设计线性相位低通特性可选择h(n)=h(N-1-n)一类,而不能选h(n)=-h(N-1-n)一类。

验算技术指标是否满足要求,为了计算数字滤波器在频域中的特性,可调用freqz子程序,如果不满足要求,可根据具体情况,调整窗函数类型或长度,直到满足要求为止[9]。

h(n)?hd(n)?(n)如果要求线性相位特性, 则h(n)还必须满足:

h(n)??h(N?1?n)

4.1.2典型的窗函数

1、矩形窗(Rectangle Window)

w(n)?RN(n)

其频率响应和幅度响应分别为:

N?1sin(N?/2)?j?2sin(N?/2)j?W(e)?e,WR(?)?

sin(?/2)sin(?/2)2、三角形窗(Bartlett Window)

N?1?2n,0?n??2 w(n)??N?12nN?1?2?,?n?N?1N?12?其频率响应为:

N?12sin(N?/4)2?j?2j?W(e)?[]e

Nsin(?/2)3、汉宁(Hanning)窗,又称升余弦窗

12n?w(n)?[1?cos()]RN(n)

2N?1其频率响应和幅度响应分别为:

22

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

?j(2?2?W(e)?{0.5WR(?)?0.25[WR(??)?WR(??)]}eN?1N?1?W(?)e?j?aj?N?1)?2

W(?)?0.5WR(?)?0.25[WR(??2?2?)?WR(??)]N?1N?14、汉明(Hamming)窗,又称改进的升余弦窗

2n?w(n)?[0.54?0.46cos()]RN(n)

N?12?2?其幅度响应为:W(?)?0.54WR(?)?0.23[WR(??)?WR(??)]

N?1N?15、布莱克曼(Blankman)窗,又称二阶升余弦窗

2n?4n?w(n)?[0.42?0.5cos()?0.08cos()]RN(n)

N?1N?1其幅度响应为:

W(?)?0.42WR(?)?0.25[WR(??2?2?)?WR(??)]N?1N?1 4?4??0.04[WR(??)?WR(??)]N?1N?16、凯泽(Kaiser)窗

I0(?1?[1?2n/(N?1)]2)w(n)?,0?n?N?1

I0(?)其中:β是一个可选参数,用来选择主瓣宽度和旁瓣衰减之间的交换关系,一般说来,β越大,过渡带越宽,阻带越小衰减也越大。I0(·)是第一类修正零阶贝塞尔函数。

若阻带最小衰减表示为As??20log10?s,β的确定可采用下述经验公式:

?0????0.5842(As?21)0.4?0.07886(As?21)?0.1102(A?8.7)s?As?2121?As?50 As?50若滤波器通带和阻带波纹相等即δp=δs时,滤波器节数可通过下式确定:

N?As?7.95?1

14.36?F其中:

???s??p?F??

2?2?在MATLAB中,实现矩形窗的函数为boxcar和rectwin,其调用格式如下: w=boxcar(N)

23

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

w=rectwin(N)

其中N是窗函数的长度,返回值w是一个N阶的向量,它的元素由窗函数的值组成。实际上,w=boxcar(N)等价于w=ones(N,1)。

在MATLAB中,实现三角窗的函数为triang,调用格式为:

w=triang(N)

在MATLAB中,实现汉宁窗的函数为hann,调用格式如下:

w=hann(N) w=hann(N, 'sflag')

Hann函数中的参数sflag为采样方式,其值可取symmetric(默认值)或periodic。当sflag=symmetric时,为对称采样;当sflag=periodic时,为周期采样,此时hann函数计算N+1个点的窗,但是仅返回前N个点。

在MATLAB中,实现海明窗的函数为hamming,调用格式分别如下:

w=hamming (N) w=hamming (N,'sflag')

其中sflag的用法同上。

在MATLAB中,实现布拉克曼窗的函数为blackman,调用格式如下:

w=blackman (N) w=blackman (N,'sflag')

在MATLAB中,实现切比雪夫窗的函数为chebwin,调用格式为:

w=chebwin (N,r)

其中r 表示切比雪夫窗函数的傅里叶变换旁瓣幅度比主瓣低rdB(其默认值为100dB),且旁瓣是等纹波的。

在MATLAB中,实现巴特里特窗的函数为bartlett,调用格式为:

w=bartlett (N)

在MATLAB中,实现凯塞窗的函数为kaiser,调用格式为:

w=kaiser (N,beta)

其中beta为窗函数的参数β。

在MATLAB中,提供了基于窗函数法的两类设计函数,即函数fir1和函数fir2。 1、 函数fir1

该函数实现加窗的线性相位FIR数字滤波器,可设计标准低通、带通、高通和带阻滤波器[5]。

其调用格式如下: a. b=fir1(n,Wn) b. b=fir1(n,Wn, 'ftype') c. b=fir1(n,Wn,window) d. b=fir1(n,Wn, 'ftype', window) e. b=fir1(?,'normalization')

n表示滤波器的阶数。'ftype'代表所设计滤波器的类型:'high'表示高通滤波器;'stop'表示带阻滤波器;?DC?-1表示多通带滤波器,第一频带为通带;?DC?-0表示多通带滤波器,第一频带为阻带;默认时代表低通或带通滤波器。?window?为窗函

24

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

数,是长度为N的列向量,默认时函数自动取Hamming窗。b=fir1(n,Wn)可得到n阶低通FIR滤波器,调用后返回维数为n+1的行向量b,它是滤波器的系数。b与FIR滤波器的系统函数有如下关系:

H(z)?b(1)?b(2)z?1???b(n?1)z?n

对于高通和带阻滤波器,n取偶数,?n为滤波器的截止频率,范围是(0,1);对于带通和带阻滤波器,?n=[?1,?2],且?1??2;对于多通带滤波器,

?n=[?1,?2,?3,?4],频段为0????1,?1????2,?2????3,?。

2、函数fir2

该函数用于设计基于窗函数的任意响应FIR滤波器,其频率响应由向量f和向量m共同决定,取值在[0,1]之间;n为滤波器阶数;b向量为返回滤波器的系数;window为窗函数,长度为n+1,默认时为hamming窗;npt为对频率响应进行内插点数,默认时为512;lap参数用于指定fir2在重复频率点附近插入的区域大小。 4.2 FIR滤波器滤波实例 1、窗函数设计低通 程序设计如下: clear;close all

[z1,fs,bits]=wavread('E:\\耿博.wav') y1=z1(1:8192); Y1=fft(y1);

fp=1000;fc=1200;As=100;Ap=1;Fs=8000; wc=2*pi*fc/Fs; wp=2*pi*fp/Fs; wdel=wc-wp; beta=0.112*(As-8.7); N=ceil((As-8)/2.285/wdel); wn= kaiser(N+1,beta); ws=(wp+wc)/2/pi; b=fir1(N,ws,wn); figure(1); freqz(b,1); x=fftfilt(b,z1); X=fft(x,8192);

25

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

figure(2);

subplot(2,2,1);plot(abs(Y1));axis([0,1000,0,1.0]); title('滤波前信号频谱');

subplot(2,2,2);plot(abs(X));axis([0,1000,0,1.0]); title('滤波后信号频谱'); subplot(2,2,3);plot(z1); title('滤波前信号波形'); subplot(2,2,4);plot(x); title('滤波后信号波形'); sound(x,fs,bits);

图形分析如图18、图19:

滤波前信号频谱100Magnitude (dB)滤波后信号频谱10.80.60.40.210.800.60.40.2-100-20000.10.20.30.40.50.60.70.8Normalized Frequency (?? rad/sample)0.9100500滤波前信号波形100000500滤波后信号波形10000Phase (degrees)10.50-0.500.10.20.30.40.50.60.70.8Normalized Frequency (?? rad/sample)0.9110.50-0.5024x 1064-2000-4000-6000-8000-1-1024x 1064图18 FIR数字低通滤波器幅频-相位特性曲线 图 19 滤波前后信号频谱和波形对比

2、窗函数设计高通 程序设计如下: clear;close all

[z1,fs,bits]=wavread('E:\\耿博.wav') y1=z1(1:8192); Y1=fft(y1);

fp=2800;fc=3000;As=100;Ap=1;Fs=8000; wc=2*pi*fc/Fs; wp=2*pi*fp/Fs; wdel=wc-wp; beta=0.112*(As-8.7); N=ceil((As-8)/2.285/wdel); wn= kaiser(N,beta); ws=(wp+wc)/2/pi;

26

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

b=fir1(N-1,ws,'high',wn); figure(1); freqz(b,1); x=fftfilt(b,z1); X=fft(x,8192); figure(2);

subplot(2,2,1);plot(abs(Y1));axis([0,1000,0,1.0]); title('滤波前信号频谱');

subplot(2,2,2);plot(abs(X));axis([0,1000,0,1.0]); title('滤波后信号频谱'); subplot(2,2,3);plot(z1); title('滤波前信号波形'); subplot(2,2,4);plot(x); title('滤波后信号波形'); sound(x,fs,bits);

图形分析如图20、图21:

滤波前信号频谱1100滤波后信号频谱10.80.60.40.20.80.60.4Magnitude (dB)0-1000.2000.10.20.30.40.50.60.70.8Normalized Frequency (?? rad/sample)0.91-2000500滤波前信号波形100000500滤波后信号波形1000200010.50-0.5-1024x 10640.040.020-0.02-0.04024x 1064Phase (degrees)0-2000-4000-6000-800000.10.20.30.40.50.60.70.8Normalized Frequency (?? rad/sample)0.91

图20 FIR数字高通滤波器幅频-相位特性曲线 图21 滤波前后信号频谱和波形对比

3、窗函数设计带通 程序设计如下: clear;close all

[z1,fs,bits]=wavread('E:\\耿博.wav') y1=z1(1:8192); Y1=fft(y1);

27

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

fp1=1200 ;fp2=3000 ;fc1=1000 ;fc2=3200 ;As=100 ;Ap=1 ;Fs=8000 ; wp1=2*pi*fp1/Fs; wc1=2*pi*fc1/Fs; wp2=2*pi*fp2/Fs; wc2=2*pi*fc2/Fs; wdel=wp1-wc1; beta=0.112*(As-8.7); N=ceil((As-8)/2.285/wdel);

ws =[(wp1+wc1)/2/pi,(wp2+wc2)/2/pi]; wn= kaiser(N+1,beta); b=fir1(N,ws,wn); figure(1); freqz(b,1) x=fftfilt(b,z1); X=fft(x,8192); figure(2);

subplot(2,2,1);plot(abs(Y1));axis([0,1000,0,1.0]); title('滤波前信号频谱');

subplot(2,2,2);plot(abs(X));axis([0,2000,0,0.0003]); title('滤波后信号频谱') subplot(2,2,3);plot(z1); title('滤波前信号波形'); subplot(2,2,4);plot(x); title('滤波后信号波形'); sound(x,fs,bits);

图形分析如图22、图23:

100滤波前信号频谱10.80.62-4x 10滤波后信号频谱

Magnitude (dB)0-1000.40.21-20000.10.20.30.40.50.60.70.8Normalized Frequency (?? rad/sample)0.9100500滤波前信号波形1000005001000150020005000滤波后信号波形0.10.050-0.0510.50-0.5Phase (degrees)0-5000-10000-1500000.10.20.30.40.50.60.70.8Normalized Frequency (?? rad/sample)0.91-1024x 1064-0.1024x 1064

图22 FIR数字带通滤波器幅频-相位特性曲线 图23 滤波前后信号频谱和波形对比

28

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

5 总 结

本文在数字滤波器的设计过程中,采用的设计方法是基于MATLAB的数字滤波器的设计。本文通过IIR数字滤波器和FIR数字滤波器的设计实例,说明如何利用MATLAB来完成数字滤波器的设计。窗函数法中相位响应有严格的线性,不存在稳定性问题, 设计简单。双线性变换不会出现由于高频部分超过折叠频率而混淆到低频部分去的现象,但会产生频率混碟现象,使数字滤波器的频响偏移模拟滤波器的频响。

在对语音信号进行滤波的时候,双线性低通滤波器的滤波效果最好,滤波后的语音信号失真比较小,高通和带通的失真都较大。由此可知人的语音信号的能量主要集中在低频部分。

在用窗函数设计滤波器的时候,由于所给指标的阻带衰减比较大,对所选窗的要求比较苛刻,而且很难达到所给指标。由公式N>=a*fc/|fp-fs|或公式N=1+(As-7.95)*fc/[14.36*|fp-fs|]可算得阶数N>=3。但实际在仿真的过程中发现要想达到阻带的最小衰减为100,阶数N得选很大,约在800左右。对如此大阶数的滤波器在实现的时候很困难。

由滤波器的频谱图和滤波前后的语音信号的频谱图对比可知本设计选用双线性变换法设计的IIR滤波器比较好。在同样的技术指标的要求下,IIR滤波器所要求的阶数N也比较小,实现起来比较容易。

事实上IIR滤波器系统函数的极点可以在单位圆内的任何位置,实现IIR滤波器的阶次较低,所用的存储单元较少,效率高,且IIR数字滤波器能够保留一些模拟滤波器的优良特性。但是这些特性是以牺牲线性相位频率特性为代价的,即用Butterworth、chelbchev滤波器的幅度频率特性,得到的滤波器往往是非线性的。在许多电子系统中,对幅度频率特性和线性相位特性都有较高的要求,所以IIR滤波器在这些系统中往往难以胜任。有限长单位冲激响应(FIR)数字滤波器具有以下优良的特点:可在设计任意幅度频率特性滤波器的同时,保证精确、严格的线性相位特性。FIR数字滤波器的单位冲激响应h(n)是有限长的,可以用一个固定的系统来实现,因而FIR数字滤波器可以做成因果稳定系统。

29

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

参考文献

[1] 丁玉美,高西全等.数字信号处理[M].西安:西安电子科技大学出版社,2000. [2] 程佩青,数字信号处理教程[M].北京:清华大学出版社,2003.

[3] 罗军辉,罗勇江等.MATLAB7.0在数字信号处理中的应用[M].北京:机械工

业出版社,2005.

[4] 郭仕剑,王宝顺等.MATLAB7.X数字信号处理[M].北京:人民邮电出版社,

2006.

[5] 门爱东,杨波等.数字信号处理.北京:人民邮电出版社,2003.

[6] 余成波,杨菁等.数字信号处理及MATLAB实现[M].北京:清华大学出版社,2005.

[7] T.W.ParksandC.S.Burrus.DigitalFilterDesign.NewYork,NewYork: Wiley-Interscience,1987.

[8]S.K.Mitra.Digital Signal Processing:A Computer-Based Approach. NewYork,NewYork:McGraw-Hill,thirded,2006.

[9] 赵琳,王淑伟,邢帆.基于 Matlab的数字滤波器设计方法研究[J].产业与科技论坛.2008(9).157~158

[10] 张亚妮.基于MATLAB的数字滤波器设计[J].辽宁工程技术大学学报.2005年10月第24卷第5期.716~718

[11] 阎晓艳,傅丰林等.FIR数字滤波器的设计及其在MATLAB中的仿真实现[J].

电子科技.2004(5).43~46

30

淮北煤炭师范学院学士学位论文 基于MATLAB的数字滤波器设计

致 谢

诚挚地感谢邹锋等老师的帮助和精心指导!本文是在邹锋老师的指导和支持下完成的。如果没有老师的指导,我想我是很难在这么短的时间内顺利完成论文的设计。邹老师的治学态度、敬业精神对我产生了很大影响。他的教导和思路给予我无尽的启迪。

最后,再次对在大学四年里关心和帮助我的老师表示衷心地感谢!

31

本文来源:https://www.bwwdw.com/article/m5so.html

Top