CN103542981A - 一种双目视觉测量转动惯量的方法 - Google Patents

Abstract

本发明一种双目视觉测量转动惯量的方法属于计算机视觉测量技术领域，涉及一种基于双目视觉的快速、精确转动惯量测量方法。该测量方法采用双目视觉系统，利用两台超高速摄像机实时采集被测物体表面标记点的运动图像，由图形工作站将两个超高速摄像机传输来的图片信息进行处理，得知标记点在空间运动的轨迹，通过对轨迹的拟合计算出被测物体的转动惯量。本发明利用双目视觉系统将被测物每一时刻的空间位置进行平面拟合与投影，使得被测物空间运动坐标在同一平面内，减小摆动运动不在一个平面给测量精度带来的影响；并通过考虑阻尼系数的被测物实时空间运动坐标曲线拟合，减小空气阻力以及机构阻力对测量系统精度的影响，提高了转动惯量测量精度。

Description

一种双目视觉测量转动惯量的方法

技术领域

本发明属于计算机视觉测量技术领域，涉及一种基于双目视觉的快速、精确转动惯量测量方法。

背景技术

转动惯量是表征转动物体惯性大小的物理量，它是研究转动物体运动规律的重要参数之一，涉及到运动物体测量以及运动动力学的问题，都需要精确知道转动物体的转动惯量。随着现代科技的不断进步，转动惯量的测量越来越受到重视，其已经成为现代测量领域不可缺少的重要组成，其在各领域中都有着十分重要的地位。在航天工业中，人造卫星、运载火箭、载人飞船等都需要对转动惯量进行，以确定产品是否符合设计要求加以修正；转动惯量直接影响飞行器飞行地稳定性，因此测量转动惯量为研究和设计飞行器提供了重要的参数，同时也可为形状、飞行姿态的优化提供信息。在汽车工业，各种车辆以及转动部件必须测量其转动惯量，通过修正偏心来提高车辆的性能和寿命，因此测定物体的转动惯量具有重要的实际意义。

专利号为CN1646971A的名为《一种转动惯量和惯性积测量方法及其装置》，该专利采用复摆装置，但是该复摆装置操作繁琐。专利号为CN102692264A的《一种用于质量、质心位置与转动惯量的测试台及测试方法》采用扭摆机构和安装在机构上的传感器来测量转动惯量，传感器需进行测前标定，且扭摆装置本身操作复杂，使得整个测量过程十分繁琐。

发明内容

本发明要解决的技术难题是消减转动惯量测量时摆动轨迹不共面与摆动时机构、空气阻力带来的转动惯量测量误差，发明一种基于双目视觉的转动惯量测量方法，该方法采用复摆装置进行双目视觉测量转动惯量。利用双目视觉系统可以得知每个时刻被测物体在空间摆动的轨迹，可将被测物每一时刻的空间位置进行平面拟合与投影，使得被测物空间运动坐标在同一平面内，减小摆动运动不在一个平面给测量精度带来的影响；并且通过考虑阻尼系数的被测物实时空间运动坐标曲线拟合，减小空气阻力以及机构阻力对测量系统精度的影响，增加了转动惯量测量精度。

本发明采用的技术方案是一种双目视觉测量转动惯量的方法，其特征是，本方法通过左、右两台超高速摄像机8、8’采集被测物6表面标标记点的位置信息，经标记点提取、匹配与重建得到序列图像的标记点圆心坐标，再进行空间曲面拟合投影，最终计算转动惯量。该方法的具体测量步骤如下：

（1）两个高速摄像机的标定

采用张氏标定方法标定出两相机的内外参数；

$s\left[\begin{array}{c}u\\ v\\ 1\end{array}\right]=\left[\begin{array}{cccc}{\mathrm{\α}}_x& 0& u_0& 0\\ 0& {\mathrm{\α}}_y& v_0& 0\\ 0& 0& 1& 0\end{array}\right]\left[\begin{array}{cc}R& t\\ 0^T& 1\end{array}\right]\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\\ 1\end{array}\right]---\left(1\right)$

其中s为比例因子，u、v为相面坐标，X_w、Y_w、Z_w是世界坐标，α_xα_x、α_y、u₀、v₀为摄像机内参数，分别代表u轴尺度因子、v轴尺度因子、相面中心u轴坐标、相面中心v轴坐标；R为旋转矩阵、t为平移矩阵，它们是摄像机的外参数。

