第 5章 IIR数字滤波器设计
5.1 数字滤波类型与指标
5.2 模拟滤波器设计
5.3 设计 IIR滤波器的脉冲响应不变法
5.4 设计 IIR滤波器的双线性变换法
5.5 设计 IIR数字滤波器频率变换法
5.6 数字陷波器设计
5.7 IIR数字滤波器的计算机辅助设计第 5章 IIR数字滤波器设计滤波的目的
为了压制输入信号的某些频率成分,从而改变信号频谱中各频率分量的相对比例
广义滤波包括对信号的检测与参量的估计
信号的检测,确定在干扰背景中信号是否存在
信号参量的估计,为识别信号而确定信号的某一个或某几个参量的估值
5.1 数字滤波类型与指标滤波技术包括,
滤波器设计,根据给定滤波器的频率特性,
求得满足该特性的传输函数,
滤波过程的实现:获得传输函数后,以何种方式达到对输入信号的进行滤波的目的
1.数字滤波器的频率特性数字滤波器
具有某种特定频率特性的线性时不变系统
广义上,任何线性时不变离散系统都是一个数字滤波器设计数字滤波器的任务
寻求一个因果稳定的线性时不变系统,使其系统函数 H(z)具有指定的频率特性

0
)()()(
n
nj
ez
j enhzHeH
j


0
)()()(
n
nj
ez
j enhzHeH
j

对因果稳定的线性时不变系统:
)(?jeH,滤波器的传输函数
H(z):系统函数 h(n),滤波器的单位脉冲响应
)()()( jj eHeH?
)]([Im)]([Re)( 22 jj eHeHH
)](R e [
)](I m [)(
j
j
eH
eHa r c t g?
)(?H,幅度响应 )(,相位响应
DF按频率特性分类可分为低通、高通、带通、带阻和全通特点为
数字频率以 周期
2
频率特性只限于 范围,依取样定理,对应于实际模拟抽样频率的一半
)( sfT频率变量以数字频率 表示其中 模拟角频率,T抽样时间间隔,f
s 抽样频率?
理想滤波器的频率响应
0
s?
p1
p?,通带波纹
s?,阻带波纹
,过渡带
:通带截止频率p?
:阻带截止频率s?
DF的性能要求(低通为例)
11
11
2?
c? st

) (?j e H
从信号不失真角度讲通常要求
相位线性
)( 为时延常数?
– 具有群恒时延特性常数 )(
)(


d
d
)(相位响应
2,IIR和 FIR数字滤波器
IIR滤波器的系统函数通常可表示成的有理分式
FIR滤波器的系统函数则可表示为的多项式设计过程一般包括以下三个基本问题:
根据实际要求确定滤波器性能指标;
用一个因果稳定的系统函数去逼近这个指标;
用一个有限精度的运算去实现这个传输函数
问题 1,3与实际的要求及实现的硬件条件有关
本章主要讨论问题 2,即系统函数的设计 (或逼近 )问题。
3.设计 IIR滤波器的几种方法
IIR数字滤波器的系统函数可表示为的有理分式
1
1
1
0
1
)(
za
zb
zH
k
N
k
k
N
k
设计 IIR滤波器的系统函数,就是要确定 H(z)的阶数 N( 通常称 N为滤波器的阶数)以及分子分母多项式的系数
jezjkk zHeHba )()(,使其、
满足指定的频率特性
(1)利用模拟滤波器的理论来设计模拟滤波器研究较早,理论已经十分成熟,有许多简单而严谨的设计公式和大量的图表可以利用,利用这些现有技术来解决数字滤波器的设计问题采用这种方法时,要先要设计一个合适的模拟滤波器,然后将它转换成满足给定指标的数字滤波器这种方法适合于设计幅频特性比较规则的滤波器,例如低通、高通、带通、带阻等当把模拟滤波器的 H(s)转换成数字滤波器的 H(z) 时,
要实现 S平面向 Z平面的映射,必须满足两个条件必须保证模拟频率映射为数字频率,且保证两者的频率特性基本一致
要求变换后代表 S平面的虚轴 jΩ应映射到 Z片面的单位圆
且数字滤波器的频率响应和模拟滤波器频率响应的形状应基本保持不变;
因果稳定的模拟滤波器系统函数 H(s)转换成数字滤波器传输函数 H(z)后,仍然是因果稳定的
要求 S平面左半平面的极点必须映射到 Z平面的单位圆内两种常用的方法
脉冲响应不变法:从时域的角度出发进行映射
双线性不变法:从频域角度出发进行映射
(2)利用最优化技术进行 CAD设计若需设计滤波器的幅频特性是任意的或者形状比较复杂,可采用计算机辅助设计 (CAD)方法进行优化设计设计思想
)(?jd eH希望滤波器的幅频响应,
)(?jeH设计 滤波器的幅频响应,
选择一种最优化的准则,例如采用最小均方误差准则
)()( jjd eHeH,设 在指定的一组离散的频率点
Mii,,2,1},{ 的均方误差
2
1
2 ])()([ jj
d
M
i
eHeH
求解 H(z)的系数,,kk ba 使 均方误差 最小
当 滤波器阶数 N 较高时,转换为 一个多变量最优化问题,需要大量的迭代运算,因此必须采用 CAD的方法。
5.2 模拟滤波器设计
IIR滤波器的设计是基于模拟滤波器的成熟技术而完成的简单介绍模拟滤波器设计的一些基本概念,并介绍两种常用的滤波器的设计方法,
巴特沃思 (Butterworth)滤波器
切比雪夫 (Chebyshev)滤波器
5,2,1模拟滤波器设计的基本概念
1,模拟滤波器的频率特性与衰减特性滤波器的频率特性主要取决于构成滤波器系统的系统函数
jssHjH |)()(
工程设计中给定的指标往往是通带和阻带的衰减,它一般用反映功率增益的幅度平方函数或称模方函数来定义
dBjHjHA |)(|lg20|)(|lg10)( 2
当要求滤波器具有线性相位特性 ( 延时 τ 为常数 ) 时滤波器的频率特性为
H j H j e j( ) | ( )| ( ) )(,
2,归一化与频率变换采用归一化参数
设计结果具有普遍性
计算方便归一化包含:
电路参数归一化:将系统中无源元件的阻抗或运算阻抗分别除以基准电阻 (系统的负载电阻值 )。
频率归一化:将所有的频率都除以基准频率 ( 滤波器的截止频率 )
计算实际电路参数时应要将归一化频率乘以截止频率,进行反归一化
频率变换:从归一化低通原型滤波器到高通,带通,带阻等其它类型的滤波器的变换方法
3,从模方函数 求模拟滤波器的系统函数 H(s)
当不含有源器件,作为一个因果稳定、物理可实现的系统函数必须满足的条件
| ( )|H j? 2
a,是一个具有实系数的 s有理函数 )(/)()( sDsNsH?
b,所有极点必须全部分布在 s的左半平面内
c,分子多项式式 N(s)的阶次必须小于或等于分母多项式
D(s)的阶次
—— 正实函数实函数的傅立叶变换存在共轭对称的性质
H j H j* ( ) ( )
)()(|)(| 2 jHjHjH )()( sHsHjs
| ( )| ( ) * ( ) ( ) ( )H j H j H j H j H j2有得平面的虚轴,解析延拓代表 sj?
从给定的模方函数求出所需要的系统函数的方法
a,解析延拓,令 s= 代入模方函数得到,
并求其零极点
j? H s H s( ) ( )?
b,取 所有在左半平面的极点作为 的极点)(sHH s H s( ) ( )?
c,按需要的相位条件 (最小相位,混合相位等 )取一半的零点构成 的零点)(sH
H s H s( ) ( )?
4,模拟滤波器的设计 -逼近问题
pA
sA
:通带衰减
:阻带衰减
,
:与通带衰减、阻带衰减有关的系数
:通带截止频率
:阻带截止频率
p?
s?
寻找一个恰当的近似函数来逼近理想特性 —— 谓逼近问题最常用的具有优良性能的滤波器:
巴特沃思 (Butterworth)滤波器
切比雪夫 (Chebyshev)滤波器
椭圆 (elliptic)函数或考尔( Cauer) 滤波器
实现线性相位的贝塞尔滤波器
)1l o g (20)1 1l o g (10 2 ppA
ssA l o g20)1
1l o g (10
2
之间的关系、与 sP AA,
5,2,2巴特沃思 Butterworth低通滤波器
1.基本性质
BW滤波器以巴特沃思函数来近似滤波器的系统函数
BW的低通模平方函数表示
,,2,1)/(1 1|)(| 22 NjjjH N
c
指定,后,带 到上式,得pAp? p
21.0
22
2 )1(10
1
1
)/(1
1|)(|
p
A
N
cp
p
pjH?

10 10 1,A p
13 =时,=当?dBA p
指定,后,带 到上式,得sAs? s
21.0
22
2 10
1
1
)/(1
1|)(|
s
A
N
cs
s
sjH?

10 10 1,A s
用 3dB截止频率 来规一化:对频率进行,下式变为?
c/ c
,,2,1
)/(1
1|)(|
2
2

N
jj
jH N
c
| ( )|
( )
H j N?
2
2
1
1
讨论:
| ( )| ( )H j N2 211
当 =0 时,=1,取最大值? | ( )|H j? 2
当 时,=0.5,取 3dB值)(1
c | ( )|H j?
2
1)(/,1/ 22 接近时,很小,通带 jHNNcc
0)(/,1/ 22 接近时,很大,通带 jHNNcc
阻带内,由于
| ( )|
( )
H j j
j c
N

2
2
1?

A H j N jjs
c
10 202l g | ( )| l g ( )
幅度随着 N的增加阻带衰减近似为 6N db/倍频程。
N越大,频带特性越接近理想矩形特性
| ( )| ( )H j N2 211
上式的台劳级数展开为,
NNjH 422 1|)(|
12,,2,1,0|)(| 02 NkjHdd
k
=0处函数对 2N— l阶导数都等于零 —— 曲线在 =0附近是最
,平坦,,—— 巴特沃思滤波器又叫做,最大平坦滤波器,

归一化巴特沃思低通滤波器的幅度特性
2.设计过程
(a)按给定指标确定阶次 N
222 )(
N
p
s
)/lg (
)]110/()110lg [ (
2
1
)/lg (
)/lg ( 1.01.0
ps
AA
ps
ps
N




)/lg (
)110lg (
)/lg (
)lg (
)/lg (
)/lg ( 1.0
cs
A
cscs
s
N





实际计算时,要对上式求得的数值取整加 1。
若给定的指标 =3dB,即通带边频 时,ε =1,可求得
cppA
(b)从模方函数求系统函数 H(s)
①求得极点
,,2,1)/(1 1|)(| 22 NjjjH N
c
带入上式,得: js
0)(1 2 N
cj
s
)12()1( kje由于
Nkes NNkjck 2,,2,12/)12(
分析讨论
c在归一化频率的情况 =1,极点均匀分布在单位圆上
Nkes NNkjk 2,,2,12/)12(
对于物理可实现系统,它的所有极点均应在 s的左半平面上
②系统函数的构成滤波器的极点求出后,可取左平面上的所有极点构成系统函数
)(
1
)(
1
i
N
i
ss
AsH
对于低通滤波器,为了保证在频率零点 =0处,=1,可取? |)(|?jH

N
i
i
N sA
1
)1(
)(
)1()(
1 i
i
N
i
N
ss
ssH

因此得例 5-2-1举例说明系统函数的构成设计一巴特沃思滤波器,使其满足以下指标,通带边频 =100k
rad/s,通带的最大衰减为 =3dB,阻带边频为 =400k rad/s,
阻带的最小衰减为 =35 dB
p?
Ap s?
As
解,由于通带边频就是 3dB 截止频率,即
cp
1110 1.0 pA? 2.56110 1.0 sA?
确定阶次 N
,9.24lg 2.56lg)/lg ( )/lg (
cs
N 3?N取求左半平面的极点,
s e kk c j k N N ( ) /,,2 1 2 1 2 3?
,3/21?jc es,2?jc es s ec j3 2 3 /
得极点:
构成巴特沃思滤波器传输函数 H(s)为
H s s s ss s s s s s s s sc
c c c
( ) ( )( )( )1 2 3
1 2 3
3
3 2 2 32 2

相对截止频率 归一化,得归一化巴特沃思滤波器传输函数?
c H sa ( )
H s s s sa ( )12 2 13 2
③一般 N阶归一化巴特沃思滤波器传输函数 表示
H s
s s a s a s a s s
a
i
N
i
N
N N( )
( )

1 1
1
1
1 2
2
1
1?
是 =1时的极点,分布在单位圆上s
i?c
分母一般称为巴特沃思多项式,其系数可通过查表求得,
见表 5-2-1
表 5-2-1 巴特沃思多项式系数
N a1 a2 a3 a4 a5 a6 a7 a8 a9
2 1.4142
3 2.0000 2.0000
4 2.6131 3.4142 2.6131
5 3.2361 5.2361 5.2361 3.2361
6 3.8637 7.4641 9.1416 7.4641 3.8637
7 4.4940 10.097 14.592 14.592 10.097 4.4940
8 5.1528 13.137 21.846 25.688 21.846 13.137 5.1528
9 5.7588 16.581 31.163 41.986 41.986 31.163 16.581 5.7588
10 6.3925 20.431 42.802 64.882 74.233 64.882 42.802 20.431 6.3925
表 5-2-2 巴特沃思多项式因式分解
N 巴 特 沃 思 多 项 式
1 s+1
2 s2+1.4142s+1
3 (s+1)(s2+s+1)
4 (s2+0.7654s+1)(s2+1.8478s+1)
5 (s+1)(s2+0.6180s+1)(s2+1.6180+1)
6 (s2+0.5176s+1)(s2 +1.412s+1)(s2 +1.9319s+1)
7 (s+1)(s2+0.4450s+1)(s2+1.2470s+1)(s2+1.8019s+1)
8 (s2+0.3092s+1)(s2+1.1111s+1)(s2+1.6629s+1)(s2+1.9616s+1)
9 (s+1)(s2+0.3473s+1)(s2+s+1)(s2+1.5321s+1)(s2+1.8794s+1)
上述归一化公式和表格是相对 3dB 截止频率给出的。由指定的技术指标 利用上述公式和表格进行设计时,最关键的 2个参数是滤波器的节数 N和 3dB 截止频率 。
N用来求巴特沃思多项式,用来反归一化,求实际滤波器的参数。
ssPp AA,,,
c
c
c
5,2,3 切比雪夫滤波器
5,3 设计 IIR滤波器的脉冲响应不变法
1.设计的基本原理和方法原理,从时域响应出发,使求得的数字滤波器的单位脉冲响应 h(n)等于模拟滤波器的单位冲激响应 h(t)的抽样值。
nTtthnh )()(如果:
][)(
,)(
1 )(则由:
已知
sHLth
sH

则可有下式求的 H( z):
)]([)]([)( nThZnhZzH
方法,将 H(s)表示为部分分式形式
k
k
N
k ss
AsH
1)(
其拉氏反变换为
得到数字滤波器的单位脉冲响应
对上式两边取 Z变换得
)()(
1
tueAth tsk
N
k
k?



N
k
nTs
k
N
k
nTs
k nueAnTueAnh
kk
11
)()()()(
1
1 1
)(?

zeAzH Ts k
N
k k
( 5-2-3)
如果模拟滤波器的系统函数是稳定的,其极点应位于左半平面
1
1 1
)(?


ze
AzH
Ts
k
N
k
k
( 5-2-3)
0][?ke sR
对 Z平面的极点有 1 Tsk kez
位于单位园内。因此 H(z)是一个稳定的离散系统函数,这说明由一个稳定的模拟滤波器得到了一个稳定的数字滤波器
2.脉冲响应不变法设计的滤波器的频率响应根据抽样定理,序列 h(n)的频谱是原模拟信号频谱的周期延拓原模拟滤波器的频率响应为,)(?jH 由于 h(n)是 h(t)的等间隔抽样
H e
T
H j
T
j
T
mj
n
( ) ( )

1 2
如果模拟滤波器的频率响应是带限于折叠频率之内,即这样数字滤波器的频率响应才能等于模拟滤波器的频率响应然而,高通和带阻滤波器不能满足 (-2-)式的要求,将会产生混叠脉冲响应不变法不适合用来设计高通和带阻数字滤波器。
H j Ts( )0 2?
H e T H j Tj( ) ( )1
3.几点修正
1)、消去 T的影响
H e T H j Tj( ) ( )1
由 上 式可见,数字频率响应与模拟频率响应的第一差别是具有一个乘法因子 (1/T)= fS,当采样频率 fS很高时,将会使滤波器的增益很大,这往往是不希望的,为此可对 下 式作修正,
nTtthnh )()(
h n Th t t nT( ) ( )

2)、直接用数字频率表示的求 H(z)的公式在实际滤波器设计中,因模拟滤波器系统函数的表格大都是归一化低通原型,其滤波器 3dB点截止频率都归一化在?
c? 1H sa ( )
原因:可将设计公式及有关参数表格化,使之更通用。我们只要知道滤波器的阶数,就可直接查出低通原型的系统函数。
H sa ( )?c
当滤波器的实际截止频率不等于 1时,须进行所谓反归一化,
以 (s/ )代替 中的 s,
即实际低通滤波器的系统函数 H(s)应为
H s H sa c( ) ( / )
H s H s A
s s
A
s sa c k
N
k
c k k
N
k c
k c
( ) ( )
( / )


