航磁探测系统和方法与流程

本发明涉及航磁技术，尤其涉及一种航磁探测系统和方法。

背景技术

在航磁探测过程中，需要获取飞机的姿态信息以及飞机当前姿态所引入的磁场干扰大小。一般而言，通过三分量磁力仪测量磁场矢量值可以获取飞机的姿态信息，通过总场磁力仪测量磁场的标量值获取此姿态的干扰大小。为了测量过程中数据精度尽可能高，三分量磁力仪和总场磁力仪的探头都需要放置在飞机机头前方，或者机尾后方，或者翼尖部位。

对于有人机而言，飞机空间大，改造容易，但是对于无人机航磁，无人机自身载重不如有人机，并且改造难度大。

技术实现要素：

在下文中给出了关于本发明的简要概述，以便提供关于本发明的某些方面的基本理解。应当理解，这个概述并不是关于本发明的穷举性概述。它并不是意图确定本发明的关键或重要部分，也不是意图限定本发明的范围。其目的仅仅是以简化的形式给出某些概念，以此作为稍后论述的更详细描述的前序。

鉴于此，本发明提供了一种航磁探测系统和方法，以至少解决由于无人机自身载重小难以安装三分量磁力仪等设备的问题。

根据本发明的一个方面，提供了一种航磁探测系统，航磁探测系统包括固定翼无人机主体、总场磁力仪、磁补偿器、触发单元和航姿仪，总场磁力仪、磁补偿器以及航姿仪分别设置在固定翼无人机主体上；航姿仪适于实时获取固定翼无人机主体的飞机姿态数据；触发单元适于当检测到外部触发脉冲后，触发总场磁力仪进行频率计数，并同时触发磁补偿器读取来自航姿仪的当前飞机姿态数据；总场磁力仪适于将通过频率计数获得的频率值转换为磁场值，作为当前磁总场数据发送给磁补偿器；磁补偿器适于对读取到的当前飞机姿态数据进行滑动均值滤波，以利用当前磁总场数据和经过滑动均值滤波后的当前飞机姿态数据进行航磁补偿。

进一步地，总场磁力仪安装在固定翼无人机主体的机尾后部，而磁补偿器设置于固定翼无人机主体的机身腹部。

进一步地，波特率设置为115200。

进一步地，磁补偿器适于按照如下过程进行航磁补偿：设任意时刻地磁场hearth与固定翼无人机主体的机体的三轴x、y和z之间的角度分别为α、β和γ，用cx、cy和cz分别表示α、β和γ各自的方向余弦，使得且设飞机干扰场为hd，将hd划分为恒定场hcon、感应场hinduce和涡流场hvortex；设恒定场系数为pi(i＝1，2，3)，则得：设感应场系数aij(i＝1,2,3；j＝1,2,3)，则得设涡流场系数为bij(i＝1,2,3；j＝1,2,3)，则得设磁力仪测得的总场为htotal，则

飞机在校准过程中是出在海拔3–4km处的高空中，假设地磁场hearth保持不变，则则：

将式中的(a11+1)hearth、(a22+1)hearth(a33+1)hearth、(a12+a21)hearth、(a13+a31)hearth、(a23+a32)hearth分别记为同时，bijhearth(i＝1,2,3；j＝1,2,3)记为则：

其中，δ表示由方向余弦及其导数构成的行向量，x为未知航磁干扰补偿系数构成的行向量，即：

根据htotal＝δ^tθ得到如下方程组，用于估计未知系数构成的行向量的值：h＝α^tθ；其中，θ包括18项未知系数，α是n×18的矩阵，而h是n×1的列向量；则得如下线性方程组：

即α＝(δ1δ2…δn)；其中，δi(i＝1,2,3…n)为第i次采样得的方向余弦及其导数构成的行向量，htotal(i)为第i次采样得到的总场磁力仪输出值，n为采集次数；记通带为0.06-1.0hz的滤波器为bpf，则对线性方程组进行带通滤波：bpf(h)＝bpf(α^t)θ；其中，bpf(h)表示利用滤波器对总场列向量进行带通滤波，bpf(α^t)表示利用滤波器对方向余弦矩阵的每一列进行带通滤波；对滤波后的线性方程组中未知航磁干扰补偿系数构成的行向量θ进行估计,利用递推最小二乘法求取未知航磁干扰补偿系数构成的行向量θ的估计值在实际航磁测量的时候，飞机产生的磁干扰的值hd根据下式计算得到：获得不含飞机磁干扰的有效值如下：he＝htotal-hd。

