上手先讲个真实场景。我前阵子给一块用低成本MCU做的采集板调信号质量输入是压力传感器的模拟量经ADC采样后毛刺特别明显高频抖动能把显示数字最后两位带得来回跳。首选的方案自然是经典二阶IIR低通但芯片没有硬件浮点单元软件浮点一次乘加就要好几条指令传感器采样率一提到10kHz以上CPU占用率直接飙到四成多别说再做显示刷新和通信协议连基本任务调度都开始抖。后来我把滤波器整个改成纯整数二阶滤波器乘法全部用固定点整数运算CPU占比一下降到5%以内滤波效果和浮点版本几乎看不出差别。这篇文章就把我从原理、系数推导、定标到最终嵌入式代码的完整过程摊开讲希望能给同样在资源受限平台上做信号处理的朋友一个可以直接抄作业的参考。1. 为什么嵌入式里要死磕纯整数滤波嵌入式工程师听到“浮点”两个字第一反应往往是“贵”。这个贵不是指芯片价格而是指运算代价。带硬件FPU的MCU一次浮点乘加可能只要一两个周期但这类芯片往往要多花不少预算。很多低成本场景用的是不带FPU的Cortex-M0或者老式8位机浮点运算完全靠软件库模拟一次乘法要几十个周期跑起来相当心疼。纯整数二阶级滤波器就能把这个问题绕开。它的所有系数预先算成整数滤波过程里没有一位小数点运算只有整数乘加和移位几个周期就能完成一次输出。尤其适合传感器数据调理、电机电流环反馈、电源环路控制这类实时性要求高、又不方便外挂DSP的场景。还有个容易忽略的优点确定性。浮点运算在不同编译器、不同优化等级下舍入行为可能有细微差异。整数运算的中间结果在所有平台上一模一样这在多机联调、离线回放复现数据时特别有用。我调试某工业项目时就是靠这份确定性在PC上复现了MCU里的每一个滤波输出点排查问题省了不少时间。说白了纯整数滤波器不是“买不起FPU的将就方案”而是一种稳定性、实时性、可移植性都更可控的工程选择。2. 二阶滤波器原理与系数设计2.1 从差分方程说起二阶无限脉冲响应滤波器的离散域差分方程长这样y[n] b0*x[n] b1*x[n-1] b2*x[n-2] - a1*y[n-1] - a2*y[n-2]式子里x是输入序列y是输出序列四个历史量分别代表前两次的输入和输出。b0、b1、b2是前馈系数a1、a2是反馈系数。滤波器的性格——是低通、高通还是带通截止频率多少有没有过冲——全部由这五个系数决定。这里有一个新手容易绕晕的地方符号。差分方程的反馈项写的是减号这跟教科书传递函数H(z)的分母写法是对应的。分母是1 a1*z^-1 a2*z^-2对应差分方程里就是- a1*y[n-1] - a2*y[n-2]。在你把系数存入结构体时最好直接按差分方程的“减”来理解否则程序里符号很容易搞混。我自己早年就吃过这个亏后来统一约定结构体里存的a1、a2已经带符号运算时一律写成减。2.2 双线性变换得到低通系数设计滤波器系数的方法很多双线性变换是嵌入式里最常用的一个。它的思路很简单先写一个模拟低通滤波器的传递函数再做一次从模拟域到数字域的映射。核心公式是s K * (1 - z^-1) / (1 z^-1)其中K tan(π * fc / fs)fc是希望的截止频率fs是采样率。注意这里用tan做预畸变。因为双线性变换在数字角频率和模拟角频率之间存在非线性压缩如果不做预畸变设计出来的实际截止频率会偏离目标值高频偏差尤其明显。以最简单的二阶巴特沃斯低通为例。模拟归一化传递函数的分母是1 2ζs s^2对二阶巴特沃斯来说阻尼比ζ0.7071对应的Q0.7071。代入双线性变换整理后可以得到一组比较干净的系数公式a0 1 2ζK K^2 a1 2*(K^2 - 1) a2 1 - 2ζK K^2 b0 K^2 b1 2*K^2 b2 K^2再把a0除到分母外面让分母首项变成1就得到差分方程里需要的那五个最终系数b0 b0 / a0 b1 b1 / a0 b2 b2 / a0 a1 a1 / a0 a2 a2 / a0顺手算一组。设fc50Hzfs1000Hzζ0.7071。Ktan(π*50/1000)0.1584K^20.02512ζK0.224。于是a01.2491a1-1.9498a20.8011b0b20.0251b10.0502。归一化得到b00.0201b10.0402b20.0201a1-1.5611a20.6414。这组数后面做整数化还会用到。如果你想设计其他响应比如巴特沃斯之外的切比雪夫、贝塞尔只需改变模拟原型里的ζ值或者引入纹波系数公式结构不变所以这套计算流程可以覆盖大多数场景。2.3 为什么二阶就够用一阶低通的滚降只有20dB/十倍频衰减能力偏弱。传感器信号一旦混入宽频噪声一阶滤波后高频部分依然残留不少。二阶巴特沃斯的滚降是40dB/十倍频比一阶强得多实现代价却只多了一个状态变量和几次乘加。三阶、四阶衰减厉害但在嵌入式定点环境下成本上升很快。系数敏感度变高量化误差更容易把极点推到单位圆外。而且很多传感器噪声频谱的主要能量集中在几百赫兹以内二阶低通配合适中的截止频率已经能把信噪比改善几个数量级。这就像盖房子一阶是个棚子二阶是个正经平房三阶四阶是楼房——但平房在大多数情况下已经完全够住。3. 整数化的核心定标、动态范围与溢出控制3.1 把小数系数变成整数上一节算出来的系数全是小数比如b00.0201。MCU处理这种数要么用浮点要么用一个技巧——定标。定标的思路很朴素既然你小数位数太多我干脆把整个差分方程乘以2的N次方把所有系数变成整数算完再右移N位把结果缩回来。差分方程变成y[n] (B0*x[n] B1*x[n-1] B2*x[n-2] - A1*y[n-1] - A2*y[n-2]) N其中B0 (int)(b0 * 2^N) B1 (int)(b1 * 2^N) B2 (int)(b2 * 2^N) A1 (int)(a1 * 2^N) A2 (int)(a2 * 2^N)取整可以选择就近取整或截断我习惯用(int)floor(x*scale 0.5)做就近取整误差更小。上面那组50Hz系数选N12B082B1165B282A1-6394A22627。差分方程就变成y[n] (82*x[n] 165*x[n-1] 82*x[n-2] 6394*y[n-1] - 2627*y[n-2]) 12计算过程里没有任何小数全是int或int64乘法。3.2 定标位数怎么选中间量的动态范围选N不是越大越好也不是越小越好关键看中间累加结果会不会溢出。还是上面那组。假设输入x是16位ADC的满幅输出范围±32768滤波器的输出y也控制在±32768量级。累加器里每一项的绝对值上限是系数绝对值乘32768。最大的是A1对应项约6394乘32768约2.1e8。五项累加最坏情况可能到8e8到1e9。这个量级已经逼近int32的极限2^31约等于2.147e9。所以用32位int做中间累加不是不行但余量不宽。只要输入信号动态范围再大一点或者系数因为截止频率靠近奈奎斯特频率而变大累加器就可能在瞬时极端输入下翻车。我的习惯是直接用int64_t做累加器。嵌入式编译器对64位乘法会生成优化代码虽然比32位慢一点但在低通滤波这种每样本只做五次乘加的场合开销完全可以接受。求稳别省这一个类型。如果实在想用32位累加可以把定标位数降到N10或N8。但系数量化误差会变大滤波器的实际幅频响应会和设计值有偏差。取舍关系如表所示定标位数N系数量化误差中间累加器安全度适用场景8较大截止频率可能偏1%-3%很安全int32随便用对频率精度不敏感的粗滤波12较小适合多数传感器int32需小心推荐int6416位信号常见选择15很小接近浮点精度必须int64高精度信号处理3.3 一个很实用的小验证直流增益定标换算做完之后怎么快速判断系数没算错有个非常实用的检查方法。一个稳定的低通滤波器在直流处的增益应当是1也就是说输入一个恒定不变的值输出最终要收敛到输入值。用整数系数验证起来很简单把五个整数系数全部代入公式dc_gain (B0 B1 B2 - A1 - A2) / 2^N代回上面的例子82165826394-26274096正好等于2^12。这个结果如果等于1乘以定标后等于2^N说明系数换算的直流增益正确。如果不是说明至少某个系数取整错了一位数或者符号写反了。这个办法我每次写滤波器都用两分钟就能排除一大半低级错误。4. 嵌入式代码实现4.1 结构体与核心滤波函数话不多说直接上代码。这里用纯C实现任何嵌入式平台都能用。#include stdint.h typedef struct { /* 定标后的整数系数a1/a2 已带符号 */ int32_t b0; int32_t b1; int32_t b2; int32_t a1; int32_t a2; /* 定标位数对应差分方程的右移量 */ uint32_t shift; /* 历史状态 */ int32_t x1, x2; int32_t y1, y2; /* 输出限幅用可选 */ int32_t min_out; int32_t max_out; } biquad_i_t;核心滤波函数长这样int32_t biquad_i_process(biquad_i_t *f, int32_t x) { int64_t acc; acc (int64_t)f-b0 * x (int64_t)f-b1 * f-x1 (int64_t)f-b2 * f-x2 - (int64_t)f-a1 * f-y1 - (int64_t)f-a2 * f-y2; /* 可选加舍入偏置再右移 */ #if 1 if (acc 0) { acc (int64_t)1 (f-shift - 1); } else { acc - (int64_t)1 (f-shift - 1); } #endif /* 状态更新 */ f-x2 f-x1; f-x1 x; f-y2 f-y1; f-y1 (int32_t)(acc f-shift); /* 可选限幅 */ if (f-y1 f-min_out) { f-y1 f-min_out; } if (f-y1 f-max_out) { f-y1 f-max_out; } return f-y1; }这里有两个细节值得展开说。第一为什么中间累加器用int64_t。看最坏情况乘法项里A1可能到六千多乘上y1的32768单这一项就超过两个亿再累加五项容易突破int32的2.1亿上限。int64_t范围大得多这是最简单可靠的办法。第二右移前为什么要加舍入偏置。acc shift是向下取整对正数来说相当于直接砍掉小数部分。这会造成一个系统性的负偏差尤其在信号接近零的位置会产生可感知的波形台阶。加上偏置后变成就近舍入误差更小。注意我对acc0和acc0分开处理避免负数的算术右移偏置方向搞反。4.2 初始化与硬编码系数初始化有两条路线。路线A初始化函数在运行时接收浮点系数内部完成定标换算。好处是同一套代码可以在不同截止频率间切换坏处是库函数里还得带浮点运算和数学库代码体积增大。路线B把定标后的整数系数直接硬编码写进结构体。这是我最推荐嵌入式场景使用的。系数在PC或主机端用公式算好写死到C文件里目标MCU上完全不需要浮点计算程序连数学库都不必链接。切换滤波器特性时改几个常量重编译即可。以fs1kHz、fc50Hz、Q0.707、定标N12为例const biquad_i_t biquad_lpf_50 { .b0 82, .b1 165, .b2 82, .a1 -6394, .a2 2627, .shift 12, .x1 0, .x2 0, .y1 0, .y2 0, .min_out -32767, .max_out 32767, };上电后将结构体清零或者在每次输入值发生变化前手动resetvoid biquad_i_reset(biquad_i_t *f) { f-x1 0; f-x2 0; f-y1 0; f-y2 0; }注意reset后前几个输出存在过渡过程。如果系统对启动瞬态敏感可以把历史状态初始化成当前输入值比如f-y1 f-y2 f-x1 f-x2 首采样值这样输出从第一个点就基本贴住输入没有爬升过程。4.3 代码再往细里抠一抠类型转换的坑值得多说一句。(int64_t)f-b0 * x为什么要显式把b0转成int64_t再乘因为b0是int32x也是int32C语言会把两者都提升为int32做乘法结果溢出后才赋值给int64到时候已经晚了。写成(int64_t)f-b0 * x强制让乘法在64位下进行是必须的。状态更新顺序也有讲究。必须先更新x2和y2再更新x1和y1顺序搞反会串位导致输出变成另一个滤波器。我自己用代码生成器批量产生滤波器时专门给这个写过静态断言。5. 实战给压力传感器信号去噪5.1 参数选型思路场景是一块采集板压力传感器输出0.5V到4.5V经12位ADC采样采样率fs1000Hz信号里主要噪声来自工频耦合和高频开关毛刺目标是把噪声压下去同时保留200Hz以内有用的脉动波形。截止频率选择要遵循一个原则截止频率要比有用信号最高频率高但又要比噪声主频低。这里有用信号脉动大概在几十赫兹工频干扰是50Hz附近但不算太强因此我把fc定在80Hz比较合适。当然这也意味着50Hz工频还是会漏过来一些想要更狠的抑制需要加陷波器那是另一个话题。Q值这里选0.7071即巴特沃斯响应通带最平坦没有明显过冲。如果后面跟的环节怕过冲可以选Q0.5也就是贝塞尔响应群延迟更平但下降更缓。对压力信号巴特沃斯是最稳的起点。5.2 系数换算全过程还是上面的公式先把fc80Hz、fs1000Hz代入Ktan(π*80/1000)tan(0.2513)0.2572。K^20.06622ζK1.4142*0.25720.3637。a010.36370.06621.4299 a12*(0.0662-1)-1.8676 a21-0.36370.06620.7025 b0b20.0662 b12*0.06620.1324归一化b00.0463 b10.0926 b20.0463 a1-1.3061 a20.4913定标N12取整B0190 B1379 B2190 A1-5350 A22012检查直流增益1903791905350-20124097相比2^124096偏差1。这个偏差来自系数取整的累积误差完全可以忽略实际直流增益为0.9998。写进C代码const biquad_i_t biquad_lpf_80 { .b0 190, .b1 379, .b2 190, .a1 -5350, .a2 2012, .shift 12, };这样整个滤波器就完成了。整个过程不需要在MCU上做任何浮点运算所有设计计算都在PC上完成。5.3 实测效果与波形表现实际跑起来先给一个1kHz采样率的伪方波叠加约满幅5%的白噪声。滤波后明显看到高频毛刺被压平上升沿还能保持锐利没有明显的高频损失。再给一个固定直流输入输出大约在几十个采样点内稳定下来稳定后波动在±1个LSB以内Q12定标下这就是16位整数量化的理论地板。对比浮点仿真结果两者在主波形上几乎不可区分最大瞬时误差在输出满幅的0.1%以内。对绝大多数传感器应用来说这个精度足够了。如果你的项目要求更高精度把定标位数提到N15同时所有中间累加保持int64量化误差还能再降一个数量级。6. 常见问题与调优经验6.1 问题速查表现象原因排查与解决输出缓慢漂移或不稳定系数符号写反a1/a2存储与差分方程不一致用直流增益验证公式检查确认B0B1B2-A1-A2近似等于2^N输出前几十点严重偏小启动时历史状态为零滤波器需要建立时间用首个采样值初始化全部状态或丢弃前N个输出样本突然出现大幅跳变中间累加器溢出通常是int32累加且输入接近满幅换成int64_t累加器或降低定标位数N固定输入时输出周期性小幅振荡定点滤波器的极限环现象提高定标位数使用舍入而非截断必要时加1LSB抖动滤波后仍有很多高频残留截止频率选太高或阶数不够降低fc检查是否真的需要更高阶滤波器输出与浮点仿真差别明显定标位数太低系数量化误差大提高N到15状态更新保持int646.2 溢出的深水区什么时候必须上64位一个容易被忽略的溢出场景是y1、y2不再局限于16位范围。比如你在做闭环控制输出作为误差量叠加到控制器y的幅值可能超过ADC量程这时候A1*y1单一项就能轻松到几十亿。另一种场景是滤波器串联。两个二阶滤波器的中间量如果都放在16位很可能截断正确做法是每个滤波器的输出仍然保留int32状态变量全部int32累加器一律int64。更隐蔽的陷阱输入信号本身如果来自32位累加器如定时器计数幅值范围可能是±2^31。这时候无论定标N取多大累加器都必须使用int64否则输入那一步就先爆了。我一般这样判断——写代码之前先估算max_acc (|B0||B1||B2||A1||A2|) * max_amplitude这个值如果明显小于2^31int32可行接近或超过直接int64。五项的绝对值之和我算过很多常见配置下大概在2e9附近刚好逼近int32的边界所以保守派直接int64。6.3 量化误差、极限环与稳定性量化误差的来源有两个一是系数的取整误差二是输出每次右移时丢掉的小数部分。系数取整误差可以靠提高N来减小但会挤占累加器的动态范围。输出量化误差则与舍入策略强相关我测试过截断与就近舍入在中低频段差别不大但在输入信号幅度较小时就近舍入的波形明显干净得多。极限环是个很经典的现象输入固定不变时理想滤波器输出应该收敛到固定值但定点实现可能出现几个LSB的周期性振荡。根因是状态变量的小数部分丢失使得系统在平衡点附近形成了一个不收敛的小环。缓解办法不多一是提高定标位数二是在状态更新时用舍入而不是截断三是在系统层面接受这个误差因为几个LSB的振荡在大多数控制回路里都在死区范围内。稳定性的检查我习惯给滤波器喂单位阶跃信号看输出是否收敛。判断条件很简单|A2| 2^N且|A1| 2^N A2时要小心整数环境下靠近边界时极容易越过。如果你的系数经过取整后触发了这个判断稳妥起见可以把极点往圆心方向挪一点即把A1、A2的绝对值各减小万分之几换来的是稳定余量。6.4 一些个人调试心得第一永远先测直流增益。每次改了系数先用恒定输入跑一遍输出如果收敛不到输入值系数换算或者定标一定有问题别急着接传感器找原因。第二做一个扫频测试。嵌入式上不方便用标准库的sin可以预先算一张正弦波表频率从1Hz到400Hz步进扫记录每个频率点输出的幅值画出来和理想幅频曲线对比。这个测试能同时发现系数错误、极限环、截止频率偏移三个问题性价比很高。第三把滤波器的中间变量用调试器全部拉出来加入一个只用于测试的峰值记录。跑一轮后检查峰值有没有逼近int64边界这比事后看输出波形更早发现隐患。第四代码里尽可能把shift做成编译期常量不要留成运行期变量。编译器遇到常量移位会生成更紧凑的指令我实测大概能省4到6个周期每样本在高速采样里这个收益不小。我自己现在做MCU上的信号调理第一反应就是列一列能不能用纯整数搞定而不是急着挂浮点库。这个思路从低通滤波器出发扩展到高通、带通、带阻也只是改系数的问题整套框架可以一直复用。如果你在项目中试出了更极端的定标方案或者踩到了我没提到的坑欢迎交流这种问题往往聊一聊就能少走一大段弯路。