（2）标记点特征的提取

对采集到的视频每一帧进行处理，每一帧只有一个有效标记点，采用梯度重心法提取被被测物（6）运动时表面的标记点，首先采用高斯一阶微分算子对图像进行卷积操作，以获得图像在各点的梯度，则对于标志图案所在子区域的梯度重心点坐标的求取方法如下：

$C(x,y)=\underset{i=-h}{\overset{h}{\mathrm{\Σ}}}\underset{j=-w}{\overset{w}{\mathrm{\Σ}}}\left[\right|G(i,j)|\·P(i,j)]/\underset{i=-h}{\overset{h}{\mathrm{\Σ}}}\underset{j=-w}{\overset{w}{\mathrm{\Σ}}}\left|G(i,j)\right|---\left(2\right)$

其中：C（x,y）为应用梯度重心法所提取的标记点中心像素坐标值，|G(i,j)|为(i,j)点的梯度幅值，w，h分别为目标图像的宽和高，P(i,j)为(i,j)点的图像坐标。

（3）标记点的匹配与重建

对采集到的视频每一帧所提取的标记点进行匹配与重建，图像标记点的匹配方法如下：先利用Longguet-Higgins提出的归一化8点算法计算出左右两高速摄像机（8、8’）的基本矩阵F，再通过左高速摄像机(8)拍摄图片与右高速摄像机(8’)拍摄图片的极限约束关系进行图像标记点的初匹配，若左高速摄像机(8)拍摄的图像标记点x_l与右高速摄像机(8’)拍摄的图像标记点x_r相匹配，即满足极限约束条件，极限约束条件如下所示：

$x_l^{T}F{x}_r=0---\left(3\right)$

其中：x_l为左高速摄像机(8)拍摄的图像标记点；x_r为右高速摄像机(8’)拍摄的与x_l匹配的图像标记点；F为左右两个高速摄像机（8、8’）间的基本矩阵；

然后将所有满足极限约束条件的左右图像标记点进行三维重建，重建出空间标记点在世界坐标系下的三维坐标值，其重建公式如下所示：

$x=\frac{z{X}_1}{f_1}$

$y=\frac{{\mathrm{zY}}_1}{f_1}---\left(4\right)$

$\frac{f_1(f_2{t}_y-y_2{t}_z)}{Y_1(r_7{X}_1+r_8{Y}_1+r_9{f}_1)f_2(r_4{X}_1+r_5{Y}_1+r_6{f}_1)}$

其中：x_l=[X₁ Y]₁，X₁,Y₁分别为左高速摄像机(8)拍摄的图像标记点x_l的横、纵坐标；x_r=[X₂ Y₂]，X₂,Y₂分别为右高速摄像机(8’)拍摄的图像标记点x_r的横、纵坐标；(x y z)为由左两图像标记点x_l和右两图像标记点x_r重建出来的空间标记点的三维坐标；f₁、f₂分别为左右高速摄像机(8、8’)的焦距；

$\left[\begin{array}{ccc}{r}_1& r_2& r_3\\ r_4& r_5& r_6\\ r_7& r_8& r_9\end{array}\right]$ 为右高速摄像机（8’）相对于左高速摄像机(8)的旋转矩阵，

[t_x t_y t_z]是右高速摄像机(8’)相对于左高速相机(8)的平移矩阵；

（4）转动惯量计算

●空间标记点平面拟合

将视频每一帧重建出的空间点坐标在同一坐标系下进行显示，并且将这些序列点进行平面拟合。平面拟合公式如下所示：

空间平面方程的通常可采用公式为：

Ax+By+Cz+1＝0(C≠0)          (5)

对于视频第i帧提取的空间坐标点，即i帧空间点坐标为：(x_i,y_i,z_i)i＝1,2,Ln(n≥3)，所要拟合的平面可表达为如下的矩阵形式：

$\left[\begin{array}{ccc}{x}_1& y_1& z_1\\ M& M& M\\ x_n& y_n& z_n\end{array}\right]\left[\begin{array}{c}A\\ B\\ C\end{array}\right]=\left[\begin{array}{c}-1\\ -1\\ -1\end{array}\right]---\left(6\right)$

