开场别被老技术这三个字骗了IIR滤波器无限脉冲响应滤波器大学数字信号处理课上最劝退的那一章也是我工作这些年里用得最多的滤波器。你可能觉得FIR才是万金油什么场合都能套一个窗函数上去但我告诉你在很多资源受限的嵌入式场景里IIR才是真正的保命方案。为什么一句话同样的滤波效果IIR用的阶数更低算得快省内存。一个3阶的巴特沃斯低通性能大致能抵得上十几阶甚至几十阶的FIR。对跑在STM32这种Cortex-M内核上的程序来说这差别不是一点点是实实在在的算力开销和RAM占用。这篇文章我不会从Z变换的严格定义开始讲那套东西教材里多得是。我按自己实际用下来的思路来先搞明白IIR到底是个什么东西它和FIR的核心差异在哪然后讲清楚SOS矩阵和直接I型这两种最常用的实现方式再落到STM32上系数怎么算、代码怎么写、怎么避开那些害死人的坑。读完你能直接上手至少不会再用错结构、算错系数。1. IIR滤波器的本质输出不仅取决于输入还取决于过去的输出1.1 从差分方程看IIR的反馈本质FIR滤波器的差分方程长这样y[n] b0·x[n] b1·x[n-1] ... bN·x[n-N]注意输出只和当前及过去的输入有关没有输出反馈。所以一个单位脉冲进去经过N个采样点之后输出就归零了脉冲响应是有限长的这就是Finite Impulse Response名字的由来。IIR滤波器的差分方程多了一项y[n] b0·x[n] b1·x[n-1] ... bM·x[n-M] - a1·y[n-1] - a2·y[n-2] - ... - aN·y[n-N]后面这一串带a系数的项是把过去的输出又加权加回来这就是反馈。因为输出会不断反馈到输入端所以一个单位脉冲进去理论上输出永远不会完全归零虽然实际上因为有限精度会逐渐衰减到零脉冲响应是无限长的这就是Infinite Impulse Response。这个反馈是IIR的灵魂也是它所有优缺点的根源。反馈让同样的滤波效果只需要很少的系数——阶数低、运算量小、内存占用小。但反馈也带来了稳定性问题如果反馈系数不合适输出可能发散直接飘到天上去。FIR是绝对稳定的只要系数是有限值IIR则必须检查极点位置。1.2 频域视角IIR的优势从哪来IIR滤波器在频域上可以做到极陡的过渡带。比如你有一个50Hz的工频干扰旁边就是你要的100Hz信号你用FIR想把这个干扰压下去可能需要60阶甚至80阶但用IIR2阶到4阶就差不多能做到。原因是IIR的系统函数可以写成H(z) B(z) / A(z)分子多项式B(z)决定零点分母多项式A(z)决定极点。FIR只有分子相当于只有零点IIR既有多项式又有分母可以实现更复杂的频率形状。极点可以看成是让滤波器在某些频率上谐振从而实现高增益的窄带特性或者急剧变化的相位特征。用个生活化的类比FIR是纯粹靠记住更多历史数据来平滑结果的算法像一个记性好但反应慢的人IIR是记得住过去的结果并据此调整当前判断的算法像一个经验丰富但偶尔会一根筋不稳定的老手。1.3 一个看得见的例子同参数下直接对比我用一个具体的例子来说明。假设采样率Fs 1000Hz截止频率Fc 50Hz的低通滤波器MATLAB/Octave里用同样的设计规格来对比。如果用切比雪夫I型IIR3阶就能做到通带波纹0.5dB、阻带衰减40dB。但如果用FIR要达到差不多的过渡带宽度和阻带衰减用窗函数法设计估算需要的阶数大约是N ≈ (Attenuation_dB - 8) / (2.285 × Δω)其中Δω是归一化过渡带宽度。算下来至少需要20多阶。在STM32F103这种72MHz的MCU上每秒钟采样1000次每个采样点要算20次乘加运算IIR只需要算6次乘加。高下立判。2. 从传递函数到实际代码IIR的三种常见实现结构这一节是实操的基石。很多新手直接拿高阶IIR系数往代码里一塞发现输出全是NaN根本不知道自己错在哪。问题往往出在实现结构上。2.1 直接I型Direct Form I最直观但也最浪费直接I型是最容易理解的实现方式就是把差分方程照抄成代码y[n] b0·x[n] b1·x[n-1] ... bM·x[n-M] - a1·y[n-1] - ... - aN·y[n-N]在代码里你需要维护两个数组一个是输入的历史缓冲x_buffer一个是输出的历史缓冲y_buffer。处理每个新样本的伪代码如下// Direct Form I 伪代码 y b0*x b1*x_buf[0] b2*x_buf[1] - a1*y_buf[0] - a2*y_buf[1]; // 更新缓冲 x_buf[1] x_buf[0]; x_buf[0] x; y_buf[1] y_buf[0]; y_buf[0] y;这样有两个历史缓冲各占M和N个位置。看起来逻辑简单但要注意两个问题第一存储浪费。虽然IIR阶数低但每个滤波器都要维护两个缓冲。如果是多通道的音频处理或者多路传感器信号内存消耗会成倍增加。第二数值灵敏度。直接I型在系数取值范围很大时中间结果可能非常大在定点DSP或单片机上很容易溢出。这是它最大的问题。但直接I型也有好处它对系数量化误差的敏感度相对较低相比直接II型来说而且移植简单、容易查错。STM32这类带FPU的MCU上用浮点运算实现直接I型性能基本不是问题出不了大乱子。2.2 直接II型Direct Form II省内存但更矫情直接II型也叫Canonical Form是一种更节省内存的结构。它首先计算一个中间变量w[n]w[n] x[n] - a1·w[n-1] - ... - aN·w[n-N]然后再计算输出y[n] b0·w[n] b1·w[n-1] ... bM·w[n-M]这样做的好处是你只需要维护一组状态变量w_buffer而不是两组内存开销减少了一半。在很多老的DSP课程里它被大大推崇但实际工程中用得反而少原因在于它对系数误差更敏感。特别是在把高阶滤波器系数直接量化成16位定点数时直接II型的极点位置偏移可能比直接I型更大稍不留神滤波器就不符合设计规格了。2.3 级联型SOS高阶滤波器的正确打开方式如果你搜索过iir滤波器 sos 矩阵你看到的SOSSecond-Order Sections就是级联型的标准格式。它的核心思想是不把高阶IIR滤波器作为一个整体去实现而是把它拆成多个二阶滤波器串行级联起来。为什么这么做因为高阶多项式在数值计算上非常脆弱。假设你设计了一个8阶的巴特沃斯滤波器分母多项式A(z)有8个系数这些系数的动态范围可能非常大小到10的负几次方大到几百。在浮点运算中可能还好但在定点处理器上量化误差会迅速放大极点的实际位置和理论位置偏差很大滤波器可能变成振荡器——输出持续抖动甚至发散。如果把8阶滤波器拆成4个二阶节每个二阶节的系数范围小得多、数值稳定性好得多依次计算每一级的输出是下一级的输入整体效果不变但数值表现优秀得多。SOS矩阵的格式通常是这样的每一行表示一个二阶节[b0, b1, b2, 1, a1, a2]举个例子一个4阶巴特沃斯低通滤波器Fs1000HzFc50Hz在MATLAB/Octave中设计后用sos函数输出的就是一个Nx6的矩阵sos 0.0201, 0.0402, 0.0201, 1.0000, -1.5606, 0.6414 1.0000, 2.0000, 1.0000, 1.0000, -1.3643, 0.5098每一行代表一个二阶节。第一行增益较低b系数小第二行增益较高b系数接近1。两级串联在一起整体的传递函数就是这两个二阶节的乘积。级联SOS是实际工程中最推荐的选择。TI的DSP库、CMSIS-DSP里的arm_biquad_cascade_df1_f32函数以及各种CMSIS滤波器库都是基于SOS级联结构实现的。3. 用SOS矩阵串联还是直接算代码实现的关键区别现在问题来了已知SOS矩阵怎么在代码里实现3.1 逐级处理的实现方法最直接的方法是按顺序处理每一级。注意每一级的输出就是下一级的输入处理完第一级得到中间信号再传给第二级。伪代码如下// 假设有numSections个二阶节sos是numSections x 6的系数矩阵 // x_in是当前采样值state是numSections x 4的状态缓冲 float filter_process(float x_in) { float y x_in; for (int i 0; i numSections; i) { y biquad_process(y, sos[i], state[i]); } return y; }其中biquad_process输入一个值输出这个二阶节的输出。它的内部实现直接I型二阶节float biquad_process(float x) { float y b0*x b1*state[0] b2*state[1] - a1*state[2] - a2*state[3]; state[1] state[0]; state[0] x; state[3] state[2]; state[2] y; return y; }这就是CMSIS-DSP中arm_biquad_cascade_df1_f32函数的思路。这种逐级处理的方法也有额外的调试好处你可以在任意一级输出处加打印或者断点确认是哪一级出了问题排错方便很多。3.2 直接实现一个塞满系数的高阶滤波器为什么危险有人会想既然我有8阶滤波器的全部系数为什么不直接把差分方程里所有a和b系数写进代码理论上可以但实际工程中我强烈不建议。原因还是之前提到的数值稳定性问题。用MATLAB/Octave设计一个8阶切比雪夫II型滤波器它的极点和零点可能非常接近单位圆任何微小的量化误差都可能把极点推出单位圆外导致滤波器不稳定。就算不推到单位圆外极点位置偏移也会让实际频响曲线和设计规格差异很大可能原来要求-40dB处实际上只有-32dB。还有一个隐患如果输入的x信号太大中间各级的输出可能很大在级联结构里你可以在每一级之间重新归一化很多库会自动做但直接实现的高阶结构没有这个灵活度。如果中间变量超出浮点数范围直接就NaN了。3.3 从MATLAB/Octave导出SOS矩阵的操作流程这段操作值得仔细看因为它是搜索iir滤波器 sos 矩阵的人最想找到的答案。在MATLAB/Octave中的标准流程是设计原型滤波器比如一个5阶的巴特沃斯低通采样率Fs1000Hz截止频率50Hz[z, p, k] butter(5, 50/(1000/2), low);注意50/(1000/2)是归一化截止频率数字信号处理中频率必须用奈奎斯特频率Fs/2归一化这是最容易搞错的地方。把零极点增益模型转换成SOS矩阵[sos, g] zp2sos(z, p, k);sos就是二级节矩阵g是全局增益。这里有一个重要技巧zp2sos可以带参数up或down指定级联顺序通常用up把最靠近单位圆的极点放在最后数值稳定性更好[sos, g] zp2sos(z, p, k, up);把全局增益分配进第一级和第二级。许多库的biquad结构自带增益系数每个节都有b0/b1/b2所以你可以把增益g乘到第一级的b系数上也可以平均分配到所有级。为了缩小每级的动态范围通常建议g乘到第一级sos(1, 1:3) sos(1, 1:3) * g;这样后面每一节都不需要额外乘增益。但如果g特别大或特别小建议在级间观察信号幅度必要时手动调整增益分配。如果要在STM32上用定点实现还需要把系数转换成Q格式。比如用16位Q15表示系数取值范围-1到1之间但要注意IIR系数的负数范围可能超过-1这需要额外处理。通常建议直接用浮点STM32F4以上带FPU或者用32位定点而不是用16位。后面我会详细讲这个问题。4. FIR和IIR到底怎么选从几个实际工程场景看4.1 一张表站在需求角度做对比先给一个实用对比表这是我做选型时必看的对比维度IIRFIR相同滤波效果的阶数低通常2~8阶高通常几十阶计算量小大内存占用小大线性相位不支持相位非线性天然支持关于中心对称稳定性需要检查极点始终稳定数值敏感性高需小心实现低适合场景资源受限、实时性要求高的嵌入式音频处理、需要无相位失真的数据采集表格里最关键的指标是线性相位。4.2 相位敏感场景FIR胜出如果你的应用是对信号做高精度测量比如振动分析、心电信号处理、音频滤波波形的时域形状很重要那么IIR的非线性相位可能是个问题。IIR滤波器在某些频率点会有较大的群延迟波动信号经过滤波器后不同频率成分的延迟不同时域波形会发生畸变。举个例子你用IIR滤波器处理心电信号它的P波、QRS波群、T波含有不同的频率成分IIR滤波器会导致这些波的相对时间关系发生偏移医生一看波形就觉得不对。你以为滤波后信号变干净了实际上它已经被扭曲了。这种情况下FIR的线性相位特性就是刚需。但如果你只是从传感器数据里滤掉一些噪声不关心波形的具体相位比如做温控的PID控制里滤掉高频抖动、电池管理系统里读取电流电压的平均值那IIR完全没有问题。4.3 资源有限场景IIR胜出这又要说回STM32了。假设你用STM32F103主频72MHz没有FPU只有单精度浮点仿真库做一个电机电流环的采样滤波。电流环的采样频率可能要到10kHz甚至20kHz。如果你用FIR滤波器阶数50每次采样要做50次乘加运算在无FPU的芯片上这已经是很重的开销了。但如果用IIR 2阶滤波器每次采样只需要大约6次乘加运算少了一个数量级。更关键的是状态缓冲只需要4个浮点数而FIR需要50个。所以在硬实时系统中IIR往往是唯一的合理选择。你可能觉得反正计算量也不大但你要考虑整个系统的预算——中断里除了滤波还有PID计算、通信协议处理、状态机判断。每多1微秒的开销在高速控制回路里都是要命的。5. STM32上实现IIR从系数获取到实测这是搜索stm32 iir滤波器直接i型 系数的人最关心的部分我直接按完整流程来。5.1 在STM32上用浮点直接I型实现2阶IIR先提供一个可以直接用的函数这是最经典的2阶直接I型实现系数从MATLAB/Octave导出后手动填入typedef struct { float b0, b1, b2; float a1, a2; float x1, x2; // 输入历史 float y1, y2; // 输出历史 } iir_biquad_t; float iir_biquad_process(iir_biquad_t *f, float x) { float y f-b0 * x f-b1 * f-x1 f-b2 * f-x2 - f-a1 * f-y1 - f-a2 * f-y2; // 更新历史 f-x2 f-x1; f-x1 x; f-y2 f-y1; f-y1 y; return y; }注意这里的符号约定在MATLAB中滤波器系数形式一般是y[n] b0*x[n] ... - a1*y[n-1] - ...所以代码里面减法要对应好。很多人就是在这里搞错了符号导致滤波器完全不对。5.2 从MATLAB/Octave获取系数并检查以STC/STM32项目中常用到的一个50Hz工频陷波器为例采样率Fs500Hz希望滤除50Hz干扰。可以用iirnotch函数Fs 500; fo 50; bw 5; % 3dB带宽 [b, a] iirnotch(fo/(Fs/2), bw/(Fs/2));得到类似这样的系数b [0.97551, -1.10060, 0.97551] a [1.00000, -1.10060, 0.95102]注意这里a的第一个元素是1代码里不需要用它但需要确认。把b0、b1、b2、a1、a2分别填到结构体的对应字段里就行。在把系数烧到单片机之前强烈建议先在电脑上做一次快速仿真验证。可以在MATLAB/Octave里用freqz(b, a, 1024, Fs)看频响曲线确认陷波点位置和带宽正确再用stepz或impz看时域响应是否收敛。这一步花5分钟能省去后面在板子上调试几小时。5.3 定点还是浮点STM32上怎么选STM32F0、F1这类不带FPU的芯片用浮点数做实时的IIR是有代价的。虽然编译器有软件浮点库但乘法运算会扩展到几十条汇编指令速度慢。这种情况下有两种方案方案一用Q15或Q31格式实现在STM32上不太推荐自己搞。STM32官方库里的arm_biquad_cascade_df1_f32是单精度浮点版本需要用带FPU的Cortex-M4/M7/M33。CMSIS-DSP里也有arm_biquad_cascade_df1_q15和arm_biquad_cascade_df1_q31定点版本可以直接用省心很多。方案二如果你必须用STM32F1系列说实话我建议先算算实际负载再决定。一个2阶IIR每秒钟跑1000次软件浮点也就需要几百微秒对于慢速采样场景其实无所谓。只有当采样率很高或多路滤波时才需要考虑软件浮点性能瓶颈。方案三用定点手动实现系数非常不建议初学者做因为会遇到饱和、舍入误差、极限环振荡等问题。我个人的偏好是如果成本允许直接上STM32G4或STM32F4带FPU用单精度浮点实现IIR开发速度和调试体验远超定点方案价格差距也就几块钱。5.4 一个简易测试方法板子上怎么确认滤波器没写错代码写完了板子跑起来了怎么确认滤波器正常工作我的做法是在滤波函数入口处注入一个已知信号观察输出是否和MATLAB对同一信号的仿真结果一致。具体操作在单片机里生成一个固定频率的正弦波数组比如放在Flash里的const数组或者直接在代码里用简单的sin函数生成。给定一个包含50Hz低频和300Hz高频混合的信号经过滤波器后观察输出波形是否只留下低频部分。把UF问题在输出端加一个调试变量通过串口或者JTAG调试器把数据点抓出来画波形。如果输出波形形状和MATLAB仿真基本一致幅值允许有微小误差滤波器就通了。这一步别跳曾经我身边不止一个人跳过这一步直接接到实际传感器上结果滤波器系数里符号弄错折腾了一天最后发现低频全被滤掉了。6. 实际踩坑记录IIR滤波器工程化的五个教训最后这部分是我自己踩坑踩出来的每一件都付出过时间成本。6.1 初始瞬态问题滤波器的启动冲击IIR滤波器有反馈所以它有一个建立时间。假设初始化时状态变量全部清零突然输入一个大的阶跃信号输出会有一个很大的过冲可能持续几毫秒甚至更长。在很多控制系统中这个过冲会导致执行机构猛地动一下这是不可接受的。我的处理办法是系统启动时先不要马上让滤波器输出控制量而是先让传感器信号稳定一段时间或者在初始化时把滤波器的输入历史状态赋值为当前采样值的合理估计让滤波器从工作点开始。比如传感器刚上电时输出接近0就把x1、x2、y1、y2都初始化成0如果预期传感器稳定在某个偏置电压你可以先采N个样本求平均然后设置初始状态。6.2 系数量化误差16位定点的灾难我第一次在STM32F103上做IIR时天真地把浮点系数直接转成Q15格式结果滤波器特性完全乱套。原因很简单IIR系数的动态范围经常超过Q15能表示的范围比如系数可能在0.998到1.002之间Q15量化后直接变成一个小整数极点的精度损失太严重。后来我改用浮点加上CMSIS-DSP的库问题就消失了。如果确实只能用定点至少要选Q31并且考虑用级联SOS结构。16位定点做IIR尤其是高阶IIR基本等于自找麻烦。6.3 采样率不稳定的影响IIR系数是基于固定采样率设计的。如果一个系统用定时器中断触发采样但是定时器优先级设置不当导致采样周期抖动滤波器实际表现会和设计差很多。最典型的例子是你用HAL_Delay或者简单的while循环来做ADC采样触发时钟不准确采样率漂移你以为是5kHz采样并计算截止频率实际可能是4.7kHz滤波器特性偏移还不算大但如果触发被其他中断打断采样间隔忽长忽短IIR的反馈会累积误差输出出现抖动。解决办法是用硬件定时器触发ADC采样或者用DMA保证采样间隔精确。这是嵌入式工程师做任何数字滤波都必须刻进DNA的规范。6.4 状态变量精度float还是doubleCortex-M4的FPU有单精度浮点速度极快但精度只有大约7位有效数字。对于2阶IIR来说这精度足够。但对于8阶以上的高阶IIR单浮点可能不够用状态变量累积误差会导致输出噪声增加。我的建议8阶以下放心用float更高阶建议用double或者保证SOS级联顺序合理。还有如果你需要极窄带滤波器比如Q值很高的带通或陷波器float可能不够因为极点离单位圆太近对精度极其敏感。我的经验阈值是如果滤波器的极点模值超过0.999用double。6.5 滤波器稳定性必须在真实运行范围内验证很多滤波器的极点位置是在设计阶段算好的但实际运行过程中信号频率和幅度的变化可能导致滤波器中间变量过大出现饱和等非线性现象使得系统不稳定。特别是定点实现中饱和溢出后状态变量变成极大值反馈恢复不过来滤波器就卡在抱死状态。所以稳定性验证不能只看理论要加到最恶劣的输入信号进行测试比如信号满幅值、接近截止频率的阶跃信号。如果出现异常首先检查各节输出是否触及数据范围上界必要时加入饱和保护逻辑。结尾一些关于IIR的个人心得IIR滤波器不是银弹但它绝对是嵌入式工程师和算法工程师工具箱里不可缺少的一把锤子。我用它做过电机电流环的陷波、心率传感器的噪声抑制、电源纹波数据的平滑、音频回声消除的前置滤波每个场景都踩过不同的坑但核心方法论是一样的先算清需求指标选对阶数和类型用MATLAB/Octave设计并用SOS级联实现最后在板子上用已知信号验证。如果你想从IIR开始入手我建议先从2阶巴特沃斯低通开始把差分方程、直接I型实现、频响验证这套链路走通再尝试高阶的切比雪夫、椭圆滤波器。在这个过程中要特别留意那个符号约定的问题——MATLAB里a系数带着负号出现在差分方程里很多初学者在这里翻车不是少数。最后分享一个小技巧如果你的系统里信号频率范围比较宽又担心IIR的相位失真影响波形可以在信号链路里先做一次正向滤波再做一次反向滤波所谓零相位滤波但这只适合离线数据处理实时系统里就别想了。实时系统的话要么接受IIR的相位特性要么踏踏实实上FIR。搞清楚自己的需求边界比纠结哪个技术更高级重要得多。