CN104374338A

CN104374338A - Single-axis rotation angle vision measurement method based on fixed camera and single target

Info

Publication number: CN104374338A
Application number: CN201410506279.0A
Authority: CN
Inventors: 全权; 董洪信; 张瑞峰; 付强; 蔡开元
Original assignee: Beihang University
Current assignee: Beihang University
Priority date: 2014-09-28
Filing date: 2014-09-28
Publication date: 2015-02-25
Anticipated expiration: 2034-09-28
Also published as: CN104374338B

Abstract

The invention provides a single-axis rotation angle vision measurement method based on a fixed camera and a single target. The method includes the four steps of firstly, calibrating internal parameters and external parameters of the camera; secondly, obtaining the initial values of the external parameters related to pivoting movement; thirdly, conducting nonlinear optimization on the internal parameter pin and the external parameter pout; fourthly, measuring the rotation angle through the calibrated internal parameters and the calibrated external parameters. The method is a method for measuring the rotation angle of the pivoting movement through the monocular camera, and solves the problem that the target is hard to install in an existing angle vision measurement method; meanwhile, the method is suitable for the conditions with blockage as well and is quite high in practicability.

Description

A Visual Measurement Method of Single-axis Rotation Angle Based on Fixed Camera and Single Target

技术领域technical field

本发明涉及一种基于固定相机和单靶标的单轴旋转角的视觉测量方法，属于视觉测量技术领域。The invention relates to a visual measurement method based on a fixed camera and a single-axis rotation angle of a single target, belonging to the technical field of visual measurement.

背景技术Background technique

角位移测量广泛应用于工业、航空航天等领域。由于在恶劣的环境中，角位移传感器的电气、机械等参数会发生变化，因此需要定期对其进行标定矫正。视觉测量技术因具有结构简单，非接触实时测量，现场安装和调试方便等优点而备受青睐。现有的视觉测角方法分为双目视觉和单目视觉的方法。双目视觉方法，因其公共视场小，安装复杂而受到很大的限制。现有的单目视觉的方法大都对靶标的安装有一定的要求，安装难度大。因此，本专利申请提出了一种基于固定相机和单靶标的单轴旋转角的视觉测量方法，该方法对靶标的安装没有任何要求。另外，本发明对有遮挡的情况同样成立，更加实用。Angular displacement measurement is widely used in industry, aerospace and other fields. Since the electrical and mechanical parameters of the angular displacement sensor will change in harsh environments, it is necessary to calibrate and correct them regularly. Vision measurement technology is popular because of its simple structure, non-contact real-time measurement, and convenient on-site installation and debugging. The existing visual angle measurement methods are divided into binocular vision and monocular vision methods. The binocular vision method is greatly limited due to its small public field of view and complicated installation. Most of the existing monocular vision methods have certain requirements for the installation of the target, and the installation is difficult. Therefore, this patent application proposes a visual measurement method based on a fixed camera and a single-axis rotation angle of a single target, which does not require any installation of the target. In addition, the present invention is also valid for the situation of occlusion, which is more practical.

发明内容Contents of the invention

本发明提出了一种基于固定相机和单靶标的单轴旋转角的视觉测量方法，它是一种利用单目相机对绕轴旋转运动的旋转角进行测量的方法。它解决了现有视觉测角方法中靶标安装困难的问题。同时，本发明对有遮挡的情况同样适用，实用性很好。The invention proposes a single-axis rotation angle visual measurement method based on a fixed camera and a single target, which is a method for measuring the rotation angle of a rotation around an axis by using a monocular camera. It solves the problem of difficult target installation in existing visual goniometric methods. Simultaneously, the present invention is equally applicable to the situation that has shelter, and practicability is very good.

该方法，通过二维靶标，对相机内参及相对于轴的外参数进行标定，进而测量旋转角，该方法对靶标安装没有任何要求，可以将靶标固定到任意位置。In this method, the internal parameters of the camera and the external parameters relative to the axis are calibrated through a two-dimensional target, and then the rotation angle is measured. This method does not require any installation of the target, and the target can be fixed at any position.

本发明采用的摄像机模型为非线性透视投影模型，模型描述如下：The camera model adopted in the present invention is a nonlinear perspective projection model, and the model is described as follows:

1.线性摄像机模型1. Linear camera model

如图(1)所示，空间任何一点P在图像中的成像位置可以用针孔成像模型近似表示，即点P在图像中的投影位置p，为光心O与点P的连线与图像平面的交点。因而世界坐标系下P点坐标(X_w,Y_w,Z_w)^T与投影点p的像素坐标(u,v)^T之间的关系如下：As shown in Figure (1), the imaging position of any point P in the image can be approximated by the pinhole imaging model, that is, the projection position p of point P in the image is the connection line between optical center O and point P and the image intersection of the planes. Therefore, the relationship between the coordinates (X _w , Y _w , Z _w ) ^T of point P in the world coordinate system and the pixel coordinates (u, v) ^T of the projected point p is as follows:

$s the s [\begin{matrix} u u \\ v v \\ 11 \end{matrix}] = = [\begin{matrix} {α α}_{x x} & 00 & {u u}_{00} & 00 \\ 00 & {α α}_{y the y} & {v v}_{00} & 00 \\ 00 & 00 & 11 & 00 \end{matrix}] [\begin{matrix} {R R}_{cw cw} & {T T}_{cw cw} \\ 00 & 11 \end{matrix}] [\begin{matrix} {X x}_{w w} \\ {Y Y}_{w w} \\ {Z Z}_{w w} \\ 11 \end{matrix}] - - - - - - ((1.1 1.1))$

其中α_x＝f/dX为u轴上的尺度因子，α_y＝f/dY为v轴上的尺度因子。α_x,α_y,u₀,v₀只与摄像机内部参数有关，称为摄像机的内参数。R_cw,T_cw分别为相机坐标系与世界坐标系的旋转矩阵和平移向量，称为摄像机的外部参数。Where α _x =f/dX is the scale factor on the u axis, and α _y =f/dY is the scale factor on the v axis. α _x , α _y , u ₀ , v ₀ are only related to the internal parameters of the camera, which are called internal parameters of the camera. R _cw , T _cw are the rotation matrix and translation vector of the camera coordinate system and the world coordinate system respectively, which are called the external parameters of the camera.

