CN103308028B - 三轴气浮台姿态角双目立体视觉测量装置及其测量方法 - Google Patents

Abstract

本发明涉及了一种三轴气浮台姿态角双目立体视觉测量装置及其测量方法，测量装置由2台摄像机、摄像机基座、测量标靶、光源、计算机等组成，在三轴气浮台的上方安装两台黑白数字CCD摄像机，摄像机安装在摄像机基座上，人工光源安装在摄像机附近，立体测量标靶安装在气浮球轴承台面上。方法是：建立坐标系，两台摄像机各采集一幅立体测量标靶的数字图像并传输到计算机。综合两台摄像机各自与测量标靶之间的投影矩阵、以及测量光标的定位信息，利用立体视觉测量理论可以计算各标志点在世界坐标系下的三维坐标，再通过坐标转换，可以分解出气浮球轴承的绝对姿态角信息。通过这种方式实现对三轴气浮台的非接触、无扰动、高频率的测量。

Description

三轴气浮台姿态角双目立体视觉测量装置及其测量方法

技术领域

本发明涉及测量技术领域，具体涉及一种三轴气浮台姿态角双目立体视觉测量装置及其测量方法。

背景技术

三轴气浮台是空间飞行器运动模拟器，用于在全物理仿真实验中检验卫星控制系统的性能，是空间飞行器研制过程中的重要手段和方法。

三轴气浮台由球面气浮轴承和轴承座组成，轴承座上有气孔，压缩气瓶不断通过气孔向球面气浮轴承与轴承座之间输送压缩空气，使球面气浮轴承浮起，球面气浮轴承台面可以进行三维转动。由于球面气浮轴承和轴承底座之间摩擦极小，所以三轴气浮台的姿态变化可以模拟卫星在太空无摩擦力学环境中的运动。

三轴气浮台在试验过程中需要通过姿态测量系统动态地给出气浮台的三维运动角度、角速度，以便完成控制闭环。由于三轴气浮台需要为卫星姿态控制系统物理仿真提供无干扰环境，以往接触式的测量装置和方法(如旋转变压器、感应同步器、光电码盘、光栅等)不适用于三轴气浮台的测量，需要考虑新的测量方法和装置。

经检索文献发现，中国发明专利申请号：200610010260.2，专利名称为：三轴气浮台姿态角测量装置及其测量方法，该专利在三轴气浮台上方安装有CCD摄像机，在气浮台台面上安装有测量LED光标系统，利用计算机视觉理论结合测量光标点间的距离信息，计算出气浮台台面相对于摄像机的相对运动参数。该测量方法对光标点的几何结构有严格要求。

中国发明专利申请号：200610010435.X，专利名称为：非接触式三轴气浮台转角测量装置及其测量方法，该专利采用彩色CCD摄像机采集由4个红光LED和微型绿色激光器在测量靶标白色底板上形成的光点构成的图像，利用计算机视觉理论并结合测量靶标及CCD摄像机安装信息，计算出气浮台台面的转动参数。该方法必须采用彩色CCD摄像机与激光器实现，成本高。

在文献“三轴气浮台视觉测量系统研究”(发表于中文核心期刊《微计算机信息》(测控自动化)2008年第24卷第4-1期，124-126)中，北京理工大学的刘伟、徐斌介绍一种利用机器视觉识别物体颜色特征的方法来检测三轴气浮台的姿态变化的方法，通过摄像头获取的图像检测到目标点颜色，以此计算求出气浮平台的姿态。本文通过采用自适应阈值、快速预测搜索算法。但该系统使用颜色信息作为测量特征，其精度受到限制，同时该文没有给出明确的角度计算方法。

在文献“五自由度气浮台姿态的计算机视觉辅助确定”(2009年第41卷第4期，220-226)设计了一种“单目视觉+两轴倾角仪+三轴陀螺仪”的组合定姿方案来确定五自由度气浮台姿态信息，其中单目视觉模块只能解算出气浮台的偏航角。

在文献“基于单目视觉的三维姿态测量方法与系统实现”(哈尔滨工业大学工程硕士论文，2011年6月)介绍了一种基于单目视觉的姿态测量系统，该系统利用摄像机标定原理实现对目标的运动姿态实时测量，在这个过程中需要进行多次坐标转换以及参数解析。

发明内容

