|
FIR滤波器设计-窗函数法（2） 来自信号与系统和数字信号处理 00:00 29:30 本文是FIR滤波器的第二个问题——窗函数法的第二篇。我们介绍各种窗函数。 二 窗函数法 1. 设计原理 上一篇以低通滤波器加矩形窗为例，详细介绍了窗函数法设计FIR滤波器的原理，辅以动画演示，形象生动。链接如下： 数字信号处理系列串讲第17篇（数字滤波器之二）——FIR滤波器（2）：窗函数法（1） 2. 各种窗函数 下面，给出几种常见窗函数的时域和幅度函数的表达式。由于其时域都满足关于(N-1)/2 偶对称，都具有线性相位特性。下面重点分析其时域表达式和频域的幅度函数。 （1）矩形窗 大家可能会有疑问：矩形窗的旁瓣与主瓣幅度之比，真的与长度N无关吗？我们详细分析一下： 矩形窗的幅度函数：sin(Nw/2)/sin(w/2)，将w=0代入得到主瓣幅度峰值为N，将w=3Π/N代入得到第一旁瓣（最高的旁瓣）峰值为1/sin(3Π/2N)，所以旁瓣与主瓣峰值之比为：1/(N×sin(3Π/2N))。我们分别以该数值以及该数值取20log10为纵轴，N为横轴，画出图形如下： 可见，当N大于30时，矩形窗谱的旁瓣与主瓣幅度之比基本保持为常数0.21（-13.5dB）。 （2）三角窗 既然矩形窗的旁瓣和主瓣幅度之比基本是一固定值，我们就想到，如果把它平方一下，这一比值会减小，0.21的平方为0.044（-27dB）。 而根据傅里叶变换的性质，频域相乘，对应时域卷积，两个宽度相同的矩形脉冲卷积为三角形。这就是三角窗。 与矩形窗相比，三角窗的旁瓣虽然显著降低了，但不幸的是，主瓣宽度也增加了一倍，为8Π/N。 （3）汉宁（Hanning）窗（又称为升余弦窗） 汉宁窗的设计思想是，将矩形窗的窗谱分别左移和右移半个主瓣宽度，使其旁瓣相互抵消，如下图所示。 蓝色实线为矩形窗谱（幅度函数），红色虚线和绿色虚线分别为其右移和左移半个主瓣，用红色的第一旁瓣和绿色的第三旁瓣去抵消蓝色的第一旁瓣，以此类推，达到降低旁瓣幅度的目的。 抵消的结果是令人惊喜的，旁瓣幅度可以降至-31dB。比三角窗还要低。而其主瓣宽度与三角窗相同，为8Π/N。 由于其时域表达式为sin的平方，所以称为升余弦窗。 （3）海明（Hamming）窗（又称为改进的升余弦窗） 观察前面的汉宁窗我们发现，旁瓣幅度下降的很快，也就是说用这种窗设计出来的滤波器，如果在阻带边沿处刚好满足性能指标，则在阻带远离边沿处，会远远高于性能要求。 而下面说的这种窗函数——海明窗则具有比较均匀的旁瓣特性。 海明窗的表达式与汉宁窗类似，只是系数有变化。所以又称为”改进的升余弦窗“。别看这一小小的调整，与汉宁窗相比，海明窗的第一旁瓣幅度降低了10个dB。其旁瓣幅度比较均匀（我们称这种特性为”等波纹性“）。 也就是说，分别用相同长度的海明窗与汉宁窗设计出来的FIR滤波器，在阻带边沿处，海明窗设计的滤波器衰减更大，但在阻带内部远离边沿处，反倒不如汉宁窗设计的滤波器衰减大。 （4）布莱克曼（Blackman）窗（又称为二阶升余弦窗） 依然沿用”左移右移，旁瓣相互抵消“的思想，用更多的旁瓣（左/右移2Π/N和4Π/N）去抵消，就是布莱克曼窗。 其旁瓣衰减可以达到57dB，但同时主瓣宽度也增加了，为12Π/N。 以上几种窗，窗的形状都是固定的，我们能调整的只有一个参数：窗的长度N。 （5）凯赛（Kaiser）窗 凯赛窗的表达式如下（不要被吓倒，不需要记，考试不会考的 大家不需要记住这些公式（包括前面的几种窗函数，也不需要记住其表达式，只需要知道其特点即可）。 大家只要知道，凯赛窗是一种形状可变的窗，窗的旁瓣衰减和主瓣宽度可调，通过参数β来调整。 β=0时，分子分母都为1，即为矩形窗。β越大，旁瓣衰减越大，但主瓣宽度也相应地增大，典型值为4~7。 下图是相同长度（N=21），不同 β 值的凯塞窗的时域波形图（包络）和幅度谱。 β 值如何确定有经验公式可以直接使用，这里就不摆出来吓唬大家了。 （6）几种窗函数的性能比较 下表为几种窗函数的性能比较。第二、三列为窗函数频谱的两个主要指标（旁瓣峰值和主瓣宽度），第四~六列为使用该窗函数得到的FIR滤波器的性能指标（过渡带宽、阻带最小衰减和通带最大衰减）。
加窗后FIR滤波器的过渡带宽取决于窗谱的主瓣宽度，但二者数值不同；同样，滤波器的阻带最小衰减取决于窗谱的旁瓣峰值，但二者数值也不同。 另外，我们也注意到，通带边沿衰减数值一般很小（最大的也只有0.815dB，因为20*log10(1-0.0895)≈-0.815dB），所以，用窗函数法设计FIR滤波器，我们一般只关注过渡带宽和阻带最小衰减这两个指标。 上面这张表（主要是第4第5列），就是我们用窗函数法设计FIR滤波器的依据。 用矩形窗设计的FIR滤波器的过渡带最窄，但阻带特性最差，而用布莱克曼窗，阻带特性最好，但过渡带最宽。因为阻带性能太差，矩形窗和三角窗在实际中应用很少。海明窗和汉宁窗在实际中应用较多。凯赛窗较为灵活，但计算β值较为复杂。 3. 设计步骤及设计实例 首先归纳总结一下用窗函数法设计FIR滤波器的步骤，然后结合具体实例来说明。 设计步骤如下图：
【解】
Matlab代码如下：
画图语句省略...... 其中，自编函数Ideal_lp代码如下：
上例中，为了展示窗函数法的原理，采用了自己编写理想低通滤波器的单位冲激响应hd(n)的方式，实际上，Matlab中提供了函数直接调用：
使用fir1函数时，需要注意三点： 第一，第一个参数N是指FIR滤波器的阶数，对于M点长（0≤n≤M-1）的FIR滤波器，阶数为M-1阶。 第二，设计高通、带阻滤波器时，长度M必须为奇数，即阶数M-1必须为偶数，若用户误设阶数为奇数，该函数会自动加1。如果此处感到疑惑，请复习以下链接： 数字信号处理系列串讲第16篇（数字滤波器之二）——FIR滤波器（1）：线性相位FIR滤波器的特点 第三，该函数中的第二个参数wc的范围是0~1，也就是“数字域频率/Π'。 另外，可能有同学会问，wc为何叫”6dB截止频率“？结合下面这个图，以低通滤波器为例解释一下：
上图中，绿色实线为理想低通滤波器的幅度函数，蓝色实线为加窗后FIR滤波器的幅度函数，wc是理想低通滤波器的截止频率，即FIR滤波器的通带和阻带截止频率的中点（也就是肩峰和谷底的中间）。 当w=wc时，FIR滤波器的幅度函数取值为1/2，转换为dB就是20*log10(1/2)≈-6dB，所以，wc称为6dB截止频率。 第四，如果指定窗函数，其长度必须等于N+1（N为第一个参数）。 下面再给一道例题，用Matlab函数fir1进行设计。
首先，需要将模拟频率转换为数字域频率。二者的关系为：数字域频率=模拟频率*2Π/采样频率。 故转换为数字域频率的带阻滤波器边沿频率分别为：
以上是分析结果，下面就是Matlab程序。
运行结果如下图，对图形局部放大，验证，在通带边沿和阻带边沿处，满足性能指标的要求。
上例中，计算出来的海明窗长度M≥41.25，也就是M取大于等于42的整数，即可满足过渡带宽的要求。但是需要注意的是，我们设计的是带阻滤波器，M必须为奇数。 4. 算法特点 窗函数法具有以下特点： （1）设计原理简单，有闭合公式可循； （2）是从时域出发的一种设计方法，需要求hd(n)，不直接； （3）不能准确控制通带及阻带的截止频率； 在前面的例2中，大家可以尝试以下，将长度M修改为43，运行程序，得到的FIR滤波器的幅度函数如下图：
我们对阻带部分进行局部放大，如下图：
发现阻带中有些频率范围内衰减没有达到50dB，不能满足性能指标要求。 这是窗函数法设计FIR滤波器的一个经常会出现的现象：需要对设计结果进行验证，有时可能需要返工。 （4）通带阻带最大波纹相等，不能分别控制。 这一点也很好理解，因为不管是通带还是阻带，其波动都是由于窗谱与理想低通幅度函数卷积时，窗谱旁瓣的起伏引起的。因此通带和阻带的波动情况是相同的。 窗函数法的内容到此为止。 下一篇预告：FIR滤波器设计的另外一种方法：频率抽样法。 |
|
|
）：