根据本发明的另一方面，还提供了一种磁探测方法，该航磁探测方法通过航磁探测系统实现，所述航磁探测系统包括固定翼无人机主体、总场磁力仪、磁补偿器、触发单元和航姿仪；所述总场磁力仪、所述磁补偿器以及所述航姿仪分别设置在所述固定翼无人机主体上；所述航姿仪适于实时获取固定翼无人机主体的飞机姿态数据；航磁探测方法包括：当触发单元检测到外部触发脉冲后，触发所述总场磁力仪进行频率计数，并同时触发所述磁补偿器读取来自所述航姿仪的当前飞机姿态数据；通过所述总场磁力仪将通过频率计数获得的频率值转换为磁场值，作为当前磁总场数据发送给所述磁补偿器；通过所述磁补偿器对读取到的当前飞机姿态数据进行滑动均值滤波，以利用所述当前磁总场数据和经过滑动均值滤波后的当前飞机姿态数据进行航磁补偿。

进一步地，所述总场磁力仪安装在所述固定翼无人机主体的机尾后部，而所述磁补偿器设置于所述固定翼无人机主体的机身腹部。

进一步地，波特率设置为115200。

进一步地，所述利用所述当前磁总场数据和经过滑动均值滤波后的当前飞机姿态数据进行航磁补偿的步骤包括：

设任意时刻地磁场hearth与固定翼无人机主体的机体的三轴x、y和z之间的角度分别为α、β和γ，用cx、cy和cz分别表示α、β和γ各自的方向余弦，使得且

设飞机干扰场为hd，将hd划分为恒定场hcon、感应场hinduce和涡流场hvortex；

设恒定场系数为pi(i＝1，2，3)，则得：

设感应场系数aij(i＝1,2,3；j＝1,2,3)，则得

设涡流场系数为bij(i＝1,2,3；j＝1,2,3)，则得

设磁力仪测得的总场为htotal，则

假设地磁场hearth保持不变，则

则：

将式中的(a11+1)hearth、(a22+1)hearth(a33+1)hearth、(a12+a21)hearth、(a13+a31)hearth、(a23+a32)hearth分别记为同时，bijhearth(i＝1,2,3；j＝1,2,3)记为则：

其中，δ表示由方向余弦及其导数构成的行向量，x为未知航磁干扰补偿系数构成的行向量，即：

根据htotal＝δ^tθ得到如下方程组，用于估计未知系数构成的行向量的值：

h＝α^tθ

其中，θ包括18项未知系数，α是n×18的矩阵，而h是n×1的列向量；则得如下线性方程组：

即

α＝(δ1δ2…δn)

其中，δi(i＝1,2,3…n)为第i次采样得的方向余弦及其导数构成的行向量，htotal(i)为第i次采样得到的总场磁力仪输出值，n为采集次数；

记通带为0.06-1.0hz的滤波器为bpf，则对线性方程组进行带通滤波：

bpf(h)＝bpf(α^t)θ

其中，bpf(h)表示利用滤波器对总场列向量进行带通滤波，bpf(α^t)表示利用滤波器对方向余弦矩阵的每一列进行带通滤波；

对滤波后的线性方程组中未知航磁干扰补偿系数构成的行向量θ进行估计,利用递推最小二乘法求取未知航磁干扰补偿系数构成的行向量θ的估计值

在实际航磁测量的时候，飞机产生的磁干扰的值hd根据下式计算得到：

获得不含飞机磁干扰的有效值如下：

he＝htotal-hd。

进一步地，该方法还包括：获取飞机在fom校准过程中的多个校准数据作为多个训练数据，其中每个训练数据包括在所述fom校准过程中由所述飞机搭载的三轴磁通门磁力仪采集的三分量数据中的沿笛卡尔坐标系的x轴和y轴的两个分量；其中，所述fom校准过程包括4个正交航向；其中，所述飞机为所述固定翼无人机主体，且该飞机搭载有所述三轴磁通门磁力仪；随机产生4个质心，作为所述4个正交航向各自初始的质心；针对每一个训练数据，对该训练数据进行分类，以获得该训练数据所属的航向类型，以及利用该训练数据重新计算其所属航向类型对应的质心；获取待处理的多个作业数据，其中，每个作业数据包括实际架次作业中由所述飞机搭载的三轴磁通门磁力仪采集的三分量数据中的沿笛卡尔坐标系的x轴和y轴的两个分量；按作业数据的采样顺序依次对各个作业数据进行分类，并将连续分到同一航向类型的多个作业数据作为该航向类型下的同一段作业数据；针对每个航向类型的每一段作业数据，利用该段作业数据中每个作业数据与该航向类型对应的质心之间的距离、按照采样顺序形成该段作业数据的距离向量，将该距离向量中出现的第一个极小值点和最后一个极小值点分别作为该段作业数据的起始点和截止点，以将该段作业数据缩短为从所述第一个极小值点至所述最后一个极小值点对应的数据段。