本发明的目的在于提供一种能够动态测量三轴气浮台三个姿态角、并且不会对气浮台产生干扰的三轴气浮台双目立体视觉测量装置及其测量方法。

本发明的目的是这样实现的：一种三轴气浮台姿态角双目立体视觉测量装置，包括三轴气浮台，安装在三轴气浮台上的气浮轴承台面能够浮起进行三维转动，还包括2台摄像机及其基座、立体测量标靶、辅助光源和计算机，在气浮轴承台面上方安装两台黑白数字CCD摄像机，摄像机安装在基座上，辅助光源安装在摄像机附近，立体测量标靶安装在气浮轴承台面上，立体测量标靶上设置有多个人工标志点，数字CCD摄像机连续采集人工标志点的图像并传输到计算机，计算机将图像信息进行分析与处理。

本发明还具有如下特征：

以上所述的数字CCD摄像机通过数据线与转换器相连，转换器与计算机相连，数字CCD摄像机的镜头为物方远心镜头。

本发明公开的一种三轴气浮台姿态角双目立体视觉测量方法，方法如下：

(1)、将立体测量标靶安装在三轴气浮台的台面上，建立世界坐标系；

(2)、立体测量标靶上有多个标志点，设置在三轴气浮台的台面上方的两台计算机对立体测量标靶的图像信息进行分析与处理，针对这些标志点完成特征提取和亚像素定位；

(3)、利用标志点的定位信息，根据计算机视觉成像原理，计算两台摄像机各自的图像像素坐标系与世界坐标系之间的投影矩阵H₁和H₂；

(4)、两台摄像机实时、连续地采集立体测量标靶上标志点的图像，并传输至计算机；

(5)、计算机对来自两台摄像机的两幅图像分别进行特征提取，利用亚像素定位方法确定上平面测量标靶的标志点分别在两幅图像中的图像像素坐标系中的坐标；

(6)、根据步骤(5)的结果，结合步骤(3)求得的投影矩阵H₁和H₂，利用双目视觉测量原理，计算三轴气浮台在世界坐标系中的三维坐标；

(7)、根据气浮球轴承的坐标转换关系，解算其三维姿态角信息。

本发明具有以下优点：

1.本发明采用非接触式视觉测量，三轴气浮台的姿态角测量过程不会对气浮平台的运动产生干扰；

2.本发明采用黑白数字摄像机通过转换器向计算机传输数据，对三轴气浮台工作的电磁环境具有抗干扰能力；

3.本发明所涉及视觉测量技术，不需要对摄像机进行内外参数标定，计算过程简洁、精度高、速度快；

4.本发明所涉及视觉测量技术使用2台摄像机采集图像，它们的安装位置保证测量标靶始终位于2台摄像机的视场中即可；

5.本发明所涉及视觉测量技术，能够得到气浮球轴承的姿态角信息。

附图说明

图1为三轴气浮台双目立体视觉姿态角测量装置组成示意图；

图2为三轴气浮台双目立体视觉姿态角测量流程图；

图3为三轴气浮台图像处理模块流程图；

图4为本发明测量系统涉及到的坐标系，世界坐标系O_wx_wy_wz_w和气浮球轴承本体坐标系O_bx_by_bz_b的位置示意图；

图5为摄像机成像面的图像像素坐标系Ouv与摄像机坐标系O_cx_cy_cz_c的关系示意图。

具体实施方式

下面结合附图举例对本发明作进一步说明。

实施例1

参见图1，本发明主要由以下部分组成：

第一数字CCD摄像机101、第二数字CCD摄像机102、第一镜头103、第二镜头104、第一摄像机基座105、第二摄像机基座106、第一人工光源201、第二人工光源202、气浮球轴承台面301、球面气浮轴承302、轴承底座303、测量标靶304、仪表平台305、转换器4、计算机5。

第一和第二CCD摄像机101.102都为黑白数字CCD摄像机，使用转换器4与计算机5相连，抗干扰。

第一和第二摄像机101.102分别安装在第一和第二基座105.106上。

第一和第二镜头103.104为物方远心镜头。

第一和第二人工光源201.202为图像采集提供恒定、可靠的照明。