根据最小二乘法法拟合空间坐标点所在平面，将矩阵左乘 ${\left[\begin{array}{ccc}{x}_1& y_1& z_1\\ M& M& M\\ x_n& y_n& z_n\end{array}\right]}^T,$ 式（6）化简为：

$\left[\begin{array}{ccc}\mathrm{\Σ}{x_i}^2& \mathrm{\Σ}{x}_i{y}_i& \mathrm{\Σ}{x}_i{z}_i\\ \mathrm{\Σ}{x}_i{y}_i& \mathrm{\Σ}{y_i}^2& \mathrm{\Σ}{y}_i{z}_i\\ \mathrm{\Σ}{x}_i{z}_i& \mathrm{\Σ}{y}_i{z}_i& \mathrm{\Σ}{z_i}^2\end{array}\right]\left[\begin{array}{c}A\\ B\\ C\end{array}\right]=\left[\begin{array}{c}-\mathrm{\Σ}{x}_i\\ -\mathrm{\Σ}{y}_i\\ -\mathrm{\Σ}{z}_i\end{array}\right]---\left(7\right)$

可推导得：

$\left[\begin{array}{c}A\\ B\\ C\end{array}\right]={\left[\begin{array}{ccc}\mathrm{\Σ}{x_i}^2& \mathrm{\Σ}{x}_i{y}_i& \mathrm{\Σ}{x}_i{z}_i\\ \mathrm{\Σ}{x}_i{y}_i& \mathrm{\Σ}{y_i}^2& \mathrm{\Σ}{y}_i{z}_i\\ \mathrm{\Σ}{x}_i{z}_i& \mathrm{\Σ}{y}_i{z}_i& \mathrm{\Σ}{y}_i{z}_i\end{array}\right]}^{-1}\left[\begin{array}{c}-\mathrm{\Σ}{x}_i\\ -\mathrm{\Σ}{y}_i\\ -\mathrm{\Σ}{z}_i\end{array}\right]---\left(8\right)$

即所得平面为拟合平面。

●空间标记点平面投影

将视频每帧提取的空间坐标点P_i＝(x_i,y_i,z_i)i=1,2,Ln(n≥3)投影到根据公式（8）拟合的平面，空间平面的单位法向量为：

$N_P=\left[\begin{array}{ccc}\frac{1}{B*C}& \frac{1}{A*C}& \frac{1}{A*B}\end{array}\right]---\left(9\right)$

且平面过P_p点(1 1(1+A+B)/(-C))。则空间坐标点在空间平面的投影点坐标为：

P＝P_i-N_p×((P_i-P_p)×N_p)      (10)

P为P_i投影到由公式（10）拟合的平面上的投影点。

●被测标记点轨迹拟合与摆动周期计算

建立P中x轴空间坐标与时间的曲线，并且利用公式（12）进行曲线拟合。

x＝a×e^-b×t×sin(c×t+d)+e           (11)

其中a,d,e为积分常数，由初始条件决定。b为阻尼因子，

ω₀为系统固有频率，然后根据公式（12）计算得出被测物的摆动周期。

$T=\frac{2\mathrm{\π}}{c}---\left(12\right)$

其中T为被测物体的摆动周期。

●被测物体转动惯量计算

$I=\frac{T^2\mathrm{mgl}}{4{\mathrm{\π}}^2}---\left(13\right)$

其中,r为模型做单摆运动的周期,I为其转动惯量,m为模型质量,g为重力加速度,l为模型质心距转轴的距离。通过对模拟目标进行视频拍摄,根据目标标志的运动周期推算其做单摆运动的周期,可由上式对模型的转动惯量进行求解。

本发明的有益效果是本发明利用双目视觉系统将被测物每一时刻的空间位置进行平面拟合与投影，使得被测物空间运动坐标在同一平面内，减小摆动运动不在一个平面给测量精度带来的影响；并通过考虑阻尼系数的被测物实时空间运动坐标曲线拟合，减小空气阻力以及机构阻力对测量系统精度的影响，提高了转动惯量测量精度。

附图说明

