基于EKF算法的三轴磁力仪一体化补偿方法*

孙 欢，杨宾峰，吴 昊，李 驰，王 润

(空军工程大学信息与导航学院， 西安 710077)

0 引言

地磁场作为一个天然的物理场，具有全地域、全天候、全天时等特点，以地磁场为基础的地磁导航，相比于传统无线电导航和卫星导航，具有隐蔽性好，自主性高，抗干扰能力强的特点；相比于地形匹配导航和天文导航，地磁导航受天气影响小，地磁场随时间变化缓慢，因此导航定位更稳定[1-2]。但是在实际测量中，在诸多因素的影响下，磁场测量结果将会受到多方面的干扰。首先，三轴磁力仪制造及在安装载体过程中存在工差；其次，载体内部设备结构将形成硬磁误差和软磁误差；再者，对于频繁做姿态变换的载体，比如战斗机，机体姿态转动会产生较强的涡流，进而产生涡流磁场，也会极大干扰磁力仪的测量[3-9]。

针对三轴磁力仪的标定和误差补偿问题，国内外学者进行了诸多研究，Geber-Egziabher提出two-step算法补偿三轴磁力仪刻度因子误差、零偏误差和硬磁误差[10]；Vasconcelos对三轴磁力仪误差和安装误差进行一体化建模，并提出用椭球面旋转、尺度变换、平移变换的方法进行估计，并设计了基于最大似然算法的自校正算法[11]；文献[12-13]对三轴磁力仪的仪表误差、安装误差和软硬磁误差分别进行建模，后者采用ANSYS软件仿真建模，并采用最小二乘法进行补偿；文献[14]对三轴磁力仪本身误差和软硬磁误差进行一体化建模，并考虑到观测向量和数据矩阵均存在误差，因此用总体最小二乘进行了补偿；针对载体姿态变换的问题，文献[15] 采用COMSOL Multiphysics仿真软件建立了载体涡流干扰场的仿真模型，并进行定性分析；文献[16-18]建立等效磁强计误差模型，并采用两步估计算法和圆约束非对准误差估计算法实现了补偿。

通过分析，以往的研究存在以下问题：1)模型中没有考虑涡流磁场的影响，对于低速、姿态稳定的载体尚可补偿，但对于高速、频繁进行姿态变换的载体，其涡流磁场的影响将会降低补偿精度；2)在以往研究所普遍采用的补偿方法中，最小二乘算法解算依赖大量先采数据，对数据量敏感度高，且不能做到实时补偿；椭圆约束法和椭球约束法只适用于以软磁硬磁为核心的简单建模方式，没有考虑仪表、安装、涡流所带来的误差。文中对三轴磁力仪的仪表误差和安装误差进行了分析，结合硬磁误差和软磁误差进行统一处理，考虑到飞行器载体进行频繁姿态变换的现状，引进涡流磁场进行一体化建模，并采用EKF(扩展卡尔曼滤波)算法进行补偿，最后通过模拟仿真和实验来检验建模方法和补偿方法的实用性和有效性。

1 地磁场观测值误差分析

在磁力仪标定过程中，三轴磁力仪的三轴x′、y′、z′和真实三轴x、y、z的不完全重合，将会造成非正交误差，从而影响磁力仪的测量精度。令z轴和z′轴重合，yoz和y′oz′共面， 记oy轴和oy′轴之间的夹角为β，ox轴和面x′oy′之间的夹角为γ，ox轴和面x′oz′之间的夹角为α，其对应关系如图1所示。

图1 三轴非正交误差

则非正交误差矩阵CNO可以表示为：

(1)

在磁力仪安装过程中，原则上传感器坐标系应该与载体坐标系严格一致，但是在安装过程中，将不可避免的出现安装误差，从而造成安装误差，设磁力仪三轴与载体坐标系三轴的误差角为εx、εy、εz，则安装误差CM为：

(2)

由于制作工艺精度的问题，三轴磁力仪三轴存在细微差异，由此产生刻度因子误差，设三轴的刻度因子误差分别为sx、sy、sz，则刻度因子误差矩阵CSF为：

CSF=diag(sx,sy,sz)

(3)

制作工艺精度和安装误差同样也会造成零刻度漂移误差w和剩磁误差bs，矩阵表达式分别为：

w=[δ1,δ2,δ3]T，bs=[bx1,by1,bz1]T

