刘 玉,陈 凤,黄建明,魏祥泉
(上海宇航系统工程研究所,上海201108)
空间站的组装建造、在轨检测、维护维修、辅助航天员出舱活动、支持空间应用[1]等都需要借助空间机械臂,完成多种复杂的舱外操作。机械臂视觉系统作为空间机械臂遥操作控制的伺服输入,会直接影响空间机械臂的控制精度,进而进一步影响机械臂的空间操作能力。
空间机械臂视觉测量技术主要包括手眼关系标定、标志器识别、相对三维位姿测量等关键技术[2-6]。目前国外已经在空间站使用的标志器各有优缺点,结合我国空间站实际任务需求,选择设计的标志器上应有容易识别的点、线以及圆等丰富的特征信息,以为后续机械臂视觉高精度测量奠定基础。本文重点对标志器识别、相对三维位姿测量等关键技术进行深入研究,并以加拿大机械臂为例提出了一种基于边缘特征的标志器识别方法和基于非迭代的位姿测量算法,并给出了仿真实验结果以及机械臂原理样机的实验结果。
空间机械臂视觉测量技术的主要工作流程如图1所示。
图1 空间机械臂视觉测量技术流程框图Fig.1 The flow chart of the vision measurement technology of Space Robotic Arm
首先,通过手眼关系标定技术确定空间机械臂末端与摄像机之间的关系;其次,通过手眼相机采集标志器图像,利用图像处理与模式识别等技术完成标志器的识别,从中提取标志器特征信息;然后,利用相对三维位姿测量技术获得目标标志器相对摄像机的位置与姿态,并通过手眼关系矩阵确定目标标志器相对机械臂末端的位置姿态;最后,将解算出的标志器相对机械臂末端的位置姿态信息作为空间机械臂视觉伺服控制的输入,反馈给机械臂伺服控制器。伺服控制器则可根据实际需要发送控制指令给视觉测量系统执行相应的功能。
假定摄像机内/外部参数已标定,那么标定手眼关系则主要通过设计合理的机械臂运动,计算出手眼关系矩阵。如图2所示,A为机械臂末端坐标系在这两个相应位置间的相对变换矩阵,B则为摄像机坐标系在这两个不同位置间的相对变换矩阵,则有:AX=XB。机械臂的手眼关系标定即转换为求解上述方程[7],具体求解方法可详见文献[8]。
图2 手眼关系标定原理图Fig.2 The principle of calibration of hand-eye relationship
2.3.1 标志器设计
为使空间机械臂接近抓捕目标并准确地对其进行定位,需要在抓捕目标上安装标志器作为特征,并且标志器的设计应以与周围物体显著不同为原则。参考国外有关技术文献,标志器设计一般有两种:第一种是被动标志,第二种是主动标志。主动标志主动地辐射光能,被动标志包括特殊设计的观测靶标和角反射器[2-6]。目前用于机械臂视觉测量的标志均采用被动标志,例如加拿大机械臂的十字靶标[2](图3)、日本实验舱暴露平台上的人工标志[4-5]等。由于被动标志中角反射器标志易受到杂光干扰,会增加目标识别算法的复杂性。因此,本文主要采用基于特殊设计的标志器,考虑到加拿大机械臂中设计的十字标志器上有容易识别的点、线以及圆等丰富的特征信息[2],基本能够满足机械臂精细化操作要求,因此,本文主要参考该标志器进行设计。
图4所示为本文设计的标志器,其上有14个有效特征点。摄像机坐标系定义和标志器坐标系定义一致。标志器坐标系的原点确定为白色圆环的圆心;x轴正方向为原点指向靶杆的负方向,y轴为较长白色线条之一的对称线,z轴正方向为较短白色线条的对称线的反方向,坐标系满足右手法则。由于标志器上的直线特征比较明显,根据较长白色线条和较短白色线条具体尺寸和位置关系,可以确定出标志器区域。
图3 加拿大十字靶标Fig.3 The cross marker of the Canadian Robotic Arm
图4 标志器上14个特征点Fig.4 The fourteen feature points in the maker
2.3.2 标志器识别
由于标志器上有明显的边缘特征,且存在互相平行的直线。因此,为充分识别标志器,关键是识别到正确的直线边缘。图5给出了本文的具体识别流程。其中,图像预处理主要采用高斯滤波对图像进行去噪处理,并采用自适应阈值法[9]对图像进行二值化处理。下面重点对边缘特征提取、边缘特征匹配以及特征点提取作说明。
1)边缘特征提取
在二值化图像中,首先,采用Freeman链码法[10]获取以点序列形式表示的边缘特征信息。然后,为避免程序处理过多的非标志物区域,对检测到的区域进行筛选,提取每个区域的面积,凸凹性,最大包围盒的长宽比信息并进行判断,初步筛选出可能的标志器区域。最后,为进一步逼近边缘信息,用更少的边缘点信息表征边缘信息,采用Douglas-Peucker算法[11]进行简化处理,并以矢量化形式表征一个边缘信息,矢量化信息主要包括边缘线段的两个端点,边缘线段所在的轮廓等。
2)边缘特征匹配
图5 标志器识别流程Fig.5 The flowchart of maker identification
根据提取到的边缘特征,需要依据标志器CAD模型提供的几何形状、尺寸信息以及直线段间的几何位置关系对候选的边缘特征进行匹配筛选,以得到标志器特征。首先,对边缘特征中的直线段特征是否平行进行判断,同时辅以直线段特征长度范围、大致比例判断等。这里采用方程xcosθ+ysinθ= ρ( ρ,θ为常数)来表达一根直线。两直线平行,这意味着两直线方程的θ参数差值小于指定阈值。通过该方法可得到平行线段组的集合。其次,在平行线段组合中,根据直线段特征是否共线、垂直,同时辅以直线段特征是否属于同一个轮廓、长度范围、大致比例等,寻找标志器中四条较长白色线条和两条较短白色线条。若两直线共线,则表明两直线方程的ρ,θ参数的差分别小于相应的指定阈值;两直线垂直,则表明两直线方程的θ参数的差值的绝对值小于指定阈值。
3)特征点提取
根据识别出的标志器可进一步得到如图4所示的P0到P13特征点。编号为P1~P12的特征点都是边缘直线段间的交点,主要通过对边缘直线分别进行拟合,然后计算两条直线的交点精确定位。而特征点P0则直接根据提取到四条较长白色直线段间的中心线和两条较短直线的中心线的交点得到。特征点P13则是标志器中间圆杆顶部圆的中心点,在图像中则表现为椭圆中心点,一般通过在特定区域中采用最小二乘法拟合椭圆[12]得到。
根据标志器识别出的14个特征点信息,利用三维位姿解算方法可得到标志器相对摄像机的位姿。目前,基于单目视觉物体三维位姿求解算法有很多,大致可分为线性算法[13]和迭代算法[14]。线性算法可直接得到解析解,计算速度很快,但该算法对噪声较为敏感。常用的迭代算法则精度较高,但需要提供较好的初始值,否则很难收敛到正确解。文献[15]提出一种RPNP算法,该方法是一种非迭代算法,通过求解一元七次方程获得最小平方误差意义上的最优解。和传统的非迭代法相比,该方法具有计算精度高,且能够得到相对稳定的位姿解的优势。因此,本文主要采用该方法实现相对三维位姿解算。RPNP算法的中心思想是将参考点划分为3点子集而得到一个四次方程组,用该四次方程组的平方和构造代价函数,然后通过求解代价函数的导数确定最优解。该算法能够稳定地处理平面情况、一般3D情况与准奇异情况,达到迭代算法精度水平,却只占用很小的计算资源;能够在较少冗余点(n≤5)的情况下得到比迭代算法更高精度的解;算法效率高,能够高效处理大量参考点。RPNP算法的流程如图6所示,总体可分解为以下四个步骤[15]:
图6 RPNP算法的流程图Fig.6 The flowchart of RPNP algorithm
1)选择一个旋转轴
假定,摄像机坐标系为OcXcYcZc,标志器坐标系为OoXoYoZo,如图7所示。一般选择边集中投影长度最大边‖PiPj‖对应的PiPj作为旋转轴Za轴。为确定Za轴方向,端点Pi0和Pj0深度值必须知道。将n个参考点划分为n-2个3点子集,即{Pi0Pj0Pk|k≠ i0,k≠ j0}。通过 P3P约束,每个子集可得到一个四次多项式,如式(1)所示所示。
其中,未知参数x为Pi0的未知深度的平方。
图7 三维参考点的二维投影Fig.7 The projection of the reference points
2)确定最小二乘误差意义上的旋转轴
通过最小二乘误差方法[16]可解算上述方程组。首先,定义一个代价方程为其极小值可以通过求解它导数 F'=得到。F'为一个一元七次方程,可通过特征值方法求解。当x与端点Pi0的深度被确定,另一个端点Pj0的深度可以通过P2P约束求得。最后可以得到坐标系OaXaYaZa的旋转轴 Za=Pi0Pj0/‖Pi0Pj0‖ 。
3)求解旋转角与平移向量
当OaXaYaZa的Za轴被确定,OaXaYaZa到摄像头坐标系OcXcYcZc的旋转矩阵可以表达如式(2)。
从三维参考点到二维归一化图像平面的投影可以表达为式(3)。
其中,(ui,vi)是图像点 pi的归一化坐标,t=是平移向量。
我们对上式各项组成2n×6齐次线性方程组即可解算出未知变量向量[c s txtytz1]T
4)计算标志器相对摄像机的位置姿态
在实际使用中,由于噪声干扰,上述齐次线性方程组的解可能不满足三角函数约束。因此,须对旋转矩阵R施加正交约束。为对旋转矩阵R进行归一化处理[17],首先使用未归一化的R和t估计参考点在摄像机坐标系中的三维坐标,然后采用标准的三维对齐方法求解旋转与平移矩阵。代价函数F为多项式的平方和,至多包含4个局部最小值。根据每个局部最小值估计标志器相对摄像机的位置,并选择反向投影误差最小的解作为算法输出的最优解。
本文算法结果均在CPU为Intel(R)Core(TM)2 Duo CPU,2.99 GHz,内存为 1.98 GB 的PC机上运行所得。操作系统为Windows XP,开发平台为VS2008环境。
结合空间机械臂需要完成的空间操控任务,提出在0.5~1.5 m之间,视觉测量的位置精度≤0.005R(R表示机械臂末端与标志器之间距离),姿态精度≤0.5°,算法处理周期≤250 ms。为充分验证本文提出的空间机械臂视觉测量技术方法的可行性,本文分别对仿真图像和机械臂原理样机采集图像进行验证分析。
图像的数字仿真主要是建立包括摄像机模型、标志器模型在内的视觉成像仿真场景。这里主要采用Creator进行机械臂末端效应器与目标适配器即标志器的三维建模,并对模型表面材料属性及纹理进行渲染;然后根据机械臂末端效应器与标志器之间相对运动关系,利用Vega视景仿真技术,结合相机参数,生成一组相对运动的图像序列,其中图像分辨率为1292*964。
图8为标志器识别过程,其中(a)为原图像,(b)为边缘初步检测结果,(c)为边缘特征提取结果,(d)为根据标志器上边缘直线之间平行关系得到的边缘特征匹配结果,显然还有许多干扰的平行直线对,(e)为根据标志器上直线对之间的关系进一步筛选得到的结果,(f)为标志器上特征点提取结果。通过这一组实验可知,标志器识别算法能够有效识别出标志器。
图8 标志器识别过程Fig.8 The processof maker identification
图9 给出了摄像机与标志器距离650.92 mm处的仿真图像边缘检测及特征点提取结果,其中(b)中红色加号表示提取出的特征点。图10则给出了不同距离下根据特征点信息解算的位姿结果误差曲线,由曲线可知,根据仿真图像解算出的位姿测量结果误差在技术指标要求的范围内。由于本文主要采用基于非迭代思想的RPNP算法[15]实现位置姿态解算,因此,算法耗时量会明显下降,通过实验统计发现,本文算法对仿真图像的处理速度能够达到15 Fps,即算法处理周期在67 ms左右,显然满足预期速度指标要求。
图9 摄像机与标志器相距650.92 mm处检测结果Fig.9 The result of maker identificationat the 650.92 mm distance
图10 基于仿真图像的相对位置姿态误差曲线Fig.10 Theerror curve of the position and pose based on the simulated images
为进一步验证提出的空间机械臂视觉测量技术有效性,又开展了其在机械臂原理样机上的集成测试与验证工作。测试过程中移动机械臂末端到不同位置,利用安装于机械臂末端效应器的手眼相机进行实时图像采集,实时识别标志器,并将相关特征点高亮动态显示(图11中红色加号表示提取出的特征点),以判断标志器是否准确识别,识别之后进行三维位姿的实时测量和实时显示。
图11 摄像机相对标志器距离770.98 mm处检测结果Fig.11 The result of maker identificationat the 770.98 mm distance
图12 基于机械臂原理样机的相对位置姿态测量误差曲线Fig.12 Theerror curve of the position and pose based on the prototype of the robotic arm
为对测量精度进行验证,在某一位置处利用激光跟踪仪测量出机械臂末端与标志器之间的精确位姿,并与视觉测量结果进行比较,结果见图12。其中,图12(a)中的红色‘-x-’线为位置理想误差范围,(b)中的红色直线为姿态的理想误差范围。由图12可以看出,随着摄像机与标志器之间相对距离越来越小,测量误差呈减小趋势,并且位置测量误差没有超出理想误差范围,姿态测量误差在相对距离1 m以内时也未超出理想误差范围,在相对距离1米以外,只有个别位置处姿态误差较大。从总体趋势来看,视觉测量算法精度基本满足技术指标要求。图13展示了摄像机在移动过程中识别的标志器效果。从中可看出,标志器识别相对准确,进一步说明本文提出算法的鲁棒性。此外,实验统计了针对实际图像的算法耗时量,由于实际环境复杂,易受外界干扰影响,目前能达到平均每秒处理10帧图像,即算法耗时量在100 ms左右,比仿真图像的处理速度要慢些,但仍能满足了预期速度指标要求。
图13 摄像机与标志器在不同距离下识别结果Fig.13 The result of maker identification at the different distances
本文主要对空间机械臂视觉测量技术进行了深入研究,梳理了空间机械臂视觉测量技术的技术流程,介绍了手眼关系标定技术、标志器识别技术以及相对三维位姿测量技术等关键技术,并以加拿大机械臂为标志器,提出了一种基于边缘特征的标志器识别方法以及基于非迭代的三维位姿解算方法。最后,根据动力学仿真图像以及机械臂原理样机对本文提出的方法进行实验验证。实验结果表明:本文提出的空间机械臂视觉测量方法合理可行,具有较强的工程应用价值。
[1]罗亚中,林鲲鹏,唐国金.空间站运营任务规划技术评述[J].载人航天,2012,18(2):7-13.
[2]Beland S,潘科炎.加拿大的空间机器人——从国际空间站的灵敏作业机器人到行星探测机器人[J].控制工程(北京),2001(2):22-29.
[3]Ravindran R,Doetsch K.Design aspect of the shuttle remote manipulator system[C]//Proc.of the AIAA Guidance and control conference,1982:456-465.
[4]Mack B.,McClure S.,Ravindran R.A ground testbed for evaluating concepts for the special purpose dexterous manipulator[C]//IEEE International Conference on Robotics and Automation,1991:884-889.
[5]Moe R.An Early In-Space Platform for Technology Development& Demonstration of Robotic Assembly & Servicing[C]//AIAA 1st Space Exploration Conference:Continuing the Voyage of Discovery,2005.
[6]Inaba N,Oda M,Asano M.Rescuing a stranded satellite in space-Experimental robotic capture of non-cooperative satellites[J].Transactions of the Japan Society for Aeronautical and Space Sciences.2006,48(162):213-220.
[7]张广军.视觉测量[M].北京:科学出版社,2007:125-132.
[8]Shiu Y C,AhmadS.Calibration of wrist-mounted robotic sensors by solving homogenous transform equations of the form AX=XB[J].IEEE Trans.Robotics and Automation,1989,5(1):16-29.
[9]Wu H M,Zheng X H.A New Thresholding Method Applied to Motion Detection[C]//Wuhan:Pacific-Asia Workshop.2008:119-122.
[10]裘镇宇,危辉.基于Freeman链码的边缘跟踪算法及直线段检测[J].微型电脑应用.2008,24(1):17-20.
[11]Douglas D H,Peucker T K.Algorithms for the reduction of the number of points required to represent a line or its caricature[J].The Canadian Cartographer,1973,10(2):112-122.
[12]陈海峰,雷华,孔燕波,等.基于最小二乘法的改进的随机椭圆检测算法[J].浙江大学学报(工学版):2008,42(8):1360-1364.
[13]张世杰,曹喜滨,李晖.交会对接航天器间相对位姿参数单目视觉测量的解析算法[J].光学技术.2005,31(1):6-10.
[14]Fiore P D.Efficient Linear Solution of Exterior Orientation[J].IEEE Trans Pattern Analysis and Machine Intelligence,2001,23(2):140-148.
[15]Li Shiqi,Chi Xu,Ming Xie.A Robust O(n)Solution to the Perspective-n-Point Problem[J].IEEE Transactions on Pattern Analysis and Machine Intelligence,2012,34(7):1444-1450.
[16]Quan L,Lan Z.Quan and Z.Lan,Linear n-Point Camera Pose Determination[J].IEEETrans.Pattern Analysis and Machine Intelligence,1999,21(8):774-780.
[17]Umeyama S.Least-squares estimation of transformation parameters between two point patterns[J].IEEE Transactions onPattern Analysis and Machine Intelligence,1991,13(4):376-380.