进一步地，对每个训练数据按照如下方式分类以获得该训练数据所属的航向类型：计算该训练数据与当前4个质心各自之间的距离，以及将该4个质心中与该训练数据之间距离最小的质心所对应的航向类型确定为该训练数据所属的航向类型。

无人机相比于有人机，由于必须具备自主飞行功能，因而无人机飞控需要知晓无人机的姿态信息。因此，在无人机航磁探测过程中，飞控模块本身就可以获取飞机的姿态信息，又在飞机上装备三分量磁力仪，这对于无人机空间是一种浪费。相比之下，本发明的航磁探测系统和方法，通过在固定翼无人机主体的机尾后部安装总场磁力仪，以及在固定翼无人机主体的机身腹部设置磁补偿器，将飞控模块的飞控姿态数据和总场磁力仪的磁总场数据发送给磁补偿器，以进行滑动均值滤波和航磁补偿，从而实现基于固定翼无人机的航磁探测，而不需安装三分量磁力仪，克服了现有技术的不足。

通过以下结合附图对本发明的最佳实施例的详细说明，本发明的这些以及其他优点将更加明显。

附图说明

本发明可以通过参考下文中结合附图所给出的描述而得到更好的理解，其中在所有附图中使用了相同或相似的附图标记来表示相同或者相似的部件。所述附图连同下面的详细说明一起包含在本说明书中并且形成本说明书的一部分，而且用来进一步举例说明本发明的优选实施例和解释本发明的原理和优点。在附图中：

图1是示意性地示出本发明的航磁探测系统的一个示例的结构框图；

图2是示意性地示出本发明的航磁探测系统的外观仰视图；

图3是示出本发明的航磁探测方法的一种示例性处理的流程图；

图4是示出机体坐标系的示意图；

图5是示出fom校准飞行的示意图；

图6是示出不同航向的(x,y)点状分布图；

图7是示出某段数据到其质心的距离的示意图。

本领域技术人员应当理解，附图中的元件仅仅是为了简单和清楚起见而示出的，而且不一定是按比例绘制的。例如，附图中某些元件的尺寸可能相对于其他元件放大了，以便有助于提高对本发明实施例的理解。

具体实施方式

在下文中将结合附图对本发明的示范性实施例进行描述。为了清楚和简明起见，在说明书中并未描述实际实施方式的所有特征。然而，应该了解，在开发任何这种实际实施例的过程中必须做出很多特定于实施方式的决定，以便实现开发人员的具体目标，例如，符合与系统及业务相关的那些限制条件，并且这些限制条件可能会随着实施方式的不同而有所改变。此外，还应该了解，虽然开发工作有可能是非常复杂和费时的，但对得益于本公开内容的本领域技术人员来说，这种开发工作仅仅是例行的任务。

在此，还需要说明的一点是，为了避免因不必要的细节而模糊了本发明，在附图中仅仅示出了与根据本发明的方案密切相关的装置结构和/或处理步骤，而省略了与本发明关系不大的其他细节。

本发明的实施例提供了一种航磁探测系统，航磁探测系统包括固定翼无人机主体、总场磁力仪、磁补偿器、触发单元和航姿仪，总场磁力仪、磁补偿器以及航姿仪分别设置在固定翼无人机主体上；航姿仪适于实时获取固定翼无人机主体的飞机姿态数据；触发单元适于当检测到外部触发脉冲后，触发总场磁力仪进行频率计数，并同时触发磁补偿器读取来自航姿仪的当前飞机姿态数据；总场磁力仪适于将通过频率计数获得的频率值转换为磁场值，作为当前磁总场数据发送给磁补偿器；磁补偿器适于对读取到的当前飞机姿态数据进行滑动均值滤波，以利用当前磁总场数据和经过滑动均值滤波后的当前飞机姿态数据进行航磁补偿。