图1为基于双目视觉转动惯量的测量装置图。其中：1转动惯量支架、2轴承、3螺母、4复摆、5夹具、6被测物、7左方形光源、7’右方形光源、8左超高速摄像机、8’右超高速摄像机、9左电控平台、9’右电控平台、10气浮平台、11图形工作站。

图2基于双目视觉的转动惯量测量方法流程图。

图3被测物体空间摆动轨迹，其中，x、y、z为空间三个坐标轴。

图4被测物体摆动轨迹在X轴坐标与时间轴的拟合曲线。

具体实施方式

以下结合技术方案和附图详细叙述本发明的具体实施方式。附图1为基于双目视觉转动惯量的测量装置图。本方法通过左、右两台超高速摄像机8、8’采集被测物6表面标标记点的位置信息，经标记点提取、匹配与重建得到序列图像的标记点圆心坐标，再进行空间曲面拟合投影，最终计算转动惯量。

先安装测量装置，左、右两个电控平台9、9’通过螺栓固定在气浮平台10上，左、右两台超高速摄像机8、8’利用螺栓分别安装在左、右两个电控平台9、9’上，左、右两个方形光源7、7’使用螺栓固定于左、右两台超高速摄像机8、8’镜头周边，用传输线连接左、右两台超高速摄像机8、8’与图形工作站11。将轴承2安装进复摆4的轴承孔内，利用螺母3把轴承2固定在转动惯量支架1上，通过夹具5把被测物6固定在复摆4上。将安装好的转动惯量支架1放在左、右两台超高速这相机8、8’前，使得被测物6上的标记点清晰可见。

实施例1，本发明采用分别带有广角镜头的两个高速摄像机8、8’拍摄物体运动情况，两个超高速摄像机型号为FASTCAM SAX摄像机，分辨率：1024×1024，CCD面积：2cm×2cm，帧频：保证分辨率1024×1024情况下最多可达10000fps，重量：12.6kg。广角镜头型号为AF-S17-35mm f/2.8D IF-ED参数如下所示，镜头焦距：f=17-35，APS焦距：25.5-52.5，最大光圈：F2.8，最小光圈：F22，镜头重量：745g，镜头尺寸：82.5×106。拍摄条件如下：高速摄像机帧频为3000fps，图片像素为1024×1024，镜头焦距为17mm，物距为750mm，视场约为800mm×800mm。

首先通过图形工作站11控制左、右两台电控平台9、9’移动，调整左、右两台超高速摄像机8、8’的测量位置，然后打开左、右两个方形光源7、7’，以提高测量空间的亮度，最后由图形工作站11进行双目图像标定、双目图像的特征点提取、滚转体特征点的匹配、滚转体位姿测量等工作。

附图2是基于双目视觉的转动惯量测量方法流程图，测量方法主要步骤为两台超高速摄像机8、8’的标定，被测物6装夹及标记点黏贴，标记点提取，

标记点匹配与重建，空间点平面拟合与投影、转动惯量求解。

（1）两个高速摄像机的标定

采用张氏标定方法标定出两相机的内外参数；

$s\left[\begin{array}{c}u\\ v\\ 1\end{array}\right]=\left[\begin{array}{cccc}{\mathrm{\α}}_x& 0& u_0& 0\\ 0& {\mathrm{\α}}_y& v_0& 0\\ 0& 0& 1& 0\end{array}\right]\left[\begin{array}{cc}R& t\\ 0^T& 1\end{array}\right]\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\\ 1\end{array}\right]---\left(1\right)$

其中s为比例因子，u、v为相面坐标，X_w、Y_w、Z_w是世界坐标，α_xα_x、α_y、u₀、v₀为摄像机内参数，分别代表u轴尺度因子、v轴尺度因子、相面中心u轴坐标、相面中心v轴坐标；R为旋转矩阵、t为平移矩阵，它们是摄像机的外参数。

（2）标记点特征的提取

对采集到的视频每一帧进行处理，每一帧只有一个有效标记点，采用梯度重心法提取被被测物（6）运动时表面的标记点，首先采用高斯一阶微分算子对图像进行卷积操作，以获得图像在各点的梯度，则对于标志图案所在子区域的梯度重心点坐标的求取方法如下：

