一种基于多传感器数据融合探测地下管网的方案

胡云斌,周熙人,陈欢欢

(中国科学技术大学 计算机学院,合肥 230026)

1 引 言

地下管网的管理在当今的城市发展中十分重要.电力、水、煤气、污水和电信信号通过地下管网进行传输.这些管道经常需要在非开挖的情况下进行准确定位,而相关的地球物理探测传感器技术已经有很多[1].因为不同传感器适用于不同的情况,所以研究者们研究了多传感器探测设备[2].同时通过程序智能分析传感器数据获得管道信息来降低人力成本和时间也是当前趋势.

已有很多关于多传感器设备的研究成果[3-6].Li 等人[4]展示了一种3G系统来可视化测量管道结果的不确定程度,并研究了探地雷达探测的误差梯度和埋藏物体的深度之间的关系.但是它没有涉及数据智能分析领域.Hafsi等人[5]展示了一种电磁场传感器和探地雷达结合的解决方案.他们设计了算法先解释探地雷达图像并给出管道可能存在的空间点,然后融合多个空间点分析出管道走势.Dou等人[6]在一个多传感器设备上设计了一种多横截面前进算法来实现地下管道的精确定位.设备包含五种传感器:探底雷达,被动电磁场探测仪,磁场梯度仪,低频电磁场探测仪和声振传感器.

本文提出了一种快速探测方案能探测地下管网且在地图上实时地绘制最新结果.该方案是针对上述已有方案的不足而提出,具体对比在下一小节给出.该方案建立在一辆带有被动电磁感应线圈和探地雷达的探测车之上,并包含了如下内容:设计扫描轨迹的形状,设计解释算法解释一段扫描轨迹上的传感器数据,还有设计数据融合方法得到管网位置和方向.本文主要的贡献如下所示:

1)提出扫描线的 “Z” 形布局方式,因为它在管道方向未知时比平行线布局更好,并且扫描路径成本不高于两倍.

2)提出了一种基于自修正和筛选的管道数据融合算法用于处理由噪声引起的虚警空间点问题.解释算法给出的空间点既有反映真实管道的,也有误报的,该算法基于管道的线性特征分析不同假设的概率并筛选出高可信的假设管道.

3)提出了一种简单的利用被动电磁感应线圈来探测地下电缆的方法.结合 GPS模块,设备只要沿着一条扫描轨迹扫描一遍,就可以分析出地下电缆的位置和方向.

2 探测问题描述

下面介绍本文研究的管网探测问题的细节.假设在地下浅层处水平埋藏着多根管道,管道连续但材料未知.有传感器能够探测地下管道位置,但工作模式比较特殊——设备沿着地上某条和地下管道在水平方向交叉的扫描轨迹前进,解释算法分析采集到的数据并给出交叉点的位置和交叉处管道的切线方向.为了得到某管道的走势,需要沿着管道的方向完成多条这样的扫描轨迹,融合多个交叉点信息分析走势.假定解释算法给出了一些空间点,空间点和真实交叉点之间存在三种情况:1)真实情况:空间点反映了交叉点位置.2)虚警情况:空间点标注的位置实际不是交叉点.3)缺失情况:某个交叉点没有被探测到,故没有空间点反应它.

图1从俯视角度说明了探测问题.图中细长条表示地下的真实管道,带箭头的线表示地上设备水平移动走过的扫描线,三角形表示真实管道和扫描线的交叉点.图中多条扫描线的布局没有规律,沿着一条扫描线采集的传感器数据会被解释得到空间点,空间点的三种情况在第1、2、3条扫描线中展示.宽长条是通过某种算法分析所有空间点得到的一个假设管道,和真实管道间有误差.

图1 探测问题示意图Fig.1 Diagram of detection problem

扫描线的布局关系到数据采集的效率和质量,Hafsi等人的方案采用网格式布局,而Dou等人采用了平行线布局.为了获得更好的效果,本文提出了一种扫描线的“Z” 形布局.

解释算法需要针对不同传感器而特殊设计.比如探地雷达返回的是二维扫描图像,一般地下管道会以双曲线的形式表示.关于从图像中自动识别和拟合双曲线的问题,已有一些发表的成果[7-9].地下管道的深度和直径能通过拟合双曲线的参数来获得[10].而对于电磁感应类的探测仪器,目前主要以人工读数并分析地下管道位置为主[11].为了实现自动化分析的过程,本文基于被动电磁感应线圈提出了一种探测地下电缆的方法并给出了解释算法.

