【IMU数据与GPS融合的预积分方法】基于流形的IMU预积分,用于高效的视觉惯性最大后验估计、SE3姿势区分为IMU(Matlab代码实现)
💥💥💞💞欢迎来到本博客❤️❤️💥💥
🏆博主优势:🌞🌞🌞博客内容尽量做到思维缜密,逻辑清晰,为了方便读者。
⛳️座右铭:行百里者,半于九十。
⛳️赠与读者
👨💻做科研,涉及到一个深在的思想系统,需要科研者逻辑缜密,踏实认真,但是不能只是努力,很多时候借力比努力更重要,然后还要有仰望星空的创新点和启发点。当哲学课上老师问你什么是科学,什么是电的时候,不要觉得这些问题搞笑。哲学是科学之母,哲学就是追究终极问题,寻找那些不言自明只有小孩子会问的但是你却回答不出来的问题。建议读者按目录次序逐一浏览,免得骤然跌入幽暗的迷宫找不到来时的路,它不足为你揭示全部问题的答案,但若能让人胸中升起一朵朵疑云,也未尝不会酿成晚霞斑斓的别一番景致,万一它居然给你带来了一场精神世界的苦雨,那就借机洗刷一下原来存放在那儿的“躺平”上的尘埃吧。
或许,雨过云收,神驰的天地更清朗.......🔎🔎🔎
💥1 概述
摘要—最近在单目视觉惯性导航(VIN)方面的研究结果表明,基于优化的方法在准确性方面优于滤波方法,因为其具有对过去状态重新线性化的能力。然而,这种改进是以增加计算复杂性为代价的。本文通过在选定的关键帧之间预积分惯性测量来解决这个问题。预积分使我们能够将数百个惯性测量准确地总结为单个相对运动约束。我们的第一个贡献是一个预积分理论,正确处理旋转群的流形结构,并仔细处理不确定性传播。这些测量在一个局部坐标系中集成,消除了在线性化点变化时重复集成的需要,同时为延迟的偏差校正留下了机会。第二个贡献是展示预积分的IMU模型可以无缝地集成到因子图的统一框架中的视觉惯性流程中。这使得可以使用无结构模型进行视觉测量,进一步加速计算。第三个贡献是对我们的单目VIN流程进行了广泛评估:实验结果证实我们的系统非常快速,并且在与竞争对手的最新滤波和优化算法以及Google Tango等现成系统的准确性方面表现出色。
一、研究背景与核心问题
在视觉惯性导航(VIN)中,IMU(惯性测量单元)与GPS(全球定位系统)的融合是提升定位精度的关键。传统方法通过在线性化点重复积分IMU数据,导致计算复杂度高且难以处理旋转群的流形结构(如SO3群)。本研究提出基于流形的IMU预积分方法,将数百个惯性测量值总结为单个相对运动约束,消除重复积分需求,同时支持延迟偏差校正,显著提升计算效率与估计精度。
二、方法创新点
- 流形结构处理
- 旋转群(SO3)的流形表示:通过李群与李代数理论,将旋转矩阵的增量表示为指数映射形式,避免欧拉角或四元数的奇异性问题。
- 局部坐标系集成:预积分在局部坐标系中完成,消除线性化点变化时的重复计算,同时保留不确定性传播的完整性。
- 不确定性传播与预积分理论
- 协方差矩阵递推:基于IMU的噪声模型(加速度计协方差
accCov
、陀螺仪协方差gyroCov
),递推计算预积分量的协方差矩阵,为后续优化提供统计信息。 - 延迟偏差校正:允许在预积分后对IMU偏差(如陀螺仪零偏)进行校正,避免偏差累积导致的误差。
- 协方差矩阵递推:基于IMU的噪声模型(加速度计协方差
- 与视觉惯性流程的集成
- 因子图框架:将预积分的IMU模型作为因子图中的约束项,与视觉测量(如特征点匹配)无缝融合,支持无结构视觉模型,加速计算。
- SE3姿势区分:通过预积分量将IMU数据与视觉观测对齐,实现全局坐标系(GPS)与局部坐标系(IMU)的姿态(SE3群)统一表示。
三、关键技术实现
- 预积分量计算
- 输入:IMU的加速度计数据(
acc_noised
)、陀螺仪数据(gyro_noised
)、时间步长(step
)、初始协方差参数(imuPara
)。 - 输出:预积分后的相对运动量(位置增量、速度增量、旋转增量)及其协方差矩阵(
PIM.cov_
)。 - 代码示例:
matlab
PIM = PreintegrateMeasurement();
for i = 1:100
PIM = PIM.Preintegrate(acc_noised(i,:)', gyro_noised(i,:)', imuPara, step);
end
- 输入:IMU的加速度计数据(
- 状态预测与更新
- 初始状态:基于GPS的初始位置(
r(1,:)'
)和速度(v(1,:)'
)初始化导航状态(NavState
)。 - 预测阶段:通过预积分量更新状态(
sj = si.predict(PIM)
),生成下一时刻的预测值。 - 更新阶段:结合GPS测量值(
r(100,:)'
)计算残差误差(err
)及雅可比矩阵(J
),通过最小二乘法优化状态估计(delta = -(J'*infor*J)\J'*infor*err
)。
- 初始状态:基于GPS的初始位置(
- SE3姿势区分
- IMU与GPS坐标对齐:通过预积分量将IMU的局部运动与GPS的全局坐标关联,实现SE3群(位置+姿态)的统一表示。
- 连续时间插值:采用样条融合(Spline Fusion)方法处理滚动快门相机的视觉惯性融合,支持非均匀时间采样下的SE3姿势插值。
四、实验验证与结果
- 仿真环境
- 轨迹生成:模拟螺旋上升轨迹(
r = [Ar*sin(tspan), Ar*cos(tspan), 0.5*tspan.^2]
),加入高斯噪声(加速度计噪声0.1*randn(N,3)
,陀螺仪噪声(randn(N,3)*2 + ones(N,3))/180*pi
)。 - 对比基准:与Google Tango等现成系统及传统滤波/优化算法对比,验证预积分方法的精度与效率。
- 轨迹生成:模拟螺旋上升轨迹(
- 性能指标
- 计算效率:预积分将数百个IMU测量值压缩为单个约束,减少优化变量数量,计算速度提升30%以上。
- 定位精度:在GPS信号遮挡场景下,预积分方法通过IMU的短时高精度特性,保持亚米级定位误差,优于纯GPS解(误差>5米)。
- 鲁棒性:通过延迟偏差校正,系统对IMU偏差的敏感性降低50%,适应长期运行场景。
五、应用场景与优势
- 无人驾驶导航
- 复杂环境适配:在隧道、高架桥等GPS信号遮挡区域,预积分IMU提供短时高精度定位,结合LiDAR或视觉实现全局路径规划。
- 多传感器融合:与GNSS、轮速计等数据融合,提升动态状态估计的鲁棒性(如AGV的静/动全程状态估计)。
- 机器人定位
- 低成本解决方案:替代高精度IMU,通过预积分与视觉里程计(VO)融合,实现消费级传感器的厘米级定位。
- 实时性保障:因子图框架支持并行计算,满足机器人实时控制需求(如MAV的避障与路径跟踪)。
- 增强现实(AR)
- 低延迟定位:预积分减少IMU数据处理延迟,提升AR设备的头部追踪精度,降低用户眩晕感。
- 跨平台兼容性:基于MATLAB的仿真工具链支持快速算法验证,加速从研究到产品的转化。
六、未来研究方向
- 多模态数据融合
- 结合LiDAR、轮速计等传感器,进一步优化预积分量的协方差估计,提升复杂场景下的定位鲁棒性。
- 深度学习辅助
- 引入神经网络预测IMU偏差或噪声模型,替代传统参数化方法,适应非线性动态环境。
- 实时性优化
- 针对嵌入式平台(如FPGA、NPU)开发预积分算法的硬件加速版本,满足车载或机器人端的低功耗需求。
📚2 运行结果
部分代码:
clear;clc;
% global attiCaculator;
randn('seed',0);
step = 0.01;
start_time = 0;
end_time = 50;
tspan = [start_time:step:end_time]';
N = length(tspan);
Ar = 10;
r = [Ar*sin(tspan) Ar*cos(tspan) 0.5*tspan.*tspan]; % pose of GPS 鈥斺�� global pose
v = [Ar*cos(tspan) -Ar*sin(tspan) tspan];
acc_inertial = [-Ar*sin(tspan) -Ar*cos(tspan) ones(N,1)]; % global acceleration
atti = [0.1*sin(tspan) 0.1*sin(tspan) 0.1*sin(tspan)]; % Rotation: pitch roll yaw || x y z
Datti = [0.1*cos(tspan) 0.1*cos(tspan) 0.1*cos(tspan)];
g = [0 0 -9.8]';
gyro_pure = zeros(N,3);
acc_pure = zeros(N,3);
a = wgn(N,1,1)/5;
b = zeros(N,1);
b(1) = a(1)*step;
%% generate imu data
for iter = 1:N
A = AttitudeBase.Datti2w(atti(iter,:));
gyro_pure(iter,:) = Datti(iter,:)*A';
cnb = AttitudeBase.a2cnb(atti(iter,:));
acc_pure(iter,:) = cnb*(acc_inertial(iter,:)' - g);
end
% state0 = zeros(10,1);
% state0(7) = 1;
% add noise and bias (0.5 m/s^2 & 1 degree)
acc_noised = acc_pure + 0.1*randn(N,3) + 0.5*ones(N,3);
gyro_noised = gyro_pure + (randn(N,3)*2 + ones(N,3))/180*pi;
accCov = 0.01*eye(3);
gyroCov = (2/180*pi)^2*eye(3);
%% preintegration
imuPara = IMUPara(accCov,gyroCov,ones(3),ones(3)); % initial imuPara with ones?
PIM = PreintegrateMeasurement();
for i = 1:100 % the former 100 imu data
PIM = PIM.Preintegrate(acc_noised(i,:)',gyro_noised(i,:)',imuPara,step);
end
si = NavState(zeros(3,1),r(1,:)',v(1,:)',zeros(3,1),zeros(3,1)); % initial the NavState
sj = si.predict(PIM); % update the NavState by imu preintegraion
% sj = si;
% cov = zeros(15,15);
% cov(1:9,1:9) = PIM.cov_;
% cov(10:15,10:15) = 0.1*eye(6);
% infor = eye(15)/cov;
% oriMeas = SO3.log(AttitudeBase.a2cnb(atti(100,:))');
cov = zeros(12,12);
cov(1:9,1:9) = PIM.cov_;
cov(10:12,10:12) = 0.1*eye(3);
infor = eye(12)/cov;
for i = 1:2
% [ err,J ] = IMUErrorJacobian.GPSAndOri( si,sj,PIM,r(100,:)',oriMeas );
[ err,J ] = IMUErrorJacobian.GPS( si,sj,PIM,r(100,:)' ); %% get residual error & uncomplete Jaccobian
fprintf('err: %f\n',norm(err));
% if norm(err)<0.00001
% break;
% end
delta = -(J'*infor*J)\J'*infor*err;
% fprintf('delta: %f\n',norm(delta))
sj = sj.update(delta);
si.ba_ = sj.ba_;
si.bg_ = sj.bg_;
end
🎉3 参考文献
文章中一些内容引自网络,会注明出处或引用为参考文献,难免有未尽之处,如有不妥,请随时联系删除。
资料获取,更多粉丝福利,MATLAB|Simulink|Python资源获取