一种抗混叠失真的IIR数字滤波器设计方法

刘文斐

(齐鲁师范学院 物理与电子工程学院,山东 济南 250200)

滤波器,顾名思义,是用于“滤波”的,即将某些频率范围的信号进行保留,抑制其他频率的信号．按照处理信号种类的不同，滤波器分为数字滤波器和模拟滤波器两大类,其中,数字滤波器是通过程序运算将数据进行处理[1],对数字信号进行滤波处理得到期望响应特性的离散时间系统．实际上,数字滤波器不仅可以处理数字信号,如果在滤波器的输入端接入模-数转换器,它也可对模拟信号进行处理．与模拟滤波器相比,数字滤波器具有高精度、高可靠性、可程控改变特性或复用、便于集成等优点,广泛应用于语音、图像、医学生物信号处理等领域．数字滤波器按照其系统中冲激响应的特点可分为无限长单位冲激响应(infinite impulse response，IIR)滤波器和有限长单位冲激响应(finite impulse response，FIR)滤波器．IIR滤波器具有结构简单、效率高、与模拟滤波器有对应关系、易于解析控制及计算机辅助设计等优点[2],与FIR滤波器相比,它的相同阶次滤波效果更好．

1 IIR数字滤波器的设计方法

IIR数字滤波器的设计方法[3]有间接法和直接法,一般采用间接法进行设计．间接法是通过模拟滤波器设计数字滤波器,这是因为模拟滤波器的设计方法较为成熟,有完整的设计公式以及数据表格可以直接使用,可使数字滤波器的设计更为简便快捷．直接法设计IIR数字滤波器,如果给定的设计要求比较复杂,需通过计算机辅助设计求解线性或非线性的联立方程组,计算量大．

间接法设计IIR数字滤波器有3种设计方案[3].第1种方案分3步:首先通过给定的技术指标设计归一化样本的模拟低通滤波器；然后通过模拟频带变换成其他模拟各种(低通、高通、带通、带阻)滤波器；最后将模拟滤波器数字化成相应的数字滤波器.该方法的实现需要经过“模拟—模拟频带变换”．第2种方案是直接由设计出的归一化样本模拟低通滤波器变换成各种通带的数字滤波器,再将频带变换与数字化结合在一起,一步变换完成.该方法虽然简化了运算步骤但实质也是经过“模拟—模拟频带变换”．第3种方案也分3步:首先通过给定的技术指标设计归一化样本的模拟低通滤波器；然后直接数字化成数字低通滤波器；最后用数字频带变换成其他各种通带的数字滤波器.该方法的实现需经过“数字—数字频带变换”．这3种设计方案的根本区别在于频带变换的信号领域不同,前2种是在模拟信号域中进行,且对模拟滤波器的运算设计有可直接使用的公式和表格,所以是较常用的方法;第3种是在数字信号域中进行频带变换,需找到不同数字频带间的函数关系．

2 冲激响应不变法产生混叠失真的解决方案

在间接法设计数字滤波器的3种设计方案中,无论用哪一种对数字滤波器进行设计,都有一个共同的关键步骤——模拟滤波器的数字化过程[3],也就是将模拟滤波器的系统函数H(s)变换成所需的数字滤波器的系统函数H(z),即把s平面映射到z平面上．数字化方法可采用冲激响应不变法、阶跃响应不变法和双线性变换法,其中较常用的是冲激响应不变法和双线性变换法．

冲激响应不变法是从滤波器的单位抽样响应出发,使数字滤波器的单位抽样响应h(n)逼近模拟滤波器的单位抽样响应h(t),频率间的变化是线性变换关系,克服了双线性变换法中非线性频率变换带来的临界频率点的频率畸变,所以是最简单、最直接的把s平面映射到z平面的映射方法．但是该方法要求模拟滤波器是严格带限于抽样频率的1/2,如果不满足该要求,数字滤波器的频率响应将产生混叠失真．

混叠失真是指在对连续信号进行等间隔采样时[4],如果不能满足奈奎斯特采样定理,即信号的最高频率成分大于采样频率的1/2,那么采样后信号的频谱会发生频谱分量的交叠现象,高于和低于采样频率的信号混杂在一起.如果发生频谱混叠失真,原信号将无法从采样信号中恢复还原．

如果通过“模拟—模拟频带变换”方法设计IIR数字高通或者带阻滤波器,冲激响应法确实会产生混叠失真现象,但如果通过“数字—数字频带变换”方法则不存在该问题,因为该方法的数字化过程是将模拟低通滤波器的系统函数映射为数字低通滤波器的系统函数,而模拟低通滤波器是严格带限于抽样频率的1/2、是抗混叠的,所以不会出现频率混叠失真现象．