最后一个问题是设计算法融合多个空间点分析管道走势.Hafsi等人的方案直接假设管道方向是已知的,所以所设计的算法实用性不强.Dou等人的方案提出了多横截面前进算法,其中扫描横截面是以竖直方向和水平扫描线线段形成的平面.它结合了扩展卡尔曼滤波(EKF)方法来自动推测管道如何从当前扫描横截面延伸到下一个横截面.该算法能处理真实空间点的位置误差问题,但是它无法处理虚警空间点问题.为此,本文设计了一种基于自修正和筛选的管道数据融合算法专门处理虚警空间点问题.

针对上述三个方面,本文的几个成果将在下文给出具体描述.其中解释算法通过真实环境实验的实验数据来说明其细节.

3 扫描线布局问题

设θbs是一根管道和一条扫描线之间的夹角.根据本文设备上的被动电磁感应线圈和探地雷达的使用经验,θbs大于 60 度时,传感器在交叉点处才有明显信号变化.定义θbs大于60度的交叉点为有效交叉点,它的空间点不会缺失.在图2中,细长条表示地下管道,总行进方向表示扫描线依次安置的方向,a和b的含义如图所示,定义布局密度d=b/a.每个图中两条地下管道是在该总行进方向下可被探测到的极限情况,方向在它们之间的管道都可被探测到,定义探测范围角θr.可以看到,在d=6时,平行线布局的管道上有效交叉点为10个/每根且θr= 60°,而“Z”字形布局的管道上有效交叉点为5个/每根(图中用小矩形框标注)且θr=60°+2arctan(1/d).虽然有效交叉点数在平行线布局上更多,但“Z”字形布局的θr更大.

图2 平行线布局和“Z”字形布局Fig.2 Parallel layout and zigzag layout

下面分析在总前进方向上探测单位长度时两种布局的扫描路径长度.首先假定一条管道和某布局的有效交叉点数k和夹角θbs,可以给出d关于k的函数:

(1)

(2)

用d可以算出在总前进方向上前进单位长度时扫描路径长度R:

(3)

Rparallel(d)=d(b+a)

(4)

如果设θbs=60°且b=1,可以做出公式(3)和公式(4)之比关于密度的曲线,可得平行线布局的成本永远低于“Z”字形布局,但是扫描路径成本后者不会高于前者两倍.

总结上面分析,“Z”字形布局θr更大,更适用于管道方向未知时的情况,同时k相同时它的扫描路径成本不会比平行线布局高两倍.

4 基于自修正和筛选的管道数据融合算法

4.1 算法的数学模型

在之前章节中已经描述了探测问题的细节,解释算法返回的空间点会有三种情况:真实的、虚警的和缺失的.本文的数据融合算法首先建立在下面的四个假设上:

1)地下管道几乎是直的.

2)虚警率>0,而缺失率为0.

3)真实的空间点在交叉点周围且包含高斯噪声.

4)虚警的空间点的出现遵循泊松过程.

第一条假设可以通过增加扫描线的密度来满足,第二条假设可以通过提高解释算法的灵敏度来满足.比如探地雷达的原始数据是二维图像,管道在图像中反映为双曲线形式,解释算法可以将尽可能多的区域判定为有双曲线,即使判据不那么明确.第三和第四条假设依据之前的研究[4]和合理推测而定.

第一条假设意味着一条地下管道和一条扫描线最多只有一个交点.选择n个相邻的扫描线并从每条扫描线中选择一个空间点,然后拟合一条直线.尝试n的值和选取的空间点的不同组合会得到很多拟合线,如果拟合线所用空间点全部是真实的,则代表了地下管道的真实情况,定义这种线是“正确的”,而其它的是“错误的”.但是空间点真实与否是未知的,这里可以设计一种筛选条件将拟合线分成不同可信度,其中低可信的拟合线所用的空间点有大概率是虚警的.

下面给出了一种合适的筛选条件.假设一条管道具有N1个与扫描线的交叉点,并且沿着其扫描线,第i个交叉点的真实空间点的坐标为Xi(Xi～(μi,σ),其中μi是交叉点的位置).对于有n个点的拟合线,通常n的数量越大且所用到的空间点距离拟合线的距离方差越小越能说明拟合线可信,筛选条件可以如下:

(5)

先用Pin来表示空间点的X满足|X-μ|

Pin=P(-kσ≤X-μ≤kσ|X～N(μ,σ))
=P(-kσ≤X≤kσ|X～N(0,σ))

(6)