图1给出了本发明的航磁探测系统的一个示例的结构框图，图2给出了该航磁探测系统从底部向上看的一个示例，即仰视图。

如图1所示，航磁探测系统包括固定翼无人机主体1、总场磁力仪2、磁补偿器3、触发单元4以及航姿仪5。

总场磁力仪2、磁补偿器3以及航姿仪5分别设置在固定翼无人机主体1上。

航姿仪5的姿态数据输出端连接磁补偿器3的姿态数据输入端，以将航姿仪5实时获得的固定翼无人机主体1的飞机姿态数据发送给磁补偿器3；总场磁力仪2的磁总场数据输出端连接磁补偿器3的磁总场数据输入端，以将总场磁力仪2的磁总场数据发送给磁补偿器3。

触发单元4适于当检测到外部触发脉冲后，触发总场磁力仪2进行频率计数，并同时触发磁补偿器3读取来自航姿仪5的当前飞机姿态数据。

总场磁力仪2适于将通过频率计数获得的频率值转换为磁场值，作为当前磁总场数据发送给磁补偿器3。

磁补偿器3适于对读取到的当前飞机姿态数据进行滑动均值滤波，以利用当前磁总场数据和经过滑动均值滤波后的当前飞机姿态数据进行航磁补偿。

根据本发明的实施例，总场磁力仪2例如可以安装在固定翼无人机主体1的机尾后部，而磁补偿器3例如可以设置于固定翼无人机主体1的机身腹部。

根据本发明的实施例，串口通信的波特率例如设置为115200。

根据本发明的实施例，磁补偿器3例如按照如下过程进行航磁补偿：

设任意时刻地磁场hearth与固定翼无人机主体(1)的机体的三轴x、y和z之间的角度分别为α、β和γ，用cx、cy和cz分别表示α、β和γ各自的方向余弦，使得且

设飞机干扰场为hd，将hd划分为恒定场hcon、感应场hinduce和涡流场hvortex；

设恒定场系数为ki(i＝1，2，3)，则得：

设感应场系数aij(i＝1,2,3；j＝1,2,3)，则得

设涡流场系数为bij(i＝1,2,3；j＝1,2,3)，则得

设磁力仪测得的总场为htotal，则

飞机在校准过程中是出在海拔3–4km处的高空中，假设地磁场hearth保持不变，则

则：

将式中的(a11+1)hearth、(a22+1)hearth(a33+1)hearth、(a12+a21)hearth、(a13+a31)hearth、(a23+a32)hearth分别记为同时，bijhearth(i＝1,2,3；j＝1,2,3)记为则：

其中，δ表示由方向余弦及其导数构成的行向量，x为未知航磁干扰补偿系数构成的行向量，即：

根据htotal＝δ^tθ得到如下方程组，用于估计未知系数构成的行向量的值：

h＝α^tθ

其中，θ包括18项未知系数，α是n×18的矩阵，而h是n×1的列向量；则得如下线性方程组：

即

α＝(δ1δ2…δn)

其中，δi(i＝1,2,3…n)为第i次采样得的方向余弦及其导数构成的行向量，htotal(i)为第i次采样得到的总场磁力仪输出值，n为采集次数；

记通带为0.06-1.0hz的滤波器为bpf，则对线性方程组进行带通滤波：

bpf(h)＝bpf(α^t)θ

其中，bpf(h)表示利用滤波器对总场列向量进行带通滤波，bpf(α^t)表示利用滤波器对方向余弦矩阵的每一列进行带通滤波；

对滤波后的线性方程组中未知航磁干扰补偿系数构成的行向量θ进行估计,利用递推最小二乘法求取未知航磁干扰补偿系数构成的行向量θ的估计值

在实际航磁测量的时候，飞机产生的磁干扰的值hd根据下式计算得到：

获得不含飞机磁干扰的有效值如下：

he＝htotal-hd。

