虽然窗函数设计法设计数字滤波器过程简单,但不能设计具有复杂频率响应特性的数字滤波器。频率采样法直接从频域出发,可方便地设计具有任意频率响应的数字滤波器。
由于希望把预期频率特性的傅里叶反变换作为滤波器系数,故在设定预期频率特性时就必须遵循下列几条基本原则:
①预期频率特性的样本点数应等于滤波器的长度N,并在单位圆上等间隔分布。
②作为复数序列的预期频率特性应具有共轭对称性,以保证其傅里叶反变换所得系数为
实序列,因此其幅频特性应为偶函数,相频特性应为奇函数。
③预期频率特性的相位特性应该与频率呈线性关系;这意味着其幅频特性及其反变换所
得序列应该具有对称或反对称的特点。
频率采样法的设计过程
给定理想低通滤波器Hd(w),先选择滤波器长度N,然后对Hd(w)在0到2π上的N个等间隔频率上采样,得到H(k),而脉冲响应由离散傅里叶反变换式得到。这种反变换在一定意义上也隐含了对样本H(k)进行内插而得到实际响应H(w)的过程。
从Z变换的定义可以得到,h(k)的Z变换可得
这是频率采样法的理论基础,其中(WN)^(-k)为e^(j2πk/N),式(3-10-7)可进一步表示为
实际滤波器的频率响应为
从上式可以看到,在采样点上,实际滤波器的幅值特性与期望特性完全一致,
在非采样点上,实际滤波器的幅值由采样点处的增益通过函数(a)加权后积分得到,非采样点
一般存在波纹,下图所示。
图(a)是预期低通滤波器的频率特性示意图,只给出正频率部分,经过频率采样得到图(b),可以看到非采样点存在波纹。
理想滤波器的过渡带宽为零,边界频率(通带和阻带交界点)附近的幅值特性会发生突变。如果直接采样,则边界处的幅值突变会造成实际滤波器在该频点出现波纹(如图(b)所示)。为了减少波纹,可适当增加边界频率处的过渡带宽度,减少通带和阻带转换时的非线性。下图是对边界频率处的增加1个或2个样本点使幅值响应进行修正示意图。与直接频率采样法相比,修正边界频率采样值后,减少了通带的波纹,阻带衰减也得到了改善,付出的代价是过渡带增宽。
频率采样法和窗函数法的结合
频率采样法可设计任意频率响应的数字滤波器,而窗函数法可通过选择窗函数调节边界频率处的滤波特性(如过渡带宽度、阻带衰减等)。可以结合两种方法的优点,设计出符合要求的数字滤波器。
实际应用中,可首先采用频率采样法设计一个阶次足够大(例如N=512)的线性相位数字滤波器h(n)。然后,根据滤波器的类型,对过渡带、阻带衰减等的要求,选择特定窗函数w(n)对h(n)进行加窗处理,修正频率采样法滤波器的边界特性,从而满足设计要求。
案例、在信号处理中对信号要乘以一个阶梯形的权函数,信号采样频率为2000Hz,权函数的要求是150~250Hz权函数值为1,350~500Hz权函数值为0.5,600~700Hz权函数值为0.25。用fir2函数设计一个权函数的滤波器满足这些要求。程序如下:
clear all; clc; close all;
Fs=2000; % 采样频率
Fs2=Fs/2; % 奈奎斯特频率
fc1=150; fc2=250; fc3=350; % 各频点的频率值
fc4=500; fc5=600; fc6=700;
fd=[0 fc1 fc2 fc3 fc4 fc5 fc6 800 Fs2]/Fs2; % 各频率点构成频率矢量
Hd=[0 1 1 0.5 0.5 0.25 0.25 0 0]; % 对应各频点的理想幅值
dw=(fc3-fc2)*pi/Fs2; % 求出过渡带宽
N=ceil(11*pi/dw); % 计算滤波器阶数
wind=blackman(N+1)'; % 布莱克曼窗函数
hn=fir2(N, fd, Hd, wind); % 用fir2函数义滤波器系数
[H, f] = freqz(hn, 1, 512, Fs); % 求滤波器的响应
% 作图
plot(f, abs(H),'k','linewidth',2),
xlabel('频率/Hz'); ylabel('幅值');
title('一个频域阶梯形权函数的幅值响应曲线')
set(gca,'XTickMode','manual','XTick',[0,150,250,350,500,600,700,800])
set(gca,'YTickMode','manual','YTick',[0.25,0.5,1])
grid on; ylim([0 1.2]);
set(gcf,'color','w');
运行结果如下:
参考文献:MATLAB数字信号处理85个实用案例精讲——入门到进阶;宋知用(编著)