1.非线性摄像机模型1. Nonlinear camera model

实际上，镜头并不是理想的透视成像，而是带有一定程度的畸变。本发明采用的是非线性畸变模型，由世界坐标(X_w,Y_w,Z_w)^T求取像素坐标(u,v)^T的过程如下：In fact, the lens is not an ideal perspective imaging, but has a certain degree of distortion. The present invention uses a nonlinear distortion model, and the process of obtaining pixel coordinates (u, v) ^T from world coordinates (Xw _, _Yw , _Zw ) ^T is as follows:

$[\begin{matrix} {X x}_{c c} \\ {Y Y}_{c c} \\ {Z Z}_{c c} \end{matrix}] = = {R R}_{cw cw} [\begin{matrix} {X x}_{w w} \\ {Y Y}_{w w} \\ {Z Z}_{w w} \end{matrix}] + + {T T}_{cw cw},, - - - - - - ((1.2 1.2))$

x＝X_c/Z_c,y＝Y_c/Z_c,r²＝x²+y², (1.3)x=X _c /Z _c , y=Y _c /Z _c , r ² =x ² +y ² , (1.3)

${x x}_{d d} = = [\begin{matrix} ((11 + + {k k}_{11} {r r}^{22} + + {k k}_{22} {r r}^{44} + + {k k}_{55} {r r}^{66})) x x + + {22 k k}_{33} xy xy + + {k k}_{44} (({r r}^{22} + + {22 x x}^{22})) \\ ((11 + + {k k}_{11} {r r}^{22} + + {k k}_{22} {r r}^{44} + + {k k}_{55} {r r}^{66})) y the y + + {22 k k}_{44} xy xy + + {k k}_{33} (({r r}_{22} + + {22 y the y}^{22})) \end{matrix}],, - - - - - - ((1.4 1.4))$

$[\begin{matrix} u u \\ v v \\ 11 \end{matrix}] = = [\begin{matrix} {α α}_{x x} & 00 & {u u}_{00} \\ 00 & {α α}_{y the y} & {v v}_{00} \\ 00 & 00 & 11 \end{matrix}] [\begin{matrix} {x x}_{d d} \\ 11 \end{matrix}],, - - - - - - ((1.5 1.5))$

其中，k₁,k₂,k₃,k₄,k₅为畸变系数。α_x,α_y,u₀,v₀,k₁,k₂,k₃,k₄,k₅构成摄像机的内参数，R_cw,T_cw为摄像机的外参数。Among them, k ₁ , k ₂ , k ₃ , k ₄ , and k ₅ are distortion coefficients. α _x , α _y , u ₀ , v ₀ , k ₁ , k ₂ , k ₃ , k ₄ , and k ₅ constitute the internal parameters of the camera, and R _cw and T _cw are the external parameters of the camera.

为描述方便，我们记u,v为p_in＝(α_x,α_y,u₀,v₀,k₁,k₂,k₃,k₄,k₅)^T,R_cw,T_cw,X_w,Y_w,Z_w的函数，即For the convenience of description, we record u and v as _pin =(α _x ,α _y ,u ₀ ,v ₀ ,k ₁ ,k ₂ ,k ₃ ,k ₄ ,k ₅ ) ^T ,R _cw ,T _cw ,X function of _w , Y _w , Z _w , namely

u＝f_u(p_in,R_cw,T_cw,X_w,Y_w,Z_w),v＝f_v(p_in,R_cw,T_cw,X_w,Y_w,Z_w) (1.6)u＝f _u (p _in , R _cw , T _cw , X _w , Y _w , Z _w ), v=f _v (p _in , R _cw , T _cw , X _w , Y _w , Z _w ) (1.6)

本发明所建立的世界坐标系如图(2).取转轴为Z_w轴，并使得在基准位置时靶标原点位于X_w轴上。假设基准位置靶标坐标系与世界坐标系之间的旋转矩阵和平移向量为并设相机坐标系与世界坐标系的旋转矩阵为R_cw，平移向量为T_cw.The world coordinate system established by the present invention is shown in figure (2). The axis of rotation is taken as the _Zw axis, and the target origin is located on the _Xw axis at the reference position. Suppose the rotation matrix and translation vector between the reference position target coordinate system and the world coordinate system are And set the rotation matrix between the camera coordinate system and the world coordinate system as R _cw , and the translation vector as T _cw .

本发明所提出的方法包括两步：内外参数信息标定和任意角度求解。内外参数信息标定，是在普通相机标定的基础上，利用绕轴旋转运动坐标系之间的特殊几何关系，求取所需要的外部参数信息的过程。t,R_cw,T_cw及各位置与基准位置之间的夹角是内外参信息标定所要估计的外参数。任意角度求解则是根据已求得的内外参数信息，求取旋转角。在任意角度求解过程中，我们提出了一种点匹配算法，从而使测角方法可以求取遮挡问题的解。The method proposed by the invention includes two steps: calibration of internal and external parameter information and solution at any angle. Calibration of internal and external parameter information is the process of obtaining the required external parameter information on the basis of ordinary camera calibration, using the special geometric relationship between the rotating coordinate system around the axis. t, R _cw , T _cw and the angle between each position and the reference position are the external parameters to be estimated in the calibration of internal and external parameter information. Arbitrary angle solution is to calculate the rotation angle based on the obtained internal and external parameter information. During the arbitrary angle solution, we propose a point-matching algorithm that enables the goniometric method to solve the occlusion problem.

综上所述，本发明一种基于固定相机和单靶标的单轴旋转角的视觉测量方法，该方法具体步骤如下：In summary, the present invention is a visual measurement method based on a fixed camera and a single-axis rotation angle of a single target. The specific steps of the method are as follows:

步骤一：相机内外参数标定。Step 1: Camera internal and external parameters calibration.

