马小陆 ,巩朝光 ,王兵 ,郑睿
(1.安徽工业大学 电气与信息工程学院,安徽 马鞍山,243002;2.安徽师范大学 物理与电子信息学院,安徽芜湖,241002)
单线激光雷达特征提取对移动机器人的目标追踪、目标提取和特征融合等具有重要作用,如何准确且可靠检测出点云特征成为研究热点。国内外学者针对单线激光雷达特征提取已经做出了大量研究。例如,Siadat[1]等在基于点间距离分割(Point-Distance-Based Segmentation,PDBS)算法基础上提出了以相邻点之间的差值提取直线的连续边沿追踪(Successive Edge Following,SEF)算法,刘朋[2]等提出使用相邻点之间的斜率差提取环境特征。Arras[3]等提出线性回归(Line-Regression,LR)算法,实现移动机器人在地图上的定位。Borges[4]、Zhang[5]等在分割—合并(Split-and-Merge,SM)算法[6-7]基础上,提出了迭代端点拟合(Iterative End Point Fit,IEPE)方法提取环境中的直线特征,该方法具有速度快和准确度高的特点,但由于该方法对阈值敏感,所以容易出现分割、合并效果不好的现象,因此满增光[8]等将IEPE 方法融合端点拟合(End Point Fit,EPF)和总体最小二乘方法,有效降低了发生欠分割或过分割的可能性,并且将该方法应用于室内AGV 地图创建。除此之外,张岩[9]等提出了随机采样一致(Random Sample Consensus,RANSAC)算法,该方法可以较为稳健的进行特征提取,但是特征提取正确率不高[10]。Fernandes[11]等、李万莉[12]等将霍夫变换(Hough Transform)应用于单线激光雷达特征提取上,邢亚蒙[13]等提出了改进的Hough 变换算法,使用激光扫描信息提取墙面特征,精度较高,但对噪声比较敏感,稳定性差。
综上,当前的特征提取方法存在着对噪声敏感、稳定性差等问题。针对上述问题,本文提出了一种基于前瞻窗口的单线激光雷达特征提取方法。
本文的单线激光雷达特征提取方法包括双边滤波、Harris 角点提取、基于前瞻窗口的区域搜索,相交直线拆分。激光点云原始数据使用极坐标表示为
本方法在笛卡尔坐标系下完成,转化到直角坐标系:
un为在笛卡尔坐标系下的扫描点,本文中N取909。
由于激光雷达在测距过程中受到传感器噪声、环境噪声以及物体对激光的反射率影响,需要对点云数据进行预处理,这有利于提升测量的准确性。预处理的方法要求满足在过滤噪声的同时又能很好的保护图像的边缘。激光雷达噪声主要变现为非线性噪声,常用滤波的方法主要有中值滤波、双边滤波以及小波变换[14]等方法,其中双边滤波[15]具有非线性、噪声抑制、数据平滑等特性以及更好的边缘保护特征,所以本文选择双边滤波。
双边滤波定义为
式(3)中权重Wp定义为
式(4)中,Ifiltered(x)表示当前滤波的扫描点深度距离值,Ω表示当前窗口区域,xi表示窗口内非窗口中心的其他扫描点位置,I(xi)为其对应的深度距离值,fr表示平滑值域中强度差异的核函数,gs表示平滑空间域差异的核函数。本文中使用高斯函数作为值域、空间域核函数
每个扫描点权重
式(7)中,σd为空间域核函数方差,σr为值域核函数方差,本文分别取3、0.1,对式(7)进行归一化后得
式(8)中Id为滤波后的扫描点,重复以上过程对扫描点进行滤波。
使用激光雷达真实扫描到的数据进行测试,结果如图1 所示。
图1 橙色散点表示滤波后的点云,蓝色点为原始点云,原始点云受到反射率或噪声影响后,物体轮廓点云出现突变、波动,经过双边滤波后波动得到了抑制,获得了平滑的点云。验证了双边滤波对激光雷达滤波的有效性。
图1 点云双边滤波效果图
环境中的特征点主要由角点、直线、圆弧、断点等组成,首先对角点进行定位,再进行特征线的提取,有利于提升检测的速度与准确性。Harris 角点检测方法通过计算扫描点的曲率与梯度来检测角点,具有良好的噪声鲁棒性、稳定性[16],下面给出计算详细过程。首先计算扫描距离点u(x,y),在X和Y方向上的2 个梯度
计算梯度乘积
对ux2、ux2、uxy使用σ取1的高斯函数加权得
式(11)中w为高斯加权函数。图像u(x,y)在点(x,y)处平移(Δx,Δy)后的自相关函数可以近似为
将自相关函数表示为椭圆方程:
椭圆的扁率与尺寸由M(x,y)的特征值决定,与物理特征对应关系如表1 所示。
表1 环境特征与特征值λ1、λ2的对应关系
将激光点云的投影点作为待处理点,这时方法沿着投影点轨迹进行计算,而非对整个二维平面进行计算(图2),提升准确度的同时,提升了检测效率。
图2 沿序列位置角点提取示意图
计算每一个投影点的角点响应值
式(15)中α为经验常数。角点位置响应值阈值TR=pr*Rmax,pr为角点响应比例系数,Rmax为遍历一帧点云后计算得到的最大角点响应值。
对角点响应值进行局部非极大值抑制,取半径为d的窗口,满足条件
取α=0.06、pr=0.000 2、d=2,获取单线激光雷达在室内扫描的点云,使用经过2.1 节滤波处理后,提取点云中的角点,结果如图3 所示。图3中红圈处为角点位置,蓝色为初始点云。
由图3 可知,该方法可以准确地提取点云中的特征角点。
图3 角点提取结果
为提取到完整的特征,使用基于前瞻窗口的搜索方法,以提取到的角点位置为起点,计算当前点Pcurrent(x0,y0)与下一点Pnext(x1,y1)组成直线的斜率ki,若|ki|小于阈值,将Pnext归为当前类,同时将Pnext更新为Pcurrent,若|ki|大于阈值,则进入前瞻窗口搜索模式。前瞻窗口搜索原理如图4 所示,图中P0~P9属于同一特征类,P5、P6为跳跃断点或噪点,假设搜索窗口宽度为3,以P0为起点搜索,P0至P4处的|ki|均小于阈值,当Pcurrent为P4,Pnext为P5时,|ki|大于阈值,进入前瞻窗口搜索模式,取进入窗口搜索前的点为Pcurrent,即将P4设为Pcurrent,在窗口内依次取点作为Pnext计算|ki|值。图中P5、P6与P4点计算的|ki|均大于阈值,将下一点设为Pnext,且P7与Pcurrent所在直线的斜率小于阈值,此时退出窗口搜索模式,并更新Pcurrent为P7,Pnext等于P8,继续根据斜率继续搜索。图4 搜索完成后,归为此类的点集为{P0,P1,P2,P3,P4,P7,P8,P9},完成了特征点聚类。
图4 前瞻窗口搜索、聚类示意图
若搜索完整个前瞻窗口,均未出现斜率小于设定阈值的点,将Pcurrent点标记为特征终止端点,示意图如5 所示,图5搜索后对应的点集为{P0,P1,P2,P3,P4,P5}。按照同种方式向前继续搜索,获得特征点云的起始点。对每一个角点进行前后方向搜索,若角点处于之前搜索得到的点云中,则跳过此角点的搜索,遍历完角点后,获得各类特征点云。
图5 特征终止位置前瞻窗口搜索示意图
完整的特征提取流程图如图6 所示。图6 中,i为点云投影点数组ranges 下标,num_max 表示点云总个数。变量ranges_last 用于暂存符合条件的距离参考点,搜索过程中对计算当前点Pcurrent(x0,y0)与下一点Pnext(x1,y1)所在直线斜率为
图6 特征线区域搜索流程图
若斜率|ki|不小于设定阈值T_k,进入前瞻窗口扫描模式,否则继续搜索,搜索过程中|ki|计算方式为
式(18)中,Ps(xs,ys)为进入窗口扫描模式前的一个扫描点,依次计算窗口内扫描点与Ps(xs,ys)点所在直线的斜率绝对值,若在前瞻窗口宽度Tw内|ki|大于T_k的点个数(cnt)等于前瞻窗口宽度,则Ps(xs,ys)为特征的边界点,若在前瞻窗口范围内|ki|小于T_k,则退出前瞻窗口搜索模式,并将i处的点坐标更新为Ps(xs,ys),且将i点聚类到当前特征点集内。对提取到的角点重复以上搜索步骤,获取点云中存在的特征点集。沿着投影点轨迹上下搜索,找出点云特征起始位置,搜索判定条件总结如表2 所示。
表2 前瞻窗口搜索判定条件
向下搜索与向上搜索原理一致,重复此过程对每一个角点进行上下区域搜索。若搜索角点处于已搜索到的特征点云中,则放弃对当前角点搜索,并选取下一角点进行处理。使用上述1.2 节提取到的角点为中心进行搜索,结果如图7 所示。
图7 中使用方框圈出区域搜索后的角点位置,经前瞻窗口搜索后,特征的边界点被准确的提取了出来,且对噪声不敏感,稳定性较好。
图7 方法角点提取结果
实际的物理环境中,存在着弧线、直线、折线等特征,需要将这些特征有效的区分开。对于折线,需要将两相交直线进行分割,拆分出两条独立的直线,从而有利于后续的拟合处理。首先创建一个长度为2m+1的窗口(文中m取3),分别在前后两半窗口中取值,并计算相对于窗口中心点的斜率
则两段线段的夹角
式(20)中,k1=(y0-ym)/(x0-xm),k2=(ym-y2m+1)/(xm-x2m+1)。按照此方法遍历完整个区域,计算窗口内线段的夹角,若夹角大于设定阈值,则窗口中心即为交点所在位置,此时的θ近似等于两条直线的夹角。以交点位置为分界点,将相交的特征进行拆分,最终提取结果如图8 所示。结果表明,该方法可以将特征有效、稳定地提取出来,且对于相交的直线也可以准确进行分割。
图8 特征拆分后提取结
为验证算法的实用性与可靠性,将算法移植到自制的机器人平台,结构示意图如图9 所示。机器人平台包括主体支撑、执行机构、处理单元、电源管理以及传感器部分。机器人平台搭载了ROS Kinect 版本系统,配备了光电码盘、基于三角测量的激光雷达[17]、RGBD 摄像机等传感器。平台前轮平衡结构保证了机器人平台的平稳,可减少地面环境对激光雷达等传感器的干扰。相关软硬件配置如表3所示。
图9 机器人结构示意图
表3 机器人软硬件配置
实验选取了实验室环境进行测试,环境中存在着桌子、墙壁等结构化特征,也存在着废旧电缆、凳子腿等非结构化特征。实验从距离测量、长度测量、提取成功次数角度对方法特征提取效果进行分析,并与准确度较高的Hough 变换方法进行对比,验证本方法在特征提取上的准确性与鲁棒性。
首先对方法的滤波效果进行验证。提取真实环境中的特征,环境如图10(a)所示,方法处理完成后,将提取到的特征以点云形式发布出去,如图10(b)所示。
图10 特征提取实验环境及结果
图10(b)中,取一段特征直线“LA-LB”,计算特征区域内雷达中心到扫描平面的距离,并在同一时间戳下计算使用相同区域内初始点云的扫描距离值。重复上述测试过程100 次,记录每次获取的距离值,统计结果如图11 所示。图11 中,蓝色线为每次使用原始点云数据获取的距离值,红色线为使用该方法处理后的测距结果。由图11 可知,使用该方法对点云数据处理后,与原始点云测距相比测距结果稳定,数据更为准确。对比2 次实验结果,统计两者的测量均值、标准差,并以均值为参考值统计绝对平均误差,如表4 所示。由表4 可知,使用该方法处理后的点云,在测距方面结果更为稳定,方法有效提升了点云距离扫描值的精确性。
表4 点云测距结果统计
图11 该方法输出点云数据与原始点云数据测距结果对比图
为验证方法的稳定性与准确性,首先对点云进行滤波,再分别使用本方法与Hough 变换方法对同一场景下的环境提取特征,结果如图12(a)~(b)所示。
图12 使用本方法与Hough 变换方法特征提取结果对比
图12(a)中,橙色圆点为提取到的特征起始端点,图12(b)中Hough变换方法提取到特征的起始点使用黑色圆点进行标记,为便于统计使用“L_A”、“L_B”、“L_C”…“L_Q”分别对两种方法提取到的特征进行编号。
将本方法与Hough 变换方法提取结果进行对比,在同一帧激光数据下,统计在图12 中编号位置2种方法提取情况,若提取到相应特征则对应特征的成功提取次数加一,并计算2种方法的特征起始点序列号均值。重复上述检测、统计过程100 次,统计结果如表5 所示。
表5 2 种方法统计结果
续表5
表5 中,将成功提取次数转换成百分比,使用柱状图进行表示,如图13 所示。
由图13和表5可知,两种方法提取到的特征范围结果相近,表明两种方法提取精确度相近。但是本方法特征提取稳定性明显高于Hough 变换方法,且对于“L_A”、“L_B”、“L_E”等非结构化特征也能够保证较高的正确率。实验表明本方法对于单线激光雷达特征提取具有较高的鲁棒性、稳定性。
图13 2 种方法特征提取成功率统计图
为对比特征提取的准确性,分别使用本方法与Hough 变换方法提取到的特征点云测量特征“L_I”的长度,根据提取到的特征起始点n1、n2查找所对应的扫描距离值r1、r2,待测量物体在扫描平面的长度
式(21)中,αi=angle_delta*ni+angle_start 其中angel_start、angel_delta 表示角度起始值与角度增量,实验中分别取-π/2、π/909。获取100 帧不同时刻点云扫描数据,使用两种方法计算每一次特征“L_I”处对应的特征长度,并与手工测量结果(746 mm)进行对比,如图14 所示。
图14 中,使用红色线表示使用本方法长度测量结果,蓝色线表示Hough 变换方法测量结果,绿色线表示手工测量结果。统计测量结果的均值、标准差,以手工测量值为参考计算绝对平均误差以及均值偏差百分比,如表6 所示。
结合图14 和表6 可知,本次长度测量实验中,本方法特征提取精度较高,数据波动较小,平均长度与真实长度仅相差1.5 mm,而使用Hough 变换方法提取结果与真实值相差5.7 mm。且由表6 可知,本方法测量结果的标准差、平均误差均低于Hough变换方法,表明本方法具有更高的稳定性、准确性。
图14 使用该方法处理后点云长度测量与手工测量对比图
表6 长度测量对比结果统计
目前的单线激光雷达特征方法,仍存在着对噪声敏感、稳定性不高的问题,导致机器人在单线激光雷达特征提取上效果不够理想。本文针对现有特征提取方法存在的问题,提出了基于前瞻窗口的特征提取方法,方法中使用了前瞻窗口搜索,有效提升了特征提取的稳定性、准确性,且对噪声不敏感。最后将方法应用在真实机器人平台,实验结果表明,特征提取稳定性、准确性较高,对环境适应性更好,满足机器人特征提取需求。