气浮球轴承302放置轴承底座303上，轴承底座303上有供气系统将高压气体通入气浮轴承中，底座和气浮轴承的接触面上设计有节流孔，可以在气浮轴承的球面和底座之间形成气膜，靠气体压力将球面气浮轴承台浮起，从而使气浮轴承台面具有姿态变换的能力。测量标靶304安装在气浮轴承台面301上，气浮轴承台面301运动时带动测量标靶进行相同的姿态角变换。

该系统的基本工作原理是：

将测量标靶放置于在气浮轴承平面上，两台摄像机各采集一幅测量标靶的数字图像并传输到计算机。立体测量标靶上有多个人工标志点，计算机对2幅立体测量标靶的图像信息进行分析与处理，针对这些标志点完成特征提取、亚像素定位。利用标志点的定位信息，根据计算机视觉成像原理，计算每一台摄像机与测量标靶之间的投影矩阵，也就是说建立每一台摄像机与测量标靶之间对应关系。

完成以上过程以后，两台摄像机连续采集测量标靶的数字图像并传输给计算机。计算机运行图像处理程序对测量光标进行特征提取、亚像素定位。综合两台摄像机各自与测量标靶之间的投影矩阵、以及测量光标的定位信息，利用立体视觉测量理论可以计算各标志点在世界坐标系下的三维坐标，再通过坐标转换，可以分解出气浮球轴承的姿态角信息。通过这种方式实现对三轴气浮台的非接触、无扰动、高频率的测量。

实施例2

本发明利用双目摄像机采集测量标靶的图像，利用立体视觉测量原理计算得到三轴气浮球轴承在本体坐标系下的姿态角。主要流程图示意图如图3所示，主要包括以下主要步骤：

1.将测量标靶安装于气浮轴承台面上，建立本发明测量系统中的世界坐标系，在这个坐标系中立体标靶上所有标志点具有唯一的三维坐标；

2.两台摄像机各采集一幅立体测量标靶的图像传输至计算机；

3.对两幅图像进行特征提取，利用亚像素定位技术确定两幅图像中所有标志点在图像像素坐标系中的坐标；相关算法参见于起峰等的著作“基于图像的精密测量与运动测量”，科学出版社出版，2002年7月第一版；

4.结合步骤2、3的结果，根据摄像机成像模型，计算摄像机1的图像像素坐标系与世界坐标系之间的投影矩阵H₁；

5.结合步骤2、3的结果，根据摄像机成像模型，计算摄像机2的图像像素坐标系与世界坐标系之间的投影矩阵H₂；

6.以上步骤只需要进行一次；

7.测量过程中，两台摄像机实时地、连续地采集平面测量标靶上标志点的图像传输至计算机；

8.图像处理程序对来自两台摄像机的两幅图像分别进行特征提取，利用亚像素定位技术确定标靶上的标志点分别在两幅图像中的图像像素坐标系中的坐标；

9.根据步骤8的结果，结合步骤4、5求得的投影矩阵H₁和H₂，利用双目视觉测量原理，计算测量标靶上各个标志点在世界坐标系中的三维坐标；

10.根据步骤9的结果，通过坐标变换，解算气浮球轴承在本体坐标系下的姿态角；

11.不断重复步骤7-10，实现对三轴气浮平台轴承的姿态角信息实时测量。

本发明所述三轴气浮台姿态角测量过程，如附图5所示，涉及到以下坐标系：

(1)世界坐标系O_wx_wy_wz_w，该坐标系由测量标靶建立，为静止坐标系；

(2)气浮球轴承的本体坐标系O_bx_by_bz_b，原点O_b位于气浮球球心，为静止坐标系；

(3)摄像机1成像面的图像像素坐标系O₁u₁v₁；

(4)摄像机2成像面的图像像素坐标系O₂u₂v₂；

(5)摄像机坐标系O_cx_cy_cz_c，O_c为光心，z_c轴为光轴。

实施例3

1、投影矩阵计算原理：

本发明所涉及摄像机的图像像素坐标系与世界坐标系之间的投影矩阵H的计算方法，如下所示：

测量标靶上第i个标志点与其像点之间的成像模型为：

$z_{\mathrm{ci}}\left[\begin{array}{c}{u}_i\\ v_i\\ 1\end{array}\right]=H\left[\begin{array}{c}{x}_{w,i}\\ y_{w,i}\\ z_{w,i}\\ 1\end{array}\right]=\left[\begin{array}{cccc}{h}_{11}& h_{12}& h_{13}& h_{14}\\ h_{21}& h_{22}& h_{23}& h_{24}\\ h_{31}& h_{32}& h_{33}& h_{34}\end{array}\right]\left[\begin{array}{c}{x}_{w,i}\\ y_{w,i}\\ z_{w,i}\\ 1\end{array}\right],---\left(1\right)$