在具体实现过程中，我们利用二维棋盘格(如图5)作为靶标。我们采用的图片序列，包含一张基准位置时的图片和m张旋转部件旋转到不同角度时的图片。假设这m个角度为θ₁,θ₂,...,θ_m。接下来，可以利用MATLAB 2014自带的棋盘格角点提取函数detectChecker--boardPoints进行角点提取，并结合MATLAB标定工具进行内参标定,求得内参数p_in和m+1个位姿的外部参数。该标定工具箱，可以在http://www.vision.caltech.edu/bouguetj/calib_doc/上下载。这些外部参数可以表示为In the specific implementation process, we use a two-dimensional checkerboard (as shown in Figure 5) as the target. The picture sequence we use includes a picture at the reference position and m pictures when the rotating part is rotated to different angles. Suppose these m angles are θ ₁ , θ ₂ ,...,θ _m . Next, you can use the checkerboard corner extraction function detectChecker--boardPoints that comes with MATLAB 2014 to extract the corner points, and combine the MATLAB calibration tool for internal parameter calibration to obtain the internal parameters p _in and external parameters of m+1 poses . The calibration toolbox can be downloaded at http://www.vision.caltech.edu/bouguetj/calib_doc/. These external parameters can be expressed as

$[[{R R}_{{cb cb}_{i i}},, {T T}_{{cb cb}_{i i}}]],, i i = = 0,1,2 0,1,2,, . . . . . .,, m m . .$

表示在旋转角为θ_i时，相机与靶标之间的旋转矩阵与平移向量，i＝0时表示基准位置时，相机与靶标之间的旋转矩阵与平移向量。 Indicates the rotation matrix and translation vector between the camera and the target when the rotation angle is θ _i , and when i=0 indicates the rotation matrix and translation vector between the camera and the target at the reference position.

步骤二：绕轴旋转运动相关的外参数初值求取Step 2: Calculate the initial value of the external parameters related to the rotation around the axis

1、我们根据步骤一中求得的外参数，由(1.7)求取在旋转角为θ_i时，靶标相对于基准时靶标的旋转矩阵和平移向量 $[R_{b_{i} b_{0}}, T_{b_{i} b_{0}}], i = 1,2, . . ., m .$ 1. According to the external parameters obtained in step 1, we use (1.7) to obtain the rotation matrix and translation vector of the target relative to the reference time when the rotation angle is θ _i $[R_{b_{i} b_{0}}, T_{b_{i} b_{0}}], i = 1,2, . . ., m .$

${R R}_{{b b}_{i i} {b b}_{00}} = = {R R}_{{cb cb}_{i i}}^{T T} {R R}_{{cb cb}_{00}},, {T T}_{{b b}_{i i} {b b}_{00}} = = {R R}_{{cb cb}_{i i}}^{T T} (({T T}_{{cb cb}_{00}} - - {T T}_{{cb cb}_{i i}})) . . - - - - - - ((1.7 1.7))$

2、求取t及θ_i初值：2. Seek Initial value of t and θ _i :

定义r₁,r₂,r₃分别为的第1,2,3列，且向量[t]_×表示如下矩阵Define r ₁ , r ₂ , r ₃ as Columns 1, 2, and 3 of , and the vector [t] _× represents the following matrix

${[\begin{matrix} {t t}_{x x} \\ {t t}_{y the y} \\ {t t}_{z z} \end{matrix}]}_{\times \times} = = [\begin{matrix} 00 & {- - t t}_{z z} & {t t}_{y the y} \\ {t t}_{z z} & 00 & {- - t t}_{x x} \\ {- - t t}_{y the y} & {t t}_{x x} & 00 \end{matrix}] . . - - - - - - ((1.8 1.8))$

则，我们可以求得Then, we can obtain

${r r}_{22} = = \frac{{(({S S}_{3,2 3,2},, {S S}_{1,3 1,3},, {S S}_{2,1 2,1}))}^{T T}}{| | | | (({S S}_{3,2 3,2},, {S S}_{1,3 1,3},, {S S}_{2,1 2,1})) | | | |},, μ μ = = t t sin sin θ θ = = \frac{11}{22 | | | | (({S S}_{3,2 3,2},, {S S}_{1,3 1,3},, {S S}_{2,1 2,1})) | | | |},, - - - - - - ((1.9 1.9))$

其中，S_i,j表示矩阵S第i行j列的元素。然后，r₁可由下式求得in, S _i,j represents the element in row i and column j of matrix S. Then, _r1 can be obtained by

$((\frac{11}{μ μ} {[[{T T}_{{b b}_{i i} {b b}_{00}}]]}_{\times \times} {R R}_{{b b}_{i i} {b b}_{00}} + + {[[{r r}_{22}]]}_{\times \times})) x x = = 00,, | | | | x x | | | | = = 11 . . - - - - - - ((1.10 1.10))$

上式可以通过对作SVD分解求得，若 $\frac{1}{μ} {[T_{b_{i} b_{0}}]}_{\times} R_{b_{i} b_{0}} + {[r_{2}]}_{\times} = {USV}^{T},$ 则x为V的第三列。然后，我们可以得到The above formula can be passed to It can be obtained by SVD decomposition, if $\frac{1}{μ} {[T_{b_{i} b_{0}}]}_{\times} R_{b_{i} b_{0}} + {[r_{2}]}_{\times} = {USV}^{T},$ Then x is the third column of V. Then, we can get

${R R}_{{b b}_{00} w w} = = [\begin{matrix} {r r}_{11} & {r r}_{22} & {r r}_{11} \times \times {r r}_{22} \end{matrix}] . . - - - - - - ((1.11 1.11))$

θ_i＝2arctan((A^TA)^-1(A^TB))， (1.12) _θi = 2arctan((A ^T A) ^-1 (A ^T B)), (1.12)

其中 $A = vec (r_{2} r_{3}^{T} + r_{3} r_{2}^{T}),$ $B = vec (\frac{1}{μ} {[T_{b_{i} b_{0}}]}_{\times} R_{b_{i} b_{0}} + {[r_{2}]}_{\times}) .$ 这里的vec(·)表示矩阵的拉直。所以，in $A = vec (r_{2} r_{3}^{T} + r_{3} r_{2}^{T}),$ $B = vec (\frac{1}{μ} {[T_{b_{i} b_{0}}]}_{\times} R_{b_{i} b_{0}} + {[r_{2}]}_{\times}) .$ Here vec( ) represents the straightening of the matrix. so,

${R R}_{{wb wb}_{00}} = = {R R}_{{b b}_{00} w w}^{T T} - - - - - - ((1.13 1.13))$

