CN112066879A - 基于计算机视觉的气浮运动模拟器位姿测量装置及方法 - Google Patents

基于计算机视觉的气浮运动模拟器位姿测量装置及方法 Download PDF

Info

Publication number
CN112066879A
CN112066879A CN202010952951.4A CN202010952951A CN112066879A CN 112066879 A CN112066879 A CN 112066879A CN 202010952951 A CN202010952951 A CN 202010952951A CN 112066879 A CN112066879 A CN 112066879A
Authority
CN
China
Prior art keywords
coordinate system
image
target
camera
points
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202010952951.4A
Other languages
English (en)
Other versions
CN112066879B (zh
Inventor
孙浩然
李莉
夏红伟
马广程
温奇咏
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Harbin Institute of Technology
Original Assignee
Harbin Institute of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Harbin Institute of Technology filed Critical Harbin Institute of Technology
Priority to CN202010952951.4A priority Critical patent/CN112066879B/zh
Publication of CN112066879A publication Critical patent/CN112066879A/zh
Application granted granted Critical
Publication of CN112066879B publication Critical patent/CN112066879B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
    • G01B11/00—Measuring arrangements characterised by the use of optical techniques
    • G01B11/002—Measuring arrangements characterised by the use of optical techniques for measuring two or more coordinates
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C21/00—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/70—Determining position or orientation of objects or cameras
    • G06T7/73—Determining position or orientation of objects or cameras using feature-based methods

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Automation & Control Theory (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Theoretical Computer Science (AREA)
  • Image Analysis (AREA)

Abstract

本发明提供了一种基于计算机视觉的气浮运动模拟器位姿测量装置及方法，属于仿真测试领域技术领域。本发明中针对单轴气浮平台，将其视为合作待测目标，首先在单轴气浮台上安装设计好的测量靶标，之后在测量图像中提取测量靶标的三维坐标，最后根据测量靶标的三维坐标反推单轴气浮平台的位置和姿态信息；针对三轴气浮台时，将其视为非合作目标，首先在三轴气浮台上预先提取一系列特征点，之后在测量图像中提取特征点的三维坐标，最后根据特征点的三维坐标重构三轴气浮台台体坐标系，并求解三轴气浮台台体坐标系与世界坐标系之间的坐标变换关系，进而得到三轴气浮台的位置和姿态信息。本发明具有实现简单、非接触、低成本、高精度的特点。

Description

基于计算机视觉的气浮运动模拟器位姿测量装置及方法
技术领域
本发明涉及一种基于计算机视觉的气浮运动模拟器位姿测量装置及方法，属于仿真测试领域技术领域。
背景技术
论文“大目标高精度位置姿态测量技术研究”(中国科学院大学硕士论文，崔艺涵，2018年6月)中提出了一种基于视觉引导的大目标高精度位置姿态测量方法。在摄像机附近固连激光测距仪，建立基于单目摄像机—激光测距仪的位姿测量系统。在待定位目标上设置合适的合作目标，利用单目摄像机获取合作目标中心位置，引导激光测距仪/激光经纬仪自动瞄准合作目标中心位置，进而通过参数解算得到待定位目标的空间位姿，实现高精度自动引导测量。但是其方法虽然能够有效地对目标进行高精度的位姿测量，但是采用激光跟踪仪价格昂贵，成本较高，因而具有一定的局限性。
论文“卫星高精度姿态确定算法及仿真技术研究”(哈尔滨工业大学硕士论文，刁雁华，2015年6月)中提出了一种基于高精度光学陀螺和两台光电自准直仪的位姿测量方法。在小角度时采用双矢量定姿的方法利用两台光电自准直仪测量数据得到位姿，在大角度时采用光学陀螺独立定姿的方法确定位姿。该论文在理论研究的基础上，完成联合定姿、EKF/UKF/PF/CKF定姿等多项试验，验证了该方法的正确性。但是其中提出的陀螺和光电自准直仪联合定姿方法虽然能够准确测量待测目标的位姿，但高精度陀螺仪具有抗过载能力差、体积大、价格昂贵等缺点，因而难以普遍使用。
论文“一种简易的单目视觉位姿测量方法研究”(期刊，谷凤伟、高宏伟、姜月秋，2018年8月)中针对特定工作条件下机器人位姿测量问题，研究和设计了一种简易的基于单目视觉的位姿测量系统，提出了一种精确的靶标检测算法与单目视觉测量方法相结合的方案，可实现机器人位姿的测量。实验证明，该方法准确有效，满足测量要求。但是其中提出的靶标检测算法与单目视觉测量方法相结合的位姿测量方法，虽然能够实现对单轴气浮平台的位姿测量，但是不能用于三轴气浮平台的位姿测量，因而具有一定的局限性。
基于上述现有技术，本专利提出一种基于计算机视觉的气浮运动模拟器位姿测量方法。采用单目相机实现对单轴气浮平台的位姿测量；采用双目相机实现对三轴气浮平台的位姿测量，具有非接触、低成本、高精度的特点。
发明内容
本发明的目的是为了解决上述现有技术存在的问题，进而提供一种基于计算机视觉的气浮运动模拟器位姿测量装置及方法。
本发明的目的是通过以下技术方案实现的：
一种基于计算机视觉的气浮运动模拟器位姿测量装置及方法，针对单轴气浮平台，其运动形式为二维平动和一维转动，为了最大限度的使用这一先验信息，将其视为合作待测目标，并设计了一种适用于单轴气浮台的平面测量靶标。由于采用了合作目标的测量方案，标志点的深度信息可以通过测量靶标本身的几何关系计算得到，在相机选择上可以使用单目相机。在对单轴气浮台进行位姿测量时，首先在单轴气浮台上安装设计好的测量靶标，之后在测量图像中提取测量靶标的三维坐标，最后根据测量靶标的三维坐标反推单轴气浮平台的位置和姿态信息，实现对单轴气浮平台的位姿测量。
针对三轴气浮台，其运动形式为二维平动和三维转动。由于转动轴较多，使用平面测量靶标会导致靶标在特定姿态下被遮挡，使用立体靶标，又会因为标志点在三维空间中的几何特征匹配步骤消耗大量的计算时间。因此，更为有效的方式是将其视为非合作目标，使用特征点描述符的匹配来代替标志点在三维空间中的几何特征匹配。由于使用了非合作目标的测量方案，特征点的深度信息无法根据先验知识获得，在相机选择上只能使用双目相机。
在对三轴气浮台进行位姿测量时，首先在三轴气浮台上预先提取一系列特征点，之后在测量图像中提取特征点的三维坐标，最后根据特征点的三维坐标重构三轴气浮台台体坐标系，并求解三轴气浮台台体坐标系与世界坐标系之间的坐标变换关系，进而得到三轴气浮台的位置和姿态信息，完成对三轴气浮台的位姿测量。
本发明的有益效果为：
本发明采用基于计算机视觉原理的位姿测量方法，只需要工业相机以及相应的靶标标志点即可实现高精度测量，具有实现简单、非接触、低成本、高精度的特点。
本发明采用基于计算机视觉原理的位姿测量方法，该方法不会因积分累计误差而导致精度漂移，在仿真实验需要长时间运行时可以输出精度稳定的位姿信息。
本发明针对单轴气浮平台采用单目相机进行测量，针对三轴气浮平台采用双目相机进行测量，因而适用范围更加广泛，具有良好的适应性。
本发明中针对单轴气浮平台的基于计算机视觉的气浮运动模拟器位姿测量装置，相机安装在相机支架上，相机位置固定，气浮机器人运动时相机不会产生细微抖动，具有更高的稳定性和测量精度，并且装置的成本更低。
本发明中针对三轴气浮台的基于计算机视觉的气浮运动模拟器位姿测量装置，系统大大简化了组成，仅仅用相机支架和双目相机即可完成相关的测量，不含有视觉导航模块和仿真模拟模块等设备，这就在保证精度及稳定性的同时很大程度上降低了成本，因而更加具有实际意义。
附图说明
图1为本发明单轴气浮平台位姿测量示意图。
图2为本发明单轴气浮平台位姿测量流程图。
图3为本发明图像的预处理流程图。
图4为本发明非合作目标的位姿测量示意图。
图5为本发明三轴气浮平台位姿测量流程图。
图6为本发明三角测距模型示意图。
图7为本发明立体校正模型示意图。
图8为本发明搜索极值点的范围示意图。
图中的附图标记，1为单轴气浮平台，2为单目相机，3为相机支架，4为三轴气浮台，5为双目相机，6为双目相机支架。
具体实施方式
下面将结合附图对本发明做进一步的详细说明：本实施例在以本发明技术方案为前提下进行实施，给出了详细的实施方式，但本发明的保护范围不限于下述实施例。
如图1至图8所示，本实施例所涉及的基于计算机视觉的气浮运动模拟器位姿测量装置及方法，包括：针对的气浮台，按照运动形式的不同，可分为单轴气浮台和三轴气浮台两类。下面分别针对两类气浮台提出相应的位姿测量方法。
实施例1
针对单轴气浮平台1，其运动形式为二维平动和一维转动，为了最大限度的使用这一先验信息，将其视为合作待测目标，并设计了一种适用于单轴气浮台1的平面测量靶标。由于采用了合作目标的测量方案，标志点的深度信息可以通过测量靶标本身的几何关系计算得到，在相机选择上可以使用单目相机2。单轴气浮平台1的位姿测量示意图如图1所示，单目相机2安装在相机支架3上，单目相机2位于单轴气浮台1的上方，在单轴气浮台1上安装测量靶标。在对单轴气浮台1进行位姿测量时，首先在单轴气浮台1上安装设计好的测量靶标，之后在测量图像中提取测量靶标的三维坐标，最后根据测量靶标的三维坐标反推单轴气浮平台的位置和姿态信息，实现对单轴气浮平台1的位姿测量。
单轴气浮平台位姿测量过程如图2所示。
首先进行相机建模与标定。
单轴气浮台单目测量采用针孔成像模型，根据针孔成像模型，可以实现由相机坐标系到世界坐标系的坐标变换。
根据相似三角形原理，可以得到
Figure BDA0002677630240000041
其中，(X,Y，1)T为空间点P所成的像P′在像平面坐标系下的齐次坐标，(XC,YC,ZC,1)T为空间点P在相机坐标系下的齐次坐标，f为相机镜头的焦距。
空间点P在世界坐标系下的齐次坐标为(XW,YW,ZW,1)T，则有式(2)所示的变换关系
Figure BDA0002677630240000051
其中R和T分别为旋转矩阵和平移向量。
空间点所成的像P′在像素坐标系下的齐次坐标为(u,v,1)T，由像平面坐标系到像素坐标系的变换关系为
Figure BDA0002677630240000052
其中，u0、v0为像平面坐标系的原点在像素坐标系下的横坐标和纵坐标，dx、dy分别为像素坐标系下X轴、Y轴方向一个像素对应的像平面坐标系下的长度。
根据上述关系式，可以得到相机的针孔成像模型。
Figure BDA0002677630240000053
令
Figure BDA0002677630240000054
Figure BDA0002677630240000055
其中，ax、ay分别为相机在X轴、Y轴方向上的尺度因子，M1、M2分别为相机的内参矩阵和外参矩阵。
对于单目相机，标定的作用主要是确定相机模型中的内参矩阵，从而建立起三维世界到相机成像平面的几何模型。采用张氏标定法，采用棋盘格形式的二维平面靶标来完成标定。由于二维平面靶标相对于三维物体会缺少部分信息，因此要多次改变平面靶标的方位。
由建立的针孔成像模型可知
Figure BDA0002677630240000061
由于采用的是二维平面靶标，通过合理选择世界坐标系的位置可以使Zw＝0，即
Figure BDA0002677630240000062
其中，r1为绕X轴转动向量，r2为绕Y轴转动向量，t为平移向量。
由上式可知，图像平面与二维靶标平面之间存在单应关系，其单应矩阵H可表示为
H＝M1[r1 r2 t] (9)
其中，H为3×3的齐次矩阵，含有8个待求参数。
对于二维靶标中的每一个标志点P(i)，都可以获得一组对应的(u(i),v(i))和
Figure BDA0002677630240000063
可以为求解H提供2个独立的方程。因此共需要至少4个标志点才能确定H中的8个待求参数。
根据标志点的坐标对应求出H后，需进一步从中提取出内参矩阵。为方便进一步的计算，将H写成列向量组合的形式
H＝[h1 h2 h3]＝M1[r1 r2 t] (10)
即
r1＝M1 -1h1 (11)
r2＝M1 -1h2 (12)
t＝M1 -1h3 (13)
由于旋转矩阵是正交矩阵，因此r1和r2可以看作是一个正交矩阵中的两列，根据正交矩阵的性质可得
r1 Tr2＝0 (14)
r1 Tr1＝r2 Tr2 (15)
将r1和r2代入上式可得
Figure BDA0002677630240000071
Figure BDA0002677630240000072
内参矩阵M1中共包含5个未知参数。为了求解这些参数，至少需要3个单应矩阵，因为3个单应矩阵在上式中的2个约束下可以产生6个方程。
因此在标定时，至少应从3个不同角度拍摄二维平面靶标的图片，来产生3个不同的单应矩阵，从而求解出内参矩阵，实现对单目相机的标定。
完成相机的建模和标定后，需要进行目标特征的提取与匹配。
在对相机拍摄的测量图像进行标志点提取和特征点提取前，需要对测量图像进行预处理，从而使得测量图像中的特征更加鲜明，同时也可以滤除大部分干扰噪声的影响。图像预处理操作包括灰度化、二值化、形态学处理等步骤。如图3所示。
灰度化操作的目的是将相机采集到的彩色图像转化为灰度图像，方便之后的二值化处理。灰度化操作的原理是将一个像素上RGB三个颜色通道上的亮度值加权求和作为这个像素的灰度值。具体计算公式如下
IGray＝IR×0.299+IG×0.587+IB×0.114 (18)
其中IGray为计算得到的像素灰度值，IR、IG和IB分别为该像素RGB三个通道上的亮度值，加权系数为查阅资料得到的经验值。
二值化操作的目的是将图像中的前景点和背景点区分出来，方便之后的轮廓提取操作。二值化操作的原理是将每个像素的灰度值与固定的阈值比较大小，根据比较的结果判断该像素为前景点或背景点。具体计算公式如下
Figure BDA0002677630240000073
其中，T为比较时使用的阈值；IBinary为二值化图像的灰度取值，可取0或1，0为白色，1为黑色。
形态学处理的目的是减少图像中的噪声点，避免对之后的质心提取操作造成干扰。主要用到两种最基本的形态学滤波操作，腐蚀和膨胀。腐蚀是对原图中的高亮部分进行腐蚀，腐蚀后图像的高亮区域比原图更小。膨胀是对图像中的高亮部分进行膨胀，膨胀后图像的高亮区域比原图更大。腐蚀或者膨胀操作的原理是将图像(或图像的一部分区域，记做A)与核(记做B,核的形状为小的长度可变的长方形或椭圆形)进行运算。可以将核视作形态学变换的模板。
腐蚀是求局部最小值的操作。核B与图形A卷积，即计算核B覆盖的区域的像素点的最小值，并把这个最小值赋值给参考点指定的像素。这样就会使图像中的含有噪声的高亮部分被腐蚀。
形态学腐蚀的计算公式如下
Figure BDA0002677630240000081
其中(x,y)为图像的像素点，(x′,y′)为核上的像素点。
膨胀是和腐蚀相反的操作，用来求图像的局部极大值。通过计算核B覆盖的区域的像素点的最大值并赋值给参考点指定的像素，使图像中的高亮区域逐渐膨胀。
形态学膨胀的计算公式如下
Figure BDA0002677630240000082
经过形态学腐蚀和形态学膨胀操作后，二值图像中的孤立噪声点会被消除，便于之后对标志点的提取。
对于合作待测目标单轴气浮平台，目标特征提取的目的是求取标志物中发光标志点的坐标。图像预处理操作得到的二值图像中，白色的前景点组成了标志点的轮廓。为了实现对标志点质心坐标的求解，首先要从二值图像中提取出标志点的轮廓，之后根据标志点的轮廓计算标志点的质心坐标。
采用区域生长的方法提取轮廓，基本方法是从一组种子点开始，将与种子预先定义的性质(如特定范围的灰度或颜色)相似的那些邻域像素添加到每个种子上来形成这些生长区域。得到的每个连通域，都对应一个标志点的轮廓。
在二值图像中提取得到标志点的轮廓后，下一步是计算标志点的质心坐标，采用基于二维高斯曲面拟合的亚像素质心定位方法。
对于标志点中包含的轮廓内点，其在图像坐标系下的坐标与灰度值满足二维高斯分布。将标志点的轮廓内点个数记为N，对于第i个轮廓内点，将其在图像坐标系下的坐标记作(ui,vi)，将其对应的像素灰度值记作pi；对于待求的标志点轮廓质心，将其在图像坐标系下的坐标记作(u0,v0)。
轮廓内点满足的二维高斯分布为
Figure BDA0002677630240000091
其中，p0为高斯分布的标准灰度值，σu、σv分别为图像坐标系中u轴、v轴方向上的标准差。
为了方便地进行后续的计算和拟合，对上式两边取自然对数，可得
Figure BDA0002677630240000092
对于标志点轮廓中的N个轮廓内点，可以得到N个如上式所示的二维高斯分布方程，将其表示为矩阵描述形式，可得
Figure BDA0002677630240000093
简记为
A＝BC (25)
矩阵A中的元素为
A＝[ai]N×1＝[pilnpi]N×1 (26)
矩阵B中的元素为
Figure BDA0002677630240000094
矩阵C中的元素为
Figure BDA0002677630240000101
由式(25)、式(26)、式(27)可知，矩阵A和矩阵B中的元素全部来自于N个轮廓内点的坐标和灰度值，为已知量；矩阵C中的元素为待求的质心坐标和二维高斯分布参数，为未知量。
因此，需要构造误差向量，描述N个轮廓内点拟合二维高斯分布的产生的误差，并且根据误差取极小值的条件建立矩阵C与矩阵A和矩阵B的约束关系，从而求得矩阵C中的元素，进而求得待求质心坐标(u0,v0)。
对于标志点轮廓中的N个轮廓内点，误差向量为
E＝A-BC (29)
考虑使用最小二乘法对N个轮廓内点进行拟合，使得均方差取最小值
Figure BDA0002677630240000102
由于标志点轮廓的轮廓内点个数N往往比较大，考虑对最小二乘法拟合的过程进行优化，从而减小计算量。
对于矩阵B，进行QR分解，得到N×N的正交矩阵Q和N×5的上三角矩阵R，对于式(29)两边同乘QT，可得
QTE＝QTA-QTBC＝QTA-RC (31)
由于Q为正交矩阵，可得
Figure BDA0002677630240000103
对矩阵QTA和矩阵R进行分块操作，可得
Figure BDA0002677630240000104
其中，S为QTA的前5行，T为QTA的后N-5行，R1为R的前5行组成的5×5上三角方阵。
将式(32)、式(33)代入式(30)可得
Figure BDA0002677630240000111
由式(34)可得，使得二维高斯分布对于全体轮廓内点的均方差最小，当且仅当
S-R1C＝0 (35)
即
Figure BDA0002677630240000112
由此可求得矩阵C中的元素
Figure BDA0002677630240000113
进一步可求得标志点的质心坐标
Figure BDA0002677630240000114
通过上述方法可求得标志点的质心在图像坐标系下的坐标(u0,v0)，完成了对合作目标的目标特征提取。
根据以上步骤可以从一帧测量图像中提取出该时刻待测目标在图像坐标系下的坐标，根据建立的针孔成像模型，可以在将其进一步变换到世界坐标下，实现对相机测量数据的坐标变换。
首先将图像坐标系转换为归一化坐标系。
图像坐标系下的齐次坐标为
Figure BDA0002677630240000115
归一化坐标系下的齐次坐标为
Figure BDA0002677630240000116
相机标定得到的内参矩阵为
Figure BDA0002677630240000121
根据相机的针孔成像模型可得
Figure BDA0002677630240000122
Figure BDA0002677630240000123
式(42)两边同乘M1 -1可得
Figure BDA0002677630240000124
接下来将归一化坐标系转换为相机坐标系。
相机坐标系下的齐次坐标为
Figure BDA0002677630240000125
根据归一化坐标的定义，可得
Figure BDA0002677630240000126
根据式(45)和式(46)可得
Figure BDA0002677630240000127
然而，在式(47)中，ZC为待求量。因此，需要根据先验信息，求出ZC，使用的先验条件为所有待测目标点共面。
最后将相机坐标系转换为世界坐标系。
世界坐标系下的齐次坐标为
Figure BDA0002677630240000131
世界坐标系Oxy平面相对于相机的外参矩阵为
Figure BDA0002677630240000132
其中，R为外参矩阵中的旋转矩阵,
Figure BDA0002677630240000133
θ,β,α分别为绕z,y,x轴旋转的角度；T为外参矩阵中的平移向量，tx,ty,tz分别为沿x,y,z轴平移的距离；
根据相机的针孔成像模型可得
Figure BDA0002677630240000134
Figure BDA0002677630240000135
式(50)两边同乘M2 -1可得
Figure BDA0002677630240000136
将式(47)代入上式可得
Figure BDA0002677630240000137
即
Figure BDA0002677630240000141
取其中的第三行，整理可得
ZW＝(r13XN+r23YN+r33)ZC-(r13tx+r23ty+r33tz) (55)
通过合理定义世界坐标系，可以使待测目标与世界坐标系的Oxy平面满足共面关系，从而保证ZW＝0。代入式(55)中，可得
(r13XN+r23YN+r33)ZC-(r13tx+r23ty+r33tz)＝0 (56)
整理可得
Figure BDA0002677630240000142
由式(57)可得，已知相机外参和归一化后的坐标，就可以计算出ZC。
得到标志点在世界坐标系下的坐标之后，需要进行位姿解算，由于待测目标的运动形式为二维平动和一维转动，且转动的转轴始终与二维平动的平面垂直，因此考虑采用如下方式对待测目标的运动进行建模。
在待测目标运动平面上，建立与大地固连的世界坐标系W，使其Oxy坐标平面位于待测目标的运动平面上，Oz坐标轴与待测目标运动平面垂直。在待测目标上建立与其自身固连的物体坐标系S，使其Oxy坐标平面与世界坐标系W的Oxy坐标平面平行，Ox坐标轴与OA方向重合(A点为与待测目标固连的标志点)，Oz坐标轴与其自身的转动轴重合。
将物体坐标系S的坐标原点O在世界坐标系下的坐标定义为待测目标的位置，将物体坐标系S的Ox轴与世界坐标系W的Ox之间的夹角定义为待测目标的姿态。
由于两个标志点A和B始终与待测目标固连，因此其在物体坐标系S下的坐标在任意时刻均相同记作
Figure BDA0002677630240000151
Figure BDA0002677630240000152
将从ti时刻下的物体坐标系Si到世界坐标系W的坐标变换矩阵记作
Figure BDA0002677630240000153
则有
Figure BDA0002677630240000154
根据坐标变换矩阵的定义，可得
Figure BDA0002677630240000155
Figure BDA0002677630240000156
将从ti+1时刻下的物体坐标系Si+1到世界坐标系W的坐标变换矩阵记作
Figure BDA0002677630240000157
则有
Figure BDA0002677630240000158
根据坐标变换矩阵的定义，可得
Figure BDA0002677630240000159
Figure BDA00026776302400001510
根据式(60)和式(63)得到坐标变换矩阵
Figure BDA00026776302400001511
和
Figure BDA00026776302400001512
后，根据坐标变换矩阵的运算关系，可进一步求得ti+1时刻下的物体坐标系Si+1到ti时刻下的物体坐标系Si的姿态更新矩阵
Figure BDA00026776302400001513
Figure BDA00026776302400001514
对于姿态更新矩阵
Figure BDA0002677630240000161
也可以将其表示为旋转矩阵和平移向量的形式，如式(67)所示
Figure BDA0002677630240000162
根据旋转矩阵
Figure BDA0002677630240000163
和平移向量
Figure BDA0002677630240000164
即可求得相邻两时刻间待测目标的二维平动位移和一维转动角度，从而实现对待测目标的位姿测量。
根据上述模型，为了求得待测目标的位姿，需要对相邻两时刻待测目标的世界坐标进行接解算，求出待测目标的二维平动位移和一维转动角度。位姿解算算法主要分为两步，分别为初始位姿确定和位姿迭代更新。
首先求解待测目标的初始位置。
初始时刻t0对应的物体坐标系S0的坐标原点O在世界坐标系下的坐标如式(68)所示
Figure BDA0002677630240000165
为了求解
Figure BDA0002677630240000166
需要使待测目标绕其转动轴进行3次旋转操作，并记录下每次旋转操作后，标志点A在世界坐标系下的坐标，分别记作
Figure BDA0002677630240000167
Figure BDA0002677630240000168
Figure BDA0002677630240000171
由标志点A与目标始终保持固连关系可知，标志点A与坐标原点O之间的距离d保持不变，即
Figure BDA0002677630240000172
根据式(72)可求解出
Figure BDA0002677630240000173
和d，根据建立的位姿测量模型中对待测目标位置的定义，
Figure BDA0002677630240000174
即为待测目标的初始位置。
接下来对目标的原始姿态进行求解。根据初始时刻t0采集得到的测量图片I0可以提取出标志点A在世界坐标系下的坐标
Figure BDA0002677630240000175
进一步可求得物体坐标系Ox轴方向上的单位向量
Figure BDA0002677630240000176
如式(73)所示。
Figure BDA0002677630240000177
世界坐标系Ox轴方向上的单位向量记作
Figure BDA0002677630240000178
如式(74)所示。
Figure BDA0002677630240000179
Figure BDA00026776302400001710
与
Figure BDA00026776302400001711
之间的夹角记作θ0，则
Figure BDA00026776302400001712
根据建立的位姿测量模型中对待测目标姿态的定义，θ0即为待测目标的初始姿态。
最后，根据求得的初始位置
Figure BDA0002677630240000181
和初始姿态θ0计算标志点A和B在物体坐标系S下的坐标PA和PB。
将从t0时刻的物体坐标系S0到世界坐标系W的坐标变换矩阵记作
Figure BDA0002677630240000182
将其写为旋转矩阵
Figure BDA0002677630240000183
和平移向量
Figure BDA0002677630240000184
组合的形式，如式(76)所示。
Figure BDA0002677630240000185
旋转矩阵
Figure BDA0002677630240000186
可由初始姿态θ0计算得到，如式(77)所示。
Figure BDA0002677630240000187
平移向量
Figure BDA0002677630240000188
可由初始位置
Figure BDA0002677630240000189
计算得到，如式(78)所示。
Figure BDA00026776302400001810
将式(77)、式(78)代入式(76)得
Figure BDA00026776302400001811
可求得从世界坐标系W到t0时刻的物体坐标系S0的坐标变换矩阵
Figure BDA00026776302400001812
Figure BDA00026776302400001813
进一步可求得PA和PB，如式(81)和式(82)所示。
Figure BDA0002677630240000191
Figure BDA0002677630240000192
其中，
Figure BDA0002677630240000193
为B点的初始位置。
求得待测目标的初始位姿后，采用迭代更新的方式求解之后每个时刻的位姿。位姿迭代更新实质上是根据
Figure BDA0002677630240000194
和
Figure BDA0002677630240000195
求解
Figure BDA0002677630240000196
联立式(64)和式(67)可得
Figure BDA0002677630240000197
联立式(65)和式(68)可得
Figure BDA0002677630240000198
由
Figure BDA0002677630240000199
和
Figure BDA00026776302400001910
的结构可知，
Figure BDA00026776302400001911
应具有如下结构形式，C,S,X,Y为四个不同的未知数，C表示cosθ，S表示sinθ，X,Y分别表示x,y；
Figure BDA00026776302400001912
将式(85)代入式(83)可得
Figure BDA00026776302400001913
将式(85)代入式(84)可得
Figure BDA0002677630240000201
根据式(86)的前两行和式(87)的前两行，可得
Figure BDA0002677630240000202
根据式(88)可求得C、S、X、Y的值，即求出了
Figure BDA0002677630240000203
中未知元素的值。求得
Figure BDA0002677630240000204
后，可以根据式(89)进一步求得
Figure BDA0002677630240000205
Figure BDA0002677630240000206
求得
Figure BDA0002677630240000207
后，再根据
Figure BDA0002677630240000208
求解ti+1时刻的位置
Figure BDA0002677630240000209
和姿态θi+1。
设
Figure BDA00026776302400002010
Figure BDA00026776302400002011
代入式(66)得
Figure BDA00026776302400002012
其中
Figure BDA0002677630240000211
Figure BDA0002677630240000212
根据式(89)和式(92)可以得到cos(θi+1-θi)和sin(θi+1-θi)的值，进一步可求得Δθ＝θi+1-θi，从而得到ti+1时刻的姿态θi+1。
根据式(89)和式(92)可以得到
Figure BDA0002677630240000213
和
Figure BDA0002677630240000214
的值，代入式(93)和式(94)中可进一步求得
Figure BDA0002677630240000215
和
Figure BDA0002677630240000216
从而得到ti+1时刻的位置
Figure BDA0002677630240000217
最后，根据
Figure BDA0002677630240000218
和θi+1求出
Figure BDA0002677630240000219
用于计算下一时刻的位姿信息。
由此便得到了单轴气浮平台的位姿信息。
实施例2
针对三轴气浮台，其运动形式为二维平动和三维转动。由于转动轴较多，使用平面测量靶标会导致靶标在特定姿态下被遮挡，使用立体靶标，又会因为标志点在三维空间中的几何特征匹配步骤消耗大量的计算时间。因此，更为有效的方式是将其视为非合作目标，使用特征点描述符的匹配来代替标志点在三维空间中的几何特征匹配。由于使用了非合作目标的测量方案，特征点的深度信息无法根据先验知识获得，在相机选择上只能使用双目相机。
在对三轴气浮台进行位姿测量时，首先在三轴气浮台上预先提取一系列特征点，之后在测量图像中提取特征点的三维坐标，最后根据特征点的三维坐标重构三轴气浮台台体坐标系，并求解三轴气浮台台体坐标系与世界坐标系之间的坐标变换关系，进而得到三轴气浮台的位置和姿态信息，完成对三轴气浮台的位姿测量。非合作目标的位姿测量示意图如图4所示，双目相机5安装在双目相机支架6上，双目相机5位于三轴气浮台4的上方，将三轴气浮台4其视为非合作目标。
三轴气浮平台位姿测量过程如图5所示：
同样地，首先进行相机的建模与标定，由于采用双目测量，因此除针孔成像模型之外，还应建立三角测距模型和立体校正模型。
双目相机的三角测距模型基于三角测量原理，根据世界坐标系下的同一点在双目相机的左右视图图像坐标下的视差信息，可以在不借助先验信息的前提下，估计深度信息。
在建立相机的三角测距模型时，为了简化模型，使用平行双目作为建模的基本假设，即双目相机中的左右视图对应的像平面位于同一平面中。对于不满足平行双目假设的实际测量系统，使用立体校正模型可以将非平行双目测量系统转化为平行双目测量系统。
如图6所示，将空间点P在双目相机的左右视图对应的图像坐标系下的坐标分别记为(ul,vl)T和(ur,vr)T，将空间点P在相机坐标系下的坐标记为(XC,YC,ZC)T，将双目相机中左右相机光心之间的距离记为基线长度B，根据相似三角形原理，可以求得下变换关系
Figure BDA0002677630240000221
其中，f为相机的焦距。
在实际测量中，由于存在机械安装过程等引入的误差，满足平行双目假设的双目测量系统往往很难获得。因此，需要立体校正模型对双目测量系统进行校正。
立体校正的实质是求取双目测量系统中左右视图对应的两个图像坐标系之间的坐标变换关系。双目相机的立体校正模型示意图如图7所示，图中Ox1和Ox2为两相机的光心，起示意作用。
由此可得双目相机的立体校正模型
Pr＝RPl+T (96)
其中，Pl和Pr分别是空间点P在左右视图对应的图像坐标系下的坐标。
建立相机模型后，对双目相机进行标定。
对于双目相机，标定的主要作用是确定双目测量系统中左右视图对应的两个图像坐标系之间的旋转矩阵R和平移向量T。
对于三维空间中的给定坐标点P，在双目测量系统中左右视图对应的两个图像坐标系中的坐标分别记为Pl和Pr，则有
Figure BDA0002677630240000222
Rl和Tl以及Rr和Tr，可以通过对双目测量系统中的左右两个相机分别进行单目相机标定来获取，对于双目相机标定算法来说视为已知参数。
对式(97)中的两个方程，消去其中的空间坐标点P，可得
Pr＝RrRl -1Pl-RrRl -1Tl (98)
对比式(98)与双目相机的立体校正模型，可得
R＝RrRl -1,T＝-RrRl -1Tl (99)
此时，双目测量系统中左右视图对应的两个图像坐标系之间的旋转矩阵R和平移向量T已经求得，完成了对双目相机的立体校正模型中未知参数的标定。
同样地，完成相机的建模和标定后，需要进行目标特征的提取与匹配。预处理部分与实施例1中单轴气浮平台部分相同。
预处理后需要进行目标的特征提取。对于非合作待测目标三轴气浮台，目标特征提取的目的是求取特征点的坐标和对应的特征描述符。采用SIFT描述符作为特征点的提取和匹配时使用的特征描述符，SIFT特征描述符具有尺度空间不变性，适用于位姿测量这类待测目标的尺度和方向会随时间发生改变的测量场景。
SIFT特征提取算法的基本思想是用不同尺度(标准差)的高斯卷积核对图像进行平滑处理，然后计算平滑后图像与平滑前图像的差异程度，差异大的像素点即为潜在的特征点。主要步骤有尺度空间极值检测、特征点精确定位、特征点主方向分配和特征点描述符的生成。
在尺度空间的建立方面，SIFT算法采用DoG尺度空间，并需要建立图像对应的高斯金字塔。首先将原图像作为最底层图像G0(高斯金字塔的第0层)，利用高斯核(5*5)对其进行卷积，然后对卷积后的图像进行下采样(去除偶数行和列)得到上一层图像G1，将此图像作为输入，重复卷积和下采样操作得到更上一层图像，反复迭代多次，形成一个金字塔形的图像数据结构，即高斯金字塔。
Figure BDA0002677630240000231
其中，0≤c≤Rn,0≤d≤Cn,Rn和Cn分别为高斯金字塔第n层的行数和列数，W(k,l)为二维可分离的长度为5的高斯卷积核，k,l分别表示卷积核的行和列,Gn(2c-k,2d-l)表示对应位置图像的像素。
尺度空间用一个变尺度的高斯函数G(x,y,σ)与图像I(x,y)卷积产生，即：
L(x,y,σ)＝G(x,y,σ)*I(x,y) (101)
其中x,y为图像空间坐标，σ为尺度坐标，G(x,y,σ)为：
Figure BDA0002677630240000241
DoG尺度空间的建立即为：
Figure BDA0002677630240000242
其中，k为两相邻尺度空间倍数的常数。
若将图像金字塔分为O组，每组有S层，则DoG中参数σ可表示为：
σ(o,s)＝σ0 2o+s/S (104)
σ(o,s)∈R,o∈Z,s∈N,σ0是基准层尺度，o为组的坐标，s为层坐标，由上式易知在DoG尺度空间中
Figure BDA0002677630240000243
如图8所示，极值点的搜索是在同一组内相邻的DoG层进行的，检测点要与跟它周围的8个点和上下两层相邻尺度的18个点比较来寻找极值点，即比较是在一个3x3的立方体内进行的。搜索过程从每组的第二层开始，假如每组只有4层，则只搜索两次，这样搜索到的点既有位置坐标(DoG的图像坐标)，又有空间尺度坐标(层坐标)。由于此处图像、层数是离散的，在离散空间里检测到的极值点不一定就是连续空间中的极值点，这里还要用到子像元插值法，求得连续空间里的极值。
特征点的精确定位是为了获得稳定的特征点。对于三维子像元插值，当一个三维点相对于插值中心点的偏移量在任何方向上大于0.5时，表明插值操作的中心点已经与其邻近点重合，这样的特征点不满足稳定条件。对于极值过小的点容易受噪声干扰，也要删除。在本课题中，通过对2x2的Hessian矩阵H求主曲率来完成判别。
Figure BDA0002677630240000244
其中Dxx(x,y),Dxy(x,y),Dyy(x,y)为对应的二阶偏导数。
令α为最大特征值，β为最小特征值，γ是α和β之间的比例，则有：
Tr(H)＝Dxx+Dyy＝α+β (106)
Det(H)＝DxxDyy-(Dxy)2＝α·β (107)
Figure BDA0002677630240000251
对于某一点检查主曲率，只需要下式成立：
Figure BDA0002677630240000252
则可以判定为稳定点，否则需要删除。
特征点主方向分配是基于图像局部的梯度方向，给每个特征点分配方向。对于已经提取出来的特征点，我们知道它的尺度值σ，便可以推导出该尺度值得高斯图像：
L(x,y,σ)＝G(x,y,σ)*I(x,y) (110)
以特征点为中心，计算图像梯度的幅角和幅值，公式如下：
Figure BDA0002677630240000253
Figure BDA0002677630240000254
其中L为关键点所在的尺度空间值。
计算完成后用直方图对范围内所有像素的方向进行统计，其中直方图将0゜～360゜分为36柱，每一柱为10゜，在统计中某一方向最多的(直方图最高的)即为该特征点的主方向，高度为最高峰的五分之四的角度为副方向，对主方向和副方向进行插值拟合，就可以得到精确的方向。
SIFT特征点描述子，描述的是SIFT特征点邻域内图像梯度的累计值。将特征点的邻域进行Bp×Bp个子区域的划分，每个子区域的尺寸为mσ个像元。然后以特征点作为旋转中心，将其邻域内图像梯度的位置和方向，进行角度为θ的旋转操作，旋转公式如下：
Figure BDA0002677630240000255
旋转后再以特征点为中心，取一个mσBp×mσBp的邻域，并将其均匀划分成Bp×Bp个子区域，每个间隔mσ个像元。与求特征点主方向时不同，此时的直方图将0゜～360゜分为8柱，计算出8个方向的梯度信息。又由于有Bp×Bp个子区域，最终生成的SIFT描述子的维度为Bp×Bp×8。
特征点提取后需要进行特征点的匹配，需要根据每个SIFT特征点的SIFT特征描述符匹配测量图像与待测目标模型中的SIFT特征点。特征点的匹配主要包括两个步骤，基于欧氏距离的特征点双向匹配和基于随机抽样一致性原则的特征点筛选。
基于欧式距离的特征点双向匹配方法，将SIFT特征描述符之间的欧式距离作为相似性度量准则，将测量图像和待测目标模型中的SIFT特征点进行双向匹配。
将测量图像中的SIFT特征点集合记作
Figure BDA0002677630240000261
其中，Lm为测量图像中SIFT特征点的数目。
将待测目标模型中的SIFT特征点集合记作
Figure BDA0002677630240000262
其中，Ln为待测目标模型中SIFT特征点的数目。
当SIFT特征向量的维度为l时，相似性度量准则欧氏距离的定义为
Figure BDA0002677630240000263
第一次匹配时，遍历测量图像的SIFT特征点集合Fm中的每一个特征点，计算该特征点与待测目标模型中SIFT特征点集合Fn中的每一个特征点的欧式距离，将其中的欧式距离最近的Fn中的特征点作为其匹配特征点。
第二次匹配时，交换测量图像的SIFT特征点集合Fm和待测目标模型中SIFT特征点集合Fn，求出与Fn中每个特征点对应的Fm中的特征点。
经过两次匹配之后，每个Fm中的特征点都有一个匹配的Fn中的特征点，每个Fn中的特征点也都有一个匹配的Fm中的特征点。若第一次匹配中，Fm中第a个特征点与Fn中第b个特征点匹配，且在第二次匹配中，Fn中第b个特征点与Fm中第a个特征点匹配，则Fm中第a个特征点与Fn中第b个特征点即为一组满足条件的特征点对。
基于随机抽样一致性原则的特征点筛选算法，用于测量图像与待测目标模型之间的误匹配特征点对筛选。
特征点筛选算法的步骤如下。
(1)随机从双向匹配得到的匹配特征点对中选取不共线的4对，计算出单应矩阵H。
(2)计算所有匹配特征点对于该单应矩阵H的投影误差，对于投影误差小于阈值的匹配特征点对，将其加入内点集I。
(3)对于第h次迭代计算，将内点集记作Ih；若Ih中的匹配特征点对的计数值大于当前最优内点集Im中的匹配特征点对的计数值，则将Im更新为Ih；
(4)重复进行迭代计算，直到迭代次数大于等于预设的最大迭代次数K。
经过上述步骤的筛选后，即可消除大部分误匹配SIFT特征点对，完成对非合作待测目标的特征匹配。
根据特征点提取算法，可以从一帧测量图像中提取出待测目标所包含的SIFT特征点在图像坐标系下的坐标。根据建立的相机的针孔成像模型和三角测距模型，可以进一步将SIFT特征点的坐标从图像坐标系变换到世界坐标下，实现对相机测量数据的坐标变换。
首先将图像坐标系转换为归一化坐标系。
对于双目相机的测量方案来说，图像坐标系到归一化坐标系的坐标变换这一步骤，与单目相机的测量方案相同，都是根据相机的针孔成像模型得到的。因此不再重复推导这一坐标变换过程，而是直接引用单目测量方案中的推导结果。
Figure BDA0002677630240000271
其中，
Figure BDA0002677630240000272
为图像坐标系下的齐次坐标，
Figure BDA0002677630240000273
为归一化坐标系下的齐次坐标，M1为相机的内参矩阵。
然后将归一化坐标系转换为相机坐标系。
相机坐标系下的齐次坐标为
Figure BDA0002677630240000274
将其转化为归一化坐标，可得
Figure BDA0002677630240000281
根据式(118)和式(119)可得
Figure BDA0002677630240000282
根据双目相机三角测距模型，可得
Figure BDA0002677630240000283
根据式(120)和式(121)可得
Figure BDA0002677630240000284
其中，B为双目相机的基线，f为镜头的焦距，d为空间点P在双目相机的左右视图中的视差。
最后将相机坐标系转换为世界坐标系。
对于双目相机的测量方案来说，相机坐标系到世界坐标系的坐标变换这一步骤，与单目相机的测量方案类似。因此直接引用单目测量方案中的推导结果。
Figure BDA0002677630240000285
其中，
Figure BDA0002677630240000286
为相机坐标系下的齐次坐标，
Figure BDA0002677630240000287
为世界坐标系下的齐次坐标，M2为世界坐标系相对于双目相机中左侧相机的外参矩阵。
得到标志点在世界坐标系下的坐标之后，需要进行位姿解算，非合作待测目标三轴气浮台的运动形式为二维平动和三维转动，因此采用如下方式对待测目标的运动进行建模。
在待测目标运动平面上，建立与大地固连的世界坐标系W，使其Oxy坐标平面位于待测目标的运动平面上，Oz坐标轴与待测目标运动平面垂直。在待测目标三轴气浮台上建立与其自身固连的物体坐标系S，坐标系S的原点位于三轴气浮台的气浮轴承球心，使其Oxy坐标平面与世界坐标系W的Oxy坐标平面平行。
将物体坐标系S的坐标原点O在世界坐标系下的坐标定义为待测目标的位置，将物体坐标系S相对于世界坐标系W的旋转矩阵定义为待测目标的姿态。
对于任意一帧测量图像T，假设从中提取到了N个SIFT特征点，则这些SIFT特征点在世界坐标系下的齐次坐标为
Figure BDA0002677630240000291
这些对应的待测目标模型中的匹配SIFT点为
Figure BDA0002677630240000292
将此时物体坐标系S到世界坐标系W的坐标变换矩阵记作
Figure BDA0002677630240000293
则有
Figure BDA0002677630240000294
根据坐标变换矩阵的定义，可得
Figure BDA0002677630240000295
得到上述位姿测量模型后，对于非合作待测目标的位姿问题即为：
根据N个匹配的SIFT特征点对在世界坐标系下的坐标，计算物体坐标系S到世界坐标系W的坐标变换矩阵
Figure BDA0002677630240000296
其中的平移向量
Figure BDA0002677630240000297
为待测目标的位置，旋转矩阵
Figure BDA0002677630240000298
为待测目标的姿态。
建立起位姿测量模型后进行位姿解算。
将N对匹配的SIFT特征点分别记为
P＝{p(1),p(2),...,p(N)} (128)
Q＝{q(1),q(2),...,q(N)} (129)
根据建立的非合作待测目标的位姿测量模型，可得
Figure BDA0002677630240000301
其中，将旋转矩阵
Figure BDA0002677630240000302
简记为R，将平移向量
Figure BDA0002677630240000303
简记为t。
之后，使用最小二乘误差值的最小值作为优化目标，求使误差平方和J达到极小的R和t
Figure BDA0002677630240000304
其中n＝N，为特征点的对数。
接下来对目标函数J进行化简。
定义两组SIFT特征点的质心
Figure BDA0002677630240000305
Figure BDA0002677630240000306
代入目标函数的表达式，化简可得
Figure BDA0002677630240000307
其中，结果中的第二项的化简结果如下
Figure BDA0002677630240000308
因此，目标函数J可化为
Figure BDA0002677630240000309
观察J中的两项
Figure BDA0002677630240000311
和
Figure BDA0002677630240000312
可以发现第一项只与旋转矩阵R有关，而第二项即含有旋转矩阵R又含有平移向量t。于是，上述问题可以分为以下三个步骤求解
(1)计算两组SIFT特征点的中心
Figure BDA00026776302400003111
和
Figure BDA00026776302400003112
然后计算每个点的去中心坐标
Figure BDA0002677630240000313
和
Figure BDA0002677630240000314
Figure BDA0002677630240000315
Figure BDA0002677630240000316
(2)根据以下优化问题计算旋转矩阵
Figure BDA0002677630240000317
(3)根据(2)中计算得到的旋转矩阵R*，计算平移向量
t*＝p-Rq (140)
求解的关键在于(2)中旋转矩阵的求解。
将(2)中的误差项展开可得
Figure BDA0002677630240000318
其中，第一项与R无关，第二项也与R无关(旋转矩阵的性质：RTR＝I)。
因此，优化的目标函数变为
Figure BDA0002677630240000319
为了求解旋转矩阵R，先定义矩阵W
Figure BDA00026776302400003110
W是一个3×3的矩阵，对W进行SVD分解，可得
W＝UΣVT (144)
其中，Σ为奇异值矩阵，U和V为正交矩阵。当W满秩时，满足上述优化条件的R为
R*＝UVT (145)
之后，根据(3)中的式(140)，可求得
Figure BDA0002677630240000321
至此，旋转矩阵
Figure BDA0002677630240000322
和平移向量
Figure BDA0002677630240000323
求解完成，位姿解算完毕。
以上所述，仅为本发明较佳的具体实施方式，这些具体实施方式都是基于本发明整体构思下的不同实现方式，而且本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求书的保护范围为准。

Claims (4)

1.一种基于计算机视觉的气浮运动模拟器位姿测量装置，其特征在于，包括：单轴气浮平台(1)、单目相机(2)和相机支架(3)；
其中，单目相机(2)安装在相机支架(3)上，单目相机(2)位于单轴气浮台(1)的上方，在单轴气浮台(1)上安装测量靶标。
2.根据权利要求1所述的基于计算机视觉的气浮运动模拟器位姿测量装置的测量方法，其特征在于，包括以下步骤：
步骤一：相机建模与标定
单轴气浮台单目测量采用针孔成像模型，根据针孔成像模型，可以实现由相机坐标系到世界坐标系的坐标变换；
根据相似三角形原理，可以得到
Figure FDA0002677630230000011
其中，(X,Y，1)T为空间点P所成的像P′在像平面坐标系下的齐次坐标，(XC,YC,ZC,1)T为空间点P在相机坐标系下的齐次坐标，f为相机镜头的焦距；
空间点P在世界坐标系下的齐次坐标为(XW,YW,ZW,1)T，则有式(2)所示的变换关系
Figure FDA0002677630230000012
其中R和T分别为旋转矩阵和平移向量；
空间点所成的像P′在像素坐标系下的齐次坐标为(u,v,1)T，由像平面坐标系到像素坐标系的变换关系为
Figure FDA0002677630230000021
其中，u0、v0为像平面坐标系的原点在像素坐标系下的横坐标和纵坐标，dx、dy分别为像素坐标系下X轴、Y轴方向一个像素对应的像平面坐标系下的长度；
根据上述关系式，可以得到相机的针孔成像模型；
Figure FDA0002677630230000022
令
Figure FDA0002677630230000023
Figure FDA0002677630230000024
其中，ax、ay分别为相机在X轴、Y轴方向上的尺度因子，M1、M2分别为相机的内参矩阵和外参矩阵；
对于单目相机，标定的作用主要是确定相机模型中的内参矩阵，从而建立起三维世界到相机成像平面的几何模型；采用张氏标定法，采用棋盘格形式的二维平面靶标来完成标定；由于二维平面靶标相对于三维物体会缺少部分信息，因此要多次改变平面靶标的方位；
由建立的针孔成像模型可知
Figure FDA0002677630230000025
由于采用的是二维平面靶标，通过合理选择世界坐标系的位置可以使Zw＝0，即
Figure FDA0002677630230000031
其中，r1为绕X轴转动向量，r2为绕Y轴转动向量，t为平移向量；
由上式可知，图像平面与二维靶标平面之间存在单应关系，其单应矩阵H可表示为
H＝M1[r1 r2 t] (9)
其中，H为3×3的齐次矩阵，含有8个待求参数；
对于二维靶标中的每一个标志点P(i)，都可以获得一组对应的(u(i),v(i))和
Figure FDA0002677630230000032
可以为求解H提供2个独立的方程；因此共需要至少4个标志点才能确定H中的8个待求参数；
根据标志点的坐标对应求出H后，需进一步从中提取出内参矩阵；为方便进一步的计算，将H写成列向量组合的形式
H＝[h1 h2 h3]＝M1[r1 r2 t] (10)
即
r1＝M1 -1h1 (11)
r2＝M1 -1h2 (12)
t＝M1 -1h3 (13)
由于旋转矩阵是正交矩阵，因此r1和r2可以看作是一个正交矩阵中的两列，根据正交矩阵的性质可得
r1 Tr2＝0 (14)
r1 Tr1＝r2 Tr2 (15)
将r1和r2代入上式可得
Figure FDA0002677630230000033
Figure FDA0002677630230000034
内参矩阵M1中共包含5个未知参数；为了求解这些参数，至少需要3个单应矩阵，因为3个单应矩阵在上式中的2个约束下可以产生6个方程；
因此在标定时，至少应从3个不同角度拍摄二维平面靶标的图片，来产生3个不同的单应矩阵，从而求解出内参矩阵，实现对单目相机的标定；
步骤二：目标的特征提取与匹配
完成相机的建模和标定后，需要进行目标特征的提取与匹配；
在对相机拍摄的测量图像进行标志点提取和特征点提取前，需要对测量图像进行预处理，从而使得测量图像中的特征更加鲜明，同时也可以滤除大部分干扰噪声的影响；图像预处理操作包括灰度化、二值化、形态学处理步骤；
灰度化操作的目的是将相机采集到的彩色图像转化为灰度图像，方便之后的二值化处理；灰度化操作的原理是将一个像素上RGB三个颜色通道上的亮度值加权求和作为这个像素的灰度值；具体计算公式如下
IGray＝IR×0.299+IG×0.587+IB×0.114 (18)
其中IGray为计算得到的像素灰度值，IR、IG和IB分别为该像素RGB三个通道上的亮度值，加权系数为查阅资料得到的经验值；
二值化操作的目的是将图像中的前景点和背景点区分出来，方便之后的轮廓提取操作；二值化操作的原理是将每个像素的灰度值与固定的阈值比较大小，根据比较的结果判断该像素为前景点或背景点；具体计算公式如下
Figure FDA0002677630230000041
其中，T为比较时使用的阈值；IBinary为二值化图像的灰度取值，可取0或1，0为白色，1为黑色；
形态学处理的目的是减少图像中的噪声点，避免对之后的质心提取操作造成干扰；主要用到两种最基本的形态学滤波操作，腐蚀和膨胀；腐蚀是对原图中的高亮部分进行腐蚀，腐蚀后图像的高亮区域比原图更小；膨胀是对图像中的高亮部分进行膨胀，膨胀后图像的高亮区域比原图更大；腐蚀或者膨胀操作的原理是将图像或图像的一部分区域与核进行运算，其中，图像或图像的一部分区域记做A，核记做B,核的形状为小的长度可变的长方形或椭圆形；
可以将核视作形态学变换的模板；
腐蚀是求局部最小值的操作；核B与图形A卷积，即计算核B覆盖的区域的像素点的最小值，并把这个最小值赋值给参考点指定的像素；这样就会使图像中的含有噪声的高亮部分被腐蚀；
形态学腐蚀的计算公式如下
Figure FDA0002677630230000051
其中(x,y)为图像的像素点，(x′,y′)为核上的像素点；
膨胀是和腐蚀相反的操作，用来求图像的局部极大值；通过计算核B覆盖的区域的像素点的最大值并赋值给参考点指定的像素，使图像中的高亮区域逐渐膨胀；
形态学膨胀的计算公式如下
Figure FDA0002677630230000052
经过形态学腐蚀和形态学膨胀操作后，二值图像中的孤立噪声点会被消除，便于之后对标志点的提取；
对于合作待测目标单轴气浮平台，目标特征提取的目的是求取标志物中发光标志点的坐标；图像预处理操作得到的二值图像中，白色的前景点组成了标志点的轮廓；为了实现对标志点质心坐标的求解，首先要从二值图像中提取出标志点的轮廓，之后根据标志点的轮廓计算标志点的质心坐标；
采用区域生长的方法提取轮廓，基本方法是从一组种子点开始，将与种子预先定义的性质，如特定范围的灰度或颜色；相似的那些邻域像素添加到每个种子上来形成这些生长区域；得到的每个连通域，都对应一个标志点的轮廓；
在二值图像中提取得到标志点的轮廓后，下一步是计算标志点的质心坐标，采用基于二维高斯曲面拟合的亚像素质心定位方法；
对于标志点中包含的轮廓内点，其在图像坐标系下的坐标与灰度值满足二维高斯分布；将标志点的轮廓内点个数记为N，对于第i个轮廓内点，将其在图像坐标系下的坐标记作(ui,vi)，将其对应的像素灰度值记作pi；对于待求的标志点轮廓质心，将其在图像坐标系下的坐标记作(u0,v0)；
轮廓内点满足的二维高斯分布为
Figure FDA0002677630230000053
其中，p0为高斯分布的标准灰度值，σu、σv分别为图像坐标系中u轴、v轴方向上的标准差；
为了方便地进行后续的计算和拟合，对上式两边取自然对数，可得
Figure FDA0002677630230000061
对于标志点轮廓中的N个轮廓内点，可以得到N个如上式所示的二维高斯分布方程，将其表示为矩阵描述形式，可得
Figure FDA0002677630230000062
简记为
A＝BC (25)
矩阵A中的元素为
A＝[ai]N×1＝[pilnpi]N×1 (26)
矩阵B中的元素为
Figure FDA0002677630230000063
矩阵C中的元素为
Figure FDA0002677630230000064
由式(25)、式(26)、式(27)可知，矩阵A和矩阵B中的元素全部来自于N个轮廓内点的坐标和灰度值，为已知量；矩阵C中的元素为待求的质心坐标和二维高斯分布参数，为未知量；
因此，需要构造误差向量，描述N个轮廓内点拟合二维高斯分布的产生的误差，并且根据误差取极小值的条件建立矩阵C与矩阵A和矩阵B的约束关系，从而求得矩阵C中的元素，进而求得待求质心坐标(u0,v0)；
对于标志点轮廓中的N个轮廓内点，误差向量为
E＝A-BC (29)
考虑使用最小二乘法对N个轮廓内点进行拟合，使得均方差取最小值
Figure FDA0002677630230000071
由于标志点轮廓的轮廓内点个数N往往比较大，考虑对最小二乘法拟合的过程进行优化，从而减小计算量；
对于矩阵B，进行QR分解，得到N×N的正交矩阵Q和N×5的上三角矩阵R，对于式(29)两边同乘QT，可得
QTE＝QTA-QTBC＝QTA-RC (31)
由于Q为正交矩阵，可得
Figure FDA0002677630230000072
对矩阵QTA和矩阵R进行分块操作，可得
Figure FDA0002677630230000073
其中，S为QTA的前5行，T为QTA的后N-5行，R1为R的前5行组成的5×5上三角方阵；
将式(32)、式(33)代入式(30)可得
Figure FDA0002677630230000074
由式(34)可得，使得二维高斯分布对于全体轮廓内点的均方差最小，当且仅当
S-R1C＝0 (35)
即
Figure FDA0002677630230000081
由此可求得矩阵C中的元素
Figure FDA0002677630230000082
进一步可求得标志点的质心坐标
Figure FDA0002677630230000083
通过上述方法可求得标志点的质心在图像坐标系下的坐标(u0,v0)，完成了对合作目标的目标特征提取；
根据以上步骤可以从一帧测量图像中提取出该时刻待测目标在图像坐标系下的坐标，根据建立的针孔成像模型，可以在将其进一步变换到世界坐标下，实现对相机测量数据的坐标变换；
首先将图像坐标系转换为归一化坐标系；
图像坐标系下的齐次坐标为
Figure FDA0002677630230000084
归一化坐标系下的齐次坐标为
Figure FDA0002677630230000085
相机标定得到的内参矩阵为
Figure FDA0002677630230000086
根据相机的针孔成像模型可得
Figure FDA0002677630230000091
Figure FDA0002677630230000092
式(42)两边同乘M1 -1可得
Figure FDA0002677630230000093
接下来将归一化坐标系转换为相机坐标系；
相机坐标系下的齐次坐标为
Figure FDA0002677630230000094
根据归一化坐标的定义，可得
Figure FDA0002677630230000095
根据式(45)和式(46)可得
Figure FDA0002677630230000096
然而，在式(47)中，ZC为待求量；因此，需要根据先验信息，求出ZC，使用的先验条件为所有待测目标点共面；
最后将相机坐标系转换为世界坐标系；
世界坐标系下的齐次坐标为
Figure FDA0002677630230000097
世界坐标系Oxy平面相对于相机的外参矩阵为
Figure FDA0002677630230000101
其中，R为外参矩阵中的旋转矩阵,
Figure FDA0002677630230000102
θ,β,α分别为绕z,y,x轴旋转的角度；T为外参矩阵中的平移向量，tx,ty,tz分别为沿x,y,z轴平移的距离；
根据相机的针孔成像模型可得
Figure FDA0002677630230000103
Figure FDA0002677630230000104
式(50)两边同乘M2 -1可得
Figure FDA0002677630230000105
将式(47)代入上式可得
Figure FDA0002677630230000106
即
Figure FDA0002677630230000107
取其中的第三行，整理可得
ZW＝(r13XN+r23YN+r33)ZC-(r13tx+r23ty+r33tz) (55)
通过合理定义世界坐标系，可以使待测目标与世界坐标系的Oxy平面满足共面关系，从而保证ZW＝0；代入式(55)中，可得
(r13XN+r23YN+r33)ZC-(r13tx+r23ty+r33tz)＝0 (56)
整理可得
Figure FDA0002677630230000111
由式(57)可得，已知相机外参和归一化后的坐标，就可以计算出ZC；
步骤三：位姿解算
得到标志点在世界坐标系下的坐标之后，需要进行位姿解算，由于待测目标的运动形式为二维平动和一维转动，且转动的转轴始终与二维平动的平面垂直，因此考虑采用如下方式对待测目标的运动进行建模；
在待测目标运动平面上，建立与大地固连的世界坐标系W，使其Oxy坐标平面位于待测目标的运动平面上，Oz坐标轴与待测目标运动平面垂直；在待测目标上建立与其自身固连的物体坐标系S，使其Oxy坐标平面与世界坐标系W的Oxy坐标平面平行，Ox坐标轴与OA方向重合，A点为与待测目标固连的标志点，Oz坐标轴与其自身的转动轴重合；
将物体坐标系S的坐标原点O在世界坐标系下的坐标定义为待测目标的位置，将物体坐标系S的Ox轴与世界坐标系W的Ox之间的夹角定义为待测目标的姿态；
由于两个标志点A和B始终与待测目标固连，因此其在物体坐标系S下的坐标在任意时刻均相同记作
Figure FDA0002677630230000112
Figure FDA0002677630230000121
将从ti时刻下的物体坐标系Si到世界坐标系W的坐标变换矩阵记作
Figure FDA0002677630230000122
则有
Figure FDA0002677630230000123
根据坐标变换矩阵的定义，可得
Figure FDA0002677630230000124
Figure FDA0002677630230000125
将从ti+1时刻下的物体坐标系Si+1到世界坐标系W的坐标变换矩阵记作
Figure FDA0002677630230000126
则有
Figure FDA0002677630230000127
根据坐标变换矩阵的定义，可得
Figure FDA0002677630230000128
Figure FDA0002677630230000129
根据式(60)和式(63)得到坐标变换矩阵
Figure FDA00026776302300001210
和
Figure FDA00026776302300001211
后，根据坐标变换矩阵的运算关系，可进一步求得ti+1时刻下的物体坐标系Si+1到ti时刻下的物体坐标系Si的姿态更新矩阵
Figure FDA00026776302300001212
Figure FDA00026776302300001213
对于姿态更新矩阵
Figure FDA00026776302300001214
也可以将其表示为旋转矩阵和平移向量的形式，如式(67)所示
Figure FDA00026776302300001215
根据旋转矩阵
Figure FDA0002677630230000131
和平移向量
Figure FDA0002677630230000132
即可求得相邻两时刻间待测目标的二维平动位移和一维转动角度，从而实现对待测目标的位姿测量；
根据上述模型，为了求得待测目标的位姿，需要对相邻两时刻待测目标的世界坐标进行接解算，求出待测目标的二维平动位移和一维转动角度；位姿解算算法主要分为两步，分别为初始位姿确定和位姿迭代更新；
首先求解待测目标的初始位置；
初始时刻t0对应的物体坐标系S0的坐标原点O在世界坐标系下的坐标如式(68)所示
Figure FDA0002677630230000133
为了求解
Figure FDA0002677630230000134
需要使待测目标绕其转动轴进行3次旋转操作，并记录下每次旋转操作后，标志点A在世界坐标系下的坐标，分别记作
Figure FDA0002677630230000135
Figure FDA0002677630230000136
Figure FDA0002677630230000137
由标志点A与目标始终保持固连关系可知，标志点A与坐标原点O之间的距离d保持不变，即
Figure FDA0002677630230000141
根据式(72)可求解出
Figure FDA0002677630230000142
和d，根据建立的位姿测量模型中对待测目标位置的定义，
Figure FDA0002677630230000143
即为待测目标的初始位置；
接下来对目标的原始姿态进行求解；根据初始时刻t0采集得到的测量图片I0可以提取出标志点A在世界坐标系下的坐标
Figure FDA0002677630230000144
进一步可求得物体坐标系Ox轴方向上的单位向量
Figure FDA0002677630230000145
如式(73)所示；
Figure FDA0002677630230000146
世界坐标系Ox轴方向上的单位向量记作
Figure FDA0002677630230000147
如式(74)所示；
Figure FDA0002677630230000148
Figure FDA0002677630230000149
与
Figure FDA00026776302300001410
之间的夹角记作θ0，则
Figure FDA00026776302300001411
根据建立的位姿测量模型中对待测目标姿态的定义，θ0即为待测目标的初始姿态；
最后，根据求得的初始位置
Figure FDA00026776302300001412
和初始姿态θ0计算标志点A和B在物体坐标系S下的坐标PA和PB；
将从t0时刻的物体坐标系S0到世界坐标系W的坐标变换矩阵记作
Figure FDA00026776302300001413
将其写为旋转矩阵
Figure FDA00026776302300001414
和平移向量
Figure FDA00026776302300001415
组合的形式，如式(76)所示；
Figure FDA0002677630230000151
旋转矩阵
Figure FDA0002677630230000152
可由初始姿态θ0计算得到，如式(77)所示；
Figure FDA0002677630230000153
平移向量
Figure FDA0002677630230000154
可由初始位置
Figure FDA0002677630230000155
计算得到，如式(78)所示；
Figure FDA0002677630230000156
将式(77)、式(78)代入式(76)得
Figure FDA0002677630230000157
可求得从世界坐标系W到t0时刻的物体坐标系S0的坐标变换矩阵
Figure FDA0002677630230000158
Figure FDA0002677630230000159
进一步可求得PA和PB，如式(81)和式(82)所示；
Figure FDA00026776302300001510
Figure FDA0002677630230000161
其中，
Figure FDA0002677630230000162
为B点的初始位置；
求得待测目标的初始位姿后，采用迭代更新的方式求解之后每个时刻的位姿；位姿迭代更新实质上是根据
Figure FDA0002677630230000163
和
Figure FDA0002677630230000164
求解
Figure FDA0002677630230000165
联立式(64)和式(67)可得
Figure FDA0002677630230000166
联立式(65)和式(68)可得
Figure FDA0002677630230000167
由
Figure FDA0002677630230000168
和
Figure FDA0002677630230000169
的结构可知，
Figure FDA00026776302300001610
应具有如下结构形式，C,S,X,Y为四个不同的未知数；
Figure FDA00026776302300001611
将式(85)代入式(83)可得
Figure FDA00026776302300001612
将式(85)代入式(84)可得
Figure FDA0002677630230000171
根据式(86)的前两行和式(87)的前两行，可得
Figure FDA0002677630230000172
根据式(88)可求得C、S、X、Y的值，即求出了
Figure FDA0002677630230000173
中未知元素的值；
求得
Figure FDA0002677630230000174
后，可以根据式(89)进一步求得
Figure FDA0002677630230000175
Figure FDA0002677630230000176
求得
Figure FDA0002677630230000177
后，再根据
Figure FDA0002677630230000178
求解ti+1时刻的位置
Figure FDA0002677630230000179
和姿态θi+1；
设
Figure FDA00026776302300001710
Figure FDA00026776302300001711
代入式(66)得
Figure FDA00026776302300001712
其中
Figure FDA00026776302300001713
Figure FDA0002677630230000181
根据式(89)和式(92)可以得到cos(θi+1-θi)和sin(θi+1-θi)的值，进一步可求得Δθ＝θi+1-θi，从而得到ti+1时刻的姿态θi+1；
根据式(89)和式(92)可以得到
Figure FDA0002677630230000182
和
Figure FDA0002677630230000183
的值，代入式(93)和式(94)中可进一步求得
Figure FDA0002677630230000184
和
Figure FDA0002677630230000185
从而得到ti+1时刻的位置
Figure FDA0002677630230000186
最后，根据
Figure FDA0002677630230000187
和θi+1求出
Figure FDA0002677630230000188
用于计算下一时刻的位姿信息；
由此便得到了单轴气浮平台的位姿信息。
3.一种基于计算机视觉的气浮运动模拟器位姿测量装置，其特征在于，包括：三轴气浮台(4)、双目相机(5)和双目相机支架(6)；
双目相机(5)安装在双目相机支架(6)上，双目相机(5)位于三轴气浮台(4)的上方，将三轴气浮台(4)其视为非合作目标。
4.根据权利要求3所述的基于计算机视觉的气浮运动模拟器位姿测量装置的测量方法，其特征在于，包括以下步骤：
步骤一：相机建模与标定
首先进行相机的建模与标定，由于采用双目测量，因此除针孔成像模型之外，还应建立三角测距模型和立体校正模型；
双目相机的三角测距模型基于三角测量原理，根据世界坐标系下的同一点在双目相机的左右视图图像坐标下的视差信息，可以在不借助先验信息的前提下，估计深度信息；
在建立相机的三角测距模型时，为了简化模型，使用平行双目作为建模的基本假设，即双目相机中的左右视图对应的像平面位于同一平面中；对于不满足平行双目假设的实际测量系统，使用立体校正模型可以将非平行双目测量系统转化为平行双目测量系统；
将空间点P在双目相机的左右视图对应的图像坐标系下的坐标分别记为(ul,vl)T和(ur,vr)T，将空间点P在相机坐标系下的坐标记为(XC,YC,ZC)T，将双目相机中左右相机光心之间的距离记为基线长度B，根据相似三角形原理，可以求得下变换关系
Figure FDA0002677630230000191
其中，f为相机的焦距；
在实际测量中，由于存在机械安装过程引入的误差，满足平行双目假设的双目测量系统往往很难获得；因此，需要立体校正模型对双目测量系统进行校正；
立体校正的实质是求取双目测量系统中左右视图对应的两个图像坐标系之间的坐标变换关系；
由此可得双目相机的立体校正模型
Pr＝RPl+T (96)
其中，Pl和Pr分别是空间点P在左右视图对应的图像坐标系下的坐标；
建立相机模型后，对双目相机进行标定；
对于双目相机，标定的主要作用是确定双目测量系统中左右视图对应的两个图像坐标系之间的旋转矩阵R和平移向量T；
对于三维空间中的给定坐标点P，在双目测量系统中左右视图对应的两个图像坐标系中的坐标分别记为Pl和Pr，则有
Figure FDA0002677630230000192
Rl和Tl以及Rr和Tr，可以通过对双目测量系统中的左右两个相机分别进行单目相机标定来获取，对于双目相机标定算法来说视为已知参数；
对式(97)中的两个方程，消去其中的空间坐标点P，可得
Pr＝RrRl -1Pl-RrRl -1Tl (98)
对比式(98)与双目相机的立体校正模型，可得
R＝RrRl -1,T＝-RrRl -1Tl (99)
此时，双目测量系统中左右视图对应的两个图像坐标系之间的旋转矩阵R和平移向量T已经求得，完成了对双目相机的立体校正模型中未知参数的标定；
步骤二：目标特征提取与匹配
完成相机的建模和标定后，需要进行目标特征的提取与匹配；预处理部分与单轴气浮平台部分相同；
预处理后需要进行目标的特征提取；对于非合作待测目标三轴气浮台，目标特征提取的目的是求取特征点的坐标和对应的特征描述符；采用SIFT描述符作为特征点的提取和匹配时使用的特征描述符，SIFT特征描述符具有尺度空间不变性，适用于位姿测量这类待测目标的尺度和方向会随时间发生改变的测量场景；
SIFT特征提取算法的基本思想是用不同尺度标准差的高斯卷积核对图像进行平滑处理，然后计算平滑后图像与平滑前图像的差异程度，差异大的像素点即为潜在的特征点；主要步骤有尺度空间极值检测、特征点精确定位、特征点主方向分配和特征点描述符的生成；
在尺度空间的建立方面，SIFT算法采用DoG尺度空间，并需要建立图像对应的高斯金字塔；首先将原图像作为最底层图像G0，即高斯金字塔的第0层，利用高斯核(5*5)对其进行卷积，然后对卷积后的图像进行下采样，去除偶数行和列，得到上一层图像G1，将此图像作为输入，重复卷积和下采样操作得到更上一层图像，反复迭代多次，形成一个金字塔形的图像数据结构，即高斯金字塔；
Figure FDA0002677630230000201
其中，0≤c≤Rn,0≤d≤Cn,Rn和Cn分别为高斯金字塔第n层的行数和列数，W(k,l)为二维可分离的长度为5的高斯卷积核，k,l分别表示卷积核的行和列，Gn(2c-k,2d-l)表示对应位置图像的像素；
尺度空间用一个变尺度的高斯函数G(x,y,σ)与图像I(x,y)卷积产生，即：
L(x,y,σ)＝G(x,y,σ)*I(x,y) (101)
其中x,y为图像空间坐标，σ为尺度坐标，G(x,y,σ)为：
Figure FDA0002677630230000202
DoG尺度空间的建立即为：
Figure FDA0002677630230000203
其中，k为两相邻尺度空间倍数的常数；
若将图像金字塔分为O组，每组有S层，则DoG中参数σ可表示为：
σ(o,s)＝σ0 2o+s/S (104)
σ(o,s)∈R,o∈Z,s∈N,σ0是基准层尺度，o为组的坐标，s为层坐标，由上式易知在DoG尺度空间中
Figure FDA0002677630230000211
极值点的搜索是在同一组内相邻的DoG层进行的，检测点要与跟它周围的8个点和上下两层相邻尺度的18个点比较来寻找极值点，即比较是在一个3x3的立方体内进行的；搜索过程从每组的第二层开始，假如每组只有4层，则只搜索两次，这样搜索到的点既有位置坐标，即DoG的图像坐标，又有空间尺度坐标，即层坐标；由于此处图像、层数是离散的，在离散空间里检测到的极值点不一定就是连续空间中的极值点，这里还要用到子像元插值法，求得连续空间里的极值；
特征点的精确定位是为了获得稳定的特征点；对于三维子像元插值，当一个三维点相对于插值中心点的偏移量在任何方向上大于0.5时，表明插值操作的中心点已经与其邻近点重合，这样的特征点不满足稳定条件；对于极值过小的点容易受噪声干扰，也要删除；在本课题中，通过对2x2的Hessian矩阵H求主曲率来完成判别；
Figure FDA0002677630230000212
其中Dxx(x,y),Dxy(x,y),Dyy(x,y)为对应的二阶偏导数；
令α为最大特征值，β为最小特征值，γ是α和β之间的比例，则有：
Tr(H)＝Dxx+Dyy＝α+β (106)
Det(H)＝DxxDyy-(Dxy)2＝α·β (107)
Figure FDA0002677630230000213
对于某一点检查主曲率，只需要下式成立：
Figure FDA0002677630230000221
则可以判定为稳定点，否则需要删除；
特征点主方向分配是基于图像局部的梯度方向，给每个特征点分配方向；对于已经提取出来的特征点，我们知道它的尺度值σ，便可以推导出该尺度值得高斯图像：
L(x,y,σ)＝G(x,y,σ)*I(x,y) (110)
以特征点为中心，计算图像梯度的幅角和幅值，公式如下：
Figure FDA0002677630230000222
Figure FDA0002677630230000223
其中L为关键点所在的尺度空间值；
计算完成后用直方图对范围内所有像素的方向进行统计，其中直方图将0゜～360゜分为36柱，每一柱为10゜，在统计中某一方向最多的，即直方图最高的，即为该特征点的主方向，高度为最高峰的五分之四的角度为副方向，对主方向和副方向进行插值拟合，就可以得到精确的方向；
SIFT特征点描述子，描述的是SIFT特征点邻域内图像梯度的累计值；将特征点的邻域进行Bp×Bp个子区域的划分，每个子区域的尺寸为mσ个像元；然后以特征点作为旋转中心，将其邻域内图像梯度的位置和方向，进行角度为θ的旋转操作，旋转公式如下：
Figure FDA0002677630230000224
旋转后再以特征点为中心，取一个mσBp×mσBp的邻域，并将其均匀划分成Bp×Bp个子区域，每个间隔mσ个像元；与求特征点主方向时不同，此时的直方图将0゜～360゜分为8柱，计算出8个方向的梯度信息；又由于有Bp×Bp个子区域，最终生成的SIFT描述子的维度为Bp×Bp×8；
特征点提取后需要进行特征点的匹配，需要根据每个SIFT特征点的SIFT特征描述符匹配测量图像与待测目标模型中的SIFT特征点；特征点的匹配主要包括两个步骤，基于欧氏距离的特征点双向匹配和基于随机抽样一致性原则的特征点筛选；
基于欧式距离的特征点双向匹配方法，将SIFT特征描述符之间的欧式距离作为相似性度量准则，将测量图像和待测目标模型中的SIFT特征点进行双向匹配；
将测量图像中的SIFT特征点集合记作
Figure FDA0002677630230000231
其中，Lm为测量图像中SIFT特征点的数目；
将待测目标模型中的SIFT特征点集合记作
Figure FDA0002677630230000232
其中，Ln为待测目标模型中SIFT特征点的数目；
当SIFT特征向量的维度为l时，相似性度量准则欧氏距离的定义为
Figure FDA0002677630230000233
第一次匹配时，遍历测量图像的SIFT特征点集合Fm中的每一个特征点，计算该特征点与待测目标模型中SIFT特征点集合Fn中的每一个特征点的欧式距离，将其中的欧式距离最近的Fn中的特征点作为其匹配特征点；
第二次匹配时，交换测量图像的SIFT特征点集合Fm和待测目标模型中SIFT特征点集合Fn，求出与Fn中每个特征点对应的Fm中的特征点；
经过两次匹配之后，每个Fm中的特征点都有一个匹配的Fn中的特征点，每个Fn中的特征点也都有一个匹配的Fm中的特征点；若第一次匹配中，Fm中第a个特征点与Fn中第b个特征点匹配，且在第二次匹配中，Fn中第b个特征点与Fm中第a个特征点匹配，则Fm中第a个特征点与Fn中第b个特征点即为一组满足条件的特征点对；
基于随机抽样一致性原则的特征点筛选算法，用于测量图像与待测目标模型之间的误匹配特征点对筛选；
特征点筛选算法的步骤如下；
(1)随机从双向匹配得到的匹配特征点对中选取不共线的4对，计算出单应矩阵H；
(2)计算所有匹配特征点对于该单应矩阵H的投影误差，对于投影误差小于阈值的匹配特征点对，将其加入内点集I；
(3)对于第h次迭代计算，将内点集记作Ih；若Ih中的匹配特征点对的计数值大于当前最优内点集Im中的匹配特征点对的计数值，则将Im更新为Ih；
(4)重复进行迭代计算，直到迭代次数大于等于预设的最大迭代次数K；
经过上述步骤的筛选后，即可消除大部分误匹配SIFT特征点对，完成对非合作待测目标的特征匹配；
根据特征点提取算法，可以从一帧测量图像中提取出待测目标所包含的SIFT特征点在图像坐标系下的坐标；根据建立的相机的针孔成像模型和三角测距模型，可以进一步将SIFT特征点的坐标从图像坐标系变换到世界坐标下，实现对相机测量数据的坐标变换；
首先将图像坐标系转换为归一化坐标系；
对于双目相机的测量方案来说，图像坐标系到归一化坐标系的坐标变换这一步骤，与单目相机的测量方案相同，都是根据相机的针孔成像模型得到的；因此不再重复推导这一坐标变换过程，而是直接引用单目测量方案中的推导结果；
Figure FDA0002677630230000241
其中，
Figure FDA0002677630230000242
为图像坐标系下的齐次坐标，
Figure FDA0002677630230000243
为归一化坐标系下的齐次坐标，M1为相机的内参矩阵；
然后将归一化坐标系转换为相机坐标系；
相机坐标系下的齐次坐标为
Figure FDA0002677630230000244
将其转化为归一化坐标，可得
Figure FDA0002677630230000245
根据式(118)和式(119)可得
Figure FDA0002677630230000251
根据双目相机三角测距模型，可得
Figure FDA0002677630230000252
根据式(120)和式(121)可得
Figure FDA0002677630230000253
其中，B为双目相机的基线，f为镜头的焦距，d为空间点P在双目相机的左右视图中的视差；
最后将相机坐标系转换为世界坐标系；
对于双目相机的测量方案来说，相机坐标系到世界坐标系的坐标变换这一步骤，与单目相机的测量方案类似；因此直接引用单目测量方案中的推导结果；
Figure FDA0002677630230000254
其中，
Figure FDA0002677630230000255
为相机坐标系下的齐次坐标，
Figure FDA0002677630230000256
为世界坐标系下的齐次坐标，M2为世界坐标系相对于双目相机中左侧相机的外参矩阵；
步骤三：位姿解算
得到标志点在世界坐标系下的坐标之后，需要进行位姿解算，非合作待测目标三轴气浮台的运动形式为二维平动和三维转动，因此采用如下方式对待测目标的运动进行建模；
在待测目标运动平面上，建立与大地固连的世界坐标系W，使其Oxy坐标平面位于待测目标的运动平面上，Oz坐标轴与待测目标运动平面垂直；在待测目标三轴气浮台上建立与其自身固连的物体坐标系S，坐标系S的原点位于三轴气浮台的气浮轴承球心，使其Oxy坐标平面与世界坐标系W的Oxy坐标平面平行；
将物体坐标系S的坐标原点O在世界坐标系下的坐标定义为待测目标的位置，将物体坐标系S相对于世界坐标系W的旋转矩阵定义为待测目标的姿态；
对于任意一帧测量图像T，假设从中提取到了N个SIFT特征点，则这些SIFT特征点在世界坐标系下的齐次坐标为
Figure FDA0002677630230000261
这些对应的待测目标模型中的匹配SIFT点为
Figure FDA0002677630230000262
将此时物体坐标系S到世界坐标系W的坐标变换矩阵记作
Figure FDA0002677630230000263
则有
Figure FDA0002677630230000264
根据坐标变换矩阵的定义，可得
Figure FDA0002677630230000265
得到上述位姿测量模型后，对于非合作待测目标的位姿问题即为：
根据N个匹配的SIFT特征点对在世界坐标系下的坐标，计算物体坐标系S到世界坐标系W的坐标变换矩阵
Figure FDA0002677630230000266
其中的平移向量
Figure FDA0002677630230000267
为待测目标的位置，旋转矩阵
Figure FDA0002677630230000268
为待测目标的姿态；
建立起位姿测量模型后进行位姿解算；
将N对匹配的SIFT特征点分别记为
P＝{p(1),p(2),...,p(N)} (128)
Q＝{q(1),q(2),...,q(N)} (129)
根据建立的非合作待测目标的位姿测量模型，可得
Figure FDA0002677630230000269
其中，将旋转矩阵
Figure FDA00026776302300002610
简记为R，将平移向量
Figure FDA00026776302300002611
简记为t；
之后，使用最小二乘误差值的最小值作为优化目标，求使误差平方和J达到极小的R和t
Figure FDA0002677630230000271
其中n＝N，为特征点的对数；
接下来对目标函数J进行化简；
定义两组SIFT特征点的质心
Figure FDA0002677630230000272
Figure FDA0002677630230000273
代入目标函数的表达式，化简可得
Figure FDA0002677630230000274
其中，结果中的第二项的化简结果如下
Figure FDA0002677630230000275
因此，目标函数J可化为
Figure FDA0002677630230000276
观察J中的两项
Figure FDA0002677630230000277
和
Figure FDA0002677630230000278
可以发现第一项只与旋转矩阵R有关，而第二项即含有旋转矩阵R又含有平移向量t；于是，上述问题可以分为以下三个步骤求解
(1)计算两组SIFT特征点的中心
Figure FDA0002677630230000279
和
Figure FDA00026776302300002710
然后计算每个点的去中心坐标
Figure FDA00026776302300002711
和
Figure FDA00026776302300002712
Figure FDA0002677630230000281
Figure FDA0002677630230000282
(2)根据以下优化问题计算旋转矩阵
Figure FDA0002677630230000283
(3)根据(2)中计算得到的旋转矩阵R*，计算平移向量
t*＝p-Rq (140)
求解的关键在于(2)中旋转矩阵的求解；
将(2)中的误差项展开可得
Figure FDA0002677630230000284
其中，第一项与R无关，第二项也与R无关(旋转矩阵的性质：RTR＝I)；
因此，优化的目标函数变为
Figure FDA0002677630230000285
为了求解旋转矩阵R，先定义矩阵W
Figure FDA0002677630230000286
W是一个3×3的矩阵，对W进行SVD分解，可得
W＝UΣVT (144)
其中，Σ为奇异值矩阵，U和V为正交矩阵；当W满秩时，满足上述优化条件的R为
R*＝UVT (145)
之后，根据(3)中的式(140)，可求得
Figure FDA0002677630230000287
至此，旋转矩阵
Figure FDA0002677630230000288
和平移向量
Figure FDA0002677630230000289
求解完成，位姿解算完毕。
CN202010952951.4A 2020-09-11 2020-09-11 基于计算机视觉的气浮运动模拟器位姿测量装置及方法 Active CN112066879B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010952951.4A CN112066879B (zh) 2020-09-11 2020-09-11 基于计算机视觉的气浮运动模拟器位姿测量装置及方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010952951.4A CN112066879B (zh) 2020-09-11 2020-09-11 基于计算机视觉的气浮运动模拟器位姿测量装置及方法

Publications (2)

Publication Number Publication Date
CN112066879A true CN112066879A (zh) 2020-12-11
CN112066879B CN112066879B (zh) 2021-11-19

Family

ID=73696192

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010952951.4A Active CN112066879B (zh) 2020-09-11 2020-09-11 基于计算机视觉的气浮运动模拟器位姿测量装置及方法

Country Status (1)

Country Link
CN (1) CN112066879B (zh)

Cited By (21)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112700501A (zh) * 2020-12-12 2021-04-23 西北工业大学 一种水下单目亚像素相对位姿估计方法
CN112927294A (zh) * 2021-01-27 2021-06-08 浙江大学 一种基于单敏感器的卫星定轨定姿方法
CN113112543A (zh) * 2021-04-08 2021-07-13 东方电气集团科学技术研究院有限公司 一种基于视觉移动目标的大视场二维实时定位系统及方法
CN113112545A (zh) * 2021-04-15 2021-07-13 西安电子科技大学 基于计算机视觉的手持移动打印装置定位方法
CN113137941A (zh) * 2021-04-29 2021-07-20 荣旗工业科技(苏州)股份有限公司 一种基于点激光测试产品平面度方法
CN113175929A (zh) * 2021-03-12 2021-07-27 南京航空航天大学 一种基于upf的空间非合作目标相对位姿估计方法
CN113436249A (zh) * 2021-06-01 2021-09-24 中国人民解放军63628部队 一种快速稳健的单目相机位姿估计算法
CN113432603A (zh) * 2021-06-25 2021-09-24 中国船舶重工集团公司第七0七研究所 一种复杂环境下母平台与子平台相对位姿测量方法
CN113884319A (zh) * 2021-08-30 2022-01-04 武汉科技大学 基于单目视觉的车辆最小转弯直径测量方法及系统
CN114062265A (zh) * 2021-11-11 2022-02-18 易思维（杭州）科技有限公司 一种用于视觉系统的支撑结构稳定性的评估方法
CN114166252A (zh) * 2022-02-10 2022-03-11 成都飞机工业（集团）有限责任公司 一种工业机器人集成系统综合定位精度测试方法
CN114199887A (zh) * 2021-12-13 2022-03-18 苏州华星光电技术有限公司 显示面板的曲面度外观检测设备
CN114777728A (zh) * 2022-04-06 2022-07-22 南京理工大学 三自由度气浮台获取自身姿态与定位的测量方法
CN114964024A (zh) * 2022-07-27 2022-08-30 承德建龙特殊钢有限公司 一种钢管测量装置
CN115761684A (zh) * 2023-01-10 2023-03-07 常熟理工学院 基于机器视觉的agv目标识别与姿态角解算方法及系统
CN116026342A (zh) * 2023-03-29 2023-04-28 中国科学院西安光学精密机械研究所 一种基于集群弹性分散的空间目标位姿测量方法
CN114926526B (zh) * 2022-05-23 2023-05-05 南京航空航天大学 一种基于变焦相机的位姿测量方法
CN116934871A (zh) * 2023-07-27 2023-10-24 湖南视比特机器人有限公司 一种基于标定物的多目系统标定方法、系统及存储介质
CN117020246A (zh) * 2023-10-09 2023-11-10 中北大学 一种大型筒体姿态调整方法
CN117315018A (zh) * 2023-08-31 2023-12-29 上海理工大学 基于改进PnP的用户面部位姿检测方法、设备、介质
CN117516485A (zh) * 2024-01-04 2024-02-06 东北大学 一种飞机发动机自动引导安装的位姿视觉测量方法

Citations (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103308028A (zh) * 2013-05-27 2013-09-18 哈尔滨工业大学 三轴气浮台姿态角双目立体视觉测量装置及其测量方法
CN104019745A (zh) * 2014-06-18 2014-09-03 福州大学 基于单目视觉间接标定方法的自由平面尺寸测量方法
CN104537671A (zh) * 2015-01-04 2015-04-22 长沙理工大学 一种基于机器视觉的香烟滤棒在线计数和质量检测方法
CN105760812A (zh) * 2016-01-15 2016-07-13 北京工业大学 一种基于Hough变换的车道线检测方法
CN106595638A (zh) * 2016-12-26 2017-04-26 哈尔滨工业大学 基于光电跟踪技术的三轴气浮台姿态测量装置及测量方法
JP2017219365A (ja) * 2016-06-06 2017-12-14 セイコーエプソン株式会社 位置姿勢算出装置、ロボット制御装置、及びロボット
JP6370038B2 (ja) * 2013-02-07 2018-08-08 キヤノン株式会社 位置姿勢計測装置及び方法
CN108645416A (zh) * 2018-03-30 2018-10-12 北京空间飞行器总体设计部 用于非合作目标相对导航仿真验证的视觉测量系统及方法
CN109087355A (zh) * 2018-06-15 2018-12-25 哈尔滨工业大学 基于迭代更新的单目相机位姿测量装置与方法
CN109099883A (zh) * 2018-06-15 2018-12-28 哈尔滨工业大学 高精度大视场机器视觉测量与标定装置及方法
KR20190051253A (ko) * 2017-11-06 2019-05-15 한국철도기술연구원 듀얼 카메라를 이용하여 카메라의 3차원 위치 오차 보정이 가능한 영상 변위 계측 방법
CN110084853A (zh) * 2019-04-22 2019-08-02 北京易达图灵科技有限公司 一种视觉定位方法及系统

Patent Citations (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP6370038B2 (ja) * 2013-02-07 2018-08-08 キヤノン株式会社 位置姿勢計測装置及び方法
CN103308028A (zh) * 2013-05-27 2013-09-18 哈尔滨工业大学 三轴气浮台姿态角双目立体视觉测量装置及其测量方法
CN104019745A (zh) * 2014-06-18 2014-09-03 福州大学 基于单目视觉间接标定方法的自由平面尺寸测量方法
CN104537671A (zh) * 2015-01-04 2015-04-22 长沙理工大学 一种基于机器视觉的香烟滤棒在线计数和质量检测方法
CN105760812A (zh) * 2016-01-15 2016-07-13 北京工业大学 一种基于Hough变换的车道线检测方法
JP2017219365A (ja) * 2016-06-06 2017-12-14 セイコーエプソン株式会社 位置姿勢算出装置、ロボット制御装置、及びロボット
CN106595638A (zh) * 2016-12-26 2017-04-26 哈尔滨工业大学 基于光电跟踪技术的三轴气浮台姿态测量装置及测量方法
KR20190051253A (ko) * 2017-11-06 2019-05-15 한국철도기술연구원 듀얼 카메라를 이용하여 카메라의 3차원 위치 오차 보정이 가능한 영상 변위 계측 방법
CN108645416A (zh) * 2018-03-30 2018-10-12 北京空间飞行器总体设计部 用于非合作目标相对导航仿真验证的视觉测量系统及方法
CN109087355A (zh) * 2018-06-15 2018-12-25 哈尔滨工业大学 基于迭代更新的单目相机位姿测量装置与方法
CN109099883A (zh) * 2018-06-15 2018-12-28 哈尔滨工业大学 高精度大视场机器视觉测量与标定装置及方法
CN110084853A (zh) * 2019-04-22 2019-08-02 北京易达图灵科技有限公司 一种视觉定位方法及系统

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
FORSYTH D A ET AL.: "Projective invariant representation using implicit algebraic curves", 《IMAGE AND VISION COMPUTE》 *
刘宇航等: "气浮实验台位姿视觉测量方法", 《红外与激光工程》 *
谷凤伟等: "一种简易的单目视觉位姿测量方法研究", 《光电技术应用》 *

Cited By (36)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112700501B (zh) * 2020-12-12 2024-03-05 西北工业大学 一种水下单目亚像素相对位姿估计方法
CN112700501A (zh) * 2020-12-12 2021-04-23 西北工业大学 一种水下单目亚像素相对位姿估计方法
CN112927294B (zh) * 2021-01-27 2022-06-10 浙江大学 一种基于单敏感器的卫星定轨定姿方法
CN112927294A (zh) * 2021-01-27 2021-06-08 浙江大学 一种基于单敏感器的卫星定轨定姿方法
CN113175929A (zh) * 2021-03-12 2021-07-27 南京航空航天大学 一种基于upf的空间非合作目标相对位姿估计方法
CN113175929B (zh) * 2021-03-12 2021-12-21 南京航空航天大学 一种基于upf的空间非合作目标相对位姿估计方法
CN113112543A (zh) * 2021-04-08 2021-07-13 东方电气集团科学技术研究院有限公司 一种基于视觉移动目标的大视场二维实时定位系统及方法
CN113112545A (zh) * 2021-04-15 2021-07-13 西安电子科技大学 基于计算机视觉的手持移动打印装置定位方法
CN113112545B (zh) * 2021-04-15 2023-03-21 西安电子科技大学 基于计算机视觉的手持移动打印装置定位方法
CN113137941A (zh) * 2021-04-29 2021-07-20 荣旗工业科技(苏州)股份有限公司 一种基于点激光测试产品平面度方法
CN113436249A (zh) * 2021-06-01 2021-09-24 中国人民解放军63628部队 一种快速稳健的单目相机位姿估计算法
CN113436249B (zh) * 2021-06-01 2023-04-11 中国人民解放军63628部队 一种快速稳健的单目相机位姿估计算法
CN113432603B (zh) * 2021-06-25 2022-09-16 中国船舶重工集团公司第七0七研究所 一种复杂环境下母平台与子平台相对位姿测量方法
CN113432603A (zh) * 2021-06-25 2021-09-24 中国船舶重工集团公司第七0七研究所 一种复杂环境下母平台与子平台相对位姿测量方法
CN113884319A (zh) * 2021-08-30 2022-01-04 武汉科技大学 基于单目视觉的车辆最小转弯直径测量方法及系统
CN113884319B (zh) * 2021-08-30 2024-04-12 武汉科技大学 基于单目视觉的车辆最小转弯直径测量方法及系统
CN114062265A (zh) * 2021-11-11 2022-02-18 易思维（杭州）科技有限公司 一种用于视觉系统的支撑结构稳定性的评估方法
CN114062265B (zh) * 2021-11-11 2023-06-30 易思维（杭州）科技有限公司 一种用于视觉系统的支撑结构稳定性的评估方法
CN114199887A (zh) * 2021-12-13 2022-03-18 苏州华星光电技术有限公司 显示面板的曲面度外观检测设备
CN114166252B (zh) * 2022-02-10 2022-05-10 成都飞机工业（集团）有限责任公司 一种工业机器人集成系统综合定位精度测试方法
CN114166252A (zh) * 2022-02-10 2022-03-11 成都飞机工业（集团）有限责任公司 一种工业机器人集成系统综合定位精度测试方法
CN114777728A (zh) * 2022-04-06 2022-07-22 南京理工大学 三自由度气浮台获取自身姿态与定位的测量方法
CN114926526B (zh) * 2022-05-23 2023-05-05 南京航空航天大学 一种基于变焦相机的位姿测量方法
CN114964024A (zh) * 2022-07-27 2022-08-30 承德建龙特殊钢有限公司 一种钢管测量装置
CN115761684A (zh) * 2023-01-10 2023-03-07 常熟理工学院 基于机器视觉的agv目标识别与姿态角解算方法及系统
CN115761684B (zh) * 2023-01-10 2023-06-16 常熟理工学院 基于机器视觉的agv目标识别与姿态角解算方法及系统
CN116026342B (zh) * 2023-03-29 2023-08-18 中国科学院西安光学精密机械研究所 一种基于集群弹性分散的空间目标位姿测量方法
CN116026342A (zh) * 2023-03-29 2023-04-28 中国科学院西安光学精密机械研究所 一种基于集群弹性分散的空间目标位姿测量方法
CN116934871A (zh) * 2023-07-27 2023-10-24 湖南视比特机器人有限公司 一种基于标定物的多目系统标定方法、系统及存储介质
CN116934871B (zh) * 2023-07-27 2024-03-26 湖南视比特机器人有限公司 一种基于标定物的多目系统标定方法、系统及存储介质
CN117315018A (zh) * 2023-08-31 2023-12-29 上海理工大学 基于改进PnP的用户面部位姿检测方法、设备、介质
CN117315018B (zh) * 2023-08-31 2024-04-26 上海理工大学 基于改进PnP的用户面部位姿检测方法、设备、介质
CN117020246A (zh) * 2023-10-09 2023-11-10 中北大学 一种大型筒体姿态调整方法
CN117020246B (zh) * 2023-10-09 2023-12-15 中北大学 一种大型筒体姿态调整方法
CN117516485A (zh) * 2024-01-04 2024-02-06 东北大学 一种飞机发动机自动引导安装的位姿视觉测量方法
CN117516485B (zh) * 2024-01-04 2024-03-22 东北大学 一种飞机发动机自动引导安装的位姿视觉测量方法

Also Published As

Publication number Publication date
CN112066879B (zh) 2021-11-19

Similar Documents

Publication Publication Date Title
CN112066879B (zh) 基于计算机视觉的气浮运动模拟器位姿测量装置及方法
CN109035320B (zh) 基于单目视觉的深度提取方法
CN107301654B (zh) 一种多传感器的高精度即时定位与建图方法
CN105913489B (zh) 一种利用平面特征的室内三维场景重构方法
Zhang et al. Vision-based pose estimation for textureless space objects by contour points matching
CN109708649B (zh) 一种遥感卫星的姿态确定方法及系统
CN105069743B (zh) 探测器拼接实时图像配准的方法
US20130156262A1 (en) Voting-Based Pose Estimation for 3D Sensors
JP2009093611A (ja) 三次元オブジェクト認識のためのシステムおよび方法
CN108596979A (zh) 一种用于激光雷达和深度相机的标定装置和方法
CN109146958B (zh) 一种基于二维图像的交通标志空间位置测量方法
Nagy et al. Online targetless end-to-end camera-LiDAR self-calibration
CN111915517A (zh) 一种适用于室内光照不利环境下rgb-d相机全局定位方法
CN111123242A (zh) 一种基于激光雷达和相机的联合标定方法及计算机可读存储介质
CN114463442A (zh) 一种非同轴相机的标定方法
Eichhardt et al. Affine correspondences between central cameras for rapid relative pose estimation
Kaufmann et al. Shadow-based matching for precise and robust absolute self-localization during lunar landings
CN112581511A (zh) 基于近似直立扫描点云快速配准的三维重建方法及系统
CN116563377A (zh) 一种基于半球投影模型的火星岩石测量方法
JP6761388B2 (ja) 推定装置及びプログラム
CN113012279B (zh) 一种非接触三维成像测量方法、系统及计算机可读存储介质
Wan et al. A performance comparison of feature detectors for planetary rover mapping and localization
CN112906573B (zh) 基于轮廓点集的行星表面导航路标匹配方法
Abratkiewicz et al. The concept of applying a SIFT algorithm and orthophotomaps in SAR-based augmented integrity navigation systems
CN114723811A (zh) 一种非结构化环境四足机器人立体视觉定位建图方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant