东南大学-数字信号处理-吴镇扬-4.ppt

第四章有限长单位脉冲响应 FIR 滤波器的设计方法 序言 4 1线性相位FIR数字滤波器的特性 4 2窗口设计法 时间窗口法 4 3频率采样法 4 4FIR数字滤波器的最优化设计 4 5IIR与FIR数字滤器的比较 序言 FIR数字滤波器的差分方程描述 对应的系统函数 因为它是一种线性时不变系统 可用卷积和形式表示 比较 得 FIR数字滤波器的特点 与IIR数字滤波器比较 优点 1 很容易获得严格的线性相位 避免被处理的信号产生相位失真 这一特点在宽频带信号处理 阵列信号处理 数据传输等系统中非常重要 2 可得到多带幅频特性 3 极点全部在原点 永远稳定 无稳定性问题 4 任何一个非因果的有限长序列 总可以通过一定的延时 转变为因果序列 所以因果性总是满足 5 无反馈运算 运算误差小 缺点 1 因为无极点 要获得好的过渡带特性 需以较高的阶数为代价 2 无法利用模拟滤波器的设计结果 一般无解析设计公式 要借助计算机辅助设计程序完成 4 1线性相位FIR数字滤波器的特性 4 1 1线性相位的条件线性相位意味着一个系统的相频特性是频率的线性函数 即 式中 为常数 此时通过这一系统的各频率分量的时延为一相同的常数 系统的群时延为 FIR滤波器的DTFT为 式中H 是正或负的实函数 等式中间和等式右边的实部与虚部应当各自相等 同样实部与虚部的比值应当相等 将上式两边交叉相乘 再将等式右边各项移到左边 应用三角函数的恒等关系 满足上式的条件是 另外一种情况是 除了上述的线性相位外 还有一附加的相位 即 利用类似的关系 可以得出新的解答为 偶对称 奇对称 图1线性相位特性 分四种情况 4 1 2线性相位FIR滤波器的幅度特性 分四种情况1 h n 偶对称 N为奇数h n h N 1 n 4 1 2线性相位FIR滤波器的幅度特性 令 则 令 则 由于偶对称 因此对这些频率也呈偶对称 2 h n 偶对称 N为偶数h n h N 1 n 令 则 或写为 由于奇对称 所以对也为奇对称 且由于时 处必有一零点 因此这种情况不能用于设计时的滤波器 如高通 带阻滤波器 3 h n 奇对称 N为奇数 h n h N 1 n 令n m N 1 2 得 所以 由于点呈奇对称 所以对这些点也奇对称 由于时 相当于H z 在处有两个零点 不能用于的滤波器设计 故不能用作低通 高通和带阻滤波器的设计 4 h n 奇对称 N为偶数 令 由于在 0 处为零 所以H 在 0 2 处为零 即H z 在z 1上有零点 并对 0 2 呈奇对称 四种线性相位FIR滤波器 四种线性相位FIRDF特性 参考表4 1第一种情况 偶 奇 四种滤波器都可设计 第二种情况 偶 偶 可设计低 带通滤波器 不能设计高通和带阻 第三种情况 奇 奇 只能设计带通滤波器 其它滤波器都不能设计 第四种情况 奇 偶 可设计高通 带通滤波器 不能设计低通和带阻 例1N 5 h 0 h 1 h 3 h 4 1 2 h 2 2 求幅度函数H 解 为奇数并且h n 满足偶对称关系a 0 h 2 2 a 1 2h 3 1 a 2 2h 4 1 H 2 cos cos2 2 cos cos2 小结 四种FIR数字滤波器的相位特性只取决于h n 的对称性 而与h n 的值无关 幅度特性取决于h n 设计FIR数字滤波器时 在保证h n 对称的条件下 只要完成幅度特性的逼近即可 注意 当H 用 H 表示时 当H 为奇对称时 其相频特性中还应加一个固定相移 4 1 3线性相位FIR滤波器的零点特性 由该式可看出 若z zi是H z 的零点 则z z 1i也一定是H z 的零点 由于h n 是实数 H z 的零点还必须共轭成对 所以z z i及z 1 z 也必是零点 所以线性相位滤波器的零点必须是互为倒数的共轭对 即成四出现 这种共轭对共有四种可能的情况 既不在单位园上 也不在实轴上 有四个互为倒数的两组共轭对 ziz i1 zi1 z i图4 2 a 在单位圆上 但不在实轴上 因倒数就是自己的共轭 所以有一对共轭零点 zi z i图4 2 b 不在单位圆上 但在实轴上 是实数 共轭就是自己 所以有一对互为倒数的零点 zi 1 zi图4 2 c 又在单位圆上 又在实轴上 共轭和倒数都合为一点 所以成单出现 只有两种可能 zi 1或zi 1图4 2 d p92我们从幅度响应的讨论中已经知道 对于第二种FIR滤波器 h n 偶对称 N为偶数 即是的零点 既在单位圆 又在实轴 所以 必有单根 同样道理 对于第三种 FIR滤波器 h n 奇对称 N为奇数 因所以z 1 z 1都是H z 的单根 对于第四种滤波器 h n 奇对称 N为偶数 H O 0 所以z 1是H z 的单根 所以 h n 奇对称 H 0 0N为偶数 H 0线性相位滤波器是FIR滤波器中最重要的一种 应用最广 实际使用时应根据需用选择其合适类型 并在设计时遵循其约束条件 4 2窗口设计法 时域 如果希望得到的滤波器的理想频率响应为 那么FIR滤波器的设计就在于寻找一个系统函数 频率响应去逼近 逼近方法有三种 窗口设计法 时域逼近 频率采样法 频域逼近 最优化设计 等波纹逼近 时间窗口设计法是从单位脉冲响应序列着手 使h n 逼近理想的单位脉冲响应序列hd n 我们知道hd n 可以从理想频响通过付氏反变换获得 但一般来说 理想频响是分段恒定 在边界频率处有突变点 所以 这样得到的理想单位脉冲响应hd n 往往都是无限长序列 而且是非因果的 但FIR的h n 是有限长的 问题是怎样用一个有限长的序列去近似无限长的hd n 最简单的办法是直接截取一段hd n 代替h n 这种截取可以形象地想象为h n 是通过一个 窗口 所看到的一段hd n 因此 h n 也可表达为h n 和一个 窗函数 的乘积 即h n w n hd n 在这里窗口函数就是矩形脉冲函数RN n 当然以后我们还可看到 为了改善设计滤波器的特性 窗函数还可以有其它的形式 相当于在矩形窗内对hd n 作一定的加权处理 设计步骤 1 由定义 3 卷积 插值 一 矩形窗口法 则 以一个截止频率为 c的线性相位理想低通滤波器为例 讨论FIR的设计问题 a 对于给定的理想低通滤波器 计算 低通滤波器的延时 理想特性的hd n 和Hd 这是一个以为中心的偶对称的无限长非因果序列 如果截取一段n 0 N 1的hd n 作为h n 则为保证所得到的是线性相位FIR滤波器 延时应为h n 长度N的一半 即 其中 b 计算 c 计算 设为窗口函数的频谱 用幅度函数和相位函数来表示 则有其线性相位部分则是表示延时一半长度 矩形窗函数及其幅度函数 见P94图4 4 对频响起作用的是它的幅度函数 理想频响也可以写成幅度函数和相位函数的表示形式Hd ej Hd e j 其中幅度函数为两个信号时域的乘积对应于频域卷积 所以有 如果也以幅度函数和相位函数来表示H ej 则实际FIR滤波器的幅度函数H 为正好是理想滤波器幅度函数与窗函数幅度函数的卷积 矩形窗的卷积过程 P95的图4 5来说明 4个特殊频率点看卷积结果 1 0时 H 0 等于 在 c c 内的积分面积 因一般 故H 0 近似为 在 内的积分面积 2 c时 一半重叠 H c 0 5H 0 3 c 2 N时 第一旁瓣 负数 在通带外 出现正肩峰 4 c 2 N时 第一旁瓣 负数 在通带内 出现负肩峰 窗口函数对理想特性的影响 改变了理想频响的边沿特性 形成过渡带 宽为 等于WR 的主瓣宽度 决定于窗长 过渡带两旁产生肩峰和余振 带内 带外起伏 取决于WR 的旁瓣 旁瓣多 余振多 旁瓣相对值大 肩峰强 与N无关 决定于窗口形状 N增加 过渡带宽减小 肩峰值不变 因主瓣附近其中x N 2 所以N的改变不能改变主瓣与旁瓣的比例关系 只能改变WR 的绝对值大小和起伏的密度 当N增加时 幅值变大 频率轴变密 而最大肩峰永远为8 95 这种现象称为吉布斯 Gibbs 效应 改变窗函数的形状 可改善滤波器的特性 窗函数有许多种 但要满足以下两点要求 窗谱主瓣宽度要窄 以获得较陡的过渡带 相对于主瓣幅度 旁瓣要尽可能小 使能量尽量集中在主瓣中 这样就可以减小肩峰和余振 以提高阻带衰减和通带平稳性 但实际上这两点不能兼得 一般总是通过增加主瓣宽度来换取对旁瓣的抑制 肩峰值的大小决定了滤波器通带内的平稳程度和阻带内的衰减 所以对滤波器的性能有很大的影响 几种常用的窗函数 1 矩形窗 上面已讲过 不再细述2 汉宁窗 升余弦窗 利用付氏变换的移位特性 汉宁窗频谱的幅度函数W 可用矩形窗的幅度函数表示为 三部分矩形窗频谱相加 使旁瓣互相抵消 能量集中在主瓣 旁瓣大大减小 主瓣宽度增加1倍 为 0 80 60 44 20 0 矩形窗 Hanning窗 dB 3 汉明窗 改进的升余弦窗 它是对汉宁窗的改进 在主瓣宽度 对应第一零点的宽度 相同的情况下 旁瓣进一步减小 可使99 96 的能量集中在窗谱的主瓣内 4 布莱克曼窗 三阶升余弦窗 增加一个二次谐波余弦分量 可进一步降低旁瓣 但主瓣宽度进一步增加 为 增加N可减少过渡带 频谱的幅度函数为 窗口函数的频谱N 51 A 20lg W W 0 四种窗函数的比较 5 凯塞窗以上四种窗函数 都是以增加主瓣宽度为代价来降低旁瓣 凯塞窗则可自由选择主瓣宽度和旁瓣衰减 I0 x 是零阶修正贝塞尔函数 参数 可自由选择 决定主瓣宽度与旁瓣衰减 越大 w n 窗越窄 其频谱的主瓣变宽 旁瓣变小 一般取4 9 5 44接近汉明 8 5接近布莱克曼 0为矩形 图2凯塞窗函数 图1零阶修正贝塞尔函数 I0 x x 0 1 而当M N时 hM n hd n 零阶贝塞尔函数 窗口设计法的主要工作是计算hd n 和w n 当较为复杂时 hd n 不容易由反傅里叶变换求得 这时一般可用离散傅里叶变换代替连续傅里叶变换 求得近似值 令 过渡带宽At阻带最小衰减 窗函数法的MATLAB实现 FIR滤波器的4个频带分别为0 w pf1pf2 w pf3pf4 w pf5pf6 w pFIR滤波器在4个频带中的幅度值为a1a2a3a4 通带取1 阻带取0 FIR滤波器在4个频带中的波动d1d2d3d4f f1f2f3f4f5f6 a a1a2a3a4 dev d1d2d3d4 M Wc beta ftype kaiserord f a dev h fir1 N Wc ftype kaiser M 1 beta 例 试用Kaiser窗设计满足下列指标的具有2个通带FIR滤波器 ws1 0 1p wp1 0 2p wp2 0 4p ws2 0 5p ws3 0 6p wp3 0 7p wp4 0 8p ws4 0 9p ds 0 01 f 0 10 20 40 50 60 70 80 9 a 0 1 0 1 0 Rs 0 01 dev Rs ones 1 length a N Wc beta ftype kaiserord f a dev h fir1 N Wc ftype kaiser N 1 beta omega linspace 0 pi 512 mag freqz h 1 omega plot omega pi 20 log10 abs mag xlabel Normalizedfrequency ylabel Gain db grid axis 01 805 Ch4课堂练习 已知N 8的FIR滤波器的DTFT为 求 h n 并问该系统是否是线性相位系统 为什么 解 计算h n N 8 并且 根据计算出来的h n 可以发现 h n 不满足线性相位的条件 h n h N 1 n 或h n h N 1 n 系统为非线性相位系统 与之对应的相位曲线如下 与之对应的幅频曲线如下 4 3频率采样法 工程上 常给定频域上的技术指标 所以采用频域设计更直接 一 基本思想使所设计的FIR数字滤波器的频率特性在某些离散频率点上的值准确地等于所需滤波器在这些频率点处的值 在其它频率处的特性则有较好的逼近 内插公式 逼近误差 由得到了H z 或 要讨论与的逼近程度 以及与H k 的关系 由 令 则 单位圆上的频响为 这是一个内插公式 式中为内插函数令则 内插公式表明 在每个采样点上 逼近误差为零 频响严格地与理想频响的采样值H k 相等 在采样点之间 频响由各采样点的内插函数延伸迭加而形成 因而有一定的逼近误差 误差大小与理想频率响应的曲线形状有关 理想特性平滑 则误差小 反之 误差大 在理想频率响应的不连续点附近 会产生肩峰和波纹 N增大 则采样点变密 逼近误差减小 二 设计方法 1 确定2 计算3 计算 并得到 三 约束条件 为了设计线性相位的FIR滤波器 采样值H k 要满足一定的约束条件 前已指出 具有线性相位的FIR滤波器 其单位脉冲响应h n 是实序列 且满足 由此得到的幅频和相频特性 就是对H k 的约束 表4 1 例如 要设计第一类线性相位FIR滤波器 即N为奇数 h n 偶对称 则幅度函数H 应具有偶对称性 令则必须满足偶对称性 而必须取为 同样 若要设计第二种线性相位FIR滤波器 N为偶数 h n 偶对称 由于幅度特性是奇对称的 因此 Hk也必须满足奇对称性 相位关系同上 其它两种线性相位FIR数字滤波器的设计 同样也要满足幅度与相位的约束条件 例 设计一个FIR数字LP滤波器 其理想特性为采样点数N 33 要求线性相位 解 根据表4 1 能设计低通线性相位数字滤波器的只有1 2两种 因N为奇数 所以只能选择第一种 即h n h N 1 n 幅频特性关于 偶对称 也即HK偶对称 利用HK的对称性 求 2 区间的频响采样值 根据指标要求 在0 2 内有33个取样点 所以第k点对应频率为而截止频率0 5 位于之间 所以 k 0 8时 取样值为1 根据对称性 故k 25 32时 取样值也为1 因k 33为下一周期 所以0 区间有9个值为1的采样点 2 区间有8个值为1的采样点 因此 N1 33 k1 0 N1 1 2 Wm1 2 pi k1 N1 Ad1 1 N1 1 2 1 Ad1 10 17 0 Hd1 Ad1 exp j 0 5 N1 1 Wm1 Hd1 Hd1conj fliplr Hd1 2 N1 1 2 h1 real ifft Hd1 w1 linspace 0 pi 0 1 1000 H1 freqz h1 1 w1 plot w1 pi 20 log10 abs H1 grid axis 01 10020 xlabel pi ylabel 幅度 db title 过渡带不设采样点 从图上可以看出 其过渡带宽为一个频率采样间隔2 33 而最小阻带衰减略小于20dB 对大多数应用场合 阻带衰减如此小的滤波器是不能令人满意的 增大阻带衰减三种方法 1 加宽过渡带宽 以牺牲过渡带换取阻带衰减的增加 例如在本例中可在k 9和k 24处各增加一个过渡带采样点H9 H24 0 5 使过渡带宽增加到二个频率采样间隔4 33 重新计算的H ej 见图 其阻带衰减增加到约 40dB N1 33 k1 0 N1 1 2 Wm1 2 pi k1 N1 Ad1 1 N1 1 2 1 Ad1 11 17 0 Ad1 10 0 5 Hd2 Ad1 exp j 0 5 N1 1 Wm1 Hd2 Hd2conj fliplr Hd2 2 N1 1 2 h2 real ifft Hd2 w1 linspace 0 pi 0 1 1000 H2 freqz h2 1 w1 plot w1 pi 20 log10 abs H2 grid axis 01 10020 xlabel pi title 过渡带设一个采样点 2 过渡带的优化设计根据H ej 的表达式 H ej 是Hk的线性函数 因此还可以利用线性最优化的方法确定过渡带采样点的值 得到要求的滤波器的最佳逼近 而不是盲目地设定一个过渡带值 例如 本例中可以用简单的梯度搜索法来选择H9 H24 使通带或阻带内的最大绝对误差最小化 要求使阻带内最大绝对误差达到最小 也即最小衰减达到最大 可计算得H9 0 3904 对应的H ej 的幅频特性 比H9 0 5时的阻带衰减大大改善 衰减接近 50dB 如果还要进一步改善阻带衰减 可以进一步加宽过渡区 添上第二个甚至第三个不等于0的频率取样值 当然也可用线性最优化求取这些取样值 红色为H9 0 3904 3 增大N如果要进一步增加阻带衰减 但又不增加过渡带宽 可增加采样点数N 例如 同样边界频率 c 0 5 以N 65采样 并在k 17和k 48插入由阻带衰减最优化计算得到的采样值H17 H48 0 5886 在k 18 47处插入经阻带衰减最优化计算获得的采样值H17 H48 0 1065 这时得到的H ej 过渡带为6 65 而阻带衰减增加了20多分贝 达 60dB以上 当然 代价是滤波器阶数增加 运算量增加 N 65 k 0 N 1 2 Wm 2 pi k N Ad 1 N 1 2 1 Ad 18 0 5886 Ad 19 0 1065 Ad 20 33 0 Hd Ad exp j 0 5 N 1 Wm Hd Hdconj fliplr Hd 2 N 1 2 h real ifft Hd w linspace 0 pi 0 1 1000 H freqz h 1 w plot w pi 20 log10 abs H grid axis 01 10020 xlabel pi title 过渡带设两个采样点 总采样点数提高一倍 小结 频率采样设计法优点 直接从频域进行设计 物理概念清楚 直观方便 适合于窄带滤波器设计 这时频率响应只有少数几个非零值 典型应用 用一串窄带滤波器组成多卜勒雷达接收机 覆盖不同的频段 多卜勒频偏可反映被测目标的运动速度 缺点 截止频率难以控制 因频率取样点都局限在2 N的整数倍点上 所以在指定通带和阻带截止频率时 这种方法受到限制 比较死板 充分加大N 可以接近任何给定的频率 但计算量和复杂性增加 4 4FIR数字滤波器的最优化设计 前面介绍了FIR数字滤波器的两种逼近设计方法 即窗口法 时域逼近法 和频率采样法 频域逼近法 用这两种方法设计出的滤波器的频率特性都是在不同意义上对给定理想频率特性Hd ej 的逼近 说到逼近 就有一个逼近得好坏的问题 对 好 坏 的恒量标准不同 也会得出不同的结论 我们前面讲过的窗口法和频率采样法都是先给出逼近方法 所需变量 然后再讨论其逼近特性 如果反过来要求在某种准则下设计滤波器各参数 以获取最优的结果 这就引出了最优化设计的概念 最优化设计一般需要大量的计算 所以一般需要依靠计算机进行辅助设计 最优化设计的前提是最优准则的确定 在FIR滤波器最优化设计中 常用的准则有 最小均方误差准则 最大误差最小化准则 1 均方误差最小化准则 若以E ej 表示逼近误差 则那么均方误差为 均方误差最小准则就是选择一组时域采样值 以使均方误差 这一方法注重的是在整个 频率区间内总误差的全局最小 但不能保证局部频率点的性能 有些频率点可能会有较大的误差 对于窗口法FIR滤波器设计 因采用有限项的h n 逼近理想的hd n 所以其逼近误差为 如果采用矩形窗则有 可以证明 这是一个最小均方误差 所以 矩形窗窗口设计法是一个最小均方误差FIR设计 根据前面的讨论 我们知道其优点是过渡带较窄 缺点是局部点误差大 或者说误差分布不均匀 2 最大误差最小化准则 也叫最佳一致逼近准则 表示为其中F是根据要求预先给定的一个频率取值范围 可以是通带 也可以是阻带 最佳一致逼近即选择N个频率采样值 或时域h n 值 在给定频带范围内使频响的最大逼近误差达到最小 也叫等波纹逼近 优点 可保证局部频率点的性能也是最优的 误差分布均匀 相同指标下 可用最少的阶数达到最佳化 例如 我们提到的频率采样最优化设计 它是从已知的采样点数N 预定的一组频率取样和已知的一组可变的频率取样 即过渡带取样 出发 利用迭代法 或解析法 得到具有最小的阻带最大逼近误差 即最大的阻带最小衰减 的FIR滤波器 但它只是通过改变过渡带的一个或几个采样值来调整滤波器特性 如果所有频率采样值 或FIR时域序列h m 都可调整 显然 滤波器的性能可得到进一步提高 低通滤波器的误差分配 切比雪夫最佳一致逼近如图 用等波纹逼近法设计滤波器需要确定五个参数 M c r 1 2按上图所示的误差容限设计低通滤波器 就是说要在通带0 p内以最大误差 1逼近1 在阻带 r 内以最大误差 2逼近零 要同时确定上述五个参数较困难 常用的两种逼近方法 1 给定M 1 2 以 c和 r为变量 缺点 边界频率不能精确确定 2 给定M c和 r 以 1和 2为变量 通过迭代运算 使逼近误差 1和 2最小 并确定h n 切比雪夫最佳一致逼近 特点 能准确地指定通带和阻带边界频率 一 误差函数定义逼近误差函数 为所设计的滤波器与理想滤波器的幅频特性在通带和阻带内的误差值 是已知的权函数 在不同频带可取不同的值 所要设计的滤波器的幅频特性理想滤波器的幅频特性 例如 希望在固定M c r的情况下逼近一个低通滤波器 这时有 对于表4 1中的第一种滤波器 于是 切比雪夫逼近问题变为 寻求一组系数使逼近误差的最大值达到最小 即 给定后等效于求最小 二 交替定理 最佳逼近定理 令F表示闭区间的任意闭子集 为了使在F上唯一最佳地逼近于 其充分必要条件是误差函数在F上至少应有 M 2 次 交替 即其中 且属于F 1 至少有M 2个极值 且极值正负相间 具有等波纹的性质 2 由于是常数 所以的极值也就是的极值 逼近方法 固定k M 和 以作为参变量 按照交替定理 如果F上的M 2个极值点频率已知 则由 1 式可得到M 2个方程 为极值点频率对应的误差函数值 注意 极值点频率必须位于和区间内 由于和固定 因而和必为这些极值频率中的一个 设 则应有求解上述方程组可得到全部系数问题 1 实际情况下 M 2个极值点频率未知 2 直接求解上述非线性方程组比较困难 雷米兹 Remez 算法给出了求解切比雪夫最佳一致逼近问题的方法 三