卡尔曼滤波,这四个字在控制、导航、信号处理圈子里几乎是个绕不开的坎。我最早正儿八经去碰它,是在做车辆纵向车速估计的时候,手里只有轮速传感器和加速度计,想要一个靠谱的车速信号,一阶低通滤波压不住噪声,直接积分又怕漂移,最后把希望全押在卡尔曼滤波上。模型搭在Simulink里,从建立系统方程、定状态转移矩阵,到调Q和R参数、对比波形,整个过程踩了不少坑,也把很多教科书上讲得模模糊糊的细节彻底理顺了。这篇文章就把整套东西重新梳理一遍:卡尔曼滤波到底在算什么、五个核心公式怎么理解,以及在Simulink里怎么把系统模型真正搭起来跑通。无论你是刚接触卡尔曼滤波的在校学生,还是要在工程里做状态估计的工程师,这套思路和模型都能直接拿过去改改用。
1. 卡尔曼滤波到底在解决什么问题
1.1 从一阶低通滤波说起:为什么它不够用
很多人第一次听到"滤波"两个字,第一反应就是低通滤波,尤其是一阶RC低通,公式简单、Simulink里一个Transfer Fcn模块就解决了。但实际做状态估计的时候,一阶低通有两个很致命的毛病。
第一个是滞后。只要是低通,就必然有相位滞后,输入信号一变,输出总要慢半拍。我之前做车速估计的时候,用一阶低通滤轮速信号,正常巡航没什么问题,但一到急加速或者重刹,滤波后的车速明显跟不上真实车速的变化趋势,这一慢,后面ESP或者扭矩控制逻辑就会误判。第二个问题是它完全不利用系统模型的信息。一阶低通只根据当前测量值做平滑,它不知道系统本身的运动规律,也不知道控制输入(比如电机扭矩、加速度)对状态的影响。换句话说,低通滤波是一个"盲"滤波器。
卡尔曼滤波的思路完全不一样。它把系统的状态方程拿进来,先根据上一时刻的状态和控制输入去预测当前状态,再用当前测量值去修正这个预测。每一步都同时利用"模型规律"和"观测数据",所以它既能压制噪声,又不会像低通那样产生明显的滞后,这就是它被称为"最优估计"的底气所在。
1.2 核心数学思想:预测与校正的闭环
卡尔曼滤波的数学思想可以浓缩成一句话:把不确定性问题拆成两步,先按物理规律往前走一步,再用观测数据拉回来一步,然后反复迭代。这个过程像极了你在一个陌生的城市里导航:手机GPS告诉你位置,但GPS有漂移;你手里有地图和你的运动模型,知道往前走一百米大概会到哪。你把模型预测的位置和GPS测到的位置按各自的置信度加权合并,得到一个比两者都更准的位置。
对应到离散时间系统里,系统的状态方程是:
x(k) = A * x(k-1) + B * u(k-1) + w(k-1)
观测方程是:
z(k) = H * x(k) + v(k)
其中w是过程噪声,v是测量噪声,它们都被假设为零均值的高斯白噪声。A是状态转移矩阵,B是控制输入矩阵,H是观测矩阵。卡尔曼滤波不直接处理原始测量信号,而是处理"状态估计值",也就是从带噪声的观测量里把真实状态一步步递推出来。
很多人一上来就被矩阵和高斯分布吓住了,其实抛开这些形式,核心就是:预测步用模型给你一个先验估计,更新步用观测给你一个后验估计,而两者之间的权重由卡尔曼增益K来动态分配。噪声大的系统,K就倾向于模型预测;模型不可靠,K就倾向于测量值。
1.3 五个核心公式一次讲透
卡尔曼滤波的计算流程可以压缩为五个公式,这也是我建议所有人在搭Simulink模型之前手推一遍的内容。
预测步:
x_pred = A * x_est + B * u
P_pred = A * P_est * A^T + Q
更新步:
K = P_pred * H^T * (H * P_pred * H^T + R)^(-1)
x_est = x_pred + K * (z - H * x_pred)
P_est = (I - K * H) * P_pred
每个量的含义要理清楚:P是状态估计的误差协方差矩阵,它衡量的是"我对当前状态估计到底有多信任";Q是过程噪声协方差矩阵,表示模型本身的不确定性;R是测量噪声协方差矩阵,表示传感器测量值的不确定性;K是卡尔曼增益,它决定了在预测值和测量值之间各信多少。P越小说明估计越自信,K就越小,这时更新步给的修正量也小;R越大说明测量越不可信,K同样越小。换句话说,卡尔曼滤波的"智能"全在这个K的自动调节上,它不是一个固定参数,而是随着每一步的P和R实时算出来的。
这五个公式的输入输出关系可以用一个很直观的比喻来理解:P_pred是"我对预测结果的把握",x_pred是"我觉得现在应该在哪",K是"在多大程度上相信外部观测"。(z - H*x_pred)叫新息,也就是"实际观测和我预测的差异有多大",差异越大,如果K也大,修正量就越大。整个滤波器就这么一圈一圈滚下去,直到P收敛到一个稳态值。
2. Simulink建模前的状态方程与参数准备
2.1 系统模型怎么定:以纵向车速估计为例
要在Simulink里搭卡尔曼滤波模型,第一步不是画模块,而是把状态方程和观测方程写清楚。拿我熟悉的车辆纵向车速估计举例,假设整车可以简化为一维运动模型,状态量取车速v,控制输入取加速度a,那么连续时间状态方程就是:
dv/dt = a + w
其中w代表未建模的加速度扰动。离散化之后,取采样周期Ts,可以得到:
v(k) = v(k-1) + Ts * a(k-1) + w(k-1)
也就是说状态转移矩阵A=1,控制输入矩阵B=Ts。观测来源用轮速传感器,设轮速直接测车速,那么H=1,观测方程就是:
z(k) = v(k) + v_k
这里v_k是传感器噪声。如果你要处理多维系统,比如同时估计车速和加速度,那么状态向量就是[v; a],状态转移矩阵就会变成2x2矩阵,A = [1 Ts; 0 1],控制输入矩阵B可以继续保留,也可以把加速度并进状态里做随机游走。总之,建模这一步的关键是把真实的物理过程离散化,并且每个矩阵的维度要逐一核对,Simulink里的矩阵乘法模块对这些是最敏感的,维度一错直接报错。
2.2 Q和R的工程整定思路
Q和R的设定是卡尔曼滤波建模里最玄学、也最影响效果的环节。很多教程喜欢说"根据经验试凑",但实际工程里是有章法的。
R代表测量噪声方差,这个最容易估计。你拿传感器静置一段时间,采集一组数据,算一下方差,基本就是R的近似值。比如轮速传感器静态采集数据,噪声方差如果是0.1 (m/s)^2,那R就定在0.1附近。Q代表过程噪声,它没有直接测量手段,通常要靠递推试验来标定。一个比较有效的初值方法是:Q的数值级从R的百分之一到十分之一开始试。如果滤波响应太迟钝,曲线像被焊死了一样老半天不动,说明Q给得偏小,过程噪声权重太低;如果滤波输出跟着测量噪声剧烈跳动,说明Q给得偏大,模型预测权重太低。
另外要留意,Q并不是越大越好,R也不是越小越好。二者之间需要一个平衡,最终看P是否收敛、滤波后的曲线是否平顺且跟踪及时。我习惯把Q、R做成Simulink模型里的Constant参数或者MATLAB脚本里的变量,每次仿真改数值不需要动模型结构,只需要在workspace里重新赋值然后运行,调试效率会提高很多。
2.3 滤波模型架构选型:纯模块、MATLAB Function还是C Function
在Simulink里实现卡尔曼滤波,通常有三条路。第一条路是用最基本的模块搭:增益、加法器、矩阵乘法、单位延迟。这种方案的好处是每个环节都可视化,特别适合教学和理解算法流程,但缺点是模型会很乱,尤其是多维系统,连线多到眼晕,而且矩阵运算模块配置起来也挺繁琐。第二条路是写成MATLAB Function或者MATLAB Function(Level-2),把五步迭代直接写在脚本里,代码直观、易于调试,是目前大多数人的首选。第三条路是C Function或S-Function,适合需要代码生成、嵌入到实际控制器里的场景,尤其是在MBD开发流程中,算法最终要生成C代码,那就得提前用C语言实现。
我个人建议:只是学习验证用,选第二条路;要做产品部署,直接从第三条路开始写。很多工程师习惯先在MATLAB Function里调通,再翻译成C,这中间反而容易因为语言差异引入新bug。倒是直接写C Function,再在Simulink里做单元测试,一步到位,后续生成代码几乎没有迁移成本。
3. Simulink系统模型的具体搭建与核心实现
3.1 一维卡尔曼滤波的模块级实现流程
为了把卡尔曼滤波五个公式和Simulink模块一一对应起来,可以先从一维系统入手。我以楼层加速度估计为例,搭建一个纯模块版本。你需要准备这些模块:Constant输入真实值并叠加Random Number模拟带噪声观测,Gain模块若干,Add模块若干,Unit Delay模块一个,以及Scope用来观察波形。
第一步,搭建预测支路。从Unit Delay中获取上一时刻的状态估计值x_est,通过Gain模块乘上A(一维时就是1),再加上B*u(这个案例里如果没有控制输入就省略),得到预测状态x_pred。同时从上一个误差协方差P_est出发,通过Gain乘上A^2再加上Q,得到预测协方差P_pred。
第二步,搭建更新支路。用预测协方差P_pred除以P_pred+R得到卡尔曼增益K。然后用测量值z减去预测状态x_pred得到新息,乘以K,再加上x_pred,得到当前时刻的状态估计值x_est。这一步的输出接回Unit Delay的输入,完成递推回路。同时计算P_est=(1-K)*P_pred,接入另一个Unit Delay,更新协方差。
这里有一个特别容易踩的坑:如果直接把状态估计信号引到前面的计算环节而不经过Unit Delay,Simulink会报代数环(Algebraic Loop)错误。解决办法就是在反馈回路上放一个单位延迟模块,让它把上一拍的值保存下来,下一拍再用。这也是卡尔曼滤波天然适合数字递推的原因所在,每个采样周期只依赖上一拍的状态。
3.2 多维状态与数组读取:Selector、Demux怎么选
当状态量从一维变成多维,问题就来了:状态向量是一个数组或者向量信号,你要在回路里取某个分量用于观测更新,或者要拆开送入不同的计算支路,这时就需要用到信号选择模块。
Simulink里处理这类问题,最常见的模块是Demux和Selector。Demux会把向量信号拆成多个标量输出,简单直接,但它对输入端口的维度是静态匹配的,状态向量维度一变,整个模型都要改,灵活性很差。Selector则要灵活得多,它可以从一个向量或矩阵信号里按索引抽取需要的行、列或者子块。
以三轴姿态估计为例,状态向量是[横滚角; 俯仰角; 横滚角速率],如果你只需要第一个分量和第三分量去做观测更新,用Selector的Index Vector模式,设置Index为[1 3],输出就是一个二维向量,正好接给后续的矩阵运算。使用Selector的时候要注意它的Port设置:如果你要从多维数组里抽数据,需要把Index Mode设为Index Vector或者Starting and Ending Indices,并正确设置输入维度。还有一个习惯性技巧:在总线信号里传卡尔曼滤波的状态估计,配合Bus Selector取总线里的各个分量,模型会清晰很多,避免了一大堆连线绕在一起。
另外,从Workspace读取测试数据时,我建议直接用From Workspace模块加上数组信号,数据格式用Structure或者Timeseries都行。如果你在信号线上看到维度标识是一个[3x1]的向量,但接收端需要[1x3],别慌,加一个Squeeze模块或者用Reshape模块调整维度即可。这类维度问题在Simulink里极其常见,尤其是卡尔曼滤波这种大量矩阵运算的模型,排查维度错误的时候,最好的办法是在关键线路上右键调出Signal Dimensions显示,一眼就能看出信号是几维的,比猜快得多。
3.3 用C Function实现卡尔曼滤波核心算法
在实际工程中,我强烈推荐在Simulink模型里嵌入C Function来实现卡尔曼滤波。这样做的直接好处是:单元测试很干净,代码生成很顺滑,而且最终的C代码可以直接给嵌入式工程师用,不需要二次移植。下面是一维卡尔曼滤波C Function的参考实现,输入是测量值z、上一个状态x_pre、上一个协方差p_pre,输出是更新后的x_est和p_est:
void kalman_1d_step(double z, double q, double r, double *x_pre, double *p_pre, double *x_est, double *p_est) { double x_pred = *x_pre; double p_pred = *p_pre + q; double k = p_pred / (p_pred + r); *x_est = x_pred + k * (z - x_pred); *p_est = (1.0 - k) * p_pred; *x_pre = *x_est; *p_pre = *p_est; }在Simulink的C Function块里,你需要定义输入端口(z、q、r)和输出端口(x_est、p_est),同时把x_pre和p_pre作为持久化的内部状态。一个比较常见的做法是使用DWork向量来保存这两个值,在每个采样周期调用函数时读入旧值、更新后写回。实现过程中有两个细节容易被忽略:一是要确认采样周期的设置,C Function块必须在离散采样模式下运行,否则系统会认为它连续执行,导致递推关系错乱;二是代码中要避免在函数内部printf或写文件,这类语句在普通仿真没问题,但一旦做代码生成或静态代码检查,很可能报错或直接被优化掉。
如果你的模型要导入外部C头文件,比如已经有现成的卡尔曼滤波库或者惯性导航算法库,可以在C Function块的Custom Code选项卡里配置头文件路径和源文件路径。配置之后,模型里直接调用库函数即可。这一步在MBD(基于模型设计)流程里非常常见,我建议所有做算法开发的人都尽早熟悉C Function的配置方式,因为它能让你基于Simulink做完整的单元测试,而不是每次都要等整套系统联调完再去查问题。
3.4 仿真对比:卡尔曼滤波与一阶低通滤波
模型搭完之后,一定要做对比验证。我建议在Simulink里建一个对比模型:信号源输出一个斜坡加正弦的复合真实信号,叠加白噪声后,分别进入一阶低通滤波模块和卡尔曼滤波子系统,然后用Scope同时观察三条曲线。这里的关键是给两个滤波器设置合理的参数,一阶低通的时间常数和卡尔曼滤波的Q/R要做到大致同一水平,否则对比不公平。
从波形上你通常会发现两个现象。第一,卡尔曼滤波输出的滞后明显小于一阶低通。这是因为它有模型预测这个前馈机制,不是一味地平滑。第二,卡尔曼滤波的稳态噪声抑制能力也更好。一阶低通在滞后和噪声抑制之间是此消彼长的关系,时间常数调小了噪声大,调大了滞后明显;卡尔曼滤波通过动态调整增益,能同时兼顾两者。
如果对比波形不理想,不要急着怀疑卡尔曼滤波本身,先检查Q和R的比例。我曾经遇到过目标跟踪场景里卡尔曼滤波输出几乎跟随噪声、没有滤波效果的情况,后来发现是把R设得太小了,测量值被滤波器当成了"高度可信",自然不肯滤波。把R调大一个数量级,波形立刻恢复正常。这是所有调参者都要经历的一课。
4. 工程落地中的常见问题与调试技巧
4.1 模型发散不收敛,先查这三件事
卡尔曼滤波模型最常见的故障就是输出直接飞了,无穷大或者NaN。遇到这种问题,我建议按顺序排查三个地方。
第一,检查状态转移矩阵A是不是写错了。尤其是多维矩阵,很容易把转置搞混。在Simulink里可以用Display模块直接看矩阵乘法的输出,如果某个元素出现NaN,基本就是矩阵维度不匹配或者相乘出错。第二,检查P矩阵的递推是否保持在正定范围。由于浮点数舍入误差,P可能在长时间运行后退化为非正定矩阵,导致K计算异常。解决办法是定期把P强制修正为对称矩阵,比如P=(P+P^T)/2,或者使用平方根卡尔曼滤波的变体来提升数值稳定性。第三,检查采样时间是否有混叠。卡尔曼滤波是离散递推算法,如果A矩阵里包含Ts,而模型步长设置和Ts不一致,那预测步的物理意义就全错了。我见过有人把Ts写成了模型步长,结果模型跑得越快,滤波效果越差,调了半天才发现是时间单位不一致。
检查代数环也是一个重点。纯模块搭建时如果反馈路径上没有延迟单元,Simulink诊断窗口会提示存在代数环,这种模型即使能仿真,结果也可能不符合预期。看到Algebraic Loop的告警,不要犹豫,直接把对应的反馈路径上加上Unit Delay,或者考虑改用函数实现。
4.2 从仿真到硬件:外部模式、C代码生成与FMU导出
卡尔曼滤波算法在Simulink里调通之后,下一步就是往工程落地走。这里有几个关键的Simulink功能,很多人不知道或者没用好。
外部模式(External Mode)是一个非常实用的功能。它允许你在一台主机上运行Simulink模型,同时通过串口或者以太网连接目标硬件,实时调整卡尔曼滤波的R、Q等参数,不用重新编译就能观测波形变化。我在调试实际传感器数据时,经常先用External Mode跑一遍,在线把噪声方差测出来,再填进模型参数里,效率非常高。不过要注意,外部模式对硬件的实时性要求比较高,如果模型本身计算量太大,或者通信链路有延迟,波形会出现毛刺,这时候可以尝试降低通信采样频率。
代码生成方面,如果最终要部署到MCU或者嵌入式Linux环境,建议用Embedded Coder生成C代码。在模型配置里选好目标硬件,设置好求解器为离散定步长,然后生成代码。生成的代码里要特别关注卡尔曼滤波函数是否被内联优化掉了,有时编译器会把看似"无用"的调试变量删掉,但如果你把这些变量标记为Output可避免这个问题。另外建议做一次静态代码检查,比如用Polyspace或者MISRA C检查,重点排查矩阵运算时数组越界的问题,卡尔曼滤波算法因为用了大量一维数组当矩阵用,索引稍微写错,编译期不会报错,运行期却会踩内存,这种问题越早发现越省事。
再一个值得提的功能是FMU导出。Simulink模型可以导出为FMU(Functional Mock-up Unit),这样就能和其他工具链做联合仿真,比如跟Carsim联合、跟Amesim联合,或者导入到别的仿真环境。卡尔曼滤波模块如果封装成FMU,最大的好处是接口固定、模型保密,可以在不泄露算法细节的前提下让别人集成。导出FMU的时候,需要在模型里设置好输入输出端口,并选择支持FMI 2.0的导出选项,整个过程并不复杂。
4.3 扩展卡尔曼滤波与惯性导航的衔接
卡尔曼滤波只适用于线性系统,但工程里大多数系统是非线性的,比如惯性导航中的姿态解算、四旋翼的飞行控制、车辆的运动学模型。这时候就要用扩展卡尔曼滤波(EKF),它的核心思想是在每一步把非线性系统在当前状态附近做一阶泰勒展开,得到雅可比矩阵,然后用这个线性化后的模型继续套用卡尔曼滤波的五步迭代。
在惯性导航场景里,状态量通常包括位置、速度、姿态角、陀螺仪零偏和加速度计零偏,维度经常达到15维以上。这时候Q矩阵和R矩阵的维度也会膨胀到15x15,手动设定已经不太现实,工程上通常用Allan方差分析陀螺和加速度计的噪声特性,再确定Q的对角元素。Simulink里搭建EKF没有新的魔法,依然是用MATLAB Function写雅可比矩阵计算,然后复用卡尔曼滤波的更新链路。
一个很常见的误区是认为EKF的雅可比矩阵只要初始给得好就行,实际上雅可比需要在每一步重新计算,因为它是状态依赖的。如果你在Simulink里偷懒,把雅可比设为常量,那模型可能在初始状态附近表现正常,一旦状态跑远,滤波输出就会迅速发散。这也是很多人在做四旋翼姿态解算时,在地面调试一切正常,一上天就出问题的原因之一。
4.4 联合仿真场景中卡尔曼滤波的位置
很多朋友会在Carsim和Simulink联合仿真里用到卡尔曼滤波。典型的场景是Carsim输出车辆动力学状态(比如纵向速度、横摆角速度),但由于传感器噪声和模型误差,直接拿到控制算法里可能不够可靠。这时候卡尔曼滤波就夹在Carsim输出的"真值+噪声"和下游控制算法之间,充当一个状态估计器。需要注意的一点是,Carsim和Simulink联合仿真的采样步长往往不一样,Carsim通常是连续动力学仿真,Simulink控制算法是离散的。卡尔曼滤波的步长必须和控制算法的离散步长一致,而不能直接采用Carsim的输出步长,否则滤波器递推的频率与状态方程离散化频率对不上,效果会大打折扣。
另外,CARSim输出的量单位通常已经换算成国际单位制,但有些版本会用km/h,这在搭建系统模型时特别容易忽略。卡尔曼滤波对量纲极其敏感,如果你把km/h当m/s代入状态方程,状态转移矩阵里的系数就会差2.78倍,整个模型看起来"滤波有效",实际输出却是错的。我的经验是每次在联合仿真里接传感器信号,先做一次单位换算并验证量级,再进入卡尔曼滤波器。
4.5 常见问题速查表
| 问题现象 | 可能原因 | 排查方向 |
|---|---|---|
| 滤波输出发散为NaN | P矩阵失去正定性或A矩阵错误 | 检查矩阵维度、强制P对称化 |
| 滤波曲线滞后严重 | Q设置太小,模型预测权重过低 | 增大Q值,或检查离散化系数 |
| 滤波曲线跟随噪声毛刺多 | R设置太小,测量值权重过高 | 增大R值,或重新计算测量噪声方差 |
| 模型报代数环错误 | 反馈回路缺少单位延迟 | 在反馈路径上插入Unit Delay |
| 仿真和代码生成结果不一致 | 采样时间或代码生成优化选项设置不一致 | 统一求解器为离散定步长,禁止优化内联关键函数 |
| 外部模式下波形卡顿 | 通信速率低或模型计算量大 | 降低通信采样频率,简化Scope采样频率 |
| 多维数组维度不匹配报错 | 信号维度设置错误 | 显示Signal Dimensions,检查Selector索引 |
最后说一个我一直沿用的习惯:每次搭完一个卡尔曼滤波模型,我都会刻意把测量信号停掉,只用预测步跑一段,看看状态估计在没有观测的情况下怎么变化。这个"开环测试"能帮你第一时间发现状态转移矩阵和离散化系数是否写对,也能直观感受Q对预测步的影响。这招帮我省了无数次在闭环里反复排查的麻烦。卡尔曼滤波的坑确实不少,但只要你把模型、参数、调试手段这三样东西理顺了,它在Simulink里绝对是你做状态估计最趁手的工具。