其中，[x_w，i y_w，i z_w，i]是第i个标志点的世界坐标系坐标，[u_i v_i]是第i个标志点的图像像素坐标系坐标。

化简为：

$\left[\begin{array}{ccccccccccc}{x}_{w,i}& y_{w,i}& z_{w,i}& 1& 0& 0& 0& 0& -x_w{u}_i& -y_w{u}_i& -z_w{u}_i\\ 0& 0& 0& 0& x_{w,i}& y_{w,i}& z_{w,i}& 1& -x_w{v}_i& -y_w{v}_i& -z_w{v}_i\end{array}\right]\left[\begin{array}{c}{h}_{11}\\ h_{12}\\ h_{13}\\ h_{14}\\ h_{21}\\ h_{22}\\ h_{23}\\ h_{24}\\ h_{31}\\ h_{32}\\ h_{33}\end{array}\right]=\left[\begin{array}{c}{u}_i{h}_{34}\\ v_i{h}_{34}\end{array}\right]---\left(2\right)$

对于每一个坐标点，都可以建立上述两个方程。因此取6个坐标点，可以得到12个方程，表示成矩阵形式为：

KX＝U    (3)

其中

$K=\left[\begin{array}{ccccccccccc}{x}_{w,1}& y_{w,1}& z_{w,1}& 1& 0& 0& 0& 0& -x_{w,1}{u}_1& -y_{w,1}{u}_1& -z_{w,1}{u}_1\\ 0& 0& 0& 0& x_{w,1}& y_{w,1}& y_{w,1}& 1& -x_{w,1}{v}_1& -y_{w,1}{v}_1& -z_{w,1}{v}_1\\ x_{w,2}& y_{w,2}& z_{w,2}& 1& 0& 0& 0& 0& -x_{w,2}{u}_2& -y_{w,2}{u}_2& -z_{w,2}{u}_2\\ 0& 0& 0& 0& x_{w,2}& y_{w,2}& y_{w,2}& 1& -x_{w,2}{v}_2& -y_{w,2}{v}_2& -z_{w,2}{v}_2\\ x_{w,3}& y_{w,3}& z_{w,3}& 1& 0& 0& 0& 0& -x_{w,3}{u}_3& -y_{w,3}{u}_3& -z_{w,3}{u}_3\\ 0& 0& 0& 0& x_{w,3}& y_{w,3}& z_{w,3}& 1& -x_{w,3}{v}_3& -y_{w,3}{v}_3& -z_{w,3}{v}_3\\ x_{w,4}& y_{w,4}& z_{w,4}& 1& 0& 0& 0& 0& -x_{w,4}{u}_4& -y_{w,4}{u}_4& -z_{w,4}{u}_4\\ 0& 0& 0& 0& x_{w,4}& y_{w,4}& z_{w,4}& 1& -x_{w,4}{v}_4& -y_{w,4}{v}_4& -z_{w,4}{v}_4\\ x_{w,5}& y_{w,5}& z_{w,5}& 1& 0& 0& 0& 0& -x_{w,5}{u}_5& -y_{w,5}{u}_5& -z_{w,5}{u}_5\\ 0& 0& 0& 0& x_{w,5}& y_{w,5}& z_{w,5}& 1& -x_{w,5}{v}_5& -y_{w,5}{v}_5& -z_{w,5}{v}_5\\ x_{w,6}& y_{w,6}& z_{w,6}& 1& 0& 0& 0& 0& -x_{w,6}{u}_6& -y_{w,6}{u}_6& -z_{w,6}{u}_6\\ 0& 0& 0& 0& x_{w,6}& y_{w,6}& z_{w,6}& 1& -x_{w,6}{v}_6& -y_{w,6}{v}_6& -z_{w,6}{v}_6\end{array}\right]$

$X=\frac{1}{h_{34}}{\left[\begin{array}{ccccccccccc}{h}_{11}& h_{12}& h_{13}& h_{14}& h_{21}& h_{22}& h_{23}& h_{24}& h_{31}& h_{32}& h_{33}\end{array}\right]}^T$