(4)

硬磁误差是载体内部硬磁材料产生的磁场所造成的误差，常见的硬磁材料有永磁体和高碳钢，硬磁材料的磁场随时间变化十分缓慢，在短时间内可以认为是时不变的，因此在3个轴上的投影为常值，记硬磁误差Bh为：

Bh=[bh1,bh2,bh3]T

(5)

软磁误差是软磁性材料被外部磁场激励所产生的误差，外部磁场包括地磁场和载体内部所产生的磁场，在载体内部磁场中，由于集成电路和磁场扰动所产生的干扰磁场远远小于硬磁磁场和地磁场，因此只考虑硬磁磁场Bh和地磁场H的激励效应。软磁材料被激励产生的激励磁场强度和外部激励强度呈线性比例关系，比例系数CSI取决于软磁材料的磁化特性，当外部激励磁场方向发生变化时，所产生的激励磁场方向也会随之变化，因此软磁误差BSI可以表示为：

BSI=CSI(H+Bh)

(6)

综上所述，地磁场的观测值B可以表示为：

B=CMCNOCSFCSI(H+Bh)+w+bs+ε

(7)

式中：ε为零均值高斯地磁观测噪声，其方差为∑。将其化简，则可以表示为：

B=CH+b+ε

(8)

式中：

C=CMCNOCSFCSI

(9)

b=CMCNOCSFCSIBh+w+bs

(10)

载体做姿态变换的过程中，由于自身做切割地球磁感线运动，根据法拉第电磁感应定律，将会在载体内部产生涡流，进而产生涡流磁场，对磁力仪的测量产生影响，涡流磁场的强度与感应涡流的强度成正比，而涡流的强度又与载体磁场变化率成正比，设该比例系数矩阵为P， 将其考虑在内，可得到地磁场测量值的数学模型为：

B=CH+b+P(ΔB/Δt)+ε

(11)

2 补偿算法研究

2.1 EKF算法

在利用卡尔曼滤波进行误差补偿的时候，如果出现了二次及二次以上的高次项，整个过程呈非线性，就不能用线性方法进行处理。对于非线性系统，理论上很难推导严格的递推滤波公式，因此常用近似的方法进行处理。非线性系统的线性化滤波就是利用近似的方法研究非线性问题的途径之一。EKF(扩展卡尔曼滤波)就是一种典型的方法。采用EKF进行处理时，先将非线性的传递函数f(·)和观测函数h(·)进行泰勒展开，并略去二阶及以上的高次项，从而得到非线性系统的线性化模型，再用标准卡尔曼滤波进行线性化处理。

由系统方程

(12)

Zk=h(Xk,k)+Vk

(13)

(14)

令

(15)

(16)

则状态方程为:

(17)

同理，对非线性函数h(·)进行处理，最终得观测方程为：

Zk=HkXk+yk+Vk

(18)

利用卡尔曼基本方程可得:

(19)

(20)

(21)

(22)

Pk=[I-KkHk]Pk,k-1

(23)

其中，滤波初值和滤波误差的协方差矩阵为：

(24)

(25)

2.2 EKF模型及校正原理

根据上文，地磁场观测值模型如下：

B=CH+b+P(ΔB/Δt)+ε

(26)

由于在进行姿态变换的过程中，B的不规则度较大，而地磁场真实值H的模值‖H‖不随姿态变换而变化，且以B作为自变量，H作为因变量的模型解算出的参数在实际问题中可以不经过变换直接使用，因此在研究中以‖H‖2作为研究对象，将上述模型进行变换处理，得地磁场校正模型为：

Zk=‖H‖2=(Q(B-b-P(ΔB/Δt))+ε)2=
hk(Xk)+vk

(27)

式中：Q=C-1，则变量X为：

(28)

hk(Xk)=(B-b-P(ΔB/Δt))TQTQ(B-b-P(ΔB/Δt))

(29)

可见校正模型Zk是关于变量X的二次系统，而并非常规线性系统，采用EKF进行降阶处理。

vk=(2Q(B-b-P(ΔB/Δt)))Tε+‖ε‖2

(30)

观测噪声vk统计模型服从高斯分布，其均值和方差分别为：

μk=E[vk]=-tr(Σ)

(31)

Σ(Q(B-b-P(ΔB/Δt)))+2trΣ2)

(32)

Σ={εεT}

