IIR滤波器这名字初听起来像教科书里的概念但我敢说几乎所有接触过数字信号处理的人最终都得跟它打交道。你手上那块STM32读进来的电压、麦克风收到的人声、传感器采集的振动波形想要去掉噪声或者提取特征频率简单粗暴的均值滤波不够用、FFT又不实时的时候IIR滤波器就是最直接的那把刀。这篇内容我从工程落地的角度把IIR滤波器从设计原理、系数提取、代码实现到SOS矩阵整理成一条完整链路重点覆盖STM32这类MCU上用“直接I型”和“SOS矩阵”实现滤波器时会遇到的坑和心得适合刚接触数字滤波的嵌入式工程师、音频开发者以及做信号采集相关项目但还没系统梳理过滤波器设计的朋友。1. 为什么选择IIR与FIR的本质差异1.1 从传递函数角度看两者的不同IIR滤波器的全称是无限脉冲响应滤波器它的核心特征可以用z域传递函数来表示。一个N阶IIR系统的传递函数长这样H(z) (b0 b1·z⁻¹ b2·z⁻² ... bM·z⁻ᴹ) / (1 a1·z⁻¹ a2·z⁻² ... aN·z⁻ᴺ)分母里那串系数a就构成了反馈回路这是IIR和FIR最根本的区别。FIR滤波器的传递函数只有分子没有分母也就没有反馈它的脉冲响应在有限个采样点后就归零了所以叫“有限脉冲响应”。而IIR因为有分母、有反馈理论上一个脉冲输入会在输出端留下一串无限延长的响应尾巴这才是“无限脉冲响应”这个名字的由来。别看这只是数学形式上的差别它带来的工程影响非常大。分母多了一项意味着IIR可以用低得多的阶数实现同样陡峭的过渡带。举个我实际对比过的例子一个采样率48kHz、截止频率10kHz的低通滤波器如果用FIR要达到60dB阻带衰减和较窄过渡带至少需要一百多阶每输入一个采样点就要做一百多次乘加运算而同样规格的IIR我用一个4阶的椭圆滤波器就搞定了单拍计算量只有FIR的零头。这个差距在音频、实时控制这类对延时和算力敏感的场景里基本就是能不能跑得动的区别。1.2 IIR真实性能与代价当然天下没有免费的午餐。IIR用低阶换来的代价是“相位非线性”和“稳定性风险”。相位非线性这一点很多人一开始不在意但真正用起来就会发现问题。你在示波器上看波形可能觉得没差别可如果处理的是心电信号、振动信号或者任何需要保持波形原始形态的场景IIR会把不同频率成分的延时变得不一样导致输出波形“变形”——这不是幅度上的变形而是时间轴上的错位。FIR因为系数对称设计天然能做到线性相位所有频率延时完全一致波形不会走样。所以选型的时候我一般先问自己一个问题这个应用对相位敏感吗如果只是去除电源纹波、平滑传感器噪声、做音频音调调整IIR完全合适如果要做带通滤波后的时域特征分析、波达时间估计那我建议要么用FIR要么接受IIR非线性相位带来的误差要么用零相位双向滤波离线场景弥补。稳定性风险则是另一个大坑。既然IIR有反馈分母那就天然存在“会不会发散”的问题。从z域看系统稳定的条件是所有极点都落在z平面单位圆内。问题在于当我们把设计好的滤波器搬进16位或32位定点芯片、用有限精度的浮点运算实现时系数的量化误差可能导致原本在单位圆内、离边界很远的极点漂移到圆外系统就变成了一个“振荡器”。这也是后面要重点聊SOS矩阵的直接原因——把高阶系统拆成多个二阶子系统级联能显著降低系数量化误差带来的极点漂移风险让滤波器在真实硬件上更稳定。1.3 选型判断什么时候该用IIR说了那么多我可以把这些年做工程选型的经验整理成一张表帮助快速决策对比项IIRFIR相同指标所需阶数低通常4~10阶高动辄几十上百阶单点计算量小大相位特征非线性可做到线性相位稳定性有极点需关注全零点无条件稳定适合场景实时控制、音频均衡、噪声抑制相位敏感测量、多速率处理、心理声学典型工具巴特沃斯、切比雪夫、椭圆窗函数法、等效波法这张表不是说要大家死记硬背我更想强调的是“实时性”这个维度。在STM32这类MCU上一个中断里可能同时要做ADC采集、滤波、PID控制、通信处理留给滤波器的预算往往只有几个微秒到几十个微秒。FIR要是上千阶光这部分的乘加运算就能把CPU吃满。IIR用个小几十次的乘加就完成了同样指标这种差距是实际测出来的不是纸面数据能体现的。结论很简单高性能计算平台、相位敏感场景选FIR嵌入式实时场景、噪声抑制与平滑场景IIR是性价比极高的选择。而且IIR配合SOS结构稳定性和可维护性都能得到保障这也是它至今在工业控制、音频、传感器处理领域占据重要位置的原因。2. 设计IIR滤波器的完整流程与工具实操2.1 设计规格的确定很多人拿到一个需求就急着打开工具箱敲代码我建议先花十分钟把设计规格写清楚这一步省了后面返工十倍的功夫。所谓设计规格核心就是五个参数采样率fs、通带边缘频率、阻带边缘频率、通带最大纹波、阻带最小衰减。采样率是这一切的基础奈奎斯特频率fs/2决定了你能处理的最高频率分量。在数字滤波器设计里所有频率参数最终都要归一化到这个奈奎斯特频率上。比如采样率1000Hz、想保留100Hz以内的信号那归一化截止频率就是100 / (1000/2) 0.2。这个归一化过程用工具的时候也会遇到只是大多数时候工具帮你做了而已。通带纹波和阻带衰减这两个参数决定了滤波器的“质量”。通带纹波1dB意味着信号通过后在幅度上会有大约±5.6%的波动这在音频应用里通常感知不到但在测量仪器上就可能影响精度。阻带衰减40dB意味着阻带信号能削弱到原来的1%这对大多数工业应用已经够了想做到60dB、80dB也不是不行只是滤波器阶数或者设计算法的复杂度会上升。设计规格不是越严越好跟机械加工的公差一样越严成本越高这里的“成本”就是阶数和运算量。2.2 用Python工具箱完成系数设计我自己最常用的设计工具是Python的SciPy库因为它既能在PC上快速验证又能把设计出的系数直接移植到C代码里。下面是一个实际可运行的低通滤波器设计示例import numpy as np from scipy.signal import ellip, butter, cheby1, sosfreqz fs 1000.0 # 采样率 1kHz fc 100.0 # 通带截止频率 100Hz order 4 # 滤波器阶数 rp 1.0 # 通带纹波 1dB rs 40.0 # 阻带衰减 40dB # 用椭圆滤波器设计输出SOS矩阵格式 sos ellip(order, rp, rs, fc/(fs/2), outputsos) print(sos) # 打印频率响应做验证 w, h sosfreqz(sos, worN2048, fsfs)这里用椭圆滤波器是因为它在同样的阶数下过渡带最窄适合在低阶条件下实现尽可能陡的衰减。如果你对通带内的纹波非常敏感、不想要椭圆滤波器的等纹波纹波特征可以把ellip换成butter巴特沃斯它的通带响应是最平坦的但过渡带会宽一些。切比雪夫I型在通带内有等纹波、阻带单调切比雪夫II型则相反。选型逻辑很简单追求最陡过渡带选椭圆追求通带平坦选巴特沃斯两种极端之间选切比雪夫。2.3 设计结果校验拿到系数之后别急着抄进代码。我习惯先打印频率响应曲线做人眼验证——确认通带内增益接近0dB、截止频率位置正确、阻带衰减达到设计目标。有一回我设计一个带通滤波器系数算出来一切正常结果画幅频图发现中心频率偏了5Hz仔细排查发现是设计时归一化频率的参考搞错了。这个错误要是直接烧进板子可能要调大半天才能发现。SciPy的sosfreqz函数可以直接对SOS格式的系数计算频率响应也可以先用tf2sos把普通的分子分母系数转成SOS。如果做的是离线处理还可以直接用sosfilt函数一次性跑完整段数据来验证滤波效果把输入信号和输出信号放在一起对比看看时域波形是不是符合预期。这个“PC端设计→PC端验证→移植到MCU”的流程我到现在还会用因为MCU上调试滤波的代价比PC上高太多。3. 直接I型滤波器在MCU/STM32上的代码实现3.1 直接I型结构的数学推导拿到系数之后怎么把它变成能在MCU上跑的代码呢有一个概念要先讲清楚差分方程。IIR滤波器的传递函数可以等价转换成时域里的差分方程代码就是逐样本地执行这个差分方程。假设我们有这样一个2阶传递函数H(z) (b0 b1·z⁻¹ b2·z⁻²) / (1 a1·z⁻¹ a2·z⁻²)对应的差分方程为y[n] b0·x[n] b1·x[n-1] b2·x[n-2] - a1·y[n-1] - a2·y[n-2]注意分母里的a1、a2在方程里变成了减法而且之前提到的归一化要求a0 1如果设计工具给出的系数里a0不是1需要先把所有系数除以a0。x[n-1]、x[n-2]是前两个输入样本y[n-1]、y[n-2]是前两个输出样本。实现这个方程的时候只需要维护四个历史变量就够了。这种结构之所以叫“直接I型”是因为它直接照着差分方程逐项实现输入历史先走分子输出历史再走分母。另外还有“直接II型”和“转置直接II型”它们的数学等价但数值特性和内存占用略有差别。在MCU上我通常用直接I型或转置直接II型因为它们的状态变量天然是信号延迟线上的数值便于初始化和调试。3.2 第一个能用单二阶节的C语言实现在MCU上用C写IIR滤波器我建议先做成一个简单的结构体加处理函数不要一开始就写很复杂的级联框架。下面这个例子就是最典型的单二阶节实现也就是俗称的biquad双二阶节typedef struct { float b0, b1, b2; // 分子系数 float a1, a2; // 分母系数 float x1, x2; // 输入历史 float y1, y2; // 输出历史 } Biquad; float biquad_process(Biquad *f, float in) { float out f-b0 * in f-b1 * f-x1 f-b2 * f-x2 - f-a1 * f-y1 - f-a2 * f-y2; // 状态更新 f-x2 f-x1; f-x1 in; f-y2 f-y1; f-y1 out; return out; }这段代码的核心就是差分方程的逐样本执行。我特意把系数a0省略了因为设计工具输出的a0通常已经是1如果遇到a0不等于1的情况一定要先做归一化再填入这个结构体。我在给板子移植代码的时候踩过一次直接从MATLAB导出的系数里a0 1.05没归一化就填进去了结果滤波器增益整体偏移响了好久的困惑才排查到问题。同时要注意结构体里的状态值x1/x2/y1/y2在滤波器运行前应该清零。这个“清状态”的步骤有时候比算法本身还关键——如果不清零上电后可能会有一段未知的瞬态输出尤其在控制回路里这种瞬态可能导致执行机构乱动一下风险不小。3.3 与STM32硬件集成的注意事项在STM32上跑这个过滤器不同的人有不同的集成路径我自己的习惯是放在ADC的中断回调或者DMA传输完成回调里每采集到一个样本就调用一次biquad_process。这样滤波是逐样本实时完成的输出可以直接给到后续的PID、FFT或者显示刷新。实时性和代码效率是两个需要同时考虑的问题。在Cortex-M4及以上内核中CMSIS-DSP库提供了arm_biquad_cascade_df1_f32函数它对多阶IIR滤波做了优化利用SIMD指令可以让多个滤波器并行执行。如果项目里只是简单的一两个二阶节手写的biquad代码就够了没必要引入整个DSP库但如果要处理多通道音频或者高阶数滤波强烈建议换到CMSIS-DSP的级联接口省下的不只是代码量更是实打实的CPU占用率。还有一点关系到ADC数据处理的细节读取ADC值之后通常需要先减去直流偏置再做滤波。要是对原始ADC码值直接滤波那个恒定的直流分量会在滤波结果里保留导致输出一直偏在一个非零基线上。对于需要判断阈值或者计算有效值的应用这个直流偏置会让所有后续判断都出错。我一般在滤波前做一次简单的直流扣除或者输入本来就是交流耦合信号就不用担心。3.4 定点化MCU上无法回避的精度问题很多STM32型号没有FPU或者FPU频率低、不适合大量浮点运算。这时候就要考虑把浮点滤波器转成定点实现最常用的格式是Q15和Q31对应16位和32位定点。定点化的核心思想是把系数和信号都乘以一个缩放因子放大成整数来做乘法累加最后再缩回去。比如Q15格式就是把浮点数乘以32768再取整。但这里有两个必须注意的点系数放大后会引入量化误差。原本0.9999的系数可能量化成0.9997这在反馈回路里可能让极点位置发生微小偏移。阶数高的时候这种微小偏移会累积甚至导致不稳定。这也是我为什么强烈推荐SOS结构的原因——每一级只有2阶极点离单位圆的敏感度低得多。中间乘加的溢出问题。biquad内部有5次乘法和4次加法如果用Q15做乘法结果需要32位来保存中间值。如果你每个数都是Q15那乘出来的积是Q30再加上另一个Q15就需要注意数据宽度。很多人在定点化的时候栽在这里要么丢精度要么溢出。如果MCU没有FPU但又不想写定点我有一条折中建议用float32试试。Cortex-M4以上的内核虽然有FPU但即使没有硬件浮点用软件浮点实现2阶IIR的耗时通常在几十微秒级别对采样率不高比如1kHz的应用完全够用。真正必须定点化的场景是那种采样率几十kHz、又要在中断里干很多活的场合。这时候再耐心做定点不要一上来就把自己绕进位宽地狱。4. SOS矩阵让高阶滤波器稳定落地的关键4.1 高阶直接型的数值稳定性隐患先做个思想实验。想象一个10阶IIR滤波器直接用传递函数分子分母那一堆系数去实现相当于一个10阶的反馈系统。任何微小的系数量化误差都像在一根长竹竿的顶端加重量——竹竿越高顶端轻轻一晃底部就产生巨大的偏差。极点分布对系数误差的敏感度跟阶数成正比阶数越高单位圆附近的极点越容易被推出圆外。这就是为什么高阶IIR滤波器直接用“直接型”结构几乎必出问题的根本原因。你可以在MATLAB里看理论频率响应画得完美无缺但把同样的系数写进32位定点MCU跑出来的可能就是自激振荡的噪声。我用MATLAB做过一个实验一个10阶巴特沃斯低通系数保留6位小数后直接实现极点位置跟原始设计差得不算大但一对共轭极点已经落到了单位圆外系统响应变成增长振荡。这个实验特别直观地说明不是设计的问题是实现结构的问题。4.2 什么是SOS矩阵二阶节的级联SOS是Second-Order Sections的缩写中文通常叫二阶节级联。它的核心思想是把一个N阶的传递函数分解成N/2个2阶系统的级联每个2阶系统用biquad结构实现每个biquad的系数单独归一化。SOS矩阵的每一行长这样[ b0, b1, b2, a0, a1, a2 ]其中a0通常等于1工具会在输出时自动归一化。比如一个4阶滤波器输出的SOS矩阵是2行6列就代表两个级联的biquad。第一个biquad的输出喂给第二个biquad的输入串联完成整个滤波。从数值稳定性角度看这个分解的意义在于每一级只承受2阶的极点敏感度即使系数有量化误差影响也限制在本级不会跨级级联放大。极点从“一根高竹竿的一端”变成了“几截短竹竿的连接”稳定裕度和抗量化能力都大幅提升。此外SOS结构还有一个工程上的好处单位“一节”就是最基本的biquad代码结构统一、调试方便每一个节的输入输出都可以单独断点查看过滤器到底哪一级出了问题一目了然。4.3 从普通系数转SOS以及实际代码结构设计工具通常都支持直接输出SOS。用SciPy的时候在ellip、butter这些函数里指定outputsos就可以了如果手里只有普通的分子分母系数b和a数组也可以用scipy.signal.tf2sos做转换from scipy.signal import butter, tf2sos, sosfilt b, a butter(4, 0.2, outputba) sos tf2sos(b, a) print(sos)在C代码里SOS级联实现其实非常直观——就是把上一节的输出作为下一节的输入逐节调用同一个biquad_process函数。我一般用一个数组保存所有SOS行用循环依次处理typedef struct { uint8_t num_sections; float coeffs[MAX_SECTIONS][6]; float state[MAX_SECTIONS][4]; // x1, x2, y1, y2 } SosFilter; float sos_process(SosFilter *f, float in) { float out in; for (int i 0; i f-num_sections; i) { out biquad_process((Biquad){...}, out); // 使用第i组系数 } return out; }需要注意的一个细节是SOS的级联顺序。级联顺序不是随便排的不同排序会影响数值精度。通常原则是把Q因子较高带宽窄、谐振峰尖锐的节放在前面把增益较大的节放在前面也有助于信噪比。SciPy的sosfilt内部会自动排序但如果你手动处理系数建议参考这个原则。另外SOS各节之间可以穿插分配滤波器总增益避免某一节固定增益特别大导致中间信号饱和——这个问题在定点实现里格外突出。5. 容易踩的坑与排查技巧实录5.1 现象输出爆炸先查极点而不是查算法我第一次把IIR滤波器烧进板子的时候内心是有点慌的——上电后输出直接变成接近满幅度的振荡用示波器一看整个波形都在疯狂抖动。当时第一反应是代码写错了来回检查了好几遍biquad函数的乘加顺序结果都没问题。后来仔细一查是设计的时候用了10阶直接型系数经过去浮点转换之后精密值流失极点跑出了单位圆。这种发散问题的排查套路现在已经固定了先用PC端MATLAB或者Python把固定系数的极点画出来看是不是都在单位圆内系数量化之后再看一次极点位置。如果极点没问题再查代码里是否忘记先归一化a0、是否初始化了状态变量、是否加了错误的负号。这个顺序基本能定位90%以上的发散问题。总结起来就是IIR滤波器只要输出不对劲永远先怀疑数值和极点不要怀疑硬件和主循环这是这个领域最值钱的一条经验。5.2 现象截止频率偏移记得双线性变换的预畸变另一个高频问题是按设计算出来的截止频率是100Hz实测出来变成95Hz或者105Hz而且偏差可能随截止频率升高而变大。在我遇到过的大多数情况下原因出在模拟原型变换到数字域时没有做“频率预畸变”。IIR滤波器通常是从模拟滤波器原型巴特沃斯、切比雪夫、椭圆出发通过双线性变换映射到数字域的。双线性变换会把模拟频率轴的无限范围压缩到数字频率轴上的0到π之间这种频率压缩是非线性的如果不做补偿最终数字滤波器的截止频率就会偏移。很多时候设计工具已经自动做了预畸变但如果你手写转换过程或者从模拟设计表直接查系数就容易踩雷。工具选择错了也会导致看起来“偏移”。比如设计采样率是48kHz结果你按44.1kHz的归一化频率设计了那实际效果当然对不上。如果排除代码问题后频率还是有偏差先检查你的设计工具用的采样率、归一化方式再检查预畸变处理大概率能解决。5.3 现象起始瞬态振荡还有一个看起来吓人但其实很常见的情况滤波器一跑起来输出的前几十个样本会出现很大的过冲感觉像“爆炸”但过了这段之后就恢复正常了。这不是滤波器不稳定而是起始瞬态响应——滤波器状态变量从零开始相当于给滤波器输入了一个阶跃信号系统自然会产生瞬态响应。在实时控制应用里这种瞬态可能让执行机构产生一次意外的运动必须处理。我的做法是在系统启动阶段先运行一段时间的滤波但丢弃这段时间的输出或者把状态初始化到输入信号的平均值附近减少阶跃感。在离线处理场景可以用scipy.signal.filtfilt做零相位滤波但它内部也会处理边缘效应原理是反向再滤波一次相位为零但计算量翻倍。5.4 现象信号经过IIR后幅值不对有时候滤波器跑得“很稳”但仔细对比输入输出某个频段的信号不是大了就是小了。这个现象往往是设计规格和实际需求不匹配造成的而不是代码问题。比如你做的是Butterworth 1dB带宽的低通但在通带边缘1kHz-3dB你却拿设计规格里那条“-1dB频点”跟“-3dB频点”混淆了自然会觉得不对。要规避这类问题最好是拿到滤波器后先把频率响应打印一下在0Hz到奈奎斯特频率范围内画一条完整的响应曲线确认不同频点的增益值符合预期。高通的低频、带通的两侧滚降尤其容易让人产生误解因为不同滤波器定义“通带频率”的方式不同有的按-3dB有的按阻带边缘不做校验就是给自己挖坑。我把这几年在IIR滤波器上遇到过的典型问题整理成一张速查表方便直接对照现象可能原因快速处置输出发散/自激振荡极点跑出单位圆、系数未归一化、定点溢出画极点图检查极点位置先改用双精度浮点验证截止频率偏移未做频率预畸变、设计工具采样率设置错误检查设计参数改用工具自动预畸变起始幅度过冲状态变量初始化为零产生的瞬态启动时丢弃前N个样本或初始化状态到稳态值输出幅度整体偏低/偏高增益分配错误、a0未归一化、通带定义混淆打印频率响应曲线验证0Hz通带增益是否为1定点实现噪声大系数量化误差大、中间累加溢出改用更高的定点位宽使用SOS级联结构5.5 排查技巧从“先画图”到“分段定点”关于IIR的调试我有一个坚持了很多年的习惯任何滤波器改动第一步一定是先画频率响应图第二步是给一个已知信号正弦波、方波、阶跃看时域输出第三步才连接到真实数据源上。三步走下来80%的问题能在接入真实系统前暴露剩下的20%再配合逐级断点、打印状态变量等手段也能快速定位。分段定点这个技巧也分享给大家。定点化之前先用浮点在PC端完整跑一遍系统把每一级SOS节点的输出最大值记录下来。这些值就是定点化的缩放依据——该级要用多大Q值不会溢出、不会浪费精度都有依据了。分段定点比整体定点效果好得多因为不同节之间信号动态范围差异可能很大用统一的Q值意味着低幅度节浪费精度、高幅度节面临溢出风险。写在最后回头看IIR滤波器这条路最让我感慨的一点是很多项目失败不是设计思路错了而是被“差一点点”的系数精度、“差一点点”的实现结构、或者“差一点点”的调试方法拖垮的。IIR本身是个非常成熟的技术从巴特沃斯到SOS每一步都有明确的数学基础和实践路径但真正把这条路走通、走稳还是需要一点一点踩坑换来的手感。我个人现在做滤波器项目的习惯是需求拿过来先写规格然后PC端设计校验确认无误再落到MCU代码全部用SOS结构浮点优先、定点按需最后一定要做一次长时间运行测试。这个流程看起来繁琐但它帮我挡掉了太多半夜调板子的尴尬。希望这篇从设计原理到工程落地的经验梳理能让你在IIR滤波器上少走几条我当年走过的弯路。如果你在STM32或者音频项目里遇到过跟IIR相关的怪问题欢迎沿着文中的思路排查一遍——很多时候答案就藏在那张速查表里。
IIR滤波器这名字初听起来像教科书里的概念但我敢说几乎所有接触过数字信号处理的人最终都得跟它打交道。你手上那块STM32读进来的电压、麦克风收到的人声、传感器采集的振动波形想要去掉噪声或者提取特征频率简单粗暴的均值滤波不够用、FFT又不实时的时候IIR滤波器就是最直接的那把刀。这篇内容我从工程落地的角度把IIR滤波器从设计原理、系数提取、代码实现到SOS矩阵整理成一条完整链路重点覆盖STM32这类MCU上用“直接I型”和“SOS矩阵”实现滤波器时会遇到的坑和心得适合刚接触数字滤波的嵌入式工程师、音频开发者以及做信号采集相关项目但还没系统梳理过滤波器设计的朋友。1. 为什么选择IIR与FIR的本质差异1.1 从传递函数角度看两者的不同IIR滤波器的全称是无限脉冲响应滤波器它的核心特征可以用z域传递函数来表示。一个N阶IIR系统的传递函数长这样H(z) (b0 b1·z⁻¹ b2·z⁻² ... bM·z⁻ᴹ) / (1 a1·z⁻¹ a2·z⁻² ... aN·z⁻ᴺ)分母里那串系数a就构成了反馈回路这是IIR和FIR最根本的区别。FIR滤波器的传递函数只有分子没有分母也就没有反馈它的脉冲响应在有限个采样点后就归零了所以叫“有限脉冲响应”。而IIR因为有分母、有反馈理论上一个脉冲输入会在输出端留下一串无限延长的响应尾巴这才是“无限脉冲响应”这个名字的由来。别看这只是数学形式上的差别它带来的工程影响非常大。分母多了一项意味着IIR可以用低得多的阶数实现同样陡峭的过渡带。举个我实际对比过的例子一个采样率48kHz、截止频率10kHz的低通滤波器如果用FIR要达到60dB阻带衰减和较窄过渡带至少需要一百多阶每输入一个采样点就要做一百多次乘加运算而同样规格的IIR我用一个4阶的椭圆滤波器就搞定了单拍计算量只有FIR的零头。这个差距在音频、实时控制这类对延时和算力敏感的场景里基本就是能不能跑得动的区别。1.2 IIR真实性能与代价当然天下没有免费的午餐。IIR用低阶换来的代价是“相位非线性”和“稳定性风险”。相位非线性这一点很多人一开始不在意但真正用起来就会发现问题。你在示波器上看波形可能觉得没差别可如果处理的是心电信号、振动信号或者任何需要保持波形原始形态的场景IIR会把不同频率成分的延时变得不一样导致输出波形“变形”——这不是幅度上的变形而是时间轴上的错位。FIR因为系数对称设计天然能做到线性相位所有频率延时完全一致波形不会走样。所以选型的时候我一般先问自己一个问题这个应用对相位敏感吗如果只是去除电源纹波、平滑传感器噪声、做音频音调调整IIR完全合适如果要做带通滤波后的时域特征分析、波达时间估计那我建议要么用FIR要么接受IIR非线性相位带来的误差要么用零相位双向滤波离线场景弥补。稳定性风险则是另一个大坑。既然IIR有反馈分母那就天然存在“会不会发散”的问题。从z域看系统稳定的条件是所有极点都落在z平面单位圆内。问题在于当我们把设计好的滤波器搬进16位或32位定点芯片、用有限精度的浮点运算实现时系数的量化误差可能导致原本在单位圆内、离边界很远的极点漂移到圆外系统就变成了一个“振荡器”。这也是后面要重点聊SOS矩阵的直接原因——把高阶系统拆成多个二阶子系统级联能显著降低系数量化误差带来的极点漂移风险让滤波器在真实硬件上更稳定。1.3 选型判断什么时候该用IIR说了那么多我可以把这些年做工程选型的经验整理成一张表帮助快速决策对比项IIRFIR相同指标所需阶数低通常4~10阶高动辄几十上百阶单点计算量小大相位特征非线性可做到线性相位稳定性有极点需关注全零点无条件稳定适合场景实时控制、音频均衡、噪声抑制相位敏感测量、多速率处理、心理声学典型工具巴特沃斯、切比雪夫、椭圆窗函数法、等效波法这张表不是说要大家死记硬背我更想强调的是“实时性”这个维度。在STM32这类MCU上一个中断里可能同时要做ADC采集、滤波、PID控制、通信处理留给滤波器的预算往往只有几个微秒到几十个微秒。FIR要是上千阶光这部分的乘加运算就能把CPU吃满。IIR用个小几十次的乘加就完成了同样指标这种差距是实际测出来的不是纸面数据能体现的。结论很简单高性能计算平台、相位敏感场景选FIR嵌入式实时场景、噪声抑制与平滑场景IIR是性价比极高的选择。而且IIR配合SOS结构稳定性和可维护性都能得到保障这也是它至今在工业控制、音频、传感器处理领域占据重要位置的原因。2. 设计IIR滤波器的完整流程与工具实操2.1 设计规格的确定很多人拿到一个需求就急着打开工具箱敲代码我建议先花十分钟把设计规格写清楚这一步省了后面返工十倍的功夫。所谓设计规格核心就是五个参数采样率fs、通带边缘频率、阻带边缘频率、通带最大纹波、阻带最小衰减。采样率是这一切的基础奈奎斯特频率fs/2决定了你能处理的最高频率分量。在数字滤波器设计里所有频率参数最终都要归一化到这个奈奎斯特频率上。比如采样率1000Hz、想保留100Hz以内的信号那归一化截止频率就是100 / (1000/2) 0.2。这个归一化过程用工具的时候也会遇到只是大多数时候工具帮你做了而已。通带纹波和阻带衰减这两个参数决定了滤波器的“质量”。通带纹波1dB意味着信号通过后在幅度上会有大约±5.6%的波动这在音频应用里通常感知不到但在测量仪器上就可能影响精度。阻带衰减40dB意味着阻带信号能削弱到原来的1%这对大多数工业应用已经够了想做到60dB、80dB也不是不行只是滤波器阶数或者设计算法的复杂度会上升。设计规格不是越严越好跟机械加工的公差一样越严成本越高这里的“成本”就是阶数和运算量。2.2 用Python工具箱完成系数设计我自己最常用的设计工具是Python的SciPy库因为它既能在PC上快速验证又能把设计出的系数直接移植到C代码里。下面是一个实际可运行的低通滤波器设计示例import numpy as np from scipy.signal import ellip, butter, cheby1, sosfreqz fs 1000.0 # 采样率 1kHz fc 100.0 # 通带截止频率 100Hz order 4 # 滤波器阶数 rp 1.0 # 通带纹波 1dB rs 40.0 # 阻带衰减 40dB # 用椭圆滤波器设计输出SOS矩阵格式 sos ellip(order, rp, rs, fc/(fs/2), outputsos) print(sos) # 打印频率响应做验证 w, h sosfreqz(sos, worN2048, fsfs)这里用椭圆滤波器是因为它在同样的阶数下过渡带最窄适合在低阶条件下实现尽可能陡的衰减。如果你对通带内的纹波非常敏感、不想要椭圆滤波器的等纹波纹波特征可以把ellip换成butter巴特沃斯它的通带响应是最平坦的但过渡带会宽一些。切比雪夫I型在通带内有等纹波、阻带单调切比雪夫II型则相反。选型逻辑很简单追求最陡过渡带选椭圆追求通带平坦选巴特沃斯两种极端之间选切比雪夫。2.3 设计结果校验拿到系数之后别急着抄进代码。我习惯先打印频率响应曲线做人眼验证——确认通带内增益接近0dB、截止频率位置正确、阻带衰减达到设计目标。有一回我设计一个带通滤波器系数算出来一切正常结果画幅频图发现中心频率偏了5Hz仔细排查发现是设计时归一化频率的参考搞错了。这个错误要是直接烧进板子可能要调大半天才能发现。SciPy的sosfreqz函数可以直接对SOS格式的系数计算频率响应也可以先用tf2sos把普通的分子分母系数转成SOS。如果做的是离线处理还可以直接用sosfilt函数一次性跑完整段数据来验证滤波效果把输入信号和输出信号放在一起对比看看时域波形是不是符合预期。这个“PC端设计→PC端验证→移植到MCU”的流程我到现在还会用因为MCU上调试滤波的代价比PC上高太多。3. 直接I型滤波器在MCU/STM32上的代码实现3.1 直接I型结构的数学推导拿到系数之后怎么把它变成能在MCU上跑的代码呢有一个概念要先讲清楚差分方程。IIR滤波器的传递函数可以等价转换成时域里的差分方程代码就是逐样本地执行这个差分方程。假设我们有这样一个2阶传递函数H(z) (b0 b1·z⁻¹ b2·z⁻²) / (1 a1·z⁻¹ a2·z⁻²)对应的差分方程为y[n] b0·x[n] b1·x[n-1] b2·x[n-2] - a1·y[n-1] - a2·y[n-2]注意分母里的a1、a2在方程里变成了减法而且之前提到的归一化要求a0 1如果设计工具给出的系数里a0不是1需要先把所有系数除以a0。x[n-1]、x[n-2]是前两个输入样本y[n-1]、y[n-2]是前两个输出样本。实现这个方程的时候只需要维护四个历史变量就够了。这种结构之所以叫“直接I型”是因为它直接照着差分方程逐项实现输入历史先走分子输出历史再走分母。另外还有“直接II型”和“转置直接II型”它们的数学等价但数值特性和内存占用略有差别。在MCU上我通常用直接I型或转置直接II型因为它们的状态变量天然是信号延迟线上的数值便于初始化和调试。3.2 第一个能用单二阶节的C语言实现在MCU上用C写IIR滤波器我建议先做成一个简单的结构体加处理函数不要一开始就写很复杂的级联框架。下面这个例子就是最典型的单二阶节实现也就是俗称的biquad双二阶节typedef struct { float b0, b1, b2; // 分子系数 float a1, a2; // 分母系数 float x1, x2; // 输入历史 float y1, y2; // 输出历史 } Biquad; float biquad_process(Biquad *f, float in) { float out f-b0 * in f-b1 * f-x1 f-b2 * f-x2 - f-a1 * f-y1 - f-a2 * f-y2; // 状态更新 f-x2 f-x1; f-x1 in; f-y2 f-y1; f-y1 out; return out; }这段代码的核心就是差分方程的逐样本执行。我特意把系数a0省略了因为设计工具输出的a0通常已经是1如果遇到a0不等于1的情况一定要先做归一化再填入这个结构体。我在给板子移植代码的时候踩过一次直接从MATLAB导出的系数里a0 1.05没归一化就填进去了结果滤波器增益整体偏移响了好久的困惑才排查到问题。同时要注意结构体里的状态值x1/x2/y1/y2在滤波器运行前应该清零。这个“清状态”的步骤有时候比算法本身还关键——如果不清零上电后可能会有一段未知的瞬态输出尤其在控制回路里这种瞬态可能导致执行机构乱动一下风险不小。3.3 与STM32硬件集成的注意事项在STM32上跑这个过滤器不同的人有不同的集成路径我自己的习惯是放在ADC的中断回调或者DMA传输完成回调里每采集到一个样本就调用一次biquad_process。这样滤波是逐样本实时完成的输出可以直接给到后续的PID、FFT或者显示刷新。实时性和代码效率是两个需要同时考虑的问题。在Cortex-M4及以上内核中CMSIS-DSP库提供了arm_biquad_cascade_df1_f32函数它对多阶IIR滤波做了优化利用SIMD指令可以让多个滤波器并行执行。如果项目里只是简单的一两个二阶节手写的biquad代码就够了没必要引入整个DSP库但如果要处理多通道音频或者高阶数滤波强烈建议换到CMSIS-DSP的级联接口省下的不只是代码量更是实打实的CPU占用率。还有一点关系到ADC数据处理的细节读取ADC值之后通常需要先减去直流偏置再做滤波。要是对原始ADC码值直接滤波那个恒定的直流分量会在滤波结果里保留导致输出一直偏在一个非零基线上。对于需要判断阈值或者计算有效值的应用这个直流偏置会让所有后续判断都出错。我一般在滤波前做一次简单的直流扣除或者输入本来就是交流耦合信号就不用担心。3.4 定点化MCU上无法回避的精度问题很多STM32型号没有FPU或者FPU频率低、不适合大量浮点运算。这时候就要考虑把浮点滤波器转成定点实现最常用的格式是Q15和Q31对应16位和32位定点。定点化的核心思想是把系数和信号都乘以一个缩放因子放大成整数来做乘法累加最后再缩回去。比如Q15格式就是把浮点数乘以32768再取整。但这里有两个必须注意的点系数放大后会引入量化误差。原本0.9999的系数可能量化成0.9997这在反馈回路里可能让极点位置发生微小偏移。阶数高的时候这种微小偏移会累积甚至导致不稳定。这也是我为什么强烈推荐SOS结构的原因——每一级只有2阶极点离单位圆的敏感度低得多。中间乘加的溢出问题。biquad内部有5次乘法和4次加法如果用Q15做乘法结果需要32位来保存中间值。如果你每个数都是Q15那乘出来的积是Q30再加上另一个Q15就需要注意数据宽度。很多人在定点化的时候栽在这里要么丢精度要么溢出。如果MCU没有FPU但又不想写定点我有一条折中建议用float32试试。Cortex-M4以上的内核虽然有FPU但即使没有硬件浮点用软件浮点实现2阶IIR的耗时通常在几十微秒级别对采样率不高比如1kHz的应用完全够用。真正必须定点化的场景是那种采样率几十kHz、又要在中断里干很多活的场合。这时候再耐心做定点不要一上来就把自己绕进位宽地狱。4. SOS矩阵让高阶滤波器稳定落地的关键4.1 高阶直接型的数值稳定性隐患先做个思想实验。想象一个10阶IIR滤波器直接用传递函数分子分母那一堆系数去实现相当于一个10阶的反馈系统。任何微小的系数量化误差都像在一根长竹竿的顶端加重量——竹竿越高顶端轻轻一晃底部就产生巨大的偏差。极点分布对系数误差的敏感度跟阶数成正比阶数越高单位圆附近的极点越容易被推出圆外。这就是为什么高阶IIR滤波器直接用“直接型”结构几乎必出问题的根本原因。你可以在MATLAB里看理论频率响应画得完美无缺但把同样的系数写进32位定点MCU跑出来的可能就是自激振荡的噪声。我用MATLAB做过一个实验一个10阶巴特沃斯低通系数保留6位小数后直接实现极点位置跟原始设计差得不算大但一对共轭极点已经落到了单位圆外系统响应变成增长振荡。这个实验特别直观地说明不是设计的问题是实现结构的问题。4.2 什么是SOS矩阵二阶节的级联SOS是Second-Order Sections的缩写中文通常叫二阶节级联。它的核心思想是把一个N阶的传递函数分解成N/2个2阶系统的级联每个2阶系统用biquad结构实现每个biquad的系数单独归一化。SOS矩阵的每一行长这样[ b0, b1, b2, a0, a1, a2 ]其中a0通常等于1工具会在输出时自动归一化。比如一个4阶滤波器输出的SOS矩阵是2行6列就代表两个级联的biquad。第一个biquad的输出喂给第二个biquad的输入串联完成整个滤波。从数值稳定性角度看这个分解的意义在于每一级只承受2阶的极点敏感度即使系数有量化误差影响也限制在本级不会跨级级联放大。极点从“一根高竹竿的一端”变成了“几截短竹竿的连接”稳定裕度和抗量化能力都大幅提升。此外SOS结构还有一个工程上的好处单位“一节”就是最基本的biquad代码结构统一、调试方便每一个节的输入输出都可以单独断点查看过滤器到底哪一级出了问题一目了然。4.3 从普通系数转SOS以及实际代码结构设计工具通常都支持直接输出SOS。用SciPy的时候在ellip、butter这些函数里指定outputsos就可以了如果手里只有普通的分子分母系数b和a数组也可以用scipy.signal.tf2sos做转换from scipy.signal import butter, tf2sos, sosfilt b, a butter(4, 0.2, outputba) sos tf2sos(b, a) print(sos)在C代码里SOS级联实现其实非常直观——就是把上一节的输出作为下一节的输入逐节调用同一个biquad_process函数。我一般用一个数组保存所有SOS行用循环依次处理typedef struct { uint8_t num_sections; float coeffs[MAX_SECTIONS][6]; float state[MAX_SECTIONS][4]; // x1, x2, y1, y2 } SosFilter; float sos_process(SosFilter *f, float in) { float out in; for (int i 0; i f-num_sections; i) { out biquad_process((Biquad){...}, out); // 使用第i组系数 } return out; }需要注意的一个细节是SOS的级联顺序。级联顺序不是随便排的不同排序会影响数值精度。通常原则是把Q因子较高带宽窄、谐振峰尖锐的节放在前面把增益较大的节放在前面也有助于信噪比。SciPy的sosfilt内部会自动排序但如果你手动处理系数建议参考这个原则。另外SOS各节之间可以穿插分配滤波器总增益避免某一节固定增益特别大导致中间信号饱和——这个问题在定点实现里格外突出。5. 容易踩的坑与排查技巧实录5.1 现象输出爆炸先查极点而不是查算法我第一次把IIR滤波器烧进板子的时候内心是有点慌的——上电后输出直接变成接近满幅度的振荡用示波器一看整个波形都在疯狂抖动。当时第一反应是代码写错了来回检查了好几遍biquad函数的乘加顺序结果都没问题。后来仔细一查是设计的时候用了10阶直接型系数经过去浮点转换之后精密值流失极点跑出了单位圆。这种发散问题的排查套路现在已经固定了先用PC端MATLAB或者Python把固定系数的极点画出来看是不是都在单位圆内系数量化之后再看一次极点位置。如果极点没问题再查代码里是否忘记先归一化a0、是否初始化了状态变量、是否加了错误的负号。这个顺序基本能定位90%以上的发散问题。总结起来就是IIR滤波器只要输出不对劲永远先怀疑数值和极点不要怀疑硬件和主循环这是这个领域最值钱的一条经验。5.2 现象截止频率偏移记得双线性变换的预畸变另一个高频问题是按设计算出来的截止频率是100Hz实测出来变成95Hz或者105Hz而且偏差可能随截止频率升高而变大。在我遇到过的大多数情况下原因出在模拟原型变换到数字域时没有做“频率预畸变”。IIR滤波器通常是从模拟滤波器原型巴特沃斯、切比雪夫、椭圆出发通过双线性变换映射到数字域的。双线性变换会把模拟频率轴的无限范围压缩到数字频率轴上的0到π之间这种频率压缩是非线性的如果不做补偿最终数字滤波器的截止频率就会偏移。很多时候设计工具已经自动做了预畸变但如果你手写转换过程或者从模拟设计表直接查系数就容易踩雷。工具选择错了也会导致看起来“偏移”。比如设计采样率是48kHz结果你按44.1kHz的归一化频率设计了那实际效果当然对不上。如果排除代码问题后频率还是有偏差先检查你的设计工具用的采样率、归一化方式再检查预畸变处理大概率能解决。5.3 现象起始瞬态振荡还有一个看起来吓人但其实很常见的情况滤波器一跑起来输出的前几十个样本会出现很大的过冲感觉像“爆炸”但过了这段之后就恢复正常了。这不是滤波器不稳定而是起始瞬态响应——滤波器状态变量从零开始相当于给滤波器输入了一个阶跃信号系统自然会产生瞬态响应。在实时控制应用里这种瞬态可能让执行机构产生一次意外的运动必须处理。我的做法是在系统启动阶段先运行一段时间的滤波但丢弃这段时间的输出或者把状态初始化到输入信号的平均值附近减少阶跃感。在离线处理场景可以用scipy.signal.filtfilt做零相位滤波但它内部也会处理边缘效应原理是反向再滤波一次相位为零但计算量翻倍。5.4 现象信号经过IIR后幅值不对有时候滤波器跑得“很稳”但仔细对比输入输出某个频段的信号不是大了就是小了。这个现象往往是设计规格和实际需求不匹配造成的而不是代码问题。比如你做的是Butterworth 1dB带宽的低通但在通带边缘1kHz-3dB你却拿设计规格里那条“-1dB频点”跟“-3dB频点”混淆了自然会觉得不对。要规避这类问题最好是拿到滤波器后先把频率响应打印一下在0Hz到奈奎斯特频率范围内画一条完整的响应曲线确认不同频点的增益值符合预期。高通的低频、带通的两侧滚降尤其容易让人产生误解因为不同滤波器定义“通带频率”的方式不同有的按-3dB有的按阻带边缘不做校验就是给自己挖坑。我把这几年在IIR滤波器上遇到过的典型问题整理成一张速查表方便直接对照现象可能原因快速处置输出发散/自激振荡极点跑出单位圆、系数未归一化、定点溢出画极点图检查极点位置先改用双精度浮点验证截止频率偏移未做频率预畸变、设计工具采样率设置错误检查设计参数改用工具自动预畸变起始幅度过冲状态变量初始化为零产生的瞬态启动时丢弃前N个样本或初始化状态到稳态值输出幅度整体偏低/偏高增益分配错误、a0未归一化、通带定义混淆打印频率响应曲线验证0Hz通带增益是否为1定点实现噪声大系数量化误差大、中间累加溢出改用更高的定点位宽使用SOS级联结构5.5 排查技巧从“先画图”到“分段定点”关于IIR的调试我有一个坚持了很多年的习惯任何滤波器改动第一步一定是先画频率响应图第二步是给一个已知信号正弦波、方波、阶跃看时域输出第三步才连接到真实数据源上。三步走下来80%的问题能在接入真实系统前暴露剩下的20%再配合逐级断点、打印状态变量等手段也能快速定位。分段定点这个技巧也分享给大家。定点化之前先用浮点在PC端完整跑一遍系统把每一级SOS节点的输出最大值记录下来。这些值就是定点化的缩放依据——该级要用多大Q值不会溢出、不会浪费精度都有依据了。分段定点比整体定点效果好得多因为不同节之间信号动态范围差异可能很大用统一的Q值意味着低幅度节浪费精度、高幅度节面临溢出风险。写在最后回头看IIR滤波器这条路最让我感慨的一点是很多项目失败不是设计思路错了而是被“差一点点”的系数精度、“差一点点”的实现结构、或者“差一点点”的调试方法拖垮的。IIR本身是个非常成熟的技术从巴特沃斯到SOS每一步都有明确的数学基础和实践路径但真正把这条路走通、走稳还是需要一点一点踩坑换来的手感。我个人现在做滤波器项目的习惯是需求拿过来先写规格然后PC端设计校验确认无误再落到MCU代码全部用SOS结构浮点优先、定点按需最后一定要做一次长时间运行测试。这个流程看起来繁琐但它帮我挡掉了太多半夜调板子的尴尬。希望这篇从设计原理到工程落地的经验梳理能让你在IIR滤波器上少走几条我当年走过的弯路。如果你在STM32或者音频项目里遇到过跟IIR相关的怪问题欢迎沿着文中的思路排查一遍——很多时候答案就藏在那张速查表里。
平方复杂度无法被常规优化消除｜TaoToken 统一 Key 通道实测)