1 1
H s H sa c( ) ( / )H s A
s sk
N
k
k
( )?
1
为模拟归一化原型系统函数的极点sk
H z Z Th t TA
e zk
N
k c
s Tk c( ) [ ( )]

1
11
有由,Tcc
H z
A
e zk
N
k c
s k c( )

1
11
例 5-3-1利用脉冲响应不变法设计一个 4阶巴特沃斯型数字低通滤波器,满足以下指标
(A) 若采样周期 T=10μ s,求实际模拟截止频率 fc,
(B) 3dB截止频率 =0.2π 弧度 。
k H zfsr a d fT cccc 102,/102010 2.0 35
有由,Tcc解:先计算模拟截止频率,
设计数字低通滤波器分三步:
4,,1,
,
)18 4 7 8.1)(17 6 5 4.0(
1
)(
8
)142(
22



kes
ssss
sH
k
j
k
a
第一步 查巴特沃斯数字低通滤波器原型表,求得系统函数第二步 部分分式分解并求 Ak
H s
A
s s
A
s s
ka
k
k
k
k
k i i k
( ),
( )
,,,?


1
4 4 1
1 4?
第三步 将 代入下式,kck sA 及,?
整理并化简求得 H(z)的实系数二次形式
H z
A
e zk
N
k c
s k c( )

1
11
H z
A
e zk
k c
s k c( )

1
4
11



10
1 84776 0 88482
1 1 31495 0 61823
1 84776 0 40981
1 1 08704 0 31317
1
1 2
1
1 2
.,
.,
.,
.,
z
z z
z
z z
5,4 设计 IIR滤波器的双线性变换法
1.设计方法
)(?jH)(?jeH
从频域响应出发,直接使数字滤波器的频域响应,逼近模拟滤波器的频域响应,
进而求出 H(z)。
脉冲响应不变法的主要缺点,
对时域的采样会造成频域的混叠效应,因而有可能使设计的数字滤波器的频域响应与原来模拟滤波器的频域响应相差很大,而且不能用来设计高通和带阻滤波器
– 原因,从 S平面到 Z平面的映射是多值的映射关系双线性变换的映射过程脉冲响应不变法的映射过程双线性变换法的改进为避免 频域的混叠,分两步完成 S平面到 Z平面的映射
TT,?将 S平面压缩到某一中介的 S1平面的一条横带域
通过标准的变换将此横带域映射到整个 Z平面上去,
实现方法:
tg T( )12
再通过 Z变换,将 Ω 1映射到 Z平面的单位圆上通过下面的正切变换,
将 S平面的 jΩ 轴压缩到 S1平面的 jΩ 1轴上的 TjTj,
Tsez 1?
将正切变换 延拓到整个 S平面,得到 S平面到 S1平面的映射关系
s th
s T e
e
s T
s T
( )
1
2
1
1
1
1
再将 S1平面按关系式 映射到 Z平面得到z e s T? 1
s z
z

1
1
1
1
z
s
s
1
1
双线性变换或双线性变换的映射关系满足关于映射关系可行性的两个条件
( 1) S平面的虚轴映射到 Z平面的单位圆上;
( 2)位于 S左半平面的极点应映射到 Z平面的单位圆内。
s
e
e
jtg j
j
j?

1
1 2
( )?
,?jez?令 s z
z

1
1
1
1
带入表达式得:
说明 S平面的虚轴 映射成了 Z平面的单位圆
s j令,带入表达式
z ss11
得:
z jj11
22
22
)1(
)1(


z
显然当 时,<1, 0 z
S平面的左半平面轴映射到了 Z平面的单位圆内,保证系统函数经映射后稳定性不变双线性变换的频率对应关系模拟频率与数字频率是一种非线性的关系
tg T( )12
模拟滤波器与数字滤波器的响应与对应的频率关系上发生了畸变,
也造成了相位的非线性变化,这是双线性变换法的主要缺点双线性变换法 除了不能用于线性相位滤波器设计外,仍然是应用最为广泛的设计 IIR数字滤波器的方法。
在?上刻度为均匀的频率点映射到?上时变成了非均匀的点,而且随频率增加越来越密
2.频率预畸变为了保证各边界频率点为预先指定的频率,在确定模拟低通滤波器系统函数之前必须按下式进行所谓频率预畸变
c ctg? ( )?2
然后将预畸变后的频率代入归一化低通原型 Ha(s) 确定
H s H sa c( ) ( / )
最后求得数字系统函数
H z H s s z
z
( ) ( )
1
1
1
1
1
1
1
11)(

z
zsa
c
sH
3.有关双线性变换公式的说明(抽样间隔 T的选取有些文献中双线性变换的关系为,
1
1
1
12

z
z
T
s s z
z

1
1
1
1
与右式有一个 2/T的系数的差别。现说明 T的取值。
上式可由用数值计算的方法求解一阶输入 -输出模拟系统的过程中推得。设一阶模拟系统函数为
bs
asH a
)(
对应的一阶微分方程为
)()()(' taxtbyty aaa
)(tya将 表示成 的积分形式)(' ty a
)()(')( 00 tydtyty att aa
1
1
1
12

z
z
T
s
采用梯形法近似计算定积分,令步长为 T,并设 t=nT,
t0=( n-1) T,用 T代替 d?,以上积分式可近似为:
))1(()])1((')('[21)( TnyTTnynTynTy aaaa
带入一阶微分方程
)])1(())1(()()([
2
))1(()(
TnaxTnbynTaxnTbyT
TnynTy
aaaa
aa


令,)()( nTxnx a? )()( nTyny a?
)]}1()([)]1()([{2))1(()( nxnxanynybTnyny
一阶差分方程化为:
b
z
z
T
a
zX
zY
zH

1
1
1
12)(
)(
)(
对上式两边进行 Z变换,得到:
得到了双线性变换的关系
1
1
1
12

z
z
T
s
其中 T为用梯形法近似计算定积分的步长,或将模拟信号离散为抽样信号时的 抽样间隔由以上的推导过程可见,T的取值可以任选,只要满足 Ts
即只要保证将模拟频率带限在之间,不会产生因多值映射产生频率混迭现象。
另一方面,当变换关系采用下式时
1
1
1
12

z
z
T
s
将有关系式
sT
sTz
)2/(1
)2/(1

)2(2?tgT

1
1
)2/t a n (
1
1
12
)2/t a n (
1
2?


z
z
z
z
T
Tsp
ccc
当由归一化的模拟变量 p变到数字变量 Z的过程中,系数 2/T已被消去。这再一次说明 T取值的大小是不重要的。
考虑到利用双线性变换法设计数字滤波器的过程
c
spa pHsH

|)()(
1
1
1
12)()(

z
z
Ts
sHzH
例:利用双线性变换法设计巴特沃斯型数字低通滤波器设计参数:通带数字截止率,通带内最大衰减阻带数字截止率,阻带内最小衰减
dBA p 5.0?
dBAs 15?
25.0?p
55.0?s
解 第一步 进行频率预畸变求
sp,
0,4 1 4 2 1 3 6t g 0,1 2 5)2( pp tg
1,1 7 0 8 4 9 6t g 0,2 7 5)2( ss tg
第二步 计算巴特沃斯型数字低通滤波器的阶数
2,6 5 8 6 9 9 7
)/l g (
)]110/()110l g [ (
2
1
)/l g (
)/l g ( 1.01.0?



ps
AA
ps
ps
N
取 N=3
第三步 通过查表求得 3阶巴特沃斯低通滤波器原型的系统函数
)1)(1(
1)(
2 ssssH a
第四步 求 3dB截止频率,并反归一化求得
c? H s H sa c( ) ( / )
pA
N
cp
pjH
1.0
2
2 10
)/(1
1|)(|


解得:,反归一化求 得0,5 8 8 1 4 8 c
)3 4 5 9 1 8.05 8 8 1 4 8.0)(5 8 8 1 4 8.0(
2 0 3 4 5 1.0)/()(
2 ssssHsH ca
第五步 求得 H(z)代入双线性变换公式,化简后求得
1
1
1
1)()(

z
zssHzH
)3917468.06762858.01)(259328.01(
)1(0662272.0)(
211
31



zzz
zzH
1
1
1
11)()(

z
zsa
c
sHzH
也可按下式将第三步与第四步可合并在一起进行,求出 H(z)
设计滤波器的幅度响应与相位响应进行频率变换有以下两种基本方法,
5,5 设计 IIR数字滤波器频率变换法第一种方法的简化的形式:
找出归一化模拟低通原型与数字高通,带通和带阻滤波器之间的从 S域到 Z域的变换关系
直接由归一化模型低通原型变换成所需的数字滤波器
5,5,1从 S域到 Z域的频率变换法
1.归一化模拟低通原型到数字高通滤波器的频率变换设归一化模拟低通原型滤波器的系统函数为 )( pH aL
p为模拟域内的拉氏变量,模拟域内从低通到高通的变换为以 p-1代替 p:
H p H paH aL( ) ( / )? 1
反归一化,即以 带入上式
csp
求得反归一化后的高通滤波器的传输函数 H(s)
H s H s H saH
c
aL
c( ) ( ) ( )
双线性变换,得数字高通滤波器的系统函数 H(z)
1
1
1
1
1
1
1
1 )()()(


z
zpaL
z
zs
c
pHsHzH
H p H paH aL( ) ( / )? 1
变换关系
1
1
1
1

z
zp
c
模拟滤波器与数字滤波器的频率以及 3dB截频 与 之间的关系,c
c?
2
c
c tg

直接由归一化低通原型变换成数字高通滤波器的由,可得
2
2 c
c c t g

cjp /
例 设计一个三阶巴特沃斯型高通数字滤波器,3dB
数字截频为 0.2π弧度,求滤波器的系统函数。
解 三阶巴特沃斯型归一化模拟低通原型的系统函数为
H p
p p paL
( )?

1
2 2 13 2
频率预畸变
/ s )0,3 2 4 9 ( r a d)t g ( 0,12 cc tg
得到数字高通滤波器的系统函数为
1
1
1
1
1
13249.023
1
1 122
1)()(


z
zp
z
zpaL ppppHzH
c
321
31
5 7 0 1 6.02 4 1 9 7 6.23 3 5 8.38 5 9 2 0 8.1
)1(



zzz
z
解 1、频率预畸变求模拟低通的截频 和阻带下边频?c
s?
c?
s?
例 5-5-2设计一巴特沃斯型高通数字滤波器,3dB的 = 0.2π弧度阻带下边频 = 0.05π弧度,阻带衰减 AS≥48dB,求滤波器的系统函数
/ s )0,3 2 4 9 ( r a d)t g ( 0,12 cc tg
a d / s )1,3 4 1 2 6 7 ( r 22 scs c t g?
2、计算低通滤波器的阶数 N
,3,8 9 7 6)/lg ( )110lg (21
1.0

cs
A s
N
取 N= 4
4,利用直接变换关系求数字滤波器的系统函数 H(z)
1
1
1
1)()(

z
zpaL
c
pHzH
4321
41
4324.04366.2346.54743.53103.2
)1(



zzzz
z
3、查表 5-2-2,得 N=4时巴特沃斯低通原型滤波器的系统函数为
H p p p p paL ( ) (,)(,)10 7654 1 1 8478 12 2
设计的高通滤波器的幅度响应与相位响应
2.归一化模拟低通原型到数字带通滤波器的频率变换直接寻求模拟低通到数字带通之间的映射关系
l?
0
u?
sl?
su?
s?
,中心频率
:下边频和上边频
:下阻带上边频
:上阻带下边频
:模拟的阻带下边频
0 0 00? z e sj
(零点 )
0 1,? z s j?
(极点 )
数字带通到模拟低通的映射应满足的映射关系满足以上关系的变换式为,
s s z e z ez z z zz
j j
( ) ( )( )( )( ) c o s
0 0
1 1
2 1
1
2
0
2
稳定变换,S平面稳定的函数变换到 Z平面也是稳定的例如设 0<z=r<1,则有
01 )c o s1(2)1(1 1c o s2 2 0
2
2
0
2
r rrrrrs
表明,Z平面单位圆内的极点变到了 S平面的左半平面
jez?其次令 代入上式,可得频率变换关系
jje ees j
jj


s i n
c o sc o s
1
1c o s2 0
2
0
2
表明,Z平面单位圆变换到了 S平面的虚轴模拟低通与数字带通的频率变换
su
su
sb?

s in
c osc os 0
1c当 时,反归一化处理
csp /
)1(
1c o s2
2
0
2


z
zzsp
cc
)
2
c o s (
)
2
c o s (
c o s 0
lu
ul


u
u
c?

s in
c osc os 0
例 5-5-3设计一个三阶巴氏数字带通滤波器对模拟信号滤波,上下边带的 3dB截止频率为 =12.5kHz,=37.5kHz,阻带边频 = 45.125
kHz,采样频率为 =100kHz。 并求阻带的最小衰减。
f2 1f 2sf
sf
解,1)先确定上下边带的 3dB数字截止频率,及上阻带下边频
r a dff s 75.010/105.372/2 5311
r a dff s 25.010/105.122/2 5322
r a dff ssu 925.010/10125.452/2 532
2)确定数字中心频率,
0
)
2
c o s (
)
2
c o s (
c o s
21
21
0


5.00?
得:
3)求模拟低通截止频率和阻带下边频
c1 0 1
1
0 5 0 75
0 75 1
c o s c o s
s i n
c o s,c o s,
s i n,


