📘 第 3 章：FIR 数字滤波器

这一篇在干嘛？

最常用的滤波器：从卷积公式讲到转置结构、线性相位、流水线与分布式算法（DA），最后落地到 FIR 设计与 IP 内核。 原书代码为 VHDL，本篇所有代码已改写为 Verilog。

• 3.1 数字滤波器概述
• 3.2 FIR理论
• 例3.1 可编程FIR滤波器
• 3.3 设计 FIR 滤波器
• 3.4 常系数FIR设计
• 例3.2 4抽头直接FIR滤波器
• 例3.3 改进的4抽头直接FIR滤波器
• 直接形式的流水线 FIR 滤波器
• 算法3.4 简化加法器图
• 例3.5 F6半带滤波器的简化加法器图
• 定理3.6 RAG等价变换
• 流水线 RAG FIR 滤波器
• 例 3.7 有符号 DA FIR 滤波器
• 采用逻辑单元的分布式算法
• 例3.8 5输入DA表
• 采用嵌入式阵列模块的 DA
• 加速 DA 滤波器
• 例3.9 DA FIR滤波器的开环形式
• BEGIN
• 例 3.10 F6 半带滤波器的 IP 生成
• 3.5 练习

---

自测一下

FIR 滤波器"线性相位"的条件是什么？为什么重要？

系数对称（h[n] = h[N-1-n]）时相位严格线性。线性相位意味着所有频率分量延迟相同，波形不会发生相位失真——这对音频/图像至关重要。

转置结构比直接结构的 FIR 好在哪？

转置结构天然形成流水线寄存器，不需要额外插入寄存器就能跑高时钟；加法链从串行长链变成并行短链，关键路径短。

分布式算法（DA）为什么不用乘法器？

DA 把”系数乘输入再累加”重组为对输入各 bit 的查表：系数组合预先存进 LUT，按 bit 位移位累加即可。FPGA 的 LUT 资源丰富，查表比乘法器更快更省。

FIR 系数怎么用 MATLAB 得到？

用 fir1/firpm（Parks-McClellan 等纹波）或 fdatool 设计，指定阶数/边界频率/窗或权重，输出量化为 FPGA 字长的定点系数，再验证量化后频响是否仍满足指标。

第3章 FIR 数字滤波器（上）：结构、设计与常系数实现

3.1 数字滤波器概述

什么是数字滤波器

数字滤波器的任务，是按照我们的意愿去修正或改变信号的特性——可以是在时域里，也可以是在频域里。最常见、最重要的一类数字滤波器是线性时不变滤波器（Linear Time-Invariant, LTI）。“线性”指输入叠加时输出也叠加，“时不变”指今天输入一个信号和明天输入同样的信号，系统行为完全一样。

LTI 滤波器与输入信号的相互作用，通过一个叫线性卷积的过程完成，记作 $y=f\ast x$：$f$ 是滤波器的脉冲响应（可以理解为滤波器的”指纹”），$x$ 是输入信号，$y$ 是输出。卷积的正式定义是：

$$y[n]=x[n]\ast f[n]=\sum _{k}x[k]f[n-k]=\sum _{k}f[k]x[n-k]\text{\hspace{1em}(3-1)}$$

直观地说：想要算出第 $n$ 个输出点，就把滤波器系数”反着套”在输入序列最近的一段时间上，逐点相乘再全部加起来。

FIR 与 IIR 的区别

LTI 数字滤波器分为两大类：

• FIR（Finite Impulse Response，有限脉冲响应）：脉冲响应只有有限个采样值，式(3-1)的无限累加退化成每个输出时刻的有限次累加；
• IIR（Infinite Impulse Response，无限脉冲响应）：脉冲响应理论上无限长，需要实现无限累加和（第 4 章讨论）。

“不用会怎样”的反问在这里很自然：正因为 FIR 的累加是有限的，它天然稳定、可以做到严格线性相位，也最容易在 FPGA 上用流水线跑出高速度——这是本章选择它的原因。

为什么用数字滤波器取代模拟滤波器

研究数字滤波器的动机在于：它们正日益成为最主要的 DSP 运算，正在迅速取代用 RLC 元件和运算放大器搭建的传统模拟滤波器。两者思路完全不同：

• 模拟滤波器用拉普拉斯变换的微分方程描述，在时域或 s（拉普拉斯）域分析，只能用于 IIR 设计；
• FIR 通常直接用计算机规范和算法设计，天然适合数字实现。

本章假定 FIR 滤波器已经设计完成，重点放在：先简要回顾 FIR 的设计过程，再讨论如何用 FPGA 高效地实现它。

3.2 FIR 理论

从卷积到传递函数

带常系数的 FIR 滤波器是一种 LTI 系统。对输入序列 $x[n]$，长度为 $L$（即 $N=L-1$ 阶）的 FIR 输出由有限卷积和给出：

$$y[n]=x[n]\ast f[n]=\sum _{k=0}^{L-1}f[k]x[n-k]\text{\hspace{1em}(3-2)}$$

其中 $f[0]\ne 0$ 到 $f[L-1]\ne 0$ 是滤波器的 $L$ 个系数，它们恰好就是 FIR 的脉冲响应。在 z 域中描述更便捷，传递函数 $F(z)$ 定义为：

$$Y(z)=F(z)X(z)\text{\hspace{1em}(3-3)}$$ $$F(z)=\sum _{k=0}^{L-1}f[k]z^{-k}\text{\hspace{1em}(3-4)}$$

抽头延迟线结构

图 3-1 直接形式的 FIR 滤波器

图 3-1 画出了长度为 $L$ 的 FIR：它由”抽头延迟线”、加法器和乘法器混合构成。输入信号沿一串寄存器逐拍向后传递，每个寄存器节点引出一个”抽头”，抽头信号乘以对应的滤波器系数（也叫”抽头权重”），所有乘积累加就是输出。过去人们也称 FIR 为”横向滤波器”，说的就是这种抽头延迟线结构。

看懂这张图，FIR 的硬件成本就一目了然：$L$ 个乘法器、$L-1$ 个加法器、$L-1$ 级延迟寄存器。后面所有优化都在压缩这三样东西。

全零点滤波器

式(3-4)中多项式 $F(z)$ 的根定义了滤波器的零点，根的个数就是阶数 $N=L-1$。FIR 只有零点、没有极点，所以有时被称作”全零点滤波器”。

有一个容易混淆的点：第 5 章会讲一类重要的 CIC 滤波器，它是递归的，但仍然是 FIR。这之所以可能，是因为递归部分产生的极点恰好被非递归部分抵消掉了，等效的零极点图只剩零点。也就是说：非递归滤波器一定是 FIR，但递归滤波器既可能是 FIR 也可能是 IIR。

图 3-2 结构和脉冲长度之间的关系

3.2.1 具有转置结构的 FIR 滤波器

直接 FIR 模型（图 3-1）有一个重要变形：转置 FIR 滤波器。构造方法只有三步：

1. 交换输入和输出；
2. 颠倒信号流的方向；
3. 加法器与分支节点互换角色。

图 3-3 转置结构的 FIR 滤波器

数学上转置结构与原结构完全等价（同样的传递函数），但硬件上转置结构通常是 FIR 的首选实现方式，原因有二：

• 输入 $x[n]$ 广播到所有乘法器，不需要为输入准备额外的移位寄存器链；
• 各乘积沿寄存器逐级累加，寄存器天然分布在加法器之间，无须为高吞吐量再给加法器树添加额外的流水线级。

例3.1 可编程 FIR 滤波器

先回顾用 PDSP 计算积之和（SOP）的位宽要求（参阅 2.8 节）：数据/系数位宽为 $B_x$、滤波长度为 $L$ 时，无符号 SOP 需要额外 ${\log }_2(L)$ 位，有符号运算需要 ${\log }_2(L)-1$ 个保护位。对 9 位有符号数据/系数、$L=4$ 的情况，加法器宽度必须是

$$9+9+{\log }_2(4)-1=19\text{ 位}$$

为什么是 19 位？ 两个 9 位有符号数相乘得 18 位；4 个这样的乘积相加，最坏情况下结果最多增大 $2^{{\log }_{2}4}=4$ 倍，需要 2 个额外位，但其中 1 位被”两个符号位相乘只贡献一个有效符号位”抵消，所以净增 ${\log }_2(4)-1=1$ 位，合计 19 位。少算一位会溢出，多算一位浪费资源。

原书 VHDL 改写为 Verilog，下面是长度为 4 的可编程滤波器（转置结构，系数可在线加载）：

// 原书 VHDL 改写为 Verilog：可编程 FIR 滤波器发生器（转置结构）
module fir_gen #(
    parameter W1 = 9,    // 输入数据/系数位宽
    parameter W2 = 18,   // 乘积位宽 = 2*W1
    parameter W3 = 19,   // 加法器位宽 = W2 + log2(L) - 1
    parameter W4 = 11,   // 输出位宽
    parameter L  = 4     // 滤波器长度
)(
    input  wire               clk,
    input  wire               reset,    // 异步复位
    input  wire               load_x,   // 0:加载系数 1:加载数据
    input  wire signed [W1-1:0] x_in,   // 系统输入
    input  wire signed [W1-1:0] c_in,   // 系数输入
    output wire signed [W4-1:0] y_out   // 系统输出
);
    reg signed [W1-1:0] x;
    reg signed [W1-1:0] c [0:L-1];      // 系数（抽头）阵列
    reg signed [W3-1:0] a [0:L-1];      // 部分和阵列
    wire signed [W2-1:0] p [0:L-1];     // 乘积阵列
    integer i;
 
    // 用 generate 例化 L 个乘法器（便于以后插入流水线级）
    genvar g;
    generate
        for (g = 0; g < L; g = g + 1) begin : mul_gen
            assign p[g] = c[g] * x;
        end
    endgenerate
 
    // 加载过程：load_x=0 时系数沿抽头延迟线逐拍移位
    always @(posedge clk or posedge reset) begin
        if (reset) begin
            x <= 0;
            for (i = 0; i < L; i = i + 1) c[i] <= 0;
        end else if (!load_x) begin
            c[L-1] <= c_in;             // 新系数从最高抽头进入
            for (i = L-2; i >= 0; i = i - 1)
                c[i] <= c[i+1];         // 系数向前移一格
        end else begin
            x <= x_in;                  // 每次装入一个数据样本
        end
    end
 
    // SOP 过程：计算积之和（转置结构加法链）
    always @(posedge clk or posedge reset) begin
        if (reset) begin
            for (i = 0; i < L; i = i + 1) a[i] <= 0;
        end else begin
            for (i = 0; i < L-1; i = i + 1)
                a[i] <= {{(W3-W2){p[i][W2-1]}}, p[i]} + a[i+1]; // 符号扩展后累加
            a[L-1] <= {{(W3-W2){p[L-1][W2-1]}}, p[L-1]};       // 第一个抽头只有寄存器
        end
    end
 
    assign y_out = a[0][W3-1 -: W4];    // 系数为小数（|f[k]|<=1.0），相当于除以 256
endmodule

这段代码的工作方式：当 load_x=0 时，系数被逐个装入抽头延迟线 c；否则数据字装入寄存器 x。第二段 always 块实现积之和：每个乘积 p[i] 先做一位符号扩展到 19 位，再与后级部分和 a[i+1] 相加，形成转置结构的加法链。所有乘法器通过 generate 循环统一例化，这为后续插入流水线级留好了口子。最后把 19 位 SOP 的高 11 位送给输出，等效于除以 256——因为假定的系数都是 $|f[k]|\le 1.0$ 的小数。该设计使用了 93 个 LE 和 4 个嵌入式乘法器，TimeQuest 缓慢 85C 模型下 $F_{max}=157.38\mathrm{MHz}$。

要仿真这个长度为 4 的滤波器，选一组有代表性的系数——Daubechies DB4 小波滤波器系数：

$$G(z)=\left[(1+\sqrt{3})+(3+\sqrt{3})z^{-1}+(3-\sqrt{3})z^{-2}+(1-\sqrt{3})z^{-3}\right]\frac{1}{4\sqrt{2}}$$ $$G(z)=0.48301+0.8365z^{-1}+0.2241z^{-2}-0.1294z^{-3}$$

FPGA 里没有小数，需要量化：把 4 个系数同时乘以 $2^8=256$ 后四舍五入取整，得到 8 位精度（外加 1 位符号）的整数模型：

$$G(z)=(124+214z^{-1}+57z^{-2}-33z^{-3})/256=\frac{124}{256}+\frac{214}{256}{z}^{-1}+\frac{57}{256}{z}^{-2}-\frac{33}{256}{z}^{-3}$$

仿真的时序是：复位释放后的前 4 个时钟节拍，把系数 $\{124,214,57,-33\}$ 依次移入抽头延迟线 c；接着向寄存器 x 装入 100，检查滤波器的脉冲响应。一个时钟周期后嵌入式乘法器的结果 p(i) 就绪，第一个有效输出 y_out 出现在 375 ns 之后。输出各拍按 $124:214:57:-33$ 的比例出现——这正是”脉冲输入的响应就是滤波器系数”这一事实的硬件验证。

图 3-4 加载了 Daubechies 滤波器系数的 4 抽头可编程 FIR 滤波器的仿真

3.2.2 FIR 滤波器的对称性

FIR 脉冲响应的中心是一个重要的对称点。为方便起见，把这一点定义成第 0 次采样时刻，这样得到的滤波器描述就是因果的。对奇数长度的 FIR，因果模型为：

$$F(z)=\sum _{k=-(L-1)/2}^{(L-1)/2}f[k]z^{-k}\text{\hspace{1em}(3-5)}$$

令 $z=e^{j\omega T}$，沿单位圆求传递函数就得到频率响应：

$$F(\omega )=F(e^{j\omega T})=\sum _{k}f[k]e^{-j\omega kT}\text{\hspace{1em}(3-6)}$$

用 $|F(\omega )|$ 表示幅频响应，用 $\varphi (\omega )$ 表示相位响应：

$$\varphi (\omega )=\arctan \left(\frac{\Im (F(\omega ))}{\Re (F(\omega ))}\right)\text{\hspace{1em}(3-7)}$$

实际工作中，人们更多用相位和幅值来描述数字滤波器，而较少直接使用 z 域传递函数或复频率响应。

3.2.3 线性相位 FIR 滤波器

为什么需要线性相位

在通信和图像处理等应用中，信号经过滤波器后波形形状不能被扭曲，这一点至关重要。例如图像的边缘、通信码元的形状，一旦不同频率分量被延迟了不同时间，波形就会失真。衡量相位线性度的标准尺度是组延迟：

$$\tau (\omega )=-\frac{d\varphi (\omega )}{d\omega }\text{\hspace{1em}(3-8)}$$

理想的线性相位滤波器，在关心频段内组延迟是一个常数——所有频率分量”手拉手”一起延迟，波形只平移、不变形。

对称性带来线性相位

从式(3-7)可以看出：只有当频率响应 $F(\omega )$ 是纯实或纯虚函数时，才可能实现固定组延迟。这就要求脉冲响应必须满足偶对称或奇对称：

$$f[n]=f[-n]\text{或}f[n]=-f[-n]\text{\hspace{1em}(3-9)}$$

以长度为奇数的偶对称 FIR 为例，把正负频率项配对求和：

$$F(\omega )=f[0]+\sum _{k>0}f[k]e^{-jk\omega T}+f[-k]e^{jk\omega T}\text{\hspace{1em}(3-10)}$$ $$=f[0]+2\sum _{k>0}f[k]\cos (k\omega T)\text{\hspace{1em}(3-11)}$$

复指数配对后变成纯余弦项，频率响应是频率的纯实函数——线性相位由此得证。对称（$f[n]=f[-n]$）、反对称（$f[n]=-f[-n]$）与奇数长度、偶数长度共有 4 种组合，表 3-1 总结了它们的零点位置。

表 3-1 4 种可能的线性相位 FIR 滤波器 $F(z)={\sum }_{k}f[k]z^{-k}$

对称性 / 长度	$f[n]=f[-n]$ 奇	$f[n]=f[-n]$ 偶	$f[n]=-f[-n]$ 奇	$f[n]=-f[-n]$ 偶
示例	—	—	n
零点位置	$\pm 120^{\circ }$	$\pm 90^{\circ },180^{\circ }$	$0^{\circ },180^{\circ }$	$0^{\circ },2\times 180^{\circ }$

对称性还能省一半乘法器

线性相位 FIR 固有的对称属性还能降低乘法器数量。观察图 3-5 的偶对称滤波器：先对称位置的两个输入相加（如 $x[n]+x[n-(L-1)]$），再乘以公共系数。因为 $f[k]x_1+f[k]x_2=f[k](x_1+x_2)$，乘法在数学上被”提到”了加法之后。

图 3-5 减少乘法器数量的线性相位滤波器

结果：每个滤波周期内乘法器数量正好是直接结构（图 3-1）的一半（从 $L$ 降到 $L/2$），而加法器数量不变，仍是 $L-1$ 个。在 FPGA 上乘法器远比加法器昂贵，这是一笔划算的交易。

3.3 设计 FIR 滤波器

现代数字 FIR 滤波器都采用计算机辅助工程（CAE）工具设计。本章使用的滤波器由 MATLAB 信号处理工具箱完成，其中”交互式低通滤波器设计”演示覆盖了多种典型设计方法：

• 等波纹（Equiripple，也称极小极大）FIR 设计：采用 Parks-McClellan 和 Remez 交换方法，设计线性相位（对称）的等波纹 FIR，还可用于设计微分器或希尔伯特变换器；
• 凯泽窗函数设计：采用凯泽窗加权的逆 DFT 方法；
• 最小二乘 FIR 方法：通带和阻带内同样存在波纹，但使最小均方误差最小化；
• IIR 方法：第 4 章讨论的巴特沃思、切比雪夫 I/II 和椭圆滤波器等 4 种设计。

本节专攻 FIR。大多数应用对滤波器的幅频响应都有明确要求。例如低通滤波器的规范通常包括：通带 $[0\dots {\omega }_p]$、过渡带 $[{\omega }_p\dots {\omega }_s]$ 和阻带 $[{\omega }_s\dots \pi ]$，其中假定采样频率归一化为 $2\pi$。接下来讨论从这些规范计算系数的两类方法。