本发明的实施例还提供了一种航磁探测方法，该航磁探测方法通过航磁探测系统实现，所述航磁探测系统包括固定翼无人机主体、总场磁力仪、磁补偿器、触发单元和航姿仪；所述总场磁力仪、所述磁补偿器以及所述航姿仪分别设置在所述固定翼无人机主体上；所述航姿仪适于实时获取固定翼无人机主体的飞机姿态数据；航磁探测方法包括：当触发单元检测到外部触发脉冲后，触发所述总场磁力仪进行频率计数，并同时触发所述磁补偿器读取来自所述航姿仪的当前飞机姿态数据；通过所述总场磁力仪将通过频率计数获得的频率值转换为磁场值，作为当前磁总场数据发送给所述磁补偿器；通过所述磁补偿器对读取到的当前飞机姿态数据进行滑动均值滤波，以利用所述当前磁总场数据和经过滑动均值滤波后的当前飞机姿态数据进行航磁补偿。

图3示意性地示出了本发明的航磁探测方法的一种示例性处理的流程图。该方法中的航磁探测系统例如可以采用与上文结合图1和图2所描述的航磁探测系统相同的结构，并能够达到相类似的功能和效果，这里不再赘述。

如图3所示，在步骤s310中，当触发单元检测到外部触发脉冲后，触发总场磁力仪进行频率计数，并同时触发磁补偿器读取来自航姿仪的当前飞机姿态数据。

接着，在步骤s320中，通过总场磁力仪将通过频率计数获得的频率值转换为磁场值，作为当前磁总场数据发送给磁补偿器。

然后，在步骤s330中，通过磁补偿器对读取到的当前飞机姿态数据进行滑动均值滤波，以利用当前磁总场数据和经过滑动均值滤波后的当前飞机姿态数据进行航磁补偿。

例如，总场磁力仪可以安装在固定翼无人机主体的机尾后部，而磁补偿器可以设置于固定翼无人机主体的机身腹部。

例如，波特率可以设置为115200。

根据本公开的实施例，利用当前磁总场数据和经过滑动均值滤波后的当前飞机姿态数据进行航磁补偿的步骤例如可以与上文所描述的磁补偿器3进行航磁补偿的过程相同，这里不再赘述。

根据一个实现方式，上述航磁探测方法还可以包括：获取飞机在fom(figureofmerit，品质因数)校准过程中的多个校准数据作为多个训练数据，其中每个训练数据包括在所述fom校准过程中由所述飞机搭载的三轴磁通门磁力仪采集的三分量数据中的沿笛卡尔坐标系的x轴和y轴的两个分量；其中，所述fom校准过程包括4个正交航向，分别作为第j个航向类型，j＝1,2,3,4；其中，所述飞机为所述固定翼无人机主体，且该飞机搭载有所述三轴磁通门磁力仪；随机产生4个质心，作为所述4个正交航向各自初始的质心；针对每一个训练数据，对该训练数据进行分类，以获得该训练数据所属的航向类型，以及利用该训练数据重新计算其所属航向类型对应的质心；获取待处理的多个作业数据，其中，每个作业数据包括实际架次作业中由所述飞机搭载的三轴磁通门磁力仪采集的三分量数据中的沿笛卡尔坐标系的x轴和y轴的两个分量；按作业数据的采样顺序依次对各个作业数据进行分类，并将连续分到同一航向类型的多个作业数据作为该航向类型下的同一段作业数据；针对每个航向类型的每一段作业数据，利用该段作业数据中每个作业数据与该航向类型对应的质心之间的距离、按照采样顺序形成该段作业数据的距离向量，将该距离向量中出现的第一个极小值点和最后一个极小值点分别作为该段作业数据的起始点和截止点，以将该段作业数据缩短为从所述第一个极小值点至所述最后一个极小值点对应的数据段。

fom校准飞行通常可以由诸如固定翼飞机、直升机或无人机等磁力仪搭载平台来，在本实施例中，例如也可以采用安装有总场磁力仪和三分量磁力仪的飞机(如无人机)作为上述磁力仪搭载平台的一个示例。

笛卡尔坐标系的x轴、y轴和z轴分别沿着飞机的横向轴、纵向轴和垂直轴。如图4所示，点o为坐标原点，安装有总场磁力仪和三分量磁力仪。笛卡尔坐标系的三轴x,y,z分别沿着飞机横向轴、纵向轴和垂直轴，n为地理北极的方向，he为地磁场的方向，飞机三轴与地磁场的夹角分别为α,β,γ。

如图5所示，fom飞行包括该飞机在4个正交航向的飞行、并在每个航向上进行俯仰、横滚、偏航三个类型的机动。一组标准的fom飞行需要完成4个正交航向的飞行，在每个航向上进行俯仰、横滚、偏航三个类型的机动，幅度分别为±5°、±5°、±10°，每种机动进行30秒，完成3组。