$C(x,y)=\underset{i=-h}{\overset{h}{\mathrm{\Σ}}}\underset{j=-w}{\overset{w}{\mathrm{\Σ}}}\left[\right|G(i,j)|\·P(i,j)]/\underset{i=-h}{\overset{h}{\mathrm{\Σ}}}\underset{j=-w}{\overset{w}{\mathrm{\Σ}}}\left|G(i,j)\right|---\left(2\right)$

其中：C（x,y）为应用梯度重心法所提取的标记点中心像素坐标值，|G(i,j)|为(i,j)点的梯度幅值，w，h分别为目标图像的宽和高，P(i,j)为(i,j)点的图像坐标。

（3）标记点的匹配与重建

对采集到的视频每一帧所提取的标记点进行匹配与重建，图像标记点的匹配方法如下：先利用Longguet-Higgins提出的归一化8点算法计算出左右两高速摄像机（8、8’）的基本矩阵F，再通过左高速摄像机（8）拍摄图片与右高速摄像机（8’）拍摄图片的极限约束关系进行图像标记点的初匹配，若左高速摄像机(8)拍摄的图像标记点x_l与右高速摄像机(8’)拍摄的图像标记点x_r相匹配，

即满足极限约束条件，极限约束条件如下所示：

$x_l^{T}F{x}_r=0---\left(3\right)$

其中：x_l为左高速摄像机(8)拍摄的图像标记点；x_r为右高速摄像机(8’)拍摄的与x_l匹配的图像标记点；F为左右两个高速摄像机（8、8’）间的基本矩阵；

然后将所有满足极限约束条件的左右图像标记点进行三维重建，重建出空间标记点在世界坐标系下的三维坐标值，其重建公式如下所示：

$x=\frac{z{X}_1}{f_1}$

$y=\frac{{\mathrm{zY}}_1}{f_1}---\left(4\right)$

$z=\frac{f_1(f_2{t}_y-y_2{t}_z)}{Y_1(r_7{X}_1+r_8{Y}_1+r_9{f}_1)f_2(r_4{X}_1+r_5{Y}_1+r_6{f}_1)}$

其中：x_l=[X₁ Y]₁，X₁,Y₁分别为左高速摄像机(8)拍摄的图像标记点x_l的横、纵坐标；x_r=[X₂ Y₂]，X₂,Y₂分别为右高速摄像机(8’)拍摄的图像标记点x_r的横、纵坐标；(x y z)为由左两图像标记点x_l和右两图像标记点x_r重建出来的空间标记点的三维坐标；f₁、f₂分别为左右高速摄像机(8、8’)的焦距；

$\left[\begin{array}{ccc}{r}_1& r_2& r_3\\ r_4& r_5& r_6\\ r_7& r_8& r_9\end{array}\right]$ 为右高速摄像机(8’)相对于左高速摄像机(8)的旋转矩阵，

[t_x t_y t_z]是右高速摄像机(8’)相对于左高速相机(8)的平移矩阵；

附图3被测物体空间摆动轨迹所显示的就是序列图像中被测物（6）表面标记点的空间坐标云图。

（4）转动惯量计算

●空间标记点平面拟合

将视频每一帧重建出的空间点坐标在同一坐标系下进行显示，并且将这些序列点进行平面拟合。平面拟合公式如下所示：

空间平面方程的通常可采用公式为：

Ax+By+Cz+1＝0(C≠0)         (5)

对于视频第i帧提取的空间坐标点，即i帧空间点坐标为：(x_i,y_i,z_i)

i＝1,2,Ln(n≥3)，所要拟合的平面可表达为如下的矩阵形式：

$\left[\begin{array}{ccc}{x}_1& y_1& z_1\\ M& M& M\\ x_n& y_n& z_n\end{array}\right]\left[\begin{array}{c}A\\ B\\ C\end{array}\right]=\left[\begin{array}{c}-1\\ -1\\ -1\end{array}\right]---\left(6\right)$

根据最小二乘法法拟合空间坐标点所在平面，将矩阵左乘 ${\left[\begin{array}{ccc}{x}_1& y_1& z_1\\ M& M& M\\ x_n& y_n& z_n\end{array}\right]}^T,$ 式（6）化简为：

