基于正则化方法的未爆弹三维立体成像方法、装置及系统

技术领域

本发明涉及探测技术领域，特别是涉及一种基于正则化方法的未爆弹三 维立体成像的方法、装置以及系统。

背景技术

目前，在我国的城市化建设中，未爆炸弹的存在会给人民群众的身体健 康和人身安全造成了严重威胁，同时也会带来巨大的生命财产损失。

为了解决上述问题，需要使用技术手段对遗弃化学弹进行探测，而遗弃 化学弹的探测，可以归为地球物理勘探的范畴。地球物理勘探简称“物探”， 而物探主要包括以下六种方法，分别是电磁勘探、磁法探测、重力勘探、地 震勘探、地热测量法、放射性测量法。

在遗弃化学弹探测工作中，磁法探测的应用成本低，速度快，而且便于 实施。因此，目前日本遗弃化学弹的探测中使用的主要方法就是磁法探测。 传统的磁法探测主要用于研究地质构造或探测矿产资源，其目标的空间尺度 很大，而且往往不需要得到目标精确的几何形状。但是在遗弃化学弹探测中， 探测目标小，而且对目标空间位置和几何形态的要求更高。因此传统的磁法 探测方法往往不能直接使用。亟需一种适应小目标，且相对更精确的成像方 法，从而提高遗弃化学弹发掘工作的工作效率。

发明内容

鉴于上述问题，本发明提供一种基于正则化方法的未爆弹三维立体成像 的方法、装置以及系统，提出了一种适应小目标，且相对更精确的成像方法 的技术方案。

第一方面，提供一种基于正则化方法的未爆弹三维立体成像的方法，所 述方法包括：

获取目标探测区域的磁场三分量数据；

根据所述磁场三分量数据，得到梯度张量数据和磁总场梯度数据；

根据所述梯度张量数据，得到第一先验信息；

将所述第一先验信息加入到正则化方法中的目标函数中；

根据第二先验信息和所述磁总场梯度数据，对所述目标函数的最优解进 行运算，得到所述正则化方法的反演结果，所述第二先验信息为已知先验信 息；

根据所述反演结果，得到所述未爆弹的三维立体成像。

可选地，根据所述梯度张量数据，得到第一先验信息，包括：

根据所述梯度张量数据，采用欧拉反褶积法进行运算，得到所述第一先 验信息。

可选地，根据所述梯度张量数据，采用欧拉反褶积法进行运算，得到所 述第一先验信息，包括：

根据所述梯度张量数据，采用欧拉反褶积法进行运算所基于的欧拉反褶 积公式为：

其中x

基于所述欧拉反褶积公式，进行运算过程中，选取多个相邻数据点组成 一个滑动窗口，使用所述滑动窗口内的数据求解所述欧拉反褶积公式组成的 线性方程组，得到异常体中心坐标；

随着所述滑动窗口的滑动，计算出多组异常体中心坐标，其中包括的异 常解使用主体异常距离准则和反演结果的聚散度准则进行筛选；

去除所述异常解后，剩余正常解的分布作为所述第一先验信息。

可选地，所述多组异常体中心坐标包括：Q

Q

其中，T表示转置。

可选地，所述目标函数的表达式为：

Φ＝||W

其中，d是一个N维矢量，表示在地表观测到的某类磁异常，包括：磁 总场、磁场分量或者是其它种类的数据，m是一个M维矢量，表示地下每个 长方体块的物性参数，包括：磁化强度、磁化率或者其它物性参数，A是一 个N×M维的矩阵，表示观测数据d与物性参数m之间的正演算子，矩阵中每 个元素的值根据正演公式进行求解，W

所述模型加权函数矩阵的表达式为：

W

所述模型加权函数矩阵的表达式中的W

其中，W

可选地，将所述第一先验信息加入到正则化方法中的目标函数中，包括：

将所述第一先验信息加入到所述模型加权函数矩阵中；

令第一表达式为：

令第二表达式为：