公式(5)中的阈值参数需要确定,这里通过两个概率的方程来确定.首先是拒绝管道的所有正确拟合线的概率P1.正确的拟合线用到的都是真的空间点,它被筛选条件拒绝等价于nkσ.由于从N1个真的空间点中取出一些的所有取法都将被枚举以创建拟合线,所以P1等价于在N1个空间点中,有不多于N0-1个满足X

(7)

如果N1=N0,可以简化P1:

P1=1-(Pin)N0

(8)

(9)

如果给定环境参数λ=0.3,k=3,σ=1/3,且P1<5%∧P2<5%,可以得到6.32

(10)

4.2 构造空间点和管道的对象结构

对于一个空间点,使用检测向量(DV)来表示其信息:

DV=

其中(x,y,z)为三维坐标,r为地下管道的半径,(dx,dy,dz)为管道的方向.管道对象的结构是:

1)DN:空间点的数量;

2)MVP:空间点到拟合线的距离的均值;

3)MDA:空间点到拟合线的角度偏差的均值.

DN=n

(11)

(12)

(13)

筛选条件在下面将被称为“边界”.一个边界可以通过结合这三种评价指标来产生.它可以是逻辑结合:

boundary=DN>N0∧MVP≤P0∧MDA≤A0

(14)

或者一个运算结合:

(15)

与数学模型中不同,这里增加了A0,但都可以以类似于公式(10)的思路设置.

4.3 算法流程

包含该算法的整体流程在图3中说明,左半部展示了探测过程中操作者的工作,而右半部展示了数据处理算法流程,其中数据融合算法以一条扫描线为周期迭代运行.算法构造的管道对象根据其可信度水平分为三个层级(按升序排列):L1,L2,L3.相应的需要三个边界来判断管道是否符合此层级.另一方面,为了说明方便,定义了管道上的二元关系:

pipe1 conflictpipe2≡
∃DV0:DV0∈pipe1.DVList∧DV0∈pipe2.DVList

一条管道是“无冲突”的表示它的所有空间点都不属于其它管道.融合算法具体分三步:

第一步是修改管道,具体工作是为已有的高可信管道对象提供新的空间点来优化其评估指标.如果一条“正确的”管道的延长线与最新扫描线相交,则在新扫描线上应有对应的真实空间点.L2和L3的管道被认为是高可信的,因此应该找到它们的对应空间点.为了估计一个空间点对管道的匹配程度,需要设计估计匹配度的方法.一般应该基于两点:空间点和管道之间的空间距离以及空间点的方向和管道方向的差异.具体操作使用一种双向选择方法:在新空间点的集合和待修改管道的集合之间,找到一对组合使空间点和管道相互最匹配.在修改该管道后从集合中删除此组合,然后重复上述操作直至一个集合为空.

图3 探测方案流程Fig.3 Procedure of detecting procedure

第二步是创建管道,目的是建立大量低可信(L1)管道作为候选者.L1边界会设置得很松,故L1管道是低可信的管道.新创建的L1管道的来源可以是一个空间点与已有L1管道的结合,也可以是单单基于一个空间点.算法会枚举所有可能来创建.L1管道集合内管道允许互相“冲突”.

第三步是通过升级产生新的L2或L3管道.要测试新创建的或修改过的管道所在层级的上一级边界,通过则升级.L1升L2要考虑两种情况.首先,L1管道允许“冲突”,但L2和L3不允许(如果一个空间点在高可信管道上,其它具有这个空间点的管道是正确的可能性很低).所以升级某L1管道后,所有与之“冲突” 的管道将被删除.这也删除了它们含有的很可能是虚警的空间点.其次,当两个或更多“冲突”的L1管道同时满足L2边界时,需要设定优先级来决定升级哪个.优先级值可以是公式(15)的值.而L2升到L3无需额外操作,因为L2和L3管道都是“无冲突” 的.

在图4中,带箭头的线表示扫描线,浅色的长条表示之前迭代得到的管道而深色的是当前迭代将生成的.第一个子图举例说明了双向选择过程.对于图中的下管道,A比B距离更近故而更匹配,但是不能组合A和下管道,因为上管道和A更匹配.使用双向选择方法后A和上管道配对,而B和下管道配对.修改后的管道表示为虚线长条.在第二个子图中,白色的空间点同时属于两个浅色管道.中间的深色管道是由浅色管道和新的空间点配对后创建.在升级步骤中,该管道会被升级,但是由于“冲突”问题,有两个管道被删除.在第三个子图中,新的空间点使得两个L1管道需要升级.直观地看,斜管道的空间点更散,所以它有低优先级.同样由于“冲突”问题,也有两个管道被删除.

图4 举例说明融合算法Fig.4 Examples for explaining the fusion method

5 实验分析

5.1 模拟实验验证数据融合算法

实验目的是验证算法在不同运行参数和环境参数下的性能.首先需要设置“Z”字形布局的密度和扫描线长度.而环境参数方面,需要设置高斯分布和泊松分布的参数,还要设定真实地下管道的位置和数目.在运行参数方面,需要设置区分管道可信度的三个边界的参数.表1展示了部分参数设置,此设置对以下所有实验都是一致的.

表1 模拟实验部分参数表
Table 1 Part of parameters of synthetic experiments

Z字形方向1(0.866,0.5,0)角度高斯(rad) (μ,σ) (0,0.5)Z字形方向2(-0.866,0.5,0)扫描线长度(m) 20 位置高斯(m)(μ,σ)(0,0.5)真实管道1的斜率10 真实管道2的斜率(如果有)-10

表2 边界1参数表
Table 2 Parameters of the 1sttype of boundaries

管道层级L1<1,2,π>L2<10,4.4,0.8>L3

在表2所示的边界参数时,针对不同真实管道条数,分别测试了不同泊松分布参数λ时的算法性能.比如真实管道为1条时,使λ=0,0.1,0.2,0.3,算法一直能找到正确的拟合管道,而错误被接受的管道数量依次为0,0,0,1.随着λ的增加,算法会无法找到正确的拟合管道,这里定义能找到时最大的λ为识别上限(DUB).这里DUB为0.6.对于真实管道有两条的情况也做了实验.在图5中(基于ArcGIS平台[13]),左边四列是1

表3 边界2参数表
Table 3 Parameters of the 2ndtype of boundaries

管道层级L1<1,2,π>L2<4,4.5,0.8>L3

条真实管道时λ=0,0.3的结果,右边四列是2条管道时λ=0,0.3的结果.同图1一样,图5是俯视图,其中锯齿状的是扫描线,扫描线上的许多线条是空间点,细虚线是真实管道的位置.在每一列中,设备沿着扫描线从图中下方到达上方,过程中算法迭代分析获得实时结果,在设备完成最上方一段扫描线后,输出的管道位置在图中用粗管道表示,最好的管道和真实管道的位置误差和其它多出来的管道条数反映了算法性能.在表3的边界设置下也开展了实验,汇总情况在表4说明.

通过上面的实验,可以发现:其他条件不变的情况下,设置L2边界的N0更大时,算法挑选管道的出错率更低,但是其DUB会相应下降.设定其他环境参数和边界相同,真实管道的数目增加时,DUB有较大幅度下降.

表4 实验结果汇总Table 4 A summary of experiment results

图5 (边界1)λ和真实管道数目不同时的实验结果Fig.5 (1st type of boundaries)Results when λ and the number of real pipes are different

通常,如果真实空间点的位置误差满足N(0,σ),P0在2σ和3σ之间最好.一般σ越小,DUB就越高.对于A0,它应该比P0更紧,因为MDA是在MVP之后使用的强化约束.如果随着传感器技术的发展使得可以得到更小的σ,则算法可以通过设置更严格的边界来适应这一进步.

5.2 真实环境实验

实验设备为一辆带有被动电磁感应线圈和探地雷达的探测车.图6展示了其实物图和俯视图.实物图中四组被动电磁感应线圈用圆圈标注,同时探地雷达天线用矩形标注;俯视图展示了传感器的布局,两个互相垂直的线圈组成一组,四组线圈标号G1-G4,G2位置上固定有一个GPS模块,正中间是探地雷达的天线.被动电磁感应线圈是2000匝的铜丝线圈,直径约10cm,可以探测空间某点的交变磁场强度.GPS模块是带有实时定位差分功能的高精度工业级GPS模块.探地雷达为GSSI SIR30型号,频率200MHz.

图6 设备的实物图和俯视图Fig.6 Practicality picture and top view of the device

图7显示了真实探测实验的结果.第一个子图包含电源线和水管的测试场地(左边的是电缆而右边的是水管).第二个子图显示了得到的空间点和扫描线,棱形点表示扫描线,浅色空间点来自探地雷达而深色的来自被动电磁感应线圈,还有一条细线表示真实管道.每条扫描线的长度接近10m.设置L2边界中N0=3和A0=1.7,而其他参数与模拟实验中的设置相同.由于两条真实管道非常接近,因此只使用一条细线表示.最后一个子图展示了算法得到的结果.