$t t = = \frac{μ μ}{sin sin {θ θ}_{i i}} . . - - - - - - ((1.14 1.14))$

3、根据求得的t及θ_i初值求解，根据(1.16)，求R_cw,T_cw的值3. According to the obtained Find the initial value of t and θ _i , according to (1.16), find the values of R _cw and T _cw

${R R}_{cw cw} = = {R R}_{{cb cb}_{00}} {R R}_{{wb wb}_{00}}^{T T},, {T T}_{cw cw} = = {((t t,, 0,0 0,0))}^{T T} - - {R R}_{cw cw} {T T}_{{wb wb}_{00}} . . - - - - - - ((1.15 1.15))$

步骤三：对内参数p_in和外参数p_out进行非线性优化。这里p_in∈□⁹,p_out∈□¹⁰表示如下Step 3: Perform nonlinear optimization on the internal parameter p _in and the external parameter p _out . Here p _in ∈ □ ⁹ , p _out ∈ □ ¹⁰ are expressed as follows

${p p}_{in in} = = {(({α α}_{x x},, {α α}_{y the y},, {u u}_{00,,} {v v}_{00},, {k k}_{11},, {k k}_{22},, {k k}_{33},, {k k}_{44},, {k k}_{55}))}^{T T},, {p p}_{out out} = = {(({v v}_{{wb wb}_{00}}^{T T},, t t,, {v v}_{cw cw}^{T T},, {T T}_{cw cw}^{T T}))}^{T T},,$

其中，分别为R_cw所对应的Rodrigues向量。一个Rodrigues向量v与其对应的旋转矩阵R间的转换关系如下in, respectively The Rodrigues vector corresponding to R _cw . The conversion relationship between a Rodrigues vector v and its corresponding rotation matrix R is as follows

$sin sin ((| | | | v v | | | |)) {[[v v]]}_{\times \times} = = \frac{R R - - {R R}^{T T}}{22},,$

$R R = = cos cos ((| | | | v v | | | |)) {I I}_{33} + + \frac{((11 - - cos cos ((| | | | v v | | | |))))}{{| | | | v v | | | |}^{22}} {vv vv}^{T T} + + sin sin ((| | | | v v | | | |)) {[[v v]]}_{\times \times},,$

其中，I₃表示3阶的单位矩阵。Among them, I ₃ represents the identity matrix of order 3.

我们对各图像点的坐标作如下标记，靶标1上，在旋转角为θ_i时，对应于第j个靶标特征点的角点坐标记为 $(u_{j}^{i}, v_{j}^{i}), i = 0,1, . . ., m, j = 1,2, . . ., n,$ 靶标上的点记为We mark the coordinates of each image point as follows. On target 1, when the rotation angle is _θi , the coordinates of the corner point corresponding to the jth target feature point are marked as $(u_{j}^{i}, v_{j}^{i}), i = 0,1, . . ., m, j = 1,2, . . ., no,$ The point on the target is marked as

X_b＝{(x_j，y_j，0)^T，j＝1，2，...，n}.X _b ={(x _j ,y _j ,0) ^T ,j=1,2,...,n}.

由于图像提取过程中是存在误差的，因此，我们选取图像的重投影误差和作为优化目标。由靶标点及p_out，获取重投影点的过程如下：Since there are errors in the process of image extraction, we choose the reprojection error sum of the image as the optimization target. Obtain the reprojection point from the target point and p _out The process is as follows:

1、旋转角为θ_i时，靶标点转化为世界坐标系下的点 1. When the rotation angle is θ _i , the target point Convert to a point in the world coordinate system

$[\begin{matrix} {x x}_{i i,, j j}^{w w} \\ {y the y}_{i i,, j j}^{w w} \\ {z z}_{i i,, j j}^{w w} \end{matrix}] = = [\begin{matrix} cos cos {θ θ}_{i i} & - - {sin sin θ θ}_{i i} & 00 \\ {sin sin θ θ}_{i i} & {cos cos θ θ}_{i i} & 00 \\ 00 & 00 & 11 \end{matrix}] (({R R}_{{wb wb}_{00}} [\begin{matrix} {x x}_{j j} \\ {y the y}_{j j} \\ 00 \end{matrix}] + + [\begin{matrix} t t \\ 00 \\ 00 \end{matrix}])) . . - - - - - - ((1.16 1.16))$

2、求取相应的重投影点 2. Find the corresponding reprojection point

重投影误差定义为The reprojection error is defined as

${E E.}_{r r} = = {Σ Σ}_{i i = = 11}^{m m} {Σ Σ}_{j j = = 11}^{n no} [[{(({u u}_{i i}^{j j} - - {g g}_{u u} (({p p}_{in in},, {p p}_{out out},, {x x}_{j j},, {y the y}_{j j},, {θ θ}_{i i}))))}^{22} + + {(({v v}_{i i}^{j j} - - {g g}_{v v} (({p p}_{in in},, {p p}_{out out},, {x x}_{j j},, {y the y}_{j j},, {θ θ}_{i i}))))}^{22}]] - - - - - - ((1.19 1.19))$

优化问题可以描述为The optimization problem can be described as

$(({p p}_{in in}^{* *},, {p p}_{out out}^{* *},, {θ θ}_{11}^{* *},, {θ θ}_{22}^{* *},, . . . . . .,, {θ θ}_{m m}^{* *})) = = arg arg \underset{{p p}_{in in},, {p p}_{out out},,,, {θ θ}_{11},, . . . . . .,, {θ θ}_{m m}}{min min} {E E.}_{r r} . . - - - - - - ((1.20 1.20))$

该优化问题可以利用稀疏的Levenberg-Marquardt算法来解决优化问题，稀疏LM算法的优化工具箱可以在http://users.ics.forth.gr/～lourakis/sparseLM/上下载。利用上述优化，我们可以获得t,R_cw,T_cw的信息，这些信息将用于步骤四的旋转角度求取。The optimization problem can be solved by using the sparse Levenberg-Marquardt algorithm. The optimization toolbox of the sparse LM algorithm can be downloaded from http://users.ics.forth.gr/~lourakis/sparseLM/. Using the above optimization, we can obtain t, R _cw , T _cw information, these information will be used in step 4 to calculate the rotation angle.