其中，σ

将所述第一表达式作为异常体中心坐标代入W

可选地，根据第二先验信息和所述磁总场梯度数据，对所述目标函数的 最优解进行运算，包括：

根据第二先验信息和所述磁总场梯度数据，利用共轭梯度法对所述目标 函数的最优解进行运算。

可选地，所述反演结果为一个三维矩阵，所述三维矩阵中每个元素的值 代表一个长方体单元的物性参数，根据所述反演结果，进行所述未爆弹的三 维立体成像，包括：

对所述三维矩阵分别沿x、y、z轴方向切片；

对切片后的三维矩阵进行二值化处理并做三维显示，得到所述未爆弹的 三维立体成像。

第二方面，提供一种基于正则化方法的未爆弹三维立体成像的装置，所 述装置包括：

获取三分量模块，用于获取目标探测区域的磁场三分量数据；

张量和总场梯度模块，用于根据所述磁场三分量数据，得到梯度张量数 据和磁总场梯度数据；

先验信息模块，用于根据所述梯度张量数据，得到第一先验信息；

加入模块，用于将所述第一先验信息加入到正则化方法中的目标函数 中；

反演模块，用于根据第二先验信息和所述磁总场梯度数据，对所述目标 函数的最优解进行运算，得到所述正则化方法的反演结果，所述第二先验信 息为已知先验信息；

成像模块，用于根据所述反演结果，得到所述未爆弹的三维立体成像。

可选地，所述先验信息模块具体用于：

根据所述梯度张量数据，采用欧拉反褶积法进行运算，得到所述第一先 验信息；

可选地，所述反演模块具体用于：

根据第二先验信息和所述磁总场梯度数据，利用共轭梯度法对所述目标 函数的最优解进行运算。

可选地，所述反演结果为一个三维矩阵，所述三维矩阵中每个元素的值 代表一个长方体单元的物性参数，所述成像模块包括：

切片单元，用于对所述三维矩阵分别沿x、y、z轴方向切片；

二值化单元，用于对切片后的三维矩阵进行二值化处理并做三维显示， 得到所述未爆弹的三维立体成像。

第三方面，提供一种基于正则化方法的未爆弹三维立体成像的系统，所 述系统包括：探测仪，所述探测仪用于执行上述第一方面任一所述的方法。

本发明提供的基于正则化方法的未爆弹三维立体成像的方法，获取目标 探测区域的磁场三分量数据；根据磁场三分量数据，得到梯度张量数据和磁 总场梯度数据；根据梯度张量数据，得到第一先验信息；将第一先验信息加 入到正则化方法中的目标函数中；根据第二先验信息和磁总场梯度数据，对 目标函数的最优解进行运算，得到正则化方法的反演结果，根据反演结果， 得到未爆弹的三维立体成像。

本发明的方法，利用梯度张量数据消除背景场的影响，使得第一先验信 息对应的边界识别效果更好。而将第一先验信息加入到正则化方法中的目标 函数中，再对该目标函数的最优值进行运算，得到的正则化方法的反演结果， 相较于目前的探测方法的精确度更高，因此最终得到的未爆弹的三维立体成 像结果也更精确，对未爆弹的识别效果自然得到了提升，从而提高了遗弃化 学弹发掘工作的工作效率。

附图说明

通过阅读下文优选实施方式的详细描述，各种其他的优点和益处对于本 领域普通技术人员将变得清楚明了。附图仅用于示出优选实施方式的目的， 而并不认为是对本发明的限制。而且在整个附图中，用相同的参考符号表示 相同的部件。在附图中：

图1是本发明实施例一种基于正则化方法的未爆弹三维立体成像的方法 的流程图；

图2是本发明实施例中一种磁异常中心水平分布俯视示意图；

图3是本发明实施例中在图2基础上，去掉不好的解后的一种磁异常中 心水平分布俯视示意图；

图4是本发明实施例中在图3基础上，显示竖直位置的统计学特性的统 计示意图；

图5是本发明实施例中遗弃化学弹的三维成像示意图；

图6是本发明实施例一种基于正则化方法的未爆弹三维立体成像的装置 的框图。

具体实施方式

为使本发明的上述目的、特征和优点能够更加明显易懂，下面结合附图 和具体实施方式对本发明作进一步详细的说明。应当理解，此处所描述的具 体实施例仅用以解释本发明，仅仅是本发明一部分实施例，而不是全部的实 施例，并不用于限定本发明。

参照图1，示出了本发明实施例一种基于正则化方法的未爆弹三维立体 成像的方法的流程图，该方法包括：

步骤101：获取目标探测区域的磁场三分量数据。

本发明实施例中，首先需要获取目标探测区域的磁场三分量数据。一般 情况下，磁探测仪对目标区域进行探测时，得到的数据即为磁场三分量数据。 当然也可以通过其它方式获取到目标探测区域的磁场三分量数据，本发明实 施例对此不作具体限定。

步骤102：根据磁场三分量数据，得到梯度张量数据和磁总场梯度数据。

本发明实施例中，获得磁场三分量数据后，进行运算，即可得到梯度张 量数据和磁总场梯度数据。一般情况下，由于实际测量得到的数据中会同时 包含背景场与磁异常，所以需要滤除数据中背景场的影响。可以对磁场三分 量数据沿不同方向求梯度，得到的梯度张量数据可以消除背景场的影响。

而对磁场的三分量数据求平方和运算，可以得到磁总场模量数据，再对 两层不同高度的磁总场模量数据求差值，即可得到磁总场梯度数据。

步骤103：根据梯度张量数据，得到第一先验信息。

本发明实施例中，运算得到梯度张量数据之后，再对梯度张量数据采用 欧拉反褶积法进行运算，可以得到第一先验信息。欧拉反褶积法属于边界识 别方法，由于梯度张量数据可以更有效地屏蔽背景场的干扰，因此可以使得 边界识别效果更好，更有利于提升正则化方法中对目标函数最优解求解的精 确度。

步骤104：将第一先验信息加入到正则化方法中的目标函数中。

本发明实施例中，得到第一先验信息，发明人创造性的将第一先验信息 加入到正则化方法中的目标函数中，这样做，是为了提高后续正则化方法中 目标函数求取最优解的精确度，具体如何将第一先验信息加入到正则化方法 中的目标函数中在下文说明，先不赘述。

步骤105：根据第二先验信息和磁总场梯度数据，对目标函数的最优解 进行运算，得到正则化方法的反演结果，第二先验信息为已知先验信息。

本发明实施例中，在将第一先验信息加入到正则化方法中的目标函数中 后，再根据第二先验信息和前面已经得到的磁总场梯度数据，对目标函数的 最优解进行运算，得到正则化方法的反演结果，所谓正则化方法的反演结果 就是对目标函数的最优解进行求解而得到的结果。

所谓第二先验信息为已知先验信息。所谓已知先验信息包括：磁化强度 等信息，例如：根据大量测试得到的已知未爆弹的最大磁化强度：M

而在正则化方法的反演中，需要提前已知异常体的磁化方向(一般以地 磁场的方向近似为异常体的磁化方向)，而梯度张量数据对磁化方向比较敏 感，因此在正则化方法的反演中，可以使用对磁化方向不敏感的磁总场梯度 数据进行反演。

本发明实施例中，根据第二先验信息和磁总场梯度数据，利用共轭梯度 法对目标函数的最优解进行运算，而第二先验信息是在求解最小值的迭代过 程中使用，求解时需要一步一步迭代，逐渐优化得到最优解，第二先验信息 就是在这个迭代过程中添加进去的。

一种优选的迭代过程为：假设已知未爆弹的最大磁化强度为M

步骤106：根据反演结果，得到未爆弹的三维立体成像。

本发明实施例中，得到正则化方法的反演结果之后，再次进行运算即可 得到未爆弹的三维立体成像。