3.3.1 直接窗函数设计方法

用逆 DFT 直接合成系数

离散傅里叶变换（DFT）在频域和时域之间建立了直接联系。既然滤波器是按频率响应定义的，那就先把理想频响采样成 $L$ 个频点，再做逆 DFT，得到一组逼近目标频响的 FIR 系数：

$$f[n]=\mathrm{IDFT}(F[k])=\sum _{k}F[k]e^{j2\pi kn/L}\text{\hspace{1em}(3-12)}$$

一个前提来自基础信号理论：实际信号的频谱是厄密共轭（Hermitian）的——实部偶对称、虚部奇对称。要想合成出的滤波器只有实系数，目标 DFT 频谱也必须是厄密共轭的，即 $F[k]=F^{\ast }[-k]$（”*” 表示共轭复数）。否则算出来的系数会带虚部，无法实现。

吉布斯振铃：绕不开的副作用

用矩形窗（即直接截取 16 个频点）设计的低通 FIR 见图 3-6(a)，其通带波纹放大后见图 3-6(b)。滤波器对理想低通给出了合理近似，但最大偏差出现在过渡带边缘，并且过渡带附近出现”振铃”。这个振铃就是著名的吉布斯（Gibbs）现象：有限项傅里叶级数在重现陡峭边缘时必然发生的振荡，振铃幅度约达滤波器阶数量程的 $\pm 7\%$。

(a) $L=16$ 的 FIR 低通滤波器的脉冲响应

(b) 传递函数 $L=16$ 的通带

关键的问题是：加长滤波器能不能消除振铃？ 答案是不能。对比图 3-6(c)、(d) 中长度为 128 的设计：滤波器长度从 16 提高到 128，通带波纹（图 3-6(d)）里边缘附近的振铃数量几乎没变——振铃”挤”得更窄了，但峰值高度依旧。

(c) $L=128$ 的 FIR 低通滤波器的脉冲响应

(d) 传递函数 $L=128$ 的通带

用数据窗抑制振铃

抑制吉布斯振铃的唯一办法是给系数乘上一个数据窗——一个在两端平滑地减小到 0 的函数。窗口覆盖 FIR 的脉冲响应，代价是过渡带变宽，换来的是更”光滑”的幅频响应。例如把凯泽窗应用到 FIR 上，吉布斯振铃被明显压低，如图 3-7(上) 所示。

(a) 传递函数

(b) 通频带的组延迟

(c) 零极点

(a) 传递函数

(b) 通频带的组延迟

(c) 零极点

图 3-7 (上) $L=59$ 的凯泽窗函数设计；(下) $L=27$ 的 Parks-McClellan 设计

已发表的经典窗函数很多，它们之间的差别在于振铃抑制与过渡带加宽之间的折中。常用窗函数 $w[n]$ 有：

• 矩形：$w[n]=1$
• 巴特莱特（三角形）：$w[n]=2n/N$
• 汉宁：$w[n]=0.5(1-\cos (2\pi n/L))$
• 汉明：$w[n]=0.54-0.46\cos (2\pi n/L)$
• 布莱克曼：$w[n]=0.42-0.5\cos (2\pi n/L)+0.08\cos (4\pi n/L)$
• 凯泽：$w[n]=I_0\left(\beta \sqrt{1-(n-L/2{)}^2/(L/2{)}^2}\right)$

表 3-2 常用窗函数的参数

名称	3dB 带宽	第一个零点	最大旁瓣	每倍频程旁瓣衰减	等价凯泽 $\beta$
矩形	0.89/T	1/T	-13dB	-6dB	0
巴特莱特	1.28/T	2/T	-27dB	-12dB	1.33
汉宁	1.44/T	2/T	-32dB	-18dB	3.86
汉明	1.33/T	2/T	-42dB	-6dB	4.86
布莱克曼	1.79/T	3/T	-74dB	-6dB	7.04
凯泽	1.44/T	2/T	-38dB	-18dB	3

读这张表的方法：3dB 带宽指传递函数从直流（DC）衰减 3dB（约 $1/\sqrt{2}$）处的带宽，它越窄过渡带越窄；“最大旁瓣”是相对于第 0 次谐波测量的，直接决定阻带抑制的下限；“每倍频程旁瓣衰减”描述窗函数旁瓣随频率远离主瓣的渐进下降速度。可以看到没有免费的午餐：布莱克曼旁瓣最低（-74dB），但 3dB 带宽也最宽（1.79/T），过渡带代价最大。

基于一阶贝塞尔函数 $I_0$ 的凯泽窗有两个特殊优点：其一，在”振铃抑制 vs 过渡带宽”这组指标上近乎最优；其二，它可以通过参数 $\beta$ 连续调谐。凯泽给出了由所需阻带衰减 $A=20\lg {\epsilon }_r$（同时也是通带波纹的 dB 数）估算 $\beta$ 的经验公式：