步骤四：利用标定的内外参数，测量旋转角度。Step 4: Use the calibrated internal and external parameters to measure the rotation angle.

首先，我们根据已经标定好的信息，对角点进行匹配。假设在θ处拍摄的照片的角点坐标为X^c＝{(u_i,v_i)^T,i＝1,2,...,N_c}，此时图片中可能存在遮挡，因此匹配并不能自动完成，可以根据如下算法：First, we match the corner points based on the calibrated information. Assuming that the corner coordinates of the photo taken at θ are X ^c ={(u _i ,v _i ) ^T ,i=1,2,...,N _c }, there may be occlusion in the picture at this time, so matching and Can not be done automatically, can be based on the following algorithm:

1、寻找可能匹配的点集。1. Find possible matching point sets.

对每个角点与靶标点对计算一个得分S_i,j，计算方法如下：For each corner point and target point pair To calculate a score S _i,j , the calculation method is as follows:

(1)求取(x_j,y_j,0)^T所对应的世界坐标 (1) Obtain the world coordinates corresponding to (x _j ,y _j ,0) ^T

$[\begin{matrix} {x x}_{j j}^{w w} \\ {y the y}_{j j}^{w w} \\ {z z}_{j j}^{w w} \end{matrix}] = = {R R}_{{wb wb}_{00}} [\begin{matrix} {x x}_{j j} \\ {y the y}_{j j} \\ 00 \end{matrix}] + + [\begin{matrix} t t \\ 00 \\ 00 \end{matrix}] . .$

(2)求取像点(u_i,v_i)^T在平面上的原像求取方法如下(2) Calculate the image point (u _i , v _i ) ^T at original image on the plane The method of obtaining is as follows

${X x}_{w w}^{j j} = = ρ ρ (({x x}_{i i}^{n no},, {y the y}_{i i}^{n no},, 11)) {r r}_{11} + + t t$

${Y Y}_{w w}^{j j} = = ρ ρ (({x x}_{i i}^{n no},, {y the y}_{i i}^{n no},, 11)) {r r}_{22},,$

${Z Z}_{w w}^{j j} = = {z z}_{w w}^{j j},,$

其中， $ρ = \frac{z_{w}^{j}}{(x_{i}^{n}, y_{i}^{n}, 1) r_{3}},$ $[r_{1}, r_{2}, r_{3}] = R_{{wb}_{0}} \cdot {(x_{i}^{n}, y_{i}^{n})}^{T}$ 为(u_i,v_i)^T对应的归一化坐标，可以由迭代求得。in, $ρ = \frac{z_{w}^{j}}{(x_{i}^{no}, {the y}_{i}^{no}, 1) r_{3}},$ $[r_{1}, r_{2}, r_{3}] = R_{{wb}_{0}} &Center Dot; {(x_{i}^{no}, {the y}_{i}^{no})}^{T}$ is the normalized coordinate corresponding to (u _i , v _i ) ^T , which can be obtained by iteration.

(3)我们定义S_i,j如下式(3) We define S _i,j as follows

${S S}_{i i,, j j} = = \sqrt{{(({X x}_{w w}^{j j}))}^{22} + + {(({Y Y}_{w w}^{j j}))}^{22}} - - \sqrt{{(({x x}_{w w}^{j j}))}^{22} + + {(({y the y}_{w w}^{j j}))}^{22}} . . - - - - - - ((1.21 1.21))$

(4)选取阈值τ_m＝0.1d,这里d是棋盘格格子的长度，然后定义可能匹配的点的集合(4) Select the threshold τ _m = 0.1d, where d is the length of the checkerboard grid, and then define the set of possible matching points

$Ω Ω = = {{{(({u u}_{i i},, {v v}_{i i}))}^{T T} &LeftRightArrow; &LeftRightArrow; {(({x x}_{j j},, {y the y}_{j j},, 00))}^{T T},, j j = = 1,2 1,2,, . . . . . .,, n no,, i i = = 1,2 1,2,, . . . . . .,, {N N}_{c c} | | {S S}_{i i,, j j} < < {τ τ}_{m m}}} . . - - - - - - ((1.22 1.22))$

2、利用直方图，计算角度。2. Use the histogram to calculate the angle.

对于Ω中的元素，由下式计算θ_i,j For elements in Ω, θ _i,j is calculated by

$[\begin{matrix} {x x}_{j j}^{w w} \\ {y the y}_{j j}^{w w} \\ {z z}_{j j}^{w w} \end{matrix}] = = [\begin{matrix} cos cos {θ θ}_{i i,, j j} & - - {sin sin θ θ}_{i i,, j j} & 00 \\ {sin sin θ θ}_{i i,, j j} & {cos cos θ θ}_{i i,, j j} & 00 \\ 00 & 00 & 11 \end{matrix}] [\begin{matrix} {X x}_{j j}^{w w} \\ {Y Y}_{j j}^{w w} \\ {Z Z}_{j j}^{w w} \end{matrix}] . . - - - - - - ((1.23 1.23))$

然后画出分布的直方图，该直方图中有360个直方块，第k个直方块中的元素为((k-181)°,(k-179)°)∩Ω中的元素个数。我们对θ的估计值为then draw The histogram of the distribution, there are 360 histograms in the histogram, and the elements in the kth histogram are the number of elements in ((k-181)°,(k-179)°)∩Ω. Our estimate of θ is

3、利用求得的角度，选择正确的匹配对。3. Use the obtained angle to select the correct matching pair.

i所对应的靶标点M(j)对应于在正确角范围内，得分最小的点，也就是The target point M(j) corresponding to i corresponds to the point with the smallest score within the correct angle range, that is,

$M (i) = \min_{j} S_{i, j},$ 满足θ_i,j∈[θ-1°,θ+1°]. $m (i) = \min_{j} S_{i, j},$ Satisfy θ _i,j ∈[θ-1°,θ+1°].

若不存在θ_i,j∈[θ-1°,θ+1°],则令M(i)＝0.If there is no θ _{i, j} ∈ [θ-1°, θ+1°], let M(i)=0.

4、非线性优化4. Nonlinear optimization

优化问题描述为The optimization problem is described as