$\left[\begin{array}{ccc}\mathrm{\Σ}{x_i}^2& \mathrm{\Σ}{x}_i{y}_i& \mathrm{\Σ}{x}_i{z}_i\\ \mathrm{\Σ}{x}_i{y}_i& \mathrm{\Σ}{y_i}^2& \mathrm{\Σ}{y}_i{z}_i\\ \mathrm{\Σ}{x}_i{z}_i& \mathrm{\Σ}{y}_i{z}_i& \mathrm{\Σ}{z_i}^2\end{array}\right]\left[\begin{array}{c}A\\ B\\ C\end{array}\right]=\left[\begin{array}{c}-\mathrm{\Σ}{x}_i\\ -\mathrm{\Σ}{y}_i\\ -\mathrm{\Σ}{z}_i\end{array}\right]---\left(7\right)$

可推导得：

$\left[\begin{array}{c}A\\ B\\ C\end{array}\right]={\left[\begin{array}{ccc}\mathrm{\Σ}{x_i}^2& \mathrm{\Σ}{x}_i{y}_i& \mathrm{\Σ}{x}_i{z}_i\\ \mathrm{\Σ}{x}_i{y}_i& \mathrm{\Σ}{y_i}^2& \mathrm{\Σ}{y}_i{z}_i\\ \mathrm{\Σ}{x}_i{z}_i& \mathrm{\Σ}{y}_i{z}_i& \mathrm{\Σ}{z_i}^2\end{array}\right]}^{-1}\left[\begin{array}{c}-\mathrm{\Σ}{x}_i\\ -\mathrm{\Σ}{y}_i\\ -\mathrm{\Σ}{z}_i\end{array}\right]---\left(8\right)$

即所得平面为拟合平面。

●空间标记点平面投影

将视频每帧提取的空间坐标点P_i＝(x_i,y_i,z_i)i=1,2,Ln(n≥3)投影到根据公式（8）拟合的平面，空间平面的单位法向量为：

$N_P=\left[\begin{array}{ccc}\frac{1}{B*C}& \frac{1}{A*C}& \frac{1}{A*B}\end{array}\right]---\left(9\right)$

且平面过P_p点(11(1+A+B)/(-C))。则空间坐标点在空间平面的投影点坐标为：

P＝P_i-N_p×((P_i-P_p)×N_p)   (10)P为P_i投影到由公式（10）拟合的平面上的投影点。

●被测标记点轨迹拟合与摆动周期计算

建立P中x轴空间坐标与时间的曲线，并且利用公式（11）进行曲线拟合。

x＝a×e^-b×t×sin(c×t+d)+e            (11)

其中a,d,e为积分常数，由初始条件决定。b为阻尼因子，

ω₀为系统固有频率。

附图4所展示的为将空间点云平面投影后的运动轨迹拟合曲线。

然后根据公式（12）计算得出被测物的摆动周期T：

$T=\frac{2\mathrm{\π}}{c}---\left(12\right)$

●计算被测物体转动惯量I

$I=\frac{T^2\mathrm{mgl}}{4{\mathrm{\π}}^2}---\left(13\right)$

其中,T为模型做单摆运动的周期,m为模型质量,g为重力加速度,l为模型质心距转轴的距离。通过对模拟目标进行视频拍摄,根据目标标志的运动周期推算其做单摆运动的周期,可由上式对模型的转动惯量进行求解。

本发明利用双目视觉系统可得知每个时刻被测物体在空间摆动的轨迹，将被测物每一时刻的空间位置进行平面拟合与投影，使得被测物空间运动坐标在同一平面内，减小摆动运动不在一个平面给测量精度带来的影响；并通过考虑阻尼系数的被测物实时空间运动坐标曲线拟合，减小空气阻力以及机构阻力对测量系统精度的影响，提高了转动惯量测量精度。

Claims (1)

1.一种基于双目视觉测量转动惯量的方法，其特征是，本方法采用双目视觉系统，通过左、右两台超高速摄像机8、8’采集被测物6表面标标记点的位置信息，经标记点提取、匹配与重建得到序列图像的标记点圆心坐标，再进行空间曲面拟合投影，最终计算转动惯量。该方法的具体测量步骤如下：

（1）左、右两个高速摄像机的标定

采用张氏标定方法标定出两相机的内外参数；