具体的，正则化方法的反演结果为一个三维矩阵，该三维矩阵中每个元 素的值代表一个长方体单元的物性参数，对三维矩阵分别沿x、y、z轴方向 切片；再对切片后的三维矩阵进行二值化处理并做三维显示，即可得到未爆 弹的三维立体成像。

通过上述实施例，本发明的方法利用梯度张量数据消除背景场的影响， 使得第一先验信息对应的边界识别效果更好。而将第一先验信息加入到正则 化方法中的目标函数中，再对该目标函数的最优值进行运算，得到的正则化 方法的反演结果，相较于目前的探测方法的精确度更高，因此最终得到的未 爆弹的三维立体成像结果也更精确，对未爆弹的识别效果自然得到了提升， 从而提高了遗弃化学弹发掘工作的工作效率。

以下，对上述第一先验信息如何加入目标函数中进行说明。

首先，基于梯度张量数据，采用欧拉反褶积法进行运算所基于的欧拉反 褶积公式为：

其中x

基于上述欧拉反褶积公式，进行运算过程中，选取多个相邻数据点组成 一个滑动窗口，使用滑动窗口内的数据求解欧拉反褶积公式组成的线性方程 组，得到异常体中心坐标；随着滑动窗口的滑动，计算出多组异常体中心坐 标，其中包括的异常解使用主体异常距离准则和反演结果的聚散度准则进行 筛选；去除异常解后，剩余正常解的分布即可作为第一先验信息。

本发明实施例中，假设多组异常体中心坐标包括：Q

Q

其中，表示转置。

本发明实施例中，正则化方法反演的目的，是对目标函数求解一个最优 解，目标函数的表达式为：

Φ＝||W

其中，d是一个N维矢量，表示在地表观测到的某类磁异常，包括：磁 总场、磁场分量或者是其它种类的数据，m是一个M维矢量，表示地下每个 长方体块的物性参数，包括：磁化强度、磁化率或者其它物性参数，A是一 个N×M维的矩阵，表示观测数据d与物性参数m之间的正演算子，矩阵中每 个元素的值根据正演公式进行求解，W

本发明实施例中，将目标函数表达式中的模型加权函数矩阵的表达式定 义为：

W

则模型加权函数矩阵的表达式中的W

其中，W

本发明实施例中，将第一先验信息加入到正则化方法中的目标函数中， 其具体包括：将第一先验信息加入到模型加权函数矩阵中。为了实现这个目 的，可以：

令第一表达式为：

令第二表达式为：

其中，σ

将第一表达式作为异常体中心坐标代入W

以下针对上述未爆弹三维立体成像方法的正确性、实用性通过模拟仿真 进行举证。

基于电磁场仿真软件CST进行仿真，用亥姆霍兹线圈模拟地磁场，对 地下40cm处的遗弃化学弹进行仿真，得到地表的磁场信息。假设遗弃化学 弹长度为1.994m，弹体切面的直径为0.183m，仿真得到的初始数据为磁场 三分量数据。再根据磁场三分量数据运算，得到梯度张量数据和磁总场梯度 数据，对梯度张量数据采用欧拉反褶积法进行运算，并且根据主体异常准则 滤除距离窗口过远的点，反演得到的磁异常中心水平分布图如图2所示，图 2为俯视图，其中横坐标y代表东西方向的坐标，纵坐标x代表南北方向的 坐标，图2中每一个圆圈点都代表一个解。前述已知，欧拉反褶积要通过滑 动窗口的方式获得很多个解，每个解都代表对异常体位置的一个估计。得到 很多个解之后，用这些解的统计学特性来表示异常体的位置的。由此可以知 晓，在原点附近有一系列反演得到的中心点呈柱状排列，可以在一定程度上 反映异常体的水平形态和位置。

