的IDTFT导出hd(n) hd(n)?12?????Hd(ejw)ejwndw (14)
由于Hd(ejw)是矩形频率特性,故hd(n)一定是无限长的序列,且是非因果的。然而FIR滤波器是有限长的,所以用有限长的h(n)来逼近无限长的hd(n),最简单的方法是截取hd(n)中最重要的一段,将无限长的hd(n)截取成长度为M的有限长序列,等效于再hd(n)上施加了一个长度为M的矩形窗口,更为一般的,可以用一个长度为M的窗口函数w(n)来截取hd(n),即
h(n)?w(n)hd(n) (15)
这一方法通常称为窗函数法,窗口函数的形状及长度M的选择是窗函数法的关键。
下面我们一低通为例,了解一下窗函数法的运用: ①提出希望频率响应函数(低通)
图3 理想低通滤波器的频响 线性相位, 具有片断特点, 即
?j????eHd(e)????0j?|Hd(ej?)|1???0.25?O0.25???|?|??c
?c?|?|??②算出
hd(n)???12????-c?H(ej?)ej?nd?
12????ce?j??ej?nd?
sin(?c(n??))(无限长)
?(n??)
0.310.20.10.50 图4 理想低通的单位脉冲响应(无限长的一部分)
③加窗w(n),长N, 得
h(n)?hd(n)w(n) (*)
要线性相位, 就要h(n)关于(N?1)/2偶对称,而hd(n)关于?偶对称, 故要求
??(N?1)/2
所以要求w(n)关于??(N?1)/2偶对称.
0.310.20.10.5000102030-0.10102030
图5 窗函数 图6 加窗后的单位脉冲响应
再回过来检验H(ej?)是否满足精度要求.
图7 图4的脉冲响应的频响 图8 理想频响与实际频响的对比
若基本满足, 则依截取的h(n), 制硬件, 编软件.
为便于选择使用, 将5种常见的窗函数基本参数如表1所示。
表1 5种常见的窗函数基本参数 窗函数过渡带加窗后滤波器的 0.501|Hd(ej?)|H(ej?)?012310.51012O???0.25?0.25?0??3类型 rectwin bartlet三角 hanning hamming blackman 的 旁瓣峰? ?13 n宽度?B 阻带最小衰减? s4?/N 8?/N 8?/N 8?/N 12?/N ?21 ?25 ?44 ?53 ?74 ?25 ?31 ?41 ?57
(2)频率取样法
窗口设计法事从时域出发,把理想的hd(n)用一定形状的窗口函数截取成有限长的h(n),以此h(n)来近似理想的hd(n),从而频率响应H(ejw)也近似于理想的频率响应Hd(ejw)。我们知道一个有限长序列可以通过其频谱的相同长度的等间隔采样值准确地恢复原有的序列。频率采样法便是从频域出发,对理想的频率响应Hd(ejw)加以等间隔采样
Hd(ejw)|2??Hd(k)
w?Mk(16)
然后以此Hd(k)作为实际FIR滤波器的频率特性的离散样本H(k),即
H(k)?Hd(k)?Hd(ejw)|w?2?kM,k?0,1,2,...M?1 (17)
由H(k)通过IDFT可求出有限长序列h(n)为
1h(n)?MM?1k?0?H(k)ej2?nkM,n?0,1,...M?1 (18)
利用M个频率的离散样本H(k)同样可求出FIR滤波器的系统函数H(z)及频率响应H(e)。H(z)??h(n)zjwM?1n?0?n1?z?M?MM?1?j2?H(k)M ,其中W?e??k?11?Wzk?0(19)
令z?ejw可得到滤波器的频率响应H(ejw)。如果设计的是线性相位的FIR数字滤波器,其采样值H(k)的相位的幅度一定要满足特定的约束
条件,这个设计时一定要注意。
(3)最优化设计法
最优化设计法事以最佳一致逼近(最大误差最小化)理论为基础,利用雷米兹算法设计的具有等波纹特性的设计方法。具体设计步骤如下:
①对设计指标进行归一化处理。
②确定remezord函数所需要的参数。包括归一化边界频率、各频带的幅度要求和波纹要求等。归一化边界频率总是从0开始到1结束,故只需递增列出中间的边界频率;频带幅度要求不含过渡区,个数是边界频率个数的一半加1;波纹要求是频带内幅度允许的波动要求,与分贝间的关系是:
Rp??20log10(1??1?),Rs??20log10(2)??20lo10g(?2) (20) 1??11??1 ③利用remezord函数确定remez所需参数。
④调用remez函数进行设计。 ⑤利用freqz函数验算技术指标是否满足要求。 3 数字滤波器的MATLAB设计 3.1 直接程序设计法
(1)IIR的直接程序设计法
例如欲设计一数字(IIR)带阻滤波器,其数字域指标为:数字阻带边缘频率分别为0.4?和0.7?,数字通带边缘频率为0.25?和0.8?,通带波动为1db最小阻带衰减为40db。 此题的MATLAB程序为:
ws=[0.4*pi 0.7*pi]; %数字阻带边缘频率 wp=[0.25*pi 0.8*pi]; %数字通带边缘频率 rp=1 ; %通带波动(db) as=40; %阻带衰减 [n,wn]=cheb2ord(wp/pi,ws/pi,rp,as);
%根据给定指标,确定滤波器的阶数N和频率缩放因子Wn [b,a]=cheby2(n,as,ws/pi,'stop');%返回的b,a分别为H(z)的分子、
分母。
[h,w]=freqz(b,a,512);%返回的h,w分别为滤波器的频率响应及
其频率
plot(w/pi,abs(h));%画出频率响应(以w/pi为横轴) grid;
xlabel('w/pi');
ylabel('幅值'); title('频率响应'); 程序运行结果为:
图9 所设计的带阻滤波器的频率响应
在设计中如果该滤波器的特性不满足要求,原有的参数必须做相应的调整,在程序中只需对参数做新的设定就可以得到所需的滤波器。接下来我们来看看此题所设计的滤波器的滤波效果:S为含有3个频率成分的信号(归一化频率(w/2?)分别为0.1、0.3、0.45),用所设计的滤波器滤除归一化频率为0.3的成分。
n=0:100;
s1=sin(pi*0.2*n); s2=sin(pi*0.6*n); s3=sin(pi*0.9*n); s4=s1+s3; s=s1+s2+s3; sf=filter(b,a,s); subplot(311)
stem(n,s);title('滤波前的信号'); subplot(312);