应当理解的是，fom飞行并不限于图5所示的航迹，也可以是其他常规航迹的fom飞行，但无论是哪一种fom飞行，都大致包括4个航向的飞行，而且这四个航向是依次相正交的。比如，在图5中，可以将a至b这段作为第1类航向，将b至c这一段作为第2类航向，将c至d这一段作为第3类航向，将d至e这一段作为第4类航向，第1类航向与第2类航向正交，第2类航向与第3类航向正交，第3类航向与第4类航向正交，第4类航向与第1类航向正交。

例如，假设当前四个质心分别为(a1，b1)、(a2，b2)、(a3，b3)和(a4，b4)，若某个训练数据进行分类后得到该训练数据所属的航向类型是1，也即，对应的当前质心是(a1，b1)，于是重新计算第1个航向类型(j＝1)的质心，可以采用该训练时候与其所属航向类型的当前质心的平均值来作为新的质心，也即，将作为第1个航向类型(j＝1)的最新的质心。

根据一个实现方式，对每个作业数据例如可以通过如下方式分类：计算该作业数据与当前4个质心各自之间的距离，将与该作业数据之间距离最小的质心所对应的航向类型确定为该作业数据所属的航向类型，完成对该作业数据的分类。

根据一个实现方式，例如还可以包括如下处理：针对每个航向类型的每一段作业数据，若该段作业数据的长度(即包含的数据数量)小于预定阈值，则舍弃该段作业数据。预定阈值例如可以根据经验值设定，比如设为10或20等。

根据一个实现方式，对于每个航向类型的每一段作业数据的距离向量，例如可以按照如下方式获得该距离向量中的极小值点：针对该距离向量中第2个元素至倒数第2个元素之间的每个元素，计算该元素的前后差分，若该元素的前差分小于0、且该元素的后差分大于0，则判定该元素为该距离向量中的一个极小值点。

根据另一个实现方式，对于每个航向类型的每一段作业数据的距离向量，例如也可以按照如下方式获得该距离向量中的极小值点：针对该距离向量中第2个元素至倒数第2个元素之间的每个元素，计算该元素的前后差分，若该元素的前差分小于0、且该元素的后差分大于0，并且该元素的值小于该距离向量中各元素的平均值，则判定该元素为该距离向量中的一个极小值点。

此外，根据本发明的实施例，针对每个航向类型的每一段作业数据，例如可以按照如下方式获得该段作业数据的距离向量：假设该段作业数据中各作业数据按照采样顺序分别到该航向类型的当前质心之间的距离依次为dk+1、dk+2、……、dk+l，k＝1,2,…,ks，ks表示该航向类型对应包括的一段或多段作业数据的总段数，k表示该段作业数据是该ks段作业数据中的第k段，dk+p表示该第k段作业数据中的第p个作业数据与该航向类型当前质心之间的距离，p＝1,2,…,l，l为该段作业数据包含的作业数据数；将向量(dk+1，dk+2，……，dk+l)作为该段作业数据的距离向量。

在一个例子中，对于每个航向类型的每一段作业数据的距离向量，可以按照如下方式获得该距离向量中的极小值点：针对k+2≤k+p≤k+l-1对应的每个元素dk+p，计算dk+p的前差分diffb^k+p＝dk+p-dk+p-1，计算dk+p的后差分difff^k+p＝dk+p+1-dk+p，若diffb^k+p＜0且difff^k+p＞0，判定该元素dk+p为该距离向量(dk+1，dk+2，……，dk+l)中的一个极小值点。

在另一个例子中，对于每个航向类型的每一段作业数据的距离向量，也可以按照如下方式获得该距离向量中的极小值点：针对k+2≤k+p≤k+l-1对应的每个元素dk+p，计算dk+p的前差分diffb^k+p＝dk+p-dk+p-1，计算dk+p的后差分difff^k+p＝dk+p+1-dk+p，若diffb^k+p＜0且difff^k+p＞0，并且dk+p≤mean(dk+1,dk+2,...,dk+l)，则判定该元素dk+p为该距离向量(dk+1，dk+2，……，dk+l)中的一个极小值点。

