FIR 数字滤波器的设计
上一章我们学习了 IIR 数字滤波器的设计方法--借助模拟低通原型(巴特沃斯、切比雪夫、椭圆),通过双线性变换把成熟的模拟设计"翻译"到数字域。IIR 滤波器可以用很低的阶次实现陡峭的过渡带,但它有一个挥之不去的痛点:相位非线性。在语音处理、图像处理、通信系统等对信号波形保真度要求高的场景中,非线性相位会造成波形失真。FIR 滤波器虽然阶次更高、计算量更大,却能实现严格的线性相位,并且天然绝对稳定。本节将系统介绍 FIR 滤波器的三种主流设计方法:窗函数法、频率采样法、等波纹(切比雪夫)最优逼近法。
前置知识回顾
如果下面这些概念不熟,先回看对应章节;本节会默认使用它们。
- DTFT 与 IDFT:序列与频谱的正反变换关系。去哪里补:DTFT 笔记 / 课件第三讲。
- 系统函数与零极点:$H(z)$ 的零极点分布与频率响应形状的对应关系。去哪里补:Z 变换笔记 / 课件第四讲。
- 理想滤波器的冲激响应:理想低通是无限长 sinc 序列。去哪里补:课件第七讲第 36 页。
- IIR 滤波器设计:双线性变换、模拟原型法的基本思想。去哪里补:IIR 滤波器设计笔记 / 课件第六讲。
FIR 滤波器(Finite Impulse Response)
一个 FIR 数字滤波器的单位冲激响应 $h(n)$ 只在有限个时刻非零:
系统函数是 $z^{-1}$ 的有限次多项式(无极点,除 $z=0$):
$N$ 是滤波器长度(抽头数),阶数为 $N-1$。$H(z)$ 在 $z\to\infty$ 处有 $N-1$ 阶零点,在 $z=0$ 处有 $N-1$ 阶极点--这些极点全部在原点,所以 FIR 滤波器天然绝对稳定。
例题：差分方程系数与系统函数的一一对应
设 LTI 系统的差分方程系数为 $b_k = [1, 2, 1]$，求系统函数 $H(z)$ 和频率响应 $H(e^{j\omega})$。
解：差分方程系数即冲激响应：
Z 变换：
令 $z = e^{j\omega}$ 得频率响应：
FIR 和 IIR 的根本区别在于冲激响应是否有限长。这个看似简单的差异,直接决定了两者在结构、稳定性、相位特性和设计方法上的全面分野:
| 维度 | IIR | FIR |
|---|---|---|
| 冲激响应 $h(n)$ | 无限长(递归) | 有限长 $N$ 点(非递归) |
| 系统函数 $H(z)$ | 有理分式,含极点 | 多项式,仅 $z=0$ 处有极点 |
| 差分方程 | $y(n)=\sum b_k x(n-k)+\sum a_k y(n-k)$ | $y(n)=\sum h(k)\,x(n-k)$ |
| 稳定性 | 需极点在单位圆内 | 始终稳定 |
| 相位特性 | 非线性,难以控制 | 可实现严格线性相位 |
| 阶次 | 低(相同指标下) | 高(通常 5~10 倍) |
| 群延迟 | 非常数 | 常数 $(N-1)/2$(线性相位时) |
| 设计方法 | 模拟原型 + 双线性变换 | 窗函数 / 频率采样 / 最优逼近 |
| 存储/计算 | 少(低阶递归) | 多(高阶卷积或 FFT) |
| 典型应用 | 音频均衡、通信信道 | 波形保真、图像处理、数据传输 |
用一个最小例子来直观感受两者的区别。输入脉冲 $x(0)=1$,其余 $x(n)=0$,观察输出如何演化。
FIR 例子:算完就忘
差分方程:
逐步计算:
- $y(0)=1+0+0=1$
- $y(1)=0+0.5\times 1+0=0.5$
- $y(2)=0+0+0.25\times 1=0.25$
- $y(3)=0+0+0=0$
- $y(4)=y(5)=\cdots=0$ ← 永远为零
脉冲进来,经过 3 个抽头,输出归零。有限长,有终点。
IIR 例子:永远记得自己说过什么
差分方程:
逐步计算:
- $y(0)=1+0=1$
- $y(1)=0+0.5\times 1=0.5$
- $y(2)=0+0.5\times 0.5=0.25$
- $y(3)=0+0.5\times 0.25=0.125$
- $y(4)=0.0625$, $y(5)=0.03125$, $\cdots$
脉冲进来,每次输出都拿上一次的结果乘 0.5 再用。无限长,永远在衰减但永远不停。
| $n$ | FIR $y(n)$ | IIR $y(n)$ |
|---|---|---|
| 0 | 1 | 1 |
| 1 | 0.5 | 0.5 |
| 2 | 0.25 | 0.25 |
| 3 | 0 | 0.125 |
| 4 | 0 | 0.0625 |
| 5 | 0 | 0.03125 |
| $\cdots$ | 0 | $(0.5)^n$ |
前 3 步输出完全一样,但第 4 步 FIR 归零了,IIR 还在继续衰减。原因是 IIR 差分方程里有 $y(n-1)$--反馈:输出反过来喂给自己。FIR 只看输入的历史 $x(n-k)$,不看自己的历史输出。
从这两个数值例子出发,可以推导出 FIR 与 IIR 在稳定性、相位和阶次上的全部差异:
从冲激响应看三个核心差异
1 稳定性。FIR 的 $H(z)$ 是有限多项式,没有分母就没有极点(除 $z=0$),不管输入什么,输出最多延续 $N$ 个采样点就归零,天然绝对稳定。IIR 的 $H(z)$ 有分母,极点 $z=0.5$ 在单位圆内所以 $h(n)=0.5^n$ 衰减;但若把系数改成 $1.1$,极点移到 $z=1.1$ 跑出单位圆,$h(n)=1.1^n$ 指数爆炸。IIR 的稳定性完全取决于极点位置,设计时必须确保所有极点在单位圆内。
2 相位。FIR 的 $h(n)$ 只有有限个非零值,可以自由安排排列--让 $h(n)$ 关于中心对称就得到严格线性相位。IIR 的 $h(n)$ 是无限长的,受极点支配呈指数衰减,无法关于某个中心对称,所以IIR 做不到严格线性相位。
3 阶次。FIR 用 3 个系数实现的衰减速率,IIR 用 1 个反馈系数就能做到类似效果,因为每个输出都在循环利用之前所有输出的累积--一个反馈系数等效于无限多个前馈系数。IIR 阶次更低、效率更高,代价是必须管住极点并接受非线性相位。
一句话:反馈是有代价的自由。FIR 没有 feedback,所以安全、线性相位、但阶次高;IIR 有 feedback,所以高效、但必须管住极点、接受非线性相位。
在深入具体设计方法之前,我们必须先回答一个根本问题:FIR 滤波器凭什么能做到线性相位?要做到这一点,它的冲激响应 $h(n)$ 必须满足什么对称条件?这些对称条件又把滤波器分成了哪几类,每类各自能做什么、不能做什么?
IIR 滤波器依靠递归结构(反馈)实现,系统函数 $H(z)$ 的极点可以靠近单位圆,从而用低阶次获得很窄的过渡带。但代价是:极点位置与相位特性深度耦合,无法独立控制。FIR 滤波器没有反馈,$H(z)$ 只有零点没有极点,这给它带来了两个天然优势:
- 绝对稳定:没有极点,不存在稳定性问题。
- 线性相位可实现:只要 $h(n)$ 满足对称条件,群延迟就是常数。
线性相位的物理意义是:信号通过滤波器后,所有频率分量被延迟的时间完全相同,波形只是整体平移,不会发生色散失真。这在图像滤波、数据通信、雷达信号处理中至关重要。
这句话的数学表述是:如果滤波器的相位响应为 $\varphi(\omega)=-\omega\tau$(关于 $\omega$ 的线性函数,斜率为常数 $\tau$),那么频率为 $\omega$ 的分量被延迟 $\tau$,频率为 $2\omega$ 的分量也被延迟 $\tau$--所有频率的延迟完全相同。结果就是波形形状不变,只是整体晚了 $\tau$ 个采样点。反过来,如果 $\varphi(\omega)$ 不是线性的,不同频率延迟不一样,波形散了,这就是相位失真。
图中对比了同一个 5 点滑动平均的两种写法:
- 上图(因果实现):输出写在 $i+4$,系数序列 $[\frac{1}{5},\frac{1}{5},\frac{1}{5},\frac{1}{5},\frac{1}{5}]$ 关于中心不对称(因果约束),相位非线性,滤波结果相对原始信号有明显位移。
- 下图(对称实现):输出写在 $i+2$(窗口中心),同样的系数但输出位置对齐了窗口中心,等效于 $h(n)=h(N{-}1{-}n)$,相位线性,滤波结果与原始信号完全重合,只有平滑效果没有位移。
这正是 $h(n)=h(N{-}1{-}n)$ 的物理后果:对称系数保证所有频率分量被延迟相同的时间 $\tau=\frac{N-1}{2}$,波形保真。
IIR vs FIR 不是简单的"谁更好"
IIR 阶次低、过渡带窄,但相位非线性且可能不稳定;FIR 相位线性且绝对稳定,但阶次高、延迟大、计算量大。工程上要根据场景取舍,而非一刀切。
FIR 滤波器的系统函数为有限项多项式:
其频率响应 $H(e^{j\omega})$ 可以写成幅度与相位的乘积。为了得到线性相位,需要相位 $\varphi(\omega)$ 与 $\omega$ 成严格的线性关系。数学上可以证明,这要求冲激响应 $h(n)$ 关于中心点 $\alpha=\frac{N-1}{2}$ 满足对称性:
线性相位的充要条件
若 FIR 滤波器具有线性相位,则其冲激响应满足
其中 "$+$" 对应偶对称($h(n)=h(N-1-n)$),"$-$" 对应奇对称($h(n)=-h(N-1-n)$)。群延迟恒为 $\tau=\frac{N-1}{2}$ 个采样点。
为什么对称能保证线性相位?直观理解:对称的 $h(n)$ 相当于把信号按中心"镜像折叠"后再加权求和,这种操作对每个频率分量的延迟完全一致,因此相位斜率是常数。
这里有一个常见的误解需要澄清。对称条件 $h(n)=h(N{-}1{-}n)$ 要求滤波器的"中心"在第 $\frac{N-1}{2}$ 个采样点,也就是说要等到输入信号的"未来"样本到达后才能计算当前输出。有人把这叫作"非因果"--但这个词容易误导。更准确的说法是非即时:
- 频域视角:线性相位意味着 $\varphi(\omega)=-\omega\cdot\tau$,即相位随频率线性变化,斜率恒为 $\tau=\frac{N-1}{2}$。这不是"相位出问题了",而是所有频率分量统一偏移了一个固定斜率。
- 时域视角:这个固定斜率 $\tau$ 就是恒定延时--输出波形整体晚 $\frac{N-1}{2}$ 个采样点出现,形状完全不变。
所以"非因果"不是真的违反因果律,而是为了对称而付出的延迟代价。因果系统完全可以实现线性相位,只要接受输出不是"即时"的,而是固定延迟 $\frac{N-1}{2}$ 个采样。在实时处理中,这个延迟是可预测、可补偿的;在离线处理中,直接把输出对齐回来就行。前面实验图里对称滑动平均的橙线和蓝线完全重合,就是因为输出写在窗口中心 $i+2$,自动补偿了这个延迟。
对称条件有"偶/奇"两种,滤波器长度 $N$ 有"奇/偶"两种,组合起来正好四类。每一类的频率响应形式和适用滤波器类型各不相同,这是 FIR 设计中最核心的约束表,做题时必须首先查表确定类型。
四类线性相位 FIR 的频率响应特征
令 $\alpha=\frac{N-1}{2}$,$M=\left\lfloor\frac{N-1}{2}\right\rfloor$,四类滤波器的频率响应可统一写成
其中 $H_g(\omega)$ 为实函数(幅度函数,可正可负),$\beta=0$(偶对称)或 $\beta=1$(奇对称)。
| 类型 | 对称性 | 长度 N | $H_g(\omega)$ 形式 | $\omega=0$ | $\omega=\pi$ | 可设计类型 |
|---|---|---|---|---|---|---|
| 第一类 | 偶对称 | 奇数 | $h(\alpha)+2\sum\limits_{n=0}^{M}h(n)\cos[(\alpha-n)\omega]$ | 任意 | 任意 | 低通、高通、带通、带阻 |
| 第二类 | 偶对称 | 偶数 | $2\sum\limits_{n=0}^{N/2-1}h(n)\cos[(\alpha-n)\omega]$ | 任意 | $H_g(\pi)=0$ | 低通、带通(不能高通/带阻) |
| 第三类 | 奇对称 | 奇数 | $2\sum\limits_{n=0}^{M}h(n)\sin[(\alpha-n)\omega]$ | $H_g(0)=0$ | $H_g(\pi)=0$ | 仅带通(微分器、希尔伯特) |
| 第四类 | 奇对称 | 偶数 | $2\sum\limits_{n=0}^{N/2-1}h(n)\sin[(\alpha-n)\omega]$ | $H_g(0)=0$ | 任意 | 高通、带通(不能低通/带阻) |
最近几张课件图还可以和这里形成一条完整的理解链：
- 理想低通图给出窗函数法的目标：矩形幅频响应 $|H(e^{j\omega})|$ 在时域对应无限长的 sinc 序列 $h_d(n)=\frac{\sin[\omega_c(n-\tau)]}{\pi(n-\tau)}$，这解释了为什么后面必须“加窗”。
- 5 点滑动平均例题说明“偶对称 + 奇数长度”就是第一类线性相位 FIR，也是最简单的线性相位低通。
- $h(n)=\{2,2,3,0,-3,-2,-2\}$ 例题说明“奇对称 + 奇数长度”属于第三类，并且天然有 $H(0)=0$，所以能滤除直流。
- 汉明窗例题说明只要构造出关于中心 $\tau=\frac{N-1}{2}$ 偶对称的窗序列，就能得到线性相位；当 $N=31$ 时，它属于第一类。
例题：判断 FIR 滤波器的线性相位类型
题目：已知系统差分方程为 $y(n)=\frac{1}{5}\sum_{i=0}^{4}x(n-i)$，判断属于哪一类线性相位滤波器。
解：
由差分方程写出冲激响应：
滤波器长度 $N=5$（$N-1=4$），系数为 $[\frac{1}{5},\frac{1}{5},\frac{1}{5},\frac{1}{5},\frac{1}{5}]$。
观察对称性：$h(0)=h(4)=\frac{1}{5}$，$h(1)=h(3)=\frac{1}{5}$，$h(2)=\frac{1}{5}$，即 $h(n)=h(N-1-n)$，偶对称。
长度 $N=5$（奇数），偶对称 → 查表得 第一类线性相位 FIR。
例题 2：判断线性相位类型与直流滤波能力
题目：已知系统单位脉冲响应 $h(n)=\{2,2,3,0,-3,-2,-2\}$，分析该滤波器属于哪类线性相位？长度是多少？能否滤除直流信号？
解：
1) 对称性与类型判断
序列长度 $N=7$，中心 $\alpha=\frac{N-1}{2}=3$。逐项检验：
- $h(0)=2$, $h(6)=-2$ → $h(0)=-h(6)$
- $h(1)=2$, $h(5)=-2$ → $h(1)=-h(5)$
- $h(2)=3$, $h(4)=-3$ → $h(2)=-h(4)$
- $h(3)=0$ → 满足 $h(3)=-h(3)$
因此 $h(n)=-h(N-1-n)$，奇对称。$N=7$（奇数）+ 奇对称 → 查表得 第三类线性相位 FIR。
部分课件将本题误判为"第二类"。正确判断：第二类 = 偶对称 + $N$ 偶数，第三类 = 奇对称 + $N$ 奇数。本题 $N=7$（奇数）且奇对称，应为第三类。
2) 滤波器阶数
$h(n)$ 中 $z$ 的最高负幂次为 $z^{-6}$，所以阶数为 $N-1=6$，即 6 阶系统。
3) 能否滤除直流
第三类的特征：$H_g(0)=0$。验证：令 $\omega=0$，
频率响应在 $\omega=0$ 处为零，可以滤除直流信号。这正是第三类"$H(0)=0$"约束的直接体现——奇对称的系数正负抵消，直流分量被完全抑制。
例题 3：由汉明窗序列判断线性相位性质
题目：已知 FIR 系统的单位脉冲响应
问：该系统属于何种类型的滤波器？是否具有线性相位特性？
解：
要判断是否具有线性相位，只需检查是否满足对称条件。令 $\tau=\frac{N-1}{2}$，考察
若要有 $h(N-1-n)=h(n)$，必须使
利用恒等式
只要取
即可成立，因此
于是 $N=31$（奇数），并且 $h(n)=h(N-1-n)$（偶对称），故该系统属于 第一类线性相位 FIR 滤波器，具有严格线性相位，群延迟为
线性相位 FIR 的零点不是随意分布的,它们具有严格的共轭倒数对称结构。若 $z_i$ 是 $H(z)$ 的零点,则 $z_i^*$、$1/z_i$、$1/z_i^*$ 也必然是零点--四个零点组成一组"四重对称"。
这种结构来源于线性相位条件的代数要求:$H(z)=\pm z^{-(N-1)}H(z^{-1})$。它强制 $H(z)$ 的零点关于单位圆镜像对称,同时关于实轴共轭对称。
- 当零点恰好落在单位圆上($|z_i|=1$)时,$z_i$ 与 $1/z_i^*$ 重合,退化为二重零点。
- 当零点恰好落在实轴上($z_i=\pm 1$)时,四个点全部重合为单点,且 $\omega=0$ 或 $\omega=\pi$ 处必有过零。
零点分布与滤波器类型的关系
第二类和第三类在 $\omega=\pi$ 处有过零,说明 $z=-1$ 处必有零点;第三类在 $\omega=0$ 处也有过零,说明 $z=1$ 处也必有零点。这是它们不能设计某些滤波器类型的根本原因--零点被对称条件"钉死"在了特定位置。
明确了线性相位的分类约束后,核心问题转化为:如何由给定的频率指标 $H_d(e^{j\omega})$,求出满足线性相位条件的有限长 $h(n)$? 三种方法对应三条不同的路径:时域截断(窗函数法)、频域采样(频率采样法)、误差最小化(等波纹最优法)。
最直观的设计思路是:先写出理想滤波器的频率响应 $H_d(e^{j\omega})$,通过 IDTFT 得到无限长的理想冲激响应 $h_d(n)$,然后用一个有限长的"窗" $w(n)$ 把它截断为 $N$ 点,得到可实现的 $h(n)=h_d(n)\cdot w(n)$。
这里有一个根本问题：为什么要加窗？理想低通滤波器的冲激响应是无限长的 sinc 函数 $h_d(n)=\frac{\sin[\omega_c(n-\alpha)]}{\pi(n-\alpha)}$，向两侧无限延伸。一个无限长的 $h(n)$ 在物理上无法实现——你不能存储无限个历史样本，也不能做无限次乘法。FIR 滤波器的定义就是冲激响应有限长：$h(n)$ 只在 $n=0,1,\ldots,N-1$ 处非零。所以必须截断。窗函数就是截断的工具：$w(n)$ 在 $0\le n\le N-1$ 内非零、外面为零，乘以 $h_d(n)$ 后就把无限长序列"剪"成了 $N$ 点。
反过来，正因为 FIR 允许我们自由选择保留哪 $N$ 个点、怎么过渡到零，所以才有了"加窗"这个设计空间。不同窗函数的区别就在于截断的方式不同：矩形窗是硬切，汉宁/汉明窗是软着陆。截断方式不同，频域的"涂抹"效果就不同——这就是窗函数法设计的全部内容。
理想低通的冲激响应
截止频率为 $\omega_c$ 的理想低通滤波器,其 IDTFT 为
这是一个以 $\alpha$ 为中心、向两侧无限延伸的 sinc 函数。$h_d(n)$ 天然关于 $\alpha$ 偶对称,因此只要窗函数 $w(n)$ 也关于 $\alpha$ 偶对称,截断后的 $h(n)$ 就保持偶对称,满足第一类线性相位。
从频域看,时域乘积对应频域卷积:
因此,$H(e^{j\omega})$ 的形状由窗函数频谱 $W(e^{j\omega})$ "涂抹"理想响应得到。窗函数的频谱特性直接决定了设计结果的三大指标:
- 过渡带宽度:由窗频谱主瓣宽度决定,近似为 $\Delta\omega\approx c\cdot\frac{2\pi}{N}$。
- 通带波纹:由窗频谱主瓣与理想响应的卷积造成。
- 阻带衰减:由窗频谱旁瓣的峰值相对幅度决定。
窗的选择是窗函数法的核心。以下五种窗从“最粗糙”到“最平滑”递进，主瓣越窄则过渡带越陡，旁瓣越低则阻带衰减越大——二者不可兼得。
从图中可以看出窗函数设计的核心权衡：矩形窗的“剪刀”最锋利（主瓣最窄，过渡带最陡），但截断太硬导致频谱旁瓣很高；越往右的窗边缘衰减越柔和，旁瓣越低，但主瓣越宽（过渡带越宽）。五种窗的关键参数对比见下表。
| 窗名称 | 时域表达式 $w(n)$ | 主瓣宽度 | 旁瓣峰值 (dB) | 阻带最小衰减 (dB) | 过渡带 $\Delta\omega$ |
|---|---|---|---|---|---|
| 矩形窗 | $1$ | $4\pi/N$ | $-13$ | $-21$ | $1.8\pi/N$ |
| 三角窗 | $1-\frac{2|n-\alpha|}{N-1}$ | $8\pi/N$ | $-25$ | $-25$ | $6.1\pi/N$ |
| 汉宁窗 | $0.5-0.5\cos\frac{2\pi n}{N-1}$ | $8\pi/N$ | $-31$ | $-44$ | $6.2\pi/N$ |
| 汉明窗 | $0.54-0.46\cos\frac{2\pi n}{N-1}$ | $8\pi/N$ | $-41$ | $-53$ | $6.6\pi/N$ |
| 布莱克曼窗 | $0.42-0.5\cos\frac{2\pi n}{N-1}+0.08\cos\frac{4\pi n}{N-1}$ | $12\pi/N$ | $-57$ | $-74$ | $11\pi/N$ |
工程上最常用的是汉明窗:它在主瓣宽度与旁瓣衰减之间取得了最佳平衡,阻带衰减约 53 dB,足以满足大多数应用。
上图展示了 $N=32$ 时矩形窗(红色)与汉明窗(蓝色)的时域波形。汉明窗在两端逐渐衰减到接近零,这种"软着陆"显著抑制了频谱旁瓣,代价是主瓣变宽一倍。
窗函数法的完整设计流程可以归纳为四步,考试时按此顺序书写,不容易遗漏。
窗函数法设计流程
- 确定滤波器类型:根据指标(低通/高通/带通/带阻)和相位要求,查四类表选定对称类型(通常第一类)。
- 由阻带衰减选定窗:根据 $A_s$ 要求选择窗函数。$A_s\leq 21$ dB 用矩形窗;$A_s\leq 53$ dB 用汉明窗;$A_s\leq 74$ dB 用布莱克曼窗。
- 由过渡带估计阶数 N:利用 $\Delta\omega\approx c\cdot\frac{2\pi}{N}$ 反解 $N$,并向上取整到最近的奇数(若用第一类)。
- 计算并加窗:求 $h_d(n)=\frac{\sin[\omega_c(n-\alpha)]}{\pi(n-\alpha)}$,乘以 $w(n)$ 得 $h(n)=h_d(n)\cdot w(n)$。
其中 $\omega_c$ 通常取通带边缘 $\omega_p$ 与阻带边缘 $\omega_{st}$ 的算术中点:$\omega_c=\frac{\omega_p+\omega_{st}}{2}$。这是因为加窗后通带和阻带各向中间"侵蚀"了一段。
例 7-1:用窗函数法设计线性相位 FIR 低通滤波器
题目:采样频率 $F_s=1.5\times 10^4$ Hz,通带截止频率 $f_p=1.5\times 10^3$ Hz,阻带截止频率 $f_{st}=3\times 10^3$ Hz,阻带衰减 $A_s\geq 50$ dB。求 $h(n)$。
目标:用汉明窗设计第一类线性相位 FIR 低通滤波器。
- 选定窗函数:阻带衰减要求 $A_s\geq 50$ dB。查表可知矩形窗($-21$ dB)和三角窗($-25$ dB)均不满足;汉宁窗($-44$ dB)仍不够;汉明窗($-53$ dB)满足要求。选汉明窗。
- 计算数字边缘频率:
$$\omega_p=2\pi\frac{f_p}{F_s}=2\pi\times\frac{1500}{15000}=0.2\pi$$$$\omega_{st}=2\pi\frac{f_{st}}{F_s}=2\pi\times\frac{3000}{15000}=0.4\pi$$
过渡带宽度 $\Delta\omega=\omega_{st}-\omega_p=0.2\pi$。
- 估计阶数 N:汉明窗过渡带近似 $6.6\pi/N$,故
$$N\approx\frac{6.6\pi}{\Delta\omega}=\frac{6.6\pi}{0.2\pi}=33$$
取 $N=33$(奇数,满足第一类)。群延迟 $\alpha=\frac{N-1}{2}=16$。
- 确定理想截止频率:取通带与阻带的几何中点或算术中点
$$\omega_c=\frac{\omega_p+\omega_{st}}{2}=\frac{0.2\pi+0.4\pi}{2}=0.3\pi$$
- 计算理想冲激响应:
$$h_d(n)=\frac{\sin[0.3\pi(n-16)]}{\pi(n-16)},\quad n=0,1,\ldots,32$$
当 $n=16$ 时,用极限 $h_d(16)=\frac{0.3\pi}{\pi}=0.3$。
- 加汉明窗:
$$w(n)=0.54-0.46\cos\frac{2\pi n}{32},\quad n=0,1,\ldots,32$$$$h(n)=h_d(n)\cdot w(n)$$
答案:$h(n)$ 为 33 点偶对称序列,中心点 $h(16)\approx 0.3\times 1=0.3$,两侧按 sinc 包络乘以汉明窗曲线递减。实际频率响应需用 MATLAB / Python 的 freqz 验证通带波纹和阻带衰减是否满足 50 dB 要求。
上图展示了 $N=33$、$\omega_c=0.3\pi$ 时理想 sinc 序列(灰色)与加汉明窗后的 $h(n)$(蓝色)。窗函数把两侧"削平",同时保持中心区域的对称性。
窗函数法走的是"时域截断"路线。另一条思路是:既然 $h(n)$ 是 $N$ 点有限长序列,它的 DTFT 在单位圆上采样 $N$ 点正好对应 DFT。那么,何不直接在频域"画"出想要的 $H_d(k)$,再反变换回时域?这就是频率采样法。
频率采样法的基本思想
对理想频率响应 $H_d(e^{j\omega})$ 在 $\omega_k=\frac{2\pi k}{N}$($k=0,1,\ldots,N-1$)处采样 $N$ 个点,得到 $H(k)=H_d(e^{j\omega_k})$,再通过 $N$ 点 IDFT 得到
实际频率响应则由内插公式重建:
其中 $\Phi(\omega)=\frac{\sin(N\omega/2)}{N\sin(\omega/2)}e^{-j\omega(N-1)/2}$ 为内插函数。
频率采样法的优点是直接:你把想要的频响"画"在 $N$ 个采样点上,IDFT 后就自动得到 $h(n)$。缺点也很明显:采样点之间的频率响应由内插函数"连接",如果采样点处跳变太剧烈(比如从通带 1 直接跳到阻带 0),内插会在中间产生较大的起伏,导致阻带衰减不够。
为保证线性相位,$H(k)$ 的相位不能随意取,必须满足与四类对称对应的约束。以最常见的第一类(偶对称、$N$ 奇数)为例,幅度采样 $H_k$ 必须为偶对称 $H_k=H_{N-k}$,相位采样则为
第二类(偶对称、$N$ 偶数)还需要额外满足 $H_{N/2}=0$,因为 $H(\pi)=0$。
频率采样法通常分为 I 型和 II 型,区别仅在于第一个采样点的位置:
- I 型:第一个采样点在 $\omega=0$ 处($k=0$),即标准的 $N$ 点 DFT 采样。
- II 型:第一个采样点在 $\omega=\frac{\pi}{N}$ 处,采样点整体偏移半个间隔。
两种类型在通带平坦度和设计灵活性上略有差异,但核心思想完全相同。
频率采样法最经典的问题出现在低通设计:通带最后一个采样点为 1,阻带第一个采样点为 0,中间没有任何过渡。内插函数在跳变边缘会产生类似吉布斯现象的振荡,阻带衰减通常只有 20 dB 左右,远不能满足工程需求。
解决方法是:在通带与阻带之间故意插入一个或几个过渡带采样点,取介于 0 和 1 之间的值。这相当于在频域"修坡",让跳变变缓。代价是过渡带变宽,但阻带衰减可以大幅提高。
过渡带采样实例
题目:用频率采样法设计低通滤波器,$N=33$,不插入过渡点时阻带衰减仅约 20 dB。若在 $k=9$ 处(即 $\omega=\frac{18\pi}{33}\approx 0.545\pi$)插入一个过渡样点 $H(9)=0.5$,其余阻带样点保持为 0,通带样点为 1。
结果:阻带衰减从约 20 dB 提升到约 40 dB。
如果进一步插入两个过渡样点($H(9)=0.5886$,$H(10)=0.1065$),阻带衰减可达 60 dB 以上。过渡样点的最优值可以通过优化算法求解。
上图示意了三种情况:理想矩形频响(蓝色实线)、无过渡带采样时的频响起伏(红色虚线,吉布斯振荡明显)、插入过渡带采样后的平滑结果(绿色)。过渡带采样"牺牲"了过渡带的陡峭度,换取了阻带内更干净的高衰减。
窗函数法和频率采样法都有一个共同缺陷:它们无法精确控制通带波纹 $\delta_1$ 和阻带衰减 $\delta_2$,设计结果是"碰运气"的--阶数 $N$ 选大了浪费计算,选小了指标不达标。有没有一种方法,能在给定阶数 $N$ 的约束下,让最大逼近误差最小化?
答案是等波纹最优逼近法,也称切比雪夫逼近法或 Parks-McClellan 算法。它的数学目标非常清晰:
切比雪夫最优逼近问题
寻找 $h(n)$,使得加权误差的最大绝对值最小化:
其中 $A$ 为通带与阻带的并集(过渡带不参与优化),$W(\omega)$ 为加权函数:
通带误差被放大了 $\delta_2/\delta_1$ 倍,使得优化结果在通带和阻带中的波纹幅度之比恰好为 $\delta_1:\delta_2$。
为什么叫"等波纹"?因为最优解有一个惊人的性质:误差函数 $E(\omega)=W(\omega)[H_d(\omega)-P(\omega)]$ 在逼近区间 $A$ 上的极值点数目达到最大,且所有极值点的绝对值相等,形成均匀分布的波纹。直观上,这就像把一块弹性膜绷紧在理想响应上,膜上的波峰波谷高度完全一致--没有任何一处"特别突出"。
等波纹最优解的存在性和唯一性由交错定理保证。
交错定理(Alternation Theorem)
设 $P(\omega)$ 为 $r$ 阶三角多项式(对应 $N$ 阶 FIR 滤波器的幅度函数),$A$ 为 $[0,\pi]$ 上若干个不相交闭区间的并集。$P(\omega)$ 是 $H_d(\omega)$ 在 $A$ 上的唯一最优切比雪夫逼近,当且仅当误差函数 $E(\omega)=W(\omega)[H_d(\omega)-P(\omega)]$ 在 $A$ 上至少有 $r+2$ 个交错点,即存在 $\omega_0<\omega_1<\cdots<\omega_{r+1}$ 使得
也就是误差在最值之间正负交替,且绝对值相等。
基于交错定理,Parks 和 McClellan 提出了高效的迭代算法--Remez 交换算法。它的核心思想是:
- 在逼近区间 $A$ 内猜测一组 $r+2$ 个极值点频率 $\{\omega_i\}$。
- 解线性方程组,求出这组极值点上的等波纹误差值 $\delta$ 和多项式系数。
- 计算整个区间上的误差函数 $E(\omega)$,找出真正的极值点。
- 用新的极值点替换旧的猜测,重复迭代直到收敛。
Remez 算法收敛非常快,通常只需 5-10 次迭代即可达到机器精度。MATLAB 的 firpm(原 remez)和 SciPy 的 scipy.signal.remez 都是这一算法的实现。
现在我们已经掌握了 FIR 滤波器的三种设计方法。窗函数法概念最直观、手算最方便;频率采样法适合需要精确控制频域某些点响应的场景;等波纹法在给定阶数下性能最优。下面从多个维度进行系统对比。
| 对比维度 | 窗函数法 | 频率采样法 | 等波纹最优法 |
|---|---|---|---|
| 设计路径 | 时域:理想响应 $h_d(n)$ + 窗截断 | 频域:$H_d(k)$ 采样 → IDFT | 优化:最小化最大加权误差 |
| 核心思想 | 用有限长窗"软化"无限长 sinc | 频域打桩 + 内插重建 | 交错定理 + Remez 迭代 |
| 能否精确控制 $\delta_1, \delta_2$ | 不能,只能靠试 $N$ | 不能,过渡带采样需试凑 | 能,直接作为优化约束 |
| 过渡带宽度 | 由窗主瓣决定,较宽 | 由采样密度决定,较宽 | 最窄(同阶下最优) |
| 阻带衰减 | 由窗旁瓣决定,固定值 | 过渡带采样可提升,但有限 | 最大(同阶下最优) |
| 通带/阻带波纹分布 | 不均匀(两端大、中间小) | 不均匀 | 等波纹(均匀分布) |
| 计算复杂度 | 低,解析公式直接求 | 低,一次 IDFT | 中,需迭代求解 |
| 手算可行性 | 高(考试常考) | 中 | 低(依赖计算机迭代) |
| 适用场景 | 快速设计、原型验证、考试 | 需要频域任意点精确控制的场景 | 阶数受限、指标苛刻的工程实现 |
工程实践的选择建议
快速原型或课堂作业 → 窗函数法(汉明窗最常用)。
需要频域陷波或特定频点精确响应 → 频率采样法(配合过渡带优化)。
产品级实现、DSP 芯片资源受限 → 等波纹法(用 MATLAB firpm 设计后量化系数量化)。
无论用哪种方法,FIR 滤波器设计都可以归纳为以下决策流程:
flowchart TD
A[给定滤波器指标] --> B{需要线性相位?}
B -->|是| C[FIR]
B -->|否| D[IIR]
C --> E{阶次受限且指标苛刻?}
E -->|是| F[等波纹最优法
Parks-McClellan]
E -->|否| G{需要手算/快速原型?}
G -->|是| H[窗函数法
汉明/汉宁/布莱克曼]
G -->|否| I[频率采样法
配合过渡带优化]
F --> J[量化系数
实现为 DSP/FPGA]
H --> J
I --> J
复习速查
- 线性相位条件:$h(n)=\pm h(N-1-n)$,群延迟 $\tau=(N-1)/2$。
- 四类 FIR:奇偶对称 × 奇偶长度。第一类万能;第二类不能高通/带阻;第三类只能带通;第四类不能低通/带阻。
- 窗函数法公式:$h_d(n)=\frac{\sin[\omega_c(n-\alpha)]}{\pi(n-\alpha)}$,$h(n)=h_d(n)\cdot w(n)$,$\omega_c=\frac{\omega_p+\omega_{st}}{2}$。
- 常用窗衰减:矩形 $-21$ dB → 汉宁 $-44$ dB → 汉明 $-53$ dB → 布莱克曼 $-74$ dB。
- 频率采样法:$H(k)$ 采样 → IDFT 得 $h(n)$。插入过渡带采样可大幅提升阻带衰减。
- 等波纹最优:最小化 $\max|W(\omega)[H_d-P]|$。交错定理保证 $r+2$ 个等幅交错极值点。Remez 算法迭代求解。
- 方法选择:考试/手算用窗函数法;频域点控用频率采样;阶次受限用等波纹最优。
参考来源
-
西安交通大学数字信号处理课程. 第七讲:FIR 数字滤波器的设计方法.
课件 PDF(114 页) -
Oppenheim, A. V. & Schafer, R. W. (2009). Discrete-Time Signal Processing (3rd ed.). Pearson.
出版社页面 -
Parks, T. W. & McClellan, J. H. (1972). Chebyshev Approximation for Nonrecursive Digital Filters with Linear Phase. IEEE Trans. on Circuit Theory, 19(2), 189-194.
IEEE Xplore