${θ θ}^{* *} = = arg arg \underset{θ θ}{min min} {Σ Σ}_{i i = = 11}^{{N N}_{c c}} δ δ ((i i)) {(((({u u}_{i i} - - {g g}_{u u} (({p p}_{in in}^{* *},, {p p}_{out out}^{* *},, {x x}_{M m ((i i))},, {y the y}_{M m ((i i))},, θ θ))))}^{22} + + (({{v v}_{i i} - - {g g}_{v v} (({p p}_{in in}^{* *},, {p p}_{out out}^{* *},, {x x}_{M m ((i i))},, {y the y}_{M m ((i i))},, θ θ))))}^{22})) - - - - - - ((1.24 1.24))$

其中， $δ (i) = \{\begin{matrix} 1, M (i) > 0 \\ 0, M (i) = 0 \end{matrix} .$ in, $δ (i) = \{\begin{matrix} 1, m (i) > 0 \\ 0, m (i) = 0 \end{matrix} .$

该优化问题可以通过MATLAB自带的lsqnonlin函数进行求解。This optimization problem can be solved by the lsqnonlin function that comes with MATLAB.

优点与功效：Advantages and effects:

本发明是一种基于固定相机和单靶标的单轴旋转角的视觉测量方法，它是一种利用单目相机和二维棋盘格对单轴旋转角进行测量的方法，其优点在于，解决了当前视觉测角问题对安装有诸多要求，不方便的问题。本发明采用二维靶标作为标定特征物，通过优化算法，可以单轴旋转角。该方法具有安装过程简单，标定精度高，成本低廉等优点，同时该方法在标定后可以解决遮挡问题。The present invention is a visual measurement method for a single-axis rotation angle based on a fixed camera and a single target. It is a method for measuring a single-axis rotation angle by using a monocular camera and a two-dimensional checkerboard. Its advantage is that it solves the problem of The current problem of visual angle measurement has many requirements for installation and is inconvenient. The present invention adopts a two-dimensional target as a calibration feature, and through an optimization algorithm, a single-axis rotation angle can be achieved. This method has the advantages of simple installation process, high calibration accuracy, low cost, etc., and at the same time, this method can solve the occlusion problem after calibration.

附图说明Description of drawings

图1：针孔成像模型示意图Figure 1: Schematic diagram of the pinhole imaging model

图2：实验设备安装示意图Figure 2: Schematic diagram of experimental equipment installation

图3a：测量误差均值随噪声变化图Figure 3a: Graph of mean measurement error versus noise

图3b：测量误差方差随噪声变化图Figure 3b: Plot of measurement error variance versus noise

图4a：角度变化时检测出的角点个数图Figure 4a: The number of corner points detected when the angle changes

图4b：角度变化时，角度误差图Figure 4b: Angle error plot when the angle changes

图5：真实实验所使用的二维棋盘格靶标图片Figure 5: Picture of the 2D checkerboard target used in the real experiment

图6：本发明流程框图Fig. 6: flow chart of the present invention

图中符号说明如下：The symbols in the figure are explained as follows:

图1中的符号说明：O表示摄像机光心,O₁表示图像坐标系原点。X,Y表示图像坐标系的坐标轴，x,y,z表示摄像机坐标系的坐标轴。P(X_w,Y_w,Z_w)表示三维点在摄像机坐标系下的坐标,p表示点P在图像上的投影。Explanation of symbols in Figure 1: O represents the optical center of the camera, and O ₁ represents the origin of the image coordinate system. X, Y represent the coordinate axes of the image coordinate system, and x, y, z represent the coordinate axes of the camera coordinate system. P(X _w , Y _w , Z _w ) represents the coordinates of a three-dimensional point in the camera coordinate system, and p represents the projection of point P on the image.

图2中的符号说明：O_w-X_wY_wZ_w为世界坐标系，Z_w轴为旋转轴，靶标在基准位置时的原点位于X_w轴上；靶标在基准位置时的靶标坐标系为该坐标系绕Z_w旋转θ_i后，靶标坐标系变为相机坐标系为O_c-X_cY_cZ_c。Explanation of the symbols in Figure 2: O _w -X _w Y _w Z _w is the world coordinate system, the Z _w axis is the rotation axis, the origin of the target is on the X _w axis when the target is at the reference position; the target coordinate system when the target is at the reference position for After the coordinate system rotates θ _i around Z _w , the target coordinate system becomes The camera coordinate system is O _c -X _c Y _c Z _c .

具体实施方式Detailed ways

见图1—图6，本发明提供了一种基于固定相机和单靶标的单轴旋转角的视觉测量方法，并进行了仿真实验和真实实验验证。As shown in Figures 1 to 6, the present invention provides a visual measurement method based on a fixed camera and a single-axis rotation angle of a single target, and has been verified by simulation experiments and real experiments.

仿真实验Simulation

仿真实验是在主频3.07GHz、内存4.00GB的计算机上，Windows XP环境下的MATLAB2014a上进行的。仿真所用的摄像机模型为二阶畸变模型，各参数如下：The simulation experiment is carried out on a computer with a main frequency of 3.07GHz and a memory of 4.00GB, and MATLAB2014a under Windows XP environment. The camera model used in the simulation is a second-order distortion model, and the parameters are as follows:

α_x＝782.5109,α_y＝782.9155,u₀＝357.3909,v₀＝264.5240α _x =782.5109, α _y =782.9155, u ₀ =357.3909, v ₀ =264.5240

k₁＝-0.4232,k₂＝0.2664,k₃＝-0.000555,k₄＝-0.0031,k₅＝0k ₁ =-0.4232, k ₂ =0.2664, k ₃ =-0.000555, k ₄ =-0.0031, k ₅ =0

靶标为6×8的棋盘格，单个棋盘格长度为39.5mm×39.5mm。仿真实验主要验证了标定算法的稳定性和高精度。仿真步骤如下：The target is a 6×8 checkerboard, and the length of a single checkerboard is 39.5mm×39.5mm. The simulation experiment mainly verifies the stability and high precision of the calibration algorithm. The simulation steps are as follows:

仿真环境设置：根据设定的内参数，及位姿，获取图像点坐标。其中转动时选取的角度为(单位是度)Simulation environment setting: According to the set internal parameters and pose, the image point coordinates are obtained. The angle selected when turning is (unit is degree)

θ_i＝2.25×i _θi = 2.25×i

步骤一：利用MATLAB标定工具箱，获取内参数的初值，结果如下(标定内参图像噪声为2像素)：Step 1: Use the MATLAB calibration toolbox to obtain the initial value of the internal parameters. The result is as follows (the image noise of the calibration internal reference is 2 pixels):

α_x＝774.34147,α_y＝773.91064,u₀＝349.58372,v₀＝269.85450α _x =774.34147, α _y =773.91064, u ₀ =349.58372, v ₀ =269.85450

k₁＝-0.43943,k₂＝0.33274,k₃＝-0.00222,k₄＝-0.00206,k₅＝0k ₁ =-0.43943, k ₂ =0.33274, k ₃ =-0.00222, k ₄ =-0.00206, k ₅ =0

步骤二：根据普通相机标定的结果，计算外参初值。Step 2: Calculate the initial value of the extrinsic parameters according to the result of ordinary camera calibration.

步骤三：根据内外参数初值，进行捆绑调整。Step 3: Perform bundle adjustment according to the initial value of internal and external parameters.

步骤四：选取角度，进行测量。Step 4: Select an angle and measure it.

为了验证算法的鲁棒性，在步骤四中我们选取100个角度，进行了测量，并将角度误差的最大值作为测量误差。在图像点中添加均值为μ＝0，标准差σ＝0□4个像素的高斯白噪声，取100次实验测量误差的均值和方差作图，结果如图3(a)、(b)所示。结果显示，即使加入了标准差为4个像素的高斯白噪声，仍能得到角度误差均值小于0.5度的结果，这证明我们的测角方法对图像噪声时鲁棒的。为了验证我们的方法对遮挡情形有效，我们对靶标在(-180°,180°)内转动时的情形，进行了仿真。图4(a)、(b)列出了仿真结果，(a)表示检测出的角点的个数，(b)表示测量的角度误差。实验结果表明，我们的方法可以有效解决遮挡问题。In order to verify the robustness of the algorithm, in step 4, we selected 100 angles for measurement, and took the maximum value of the angle error as the measurement error. Add Gaussian white noise with mean value μ=0 and standard deviation σ=0□4 pixels to the image points, and plot the mean value and variance of 100 experimental measurement errors, the results are shown in Figure 3(a) and (b) Show. The results show that even if Gaussian white noise with a standard deviation of 4 pixels is added, the mean value of the angle error is still less than 0.5 degrees, which proves that our angle measurement method is robust to image noise. In order to verify that our method is effective for the occlusion situation, we simulated the situation when the target rotates within (-180°, 180°). Figure 4(a), (b) lists the simulation results, (a) indicates the number of detected corner points, and (b) indicates the measured angle error. Experimental results show that our method can effectively solve the occlusion problem.

真实实验real experiment

为了进一步验证本发明的可行性，我们进行了真实实验验证。我们采用的靶标为打印的棋盘格(如图5所示)，规格为6×9，每个格子的长度为29.5mm，；所采用的相机为Basler相机，分辨率为658×492像素，焦距为800像素左右。角度验证采用光电跟踪系统，该系统角度测量精度为0.01度。结果如表1，误差表明在该实验条件下，角度误差的最大值小于0.15度。In order to further verify the feasibility of the present invention, we have carried out real experiment verification. The target we used is a printed checkerboard (as shown in Figure 5), the specification is 6×9, and the length of each grid is 29.5mm; the camera we use is a Basler camera with a resolution of 658×492 pixels and a focal length of It is about 800 pixels. The angle verification adopts the photoelectric tracking system, and the angle measurement accuracy of the system is 0.01 degree. The results are shown in Table 1. The error shows that under the experimental conditions, the maximum value of the angle error is less than 0.15 degrees.

表1.真实实验旋转角度估计结果Table 1. Real experimental rotation angle estimation results

11 22 33 44 55 66 77 88 99 1010 转台值turntable value -54.00-54.00 -36.00-36.00 -27.00-27.00 -10.80-10.80 -1.8-1.8 5.405.40 16.2016.20 28.8028.80 36.0036.00 54.0054.00 真值true value -53.985-53.985 -35.978-35.978 -26.966-26.966 -10.749-10.749 -0.036-0.036 5.4775.477 0.0700.070 28.89428.894 36.13536.135 54.14754.147 误差error -0.015-0.015 -0.022-0.022 -0.034-0.034 -0.051-0.051 -1.764-1.764 0.0770.077 16.27016.270 0.0940.094 0.1350.135 0.1470.147

Claims

1. A vision measurement method based on a single-axis rotation angle of a fixed camera and a single target is characterized in that: the method comprises the following specific steps:

the method comprises the following steps: calibrating internal and external parameters of a camera;

in the specific implementation process, a two-dimensional checkerboard is used as a target; the adopted picture sequence comprises a picture at a reference position and pictures when the m rotating components rotate to different angles; suppose that the m angles are θ₁,θ₂,...,θ_mNext, the checkerboard corner points of MATLAB2014 are utilized to extractTaking a function detectCheckerboardPoints for corner extraction, and combining an MATLAB calibration tool for internal parameter calibration to obtain an internal parameter p_inAnd external parameters of m +1 poses; the calibration tool box is downloaded on http:// www.vision.caltech.edu/bouguetj/caliib _ doc/; these extrinsic parameters are expressed as

[R_{{cb}_{i}}, T_{{cb}_{i}}], i = 0,1,2, . . ., m .

Expressed at a rotation angle of theta_iWhen the i is 0, the rotation matrix and the translation vector between the camera and the target are represented as the reference position;

step two: initial value of external parameter related to pivoting motion

(1) Determining the angle of rotation theta from (1.7) based on the external parameters determined in step one_iThe rotational and translational vectors of the target relative to the reference

[R_{b_{i} b_{0}}, T_{b_{i} b_{0}}], i = 1,2, . . ., m .

R_{b_{i} b_{0}} = R_{{cb}_{i}}^{T} R_{{cb}_{0}}, T_{b_{i} b_{0}} = R_{{cb}_{i}}^{T} (T_{{cb}_{0}} - T_{{cb}_{i}}) . - - - (0.1)

(2) And obtainingt and theta_iInitial value:

definition of r₁,r₂,r₃Are respectively as1,2,3, and vector [ t ]]_×Is expressed as the following matrix

Then, obtain

Wherein,S_i,jelements representing the ith row and j columns of the matrix S, r₁Obtained from the following formula

The above type passing pairObtained by SVD decomposition if

Then x is the third column of V, then, get

θ_i＝2arctan((A^TA)^-1(A^TB)), (0.6)

Wherein

A = vec (r_{2} r_{3}^{T} + r_{3} r_{2}^{T}),

Vec (-) here represents the straightening of the matrix, so,

R_{{wb}_{0}} = R_{b_{0} w}^{T} - - - (0.7)

(3) according to the obtainedt and theta_iSolving the initial value, and solving R according to (1.16)_cw,T_cwValue of (A)

R_{cw} = R_{{cb}_{0}} R_{{wb}_{0}}^{T}, T_{cw} = {(t, 0,0)}^{T} - R_{cw} T_{{wb}_{0}} . - - - (0.9)

Step three: intra-pair parameter p_inAnd an extrinsic parameter p_outCarrying out nonlinear optimization; where p is_in∈□⁹,p_out∈□¹⁰Is shown below

Wherein,are respectively asR_cwThe corresponding Rodrigues vector; the conversion relationship between a Rodrigues vector v and its corresponding rotation matrix R is as follows

Wherein, I₃An identity matrix of order 3;

the coordinates of the individual image points are marked on the target 1 at a rotation angle θ_iThen, the corner point coordinate corresponding to the jth target feature point is recorded as

(u_{j}^{i}, v_{j}^{i}), i = 0,1, . . ., m, j = 1,2, . . ., n,

Points on the target are noted

X_b＝{(x_j,y_j,0)^T,j＝1,2,...,n}.

Because errors exist in the image extraction process, the reprojection error sum of the images is selected as an optimization target, and the target point and the p_outObtaining the reprojection pointsThe process of (2) is as follows:

<1>angle of rotation theta_iTime, target point (x)_j,y_j,0)^TConversion to points in world coordinate system