例如，假设在某段作业数据(如第k段作业数据的l个作业数据)的距离向量(dk+1，dk+2，dk+3，dk+4，……，dk+l-2，dk+l-1，dk+l)中，判定dk+2、dk+4、…、dk+l-1为按照采样顺序依次出现的极小值点，并假设dk+1至dk+l的平均值为dmean，并且上述极小值点dk+2、dk+4、…、dk+l-1均小于dmean，则将该段作业数据缩短至从第一个极小值点dk+2至最后一个极小值点dk+l-1对应的数据段，也即，缩短为：换句话说，在第k段作业数据中，去除了无用数据和

根据一个实现方式，对每个训练数据可以按照如下方式分类以获得该训练数据所属的航向类型：计算该训练数据与当前4个质心各自之间的距离，以及将该4个质心中与该训练数据之间距离最小的质心所对应的航向类型确定为该训练数据所属的航向类型。

在一个了例子中，可以首先获得fom校准过程数据中的三分量数据，记录其中x和y分量的数据为和将其作为训练数据，并训练每个航向的质心。

其中，训练每个航向的质心的过程例如可以按照如下方式进行：

1.随机产生4个二元组(aj,bj)，j＝1,2,3,4作为初始质心；

2.对每一个训练数据中的二元组依据(1)式计算其与每一个质心的距离di^j；

3.最小的距离所对应的j即为其分类；

4.重新计算每一类中的元素均值，作为新的质心；

5.重复2-4，直至读取完所有训练数据。

然后，读取某架次所有作业数据中的三分量数据，记录其中x和y分量的数据为和

对于作业数据中每一个采样点(即每个作业数据)的二元组计算其与当前4个质心中每一个质心的距离di’^j，j＝1,2,3,4，取其最小的距离所对应的j作为其分类，获取对作业数据中每一个采样点的分类。

将连续分到某一类的采样点作为一段，若该段长度小于阈值th(可以根据经验值设定)，则舍弃该段。

然后，去除每一个所取得的数据段的无用数据，可以按照如下方式进行：

1.获得某段二元组到其对应分类的质心的距离dk+1,dk+2,...,dk+l；

2.对每一个k+2≤i'≤k+l-1，依据(2)、(3)式计算其前后差分diffb^i'、difff^i'；

3.若diffb^i'＜0且difff^i'＞0，同时di'≤mean(dk+1,dk+2,...,dk+l)，记录该i'为极小值点；

4.重复2，3获得整段数据的极小值点序列，记第一个极小值点为ifirst，最后一个极小值点为ilast；

5.将该航向数据段缩短为采样点号由ifirst至ilast的数据；

6.对所有数据段重复1-5步。

这样，取得该架次所有非转向的作业数据段，分段完毕，标记绘图。

在另一个例子中，

针对每一个训练数据根据计算该训练数据与当前4个质心各自之间的距离，将该4个质心中与该训练数据之间距离最小的质心所对应的航向类别确定为该训练数据所属的航向类别，并利用该训练数据重新计算其所属航向类别对应的质心；其中，di^j表示训练数据与第j类航向对应的质心之间的距离。

接着，获取待处理的作业数据，其中，所述作业数据包括多个二元组n'为所述作业数据总数。

然后，确定每个作业数据的航向类型，以获得每个航向类型的一段或多段作业数据。

接着，针对每个航向类型的每段作业数据，分别执行如下处理：利用该段作业数据中每个作业数据与该航向类别对应的质心之间的距离、按照采样顺序形成该类航向的距离向量；然后，对该段作业数据，将该距离向量中出现的第一个极小值点和最后一个极小值点分别作为该段作业数据的起始点和截止点，以将该段作业数据缩短为从所述第一个极小值点和所述最后一个极小值点对应的数据段。

其中，数据质心训练依赖航磁补偿系数校准阶段数据完成。在进行航磁探测过程中，需要对磁力仪搭载平台所产生的磁干扰进行补偿。补偿通过一组补偿系数与由所搭载的三轴磁通门产生的数据内积来实现。而求解该组补偿系数的过程被称作校准过程。在一次航磁探测开始的阶段，均需要对补偿系数进行校准。

设一组校准数据中由三轴磁通门所采集的三个分量分别为x、y、z，则每一次采样均产生一组三元组(xi,yi,zi)，i＝1,2,...,n，n为采样点数。如图4所示，磁通门所采集的三个分量可以认为是地磁场在飞机坐标系三个坐标轴上的投影，因此x和y将跟随飞机航向进行相应的变化，而z则与航向变化无关，因此，在计算时仅仅考虑一对二元组(xi,yi)即可，也即，校准数据中的二元组(xi,yi)作为训练数据而n即为训练数据总数。

