4、数字滤波器设计：FIR与IIR滤波器、窗函数法、双线性变换、滤波器系数量化
好，咱们进入第四讲。这一讲可以说是整个降噪算法的“地基工程”。你想想看，不管算法理论多漂亮，最终在DSP芯片上跑起来的，就是一堆滤波器系数。滤波器设计得好不好，直接决定了降噪效果和功耗。
我个人习惯把滤波器设计分成两条路：FIR和IIR。这两兄弟性格完全不同，咱们得摸透它们的脾气。
4.1 FIR滤波器：稳定，但代价不小
FIR滤波器，全称有限长单位冲激响应滤波器。说白了，它的输出只取决于当前和过去的输入信号，没有反馈回路。这个特性带来了一个巨大的好处：绝对稳定。
我在项目中遇到过好几次，IIR滤波器因为系数量化不当导致自激振荡，输出直接“飞”了。但FIR从来不会出这种问题。所以，如果你对稳定性要求极高，比如在医疗级助听器里，FIR是首选。
但FIR有个“硬伤”：要达到同样的滤波效果，它需要的阶数比IIR高得多。阶数高，意味着计算量大，功耗就上去了。在助听器这种电池供电的设备里，功耗是命根子。
核心要点：FIR是线性相位，不会造成相位失真，这对语音清晰度至关重要。但代价是更高的计算量和延迟。
4.2 IIR滤波器：高效，但需要小心伺候
IIR滤波器，无限长单位冲激响应滤波器。它引入了反馈，所以可以用更低的阶数实现陡峭的过渡带。说白了，就是“性价比”高。
我曾经在一个低功耗项目中，用4阶IIR实现了原本需要32阶FIR才能达到的带通效果。计算量直接降了一个数量级，电池续航从8小时延长到了12小时。嗯，这就是IIR的魅力。
但IIR有两个“坑”：
- 稳定性问题：极点必须在单位圆内。一旦系数量化后极点跑出去了，滤波器就振荡了。
- 非线性相位：会造成相位失真，对语音信号来说，可能会影响可懂度。
避坑指南：我曾经在定点DSP上实现IIR滤波器时，因为系数截断误差，导致低频段出现了自激振荡。后来我改用“直接II型”结构，并做了系数缩放，才解决了问题。记住，IIR的系数对量化误差非常敏感。
4.3 窗函数法：给FIR“剪裁”出理想响应
设计FIR滤波器，最直观的方法就是窗函数法。原理很简单：先算出理想滤波器的无限长冲激响应，然后用一个窗函数把它“截断”成有限长。
窗函数的选择很关键。不同的窗函数，对通带波纹、阻带衰减、过渡带宽的影响完全不同。
| 窗函数类型 | 主瓣宽度 | 旁瓣衰减 | 适用场景 |
|---|---|---|---|
| 矩形窗 | 最窄 | -13dB | 对阻带衰减要求不高时 |
| 汉宁窗 | 较宽 | -31dB | 语音信号处理，平衡性好 |
| 海明窗 | 较宽 | -41dB | 需要较高阻带衰减时 |
| 布莱克曼窗 | 最宽 | -57dB | 对阻带衰减要求极高时 |
我个人习惯在助听器降噪中优先使用海明窗。为什么？因为它在阻带衰减和过渡带宽之间取得了很好的平衡。你想想看，助听器里既要滤除噪声，又不能把语音的细节给“切”没了，海明窗刚好合适。
下面是一个用海明窗设计低通FIR滤波器的示例代码：
// 假设采样率Fs=16kHz，截止频率Fc=4kHz，阶数N=64
// 用海明窗设计低通FIR滤波器
float hamming[N+1];
float h[N+1]; // 滤波器系数
// 1. 生成海明窗
for (int n = 0; n <= N; n++) {
hamming[n] = 0.54 - 0.46 * cos(2 * M_PI * n / N);
}
// 2. 生成理想低通冲激响应
for (int n = 0; n <= N; n++) {
if (n == N/2) {
h[n] = 2 * Fc / Fs; // 中心点
} else {
h[n] = sin(2 * M_PI * Fc * (n - N/2) / Fs) / (M_PI * (n - N/2));
}
}
// 3. 加窗
for (int n = 0; n <= N; n++) {
h[n] *= hamming[n];
}
小技巧：实际工程中，我一般不会直接用窗函数法设计高阶FIR。因为计算量太大。我会先用MATLAB的fdatool工具设计好系数，然后直接导入到C代码里。但理解窗函数法的原理，能帮你更好地理解滤波器特性。
4.4 双线性变换：把模拟滤波器“搬”到数字域
IIR滤波器的设计，通常是从模拟滤波器出发，然后通过双线性变换映射到数字域。为什么用双线性变换？因为它能避免频率混叠。
你想想看，模拟滤波器的频率响应是无限宽的，而数字滤波器的频率响应只能到奈奎斯特频率。如果直接映射，高频部分会混叠到低频，造成失真。双线性变换通过一个非线性映射，把整个模拟频率轴压缩到了数字频率轴内，完美解决了这个问题。
双线性变换的公式很简单：
s = (2/T) * (1 - z^(-1)) / (1 + z^(-1))
其中T是采样周期。但要注意，这个映射是非线性的，会导致频率轴发生“弯曲”。所以，在变换之前，需要对模拟滤波器的截止频率进行“预畸变”校正。
我曾经犯过一个错误：直接用双线性变换设计了一个带通滤波器，结果中心频率偏了200Hz。后来才发现是忘了做预畸变。嗯，这个坑大家一定要记住。
预畸变的公式：
Ω = (2/T) * tan(ω * T / 2)
其中Ω是模拟频率，ω是期望的数字频率。
4.5 滤波器系数量化：从浮点到定点
这是最容易被忽视，但也是最容易出问题的一步。DSP芯片通常用定点数运算，而我们在MATLAB或PC上设计出来的系数是浮点数。把浮点数转成定点数，就是量化。
量化的本质就是“舍入”或“截断”。这会导致滤波器实际响应和理论响应有偏差。偏差大了，滤波器可能就不工作了。
我给大家一个经验值：对于16位定点DSP，IIR滤波器的系数至少要用Q1.15格式（1位符号位，15位小数位）。对于FIR滤波器，Q1.15通常也够用，但如果阶数很高，建议用Q2.14格式，留出更多整数位防止溢出。
量化步骤：
- 确定系数动态范围：找到所有系数的最大值和最小值。
- 选择Q格式：确保所有系数都能在所选格式内表示，且精度足够。
- 量化：将浮点系数乘以2^Q，然后四舍五入取整。
- 验证：用量化后的系数重新计算频率响应，检查是否满足设计要求。
下面是一个系数量化的示例：
// 假设浮点系数在[-1, 1]之间，使用Q1.15格式
#define Q 15
int16_t quantize_coefficient(float coeff) {
// 乘以2^Q，然后四舍五入
return (int16_t)(coeff * (1 << Q) + 0.5f);
}
// 使用示例
float b0 = 0.123456f;
int16_t b0_q = quantize_coefficient(b0); // 结果约为4044
警告：量化后的系数一定要做仿真验证。我曾经在一个项目中，量化后的IIR滤波器系数导致极点移到了单位圆外，滤波器直接振荡了。后来我改用“级联二阶节”结构，每个二阶节的系数单独量化，才解决了问题。
好了，这一讲的内容就到这里。滤波器设计是助听器降噪算法的核心，也是最能体现工程师功底的地方。下一讲，我们会深入讨论自适应滤波算法，敬请期待。