1 6 5 3.4)925.0s i n ( )925.0c o s ()5.0c o s (s i n c o sc o s 0
u
u
s
4)求,当 =1时,三阶巴氏低通滤波器的传输函数为,H sa ( )?c
H s s s sa ( )12 2 13 2
5)求 H(z),将 代入上式得s z z
z
z
z?


2
0
2
2
2
2 1
1
1
1
cos?
H z H s z z zza s z
z
( ) ( )

2
2
1
1
2 4 6
4
1
2
1 3 3
3
fs26)求相对 =45.125kHz处的阻带的最小衰减。由下式
sA
N
cs
sjH
1.0
2
2 10
)/(1
1|)(|

得:
])/(1 1[l o g10 210 N
cs
sA
代入 和 及 N=3,求得?c s?
dBA s 18.37?
例 5-5-3设计的带通滤波器的幅度响应与相位响应直接由归一化模拟低通原型求数字滤波器的频率变换公式表表 5-5-1
5,5,2 数字域频率变换法在已知数字低通滤波器时,通过在 Z域内的数字频域变换得到所需类型的数字滤波器假定变换前的平面为 z平面,变换后的为 Z平面,这种变换就是要将 z平面的系统函数 映射成 Z平面的系统函数 。
假设从 z到 Z的映射关系为
)(zHL )(ZHd
)( 11 ZGz
)( 11|)()( ZGzLd zHZH
则有:
变换函数 必须满足以下三个条件:)( 1?ZG
( 1)为满足一定的频率响应要求,z域的频率必须变换成 Z域的频率,也就是说平面的单位圆必须映射到 Z平面的单位圆上;
( 2)为保证由因果稳定的系统变换得到的系统也是因果稳定的,要求平面的单位圆内部必须映射到 Z平面的单位圆内部;
( 3)由于 是 z-1的有理函数,为了保证变换后的 也是 Z-1的有理函数,要求变换函数 必须是 Z-1的有理函数。
)(zHL )(ZHd
)( 1?ZG
设 和 分别为 z 平面和 Z 平面的数字频率变量,即
jezjeZ?
得带入 )( 11 ZGz
)](a r g [|)(|)( jeGjjj eeGeGe
有:
1|)(|jeG
)](a r g [ jeG
表明函数 在单位圆上的幅度必须恒等于 1—— 全通函数,任何全通函数都可以表示为
N
i i
i
Z
ZZG
1
1
*1
1
1
)(
)( 1?ZG
全通函数的特点保证极点在单位圆内,即,以保证变换后稳定性不变1||?i?
N
i i
i
Z
ZZG
1
1
*1
1
1
)(
)( 1?ZG 的所有零点都是其极点的共轭倒数,)/1( *i?
当 由 0变到 时,全通函数相角 的变化量为 。
选择合适的 N和,则可得到各类变换
)](a r g [?jeGN
ia
1.数字低通 — 数字低通)(zHL )(ZHd
和 都是低通系统函数,只是截止频率不同)(zH
L )(ZHd
从 0变到 时,也应从 0变到,由全通函数相角变化为,
可确定此处阶数应为,故变换函数为一阶全通函数
N
1?N
1
1
11
1)(?


Z
ZZGz
jj eZez,带入,注意?为实数。 得频率间关系为,
j
j
j
e
ee

1









c o s1
s i na r c t a n2
c o s)1(2
s i n)1(a r c t a n
2
2
由此求得 与 之间的关系:
,线性=时, 0?
时,频率压缩0
时,频率扩张0
对于幅度响应为分段常数的滤波器,变换后仍可得类似的频率响应设低通原型的截止频率为,变换后的对应截止频率为
c? c?
带入下式得:
j
j
j
e
ee

1?


2
2
s in
cc
cc


由截频为 的 得到截频为 的 变换关系:c? c?)(zHL )(ZHd
1
1
1
1
|)()(

Z
zzLd zHZH
2.数字低通 — 数字高通低通变成高通,只需将频率响应旋转 180°,即将 Z变换成 -Z
即可 —— 旋转变换
1
1
11
1)(?


Z
ZZGz
将下式中的 用 代替,1?Z 1Z
得到低通到高通的变换



1
1
1
1 Z
Zz




2
c o s
2
c o s cccc
3.数字低通 — 数字带通数字带通映射到数字低通必须满足的映射关系:
数字带通 数字低通 数字带通 数字低通
10 00 zeZ j
11,0 zZ
从 0变到 时,也应从 0变到,
因而全通函数的阶数应为
2
2?N
则有:
1
1
1
1
11
*11
*)(



Z
Z
Z
ZZGz
11122
2
1
1
2




ZZ
ZZ


变换函数参数的确定列于表 5-5-4
4.数字低通 — 数字带阻数字带阻映射到数字低通必须满足的映射关系:
100 zeZ j
110,0 zZ
从 0变到 时,也从 0变到,因而全通函数阶数为2 2?N
时,当 11Z 1)1(1 Gz,得变换关系:
1
1
1
1
11
1
*
*1)(?



Z
Z
Z
ZZGz
11122
2
1
1
2




ZZ
ZZ


变换函数参数的确定列于表 5-5-4
截止频率为 的低通数字滤波器变换成各型数字滤波器
c?
变换公式见表 5-5-4,P
5.6 数字陷波器设计也被称为点阻滤波器
一种特殊的带阻滤波器,其阻带在理想情况下只有一个频率点主要用于消除某个特殊的频率干扰
在各种测量仪器和数据采集系统中用于消除电源干扰的工频陷波器一般为 IIR滤波器设计方法
双线性变换进行设计
零极点累试的方法
1、陷波器的频率特性设信号为,为干扰频率,需要将其滤掉x t s t A t( ) ( ) s i n 0?0
则要求滤波器的传输函数为,
| ( )|H j
1
0
0
0
其数字域的传输函数为,
| ( )|H e j
1
0
0
0
/ f Ts其中 数字陷波器的频率响应
2,利用双线性变换设计数字陷波器二阶模拟陷波器传输函数为
2
0
2
2
0
2
)( Bss ssH
其幅度响应为
2222
0
22
0
)(
|)(|


B
jH
3dB带宽为 分别为在陷波频率的上边频和下边频幅度增益为 -3dB的频率。
1212,,B
采用双线性变换,可求得二阶数字陷波器传输函数如下
1
1
1
1|)()(

z
zssHzH 22
0
12
0
2
0
22
0
12
0
2
0
)1()1(2)1(
)1()1(2)1(




zBzB
zz
21
21
)1(1
)1()1(2)1(
2
1)(




zz
zzzH


B
B


2
0
2
0
1
1?
2
0
2
0
1
1


)2/t a n (1
)2/t a n (1

B
B

0c o s
采用双线性变换,可求得二阶数字陷波器传输函数如下可改写成其中:
陷波频率,数字 3dB带宽 与 之间的关系:0? 12B,
例 5-6-1 设计一个二阶数字陷波器,其模拟陷波频率为 60Hz,3dB带宽为 6Hz,采样频率为 400Hz
解,先求数字陷波频率和数字 3dB带宽:
3.0)4 0 0/60(20
03.0)4 0 0/6(2B
0,8 8 1 6 1 8 6)2/t a n (1 )2/t a n (1

B
B
0,5 8 7 7 8 5c o s 0
数字传输函数为
21
21
8 8 1 6 1 8.01 0 5 9 8 7.11
9 4 0 8 0 9.01 0 5 9 8 7.19 4 0 8 0 9.0




zz
zz
21
21
)1(1
)1()1(2)1(
2
1)(




zz
zzzH


3、零极点累试的设计方法
1)、基本方法
IIR滤波器的系统函数可表示为
H z
b z
a z
c z
d z
k
N
k
k
N
k
r
M
r
k
N
k
( )
( )
( )
0
1
1
1
1
1
1
11
1
1
由于零极点的共轭成对,先讨论二阶的情况,
H z z e z ez p z pi
j j
i i
( ) ( )( )( )( )*


1 1
1 1
1 1
1 1
0 0
传输函数在 处具有无限大的衰减,为了保证在 时,
信号无衰减,必须在 附近配置上两个极点,使
0 0
0 | ( )|H z? 1
若一个二阶的系统不能满足要求,可采用多个二阶系统的级联三个二阶系统的 级联二阶系统的零极点 分布,3个零点在,极点分布在以 Δ为半径,以 为圆心的圆弧上,以保证
0?jez?
0?jez? | ( )|H z? 1
极点坐标为,)( 00 ijji eep


3
20
1
1
1
i
i
i
i
其中总的传输函数为
H z H z H z H z( ) ( ) ( ) ( )? 1 2 3
设计的问题就是如何选取 使 时,使即如何选取参数 Δ及
pi z z e j0 0? | ( )|H z? 1
i
3)、参数 Δ及 的选择?1
将 和 代下式:pi?1
H z z e z ez p z pi
j j
i i
( ) ( )( )( )( )*