$s\left[\begin{array}{c}u\\ v\\ 1\end{array}\right]=\left[\begin{array}{cccc}{\mathrm{\α}}_x& 0& u_0& 0\\ 0& {\mathrm{\α}}_y& v_0& 0\\ 0& 0& 1& 0\end{array}\right]\left[\begin{array}{cc}R& t\\ 0^T& 1\end{array}\right]\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\\ 1\end{array}\right]---\left(1\right)$ 其中s为比例因子，u、v为相面坐标，X_w、Y_w、Z_w是世界坐标，α_xα_x、α_y、u₀、v₀为摄像机内参数，分别代表u轴尺度因子、v轴尺度因子、相面中心u轴坐标、相面中心v轴坐标；R为旋转矩阵、t为平移矩阵，它们是摄像机的外参数。

（2）标记点特征的提取

对采集到的视频每一帧进行处理，每一帧只有一个有效标记点，采用梯度重心法提取被被测物（6）运动时表面的标记点，首先采用高斯一阶微分算子对图像进行卷积操作，以获得图像在各点的梯度，则对于标志图案所在子区域的梯度重心点坐标的求取方法如下：

$C(x,y)=\underset{i=-h}{\overset{h}{\mathrm{\Σ}}}\underset{j=-w}{\overset{w}{\mathrm{\Σ}}}\left[\right|G(i,j)|\·P(i,j)]/\underset{i=-h}{\overset{h}{\mathrm{\Σ}}}\underset{j=-w}{\overset{w}{\mathrm{\Σ}}}\left|G(i,j)\right|---\left(2\right)$ 其中：C（x,y）为应用梯度重心法所提取的标记点中心像素坐标值，|G(i,j)|为(i,j)点的梯度幅值，w，h分别为目标图像的宽和高，P(i,j)为(i,j)点的图像坐标。

（3）标记点的匹配与重建

对采集到的视频每一帧所提取的标记点进行匹配与重建，图像标记点的匹配方法如下：先利用Longguet-Higgins提出的归一化8点算法计算出左右两高速摄像机（8、8’）的基本矩阵F，再通过左高速摄像机拍摄图片与右高速摄像机拍摄图片的极限约束关系进行图像标记点的初匹配，若左高速摄像机拍摄的图像标记点x_l与右高速摄像机拍摄的图像标记点x_r相匹配，即满足极限约束条件，极限约束条件如下所示：

$x_l^{T}F{x}_r=0---\left(3\right)$

其中：x_l为左高速摄像机拍摄的图像标记点；x_r为右高速摄像机拍摄的与x_l匹配的图像标记点；F为左右两个高速摄像机（8、8’）间的基本矩阵；

然后将所有满足极限约束条件的左右图像标记点进行三维重建，重建出空间标记点在世界坐标系下的三维坐标值，其重建公式如下所示：

$x=\frac{z{X}_1}{f_1}$

$y=\frac{{\mathrm{zY}}_1}{f_1}---\left(4\right)$

$z=\frac{f_1(f_2{t}_y-y_2{t}_z)}{Y_1(r_7{X}_1+r_8{Y}_1+r_9{f}_1)f_2(r_4{X}_1+r_5{Y}_1+r_6{f}_1)}$

其中：x_l=[X₁ Y]₁，X₁,Y₁分别为左高速摄像机（8）拍摄的图像标记点x_l的横、纵坐标；x_r=[X₂ Y₂]，X₂,Y₂分别为右高速摄像机（8’）拍摄的图像标记点x_r的横、纵坐标；(x y z)为由左两图像标记点x_l和右两图像标记点x_r重建出来的空间标记点的三维坐标；f₁、f₂分别为左右高速摄像机(8、8’)的焦距； $\left[\begin{array}{ccc}{r}_1& r_2& r_3\\ r_4& r_5& r_6\\ r_7& r_8& r_9\end{array}\right]$ 为右高速摄像机(8’)相对于左高速摄像机(8)的旋转矩阵，

[t_x t_y t_z]是右高速摄像机(8’)相对于左高速相机(8)的平移矩阵；

（4）转动惯量计算

●空间标记点平面拟合

将视频每一帧重建出的空间点坐标在同一坐标系下进行显示，并且将这些序列点进行平面拟合。平面拟合公式如下所示：