(33)

由于在此模型中，不含过程噪声，只有观测噪声，因此EKF模型可以简化为：

(34)

(35)

(36)

3 模拟仿真分析

在仿真实验中，根据上文地磁场观测值的模型

B=CH+b+P(ΔB/Δt)+ε

其中ε为均值为0、方差为30的高斯观测噪声。设地磁场真实值H总量为50 000 nT，为了得到各个方向的磁场分量，设参数p=(-180:0.3:510)×π/180，q=(250:0.5:1400)×π/180，则地磁场分量值为：

根据上文校正模型利用EKF算法进行校正， 为保证较好的权值收敛性，令所求参数初始值X0全为0，初始协方差矩阵P0为：

P0=107×I21×21

表1 参数估计值

各个参数均能收敛，由于参数较多，选取部分参数收敛过程，如图2～图7所示。

图2 q11收敛过程

图3 q12收敛过程

图4 q13收敛过程

图5 p31收敛过程

图6 p32收敛过程

图7 p33收敛过程

经过补偿，从数据值随机抽取10组数据进行比较，结果如表2所示。

由表2中数据可以看出，对地磁场观测值进行补偿之后，误差降低两个数量级，测量精度得到明显提高。经过补偿，数据的标准差明显减小且趋于稳定，如图8所示，表明用EKF算法补偿有效且稳定。

采用不考虑涡流磁场的模型，同样用EKF算法进行补偿，将所得结果与文中所得结果进行比较，如图9所示，可以看到，不考虑涡流磁场模型可以对地磁场实现补偿，但效果欠佳，存在400 nT左右的误差，并且会上下波动，考虑涡流磁场的一体化模型将误差降至80 nT以下，从而实现更高精度的补偿。

表2 补偿前后结果比较

图8 补偿前后标准差比较

图9 不同模型补偿结果比较

4 实验结果

在该实验中，采用Bartington公司的Mag-03MSL70三轴磁通门传感器进行补偿，该磁通门传感器能够精确测量静态或交变磁场，采用National Instruments软件进行数据采集。为了模拟飞行载体在飞行过程中的不同姿态，该实验采用3FHT30C三轴手动无磁转台来提供姿态变换信息。后期用Matlab进行数据处理。采用标准磁图确定地磁场真实值。

选取一处磁场稳定的郊外地点作为试验点，将磁通门传感器置于无磁转台，快速转动并进行全方位姿态信息测量。通过对比磁图得到所处位置地磁场真实值总量为42 587 nT，将NI软件所采数据套用文中所提出的模型，并采用EKF算法进行误差补偿，补偿结果如图11所示，在补偿前，由于各种因素的影响，误差幅值达到1 800 nT左右，在用不考虑涡流磁场的模型补偿后，误差幅值明显减小，但峰值误差仍能达到500 nT，并且补偿结果不稳定，存在波动趋势，在用考虑涡流磁场的一体化模型补偿之后，补偿结果明显趋于稳定，峰值误差减小到120 nT，证明了该补偿模型和补偿方法在实际应用中的有效性。

图10 三轴磁通门传感器和无磁转台

图11 三轴磁通门传感器校正结果

5 结论

文中采用EKF算法对地磁场观测值进行补偿，在建模过程中，在考虑磁力仪测量误差和软硬磁误差的基础上，结合飞行载体频繁进行姿态变换的实际情况，加入涡流磁场误差进行一体化建模，相比于以往研究，提出的方法具有以下优点：1)一体化建模方式对于进行频繁姿态变换的载体，能够实现更高精度的补偿；2)卡尔曼滤波算法相比于最常用的最小二乘算法，无须事先采集大量数据进行解算，对数据量敏感度低，具有较强的鲁棒性，并且可以实现即时补偿；3)EKF相比于其他非线性算法，原理简单，易于实现，降低了解算难度，有效缩短了补偿时间。

在实验中，相比于不考虑涡流磁场的模型，一体化建模的模型在有效减小误差的基础上，滤除了姿态变换所引起的涡流磁场对于补偿结果的波动干扰，极大增强了补偿结果的稳定性。补偿精度的提升和稳定性的增加使得地磁导航在飞行载体的导航中有了更加广阔的前景。考虑到EKF算法不适用于强非线性问题，下一步将针对于在其线性化过程中忽略的二次及更高项所引起的误差进行研究。