如图5所示，标准的校准飞行被称作fom飞行，包含4个正交航向上的3种机动，因此该组数据可以被分成4段。每一段的二元组分布位置将具有明显的区别，因为飞机变换到正交航向时坐标系也随之旋转90°，x与y的取值将发生巨大变化。

图6作为例证展示了一组真实数据的(xi,yi)(相当于)分布点状图。不同航向上产生的(x,y)二元组被较为明显的区分在不同位置，该性质可以用来对航向进行区分。

可以看出，每个航向的点将集中分布与某一区域，而与其他区域相隔较远。因此若计算出每个区域的质心，则可借助于一个半径将所有该航向数据涵盖，从而对航向进行区分。

质心计算例如可以通过如下处理实现：

1.随机产生4个质心，分别为(aj,bj)，j＝1,2,3,4；(aj,bj)表示4个正交航向中的第j个航向类型对应的质心；

2.重复下述步骤：

a.对于每一个训练数据计算其与每一个质心的距离：

b.最小的距离所对应的j即为其分类；

c.重新计算每一类中的元素均值，作为新的质心。

最终更新后的4个质心(aj,bj)即为训练所得。

经过上述训练，获得了四个正交航向所对应的两个分量构成的二元对的质心。接下来在测量数据(相当于上文所述的作业数据)中，若飞机沿其中一个航向飞行，则所得(x,y)(即上文所述的)将会距离其对应质心很近，而距离其他质心很远。但飞机在转向过程中的一部分数据也会被分到某一个航向的平飞段中去，这部分数据需要去除。

本发明所采用的方法是，利用距离向量，取得距离向量的第一个最小值和最后一个最小值，作为该航向的起始点和截止点。

设为一段属于某一类的二元组，共有l个采样点(即l个作业数据点)，k表示第k段数据，例如k＝1,2,3,…，假设这段数据对应的航向类别是第j0类航向(j0可以是1、2、3或4)，并假设第j0类航向的质心是同时，设该段的距离向量为dk+1,dk+2,...,dk+l，其中，dk+1表示与其所述航向质心之间的距离，dk+2表示与其所述航向质心之间的距离，……，dk+l表示与其所述航向质心之间的距离。由于飞行时，转入该航向和转出该航向的过程中所产生的二元组距离该航向对应质心较近，会被分入该航向。因此，距离向量d将会由远及近，如附图7所示的某段数据到其质心的距离符合先减小后增大的规律，通过选择极小值点可以去除掉转向中的数据段，但需要注意圆圈中的一个极小值所造成的影响。

d的第一个极小值点和最后一个极小值点的位置应属于该航向数据段，而非转弯过程。但是，有时由于气流影响，飞机在飞行过程中所产生的抖动可能会使得转向段的距离向量并非平滑向下，而是产生一个如图7圆圈中的极小值点。为了避免该情况对距离分段产生影响，在计算极小值点时应当选择极小值小于该段距离均值的点位。

具体做法描述如下：

1.对每一个k+2≤i’≤k+l-1计算前后差分

diffb^i’＝di’-di’-1(2)

difff^i’＝di’+1-di’(3)

2.若diffb^i’＜0且difff^i’＞0，同时di’≤mean(dk+1,dk+2,...,dk+l)，记录该i'为极小值点；

3.重复1，2获得整段数据的极小值点序列，记第一个极小值点为ifirst，最后一个极小值点为ilast；

4.将该航向数据段缩短为采样点号由ifirst至ilast的数据。

通过以上描述可知，本发明所提出的方法能够自适应的根据飞行航向对飞行数据进行分段，能够大大简化航磁测量数据后处理过程，节省人工成本和处理时间，大幅提高后处理效率。

尽管根据有限数量的实施例描述了本发明，但是受益于上面的描述，本技术领域内的技术人员明白，在由此描述的本发明的范围内，可以设想其它实施例。此外，应当注意，本说明书中使用的语言主要是为了可读性和教导的目的而选择的，而不是为了解释或者限定本发明的主题而选择的。因此，在不偏离所附权利要求书的范围和精神的情况下，对于本技术领域的普通技术人员来说许多修改和变更都是显而易见的。对于本发明的范围，对本发明所做的公开是说明性的，而非限制性的，本发明的范围由所附权利要求书限定。