<2>And calculating the corresponding reprojection point

The reprojection error is defined as

The optimization problem is described as

The optimization problem is solved by utilizing a sparse Levenberg-Marquardt algorithm, and an optimization tool box of the sparse LM algorithm is downloaded in http:// users. ics.forth.gr/. Loursakis/spark LM; with the above optimization, obtaint,R_cw,T_cwThe information of (1), which is used for the rotation angle calculation of the step four;

step four: measuring the rotation angle by utilizing the calibrated internal and external parameters;

firstly, matching angle points according to calibrated information; let the corner coordinate of the picture taken at θ be X^c＝{(u_i,v_i)^T,i＝1,2,...,N_cAt this time, occlusion may exist in the picture, so matching cannot be automatically completed according to the following algorithm:

scattering and searching a point set which is possibly matched;

for each corner point and target point pairCalculating a score S_i,jThe calculation method is as follows:

(1) obtaining (x)_j,y_j,0)^TCorresponding world coordinate

[\begin{matrix} x_{j}^{w} \\ y_{j}^{w} \\ z_{j}^{w} \end{matrix}] = R_{{wb}_{0}} [\begin{matrix} x_{j} \\ y_{j} \\ 0 \end{matrix}] + [\begin{matrix} t \\ 0 \\ 0 \end{matrix}] .