3 抗混叠滤波器的设计

3.1 IIR数字高通滤波器的设计

设计一个切贝雪夫Ⅰ型数字高通滤波器,要求通带的截止频率为4.5 kHz,阻带上限截止频率为3 kHz,通带内衰减不大于2 dB,阻带内衰减不小于15 dB,抽样频率为20 kHz．

如果用“模拟—模拟频带变换”进行设计,其设计流程如图1所示.用MATLAB辅助设计、仿真设计结果,主要的设计步骤如图2所示.

图1 用“模拟—模拟频带变换”进行设计的流程

图2 用MATLAB辅助设计的仿真结果

用“模拟—模拟频带变换”方法、通过冲激响应不变法完成模拟滤波器的数字化过程设计的高通滤波器,确实存在频谱混叠失真,不符合设计参数,达不到设计要求(图3)．为消除频谱失真，用“数字—数字频带变换”方法进行设计(图4).

图3 发生频谱混叠失真的数字高通滤波器频谱

由上述可知,“数字—数字频带变换”的实质就是从数字低通滤波器的Z平面映射到另一个待求所需类型数字滤波器的z平面的变化关系,关键点是找到Z到z的映射函数[3]Z-1=G(z-1),则有

设Z=ejθ,z=ejω,则数字低通映射到数字高通的映射函数为

若低通滤波器的系统函数Hl(z)的分子分母皆为Z的多项式,即

Hl(Z)=

由此可得

用该表达式的形式在MATLAB辅助设计中自定义映射函数,将数字低通滤波器系统函数Hl(z)映射为数字高通系统函数H(z)的形式(图4),用MATLAB辅助设计的设计步骤如图5所示.

图4 用“数字—数字频带变换”进行设计

图5 用MATLAB辅助设计的仿真结果

关键程序:

%自定义mapping函数:

function[bz,az]=mapping(bZ,aZ,Nz,Dz)

bzord=(length(bZ)-1); %原系统函数分子阶数

bzord1=(length(bZ)-1)*(length(Nz)-1); %变换后系统函数分子阶数

azord=(length(aZ)-1); %原系统函数分母阶数

azord1=(length(aZ)-1)*(length(Dz)-1); %变换后系统函数分母阶数

bz=zeros(1,bzord1+1); %分子系数向量初始化为零

for k=0:bzord %求各多项式乘积结果

pln=[1];

for l=0:k-1

pln=conv(pln,Nz); %求N(z)的k次乘积

end

pld=[1];

for l=0:bzord-k-1

pld=conv(pld,Dz); %求N(z)的bzord-k次乘积

end

bZZ=bZ(k+1)*conv(pln,pld);

bz=bz+bZ(k+1)*conv(pln,pld); %分子多项式系数向量求和

end

az=zeros(1,azord1+1)

for k=0:azord

pln=[1];

for l=0:k-1

pln=conv(pln,Nz);

end

pld=[1];

for l=0:azord-k-1

pld=conv(pld,Dz);

end

aZZ=aZ(k+1)*conv(pln,pld);

az=az+aZ(k+1)*conv(pln,pld); %分母多项式系数向量求和

end

az1=az(1);

az=az/az1;

bz=bz/az1;

%主程序

Fs=20000; %设置采样频率

rp=2; %设置通带衰减

rs=15; %设置阻带衰减

wp=2*pi*3000/Fs; %设置阻带截止频率

ws=2*pi*4500/Fs; %设置通带截止频率

OmegaP=2*Fs*tan(wp/2); %计算模拟低通通带截止频率

OmegaS=2*Fs*tan(ws/2); %计算模拟低通阻带截止频率

[N,OmegaC]=cheb1ord(wp,ws,rp,rs,'s'); %得到切贝雪夫Ⅰ型模拟低通滤波器阶数和通带截止频率

[b,a]=cheby1(N,rp,OmegaC,'s'); %得到模拟低通滤波器系统函数

[bZ,aZ]=impinvar(b,a,Fs); %冲激响应不变法得到数字低通滤波器系统函数

[bhp,ahp]=mapping(bZ,aZ,Nz,Dz); %调用mapping函数得到数字高通滤波器系统函数