对图2得到的结果用聚散度准则进行筛选，去除孤立存在的点，得到图 3所示的示意图。图3为俯视图，其中横坐标y代表东西方向的坐标，纵坐 标x代表南北方向的坐标。实质上，图3所示的解就相当于是把图2所示的 解中不好的解去掉，图3可以准确反映遗弃化学弹的水平位置和形态，该图 反映出第一先验信息，后续可以将第一先验信息添加到模型加权函数中。

另外，对欧拉反褶积的结果进行深度统计，得到图4所示的统计示意图， 图4是为了显示竖直位置的统计学特性，图4中横坐标表示深度，纵坐标表 示处在这个深度区间内的解的数量，由此可以知晓，异常体中心坐标的深度 主要分布在0.4m附近，与原始模型一致。

将第一先验信息添加到模型加权函数中后，根据第二先验信息和磁总场 梯度数据，对目标函数的最优解进行运算，得到正则化方法的反演结果，对 反演结果进行二值化之后进行三维显示，结果如图5所示，图5中，以x、y、 z三轴坐标显示遗弃化学弹的三维成像，结果表明本发明的未爆弹三维立体 成像方法的正确性、实用性。

基于上述未爆弹三维立体成像的方法，本发明实施例还提出一种基于正 则化方法的未爆弹三维立体成像的装置，参照图6，示出了本发明实施例一 种未爆弹三维立体成像的装置的框图，所述装置包括：

获取三分量模块610，用于获取目标探测区域的磁场三分量数据；

张量和总场梯度模块620，用于根据所述磁场三分量数据，得到梯度张 量数据和磁总场梯度数据；

先验信息模块630，用于根据所述梯度张量数据，得到第一先验信息；

加入模块640，用于将所述第一先验信息加入到正则化方法中的目标函 数中；

反演模块650，用于根据第二先验信息和所述磁总场梯度数据，对所述 目标函数的最优解进行运算，得到所述正则化方法的反演结果，所述第二先 验信息为已知先验信息；

成像模块660，用于根据所述反演结果，得到所述未爆弹的三维立体成 像。

可选地，所述先验信息模块630具体用于：

根据所述梯度张量数据，采用欧拉反褶积法进行运算，得到所述第一先 验信息；

可选地，所述反演模块650具体用于：

根据第二先验信息和所述磁总场梯度数据，利用共轭梯度法对所述目标 函数的最优解进行运算。

可选地，所述反演结果为一个三维矩阵，所述三维矩阵中每个元素的值 代表一个长方体单元的物性参数，所述成像模块660包括：

切片单元，用于对所述三维矩阵分别沿x、y、z轴方向切片；

二值化单元，用于对切片后的三维矩阵进行二值化处理并做三维显示， 得到所述未爆弹的三维立体成像。

基于上述未爆弹三维立体成像的方法，本发明实施例还提供一种基于正 则化方法的未爆弹三维立体成像的系统，所述系统包括：探测仪，所述探测 仪用于执行上述步骤101～步骤106任一所述的方法。

综上所述，本发明利用梯度张量数据消除背景场的影响，使得第一先验 信息对应的边界识别效果更好。而将第一先验信息加入到正则化方法中的目 标函数中，再对该目标函数的最优值进行运算，得到的正则化方法的反演结 果，相较于目前的探测方法的精确度更高，因此最终得到的未爆弹的三维立 体成像结果也更精确，对未爆弹的识别效果自然得到了提升，从而提高了遗 弃化学弹发掘工作的工作效率。

还需要说明的是，在本文中，诸如第一和第二等之类的关系术语仅仅用 来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗 示这些实体或操作之间存在任何这种实际的关系或者顺序。而且，术语“包 括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括 一系列要素的过程、方法不仅包括那些要素，而且还包括没有明确列出的其 他要素，或者是还包括为这种过程、方法所固有的要素。

上面结合附图对本发明的实施例进行了描述，本文中应用了具体个例对 本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解 本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发 明的思想，在具体实施方式及应用范围上均会有改变之处，综上所述，本说 明书内容不应理解为对本发明的限制。