U＝[u₁ v₁ u₂ v₂ u₃ v₃ u₄ v₄ u₅ v₅ u₆ v₆]^T。

采用最小二乘法可得：

X＝(K^TK)^-1K^TU    (4)

求出X＝[X₁ X₂ X₃ X₄ X₅ X₆ X₇ X₈ X₉ X₁₀ X₁₁]后，可得

$h_{34}=\frac{1}{\left|\begin{array}{ccc}{h}_{31}& h_{32}& h_{33}\end{array}\right|},$

h₁₁＝h₃₄X₁，h₁₂＝h₃₄X₂，h₁₃＝h₃₄X₃，h₁₄＝h₃₄X₄，  (5)

h₂₁＝h₃₄X₅，h₂₂＝h₃₄X₆，h₂₃＝h₃₄X₇，h₂₄＝h₃₄X₈，

h₃₁＝h₃₄X₉，h₃₂＝h₃₄X₁₀，h₃₃＝h₃₄X₁₁。

至此，解得投影矩阵H。

2、双目立体视觉测量原理：

假定空间任意点P在两台摄像机101和102上所成的像点p₁、p₂已经提取出来，假设摄像机101和摄像机102的投影矩阵已知，分别为H₁与H₂，则有：

$z_1\left[\begin{array}{c}{u}_1\\ v_1\\ 1\end{array}\right]=H_1\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\\ 1\end{array}\right]=\left[\begin{array}{cccc}{h}_{11}^1& h_{12}^1& h_{13}^1& h_{14}^1\\ h_{21}^1& h_{22}^1& h_{23}^1& h_{24}^1\\ h_{31}^1& h_{32}^1& h_{33}^1& h_{34}^1\end{array}\right]\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\\ 1\end{array}\right]---\left(6\right)$

$z_2\left[\begin{array}{c}{u}_2\\ v_2\\ 1\end{array}\right]=H_2\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\\ 1\end{array}\right]=\left[\begin{array}{cccc}{h}_{11}^2& h_{12}^2& h_{13}^2& h_{14}^2\\ h_{21}^2& h_{22}^2& h_{23}^2& h_{24}^2\\ h_{31}^2& h_{32}^2& h_{33}^2& h_{34}^2\end{array}\right]\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\\ 1\end{array}\right]---\left(7\right)$

式中，(u₁，v₁，1)与(u₂，v₂，1)分别为像点p₁、p₂在成像面中的像素坐标，(X_w，Y_w，Z_w，1)^T为P点在世界坐标下的齐次坐标；为H_k(k＝1，2)的第i行第j列元素。简化公式(6)和(7)，消去z₁、z₂后，可以得到关于X_w，Y_w，Z_w的四个线性方程，即

$\begin{array}{c}(u_1{h}_{31}^1-h_{11}^1)X_w+(u_1{h}_{32}^1-h_{12}^1)Y_w+(u_1{h}_{33}^1-h_{13}^1)Z_w=h_{14}^1-u_1{h}_{34}^1\\ (v_1{h}_{31}^1-h_{21}^1)X_w+(v_1{h}_{32}^1-h_{22}^1)Y_w+(v_1{h}_{33}^1-h_{23}^1)Z_w=h_{24}^1-v_1{h}_{34}^1\\ (u_2{h}_{31}^2-h_{11}^2)X_w+(u_2{h}_{32}^2-h_{12}^2)Y_w+(u_2{h}_{33}^2-h_{13}^2)Z_w=h_{14}^2-u_2{h}_{34}^2\\ (v_2{h}_{31}^2-h_{21}^2)X_w+(v_2{h}_{32}^2-h_{22}^2)Y_w+(v_2{h}_{33}^2-h_{23}^2)Z_w=h_{24}^2-v_2{h}_{34}^2\end{array}---\left(8\right)$

将两组方程联系可以求出点P坐标X_w，Y_w，Z_w。事实上，包含X_w，Y_w，Z_w三个变量的四个线性方程，可以用最小二乘法求解惟一值X_w，Y_w，Z_w。

VX′＝D    (9)

其中，

