四川网站建设培训班,网络教育平台,如何夸奖客户网站做的好,响应式网站实例重要说明#xff1a;本文从网上资料整理而来#xff0c;仅记录博主学习相关知识点的过程#xff0c;侵删。
一、参考资料
新手入门系列3——Allan方差分析方法的直观理解 惯性测量单元Allan方差分析详解 IMU误差测量模型 IMU标定之—Allan方差 IMU误差模型简介及VINS…重要说明本文从网上资料整理而来仅记录博主学习相关知识点的过程侵删。
一、参考资料
新手入门系列3——Allan方差分析方法的直观理解 惯性测量单元Allan方差分析详解 IMU误差测量模型 IMU标定之—Allan方差 IMU误差模型简介及VINS使用说明
参考文献
[1] An overview of the Allan variance method of IFOG noise analysis [2] Handbook of Frequency Stability Analysis [3] Allan Variance: Noise Analysis for Gyroscopes
i2Nav
武汉大学多源智能导航实验室 i2Nav b站 b站视频惯性导航原理 b站视频组合导航原理
二、相关介绍
1. 白噪声
1.1 白噪声的概念
白噪声white noise是指功率谱密度在整个频域内是常数的噪声。 所有频率具有相同能量密度的随机噪声称为白噪声。
白噪声是一种功率谱密度为常数的随机信号。换句话说此信号在各个频段上的功率谱密度是一样的由于白光是由各种频率颜色的单色光混合而成因而此信号的这种具有平坦功率谱的性质被称作是“白色的”此信号也因此被称作白噪声。
理想的白噪声具有无限带宽因而其能量是无限大但这在现实世界是不可能存在的。实际上我们常常将有限带宽的平整信号视为白噪声因为这让我们在数学分析上更加方便。白噪声在数学处理上比较方便因此它是系统分析的有力工具。一般只要一个噪声过程所具有的频谱宽度远远大于它所作用系统的带宽并且在该带宽中其频谱密度基本上可以作为常数来考虑就可以把它作为白噪声来处理。例如热噪声和散弹噪声在很宽的频率范围内具有均匀的功率谱密度通常可以认为它们是白噪声。
1.2 高斯白噪声
如果一个噪声它的幅度分布服从高斯分布而它的功率谱密度又是均匀分布的则称它为高斯白噪声。即白噪声未必是高斯白噪声。
热噪声和散粒噪声是高斯白噪声。
高斯白噪声是指信号中包含从负无穷到正无穷之间的所有频率分量且各频率分量在信号中的权值相同。白光包含各个频率成分的光白噪声这个名称是由此而来的。它在任意时刻的幅度是随机的但在整体上满足高斯分布函数。
2. MEMS技术
MEMS惯性传感器的研究背景与发展现状
MEMS惯性传感器
MEMSMicro-electromechanical Systems简称“微机电系统”MEMS集机械元件、微传感器、信号处理与控制电路、接口电路、通信和电源为一体一般认为MEMS是由微型机械传感器、执行机构和微电子电路组成的微机电系统。
2.1 MEMS惯性传感器简介
MEMS惯性传感器是利用MEMS微加工技术和半导体集成电路制造技术制作而成的惯性传感器由质量块、弹性元件、敏感元件等机械结构件传感器芯片和信号处理电路ASIC等部分组成。
MEMS惯性传感器行业分类可以根据根据被测物理量将MEMS惯性传感器分成三类微机械加速度计MEMS加速度计、微机械陀螺仪MEMS陀螺仪、微惯性测量组合MIMU。
MEMS惯性传感器工作原理MEMS惯性传感器的工作原理是经典力学中的牛顿定律。它的作用是测量运动物体(如车辆、飞机、导弹、舰船、卫星等)的质心运动和姿态运动进而对运动物体进行控制和导航。
由于MEMS惯性传感器基本的测量功能如加速度测量倾斜测量振动测量甚至旋转测量消费电子领域有待探索的应用将不断涌现。MEMS惯性传感器可以组成低成本的GNSS/INS组合导航系统。
2.2 MEMS加速度计
MEMS加速度计是利用感测质量的惯性力进行测量的传感器一般由标准质量块(感测元件)和检测电路组成。根据传感原理的不同MEMS加速度计可以分为压阻式、热流式、谐振式和电容式等。压阻式MEMS 加速度计容易受到压阻材料影响温度效应严重、灵敏度低横向灵敏度大精度不高。热流式加速度计受传热介质本身的特性限制器件频率响应慢、线性度差、容易受外界温度影响。因此热流式和压阻式加速度计主要用于对精度要求不高的民用领域或军事领域中的高 g 值测量。谐振式微加速度计理论上可以达到导航级的精度但目前技术状态还达不到实用化。而电容式硅微加速度计由于精度较高、技术成熟、且环境适应性强是目前技术最为成熟、应用最为广泛的 MEMS 加速度计。
2.3 MEMS微机械陀螺仪
MEMS微机械陀螺仪主要有线振动型陀螺仪和谐振环型陀螺仪线振动型陀螺仪工艺简单利于大批量、低成本生产谐振环型陀螺仪具有更高的理论精度但结构及原理更为复杂。
2.4 MIMU微惯性测量组合
微惯性测量组合Micro Inertial Messurement UnitMIMU是基于 MEMS技术的新型惯性测量器件由三个MEMS加速度计和三个MEMS陀螺仪组成。MIMU用来测量物体的三轴角速度和三轴加速度信息是实现微小型无人机、交通工具等导航制导的核心部件。
MIMU作为一种运动测量传感器IMU的数据采集直接决定了INS的精度和动态响应性能上限。
目前国内研制的微型惯性测量单元主要技术路线还是采用三个单轴加速度计和三个陀螺仪立体组装方式集成虽然产品指标能够满足一些现有领域的使用要求但还存在体积偏大、安装精度差、成本较高等问题。
2.5 惯性微系统
惯性微系统是利用3D异构集成技术将MEMS 陀螺、MEMS 加速度计、压力传感器、磁传感器和信号处理电路等在硅基片上进行集成并内置导航定位算法实现芯片级精确制导、导航、定位等功能同时能够与卫星导航共同组成组合导航系统是装备制导导航和定位的核心部件。
3. IMU组成部分
惯性测量单元Inertial Measuremnt UnitIMU是MEMS技术结合的微惯性测量组合因此在很多地方被称为MIMU。它主要由三个MEMS加速度传感器、三个陀螺仪和求解电路组成。
通常IMU由以下组件组成
3个加速度计3个陀螺仪数字信号处理硬件/软件功率调节器通信硬件/软件外壳。
三个加速度计相互成直角安装这样可以在三个轴上独立测量加速度X、Y和Z。三个陀螺仪也彼此成直角因此可以围绕每个加速度轴测量角速率。
典型的六轴IMU输出三个轴的加速度和三个轴的角速度。对加速度做积分可以得到速度再做一次积分可以得到位移对角速度积分可以得到角度。
4. INS惯性导航系统
惯性导航系统(Inertial Navigation SystemINS)简称惯导利用惯性测量单元 (IMU)的输出将加速度和旋转信息相结合带有位置、速度和姿态的初始信息。
5. IMU测量模型
5.1 陀螺仪测量模型 5.2 加速度计测量模型 6. IMU误差
IMU Errors and Their Effects
电子器件输出都会有误差又叫噪声。在IMU采集数据时会产生两种IMU误差确定性误差和随机性误差。为获得精确的数据需要对确定性误差和随机性误差这两种误差进行标定。Allan方差主要用于标定随机误差。
确定性误差包括交轴耦合误差Axis-misalignment比例因子误差Scale Factor零偏Bias。其中零偏为加性误差交轴耦合和比例因子误差为乘性误差。随机误差主要是高斯白噪声和bias随机游走。
IMU误差主要有以下几种
6.1 确定性误差
确定性误差又称为系统性误差。
6.1.1 交轴耦合误差Axis-misalignment
交轴耦合误差是由于IMU的三个轴非正交导致的其衡量了非正交的程度图中绿色的是理想的正交系定义为B系而红色的为传感器的非正交系定义为S系因此需要求解出S系到B系的映射关系才能得到每个轴较为准确的输出。
6.1.2 比例因子误差Scale Factor
i2Nav
比例因子误差指的是传感器输出值和实际值之间的比例误差。 6.1.3 零偏Bias
零偏指的是传感器输出与实际值之间的常值误差。
6.2 随机误差
Gyro Noise and Allan Deviation IMU Example
6.2.1 角度随机游走ARW
角度随机游走Angle Random Walk, ARW是由陀螺仪产生的噪声其反映了角速度信号中白噪声的特性积分后表现为角度随机游走。
在实际情况下当陀螺仪静止时由于热噪声等影响其角速度并不是一成不变的假设其服从白噪声模型那么其功率谱密度在整个频域内服从均匀分布。因此对角速度积分以后可以得到角度随时间的变化关系也就是角度随机游走。
6.2.2 速度随机游走VRW
速度随机游走Velocity Random Walk, VRM是由加速度计产生的噪声速度随机游走为是加速度白噪声的积分结果。
7. 陀螺零偏
新手入门系列1——如何区分惯性器件的零偏误差
如果只用一个指标来代表一款IMU的精度的话那毫无疑问是陀螺零偏。这是因为
惯导系统的精度主要取决于IMU中的陀螺器件精度而不是加速度计精度陀螺的精度指标中最重要的又是零偏误差它基本上决定了该惯导长时间独立工作时的误差发散速度。
但是这里需要特别注意的是陀螺零偏有好几种看产品指标时一定要弄清楚是哪一种陀螺零偏指标。
7.1 常值零偏
常值零偏就是这只陀螺生产出来后就一直固定不变的零偏值。对于传统的高性能惯性器件来说该误差在出厂标定时往往就被补偿干净了因此不会标注这个指标但对于低端MEMS IMU芯片来说例如用在手机中的价格不到10块钱的芯片则不可能做逐个的标定和补偿因此常会存在deg/s也就是几千deg/h量级的常值零偏。这看上去巨大无比但我们在实际使用中很容易对付例如在初始启动过程中利用几秒钟的静态数据求平均即可扣掉大部分。
7.2 全温零偏误差
全温零偏误差bias error over temperature又称为零偏的温度敏感系数或称温漂系数反映器件参数的温度敏感性。就是指陀螺零偏在其额定工作温度范围内相对于室温零偏值的变化量。对于传统惯性器件一般是逐个做了温漂补偿的因此这个全温零偏误差就是温补后的残差而对于低端MEMS芯片不可能逐个做温漂标定和补偿如果给出这个全温零偏误差可能会很吓人因此厂家往往给个零偏的温度敏感系数例如0.01 deg/s/℃。
零偏重复性
零偏重复性全称零偏逐次上电重复性run-to-run repeatability这是传统惯性器件的经典指标是指惯性器件不同次上电运行时的零偏的不重复程度。具体测量方法是在常温下将器件多次上电测量和记录每次上电的零偏数值然后统计其差异。很多新手会好奇一个传感器的零偏值每次上电的变化能有多大这么细微的误差因素都需要考虑这是因为传统惯性器件在出厂环节就把常值误差和温漂误差等主要误差都已仔细补偿掉了因此逐次上电重复性这个次要误差才成为了一个不可忽视的主要误差。
7.3 零偏稳定性
严格来说应该称为零偏不稳定性in-run instability反映器件上电稳定后其零偏随时间变化的情况。根据具体测算方法又分为两种 我国的国军标定义的零偏不稳定性采集几个小时的静态数据每10秒或100秒求平均以便抑制器件白噪声的影响然后统计这些平均值的标准差。 Allan方差给出的零偏不稳定性采集足够长时间的静态数据一般大于10小时越高等级的器件所需时间越长画Allan方差曲线取其谷底值。
前者对惯导的实际表现有比较直接的影响有现实指导意义而后者则只是反映器件在极端理想条件下的性能极限缺乏现实意义。从具体数值来看前者也比后者大几倍甚至高一个量级。
7.4 零偏的加速度敏感性
零偏的加速度敏感性英文是g-sensitivity或linear acceleration effect。陀螺的输出本来应该对加速度完全不敏感但由于其敏感结构的加工误差等因素多少还是会受到线加速度影响的我们就用零偏的加速度敏感性来描述。显然这种零偏误差只有在强动态载体上才会造成显著影响而对于常见车载、船载低动态载体往往可以忽略。需要注意的是越是高灵敏度的MEMS陀螺其微机械结构的敏感质量越大因此其加速度敏感性往往会比较大例如某高端MEMS陀螺的加速度敏感性为18 deg/h/g。
7.5 总结
基于上述这么多种零偏误差类型如果厂家只向用户提供一个零偏指标那么就应该是所有零偏误差的总和或者是所有零偏指标中最大的那个。但毫不意外的厂家最喜欢提Allan方差零偏尤其是MEMS厂商因为它的数值最小听上去最棒。而某些不负责任的商家甚至都不明确标注这是Allan方差给出的零偏以此误导用户。因此大家一定要擦亮眼睛认真理解具体指标及其测算方法尤其是对MEMS器件。
**我一般会重点关注MEMS零偏的温度敏感系数这个指标因为MEMS器件的温漂是其内在的关键误差对温度变化不敏感的话其它指标通常都不会太差。**然后再看零偏逐次上电重复性和零偏不稳定性国军标这两个指标。但MEMS器件手册往往不给这两个指标只给Allan方差零偏不稳定性此时无奈之下我会把它放大5~10倍来当作国军标零偏不稳定性来用。
还有一个最容易理解但很容易被忽视的陷阱就是给出的器件误差的具体数值含义到底是峰峰值p-p、最大值max、RMS值或是1σ假设某个随机误差是符合零均值正态分布的那么会有
峰峰值 2倍最大值 6倍RMS或σ值 所以使用不同的统计值来给出误差指标其具体数值会有天壤之别。大家在读指标时要特别注意毫不意外的厂家都喜欢用1σ值因为它最小。
因此我们调研选型惯性器件时一定要确保是把相同零偏指标的相同统计值来进行对比千万不要因为自己的马虎而被商家忽悠。
最后想跟大家说的是不要因为MEMS器件巨大的常值零偏和全温零偏误差而轻易否定它、拒绝它。反正MEMS惯导不会长时间单独工作一般都是与GNSS构成组合导航系统那么常值零偏和缓慢零偏变化包括温漂以及零偏重复性一般都能被组合导航算法例如增广Kalman滤波有效地进行在线估计和补偿此时我们真正在乎和需要关心的是零偏不稳定性中的中短期变化成份例如10s~1000s时间尺度上的波动变化而MEMS器件的这种中短期零偏不稳定性与传统惯性器件相比并不会差得太远。千万别因为我们的偏执而错过了一个物美价廉的好器件。
8. 时间序列分析方法
8.1 功率谱密度分析方法PSD
功率谱密度分析方法PSD描述误差序列在不同频率尺度上的成份属于频域分析。PSD擅长于分析中短期误差即中高频成份。
8.2 自相关分析方法
自相关分析擅长于分析中短期误差即中高频成份。
8.3 通用的样本方差分析方法
[2]
通用的样本方差分析方法与Allan方差不同之处是在分块、取平均之后不是采取相邻块求差而是截取连续M块求其标准方差然后统计这些标准方差的样本均值作为反映1~M倍块长度这个时间尺度区间的误差指标。
容易推导Allan方差只是上述通用样本方差分析在M2时的一个特例如下式
由于Allan方差取M2这一最小值因此它具有最高的时间尺度分辨率对不微附件同时间尺度上的误差成份的刻画最细腻因此被大家用得最多。
三、GNSS/INS组合导航系统
新手入门系列4——MEMS IMU原始数据采集和时间同步的那些坑
在进行GNSSGlobal Navigation Statellite System和INS数据融合时两者的时间同步的精度也直接影响着组合导航系统的最终导航定位精度。
GNSS/INS组合导航系统中有两个关键点数据采集和时间同步。
1. IMU数据采集
1.1 工业级MEMS IMU模块
一般来说消费级MEMS芯片仅能输出速率形式的测量值对应陀螺输出的是角速度常用单位为deg/s加速度计输出的是加速度严格说应该是比力常用单位为m/s2或g。
工业级MEMS IMU模块本文以ADI的ADIS16465为例。
1.2 消费级MEMS IMU芯片
对于工业级MEMS模块而言除了能够输出速率形式的测量值以外有些还能够输出增量形式的测量值对应陀螺输出的是角度增量deg加速度计输出的是速度增量m/s。
事实上离散化后的INS机械编排解算需要的正是增量形式的测量值。因此如果IMU提供增量形式的测量值我们推荐直接使用否则我们需要人为地将速率形式的测量值转为增量形式。一般地在采样率足够的情况下可以使用前后两个IMU历元的速率测量值的平均值乘上时间间隔梯形积分即可得到对应区间的增量测量值。
消费级MEMS IMU芯片本文以InvenSense的ICM20602为例。
2. IMU时间同步
IMU时间同步是指不同传感器的数据都打上共同的时标即采用统一的时间系统而不是指所有传感器都同步采样即采样时刻对齐。传感器同步采样一方面很难落实另一方面也没有必要因为多传感器的不同步采样可以通过组合导航算法中惯导的时间传递即状态转移作用来解决。
IMU时间同步对GNSS接收机和IMU之间进行时间同步。GNSS接收机本身就具备高精度授时的功能GNSS数据都天生具有精准的GNSS时间因此最理想的时间同步是将IMU数据也赋予GNSS时间。一般地通过GNSS接收机输出的1PPS信号1Hz的秒脉冲信号也可以是其它频率结合本地时钟或者IMU内部时钟即可实现高精度的硬件时间同步更确切地说是“时标统一”精度在微秒量级。然而很多应用系统不具备硬件同步的条件即无法连接1PPS脉冲信号或IMU采样信号那么我们只能被动地采用软件同步的方式精度只有几十毫秒量级。
2.1 硬件同步
一般地为了保证GNSS/INS组合导航系统的精度都采用硬件时间同步的方式实现GNSS信号采样时标和IMU数据采集时标的严格统一。对于GNSS接收机而言其1PPS信号严格对应其在GNSS整秒的GNSS信号采样时刻而对于MEMS IMU一般都具备“数据采集状态”信号接口具体形式为脉冲输出信号。ICM20602有中断输出引脚INT可配置为采样完成Data Ready中断信号输出ADIS16465则直接有一个“采样完成”信号的输出接口。
下图给出IMU硬件时间同步的示意图嵌入式MCU单片机同时接收IMU的“采样完成”脉冲信号和GNSS接收机的1PPS信号各自触发中断并记录当前的本地晶振时间。利用1PPS信号所对应的本地晶振时间例如间隔1us的计数器和GNSS整秒时间即可计算出本地晶振时间相对于GNSS绝对时间的偏移值再将这个偏移值补偿给每个IMU“采样完成”信号所对应的本地晶振时间即可将IMU的本地时标转化为GNSS绝对时标。IMU硬件时间同步的过程如公式(1)~(3)所示。最终实现IMU数据和GNSS数据的时标统一。1PPS信号的误差一般为几十纳秒可忽略不计本地晶振时间一般采用微秒级别的计数器那么最终GNSS数据和IMU数据的时间同步精度即为本地时间的精度即微秒级。 Δ t 2 ( C 3 − C 1 ) × 1.0 e − 6 b t 1.0 ( 1 ) Δt_2 (C_3 − C_1)× 1.0e-6 b_t 1.0 \quad (1) Δt2(C3−C1)×1.0e−6bt1.0(1) 可以推出 b t 1.0 − ( C 3 − C 1 ) × 1.0 e − 6 ( 2 ) b_t1.0-(C_3-C_1)\times1.0\mathrm{e}^{-6}\quad(2) bt1.0−(C3−C1)×1.0e−6(2) 那么IMU采样时刻的绝对时标 t 2 t 1 ( C 2 − C 1 ) × 1.0 e − 6 b t × ( C 2 − C 1 ) / ( C 3 − C 1 ) ( 3 ) \begin{aligned}t_2t_1(C_2-C_1)\times1.0e^{-6}\\ b_t\times(C_2-C_1)/(C_3-C_1)\quad(3)\end{aligned} t2t1(C2−C1)×1.0e−6bt×(C2−C1)/(C3−C1)(3) 其中我们假设本地时钟的额定频率为1MHz即计数间隔为1.0e-6秒 Δ t 2 Δt_2 Δt2 两个相邻1PPS之间的物理时间长度即1秒 C 1 C_1 C1某1PPS信号对应时刻的本地时钟计数值 C 2 C_2 C2某IMU采样完成信号对应时刻的本地时钟计数值 C 3 C_3 C3下一个1PPS信号对应时刻的本地时钟计数值 b t b_t bt本地时钟在1秒内即两个相邻1PPS之间的漂移值 t 1 t_1 t1对应于C1即1PPS时刻的GNSS绝对时间是由GNSS接收机提供的已知量 t 2 t_2 t2对应于C2即IMU采样完成时刻的GNSS绝对时间是待求量。
注如果本地晶振的品质较高例如优于10个ppm的TCXO那么其在1秒内的漂移值bt可忽略不计。
2.2 软件同步
许多GNSS/INS系统在开发过程中没有条件连接GNSS接收机的1PPS信号或者IMU根本不具备数据采集相关的硬件信号诸如“触发采样”、“采样完成”等信号管脚此时就不得不退而求其次采用软件同步。所谓的软件同步就是把数据处理器例如单片机、工控机等接收到IMU数据和GNSS数据的时刻“当作”这些数据的采样时刻记录下本地晶振时间作为其统一的时标。当然也可以进一步通过将GNSS数据接收时刻的本地时间与该GNSS数据中的GNSS绝对时标对比来获得两者的偏差量进而将IMU的本地时标都转化为GNSS绝对时标。
软件时间同步与硬件时间同步的道理是完全一样的只是将精准的数据采集时刻通过1PPS信号和“采样完成”脉冲来体现替换为数据接收时刻有几毫秒到几十毫秒的不确定延迟。其中GNSS数据接收时刻会受到若干因素的影响而造成延迟包括GNSS接收机内部的信号处理时间、接收机到嵌入式处理器的接口传输时间和嵌入式处理器的解码时间尤其是当处理器处于忙碌状态时就更加不可控而IMU数据接收时刻的延迟也类似只是因为每条IMU数据量很少而使得延迟较小。总之通过软件的方式实现的时间同步具有诸多不可控因素导致实际上通过软件时间同步的精度低且稳定性差。
软件同步是在确实无法实现硬件时间同步的情况下的无奈选择对于GNSS精密定位模式下的中高速载体的组合导航例如RTK/INS车辆导航几十毫秒的同步误差会深刻地伤害其组合导航精度因而是无法接受的。
2.3 多源导航传感器的时间同步
在多源导航定位系统中GNSS接收机已经成为必不可少的一部分无论是低至几块钱的导航型GNSS接收机模块还是价值上千元的测量型板卡都能够提供高精度的绝对授时误差一般为几十纳秒。因而GNSS接收机自然就成为了多源导航定位系统中的高精度时间基准。利用GNSS接收机的1PPS信号还可以对精度较低的MCU本地时间系统进行校正从而得到与GNSS时间同步的高精度MCU本地时间系统唐海亮等2019进而对多源导航传感器进行时间同步即给各传感器都打上GNSS绝对时标。
2.3.1 相机的时间同步
对于相机而言其光学信号采样对应着相机的曝光因此记录相机的曝光时间就获得了相机图像的采样时刻。一般的相机都具有硬件触发信号接口即通过外部信号触发相机曝光采样相机同时输出一个反映曝光过程的事件信号。那么我们可以利用经过GNSS时间校准后的本地时间系统实现相机的时间同步。例如每过一个固定的时间间隔MCU输出一个相机触发信号相机收到触发信号后经过内部很小的处理延时后开始曝光图像同时输出曝光信号曝光结束则曝光信号结束MCU则记录曝光开始和曝光结束的两个时间并取均值作为最终的相机图像采样时间这就实现了相机图像数据的时间同步。
2.3.2 LiDAR激光雷达的时间同步
对于激光雷达而言其时间同步方式相对简单因为激光雷达厂商已经将大部分工作集成到激光雷达内部为用户提供了现成的时间同步接口。一般地我们只要给激光雷达提供GNSS接收机的1PPS信号及其对应的GNSS绝对时间消息即可实现激光雷达点云数据的时间同步。例如Velodyne的VLP-16激光雷达对外提供1PPS信号的输入接口和一路消息输入串口分别连接GNSS接收机的1PPS信号和定位输出接口输出GPRMC或GPGGA格式的NMEA消息。而部分Livox激光雷达硬件上仅有1PPS信号输入接口需要在上位机将额外的GPRMC消息发送给激光雷达。激光雷达通过接收的1PPS信号及其对应的GNSS定位消息中的GNSS绝对时标再配合其内部的时间系统即可实现激光点云数据的时间同步。
四、Allan方差
Allan方差适用于分析惯导和时钟看重长期稳定性的传感器Allan方差曲线对中长期时间尺度上的误差特性有很强的表现力。Allan方差相比于普通的方差具有更好的描述长时间误差的优势因此Allan方差被用来当作IMU噪声标定的标准方法。
Allan方差是一个相对而言直观形象、朴实无华的概念它当初是在实践中针对高精度时钟稳定性分析的需求而提出的。
1. 问题引入
传统的误差指标往往是采用均值误差反映整个误差序列有无宏观偏置、标准差反映整个误差序列的波动情况以及均方根RMS可以认为是宏观偏置与波动情况的综合。如下图所示它们都是反映误差序列的整体情况的指标其中含有短时间快速变化和长时间缓变的各个成份无法细分出不同时间尺度上的误差波动情况。
2. Allan方差的概念
Allan方差又称为阿伦方差或阿兰方差Allan方差法是20世纪60年代由美国国家标准局的David Allan提出的基于时域的分析方法。
Allan方差是时频分析和惯性导航领域常用的一种误差分析方法它有效地刻画了待研究的误差时间序列在不同时间尺度上的波动水平不稳定性并可根据不同时间尺度上的Allan方差值所构成的曲线的形状特征来辨识其中包含的随机过程模型。Allan方差分析方法对中长期的随机波动具有很强的表现力它完全可以作为一个通用的时间序列分析工具来推广到其它应用领域。
在我们对惯性导航器件陀螺和加速度计进行误差分析时常采用Allan方差分析方法。这是当初从时频领域高精度时钟的频率稳定性借鉴来的一种时间序列分析方法 非常适合于对中长期波动不稳定性进行定量描述和分析因此正是我们惯性器件误差分析所需要的。
3. Allan方差计算过程
Allan方差是将误差序列在某个指定的时间尺度上的波动情况进行了精确提取其具体计算步骤如下 保持IMU静止 采集N个数据点。将采样数据划分为包含不同数量采样点的子集 : T ( n ) n T s n 1 , 2 … N / 2 T(n) nT_sn1,2…N/2 T(n)nTsn1,2…N/2 T ( n ) T(n) T(n) 为该子集的平均时间。 将N个点定长分组常用的交叠式Allen方差如下(n3) 对于每个子集计算均值后计算相邻两组的方差。这样就会得到 n n n个在每一平均时间 T ( n ) T(n) T(n)对应的方差 σ ( n ) σ(n) σ(n)。根据这 n n n个点作图可得到Allan标准差 σ σ σ随平均时间 T T T变化的双对数曲线。Allan方差有效反映了相邻两个采样段内平均频率差的起伏。
4. Allan方差图读取误差系数
Allan方差法可用于5种随机误差的标定
量化噪声误差系数为 Q Q QAllan方差双对数曲线上斜率为-1的线段延长线与t1的交点的纵坐标读数为 3 T Q \frac{\sqrt{3}}{T}Q T3 Q角度随机游走其误差系数 N N NAllan方差双对数曲线上斜率为 − 1 / 2 -1/2 −1/2的线段延长线与 t 1 t1 t1交点的纵坐标读数即为 N T \frac{N}{\sqrt{T}} T N零偏不稳定性其误差系数 B B BAllan方差双对数曲线上斜率为0的线段延长线与t1交点的纵坐标读数为 2 ln 2 π B \sqrt{\frac{2\ln2}{\pi}}B π2ln2 B 一般常取底部平坦区的最小值或取 t 1 0 1 t10^1 t101或 t 1 0 2 t10^2 t102 处的值角速率随机游走其误差系数 K K K斜率为1/2的线段延长线与 t 1 t 1 t1交点的纵坐标读数为 K T / 3 \frac{K}{\sqrt{T/3}} T/3 K角速率斜坡其误差系数 R R R斜率为1的线段延长线与 t 1 t1 t1交点的纵坐标读数为 R T 2 \frac{RT}{\sqrt{2}} 2 RT。
假设各种误差源统计独立那总的艾伦方差为各种误差源之和即将量化噪声的平方 σ Q σ_Q σQ、角度随机游走的平方 σ R A W σ_{RAW} σRAW、零偏不稳定性的平方 σ b i a s σ_{bias} σbias、角速率随机游走的平方 σ R R W σ_{RRW} σRRW、角速率斜坡的平方 σ R R σ_{RR} σRR的总和。
5. Allan方差分析
如果我们只对某一时间尺度上的误差即误差波动感兴趣的话那么比这个时间尺度更小的细节变化短时间快速跳动和比这个时间尺度更大的宏观变化长时间缓慢漂移就都不关心了希望在我们的误差指标中都被消除掉。而Allan方差是这样做到的 通过分块确定所要考察的时间块长度 利用块内求平均的办法把短于块长度的那些快速变化成份细节都抹掉 再利用相邻块求差的办法把长于两块长度的那些缓慢变化成份宏观都抹掉 最后对差值序列统计其方值这是处理任何随机样本的标准操作这样统计出来的就是介于1倍块长度和2倍块长度这样一个很窄的时间尺度范围内的误差波动情况。
6. Allan方差曲线
如果我们对误差序列的各个时间尺度上的成份都感兴趣的话可以将块长度由短到长“扫描”一遍得出一组Allan方差值然后画个“Allan方差 vs. 块长度”的曲线这样就能全面反映被研究的误差序列的特性了。具体的实际上是“Allan方差的开方Allan Deviation vs. 块长度”的双对数曲线以便在更大的范围上有更强的表现力。以下是一张经典的高精度陀螺的Allan方差示意图。
从Allan方差曲线上我们可以根据曲线的形状特征来识别出不同类型的随机误差即随机过程模型并提取其参数。例如斜率为-1/2的直线段代表白噪声。具体可参见IEEE制定的一个关于陀螺测试的技术标准中的附录C内容[1]。这个附录中还给出了一个从Allan方差曲线中分析陀螺的多种随机误差的案例如上图。但我想提醒大家实际的Allan方差曲线往往只能表现出少数两三个主要误差类型因为其它误差都被这几个主要误差给淹没了。而由于我们在工程上关心的恰恰也就是主要误差源因此我们并不受此困扰。更多的使用Allan方差分析过程中的注意事项可参考西工大严恭敏老师的博文《Allan方差分析的使用要点》写得非常精准务实。
7. Allan方差的应用
Allan方差分析方法并不局限于对惯性器件和时钟误差的分析也可推广到其它误差序列分析而且它也不局限于是对误差序列进行分析也可应用于任何时间序列的分析例如建筑形变、地质演化等。
8. Allan方差与ROS
下面是读取bag包并获取Allen方差的图片会生成对应的Allen方差图。这里主要参考了allan_variance_ros 对应的C文件。
sudo pip install allantools#!/usr/bin/env python
import rospy
import sys
import allantools
import rosbag
import numpy as np
import csv
import rospkg
import os
import matplotlib.pyplot as plt # only for plotting, not required for calculations
import mathdef getRandomWalkSegment(tau,sigma):m -0.5 # slope of random walk Find point where slope -0.5 randomWalk Nonei 1idx 1mindiff 999logTau -999while (logTau0):logTau math.log(tau[i],10) logSigma math.log(sigma[i],10)prevLogTau math.log(tau[i-1],10)prevLogSigma math.log(sigma[i-1],10)slope (logSigma-prevLogSigma)/(logTau-prevLogTau)# 随机漫步的斜率diff abs(slope-m)# 当前斜率与目标斜率的差值if (diffmindiff):mindiff diff# 更新最小差值idx ii i 1 Project line to tau 10^0 x1 math.log(tau[idx],10)# 当前点的横坐标y1 math.log(sigma[idx],10)# 当前点的纵坐标x2 0y2 m*(x2-x1)y1return (pow(10,x1),pow(10,y1),pow(10,x2),pow(10,y2))def getBiasInstabilityPoint(tau,sigma):i 1while (itau.size):if (tau[i]1) and ((sigma[i]-sigma[i-1])0): # only check for tau 10^0breaki i 1return (tau[i],sigma[i])def main(args):rospy.init_node(allan_variance_node)t0 rospy.get_time() Parameters bagfile rospy.get_param(~bagfile_path,~/data/static.bag)# 输入的bag文件路径topic rospy.get_param(~imu_topic_name,/imu)# 输入的imu topic名称axis rospy.get_param(~axis,0)# 输入的轴0为x轴1为y轴2为z轴sampleRate rospy.get_param(~sample_rate,100)# 输入的采样频率isDeltaType rospy.get_param(~delta_measurement,False)# 是否是delta measurementnumTau rospy.get_param(~number_of_lags,1000)# 输入的tau数目resultsPath rospy.get_param(~results_directory_path,None)# 输出的结果路径 Results Directory Path if resultsPath is None:paths rospkg.get_ros_paths()path paths[1] # path to workspaces develidx path.find(ws/)# 找到路径workspacePath path[0:(idx3)]# 取到workspace的路径resultsPath workspacePath av_results/# 结果输出路径if not os.path.isdir(resultsPath):os.mkdir(resultsPath)print(\nResults will be save in the following directory: \n\n\t %s\n%resultsPath) Form Tau Array taus [None]*numTau# 初始化tau数组cnt 0;for i in np.linspace(-2.0, 5.0, numnumTau): # lags will span from 10^-2 to 10^5, log spacedtaus[cnt] pow(10,i)# 将tau数组赋值,维度在10^-2 到 10^5cnt cnt 1 Parse Bagfile bag rosbag.Bag(bagfile)N bag.get_message_count(topic) # 在bag文件中找到该topic的消息数量data np.zeros( (6,N) ) # 初始化数据矩阵维度为6*Nif isDeltaType:scale sampleRateelse:scale 1.0cnt 0for topic, msg, t in bag.read_messages(topics[topic]):# 遍历bag文件中的消息data[0,cnt] msg.linear_acceleration.x * scaledata[1,cnt] msg.linear_acceleration.y * scaledata[2,cnt] msg.linear_acceleration.z * scaledata[3,cnt] msg.angular_velocity.x * scaledata[4,cnt] msg.angular_velocity.y * scaledata[5,cnt] msg.angular_velocity.z * scalecnt cnt 1bag.close()print ([%0.2f seconds] Bagfile parsed\n%(rospy.get_time()-t0)) Allan Variance if axis is 0:currentAxis 1 # 循环通过所有轴1-6else:currentAxis axis # 只需循环一次然后中断while (currentAxis 6):# taus_used 对应了是否可以使用taus的数组adev是allan deviation degree的缩写为allan偏差adev_err是allan偏差的误差adev_norm是allan偏差的标准化;(taus_used, adev, adev_err, adev_n) allantools.oadev(data[currentAxis-1], data_typefreq, ratefloat(sampleRate), tausnp.array(taus) )# 计算allan variancerandomWalkSegment getRandomWalkSegment(taus_used,adev)# 计算random walk segmentbiasInstabilityPoint getBiasInstabilityPoint(taus_used,adev)# 计算bias instability pointrandomWalk randomWalkSegment[3]# 获取random walk segment的纵坐标biasInstability biasInstabilityPoint[1]# 获取bias instability point的纵坐标 Save as CSV if (currentAxis1):fname allan_accel_xtitle Allan Deviation: Accelerometer Xelif (currentAxis2):fname allan_accel_ytitle Allan Deviation: Accelerometer Yelif (currentAxis3):fname allan_accel_ztitle Allan Deviation: Accelerometer Zelif (currentAxis4):fname allan_gyro_xtitle Allan Deviation: Gyroscope Xelif (currentAxis5):fname allan_gyro_ytitle Allan Deviation: Gyroscope Yelif (currentAxis6):fname allan_gyro_ztitle Allan Deviation: Gyroscope Zprint ([%0.2f seconds] Finished calculating allan variance - writing results to %s%(rospy.get_time()-t0,fname))f open(resultsPath fname .csv, wt)try:writer csv.writer(f)writer.writerow( (Random Walk, Bias Instability) )writer.writerow( (randomWalk, biasInstability) )writer.writerow( (Tau, AllanDev, AllanDevError, AllanDevN) )for i in range(taus_used.size):writer.writerow( (taus_used[i],adev[i],adev_err[i],adev_n[i]) )finally:f.close() Plot Result plt.figure(figsize(12,8))ax plt.gca()ax.set_yscale(log)ax.set_xscale(log)plt.plot(taus_used,adev)plt.plot([randomWalkSegment[0],randomWalkSegment[2]],[randomWalkSegment[1],randomWalkSegment[3]],k--)plt.plot(1,randomWalk,rx,markeredgewidth2.5,markersize14.0)plt.plot(biasInstabilityPoint[0],biasInstabilityPoint[1],ro)plt.grid(True, whichboth)plt.title(title)plt.xlabel(Tau (s))plt.ylabel(ADEV)for item in ([ax.title, ax.xaxis.label, ax.yaxis.label] ax.get_xticklabels() ax.get_yticklabels()):item.set_fontsize(20)plt.show(blockFalse)plt.savefig(resultsPath fname)currentAxis currentAxis 1 axis*6 # increment currentAxis also break if axis is not 0inpinput(Press Enter key to close figures and end program\n)if __name__ __main__:main(sys.argv)五、IMU标定工具
手写VIO-使用Allen方差工具标定IMU
多传感器融合定位3-惯性导航原理及误差分析6-Allan方差 进行随机误差分析
1. imu_utils工具
imu_utils
使用imu_utils工具生成IMU的Allan方差标定曲线
VIO 中 IMU 的标定流程 (1/3) - imu_utils 使用备忘
imu_utils是用于求取随机误差的开源工具。
2. kalibr工具
kalibr
VIO 中 IMU 的标定流程 (2/3) - kalibr_allan 使用备忘
Kalibr标定工具箱使用详细过程
Kalibr由ETHZ视觉组开源。
Kalibr可以解决以下校准问题
多摄像机校准具有非全局共享的重叠视场的摄像机系统的内部和外部校准。视觉惯性校准Camera-IMU带相机系统的IMU的时空校准。卷帘相机校准卷帘相机的完整固有校准投影失真和快门参数。
Kalibr工程有两种一种是已经编译好的包叫CDE package,另一种是未经编译的源文件。
前者安装简单不需要依赖ROS但是功能不全。
后者安装稍麻烦但功能全建议安装这种基于ROS。
3. imu_tk工具
『imu_tk』工具标定IMU的基本原理
4. gnss-ins-sim
GNSS-INS-SIM
gnss-ins-sim是一个由Aceinna新纳公司开发的GNSS/INS仿真项目可生成参考轨迹、IMU传感器输出、GPS输出、里程表输出和磁强计输出用户选择/设置传感器模型定义路径点并提供算法以及gnss-ins- sim可以生成算法所需的数据运行算法绘制 仿真结果保存仿真结果并将结果可视化。