(2) Finding an image point (u)_i,v_i)^TIn thatPrimary image on planeThe obtaining method is as follows

Z_{w}^{j} = z_{w}^{j},

Wherein,

is (u)_i，v_i)^TThe corresponding normalized coordinates are obtained by iteration;

(3) definition of S_i,jThe following formula

S_{i, j} = \sqrt{{(X_{w}^{j})}^{2} + {(Y_{w}^{j})}^{2}} - \sqrt{{(x_{w}^{j})}^{2} + {(y_{w}^{j})}^{2}} . - - - (0.15)

(4) Selecting a threshold τ_mWhere d is the length of the checkerboard grid, then define the set of points that may match

<math> <mrow> <mi>Ω</mi> <mo>=</mo> <mo>{</mo> <msup> <mrow> <mo>(</mo> <msub> <mi>u</mi> <mi>i</mi> </msub> <mo>,</mo> <msub> <mi>v</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> <mi>T</mi> </msup> <mo>&LeftRightArrow;</mo> <msup> <mrow> <mo>(</mo> <msub> <mi>x</mi> <mi>j</mi> </msub> <mo>,</mo> <msub> <mi>y</mi> <mi>j</mi> </msub> <mo>,</mo> <mn>0</mn> <mo>)</mo> </mrow> <mi>T</mi> </msup> <mo>,</mo> <mi>j</mi> <mo>=</mo> <mn>1,2</mn> <mo>,</mo> <mo>.</mo> <mo>.</mo> <mo>.</mo> <mo>,</mo> <mi>n</mi> <mo>,</mo> <mi>i</mi> <mo>=</mo> <mn>1,2</mn> <mo>,</mo> <mo>.</mo> <mo>.</mo> <mo>.</mo> <mo>,</mo> <msub> <mi>N</mi> <mi>c</mi> </msub> <mo>|</mo> <msub> <mi>S</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> </msub> <mo><</mo> <msub> <mi>τ</mi> <mi>m</mi> </msub> <mo>}</mo> <mo>.</mo> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>0.16</mn> <mo>)</mo> </mrow> </mrow> </math>

Utilizing the histogram to calculate an angle;

for the elements in Ω, θ is calculated by the following equation_i,j

Then draw outA histogram of distribution having 360 squares, the number of elements in the kth square being (k-181) °, (k-179) ° n Ω; an estimate of theta of

Selecting a correct matching pair by using the obtained angle;

i corresponds to the point M (j) with the smallest score in the correct angle range, i.e. the point M (j)

M (i) = \min_{j} S_{i, j},

Satisfies theta_i,j∈[θ-1°,θ+1°].

If theta is not present_i,j∈[θ-1°,θ+1°]Then, let m (i) be 0;

good and non-linear optimization

The optimization problem is described as

Wherein,

the optimization problem is solved by the lsqnolin function of the MATLAB itself.