Multi-sensor fusion positioning of detection robot for tailings pond flood discharge tunnel
-
摘要: 尾矿库排洪隧洞内环境复杂,现有探测机器人定位算法存在弱纹理室内场景定位失效的问题,难以适用于该类环境下探测机器人的精确定位和累计误差消除。针对上述问题,提出了一种适用于尾矿库排洪隧洞环境的探测机器人多传感器融合定位算法。该算法基于测程法和图论,利用ArUco码将排洪隧洞环境的远距离定位简化为多段短距离定位。首先利用UMBmark算法进行里程计标定,有效消除探测机器人轮径和轴径两类系统误差;然后利用扩展卡尔曼滤波(EKF)算法融合里程计和惯导模块(IMU)信息,利用测程法实现运动过程中探测机器人位置和姿态信息计算;最后利用ArUco码作为路标并固定于隧洞内部,探测机器人携带标定后的相机对ArUco码信息进行识别与处理,利用各传感器量测值形成约束,结合图优化方法实现位姿优化,并根据ArUco码携带的信息消除累计误差,从而实现狭长空间弱纹理场景下探测机器人的远距离高精度定位。在尾矿库排洪隧洞实际场景中运行和关闭多传感器融合定位算法,分别进行10组行进40 m的重复定位实验,结果表明:多传感器融合定位算法具有较高的稳定性和精确性,能有效校正累计误差,实现探测机器人在排洪隧洞环境中的精确定位,行进20 m的平均定位误差为19.77 cm,40 m的平均定位误差为21.23 cm,且行进20 m校正后误差均值为4.2 mm。Abstract: The environment in the flood discharge tunnel of tailings pond is complex, and the existing positioning algorithm of detection robot has the problem of positioning failure in weak texture indoor scene, which is difficult to be applied to the precise positioning and accumulative error elimination of detection robot in this kind of environment. In order to solve the above problems, a multi-sensor fusion positioning algorithm of detection robot for tailings flood discharge tunnel is proposed. Based on the odometry method and graph theory, the algorithm simplifies the long-distance positioning in the flood discharge tunnel environment to multi-segment short-distance positioning by using the ArUco code. Firstly, the odometer is calibrated by UMBmark algorithm, which effectively eliminates two types of system errors of wheel diameter and axle diameter. Secondly, the extended Kalman filter (EKF) algorithm is used to fuse the information of the odometer and the inertial measurement unit (IMU), and the odometry method is used to realize the calculation of the position and attitude information of the detection robot during the movement process. Finally, the ArUco code is used as a road sign and fixed inside the tunnel. The detection robot carries the calibrated camera to identify and process the ArUco code information. The robot uses the measurement values of each sensor to form constraints, and combines constraints with the graph optimization method to achieve position and attitude optimization. And according to the information carried by the ArUco code, the accumulative error is eliminated so as to realize the long-distance high-precision positioning of the detection robot in the narrow and long space and weak texture scene. The multi-sensor fusion positioning algorithm is operated and closed in the actual scene of the tailings pond flood discharge tunnel, and 10 groups of repeated positioning experiments traveling 40 m are carried out respectively. The results show that the multi-sensor fusion positioning algorithm has high stability and precision, can correct the accumulative error effectively and realize the precise positioning of the detection robot in the flood discharge tunnel environment. The average positioning error of traveling 20 m is 19.77 cm, the average positioning error of traveling 40 m is 21.23 cm, and the average corrected error of traveling 20 m is 4.2 mm.
-
0. 引言
尾矿库作为矿山安全生产的重大危险源之一,存在极大的潜在风险 [1-2]。对此,国家出台的《尾矿库安全监督管理规定》[3]和GB 39496−2020《尾矿库安全规程》[4]中要求定期组织尾矿库专项检查,对发现的事故隐患及时进行治理。尾矿库排洪隧洞主要用于泄洪和回收尾矿澄清水,其排洪能力和结构性能直接关系到矿山安全,因此必须对排洪隧洞进行安全巡检。目前排洪隧洞巡检以人工为主,因巡检线路较长,巡检过程存在安全风险,提出以探测机器人代替巡检人员完成巡检任务。探测机器人巡检需要在长度超过1 km的排洪隧洞内部精准定位隐患,所以有必要对探测机器人在巡检过程中的自主定位进行研究。
近年来,学者们对探测机器人自主定位进行了深入研究,提出了多种算法提高探测机器人实时位姿计算精度。文献[5-6]使用测程法进行定位,利用UMBmark算法修正探测机器人轮径和轴径两类系统误差对里程计定位造成的影响,提高了测程法定位精度。测程法易于实现且适用性广,但其鲁棒性和稳定性较差,且存在累计误差。文献[7]提出了一种基于概率估计的里程计模型。文献[8]使用有限记忆卡尔曼方法,以逐级滤波的方法融合多传感器信息,推算矿井机器人的航向和运动轨迹,减小了单一传感器测量误差的影响,但仍未解决累计误差问题,不适合长距离定位。文献[9]使用惯性模块(Inertial Measurement Unit, IMU)、里程计、全球导航卫星系统(Global Navigation Satellite System,GNSS)等传感器进行机器人定位,通过扩展卡尔曼滤波(Extended Kalman Filter, EKF)进行数据滤波,实现多传感器融合,提升了系统鲁棒性和稳定性,但累计误差的消除需依赖卫星提供的绝对位置信息,当卫星信号微弱时该类算法易因校正失败出现累计误差,受使用环境限制。文献[10]利用Fast−SLAM及包含激光雷达在内的多种传感器,利用基于粒子滤波的自适应蒙特卡洛算法构建环境栅格地图,实现了探测机器人在矿洞室内环境的精确定位,但Fast−SLAM算法在弱纹理场景中易出现粒子收敛缓慢情况,影响定位效果。文献[11-13]均利用图优化算法实现了探测机器人多传感器融合定位,并利用子图匹配实现回环检测,有效消除了定位过程的累计误差,但所用算法在弱纹理场景中匹配错误率较高,未考虑弱纹理场景的特殊性。文献[14-15]利用二维码阵列进行室内定位与导航,使得机器人能够在缺少纹理的室内环境实现精确定位,但均需要提前布置密集的二维码,且维护工作量极大,不适合大场景定位。文献[16]利用实验对比分析了目前常用标签的检测精度和效率,结果表明ArUco码的检测速度有明显优势,适用于对标签检测速度要求高的场合。
尾矿库排洪隧洞环境属于狭长室内环境,且隧洞内卫星信号弱,特征点少,现有探测机器人定位算法存在弱纹理室内场景定位失效的问题,难以实现此类环境中的精确定位和累计误差消除。本文提出了一种适用于尾矿库排洪隧洞环境的探测机器人多传感器融合定位算法。该算法将排洪隧洞环境的远距离定位以ArUco码为边界点,简化为多段短距离定位,利用EKF算法融合里程计和IMU信息实现短距离精确定位,利用预设的ArUco码作为路标消除EKF定位带来的累计误差,并实现位姿优化,最终实现狭长空间弱纹理场景下探测机器人的远距离高精度定位。该算法有效利用了多种传感器数据,避免了单一传感器定位误差,提升了定位鲁棒性。为验证探测机器人多传感器融合定位算法效果,设计了尾矿库排洪隧洞探测机器人系统实验平台,并在尾矿库排洪隧洞实际场景中运行和关闭多传感器融合定位算法,分别进行10组行进40 m的重复定位实验,结果证明了该算法具有较高的稳定性和精确性。
1. 多传感器融合定位算法
1.1 融合定位算法流程
探测机器人多传感器融合定位算法如图1所示。首先选取UMBmark算法进行里程计标定,有效消除探测机器人轮径和轴径两类系统误差。然后利用EKF算法融合IMU数据和里程计数据,并利用测程法实现运动过程中探测机器人位置和姿态信息计算。最后选用ArUco码作为路标并固定于隧洞内部,探测机器人携带经张正友标定法[17]标定后的相机对其进行识别与处理,利用图优化方式实现位姿优化,并根据ArUco码携带的信息实现累计误差消除,从而实现狭长空间弱纹理场景下探测机器人的高精度定位。
探测机器人车身上安装的里程计、IMU、单目相机数据作为算法输入,其中里程计需利用UMBmark算法进行标定,以消除系统误差。
工控机通过串口获取到传感器数据后,将里程计与IMU数据进行时间戳对齐,并且判断是否检测到ArUco码。未检测到ArUco码时,将里程计和IMU数据通过EKF算法实现融合定位。当检测到ArUco码时,对当前探测机器人位姿状态利用图优化方法进行局部优化。随后将优化结果和参考位姿信息结合,修正融合定位结果,消除累计误差影响,最后将定位结果输出。
1.2 里程计与IMU数据融合
考虑单一传感器容易出现测量错误,且依靠里程计推算得到的偏航角容易出现较大误差,因此将里程计数据和IMU数据利用EKF算法进行融合,避免单一传感器测量误差影响定位结果。
EKF算法的实现包括预测和更新2个步骤。设
$\gamma $ 时刻状态变量为$ {{\boldsymbol{x}}_\gamma } = {\left[ {{X_\gamma }\;{Y_\gamma }\;{\theta _\gamma }} \right]^{\rm{T}}} $ ,$ {X_\gamma },{Y_\gamma } $ 分别为$\gamma $ 时刻探测机器人二维坐标系的横纵坐标,$ {\theta _\gamma } $ 为$\gamma $ 时刻探测机器人偏航角。根据轮式里程计圆弧运动模型进行状态预测:$$ {{\mathop {\boldsymbol{x}}\limits^ \vee} _{{\gamma}} } = {{\mathop {\boldsymbol{x}}\limits^ \wedge } _{{{\gamma}} - {{1}}}} + \left[ {\begin{array}{*{20}{c}} {\Delta {d_{\gamma - 1,\gamma }}\cos ({\theta _{\gamma - 1}} + \Delta {\theta _{\gamma - 1,\gamma }}/2)} \\ {\Delta {d_{\gamma - 1,\gamma }}\sin ({\theta _{\gamma - 1}} + \Delta {\theta _{\gamma - 1,\gamma }}/2)} \\ {{\theta _{\gamma - 1}} + \Delta {\theta _{\gamma - 1,\gamma }}} \end{array}} \right] $$ (1) 式中:
$ {{\mathop {\boldsymbol{x}}\limits^ \vee} _{{\gamma}} } $ 为$\gamma $ 时刻状态变量的先验估计;$ {{\mathop {\boldsymbol{x}}\limits^ \wedge} _{{{\gamma}} -{{ 1}}}} $ 为$\gamma - 1$ 时刻状态变量的后验估计;$ \Delta {d_{\gamma - 1,\gamma }} $ 为$\gamma - 1$ 时刻到$\gamma $ 时刻的位移测量值;$ {\theta _{\gamma - 1}} $ 为$\gamma - 1$ 时刻探测机器人偏航角;$ \Delta {\theta _{\gamma - 1,\gamma }} $ 为$\gamma - 1$ 时刻到$\gamma $ 时刻探测机器人偏航角变化。状态变量的协方差预测方程为
$$ {{\mathop {\boldsymbol{P}}\limits^ \vee } _{{\gamma}} } = {{\boldsymbol{F}}_{{\gamma}} }{{\mathop {\boldsymbol{P}}\limits^ \wedge } _{{{\gamma}} - {{1}}}}{\boldsymbol{F}}_\gamma^{\rm{T}} + {\boldsymbol{R}} $$ (2) 式中:
$ {{\mathop {\boldsymbol{P}}\limits^ \vee} _{{\gamma}} } $ 和$ {{\mathop {\boldsymbol{P}}\limits^ \wedge } _{{{\gamma}} - {{ 1}}}} $ 分别为$\gamma $ 时刻状态变量的先验协方差矩阵和$\gamma - 1$ 时刻后验协方差矩阵;$ {{\boldsymbol{F}}_{{\gamma}} } $ 为$\gamma $ 时刻运动方程对状态变量求导得到的协方差矩阵;$ {\boldsymbol{R}} $ 为系统噪声协方差矩阵。$$ {{\boldsymbol{F}}_{{\gamma}} } = \left[ {\begin{array}{*{20}{c}} 1&0&{ - \Delta {d_{\gamma - 1,\gamma }}\sin ({\theta _{\gamma - 1}} + \Delta {\theta _{\gamma - 1,\gamma }}/2)} \\ 0&1&{\Delta {d_{\gamma - 1,\gamma }}\cos ({\theta _{\gamma - 1}} + \Delta {\theta _{\gamma - 1,\gamma }}/2)} \\ 0&0&1 \end{array}} \right] $$ (3) 将IMU数据积分得到的位置信息
$({X_{\gamma }},{Y_{\gamma }})$ 、偏航角信息${\theta _{\gamma }}$ 作为$\gamma $ 时刻观测值,则有如下观测方程:$$ {{\boldsymbol{z}}_{{\gamma}} } = \left[ {\begin{array}{*{20}{c}} {{X_{\gamma }}} \\ {{Y_{\gamma }}} \\ {{\theta _{\gamma }}} \end{array}} \right] + {{\boldsymbol{V}}_\gamma } $$ (4) 式中
${{\boldsymbol{V}}_\gamma }$ 为观测噪声。根据式(3)和式(4)进行状态更新,卡尔曼增益计算公式为
$$ {{\boldsymbol{\kappa}} _{{\gamma}} } = {{\mathop {\boldsymbol{P}}\limits^ \vee} _{{\gamma}} }{\kern 1pt} {{\boldsymbol{\varTheta}}^{\bf{T}}}{\left( {{\boldsymbol{\varTheta}}{\kern 1pt} {{\mathop {\boldsymbol{P}}\limits^ \vee }_{{\gamma}} }{{\boldsymbol{\varTheta}}^{\bf{T}}} + {{\boldsymbol{Q}}_{{\gamma}} }} \right)^{ - 1}} $$ (5) 式中:
$ {{\boldsymbol{\kappa}} _{\boldsymbol{\gamma}} } $ 为$\gamma $ 时刻计算得到的卡尔曼增益;${\boldsymbol{\varTheta}}$ 为观测方程对观测值求偏导得到的雅克比矩阵;$ {{\boldsymbol{Q}}_{{\gamma}} } $ 为观测噪声的协方差矩阵。根据式(5)可得状态变量的后验估计:
$$ {{\mathop {\boldsymbol{x}}\limits^ \wedge} _{\boldsymbol{\gamma}} } = {{\mathop {\boldsymbol{x}}\limits^ \vee } _{\boldsymbol{\gamma}} } + {{\boldsymbol{\kappa}} _{\boldsymbol{\gamma}} }\left[ {\begin{array}{*{20}{c}} {{X_{\gamma }} - {{\mathop {X}\limits^ \vee }_\gamma }} \\ {{Y_{\gamma }} - {{\mathop {Y}\limits^ \vee }_\gamma }} \\ {{\theta _{\gamma }} - {{\mathop {{\theta}} \limits^ \vee }_\gamma }} \end{array}} \right] $$ (6) 式中
${{\mathop {{X}}\limits^ \vee} _\gamma }$ ,${{\mathop {{Y}}\limits^ \vee } _\gamma }$ ,${{\mathop {{\theta}} \limits^ \vee} _\gamma }$ 分别为$\gamma $ 时刻状态变量的先验估计值。状态变量的后验协方差矩阵为
$$ {{\mathop {\boldsymbol{P}}\limits^ \wedge } _{{\gamma}} } = \left( {{\boldsymbol{I}} - {{\boldsymbol{\kappa}} _{{\gamma}} }{\boldsymbol{\varTheta }}} \right){{\mathop {\boldsymbol{P}}\limits^ \vee } _{{\gamma}} } $$ (7) 式中
$ {\boldsymbol{I}} $ 为单位矩阵。1.3 相对位姿关系
ArUco码到探测机器人坐标系需要经过多次坐标变换,假设以洞口位置探测机器人的车体中心为原点的三维世界坐标系下某点p的齐次坐标为
${{\boldsymbol{p}}_{{{{\rm{wq}}}}}} = $ $ {\left[ {{X_{\rm{w}}}\;\;{Y_{\rm{w}}}\;\;{Z_{\rm{w}}}\;1} \right]^{\rm{T}}}$ ,且其对应的相机坐标系下坐标为$ {{\boldsymbol{p}}_{{{\rm{c}}}}} = $ $ {[{X_{\text{c}}}\;\;\;{Y_{\text{c}}}\;\;\;{Z_{\text{c}}}]^{\rm{T}}} $ ,则根据针孔相机模型可得到$ {{\boldsymbol{p}}_{{{\rm{c}}}}} $ 到像素坐标系下坐标$ {\boldsymbol{p}} = {[{m}\;\;\;{n}]^{\rm{T}}} $ 的变换关系:$$ \left[ {\begin{array}{*{20}{c}} {{m_1}} \\ {{n_1}} \\ 1 \end{array}} \right] = \frac{1}{{{Z_{\rm{c}}}}}{\boldsymbol{K}}{{\boldsymbol{p}}_{{{\rm{c}}}}} = \frac{1}{{{Z_{\rm{c}}}}}\left[ {\begin{array}{*{20}{c}} {{f_{{m}}}}&0&{{c_{{m}}}} \\ 0&{{f_{{n}}}}&{{c_{{n}}}} \\ 0&0&1 \end{array}} \right]{\boldsymbol{A}}\left[ {\begin{array}{*{20}{c}} {{X_{\rm{w}}}} \\ {{Y_{\rm{w}}}} \\ {{Z_{\rm{w}}}} \\ 1 \end{array}} \right] $$ (8) 式中:
${\boldsymbol{K}}$ 为相机内参数矩阵,一般由标定得到;${f_{{m}}}$ ,${f_{{n}}}$ 分别为像素平面m轴和n轴上的尺度因子;${c_{{m}}}$ ,${c_{{n}}}$ 为相机光心偏移量;A为变换矩阵,$ {\boldsymbol{A}} = [{\boldsymbol{r}}\;\;\;{\boldsymbol{\tau}} ] $ ,$ {\boldsymbol{r}} $ ,${\boldsymbol{\tau}}$ 分别为相机当前位姿的旋转矩阵和平移矩阵。实际安装时相机并不处于探测机器人坐标系中心,故某点p在探测机器人坐标系下的坐标需进行坐标变换得到。
$$ {{\boldsymbol{p}}_{{{{\rm{rq}}}}}} = {}_{{{\rm{c}}}}^{{{\rm{r}}}}{\boldsymbol{T}}{{\boldsymbol{p}}_{{{{\rm{cq}}}}}} = {}_{{{\rm{c}}}}^{{{\rm{r}}}}{\boldsymbol{T}}{}_{{{\rm{w}}}}^{{{\rm{c}}}}{\boldsymbol{T}}{{\boldsymbol{p}}_{{{{\rm{wq}}}}}} = {}_{{{\rm{c}}}}^{{{\rm{r}}}}{\boldsymbol{T}}\left[ {\begin{array}{*{20}{c}} {\boldsymbol{A}} \\ {\begin{array}{*{20}{c}} {{{\bf{0}}^{{{\rm{T}}}}}}&{{1}} \end{array}} \end{array}} \right]{{\boldsymbol{p}}_{{{{\rm{wq}}}}}} $$ (9) 式中:
${{\boldsymbol{p}}_{{\rm{rq}}}}$ ,${{\boldsymbol{p}}_{{\rm{cq}}}}$ 分别为某点p在探测机器人坐标系、相机坐标系下的齐次坐标;${}_{\rm{c}}^{\rm{r}}{\boldsymbol{T}}$ 为探测机器人坐标系和相机坐标系之间的变换矩阵,由于相机固定在探测机器人车身上,所以其中参数均为常数;${}_{\rm{w}}^{\rm{c}}{\boldsymbol{T}}$ 为世界坐标系到相机坐标系的变换矩阵。根据式( 8)和式( 9)可由像素坐标推算得到探测机器人坐标系与世界坐标系下ArUco码之间的位姿变换关系。
1.4 ArUco码识别与解算
由于相机固定在探测机器人车身上,所以求得像素坐标系相对于ArUco码所在的世界坐标系的位姿变换矩阵及像素坐标系与相机坐标系之间的位姿变换矩阵,即可得到相机和ArUco码的位置变换关系。
首先使用张正友标定法对相机进行标定,得到相机的内参矩阵和畸变参数。然后分别通过自适应滤波对图像进行二值化,利用Suzuki算法进行轮廓提取并再次进行滤波,采用Douglas−Peucker算法进行四边形轮廓拟合并且得到4个角点。因ArUco码外侧边界均为黑色,可据此进行筛选得到ArUco码的候选轮廓,从而得到其角点在世界坐标系的位置。最后通过求解单应矩阵即可得到2帧图像间探测机器人的运动估计。
考虑角点位于ArUco所在平面,可设其为
$ {{\boldsymbol{C}}^{\rm{T}}}{{p}} + $ $ {{D}} = 0 $ ,其中${\boldsymbol{C}}$ 为系数矩阵,$ {{D}} $ 为常数。则可推导出2帧图像的像素坐标系中一对匹配的角点对$ {{\boldsymbol{p}}_{{1}}} = $ $ {[{m_1}\;\;\;{n_1}]^{\rm{T}}} $ 和$ {{\boldsymbol{p}}_{{2}}} = {[{m_2}\;\;\;{n_2}]^{\rm{T}}} $ 的转换关系:$$ s{{\boldsymbol{p}}_{{2}}} = {\boldsymbol{K}}({\boldsymbol{r}}{{\boldsymbol{p}}_{{1}}} + {\boldsymbol{\tau }}) = {\boldsymbol{K}}({\boldsymbol{r}} - \frac{{{\boldsymbol{\tau }}{{\boldsymbol{C}}^{\mathbf{T}}}}}{{{D}}}){{\boldsymbol{K}}^{ - 1}}{{\boldsymbol{p}}_{{1}}} = {{\boldsymbol{H}}_{{{3}} \times {{3}}}}{{\boldsymbol{p}}_{{1}}} $$ (10) 式中s为尺度因子。
根据匹配点参数即可求解出
${{\boldsymbol{H}}_{{{3}} \times {{3}}}}$ 矩阵:$$ s\left[ {\begin{array}{*{20}{c}} {{m_2}} \\ {{n_2}} \\ 1 \end{array}} \right] = {{\boldsymbol{H}}_{{{3}} \times {{3}}}}\left[ {\begin{array}{*{20}{c}} {{m_1}} \\ {{n_1}} \\ 1 \end{array}} \right] $$ (11) 展开式(11)得
$$ \left\{ {\begin{array}{*{20}{c}} {{\rm{s}}{m_2} = {m_1}{h_{11}} + {n_1}{h_{12}} + {h_{13}}} \\ {{\rm{s}}{n_2} = {m_1}{h_{21}} + {n_1}{h_{22}} + {h_{23}}} \\ {{\rm{s}}{{_{}^{}}_{}} = {m_1}{h_{31}} + {n_1}{h_{32}} + {h_{33}}} \end{array}} \right. $$ (12) 式中
${h_{\alpha \beta }}$ $(\alpha ,\beta=1,2,3) $ 为${{\boldsymbol{H}}_{{\mathbf{3}} \times {\mathbf{3}}}}$ 矩阵中第$\alpha $ 行第$\beta $ 列元素。进一步消去尺度因子
$ s $ ,且当$ {h_{33}} $ 非0时,设$ {h_{33}} = 1 $ ,将${{\boldsymbol{H}}_{{\mathbf{3}} \times {\mathbf{3}}}}$ 矩阵看作列向量进行化简:$$ \left[ {\begin{array}{*{20}{c}} {{m_1}}&{{n_1}}&1&0&0&0&{ - {m_1}{m_2}}&{ - {n_1}{m_2}}\\ 0&0&0&{{m_1}}&{{n_1}}&1&{ - {m_1}{n_2}}&{ - {n_1}{n_2}} \end{array}} \right] \left[ {\begin{array}{*{20}{c}} {{h_{11}}} \\ {{h_{12}}} \\ {{h_{13}}} \\ {\begin{array}{*{20}{c}} {{h_{21}}} \\ {{h_{22}}} \\ {{h_{23}}} \\ {\begin{array}{*{20}{c}} {{h_{31}}} \\ {{h_{32}}} \end{array}} \end{array}} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {{m_2}} \\ {{n_2}} \end{array}} \right] $$ (13) 式中有
$ {h_{11}} $ ,$ {h_{12}} $ ,$ {h_{13}} $ ,$ {h_{21}} $ ,$ {h_{22}} $ ,$ {h_{23}} $ ,$ {h_{31}} $ ,$ {h_{32}} $ 8个未知量,故利用4个角点的匹配点对即可建立方程组,完成相邻2帧图像间探测机器人位姿变换的初步计算。此外,根据提取的轮廓可进行二进制编码解算,由于每个ArUco码表示唯一,只有原始角点顺序下的ArUco码可在字典中找到编号,所以对提取的图像进行4次旋转,分别计算其编码值,即可得到一个唯一确定的编号。由于在首次完成ArUco码放置后,每个ArUco码在世界坐标系下的位姿均为固定值,根据其编号可确定ArUco码表示的位姿信息。1.5 姿态和位置信息优化
探测机器人携带单目相机进入隧洞,当检测到ArUco码时,按照1.3和1.4节所述方法计算得到探测机器人和ArUco码之间的相对位姿关系及ArUco码代表的绝对位置信息。考虑到探测机器人在运动过程中传感器存在误差,因此,使用图优化方式进行姿态和位置信息的优化。将全过程待优化的状态变量记为
$\chi $ ,则有$$ {\boldsymbol{\chi}} = [{{\boldsymbol{\eta}} _0}\;\;\;{{\boldsymbol{\eta}} _1}\;\cdots\;\;\;{{\boldsymbol{\eta}} _i}\;\;\;{\lambda _0}\;\;\;{\lambda _1}\;\cdots\;{\lambda _j}] $$ (14) 第k帧图像采集时刻的位姿为
$$ {{\boldsymbol{\eta}} _{{k}}} = [{w_k}\;\;\;{q_k}\;\;\;{v_k}\;\;\;b_k^{\rm{a}}\;\;\;b_k^{\rm{g}}] $$ (15) 式中:
$ {{\boldsymbol{\eta}} _i} $ 为第$i$ 个观测点处的本体位姿;${\lambda _j}$ 为观测到的第$j$ 个路标位置;${w_k}$ 为第k帧图像采集时刻的IMU位置;${q_k}$ 为第k帧图像采集时刻的IMU姿态;${v_k}$ 为第k帧图像采集时刻的IMU速度;$b_k^{\rm{a}}$ ,$b_k^{\rm{g}}$ 分别为陀螺仪和加速度计的偏移。在实际环境下,探测机器人在运动过程中存在噪声,基于高斯分布假设和最大似然估计,可建立如下目标函数:
$$ {J^*} = \arg \max \phi(z,u|{\boldsymbol{\eta}} ,\lambda ) $$ (16) 式中:
$ z,u $ 分别为观测数据和输入数据;$ {\boldsymbol{\eta}} ,\lambda $ 分别为待估计的位姿和路标位置。假设各个时刻的输入数据和观测数据之间相互独立,则联合分布可分解为
$$ \phi(z,u|{\boldsymbol{\eta}} ,\lambda ) = \prod\limits_k {\phi({u_k}|{{\boldsymbol{\eta}} _{k - 1}},{{\boldsymbol{\eta}} _k})\prod\limits_{k,j} {\phi({z_{k,j}}|{{\boldsymbol{\eta}} _k},{\lambda _j})} } $$ (17) 式中:
$ {u_k} $ 为第k帧的输入数据;$ {{\boldsymbol{\eta}} _{k - 1}} $ ,${{\boldsymbol{\eta}} _k}$ 分别为第$k-1 $ 帧和第k帧的位姿估计值;$ {z_{k,j}} $ 为第k帧探测机器人对第$j$ 个路标位置的观测值。设运动模型噪声
$ {{\boldsymbol{e}}_{u,k}} $ 服从正态分布$ N(0,{{\boldsymbol{\sigma}} _{u,k}}) $ ,观测模型噪声${{\boldsymbol{e}}_{z,k,j}}$ 服从正态分布$N(0,{{\boldsymbol{\sigma}} _{z,k,j}})$ ,$ {{\boldsymbol{\sigma}} _{u,k}},{{\boldsymbol{\sigma}} _{z,k,j}} $ 分别为2类噪声的协方差矩阵。那么单独考虑相机观测方程$ {z_{k,j}} = g({{\boldsymbol{\eta}} _k},{\lambda _j}) $ 并联立式(16)、式(17),同时对两侧取负对数并化简得$$ \begin{split} & {\rm{min}}[ - {\rm{ln}}\phi(z|{\boldsymbol{\eta}} ,\lambda )] = G + {\rm{min}}\left[ {\frac{1}{2}{{(z - g({\boldsymbol{\eta}} ,\lambda ))}^{\rm{T}}}{\boldsymbol{\sigma}} _{z,k,j}^{ - 1}(z - g({\boldsymbol{\eta}} ,\lambda ))} \right] =\hfill \\ &\quad\quad\quad\quad\quad\quad\quad\;\; G + \sum\limits_{{{k}} = 1}^{{L}} {\sum\limits_{{{j}} = 1}^{{E}} {\frac{1}{2}{\boldsymbol{e}}_{z,k,j}^{\rm{T}}} } {\boldsymbol{\sigma}} _{z,k,j}^{ - 1}{{\boldsymbol{e}}_{z,k,j}} \hfill \\[-18pt] \end{split} $$ (18) 式中:
$ G $ 为常数;L为进行优化的位姿观测点数;E为观测到的路标数。同理,单独考虑输入情况方程
$l({{\boldsymbol{\eta }}_{k - 1}},{{\boldsymbol{\eta}} _k})$ ,有$$ {J^*} = {\rm{argmin}}(\sum\limits_{{{k}} = 1}^{{L}} {\frac{1}{2}{\boldsymbol{e}}_{u,k}^{\rm{T}}{\boldsymbol{\sigma}} _{u,k}^{ - 1}{{\boldsymbol{e}}_{u,k}}} + \sum\limits_{{{k}} = 1}^{{L}} {\sum\limits_{{{j}} = 1}^{{E}} {\frac{1}{2}{\boldsymbol{e}}_{z,k,j}^{\rm{T}}} } {\boldsymbol{\sigma}} _{z,k,j}^{ - 1}{{\boldsymbol{e}}_{z,k,j}}) $$ (19) 待求解优化函数可表示为
$$ \begin{split} {J^*} =& \mathop {\arg \min }\limits_{\boldsymbol{\chi}} \Bigg\{ \sum\limits_{k \in B} {\mu \rho (||{M_{\rm{IMU}}}(\tilde Z_{{\rm{IMU}}k + 1}^{{\rm{IMU}}k},{\boldsymbol{\chi}} )||^2)} + \Bigg.\\ & \qquad\quad \Bigg. \sum\limits_{(i,j) \in W} {(1 - \mu )\rho (||{M_{\rm{f}}}(\tilde Z_{{f_j}}^{{c_i}},{\boldsymbol{\chi}} )||^2)} \Bigg\} = \hfill \\ &\qquad\quad \mathop {\arg \min }\limits_{\boldsymbol{\chi}} \left\{ {\mu {U_{{\rm{IMU}}}} + (1 - \mu ){U_{{\rm{cam}}}}} \right\} \hfill \end{split} $$ (20) 式中:
$B$ 为第k帧图像采集时刻IMU测量数据集合;$\mu $ 为比例因子,2帧图像之间探测机器人在ArUco码范围内识别的特征点信息越多且运动速度越慢,则$\mu $ 越小,表明该情况下依靠相机的位姿估计准确性更高;$\rho $ 为huber核函数,当匹配偏差超过阈值使其保持线性增长,从而降低错误匹配带来的影响;${M_{\rm{IMU}}}(\tilde Z_{{\rm{IMU}}k + 1}^{{\rm{IMU}}k},{\boldsymbol{\chi}} )$ 为IMU数据的残差项($\tilde Z_{{\rm{IMU}}k + 1}^{{\rm{IMU}}k} $ 为IMU的观测值),由IMU测量预积分形成的观测值和实际值比较形成,用于约束k和k+1连续2帧;$W$ 为匹配良好的特征点对集合;${M_{\rm{f}}}(\tilde Z_{{f_j}}^{{c_i}},{\boldsymbol{\chi}} )$ 为相机残差项($\tilde Z_{{f_j}}^{{c_i}} $ 为相机观测值),通过将特征点$ {c_i} $ 投影到第${f_j}$ 帧并与测量值计算重投影残差得到;${U_{{\rm{IMU}}}}$ 和${U_{{\rm{cam}}}}$ 分别为IMU和相机约束。利用里程计信息、IMU信息和相机信息,根据图优化的方式建立如图2所示的约束因子,图中
${t_k}$ 为时间戳对齐至第k帧的时间,v'k为里程计测量速度。在ArUco码附近范围实现探测机器人位姿优化,其中里程计约束提供的是相对约束,转换到IMU坐标系中进行速度和位姿的限制。结合绝对位置信息和局部位姿优化结果更新探测机器人位姿,同时完成对IMU和里程计累计误差的修正。
2. 隧洞探测机器人系统设计
为验证多传感器融合定位算法在尾矿库排洪隧洞实际环境中的运行效果,研究并设计了隧洞探测机器人系统,如图3所示。
尾矿库排洪隧洞探测机器人系统以履带式探测机器人实验平台为本体,基于机器人操作系统(Robot Operating System,ROS)进行软件设计,采用分布式系统,将监视端工控机布置于隧洞入口,而运行多传感器融合定位算法的工控机安装在探测机器人本体上,上下位机之间利用无线网桥进行组网通信。为实现多传感器融合定位算法,隧洞探测机器人实验平台和洞口的监控平台需要进行相应的软硬件设计,且根据系统功能将隧洞探测机器人系统划分为4个子系统:驱动控制系统、主控系统、传感通信系统、人机交互系统。驱动控制系统对接主控系统,进行电动机控制和电动机相关数据传输,保证探测机器人的运动能力。主控系统对接驱动控制系统和传感通信系统,运行多传感器融合定位、定位校正等核心算法,实现探测机器人自主定位和定位校正等功能,并解析执行控制端操作人员指令,通过串口转RS232实现与驱动控制系统的通信,通过USB获取传感器数据,在进行数据分析后将关键数据和音视频通过传感通信系统传输至人机交互系统,并实现人机交互系统的各类指令下传,以实现对探测机器人本体的控制。传感通信系统对接主控系统和人机交互系统,负责探测机器人环境感知和通信,使洞口操作人员可实时查看探测机器人实验平台的运行状态。人机交互系统对接传感通信系统,负责解析操作人员指令并通过传感通信系统将指令发送到主控系统,将主控系统数据实时显示在监视端工控机界面中。
3. 实验验证和分析
为探究搭载多传感器融合定位算法的尾矿库排洪隧洞探测机器人实验平台在实际环境中运行的稳定性,并验证多传感器融合定位算法的有效性,在尾矿库排洪隧洞实际场景中运行和关闭多传感器融合定位算法,分别进行10组行进40 m的重复定位实验,且每轮实验均记录探测机器人在20 m和40 m处的定位误差。运行多传感器融合定位算法需要在20 m处额外记录校正后误差。
在实验开始前需要标定和校正系统误差,探测机器人平均车轮直径为0.21 m,理论轮距为 0.88 m,使用UMBmark算法进行标定,根据实验数据求得校正后左轮直径为0.209 5 m,右轮直径为0.210 5 m,轮间距为0.84 m。
修改里程计参数后,在距起点20 m处的排洪隧洞墙壁上布设ArUco码。在探测机器人系统中利用ROS编写程序,对校正后的里程计和IMU数据进行EKF,利用OpenCV库进行ArUco码识别与相应图像处理操作,并利用Graph SLAM Tutorial实现图优化算法对探测机器人位姿进行优化并校正定位。在尾矿库排洪隧洞中进行重复定位实验,实验环境和实验平台实物如图4所示。
根据实验设计运行和关闭多传感器融合定位算法,在相同环境中各进行10组40 m行进实验,并记录实验数据,各组实验20 m和40 m处定位误差如图5所示。根据图5数据统计定位误差,结果见表1。
表 1 探测机器人定位误差Table 1. Positioning error of detection robot定位位置 平均定位误差/cm 关闭算法实验组 运行算法实验组 行进20 m 22.51 19.77 行进40 m 58.11 21.23 运行多传感器融合定位算法实验组利用行进20 m时的ArUco码进行误差校正。使用米尺测量ArUco码中心和探测机器人车身中心的
$X$ ,$Y$ 方向数据,与位姿估计结果比较并计算距离误差,得误差均值为4.2 mm。从图5、表1可看出,使用多传感器融合定位算法的探测机器人在行进20 m时平均定位误差为19.77 cm,在行进40 m时的平均定位误差为21.23 cm。同样条件下,关闭多传感器融合定位算法实验组在行进20 m时的平均定位误差为22.51 cm,在行进40 m时的平均定位误差为58.11 cm。可见使用多传感器融合定位算法的探测机器人在运行过程中可保持更高定位精度,且在20 m的ArUco码放置处多传感器融合定位算法可有效修正累计误差,使校正后误差均值为4.2 mm。
4. 结论
(1) 针对现有定位算法不适用于网络信号弱、场景重复率高且属于狭长室内环境的尾矿库排洪隧洞探测机器人定位的问题,提出了一种探测机器人多传感器融合定位算法。将排洪隧洞环境的远距离定位以ArUco码为边界点,简化为多段短距离定位。未检测到ArUco码时,利用EKF算法融合里程计和IMU信息,可实现位姿信息的快速更新;检测到ArUco码时,则利用图优化算法将里程计、IMU、相机信息组成约束因子图,并结合ArUco码的绝对位置实现探测机器人局部位姿优化和定位校正,避免误差累计,最终实现狭长空间内探测机器人的远距离高精度定位。
(2) 为验证多传感器融合定位算法的有效性,设计了隧洞探测机器人系统,在尾矿库排洪隧洞实际场景中运行和关闭多传感器融合定位算法,分别进行10组行进40 m的重复定位实验,实验结果表明:多传感器融合定位算法具有较高的稳定性和精确性,能够实现隧洞探测机器人在排洪隧洞环境中的高精度定位,行进20 m的平均定位误差为19.77 cm,40 m的平均定位误差为21.23 cm,且20 m处校正后误差均值为4.2 mm。
致谢:本文得到了北京科技大学、福建马坑矿业股份有限公司、中国安全生产科学研究院和中国计量大学人工智能与机器人团队的支持,对此表示感谢!
-
表 1 探测机器人定位误差
Table 1 Positioning error of detection robot
定位位置 平均定位误差/cm 关闭算法实验组 运行算法实验组 行进20 m 22.51 19.77 行进40 m 58.11 21.23 -
[1] 张彪. 尾矿库在线监测预警系统的设计与实现[D]. 北京: 中国地质大学(北京), 2017. ZHANG Biao. Design and implementation of tailings pond online safety monitoring and pre-warning system[D]. Beijing: China University of Geosciences, 2017.
[2] KOMLJENOVIC D,STOJANOVIC L,MALBAI V,et al. A resilience-based approach in managing the closure and abandonment of large mine tailing ponds[J]. International Journal of Mining Science and Technology,2020,30(5):737-746. DOI: 10.1016/j.ijmst.2020.05.007
[3] 国家安全生产监督管理总局. 尾矿库安全监督管理规定[2011]38号[Z]. 北京: 国家安全生产监督管理总局, 2011. State Administration of Work Safety. Regulations on safety supervision and management of tailings pond [2011]No.38[Z]. Beijing: State Administration of Work Safety, 2011.
[4] GB 39496—2020 尾矿库安全规程[S]. GB 39496-2020 Safety regulation for tailings pond[S].
[5] JUNG C,MOON C B,JUNG D,et al. Design of test track for accurate calibration of two wheel differential mobile robots[J]. International Journal of Precision Engineering & Manufacturing,2014,15(1):53-61.
[6] 张胜宾,赵祚喜. 基于UMBmark算法的移动机器人定位试验研究[J]. 现代电子技术,2018,41(3):80-83. ZHANG Shengbin,ZHAO Zuoxi. Research on mobile robot localization experiment based on UMBmark algorithm[J]. Modern Electronics Technique,2018,41(3):80-83.
[7] DANIEL D,JAN K,KAREL H,et al. Odometer module for mobile robot with position error estimation[J]. IFAC Papers OnLine,2016,49(25):346-351. DOI: 10.1016/j.ifacol.2016.12.063
[8] 蔡李花,方海峰,李允旺,等. 矿井救灾机器人自主定位方法研究[J]. 工矿自动化,2015,41( 7):62-67. CAI Lihua,FANG Haifeng,LI Yunwang,et al. Research of self-localization method of mine rescue robot[J]. Industry and Mine Automation,2015,41( 7):62-67.
[9] ZHANG Chao, ZHAN Quanzhong, WANG Qi, et al. Autonomous dam surveillance robot system based on multi-sensor fusion[J]. Sensors, 2020, 20( 4): 1097. DOI: 10.3390/s20041097.
[10] 杨林,马宏伟,王岩,等. 煤矿巡检机器人同步定位与地图构建方法研究[J]. 工矿自动化,2019,45( 9):18-24. YANG Lin,MA Hongwei,WANG Yan,et al. Research on method of simultaneous localization and mapping of coal mine inspection robot[J]. Industry and Mine Automation,2019,45( 9):18-24.
[11] WEN Jingren, QIAN Chuang, TANG Jian, et al. 2D LiDAR SLAM back-end optimization with control network constraint for mobile mapping[J]. Sensors, 2018, 18(11): 3668. DOI: 10.3390/s18113668.
[12] HESS W, KOHLER D, RAPP H, et al. Real-time loop closure in 2D LIDAR SLAM[C]// 2016 IEEE International Conference on Robotics and Automation (ICRA), 2016.
[13] 李猛钢. 煤矿救援机器人导航系统研究[D]. 徐州: 中国矿业大学, 2017. LI Menggang. Research on navigation system of coal mine rescue robot [D]. Xuzhou: China University of Mining and Technology, 2017.
[14] 邢伯阳,潘峰,冯肖雪. ArUco−SLAM: 基于ArUco二维码阵列的单目实时建图定位系统[J]. 北京理工大学学报,2020,40(4):427-433. XING Boyang,PAN Feng,FENG Xiaoxue. ArUco-SLAM: a monocular SLAM system based on ArUco landmark array[J]. Transactions of Beijing Institute of Technology,2020,40(4):427-433.
[15] 肖献强,程亚兵,王家恩. 基于惯性和视觉复合导航的自动导引小车研究与设计[J]. 中国机械工程,2019,30(22):2734-2740. DOI: 10.3969/j.issn.1004-132X.2019.22.013 XIAO Xianqiang,CHENG Yabing,WANG Jia'en. Research and design of AGVs based on inertial and visual composite navigation[J]. China Mechanical Engineering,2019,30(22):2734-2740. DOI: 10.3969/j.issn.1004-132X.2019.22.013
[16] KALLWIES J, FORKEL B, WUENSCHE H J. Determining and improving the localization accuracy of apriltag detection[C]// 2020 IEEE International Conference on Robotics and Automation (ICRA), 2020.
[17] ZHANG Z. A flexible new technique for camera calibration[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence,2000,22(11):1330-1334. DOI: 10.1109/34.888718