空间平面方程的通常可采用公式为：

Ax+By+Cz+1＝0(C≠0)        (5)

对于视频第i帧提取的空间坐标点，即i帧空间点坐标为：(x_i,y_i,z_i)i＝1,2,Ln(n≥3)，所要拟合的平面可表达为如下的矩阵形式：

$\left[\begin{array}{ccc}{x}_1& y_1& z_1\\ M& M& M\\ x_n& y_n& z_n\end{array}\right]\left[\begin{array}{c}A\\ B\\ C\end{array}\right]=\left[\begin{array}{c}-1\\ -1\\ -1\end{array}\right]---\left(6\right)$

根据最小二乘法法拟合空间坐标点所在平面，将矩阵左乘 ${\left[\begin{array}{ccc}{x}_1& y_1& z_1\\ M& M& M\\ x_n& y_n& z_n\end{array}\right]}^T,$ 化简为：

$\left[\begin{array}{ccc}\mathrm{\Σ}{x_i}^2& \mathrm{\Σ}{x}_i{y}_i& \mathrm{\Σ}{x}_i{z}_i\\ \mathrm{\Σ}{x}_i{y}_i& \mathrm{\Σ}{y_i}^2& \mathrm{\Σ}{y}_i{z}_i\\ \mathrm{\Σ}{x}_i{z}_i& \mathrm{\Σ}{y}_i{z}_i& \mathrm{\Σ}{z_i}^2\end{array}\right]\left[\begin{array}{c}A\\ B\\ C\end{array}\right]=\left[\begin{array}{c}-\mathrm{\Σ}{x}_i\\ -\mathrm{\Σ}{y}_i\\ -\mathrm{\Σ}{z}_i\end{array}\right]---\left(7\right)$ 可推导得：

$\left[\begin{array}{c}A\\ B\\ C\end{array}\right]={\left[\begin{array}{ccc}\mathrm{\Σ}{x_i}^2& \mathrm{\Σ}{x}_i{y}_i& \mathrm{\Σ}{x}_i{z}_i\\ \mathrm{\Σ}{x}_i{y}_i& \mathrm{\Σ}{y_i}^2& \mathrm{\Σ}{y}_i{z}_i\\ \mathrm{\Σ}{x}_i{z}_i& \mathrm{\Σ}{y}_i{z}_i& \mathrm{\Σ}{z_i}^2\end{array}\right]}^{-1}\left[\begin{array}{c}-\mathrm{\Σ}{x}_i\\ -\mathrm{\Σ}{y}_i\\ -\mathrm{\Σ}{z}_i\end{array}\right]---\left(8\right)$

即所得平面为拟合平面。

●空间标记点平面投影

将视频每帧提取的空间坐标点P_i＝(x_i,y_i,z_i)i=1,2,Ln(n≥3)投影到根据公式（8）拟合的平面，空间平面的单位法向量为：

$N_P=\left[\begin{array}{ccc}\frac{1}{B*C}& \frac{1}{A*C}& \frac{1}{A*B}\end{array}\right]---\left(9\right)$

且平面过P_p点(1 1(1+A+B)/(-C))。则空间坐标点在空间平面的投影点坐标为：

P＝P_i-N_p×((P_i-P_p)×N_p)    (10)

P为P_i投影到由公式（10）拟合的平面上的投影点。

●被测标记点轨迹拟合与摆动周期计算

建立P中x轴空间坐标与时间的曲线，并且利用公式（11）进行曲线拟合。

x＝a×e^-b×t×sin(c×t+d)+e    (11)

其中a,d,e为积分常数，由初始条件决定。b为阻尼因子，

ω₀为系统固有频率，然后根据公式（12）计算得出被测物的摆动周期T。

$T=\frac{2\mathrm{\π}}{c}---\left(12\right)$

●被测物体转动惯量计算

$I=\frac{T^2\mathrm{mgl}}{4{\mathrm{\π}}^2}---\left(13\right)$

其中,T为模型做单摆运动的周期,m为模型质量,g为重力加速度,l为模型质心距转轴的距离。通过对模拟目标进行视频拍摄,根据目标标志的运动周期推算其做单摆运动的周期,可由上式对模型的转动惯量进行求解。

Images