$V=\left[\begin{array}{ccc}{u}_1{h}_{31}^1-h_{11}^1& u_1{h}_{32}^1-h_{12}^1& u_1{h}_{33}^1-h_{13}^1\\ v_1{h}_{31}^1-h_{21}^1& v_1{h}_{32}^1-h_{22}^1& v_1{h}_{33}^1-h_{23}^1\\ u_2{h}_{31}^2-h_{11}^2& u_2{h}_{32}^2-h_{12}^2& u_2{h}_{33}^2-h_{13}^2\\ v_2{h}_{31}^2-h_{21}^2& v_2{h}_{32}^2-h_{22}^2& v_2{h}_{33}^2-h_{23}^2\end{array}\right],D=\left[\begin{array}{c}{h}_{14}^1-u_1{h}_{34}^1\\ h_{24}^1-v_1{h}_{34}^1\\ h_{14}^2-u_2{h}_{34}^2\\ h_{24}^2-v_2{h}_{34}^2\end{array}\right]$ 为已知量

求解得：

X′＝(V^TV)^-1V^TD    (10)

3、气浮球轴承在本体坐标系下的旋转矩阵计算原理：

世界坐标系与气浮轴承的本体坐标系之间的几何变换关系为：

$\left[\begin{array}{c}{X}_b\\ Y_b\\ Z_b\end{array}\right]=R_0\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\end{array}\right]+T_0---\left(11\right)$

其中[X_b Y_b Z_b]^T为气浮球轴承的本体坐标系，R₀为旋转矩阵，T₀为平移向量，这两个量可以预先测得。

由公式(11)可知：第k时刻，标志点在气浮轴承本体坐标系中的坐标可以表示为

$\left[\begin{array}{c}{X}_{b,k}\\ Y_{b,k}\\ Z_{b,k}\end{array}\right]=R_0\left[\begin{array}{c}{X}_{w,k}\\ Y_{w,k}\\ Z_{w,k}\end{array}\right]+T_0---\left(12\right)$

第k-1时刻，标志点在气浮轴承本体坐标系中的坐标可以表示为

$\left[\begin{array}{c}{X}_{b,k-1}\\ Y_{b,k-1}\\ Z_{b,k-1}\end{array}\right]=R_0\left[\begin{array}{c}{X}_{w,k-1}\\ Y_{w,k-1}\\ Z_{w,k-1}\end{array}\right]+T_0---\left(13\right)$

其中，[X_w，k Y_w，k Z_w，k]^T为第k时刻标志点在世界坐标系中的坐标，

[X_w，k-1 Y_w，k-1 Z_w，k-1]^T为第k-1时刻标志点在世界坐标系中的坐标；

[X_b，k Y_b，k Z_b，k]^T为第k时刻标志点在本体坐标系中的坐标，

[X_b，k-1 Y_b，k-1 Z_b，k-1]^T为第k-1时刻标志点在本体坐标系中的坐标。

在利用双目视觉测量原理，计算出标志点在第[X_b，K Y_b，K Z_b，K]^T和

[X_w，K-1 Y_W，K-1 Z_w，K-1]^T之后，[X_b，K Y_b，K Z_b，K]^T与[X_b，K-1 Y_b，K-1 Z_b，K-1]^T就可以由公式(12)、(13)随之计算出来。

而[X_b，K Y_b，K Z_b，K]^T与[X_b，K-1 Y_b，K-1 Z_b，K-1]^T之间存在以下关系：

$\left[\begin{array}{c}{X}_{b,K}\\ Y_{b,K}\\ Z_{b,K}\end{array}\right]=R\left[\begin{array}{c}{X}_{b,K-1}\\ Y_{b,K-1}\\ Z_{b,K-1}\end{array}\right],R=\left[\begin{array}{ccc}{r}_{11}& r_{12}& r_{13}\\ r_{21}& r_{22}& r_{23}\\ r_{31}& r_{32}& r_{33}\end{array}\right]---\left(14\right)$

其中，R就是三轴气浮球轴承在其本体坐标系下的旋转矩阵，也就是我们求解的目标。

要求解R，首先将式(14)变换如下式：

X_b，k-1r₁₁+Y_b，k-1r₁₂+Z_b，k-1r₁₃＝X_b，k

X_b，k-1r₂₁+Y_b，k-1r₂₂+Z_b，k-1r₂₃＝Y_b，k     (15)

X_b，k-1r₃₁+Y_b，k-1r₃₂+Z_b，k-1r₃₃＝Z_b，k