$$\beta =\{\begin{array}{ll}0.1102(A-8.7)& A>50\\ 0.5842(A-21{)}^{0.4}+0.07886(A-21)& 21\le A\le 50\\ 0& A<21\end{array}\text{\hspace{1em}(3-13)}$$

给定衰减等级后，滤波器长度可按下式估算（通常有 $\pm 2$ 个抽头的误差）：

$$L=\frac{A-8}{2.285({\omega }_s-{\omega }_p)}+1\text{\hspace{1em}(3-14)}$$

算一个例子：要求阻带衰减 $A=60\mathrm{dB}$、过渡带宽度 ${\omega }_s-{\omega }_p=0.2\pi$。查式(3-13)第一分支：$\beta =0.1102\times (60-8.7)\approx 5.65$；查式(3-14)：$L=(60-8)/(2.285\times 0.628)+1\approx 52/1.435+1\approx 37$。也就是说，约 37 个抽头就能达成指标。

3.3.2 等波纹设计方法

规范更完整的设计目标

窗函数法只规定了频带边界，而实际规范通常还应包括与期望传递函数之间的允许偏差（波纹）。等波纹 FIR 滤波器正是一类满足这类规范的极有效的设计：它的设计准则是最小化与理想传递函数之间的最大偏差——把误差”摊平”，让通带、阻带内的波纹等幅起伏，而不是某处特别差。

等波纹算法适用于许多 FIR 设计实例，最流行的包括：

• 低通滤波器设计（MATLAB 中 firpm(N, F, A, W)），公差方案如图 3-8(a)；
• 希尔伯特滤波器：对通带内所有频率产生 $90^{\circ }$ 相移的单位幅值滤波器（firpm(N, F, A, 'Hilbert')）；
• 微分滤波器：幅值随 $\omega$ 线性增加（firpm(N, F, A, 'differentiator')）。

(a) 公差设计方案

(b) 满足公差设计方案的示例函数

图 3-8 滤波器设计的参数

Parks-McClellan 迭代

等波纹（极小极大）算法通常采用 Parks-McClellan 迭代方法实现，其数学根基是”交错定理”：对给定的公差方案，存在唯一切比雪夫多项式（具有最小长度）满足它。图 3-8(a) 是公差方案，图 3-8(b) 是满足方案的多项式。对低通滤波器，所需滤波器长度可由下式估算：

$$L=\frac{-10\lg ({\epsilon }_p{\epsilon }_s)-13}{2.324({\omega }_s-{\omega }_p)}+1\text{\hspace{1em}(3-15)}$$

其中 ${\epsilon }_p$ 是通带波纹、${\epsilon }_s$ 是阻带波纹。

接着上面的例子算：仍取过渡带 $0.2\pi$，设通带波纹 ${\epsilon }_p=0.01$、阻带波纹 ${\epsilon }_s=0.001$，则 $-10\lg (10^{-5})-13=50-13=37$，$L=37/(2.324\times 0.628)+1\approx 26.3$，约 27 个抽头。

算法每次迭代寻找误差曲线的局部最大值位置，并压低最大误差，直到所有偏差相等；Remez 交换则通过挑选两次迭代之间误差峰值的频率集合来更新频率点。这就是 MATLAB 里的等波纹函数过去叫 remez、后来因 Parks-McClellan 方法改名为 firpm 的缘故。

与窗函数法的对比

与直接频率法（无论是否加窗）相比，等波纹法的核心优势是通带偏差和阻带偏差可以分别指定不同的值。这在音频领域特别有用：人耳只能察觉大于 3dB 的差别，通带波纹可以放宽，把”预算”全部花在阻带抑制上。回看图 3-7：同样的公差要求下，凯泽窗设计需要 59 阶，而等波纹设计只要 27 阶——阶数几乎减半。

3.4 常系数 FIR 设计

从”可编程”到”定制”

只有极少数场合（如自适应滤波器）才需要例 3.1 那种通用可编程滤波器。大多数应用中滤波器是 LTI 的，系数不随时间变化。此时硬件工作简化为：为这组固定系数开发最合适的乘法器与加法器（树）结构。难点不在生成系数（现成软件一键完成），而在把 FIR 设计映射到合适的体系结构：直接形式还是转置形式，按最大速度和最低资源占用取舍。网格滤波器主要用于自适应场合且只适合 PDSP，FPGA 上很少采用，本节集中讨论直接形式、转置形式，最后介绍分布式算法这一完全不同的思路。

3.4.1 直接 FIR 设计

直接 FIR（图 3-1）可以用行为级描述（一个时钟进程写完）实现，也可以用加法器/乘法器的元件例化手工搭建。行为级写法给综合器更多自由，例化则给设计者完全控制。下面用一个长度为 4 的 FIR 说明——虽然 4 阶对实际应用太短，但它便于扩展到高阶且编译快。

假定线性相位（对称）FIR 的脉冲响应为：

$$f[k]=\{-1.0,3.75,3.75,-1.0\}\text{\hspace{1em}(3-16)}$$

这些系数可以直接编码成 5 位小数：$3.75_{10}$ 的 5 位二进制是 $011.11_2$，”.” 是二进制小数点位置。注意：通常实现正系数更有效（正系数的非零项更少，尤其 CSD 编码下），负号可以在输出加法器处统一处理。

实际项目中系数来自设计工具，以浮点数交给设计者；定点化之后必须通过仿真或代数分析验证性能仍满足规范。本例中 3.75 和 1.0 都能被定点数精确表示，验证可以跳过——但这只是特例，不要养成跳过验证的习惯。

定点化必须回答的问题：会不会溢出？

L 阶 FIR 的动态范围最坏增长率 $G$ 很容易计算：

$$G\le {\log }_2\left(\sum _{k=0}^{L-1}|f[k]|\right)\text{\hspace{1em}(3-17)}$$

总位宽 = 输入位宽 + $G$ 的整数部分。对本例，$G={\log }_2(1.0+3.75+3.75+1.0)={\log }_2(9.5)<4$，说明内部数据寄存器至少要比输入多 4 个整数位才不会溢出。如果内部运算用 8 位，输入数据就应限制在 $\pm 128/9.5\approx \pm 13$ 以内。跳过这一步的后果是：一个正常信号突然让累加器溢出绕回，输出出现灾难性的大幅跳变。

例3.2 4 抽头直接 FIR 滤波器

原书 VHDL 改写为 Verilog，系数为 $\{-1,3.75,3.75,-1\}$ 的直接型设计如下：

// 原书 VHDL 改写为 Verilog：4 抽头直接 FIR 滤波器
// 系数 [-1, 3.75, 3.75, -1]，除以 2 的幂用算术右移实现
module fir_srg (
    input  wire               clk,
    input  wire               reset,
    input  wire signed [7:0]  x_in,   // 8 位有符号输入
    output reg  signed [7:0]  y_out   // 8 位有符号输出
);
    reg signed [7:0] tap [0:3];       // 抽头延迟线
    integer i;
 
    always @(posedge clk or posedge reset) begin
        if (reset) begin
            for (i = 0; i < 4; i = i + 1) tap[i] <= 0;
            y_out <= 0;
        end else begin
            // 系数加权求和：3.75 = 2 + 1 + 1/2 + 1/4
            y_out <= 2*tap[1] + tap[1] + (tap[1] >>> 1) + (tap[1] >>> 2)
                   + 2*tap[2] + tap[2] + (tap[2] >>> 1) + (tap[2] >>> 2)
                   - tap[3] - tap[0];
            // 抽头延迟线：逐拍移位
            tap[3] <= tap[2];
            tap[2] <= tap[1];
            tap[1] <= tap[0];
            tap[0] <= x_in;
        end
    end
endmodule

这是图 3-1 直接结构的逐行”翻译”：每个抽头输出乘以相应的 2 的幂组合（$3.75=2+1+0.5+0.25$），再全部相加；乘 $2^k$ 与除 $2^k$ 分别用左移和算术右移实现，不消耗任何乘法器。对称与非对称滤波器都适用这种写法。脉冲输入 10 时，输出的仿真波形见图 3-9：依次出现 $-10,37.5,37.5,-10$ 对应的定点值，脉冲响应与系数完全一致。

图 3-9 脉冲输入为 10 时 FIR 滤波器的 VHDL 仿真结果

三个显而易见的改进措施

1. CSD 编码：用优化的正则有符号数编码（参阅第 2 章例 2.1）实现每个系数，非零项更少；
2. 流水线：给乘法器和输出加法器插入流水线级，加法器排列成平衡树。若系数都是 2 的幂，流水线乘法器和加法器树可以合并；由于 LE 中寄存器通常闲置，流水线开销很低；
3. 利用对称性：若系数对称，乘法复杂度可降到图 3-5 的程度（乘法器减半）。

前两条适用于一切 FIR，第三条只适用于线性相位（对称）滤波器。原书用片段展示了前两条改进后的核心代码，改写为 Verilog 如下：

// 原书 VHDL 改写为 Verilog：利用对称性 + CSD 编码
// 3.75 = 4 - 0.25（CSD 形式）
wire signed [7:0] t1 = tap[1] + tap[2];  // 对称抽头预相加
wire signed [7:0] t2 = tap[0] + tap[3];
always @(posedge clk)
    y_out <= 4*t1 - (t1 >>> 2) - t2;     // 应用 CSD 编码后求和

例3.3 改进的 4 抽头直接 FIR 滤波器

例 3.2 的系数 $3.75=2^2-2^{-2}$ 正是 CSD 编码形式。再把对称性和流水线都用上，各种组合的实测性能见表 3-3。

表 3-3 改进的 FIR 滤波器

对称性	否	是	否	否	是	是
CSD	否	否	是	否	是	是
树	否	否	否	是	否	是
速度/MHz	86.31	154.66	106.62	234.85	149.86	243.01
规模/LE	117	98	64	146	57	82

这张表信息量很大，值得逐列读：

• 对称性既提速又省资源：86.31→154.66 MHz 的同时 LE 从 117 降到 98——乘法器减半的直接收益；
• CSD 省资源最狠（117→64 LE），因为它把每个系数的加法器砍掉了；
• 加法器树提速最猛（86.31→234.85 MHz），代价是流水线寄存器带来 LE 增加（117→146）；
• 三者全用得到最快设计：243.01 MHz，只比只用树多一个 LE，却比基线快近 3 倍。

最紧凑的设计是”对称 + CSD、不用加法器树”（57 LE）；最快的 design 则三种措施齐上。全流水线版本的核心代码（原书 VHDL 改写为 Verilog）：

// 原书 VHDL 改写为 Verilog：全流水线版本，所有操作都打一拍
always @(posedge clk) begin
    t1 <= tap[1] + tap[2];   // 利用系数对称性
    t2 <= tap[0] + tap[3];   // 流水线化的预加法器
    t3 <= 4*t1 - (t1 >>> 2); // 流水线化 CSD 乘法器
    t4 <= -t2;               // 构成二叉树，对齐时延
    y  <= t3 + t4;
end

注意每个中间量 t1~t4 都经过寄存器——这正是”流水线”的含义：把长组合路径切成短段，让时钟频率得以提高，代价是输出比输入晚几拍（延迟增加），但吞吐率（每拍一个结果）不变。

直接形式的流水线 FIR 滤波器

有时个别系数的流水线延迟比其他系数多，可以用 $f[n]z^{-d}$ 建模这个额外延迟。若在前面补一个正延迟：

$$f[n]=z^{d}f[n]z^{-d}\text{\hspace{1em}(3-18)}$$

两个延迟互相抵消，滤波器功能不变。落到硬件上意味着：直接形式的 FIR 中，可以改用寄存器链上提前 $d$ 拍的抽头输出去喂那个”慢”的乘法器，从而让所有路径的时延对齐。这一原则见图 3-10(a)；图 3-10(b) 给出了一个带两级内部延迟的重相位流水线乘法器示例。

(b) 重相位乘法器

图 3-10 重相位 FIR 滤波器

3.4.2 转置结构与 RAG 算法

转置结构的两大改进

3.2.1 节讨论过转置滤波器。对常系数滤波器，转置结构相对直接结构还有两个额外武器：

• 重复系数的多重利用：采用简化加法器图（Reduced Adder Graph, RAG）算法；
• 流水线加法器：采用进位保存加法器（原理见第 2 章）。

流水线加法器提高速度，但要花额外的加法器和寄存器；RAG 原理则降低规模（LE 数），有时也提速。本小节重点研究 RAG。

为什么常系数乘法可以”共享”

第 2 章提过：与逐个实现 CSD 编码相比，常系数的实现可以共享子表达式。例如 CSD 实现系数 93 需要 3 个加法器，而 $93=3\times 31$ 只需 2 个。对转置 FIR，多个系数共享公共因子的概率很高：比如系数 9 和 11，用 $8+1$ 构成 9，再用 $9+2$ 构成 11，总工作量就降到了一个加法器。但找到最优的简化加法器图是 NP 难题，只能依靠启发式算法。Dempster 和 Macleod 提出的 RAG 算法是最经典的方案。

算法 3.4 简化加法器图（RAG）

优化部分（确定性）：

1. 把输入集合中所有系数简化为正奇基数（Odd Fundamental, OF），即删去 $2^k$ 因子；
2. 用 MAG 表（第 2 章表 2-3）查每个系数的单系数加法器成本；
3. 删除集合中所有 2 的幂和重复的基数；
4. 建立能用一个加法器构造的所有系数的图集，并从输入集合中删除这些系数；
5. 检查图中是否有一对基数，能用一个加法器生成输入集合中的某个系数；
6. 重复步骤(5)，直到没有新系数进入图集。

启发部分：

7. 为输入集合中成本最小的系数添加两个加法器（必要时），选用图中基数构成的最小非输出基数（Non-Output Fundamental, NOF，即辅助系数）；
8. 步骤(7)新增的两个基数可能反过来帮助构造其他系数，返回步骤(5)；
9. 把加法器成本为 3 或更高的 OF 加入图集，并为其使用最小 NOF 累加和；
10. 返回步骤(5)，直到所有系数合成完毕。

几个乍看不明显的步骤，解释一下：

• 步骤(1)为何删 $2^k$？ 乘以 2 的幂只是移位，不花硬件；把 2 的幂因子剥掉能让尽可能多的系数共享同一个基数。系数的负号则在输出加法器处统一处理，于是 $-7$ 和 $28=4\times 7$ 可以合并成同一个基数 7。
• 步骤(5)为何允许移位？ 两个基数组合时允许做一次除法，即 $g=(2^u{f}_1\pm 2^v{f}_2)/2^w$，因为乘除 2 就是左右移位，零成本。例如系数集 $\{7,105,53\}$ 逐个 MAG 编码需要 $1+2+3=6$ 个加法器，而 RAG 合成只需 3 个：$7=8-1$；$105=7\times 15$；$53=(105+1)/2$（最后一次除法是右移，免费）。
• 步骤(7)为何选最小 NOF？ 新增的小 NOF 比大的产生更多”免费副产品”。例如需要系数 45，候选 NOF 有 3、5、9、15：选 3，则 $3,6,12,24,48,\dots$ 全部顺手可得；选 15，只多出 $15,30,45,\dots$。副产品越多，后续系数越可能免费搭车。

例3.5 F6 半带滤波器的简化加法器图

Goodman 和 Carey 定义的 F6 半带 FIR 滤波器有 4 个非零系数：$f[0]=346$、$f[1]=208$、$f[3]=-44$、$f[5]=9$。先做成本估算——把十进制转二进制，查 MAG 成本：

$$\begin{array}{lllll}f[0]=346_{10}& =2\times 173& =& 101011010_2& \text{成本 }4\\ f[1]=208_{10}& =2^4\times 13& =& 11010000_2& \text{成本 }2\\ f[3]=-44_{10}& =-2^2\times 11& =& -101100_2& \text{成本 }2\\ f[5]=9_{10}& =3^2& =& 1001_2& \text{成本 }1\\ \text{总和}& & & & 9\end{array}$$

直接 CSD 编码实现共需 9 个加法器。RAG 算法的处理过程：

步骤	要实现的	已实现的	操作
(0)	{346,208,-44,9}	{-}	初始化
(1a)	{346,208,44,9}	{-}	取非负系数
(1b)	{173,13,11,9}	{-}	删除 $2^k$ 因子
(2)	{173,13,11,9}	{-}	查系数成本：{3,2,2,1}
(3)	{173,13,11,9}	{-}	从集合中删除成本为 0 的系数
(4)	{173,13,11}	{9}	实现成本为 1 的系数：$9=8+1$
(5)	{173,13,11}	{9,11,13}	构造 $11=9+2$ 和 $13=9+4$

对剩余系数做启发处理，从最低成本、最小值的系数入手：

步骤	实现	已实现的	寻找表达式的措施
(7)	{-}	{9,11,13,173}	添加 NOF 3：$173=11\times 16-3$

图 3-11 应用 RAG 算法的 F6 的实现

结果（图 3-11）：加法器从 9 个降到 5 个，加法器路径延迟也从 4 降到 3——规模和速度同时受益。

RAG-95 到 RAG-05

实现该算法的程序 ragopt.exe 可在本书学习资料的 util 目录找到。相比最初版本（按发表年份称 RAG-95），多年来的改进包括：

• MAG 查找表扩展到 14 位（Gustafsson 等甚至扩到 19 位）；计算最小 NOF 累加和时考虑所有 MAG 加法器成本为 4 的图。14 位内只有两个系数（14709 和 15573）成本为 5，只要不用到它们，结果就是 RAG-95 系列中最优的；
• 步骤(7)考虑所有加法器成本为 2 的图（共三种形态：单基数加成本 2 的因子、两个基数之和、成本 1 的因子或三个基数之和）；
• 最关键的改进：当需要实现多个成本为 2 的系数时，RAG-95 一律选最小 NOF，可能得到次优结果。例如系数集 $\{13,59,479\}$，RAG-95 选 NOF $\{3,5,7\}$（$13=4\times 3+1$；$59=64-5$；$479=59\times 8+7$），需要 6 个加法器；而选 NOF $\{15\}$（$13=15-2$；$59=15\times 4-1$；$479=15\times 32-1$）只需 4 个加法器，节省 30%。所以与其为每个成本 2 的系数单独选最小 NOF，不如为所有成本 2 的系数联合寻找最佳 NOF。

这些修正合称 RAG-05 算法。

定理 3.6 RAG 等价变换与基准

如何评价不同 RAG 实现的优劣？随机生成滤波器做基准既无法第三方验证，也无实际相关性。更可靠的做法基于如下等价变换：

定理 3.6：设 $S_1$ 为 RAG 算法合成的系数集，$F_1$ 为输出基数集，$N_1$ 为非输出基数集。若所用系数集 $S_2=F_1\cup N_1$（包含第一个集合的全部输出基数和非输出基数），那么合成的 RAG 是一致的。

证明：$S_2$ 中所有基数都在算法优化部分合成，各需一个加法器，故最小加法器数为 $C_2=\#F_1+\#N_1$；若 $S_1$ 合成时生成同样的基数集，其 RAG 同样使用最小数量加法器。两者加法器数相同，必然等价。∎

一个推论是图可分为优化图（NOF 不超过一个）和启发式图（多个 NOF）；一个唯一的 OF 图只需一个 NOF 列表即可描述。表 3-4 给出了若干经典低通滤波器（原型来自 Goodman-Carey、Samueli、Lim-Parker）的对比：RAG-05 所需加法器（等于输出基数数量）显著少于 CSD 和公用子表达式（CSE）编码，且滤波器越长优势越大（如 L1：CSD 145 个 vs RAG 52 个）。注意 CSD 数据已利用了系数对称性 $f(k)=f(L-k)$。

表 3-4 低通滤波器的 CSD、CSE 和 RAG 算法所需加法器的数量

滤波器名称	L	B	CSD 加法器	CSE 加法器	of	nof	RAG-05 加法器	NOF 值
F5	11	8	6	—	3	0	3	—
F6	11	9	9	—	4	1	5	3
F7	11	9	7	—	3	1	4	23
F8	15	10	10	—	5	2	7	11、17
F9	19	13	14	—	5	2	7	13、1261
S1	25	9	11	6	6	0	6	—
S2	60	14	57	29	26	0	26	—
L1	121	17	145	57	51	1	52	49
L2	63	13	49	23	22	0	22	—
L3	36	11	16	5	5	0	5	—

Goodman-Carey 的半带滤波器 F5~F9 用 RAG-05 都取得了优于 RAG-95 的结果。Samueli 及 Lim-Parker 的基准数据之所以非常适合 RAG，是因为这些低通滤波器的系数在两端光滑地趋于零，提高了输出成本为 1 的基数出现的可能性。第 6 章将讨论更复杂的 DFT 系数的 RAG 设计。

流水线 RAG FIR 滤波器

RAG 图中加法器逐级串联，即使小规模图也会拉低寄存器性能。改进办法是利用 LE 中闲置的寄存器：在加法器输出端放一个寄存器不花额外逻辑资源。只做这一步（一级流水线）就比无流水线快约 50%。

要做到全流水线，图中每个加法器的输入路径必须有相同时延。以 F6 为例，各系数的构造式与时延为：

$$x9\le 8x+x;\text{时延 }1$$ $$x11\le x9+2x{z}^{-1};\text{时延 }2$$ $$x13\le x9+4x{z}^{-1};\text{时延 }2$$ $$x3\le x{z}^{-1}+2x{z}^{-1};\text{时延 }2$$ $$x173\le 16x11-x3;\text{时延 }3$$

可见输入 x 需要一个额外的流水线寄存器，图的最长时延是 3 级。现在出现新问题：各系数的时延不同（x9 是 1 拍，x173 是 3 拍），直接接到抽头延迟线上会打乱滤波器的时间对齐。两个解决方案：

1. 给所有系数补齐到相同时延（F6 补到 3 拍），输出抽头延迟线结构无需改动；
2. 流水线重定时：根据各乘法器的流水线级数，相应调整它在抽头延迟线中的接入位置——原理与图 3-10 的直接形式重相位相同，如图 3-12 所示（为了只用二输入加法器，x13 系数还需一个额外寄存器做延迟）。

图 3-12 通过流水线重定时的 F6 RAG 滤波器

表 3-5 RAG 算法的 F6 流水线解决方案

流水线级	LE	$F_{max}$ (MHz)	成本 ($LE/F_{max}$)
0	224	152.7	1.47
1	237	199.76	1.19
最大值	254	319.49	0.79
增益 (0/最大值)	-13%	109%	86%

表 3-5 的数据说明：流水线重定时让这个半带滤波器的速度几乎翻倍（152.7 → 319.49 MHz），LE 只多 13%；按 $LE/F_{max}$ 衡量的总成本改善 86%。结论明确：全流水线设计是值得的。

3.4.3 采用分布式算法的 FIR 滤波器

最后介绍一种思路完全不同的 FIR 体系结构，它建立在 2.8.1 节的分布式算法（Distributed Arithmetic, DA）基础上。与传统”乘积-求和”结构逐个样本逐个系数地算不同，DA 每一步同时处理所有系数的同一个二进制位 b：把 $L$ 个输入位的所有 $2^L$ 种组合预先存进一个小查找表，每拍查一次表、移位累加一次，$B$ 位数据 $B$ 拍算完一个输出。整个滤波器只需要一个小表和一个带移位器的累加器；有符号 DA 滤波器则需要带符号的累加器。以第 2 章例 2.25 中系数为 $\{1,-3,7\}$ 的 3 系数 FIR 为例即可看清这一计算过程。

与乘加结构相比，DA 的妙处在于：它把乘法彻底消灭了——乘法被查找表取代，特别适合乘法器资源紧张而 LUT 富余的场合。这为常系数 FIR 的实现提供了与 RAG 完全互补的另一条路。

例 3.7 有符号 DA FIR 滤波器

前面讲过的无符号 DA 滤波器有一个前提：所有输入都被当作普通的无符号二进制数。但实际的 DSP 系统里，信号几乎都是有符号的（二进制补码）。如果直接把补码当成无符号数去查表，最高位（符号位）的权重就会算错，输出全盘皆错。因此有符号 DA 滤波器必须”额外多一个状态”——专门处理符号位的那一步，把符号位的贡献从”加”改成”减”。

回忆一下二进制补码的定义：一个 4 位补码 $x$ 的值是

$$x=-8\cdot x_3+4\cdot x_2+2\cdot x_1+1\cdot x_0$$

也就是说，低 3 位的权重是正的 $2^k$，而符号位 $x_3$ 的权重是 $-2^{B-1}=-8$。把它代入 DA 的累加式中，就得到处理办法：对第 0、1、2 位查表后正常累加，而对符号位（最后一步）查表后要减去 $8\times$表值。这就是代码里变量 count 的作用——它记录当前处理到第几位，当 count = 3（处理符号位）时切换为减法。

原书 VHDL 改写为 Verilog，下面是有符号串行 DA 滤波器的顶层模块（3 抽头，系数 $\{h_0,h_1,h_2\}=\{-2,3,1\}$，4 位输入）：

// 原书 VHDL 改写为 Verilog：3 抽头有符号 DA FIR 滤波器（串行形式）
module da_fir3_serial (
    input  wire        clk,      // 系统时钟
    input  wire        reset,    // 异步复位
    input  wire [3:0]  x0_in,    // 第 1 路输入
    input  wire [3:0]  x1_in,    // 第 2 路输入
    input  wire [3:0]  x2_in,    // 第 3 路输入
    output wire [2:0]  lut,      // DA 查找表输出（测试信号）
    output reg  signed [6:0] y   // 系统输出（-64..63）
);
    localparam INI = 1'b0, RUN = 1'b1;
 
    reg         state;             // 状态机：ini / run
    reg  [3:0]  x0, x1, x2;        // 输入移位寄存器
    reg  [2:0]  count;             // 移位计数器（处理符号位的关键）
    reg  signed [6:0] p;           // 临时累加寄存器
 
    // 取三个输入字的当前最低位拼成 LUT 地址
    wire [2:0] table_in = {x2[0], x1[0], x0[0]};
    wire signed [2:0] table_out;  // 表值范围 -2..4
 
    case3s u_table (.table_in(table_in), .table_out(table_out));
    assign lut = table_out;
 
    always @(posedge clk or posedge reset) begin
        if (reset) begin                    // 异步复位
            state <= INI;
            x0 <= 4'b0; x1 <= 4'b0; x2 <= 4'b0;
            p <= 0; y <= 0; count <= 0;
        end else begin
            case (state)
                INI: begin                  // 初始化：装入新输入
                    state <= RUN;
                    count <= 0;
                    p     <= 0;
                    x0 <= x0_in; x1 <= x1_in; x2 <= x2_in;
                end
                RUN: begin                  // 处理：移位/累加
                    if (count == 3'd4) begin
                        y <= p;             // 内积完成，输出结果
                        state <= INI;       // 开始下一个内积
                    end else begin
                        if (count == 3'd3)
                            p <= p/2 - table_out*8;  // 符号位：减
                        else
                            p <= p/2 + table_out*8;  // 其余位：加
                        // 三个寄存器同时右移 1 位
                        x0 <= {1'b0, x0[3:1]};
                        x1 <= {1'b0, x1[3:1]};
                        x2 <= {1'b0, x2[3:1]};
                        count <= count + 1;
                        state <= RUN;
                    end
                end
            endcase
        end
    end
endmodule

和第 2 章的做法一样，这里用的是”移位/累加器”结构：每一步只把累加器右移 1 位（p/2），而不是每步把表值左移 $k$ 位。这样做的好处是加法器始终只做”表值 × 8”这种固定权重的运算，硬件规整、成本低。系数表由 case3s 模块给出（原书用 dagen.exe 工具自动生成），下面是它的 Verilog 版本：

// 原书 VHDL 改写为 Verilog：系数 {-2, 3, 1} 的 DA 查找表
module case3s (
    input  wire        [2:0] table_in,
    output reg  signed [2:0] table_out   // 范围 -2..4
);
    always @(*) begin
        case (table_in)
            3'b000: table_out =  0;   // 0*(-2)+0*3+0*1
            3'b001: table_out = -2;   // 1*(-2)
            3'b010: table_out =  3;   // 1*3
            3'b011: table_out =  1;   // (-2)+3
            3'b100: table_out =  1;   // 1*1
            3'b101: table_out = -1;   // (-2)+1
            3'b110: table_out =  4;   // 3+1
            3'b111: table_out =  2;   // (-2)+3+1
            default: table_out = 0;
        endcase
    end
endmodule

表里的每一个值都是”地址每一位分别乘以对应系数再求和”的预计算结果，例如地址 011（即 $x_0=1,x_1=1,x_2=0$）对应 $1\times (-2)+1\times 3+0\times 1=1$。

图 3-13 给出了输入序列 $\{1,-3,7\}$ 的仿真结果，仿真中同时显示了 clk、reset、state、count、4 个输入信号，以及为寻址 DA LUT 而从输入字中选出的 3 个位。LUT 依次输出 $\{2,1,4,3\}$，按上面讲的移位累加规则：

$$y=2+1\times 2+4\times 4-3\times 8=-4$$

注意最后一项的符号是负的，这正是符号位那一”减”的体现。这个设计消耗 52 个 LE，没有用嵌入式乘法器，时序性能 $F_{max}=258.4\mathrm{MHz}$（TimeQuest 缓慢 85C 模型）。

图 3-13 3 抽头有符号 FIR 滤波器在输入 $\{1,-3,7\}$ 时的仿真结果

用 CASE 语句定义 DA 表的好处是：综合器会直接用逻辑单元（LE）来实现这个小表，表很小时设计又快又省。但表一旦变大，CASE 方式就不行了，必须另找办法。一种替代是使用 9-kbit 嵌入式存储模块（Embedded Memory Block，M9K，第 1 章介绍过），可以把一块 M9K 配置成 $2^{10}\times 9$、$2^{11}\times 4$、$2^{12}\times 2$ 或 $2^{13}\times 1$ 的表。下面分别讨论这两条路线。

采用逻辑单元的分布式算法

对于低阶滤波器（抽头数 $L\le 4$），LUT 的地址空间很小，用 LE 实现的 DA 非常有吸引力。这里有一个重要观察可以利用：FIR 滤波器是线性滤波器，所以几个低阶滤波器的输出相加，就等价于一个高阶 FIR 滤波器的输出（见图 2-40）。也就是说，一个 8 抽头的 DA 滤波器可以拆成两个 4 抽头 DA 滤波器再做一次加法，而不需要 8 输入的大表。

Cyclone IV 器件的一个 LE 本质上是一个 $2^4\times 1$ 位的表，因此一个 4 系数的 DA 表用一组 LE 就能实现。问题在于：LE 的数量随阶数呈指数增长。整个 FPGA 里 LE 的数量远多于 M9K——例如 EP2C115 有 115K 个 LE，却只有 432 个 M9K，而且 M9K 还要留给 RAM、FIFO 等”更值钱”的用途，所以有时要节约 M9K。反过来，如果硬要用一个巨大的 CASE 语句实现 $2^b\times b$ 的表（例如 $2^9\times 9$ 的流水线表），综合结果会超过 100 个 LE，得不偿失。图 3-14 给出了用 dagen.exe 生成的 CASE 语句表（输入/输出 3~9 位）综合后所需 LE 数量的对比。

图 3-14 不同编码方式下采用带有 b 位输入和输出 CASE 语句的合成结果的规模比较

另一条出路是：仍然只用 4 输入的 CASE 表（贴合 LE 的物理结构），多于 4 个输入的部分用二叉树形的 2→1 多路复用器拼起来。这种结构还有一个附带好处：很容易插入额外的流水线寄存器——为了跑出最高速度，可以在每个 LUT 之后和每个 2→1 多路复用器之后都放一级寄存器。代价是 LE 总数可能比”一个大表最小化”的方案更多。下面的例子演示了一个 5 输入表的结构。

例 3.8 5 输入 DA 表

5 个输入超出了一个 4 输入 LE 表的地址宽度，于是拆成：低 4 位一个 CASE 表、最高位再一个 CASE 表，最后用一个 2→1 总线多路复用器选择。dagen.exe 接受滤波器长度和系数 $\{1,3,5,7,9\}$（表值的最大值为 $1+3+5+7+9=25$），自动生成代码。原书 VHDL 改写为 Verilog 后如下：

// 原书 VHDL 改写为 Verilog：系数 {1,3,5,7,9} 的 5 输入流水线 DA 表
module case5p (
    input  wire       clk,
    input  wire [4:0] table_in,
    output reg  [4:0] table_out        // 范围 0..25
);
    reg [3:0] lsbs;                     // 低 4 位（打一拍）
    reg [1:0] msbs0;                    // 最高位（打两拍，与表输出对齐）
    reg [4:0] table0out00, table0out01; // 两个 4 输入子表
 
    // 输入寄存级
    always @(posedge clk) begin
        lsbs     <= table_in[3:0];
        msbs0[0] <= table_in[4];
        msbs0[1] <= msbs0[0];
    end
 
    // 子表 00（最高位 = 0）：只含系数 1,3,5,7 的部分和
    always @(posedge clk) begin
        case (lsbs)
            4'd0:  table0out00 <=  0;  4'd1:  table0out00 <=  1;
            4'd2:  table0out00 <=  3;  4'd3:  table0out00 <=  4;
            4'd4:  table0out00 <=  5;  4'd5:  table0out00 <=  6;
            4'd6:  table0out00 <=  8;  4'd7:  table0out00 <=  9;
            4'd8:  table0out00 <=  7;  4'd9:  table0out00 <=  8;
            4'd10: table0out00 <= 10;  4'd11: table0out00 <= 11;
            4'd12: table0out00 <= 12;  4'd13: table0out00 <= 13;
            4'd14: table0out00 <= 15;  4'd15: table0out00 <= 16;
            default: table0out00 <= 0;
        endcase
    end
 
    // 子表 01（最高位 = 1）：上表再加 9
    always @(posedge clk) begin
        case (lsbs)
            4'd0:  table0out01 <=  9;  4'd1:  table0out01 <= 10;
            4'd2:  table0out01 <= 12;  4'd3:  table0out01 <= 13;
            4'd4:  table0out01 <= 14;  4'd5:  table0out01 <= 15;
            4'd6:  table0out01 <= 17;  4'd7:  table0out01 <= 18;
            4'd8:  table0out01 <= 16;  4'd9:  table0out01 <= 17;
            4'd10: table0out01 <= 19;  4'd11: table0out01 <= 20;
            4'd12: table0out01 <= 21;  4'd13: table0out01 <= 22;
            4'd14: table0out01 <= 24;  4'd15: table0out01 <= 25;
            default: table0out01 <= 0;
        endcase
    end
 
    // 最后一级 2→1 多路复用器（带流水线寄存器）
    always @(posedge clk)
        table_out <= msbs0[1] ? table0out01 : table0out00;
endmodule

5 个输入最终生成了 2 个 CASE 表和 1 个 2→1 总线多路复用器。这个多路复用器也可以不自己写，而是通过例化 LPM 库中的 busmux 函数来实现。dagen.exe 会输出名为 caseX.vhd 的文件（X 是滤波器长度，也是输入位宽）；caseXp.vhd 与之相同，只是额外带了流水线寄存器。这种元件可以直接嵌入状态机式（串行）滤波器或开环滤波器结构中使用。

参照图 3-14 可以看到，这种结构化的写法确实改善了 LE 用量。图 3-15 则从速度角度比较了不同设计方法：使用 busmux 函数生成的代码，所有流水线设计的速度可以达到 M9K 方案的近两倍；如果没有流水线级，综合工具虽然能把 LE 数量压得更低，但时序性能也会下降。综合工具无法对一个大 CASE 语句做同样程度的优化——即便用 8 级流水线实现 $2^9\times 9$ 的表能获得很高的 $F_{max}$，设计规模也大到不适合很多应用。此时可以考虑图 2-39 的分区技术（把系数分成两组分别查表再相加），或者用下面讨论的 M9K 来实现。

图 3-15 采用 CASE 语句的不同编码类型的速度比较

采用嵌入式阵列模块的 DA

用 M9K 实现短 FIR 滤波器并不划算：一方面 M9K 数量有限，另一方面 M9K 的最高时序速度是 305 MHz，而 LE 表实现往往可以更快。更别扭的是，M9K 只有一个地址译码器——如果想实现一个 $2^3\times 3$ 的小表，也得独占一整块 M9K，剩下的容量完全浪费掉，不能另作他用。

但是对于更长的滤波器，M9K 就很有吸引力了，原因有二：

• M9K 的时序吞吐量恒定为 305 MHz，不随表变大而下降；
• 表放进存储块后，LE 层面的布线工作量显著降低（布线拥塞常常才是长滤波器设计的真正瓶颈）。

加速 DA 滤波器

前面例 3.7 的串行 DA 用一个表”分时复用”，一个内积要花 $B$ 个时钟（$B$ 为输入位宽）。想要更快，可以采用开环（并行）结构：输入按字逐次采样（每次一个完整的字），位并行送入。此时输入的每一位都需要一个单独的表。有趣的是，虽然表的个数变了，每个表的内容却完全不变——因为每个表做的都是同一件事：“把这一位上各抽头的贡献加起来”。若用元件方式定义 LE 表，只需例化同一个 case3s 模块多次即可，代码规模几乎不增加。下面把前面 3 系数、4 位输入的例子改写成开环形式。

例 3.9 DA FIR 滤波器的开环形式

典型 FIR 应用中，输入是按字并行处理的（见图 3-16）：4 位输入字 $x_{in}$ 的每一位 $x_{in}(k)$ 各自进入一个 3 级移位寄存器，保存最近三个采样字的第 $k$ 位；每个移位寄存器的 3 位输出去查一个 $2^3\times 4$ 的表（内容与例 3.7 完全相同），得到该位权重的部分和 $h(k)$；最后按二进制权重把 4 个部分和加权求和，其中最高位是符号位，权重取负：

$$y=h(0)+2\cdot h(1)+4\cdot h(2)-8\cdot h(3)$$

原书 VHDL 改写为 Verilog，开环（位并行）DA 滤波器如下：

// 原书 VHDL 改写为 Verilog：并行（开环）DA FIR 滤波器
module da_fir_parallel (
    input  wire       clk,
    input  wire       reset,
    input  wire [3:0] x_in,           // 系统输入（一个完整字）
    output reg signed [6:0] y         // 系统输出（-46..44）
);
    reg [2:0] x [0:3];                // 4 个 3 位移位寄存器
    wire signed [2:0] h [0:3];        // 4 个表的输出
 
    // 每个输入位一个表，共 4 个，内容同例 3.7
    genvar k;
    generate
        for (k = 0; k < 4; k = k + 1) begin : lc_tables
            case3s u_table (.table_in(x[k]), .table_out(h[k]));
        end
    endgenerate
 
    integer i;
    always @(posedge clk or posedge reset) begin
        if (reset) begin
            for (i = 0; i < 4; i = i + 1) x[i] <= 3'b0;
            y <= 0;
        end else begin
            for (i = 0; i < 4; i = i + 1)
                x[i] <= {x_in[i], x[i][2:1]};  // 装入新字的第 i 位并右移
            // 加法器树（未加流水线寄存器的版本）
            y <= h[0] + 2*h[1] + 4*h[2] - 8*h[3];
        end
    end
endmodule

这一设计用了 4 个 $2^3\times 4$ 的表，表内容与例 3.7 相同。图 3-17 是输入序列 $\{1,-3,7\}$ 的仿真结果：因为输入是按字串行、按位并行的，每个时钟就能推进一个采样，期望结果 $-4_{10}=1111100_{2C}$（补码）在 $400\mathrm{ns}$ 内就计算完成——比串行版快得多。

图 3-16 分布式算法 FIR 滤波器的并行实现

图 3-17 并行分布式算法 FIR 滤波器的仿真结果

上面的设计消耗 39 个 LE，不用嵌入式乘法器和 M9K，运行速度 205.17 MHz。与通用 MAC（乘累加）设计相比，DA 概念的一大优点是流水线几乎零成本：在表输出和加法器树输出上都可以加额外的流水线寄存器，不需要任何额外硬件代价。

第一步，先把加法器树流水线化，在 PROCESS 中引入信号 $s_0,s_1$：

$$\begin{array}{rl}{s}_0& <=h(0)+2\cdot h(1);\\ s_1& <=h(2)-2\cdot h(3);\\ y& <=s_0+4\cdot s_1;\end{array}$$

时序性能提升到 368.60 MHz，LE 数量不变。要做成完全流水线版本，还需把每个 CASE 表的输出也寄存一拍（引入信号 $t_0\sim t_3$），Verilog 的运算部分变成：

$$\begin{array}{rl}{t}_0& <=h(0);t_1<=h(1);t_2<=h(2);t_3<=h(3);\\ s_0& <=t_0+2\cdot t_1;s_1<=t_2-2\cdot t_3;\\ y& <=s_0+4\cdot s_1;\end{array}$$

规模增加到 47 个 LE——因为保存表输出的寄存器不能再复用作 $x$ 输入的移位寄存器了——但时序性能从 205.17 MHz 提升到了 420 MHz。一句话总结：流水线级数越多越快，代价只是少量寄存器和几拍延迟。

IP 内核 FIR 滤波器设计

Altera 和 Xilinx 都随订阅服务提供 FIR 滤波器生成器，因为 FIR 滤波器是最常用的知识产权（Intellectual Property, IP）模块之一。FPGA 供应商通常更青睐基于分布式算法（DA）的 FIR 滤波器生成器，因为这类设计有四个突出优点：

• 完全流水线的体系结构；
• 编译时间短；
• 资源估计准确；
• 与 RAG 算法相比，面积结果与系数值无关。

最后一点尤其省心：基于 DA 的设计不需要做系数优化，也不用计算 RAG 图（系数集很大时 RAG 计算可能非常耗时）。基于 DA 生成的代码包含全部 HDL 代码和测试平台，用供应商的 FIR 编译器几秒钟就能搭好测试平台。

下面把例 3.5 中 Goodman 和 Carey 的 F6 FIR 滤波器用 Altera 的 FIR 编译器重新生成一遍。Altera FIR 编译器 MegaCore 函数生成的 FIR 滤波器针对 Altera 器件优化，支持 Stratix、Arria 和 Cyclone IV 系列，不支持 APEX 或 Flex 等老器件。通过 IP 工具台 MegaWizard 设计环境可以设定各种滤波器结构：固定系数、多周期变量和多速率滤波器（图 3-18(a)）。FIR 编译器内置系数生成器，也允许加载来自文件（如 MATLAB 计算结果）的预定义系数。

例 3.10 F6 半带滤波器的 IP 生成

操作流程如下：在 Tools 菜单下启动 MegaWizard Plug-In Manager，在 DSP|Filters 下找到 FIR 编译器；先为内核指定设计名称，再进入 ToolBench 工具。因为我们打算用 F6 系数，所以选 Edit Coefficient Set 并勾选 Imported Coefficient Set 加载系数文件——它是一个简单的文本文件，每行一个系数，从第一行开始；系数可以给整数也可以给浮点数（浮点数会被工具量化，因为 FIR 编译器只能生成整数系数滤波器）。图 3-18(b) 的脉冲响应窗口会显示这些系数，可按需修改。

图 3-18(a) FIR 滤波器的 IP 设计：IP 工具台

加载系数后，把 Structure 选为 Distributed Arithmetic: Fully Parallel Filter（DA 完全并行），输入系数宽度设为 8 位，用工具的 Actual Coefficients 方法计算输出位宽。由于整数系数不再量化，Coefficient Scaling 选 None，这样整数与浮点形式的传递函数应当逐行吻合（图 3-19）。FIR 编译器估算规模为 498 个 LE。在 Step 2: Set Up Simulation 中勾选全部选项，让工具生成所有可能的仿真文件；第 3 步生成 HDL 代码及全部支持文件（表 3-6 列出了部分关键文件）：不仅有带元件声明的 HDL 文件，还有 ModelSim 仿真（RTL 与门级）脚本、MATLAB（位精度）和 Quartus II（周期精度）测试向量，验证起来非常方便。

图 3-19 根据 F6 的例 3.5 设定 FIR 内核的 IP 参数化

表 3-6 为 FIR 内核生成的部分 IP 文件

文件	说明
f6.vhd	MegaCore 函数变量文件，顶层说明
f6.cmp	MegaCore 函数的元件声明
f6.bsf	Quartus II 模块图编辑器用的符号文件
f6_ast.vhd	Avalon 流接口的封装
f6_mlab.m / f6_model.m	MATLAB 仿真模型 / 位精度模型
f6.vec	Quartus II 仿真测试向量
tb_f6.vhd / f6_msim.tcl	VHDL 测试平台 / 测试平台脚本
f6.vho / f6.qip / f6.html	仿真模型 / 工程信息 / 报告文件

编译滤波器的 HDL 代码并启用时序仿真后，可以得到准确的资源数据。图 3-20 是 F6 滤波器的 ModelSim 仿真结果：先测试脉冲响应，再做随机数据测试。注意波形的数据格式与习惯不同，至少要把数据类型从 Binary（二进制）切换为 Decimal（十进制）。窗口里最关键的两个信号是 ast_sink_data（滤波器输入）和 ast_source_data（滤波器输出），sink/source 是与滤波器内核交互的”Avalon 流接口”的命名。150 ns 处的脉冲测试在大约 300 ns 之后产生了滤波器系数序列——这段初始延迟正是 DA 设计完全流水线化的直接证据。此外还合成了几个并未要求的额外控制信号。

图 3-20 FIR 内核时序仿真结果

例 3.10 的设计实际消耗 490 个 LE（比估算的 498 略少），运行速度 355.37 MHz。以 $\mathrm{LE}/F_{max}$ 衡量的总成本度量为 1.4，优于无流水线的 RAG 设计。由于 DA 是完全流水线的，从脉冲响应的较大初始时延就能看出来。与完全流水线的 RAG 设计相比：基于 DA 的设计成本更高（LE 更多），但时序性能略高。

基于 DA 和基于 RAG 的 FIR 滤波器的比较

前面 F6 半带滤波器的案例研究究竟是个特例，还是在速度/规模/成本上具有普遍规律？为了回答这个问题，研究者用 HDL 实现了一系列较大的完全流水线 RAG 滤波器（用遗传算法优化），并与 Xilinx FIR 内核编译器生成的完全并行 DA 滤波器、以及来自文献的共享子表达式（Common Sub-Expression, CSE）方法做了对比。由于文献中很多早期结果基于 Xilinx ISE 和 Virtex 4 器件，表 3-7 给出了在这套工具链和器件下的合成结果：ISE FIR 内核生成器 5.0 使用并行 DA；CSE 比较采用公开的 Verilog 基准测试代码；器件选 Virtex 4 XC4VSX25-10FF680（资源足够容纳更大的滤波器）。表中是 16 位输入数据、最小加法器深度（完全流水线需三级流水线）下的结果。

表 3-7 Virtex 4 上 CSE、DA 和 RAG 算法的规模、速度和成本比较

滤波器长度	CSE: LE / Fmax / 成本	DA: LE / Fmax / 成本	RAG: LE / Fmax / 成本
6	265 / 328.7 / 0.81	368 / 313.4 / 1.17	169 / 323.7 / 0.52
10	455 / 285.0 / 1.60	506 / 250.4 / 2.02	246 / 301.8 / 0.81
13	385 / 322.6 / 1.19	699 / 305.1 / 2.29	355 / 303.0 / 1.17
20	373 / 386.5 / 0.96	1013 / 260.4 / 3.89	442 / 292.4 / 1.51
28	1245 / 264.9 / 4.70	1531 / 284.5 / 5.38	628 / 256.0 / 2.45
41	1804 / 257.8 / 7.00	2135 / 256.0 / 8.34	893 / 262.6 / 3.40
61	2687 / 136.5 / 11.4	3159 / 271.4 / 11.6	1267 / 242.4 / 5.23
119	5058 / 235.2 / 21.5	5991 / 312.3 / 19.2	2308 / 235.1 / 9.81
151	6468 / 222.8 / 29.0	7591 / 257.1 / 29.5	2931 / 221.0 / 13.3
平均值	2082 / 282.2 / 8.68	2554 / 279.0 / 9.27	1026 / 270.9 / 4.18

从表 3-7 可以得出三条结论：

• 流水线 RAG 滤波器的规模与基于 CSE 的设计相比平均降低 103%，与基于 DA 的设计相比平均降低 149%；
• 基于 CSE 和基于 DA 的 FIR 滤波器的时序性能比流水线 RAG 设计平均分别只高 4% 和 3%；
• 以 $\mathrm{LE}/F_{max}$ 衡量总成本，基于 RAG 的设计平均优于基于 DA 的设计 122%，平均优于 CSE 方法 108%。

也就是说：RAG 用较少的资源达到几乎相同的速度，综合成本大约只有 DA/CSE 方案的一半；而 DA 的优势在于设计自动化程度高、系数任意变化时结果可预期。选择哪种方案，取决于你的设计更在意资源还是更在意开发效率。

配套练习中的两个参考设计

本章练习 3.10 和 3.11 要求分别用转置形式 CSD FIR 滤波器和分布式算法实现 GC4114 通信 IC 的 31 抽头对称 CFIR 补偿滤波器，它们的测试平台仿真波形如图 3-21 和图 3-22 所示，可作为动手实践时的结果对照。

图 3-21 练习 3.10 中 CSD FIR 滤波器的测试平台

图 3-22 练习 3.11 中基于 DA 的 FIR 滤波器的测试平台