图7 真实探测实验Fig.7 An experiment in the real environment

根据实验数据,解释算法的过程介绍如下.

为了说明被动电磁感应线圈的解释算法,用图7中间图椭圆框里扫描线的数据来举例.共八个线圈在一段扫描线上得到八条经过滤波和平滑处理后的时间-电势曲线.每两个互相正交的线圈为一组做分别处理.比如1,2号线圈组成G1组,这里设计了一种极大极小值算法处理1,2号曲线,具体是对垂直线圈分析向上梯度而对水平线圈分析向下梯度得到极大极小序列,并用极大极小序列生成评价指标评价某时刻处于电缆正上方的可能性,找到最可能的时刻.

图8展示了G1-G4四组线圈分别处于电缆正上方的时刻.结合GPS记录的空间轨迹,可以把四个时刻换成四个位置.在图9中斜率较小的线表示GPS的空间轨迹,四个实心点表示换算出来的位置.而GPS和其中一组线圈重叠,对于不在同一位置的其他三组线圈来说,其真实位置需要校正.图中每个矩形框的四个顶点表示同一时刻四组线圈的位置,空心点为校正后位置,斜率较大的线表示拟合出来的管道.利用拟合线的参数就能得到管道位置和方向.

图8 四组线圈在电缆正上方的时刻Fig.8 Times of four groups of coils being on the cable

雷达返回的数据是二维图像,已有很多算法能处理其图像[14],通常探地雷达的图像会在有电缆或管道处将显示开口向下的双曲线特征.本文的探地雷达双曲线识别和拟合算法来自合作作者周熙人的论文[15].空间点的位置结合辨别的双曲线顶点的位置和GPS的轨迹计算得到.空间点处管道的半径从双曲线方程得出.而方向被设置成与扫描线垂直,因为探地雷达无法推测方向.四个阶段分别是:原始扫描图像;二值化处理;聚类和模式识别算法处理估计双曲线参数;生成空间点.

图9 结合GPS获得管道位置和方向Fig.9 Obtaining the location and direction of the pipe with GPS

图7右边子图展示了算法输出的五条管道.最长的深色管道比较贴合真实情况,其斜率为0.047(真实管道的斜率为0.02),而其长度接近70 m,而距离真实管道的最大偏差约为0.96m.深色管道中的所有空间点都来自探地雷达,因此反映的较为可能是水管.

6 总 结

复杂的地下管道现状催生了对高效的管网探测技术的研究.为了节约人工成本,采用非开挖的探测技术是一种趋势;为了节约时间,程序智能分析传感器数据逐渐代替了专家分析.本文研究了实时的地下管道探测的情景,基于现有的探测技术,研究了一套数据处理算法从而完成了一套可行方案.设备开启后,以“Z”字形扫描方式匀速前进,对不同传感器会有不同的解释算法来将原始数据转化为统一的空间点形式,然后利用基于自修正和筛选的管道数据融合算法融合多空间点获得地下管网走势.

在介绍融合算法之前,引入了一个数学模型.基于地下管道是直的并且对不同情况的空间点建立合理的概率分布假设,然后据此设定筛选条件的参数.融合算法是一个迭代函数且标记具有不同可信度标签的假设管道,新的一批空间点的到来会触发一次新的迭代来修改现有的假设管道的属性和可信度标签.算法通过删除可信度低的管道来删除虚警空间点.但是融合算法还是有一些不足,下面给出了几个方向:

1)目前算法默认行进方向是单向的,这意味着设备不能来回探测多次采集一个地方的数据以提高精度.

2)实验中DUB不那么高,更高的DUB可以使算法具有更好的噪声容纳能力.

3)构建管道的模型,其操作涉及到寻找新点、评价分数两点.寻找新点是指完成一条扫描线并获得一些空间点后,对于一个管道寻找其中最优的空间点并延长管道本身;评价分数是指需要有一种机制评价管道的真实程度和可信程度.最简单的就是直线模型,也就是对于直线管道,其寻找新点按最小距离标准寻找,而评价分数可以用拟合管道所用的空间点关于这条管道的方差均值做评价.这也是本文所用的模型.可以修改直线模型而采用更好更复杂的模型比如可以用高阶曲线,或者滤波算法等等.也可以修改评价方式,用其它更好的方式来评价管道.

4)在一条扫描线上不同传感器给出的真实空间点可能来自同一个交叉点,所以可以增加一个环节先融合一些可能来自同一个交叉点的空间点为一个.一个思路是类比给假设管道划分可信度的方法,也可以对空间点也设定可信度划分.