再将公式(15)转换为矩阵形式，如：

$\left[\begin{array}{ccccccccc}{X}_{b,k-1}& Y_{b,k-1}& Z_{b,k-1}& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& X_{b,k-1}& Y_{b,k-1}& Z_{b,k-1}& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& X_{b,k-1}& Y_{b,k-1}& Z_{b,k-1}\end{array}\right]\left[\begin{array}{c}{r}_{11}\\ r_{12}\\ r_{13}\\ r_{21}\\ r_{22}\\ r_{23}\\ r_{31}\\ r_{32}\\ r_{33}\end{array}\right]=\left[\begin{array}{c}{X}_{b,k}\\ Y_{b,k}\\ Z_{b,k}\end{array}\right]---\left(16\right)$

则由测量标靶上每一个标志点可以建立一组上述方程，当有3个标志点时，可建立9个方程式，如：

QR′＝B    (17)

其中，

$Q=\left[\begin{array}{ccccccccc}{X^1}_{b,k-1}& {Y^1}_{b,k-1}& {Z^1}_{b,k-1}& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& {X^1}_{b,k-1}& {Y^1}_{b,k-1}& {Z^1}_{b,k-1}& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& {X^1}_{b,k-1}& {Y^1}_{b,k-1}& {Z^1}_{b,k-1}\\ {X^2}_{b,k-1}& {Y^2}_{b,k-1}& {Z^2}_{b,k-1}& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& {X^2}_{b,k-1}& {Y^2}_{b,k-1}& {Z^2}_{b,k-1}& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& {X^2}_{b,k-1}& {Y^2}_{b,k-1}& {Z^2}_{b,k-1}\\ {X^3}_{b,k-1}& {Y^3}_{b,k-1}& {Z^3}_{b,k-1}& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& {X^3}_{b,k-1}& {Y^3}_{b,k-1}& {Z^3}_{b,k-1}& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& {X^3}_{b,k-1}& {Y^3}_{b,k-1}& {Z^3}_{b,k-1}\end{array}\right]$

R′＝[r₁₁ r₁₂ r₁₃ r₂₁ r₂₂ r₂₃ r₃₁ r₃₂ r₃₃]^T

B＝[X_b，k T_b，k Z_b，k]^T

根据最小二乘法思想，可求得矩阵R′＝(Q^TQ)^-1Q^TB，从而解得旋转矩阵R。

4、气浮球轴承在本体坐标系下的三维姿态角解算原理

绕本体坐标系中的X_b轴旋转某一角度α的变换矩阵是：

$R_X\left(\mathrm{\α}\right)=\left[\begin{array}{ccc}1& 0& 0\\ 0& \cos \mathrm{\α}& -\sin \mathrm{\α}\\ 0& \sin \mathrm{\α}& \cos \mathrm{\α}\end{array}\right]---\left(18\right)$

绕本体坐标系中的Y_b轴旋转某一角度β的变换矩阵是：

$R_Y\left(\mathrm{\β}\right)=\left[\begin{array}{ccc}\cos \mathrm{\β}& 0& \sin \mathrm{\β}\\ 0& 1& 0\\ -\sin \mathrm{\β}& & \cos \mathrm{\β}\end{array}\right]---\left(19\right)$

绕本体坐标系中的Z_b轴旋转某一角度γ的变换矩阵是：

$R_Z\left(\mathrm{\γ}\right)=\left[\begin{array}{ccc}\cos \mathrm{\γ}& -\sin \mathrm{\γ}& 0\\ \sin \mathrm{\γ}& \cos \mathrm{\γ}& 0\\ 0& 0& 1\end{array}\right]---\left(20\right)$

采用Z_b→Y_b→X_b轴向的旋转顺序(约定使用右手坐标系，逆时针转动得到的角度为正)。则对应的旋转矩阵以及由相应的旋转矩阵求解欧拉角的公式如下：