设计结果如图6所示．由图6可知,该数字高通滤波器通带的截止频率约为4.5 kHz，通带内衰减小于5 dB,阻带的截止频率为3 kHz，阻带衰减大于15 dB,设计指标基本达到设计要求．所以用“数字—数字频带变换”能够解决冲激响应不变法带来的频谱混叠失真问题,实现高通滤波器的设计．

图6 消除频谱混叠失真的数字高通滤波器频谱

3.2 IIR数字带阻滤波器的设计

设计一个切贝雪夫Ⅰ型数字带阻滤波器,其技术指标:通带的截止频率分别为6,14 kHz,阻带截止频率分别为8,12 kHz,通带内衰减不大于2 dB,阻带内衰减不小于20 dB,抽样频率为400 kHz．

用“模拟—模拟频带变换”方法、通过冲激响应不变法,完成模拟滤波器的数字化过程设计带阻滤波器,频谱混叠现象严重,无法达到设计要求．如果用“数字—数字频带变换”方法,则数字低通映射为数字带阻的映射函数为[3]

θp为数字低通滤波器通带的截止频率;ωp1，ωp2为数字带阻滤波器通带的截止频率．

用MATLAB辅助设计,其设计步骤与数字高通滤波器的设计相似:

1)确定模拟低通滤波器的阶次及截止频率．

%设置带阻滤波器参数

Fs=40000; %设置采样频率

S1=2*pi*8000/Fs; %设置阻带截止频率

S2=2*pi*12000/Fs;

P1=2*pi*6000/Fs; %设置阻带截止频率

P2=2*pi*14000/Fs;

rp=2; %设置通带衰减

rs=20; %设置阻带衰减

[N,Wc]=cheb1ord(P1,S1,rp,rs,'s'); %计算模拟低通滤波器的最小阶数及通带截止频率

运行得到模拟低通滤波器阶数(N=6),通带截止频率为7 kHz.

2)创建模拟低通滤波器原型.

[b,a]=cheby1(N,rp,Wc,'s'); %由阶数和通带截止频率得到模拟低通滤波器系统函数

3)冲激响应不变法将模拟低通转化为数字低通滤波器.

[bZ,aZ]=impinvar(b,a,Fs); %得到数字低通滤波器系统函数

4)调用函数mapping()实现数字带阻滤波器系统函数.

alpha=cos(p0); %设置参数α

klpha=(tan((P2-P1)/2))*(tan(Wc/2)); %设置参数k

d1=-2*alpha/(1+klpha);

d2=(1-klpha)/(1+klpha);

Nz=[d2,d1,1];

Dz=[1,d1,d2];

[bhp,ahp]=mapping(bZ,aZ,Nz,Dz);

设计结果如图7所示．由图7可知,频率小于6 kHz和大于14 kHz的范围属于通带范围,通带衰减小于5 dB，接近指标值;频率为8～12 kHz,属于阻带,衰减大于20 dB,基本满足设计要求．

图7 消除频谱混叠失真的数字带阻滤波器频谱

将设计实现的带阻滤波器应用于杂波滤波,验证其滤波效果．

实现程序:

%引入杂波信号

t=(1:Fs)/Fs;

y=sin(2*pi*4000*t)+sin(2*pi*7000*t)+

sin(2*pi*10000*t)+

sin(2*pi*13000*t)+

sin(2*pi*15000*t);

%对杂波信号滤波,其中bZ,aZ为设计实现带阻滤波器系统函数分式中分子、分母的系数

x=filter(bZ,aZ,y).

杂波波形如图8所示.由图8可知,杂波信号中包含5个频谱分量,频率中心分别为4,9,10,13,15 kHz．进入带阻滤波器进行滤波,由于带阻滤波器的通带范围为小于6 kHz和高于14 kHz,阻带为8～12 kHz．所以经过滤波后,杂波信号中10 kHz频率处的频谱分量被滤掉,9,13 kHz频率的信号由于处于带阻滤波器的过渡带范围内,频谱被衰减,处于通带范围内的信号没有受到影响.滤波后信号如图9所示.由图8～9可知,通过“数字—数字频带变换”进行设计实现的滤波器能够达到设计要求,可实现滤波的功能．

图8 杂波信号

图9 滤波后信号

4 结论

提出用“数字—数字频带变换”设计数字滤波器的方法.该方法在保证实现设计指标的同时,消除由冲激响应不变法在将模拟滤波器系统函数映射为数字滤波器系统函数时带来的频谱混叠失真问题．通过数字高通和数字带阻滤波器的设计实例,验证该方法的有效性．用“数字—数字频带变换”方法设计数字滤波器,可满足设计要求,是一种有效的抗混叠的设计方法．