1 1
1 1
1 1
1 1
0 0
可得:
2**1
21
0
)(1
c o s21)(




zppppz
zzzH
iiii
i
其中,p pi i i* [ c o s c o s ( )]2 0 0
p pi i i* c o s1 22
设实际实现滤波运算时存在的有限字长差为,? 2 0cos?
实际为 2 0cos,实际的零点将从 移到?00
表示由于计算误差而产生的陷波频率的偏移?
从而有:
2 20 0c o s c o s ( )
或者 c o s ( ) c o s
0 0
2

当 很小时,
2 0sin
表明实际的陷波点移动了
2 0sin
为使滤波器仍然对 陷波?0
要求移动远小于 Δ,由此得出了确定 Δ的准则,
|s in2|
0?

同时,为了保证在 时,,又要求 0 | ( )|H z? 1
1
Δ越小,越接近理想的幅度特性,但 Δ太小,计算舍入误差将十分敏感,有时会发散,因而对参数 Δ的选取要根据实际情况反复调整试验
3),滤波过程的实现滤波过程就是将输入信号代入该式进行迭代求解一般不宜采用高阶迭代方程,因为高阶迭代方程计算误差影响将十分严重采用逐级求解二阶的方程的方法,求完第一节后,再将第一节的输出输入给第二节这种方法也适于用双线性变换法求得的传输函数参数 Δ及 选定后,得到下述二阶差分方程?1
y n x n x n x n p p y n p p y ni i i i i i i i i i( ) ( ) c o s ( ) ( ) ( ) ( ) ( )* *2 1 2 1 20?
4)在 点的实际衰减?0
由于有限字长的误差,使传输函数的零极点到发生了偏移,
在 点的实际衰减将变成了有限值,用矢量的方法,可求得在 点的幅度值为
0
0
3
0
|s in2||)(| 0jeH
其用 dB来表示的衰减为
dBeH j |s i n|lg6018|)(|l o g20 00
5.7 IIR数字滤波器的计算机辅助设计脉冲响应不变法和双线性变换法
设计时先设计模拟滤波器,再经过变换得到需要的数字滤波器
只适合于设计幅度特性分段恒定的滤波器,例如低通、高通、带通和带阻等滤波器,不能设计任意幅度特性或多带的滤波器
设计所得的滤波器在性能上也不一定是最优的计算机辅助设计
按照指定的滤波器特性,根据某种最优化准则,利用计算机进行反复迭代运算,从而求出所要求的滤波器
既可在时域特性进行,也可在频域特性进行,
指定滤波器的单位脉冲响应或指定滤波后的输出波形
指定滤波器的频率特性本节内容
介绍一种从频域出发,基于最小均方误差的最优化准则,求滤波器传输函数的计算机辅助设计方法
讨论最小平方逆滤波的问题
从时域出发,按照均方误差最小的准则来求取滤波器的单位脉冲响应
5.7.1 IIR数字滤波器频域最小均方误差设计给定的幅度特性形状是任意的时候,可以采用本节介绍的最小均方误差法进行设计。
最优化准则是使实际设计的滤波器频响与希望设计的滤波器频响在指定的频率点均方误差最小,从而求得传输函数中分子和分母多项式的系数设希望设计滤波器频响为,实际设计滤波器频响为其均方误差 E的表达式为
)(?jeH)(?jd eH
2
1
)()(?

i
j
d
j ii eHeHE
设 IIR滤波器是由 k个二阶网络级联而成,H(Z)可表示为:




k
i ii
ii
zczc
zbzazH
1
21
21
1
1
记:
,?EE
kkkk dcbadcba,,,...,,,1111
令:
A
eHH ij
i
)( )( i
jdd eHH
则均方误差表示成

2
1
,?

N
i
di HHAAE?
首先选择 A使 E最小,令
0,?
A
AE 022
1
2
N
i
dii HHHA
gN
i
i
N
i
di
A
H
HH
A
1
2
1
得推导 对系数的偏导数方程
gAE,?



N
i k
i
digg
k
g HHHAAAE
1
2,

k
i
i
k
i
i
ii
k
i HHHH
HH
H

*
2/1*
2
1||



k
i
ie
i
HHR
H?22
1



k
i
iei
HHRH
1
写成对系数的偏导
ij
i ez
ikik
i
ei
k
i
zbza
zRH
a
H



21
1
1
||
ij
i ez
ikii
i
ei
k
i
zbza
z
RH
b
H



21
2
1
ij
i ez
ikik
i
ei
k
i
zdzc
z
RH
c
H




21
2
1
ij
i ez
ikik
i
ei
k
i
zdza
z
RH
d
H




21
2
1
有时采用所谓最小 p误差准则进行设计,误差函数为:
pidjN
i
j
p iii eHeHeWE



1
5,7,2 最小平方逆滤波
1、问题的提出正滤波系统逆滤波系统
x(n)= s(n)* p(n)
X(z)=S(z)P(z)
逆滤波的目的在于求得一个滤波器 h(n),使
H z P z( ) ( )? 1
H(z)X(z)= H(z)P(z)S(z)=S(z)

h(n)* x(n)= h(n)*[ p(n)* s(n)]
= h(n)* p(n)]* s(n)
=δ(n) s(n)= s(n)
应满足将原滤波器的响应作为输入时其输出为 δ(n)
2、最小平方逆滤波原理通常的情况是 p(n)并不是已知的,仅仅 x(n)是可以测量的。
在地震勘探中,x(n)为地震仪器检测到的信号,p(n)为将地层作为滤波器时,地层对地震激励信号的响应,通常不是 已知的。
希望根据 x(n)= s(n)* p(n)的某些特征来寻求 p(n)的估计值,从而求得 h(n)。 最小平方逆滤波则是从使平方误差最小的准则出发来求 h(n)。
设 h(n)为有限长 { h(n),n=0,1…,N},要求在 p(n)输入时,输出为 δ(n),现在假定 p(n)输入时,其输出为
y n h n p n h m p n m
m
( ) ( ) ( ) ( ) ( )
按最小平方准则有,
为最小,即22)]()([ nny
n

n m
h m p n m n( ) ( ) ( )2 2
2为了求出使 最小的 h(n),
对 求 h(n)的偏导数并令其为零,?2
2 2h j h j h m p n m n
n m( ) ( )
( ) ( ) ( )
0)()()()(2 jnpnmnpmh
mn
或,
n m n
h m p n m p n j n p n j( ) ( ) ( ) ( ) ( )? 0
于是对 j=0,1,…,N 有
m
N
n n
h m p n m p n j n p n j p j

0
( ) ( ) ( ) ( ) ( ) ( )?
令:
r m j p n m p n j r j mp
n
p( ) ( ) ( ) ( )
N,0,1,j),()()(
0

jpjmrmh pN
m
有此方程组为
r r
r r
r N r N
r N
r N
r
h
h
h N
pp p
p p
p p
p
p
p
( ) ( )
( ) ( )
( ) ( )
( )
( )
( )
( )
( )
( )
( )0 1
1 0
1
1
0
0
1
0
0
0


假定 s(n)是不相关的序列,
r j s n s n j E js
n
s( ) ( ) ( ) ( )
x(n)的自相关函数为:
r j x n x n jx
n
( ) ( ) ( )
)()()()( kjnpkpmnsmp
kmn

)()()()( kjnsmnskpmp
nmk

)()()( kjmrkpmp s
mk

)()()( kjmEkpmp s
mk

)()()( kjmmpkpE
mk
s
)()()( jrEkjpkpE ps
k
s
可通过求解下述方程组来得到 h(n)
r r
r r
r N r N
r N
r N
r
h
h
h N
E
x x
x x
x x
x
x
x
s
( ) ( )
( ) ( )
( ) ( )
( )
( )
( )
( )
( )
( )
/0 1
1 0
1
1
0
0
1
1
0
0


上述方程组中矩阵的主对角线全为,由于有,
所以该矩阵为正定矩阵,该方程组一定有唯一的一组解,上述矩阵常称为 Toplize 矩阵。
rx( )0 r r jx x( ) ( )0?
89
89