$\begin{array}{c}R(\mathrm{\γ},\mathrm{\β},\mathrm{\α})=R_Z\left(\mathrm{\γ}\right)R_Y\left(\mathrm{\β}\right)R_X\left(\mathrm{\α}\right)\\ =\left[\begin{array}{ccc}\cos \mathrm{\γ}& -\sin \mathrm{\γ}& 0\\ \sin \mathrm{\γ}& \cos \mathrm{\γ}& 0\\ 0& 0& 1\end{array}\right]\left[\begin{array}{ccc}\cos \mathrm{\β}& 0& \sin \mathrm{\β}\\ 0& 1& 0\\ -\sin \mathrm{\β}& 0& \cos \mathrm{\β}\end{array}\right]\left[\begin{array}{ccc}1& 0& 0\\ 0& \cos \mathrm{\α}& -\sin \mathrm{\α}\\ 0& \sin \mathrm{\α}& \cos \mathrm{\α}\end{array}\right]\\ =\left[\begin{array}{ccc}\cos \mathrm{\γ}\cos \mathrm{\β}& \cos \mathrm{\γ}\sin \mathrm{\β}\sin \mathrm{\α}-\sin \mathrm{\γ}\cos \mathrm{\α}& \cos \mathrm{\γ}\sin \mathrm{\β}\cos \mathrm{\α}+\sin \mathrm{\γ}\sin \mathrm{\α}\\ \sin \mathrm{\γ}\cos \mathrm{\β}& \sin \mathrm{\γ}\sin \mathrm{\β}\sin \mathrm{\α}+\cos \mathrm{\γ}\cos \mathrm{\α}& \sin \mathrm{\γ}\sin \mathrm{\β}\cos \mathrm{\α}-\cos \mathrm{\γ}\sin \mathrm{\α}\\ -\sin \mathrm{\β}& \cos \mathrm{\β}\sin \mathrm{\α}& \cos \mathrm{\β}\cos \mathrm{\α}\end{array}\right]\end{array}---\left(21\right)$

由式(17)和(21)可以算出Z_b→Y_b→X_b欧拉角旋转矩阵为：

$R(\mathrm{\γ},\mathrm{\β},\mathrm{\α})=\left[\begin{array}{ccc}{r}_{11}& r_{12}& r_{13}\\ r_{21}& r_{22}& r_{23}\\ r_{31}& r_{32}& r_{33}\end{array}\right]---\left(22\right)$

可得各欧拉角的计算式如下：

$\left\{\begin{array}{c}\mathrm{\β}=-\arctan \frac{r_{31}}{\sqrt{{r_{11}}^2+{r_{21}}^2}}\\ \mathrm{\α}=\arctan \frac{r_{32}}{r_{33}}\\ \mathrm{\γ}=\arctan \frac{r_{21}}{r_{11}}\end{array}\right.---\left(23\right)$

Claims (1)

1.一种三轴气浮台姿态角双目立体视觉测量方法，该方法采用一种三轴气浮台姿态角双目立体视觉测量装置，包括三轴气浮台，安装在三轴气浮台上的气浮轴承台面能够浮起进行三维转动，并且还包括2台摄像机及其基座、立体测量标靶、辅助光源和计算机，在气浮轴承台面上方安装两台黑白数字CCD摄像机，摄像机安装在基座上，辅助光源安装在摄像机附近，立体测量标靶安装在气浮轴承台面上，立体测量标靶上设置有多个人工标志点，数字CCD摄像机连续采集人工标志点的图像并传输到计算机，计算机将图像信息进行分析与处理；所述的数字CCD摄像机通过数据线与转换器相连，转换器与计算机相连，数字CCD摄像机的镜头为物方远心镜头，其特征在于，方法如下：

(1)、将立体测量标靶安装在气浮轴承台面上，建立世界坐标系；

(2)、立体测量标靶上有多个标志点，设置在气浮轴承台面上方的两台计算机对立体测量标靶的图像信息进行分析与处理，针对这些标志点完成特征提取和亚像素定位；

(3)、利用标志点的定位信息，根据计算机视觉成像原理，计算两台摄像机各自的图像像素坐标系与世界坐标系之间的投影矩阵H₁和H₂；

(4)、两台摄像机实时、连续地采集立体测量标靶上标志点的图像，并传输至计算机；

(5)、计算机对来自两台摄像机的两幅图像分别进行特征提取，利用亚像素定位方法确定立体测量标靶的标志点分别在两幅图像中的图像像素坐标系中的坐标；

(6)、根据步骤(5)的结果，结合步骤(3)求得的投影矩阵H₁和H₂，利用双目视觉测量原理，计算气浮轴承在世界坐标系中的三维坐标；

(7)、根据气浮球轴承的坐标转换关系，解算其三维姿态角信息。