基于双目显微立体视觉的零件精确定位方法

著录项
  • CN201310182221.0
  • 20130516
  • CN103247053A
  • 20130814
  • 大连理工大学
  • 刘巍;贾振元;屠先明;王福吉;王文强;赵凯
  • G06T7/00
  • G06T7/00

  • 辽宁省大连市凌工路2号
  • 中国,CN,辽宁(21)
  • 大连理工大学专利中心
  • 关慧贞
摘要
本发明基于双目显微立体视觉的零件精确定位方法属于计算机视觉测量技术领域,涉及一种基于双目显微立体视觉的精密零件的精确定位方法。采用双目显微立体视觉系统,利用两个CCD摄像机采集被测零件的图像,立体显微镜将被测零件上的待测区域的图像信息放大;采用棋盘格标定板对两个CCD摄像机进行标定;采用Harris角点检测算法及亚像素提取算法进行特征点的提取。将提取后的特征点进行初匹配和匹配点对的校正,将特征点图像坐标输入到标定好的系统中得到特征点的空间实际坐标。本发明解决了由于待测目标区域小、定位精度要求高、非接触等问题所产生的测量难题。采用基于双目显微立体视觉的非接触式测量方法,很好的完成了精密零件的精确定位。
权利要求

1.一种基于双目显微立体视觉的零件精确定位方法,其特征是,采用双目显微 立体视觉系统,利用左右两个CCD摄像机(2、2’)采集被测零件(4)的图像 位置信息,立体显微镜(3)将被测零件(4)上的待测区域(6)的图像位置信 息进行放大;放大后的图像经过图像采集卡传到工业控制计算机(9)内,采用 精度较高的棋盘格标定板对左右两个CCD摄像机进行标定;采用Harris角点检 测算法及亚像素提取算法进行待测特征点的提取;将提取后的待测特征点进行 初匹配和匹配点对的校正;将已经匹配好的特征点图像坐标输入到标定好的系 统中得到特征点的空间实际坐标;测量方法的具体步骤如下:

(1)左右两个CCD摄像机的标定

左右两个CCD摄像机的标定包括摄像机内参数和摄像机外参数;摄像机内 参数包括尺度因子α、β,主点坐标u 0,v 0,以及垂直因子γ;在标定过程中, 需要先得到摄像机的五个内参数,在求得内参数矩阵的基础上,求解摄像机外 参数矩阵;尺度因子是空间点在经过平移旋转变换之后,它在摄像机坐标系中 的坐标与它在图像坐标系中的坐标的尺度关系,本发明对尺度因子的标定采用 张氏摄像机标定法;根据摄像机的针孔模型可建立图像坐标系与世界坐标系的 关系公式为:

Z c u v 1 = α γ u 0 0 β v 0 0 0 1 R t X w Y w Z w 1 - - - ( 1 )

其中:α和β就是需要标定的尺度因子,X w,Y w,Z w是空间中一点P的三维空间 坐标,u,v是P点在图像上的图像坐标,R,t代表了摄像机坐标系相对于世界坐 标系的旋转和平移矩阵;若令X′ w=(X w,Y w,Z w) T, 则

x w = ( u , v ) T , Z c = s ,

其中,H也称为单应性矩阵,且有

H=K[r 1,r 2,t]                 (3)

每一个位置的标定板图像都独立对应一个单应性矩阵,令H=[h 1,h 2,h 3]则可由 上式推得:

[h 1  h 2  h 3]=λK[r 1  r 2  t]           (4)

其中,λ为任意的比例因子,K摄像机内参数矩阵,因r 1和r 2是单位正交向量, 即r 1 Tr 1=r 2 Tr 2=1且r 1 Tr 2=0则:

h 1 T K - T K - 1 h 1 = h 2 T K - T K - 1 h 2 = 1 - - - ( 5 )

h 1K ‑TK ‑1h 2=0                   (6)

根据公式(5)公式(6),计算出立体显微镜左右光路的摄像机尺度因子;针对主点 坐标,采用变倍率法对摄像机的主点进行标定;假设空间中任意一点P在摄像 机坐标系中的坐标为x 1,y 1,z 1,在不考虑非线性畸变以及图像坐标垂直度的情况 下,该点的图像平面投影的坐标为:

u 1 = rx 1 + u 0 v 1 = ry 1 + v 0 - - - ( 7 )

其中u 1,v 1是该点的图像坐标,r是任一放大倍率,u 0,v 0是主点坐标;将式中的r 消去,可得直线方程:

u 1 - u 0 x 1 = v 1 - v 0 y 1 - - - ( 8 )

即在任何放大倍率下,点x 1,y 1,z 1的图像坐标都在同一直线上,且该直线一定经 过摄像机主点u 0,v 0;为了求取摄像机主点u 0,v 0,选定12个在图像投影平面上 不重合的点,令这12个点分别在不同放大倍率下0.7、1、2、2.4、3、5、5.8、 8进行投影,利用最小二乘法对获得的12条直线进行拟合,并用最小二乘法求 取12条直线的交点坐标即摄像机主点u 0,v 0;由于垂直因子γ对成像精度影响 不大,因此在初始参数标定时,将γ设为0,通过优化手段得到γ具体值;针对 摄像机外参数旋转矩阵R和平移矩阵t;由公式(2)知,在标定完内参数矩阵K 和单应性矩阵H后可以求得摄像机旋转矩阵R=[r 1,r 2,r 3]和平移向量t如下:

r 1=λK ‑1h 1

r 2=λK ‑1h 2                     (9)

r 3=r 1×r 2

t=λK ‑1h 3

在初步标定完左右两个CCD摄像机的主参数后,下一步是对主参数和畸变参数 的优化;为了简化优化参数,本发明采用四元数法将旋转矩阵 R = [ r 1 , r 2 , r 3 ] = r 11 r 21 r 31 r 12 r 22 r 32 r 13 r 23 r 33 中的九个未知数缩减为四个,公式如下:

q 1 = 1 + r 11 + r 22 + r 33 2 , q 2 = r 23 - r 32 4 q 1 , q 3 = r 31 - r 13 4 q 1 , q 4 = r 12 - r 21 4 q 1

q 2 = 1 + r 11 - r 22 - r 33 2 , q 1 = r 23 - r 32 4 q 2 , q 3 = r 12 + r 21 4 q 2 , q 4 r 31 - r 13 4 q 2        (10)

q 3 = 1 - r 11 + r 22 - r 33 2 , q 1 = r 31 - r 13 4 q 3 , q 2 = r 12 - r 21 4 q 3 , q 4 = r 23 - r 32 4 q 3

q 4 = 1 - r 11 - r 22 + r 33 2 , q 1 = r 12 - r 21 4 q 4 , q 2 = r 31 + r 13 4 q 4 , q 3 = r 23 - r 32 4 q 4

从公式(10)计算得到新的旋转矩阵:

R = - 1 + 2 q 2 2 + 2 q 1 2 2 q 2 q 3 - 2 q 1 q 4 2 q 1 q 3 + 2 q 2 q 4 2 q 2 q 3 + 2 q 4 q 1 - 1 + 2 q 1 2 + 2 q 3 2 - 2 q 1 q 2 + 2 q 3 q 4 - 2 q 1 q 3 + 2 q 4 q 2 2 q 1 q 2 + 2 q 3 q 4 - 1 + 2 q 1 2 + 2 q 4 3 - - - ( 11 )

由于四元数q 1,q 2,q 3,q 4的平方和等于一,因此旋转矩阵的待优化参数简化为三 个;将已求得的左右两个CCD摄像机标定参数值作为优化初始值,利用优化算 法对摄像机主参数和畸变参数进行整体优化;对左右两个CCD摄像机主参数和 畸变参数的优化是基于最大似然估计准则,对于给定的n个标定板图片提供的 m×n个标定点坐标,左右两个CCD摄像机参数的最优化问题通过下式的最小化 问题来表达:

S ( θ ) = Σ i = 1 n Σ i = 1 m ϵ i 2 = Σ i = 1 n Σ i = 1 m [ y 1 - m i ( C j , p j , X i ) ] 2 - - - ( 12 )

其中j是指参与计算的第j个摄像机,i是指第j个摄像机获取的第i个点,X i是 输入的空间点坐标,y i是空间第i点的图像坐标,C j是固定不变的摄像机参数向 量,其长度为n 0,p j是需要调整的摄像机参数向量,其长度为n 1,n 0+n 1即为摄 像机所有参数向量长度,m i(C j,p j,X i)是摄像机的成像方程;本发明采用基于 LM算法的光束平差法来解决上式的最小化优化问题;利用光束平差法的编写优 化程序时除了要提供左右两个CCD摄像机的标定参数初值外还要提供一定数量 的和形式的三维空间坐标点以及这些三维坐标点的相应图像坐标,因此为了精 确得到整个视场空间内高精密标定板的角点坐标,采用数控机床z轴进给带动 显微镜纵向移动的方法;并采用激光干涉仪(测量精度可达到0.1um)对z轴进 行测量得到z轴移动的实际距离,保证所采集角点坐标点的z向实际精度;在光 束平差法中,优化算法的迭代因子是由迭代参数的一阶导数组成的雅克比矩阵 构成的,因此需要在计算过程中引入待成像方程中调整参数的雅克比矩阵;

J = [ df dp ( 1 ) , df dp ( 2 ) , . . . , df dp ( n 1 ) ] - - - ( 13 )

光束平差法以其待调整参数的雅克比矩阵组成迭代因子,对模型参数进行调整, 得到使成像误差S(θ)达到最小的模型参数最优解;

(2)待测区域特征点的提取

本发明针对待测区域特征点的提取采用的是Harris角点检测算法;Harris 提出了用高斯函数代替方形区域计算梯度的算法,若目标点的像素坐标为x,y, 其在x方向和y方向的位移分别为u和v,则该目标点的灰度变化由如下公式表 示:

E ( x , y ) = Σ u , v ω ( u , v ) [ I ( x + u , y + v ) - I ( x , y ) ] 2

Σ u , v ω ( u , v ) ( u , v ) I x 2 I x I y I x I y I y 2 u v - - - ( 14 )

M ( x , y ) = Σ u , v ω ( u , v ) I x 2 I x I y I x I y I y 2 则:

E ( x , y ) ( u , v ) M ( x , y ) u v - - - ( 15 )

计算矩阵M的特征值λ 1,λ 2,如果两个特征值都比较大,说明目标点的图像灰度 自相关函数的两个正交方向上的值均较大,则该点即为特征点;通过Harris算 法对特征点进行初提取后;采用亚像素角点提取算法提取出精度更高的角点坐 标;对于理想角点,它附近的像素点的灰度梯度方向均垂直于该点与理想角点 的连线;这一特点可以用公式表达为:

H ( α β ) = 0 - - - ( 16 )

其中 是理想角点的灰度梯度方向, 向量为图像原点指向理想角点的坐标, 向量为图像原点指向理想角点附近任一边缘点的坐标;实际上,由于受到图像 噪声的影响,公式(10)通常都不为零,可以将其视为误差e,即

e = H ( α - β ) - - - ( 17 )

在以角点为中心的邻域内,对所有点按上式计算,误差和为E,则

E = Σ i e = Σ i H ( α - β ) - - - ( 18 )

其中,i是邻域内第i个点;这样求使误差和E最小的点即是亚像素级别的角点 坐标;

(3)特征点的立体匹配

采用归一化交叉算法(NCC)实现特征点的初匹配,以CCD摄像机2拍得 的图片作为基准图片,并以特征点p l为中心,构造一个N*N(代表图像像素)的 局部图像块作为模版图T,并令模版图T在待匹配图像S的包含有特征点的一定 范围内进行遍历,搜索模版所覆盖的子图记做S i,j,i,j为子图中心点在匹配图 像S中的像素坐标,通过使用归一化交叉算法(NCC)来计算模版图T和搜索 子图S i,j之间的归一化互相关系数,算法如下:

R ( i , j ) = Σ m , n | S i , j ( m , n ) - E ( S i , j ) | | T ( m , n ) - E ( T ) | Σ m , n ( S i , j ( m , n ) - E ( S i , j ) ) 2 Σ m , n ( T ( m , n ) - E ( T ) ) 2 - - - ( 19 )

其中m,n代表各个像素的坐标,E(S i,j)和E(T)分别是搜索子图S i,j与模版图T 的灰度平均值;互相关系数R(i,j)的值越大,则搜索子图与模板图的匹配程度 越高;设定互相关系数R(i,j)的阈值后,选出特征点p l的候选匹配点;使用NCC 算法进行初匹配后还需要利用外极限约束和距离约束来对匹配的特征点对进行 校正;通过使用Longguet‑Higgins提出的归一化8点算法在摄像机标定时计算左 右两个CCD摄像机的基本矩阵F,根据外极线理论,对基准图像特征点p l,在 待匹配图像S中对应的外极线l r可以表示为:

l r=F·p l                       (20)

初匹配的候选特征点若在外极线l r附近则可进一步确定该候选特征点为匹配点; 最后统计每一候选特征点到其他特征点的距离和基准图像中特征点到其他特征 点点距离相等的个数num;num最大值时对应的候选特征点则为基准图像中特 征点的正确匹配点;

(4)三维坐标的求取

首先计算出双目图像中特征点经亚像素算法提取出的图像坐标值,并大致 估算出特征点在世界坐标系中的估计值;在优化程序中将左右两个CCD摄像机 参数对应雅克比矩阵设为零,同时在优化程序计算迭代因子部分中嵌入关于空 间坐标x、y、z的雅克比矩阵方程;将已经优化好的左右两个CCD摄像机参数 以及特征点点图像坐标值和世界坐标估计值输入光束平差法优化程序,就可以 对空间点的三维坐标进行优化,得到特征点的三维空间测量值,从而完成三维 坐标点的求取并最终实现特征点的定位。

说明书
技术领域

本发明属于计算机视觉测量技术领域,特别涉及一种基于双目显微立体视 觉的精密零件的精确定位方法。

双目立体视觉技术通过模拟人类双眼的处理方式,具有获取待测物体深度 信息的能力,进而能够获得被测物体的空间位置信息,同时又具有无损测量及 实时测量的优点,在各行各业得到了广泛的应用。双目显微立体视觉技术是建 立在双目立体视觉技术和体式显微镜技术之上的,通过左右两个CCD摄像机分 别捕捉立体显微镜两条光路的图像,通过捕捉的图像存在同双目视觉测量时类 似的角度差,能够实现对目标的三维测量。双目显微立体视觉技术在保证非接 触性、实时性等优点的同时能够精确地获得目标点的空间位置信息。其可应用 于三维高精度伺服控制、高精度空间定位以及高精度三维测量等方面。零件的 加工定位精度直接影响着零件的加工精度,而对于精密零件定位精度有可能直 接影响加工零件是否可用,同时针对精密超精密零件的加工定位往往要求非接 触的特点,基于双目显微立体视觉的精密零件的精确定位方法应运而生。

重庆大学的汤宝平等人申请的发明专利CN102567989A“基于双目立体视觉 的空间定位方法”针对空间定位存在的摄像机标定过程繁琐、定位精度不高问题 提出了一种新的基于双目立体视觉的空间定位方法,该方法通过对摄像机进行 标定求得摄像机的单应性矩阵和畸变系数,再通过标定结果获得目标点的理论 坐标,最后结合目标点的实际坐标计算出目标点的空间坐标。然而,由于双目 显微立体视觉系统具有更加复杂的光路,焦深小,畸变因素较多,视场窄,而 目前广泛采用的摄像机标定方法都是建立在大视场下实现的,因此上述方法并 不能解决基于双目显微立体视觉的精密零件的精确定位。

本发明要解决的技术难题是克服现有技术的缺陷,发明一种基于双目显微 立体视觉的精密零件的精确定位方法,解决由于待测目标区域小、定位精度要 求高、非接触等问题所产生的测量难题。采用基于双目显微立体视觉来实现精 密零件的精确定位,拓宽了传统的基于双目立体视觉的空间定位方法的应用范 围,并提高了测量的精度,解决了待测零件区域小,定位精度要求高,难以测 量的问题。

本发明所采用的技术方案是基于双目显微立体视觉的精密零件的精确定位 方法,其特征是,采用双目显微立体视觉系统,利用左右两个CCD摄像机2、2’ 采集被测零件4的图像位置信息,立体显微镜3将被测零件4上的待测区域6 的图像位置信息进行放大;放大后的图像经过图像采集卡传到工业控制计算机9 内,采用精度较高的棋盘格标定板对左右两个CCD摄像机进行标定;采用Harris 角点检测算法及亚像素提取算法进行待测特征点的提取;将提取后的待测特征 点进行初匹配和匹配点对的校正;将已经匹配好的特征点图像坐标输入到标定 好的系统中得到特征点的空间实际坐标;测量方法的具体步骤如下:

(1)左右两个CCD摄像机的标定

左右两个CCD摄像机的标定包括摄像机内参数和摄像机外参数;摄像机内 参数包括尺度因子α、β,主点坐标u0,v0,以及垂直因子γ;在标定过程中, 需要先得到摄像机的五个内参数,在求得内参数矩阵的基础上,求解摄像机外 参数矩阵;尺度因子是空间点在经过平移旋转变换之后,它在摄像机坐标系中 的坐标与它在图像坐标系中的坐标的尺度关系,本发明对尺度因子的标定采用 张氏摄像机标定法;根据摄像机的针孔模型可建立图像坐标系与世界坐标系的 关系公式为:

<mrow> <msub> <mi>Z</mi> <mi>c</mi> </msub> <mfenced open='[' close=']'> <mtable> <mtr> <mtd> <mi>u</mi> </mtd> </mtr> <mtr> <mtd> <mi>v</mi> </mtd> </mtr> <mtr> <mtd> <mn>1</mn> </mtd> </mtr> </mtable> </mfenced> <mo>=</mo> <mfenced open='[' close=']'> <mtable> <mtr> <mtd> <mi>&alpha;</mi> </mtd> <mtd> <mi>&gamma;</mi> </mtd> <mtd> <msub> <mi>u</mi> <mn>0</mn> </msub> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mi>&beta;</mi> </mtd> <mtd> <msub> <mi>v</mi> <mn>0</mn> </msub> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> </mtr> </mtable> </mfenced> <mfenced open='[' close=']'> <mtable> <mtr> <mtd> <mi>R</mi> </mtd> <mtd> <mi>t</mi> </mtd> </mtr> </mtable> </mfenced> <mfenced open='[' close=']'> <mtable> <mtr> <mtd> <msub> <mi>X</mi> <mi>w</mi> </msub> </mtd> </mtr> <mtr> <mtd> <msub> <mi>Y</mi> <mi>w</mi> </msub> </mtd> </mtr> <mtr> <mtd> <msub> <mi>Z</mi> <mi>w</mi> </msub> </mtd> </mtr> <mtr> <mtd> <mn>1</mn> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </mrow>

其中:α和β就是需要标定的尺度因子,Xw,Yw,Zw是空间中一点P的三维空间
坐标,u,v是P点在图像上的图像坐标,R,t代表了摄像机坐标系相对于世界坐
标系的旋转和平移矩阵;若令X′w=(Xw,Yw,Zw)T,Zc=s,则

<mrow> <mi>s</mi> <msub> <mover> <mi>x</mi> <mo>&OverBar;</mo> </mover> <mi>w</mi> </msub> <mo>=</mo> <mi>H</mi> <msubsup> <mi>X</mi> <mi>w</mi> <mo>&prime;</mo> </msubsup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </mrow>

其中,H也称为单应性矩阵,且有

H=K[r1,r2,t]          (3)

每一个位置的标定板图像都独立对应一个单应性矩阵,令H=[h1,h2,h3]则可由 上式推得:

[h1  h2  h3]=λK[r1  r2  t]      (4)

其中,λ为任意的比例因子,K摄像机内参数矩阵,因r1和r2是单位正交向量, 即r1Tr1=r2Tr2=1且r1Tr2=0则:

<mrow> <msubsup> <mi>h</mi> <mn>1</mn> <mi>T</mi> </msubsup> <msup> <mi>K</mi> <mrow> <mo>-</mo> <mi>T</mi> </mrow> </msup> <msup> <mi>K</mi> <mrow> <mo>-</mo> <mn>1</mn> </mrow> </msup> <msub> <mi>h</mi> <mn>1</mn> </msub> <mo>=</mo> <msubsup> <mi>h</mi> <mn>2</mn> <mi>T</mi> </msubsup> <msup> <mi>K</mi> <mrow> <mo>-</mo> <mi>T</mi> </mrow> </msup> <msup> <mi>K</mi> <mrow> <mo>-</mo> <mn>1</mn> </mrow> </msup> <msub> <mi>h</mi> <mn>2</mn> </msub> <mo>=</mo> <mn>1</mn> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>5</mn> <mo>)</mo> </mrow> </mrow>

h1K‑TK‑1h2=0        (6)

根据公式5公式6,计算出立体显微镜左右光路的摄像机尺度因子;针对主点坐 标,采用变倍率法对摄像机的主点进行标定;假设空间中任意一点P在摄像机 坐标系中的坐标为x1,y1,z1,在不考虑非线性畸变以及图像坐标垂直度的情况 下,该点的图像平面投影的坐标为:

<mfenced open='' close='' separators=''> <mtable> <mtr> <mtd> <msub> <mi>u</mi> <mn>1</mn> </msub> <mo>=</mo> <msub> <mi>rx</mi> <mn>1</mn> </msub> <mo>+</mo> <msub> <mi>u</mi> <mn>0</mn> </msub> </mtd> </mtr> <mtr> <mtd> <msub> <mi>v</mi> <mn>1</mn> </msub> <mo>=</mo> <mi>r</mi> <msub> <mi>y</mi> <mn>1</mn> </msub> <mo>+</mo> <msub> <mi>v</mi> <mn>0</mn> </msub> </mtd> </mtr> </mtable> <mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>7</mn> <mo>)</mo> </mrow> </mrow> </mfenced>

其中u1,v1是该点的图像坐标,r是任一放大倍率,u0,v0是主点坐标;将式中的r 消去,可得直线方程:

<mrow> <mfrac> <mrow> <msub> <mi>u</mi> <mn>1</mn> </msub> <mo>-</mo> <msub> <mi>u</mi> <mn>0</mn> </msub> </mrow> <msub> <mi>x</mi> <mn>1</mn> </msub> </mfrac> <mo>=</mo> <mfrac> <mrow> <msub> <mi>v</mi> <mn>1</mn> </msub> <mo>-</mo> <msub> <mi>v</mi> <mn>0</mn> </msub> </mrow> <msub> <mi>y</mi> <mn>1</mn> </msub> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>8</mn> <mo>)</mo> </mrow> </mrow>

即在任何放大倍率下,点x1,y1,z1的图像坐标都在同一直线上,且该直线一定经 过摄像机主点u0,v0;为了求取摄像机主点u0,v0,选定12个在图像投影平面上 不重合的点,令这12个点分别在不同放大倍率下0.7、1、2、2.4、3、5、5.8、 8进行投影,利用最小二乘法对获得的12条直线进行拟合,并用最小二乘法求 取12条直线的交点坐标即摄像机主点u0,v0;由于垂直因子γ对成像精度影响 不大,因此在初始参数标定时,将γ设为0,通过优化手段得到γ具体值;针对 摄像机外参数旋转矩阵R和平移矩阵t;由公式2知,在标定完内参数矩阵K和 单应性矩阵H后可以求得摄像机旋转矩阵R=[r1,r2,r3]和平移向量t如下:

r1=λK‑1h1

r2=λK‑1h2                              (9)

r3=r1×r2

t=λK‑1h3

在初步标定完左右两个CCD摄像机的主参数后,下一步是对主参数和畸变参数
的优化;为了简化优化参数,本发明采用四元数法将旋转矩阵
<mrow> <mi>R</mi> <mo>=</mo> <mo>[</mo> <msub> <mi>r</mi> <mn>1</mn> </msub> <mo>,</mo> <msub> <mi>r</mi> <mn>2</mn> </msub> <mo>,</mo> <msub> <mi>r</mi> <mn>3</mn> </msub> <mi></mi> <mo>]</mo> <mo>=</mo> <mfenced open='[' close=']'> <mtable> <mtr> <mtd> <msub> <mi>r</mi> <mn>11</mn> </msub> </mtd> <mtd> <msub> <mi>r</mi> <mn>21</mn> </msub> </mtd> <mtd> <msub> <mi>r</mi> <mn>31</mn> </msub> </mtd> </mtr> <mtr> <mtd> <msub> <mi>r</mi> <mn>12</mn> </msub> </mtd> <mtd> <msub> <mi>r</mi> <mn>22</mn> </msub> </mtd> <mtd> <msub> <mi>r</mi> <mn>32</mn> </msub> </mtd> </mtr> <mtr> <mtd> <msub> <mi>r</mi> <mn>13</mn> </msub> </mtd> <mtd> <msub> <mi>r</mi> <mn>23</mn> </msub> </mtd> <mtd> <msub> <mi>r</mi> <mn>33</mn> </msub> </mtd> </mtr> </mtable> </mfenced> </mrow> 中的九个未知数缩减为四个,公式如下:

<mrow> <msub> <mi>q</mi> <mn>1</mn> </msub> <mo>=</mo> <mfrac> <msqrt> <mn>1</mn> <mo>+</mo> <msub> <mi>r</mi> <mn>11</mn> </msub> <mo>+</mo> <msub> <mi>r</mi> <mn>22</mn> </msub> <mo>+</mo> <msub> <mi>r</mi> <mn>33</mn> </msub> </msqrt> <mn>2</mn> </mfrac> <mo>,</mo> <msub> <mi>q</mi> <mn>2</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>r</mi> <mn>23</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>32</mn> </msub> </mrow> <msub> <mrow> <mn>4</mn> <mi>q</mi> </mrow> <mn>1</mn> </msub> </mfrac> <mo>,</mo> <msub> <mi>q</mi> <mn>3</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>r</mi> <mn>31</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>13</mn> </msub> </mrow> <msub> <mrow> <mn>4</mn> <mi>q</mi> </mrow> <mn>1</mn> </msub> </mfrac> <mo>,</mo> <msub> <mi>q</mi> <mn>4</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>r</mi> <mn>12</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>21</mn> </msub> </mrow> <msub> <mrow> <mn>4</mn> <mi>q</mi> </mrow> <mn>1</mn> </msub> </mfrac> </mrow>

<mrow> <msub> <mi>q</mi> <mn>2</mn> </msub> <mo>=</mo> <mfrac> <msqrt> <mn>1</mn> <mo>+</mo> <msub> <mi>r</mi> <mn>11</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>22</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>33</mn> </msub> </msqrt> <mn>2</mn> </mfrac> <mo>,</mo> <msub> <mi>q</mi> <mn>1</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>r</mi> <mn>23</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>32</mn> </msub> </mrow> <msub> <mrow> <mn>4</mn> <mi>q</mi> </mrow> <mn>2</mn> </msub> </mfrac> <mo>,</mo> <msub> <mi>q</mi> <mn>3</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>r</mi> <mn>12</mn> </msub> <mo>+</mo> <msub> <mi>r</mi> <mn>21</mn> </msub> </mrow> <msub> <mrow> <mn>4</mn> <mi>q</mi> </mrow> <mn>2</mn> </msub> </mfrac> <mo>,</mo> <msub> <mi>q</mi> <mn>4</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>r</mi> <mn>31</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>13</mn> </msub> </mrow> <msub> <mrow> <mn>4</mn> <mi>q</mi> </mrow> <mn>2</mn> </msub> </mfrac> </mrow>      (10)

<mrow> <msub> <mi>q</mi> <mn>3</mn> </msub> <mo>=</mo> <mfrac> <msqrt> <mn>1</mn> <mo>-</mo> <msub> <mi>r</mi> <mn>11</mn> </msub> <mo>+</mo> <msub> <mi>r</mi> <mn>22</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>33</mn> </msub> </msqrt> <mn>2</mn> </mfrac> <mo>,</mo> <msub> <mi>q</mi> <mn>1</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>r</mi> <mn>31</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>13</mn> </msub> </mrow> <msub> <mrow> <mn>4</mn> <mi>q</mi> </mrow> <mn>3</mn> </msub> </mfrac> <mo>,</mo> <msub> <mi>q</mi> <mn>2</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>r</mi> <mn>12</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>21</mn> </msub> </mrow> <msub> <mrow> <mn>4</mn> <mi>q</mi> </mrow> <mn>3</mn> </msub> </mfrac> <mo>,</mo> <msub> <mi>q</mi> <mn>4</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>r</mi> <mn>23</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>32</mn> </msub> </mrow> <msub> <mrow> <mn>4</mn> <mi>q</mi> </mrow> <mn>3</mn> </msub> </mfrac> </mrow>

<mrow> <msub> <mi>q</mi> <mn>4</mn> </msub> <mo>=</mo> <mfrac> <msqrt> <mn>1</mn> <mo>-</mo> <msub> <mi>r</mi> <mn>11</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>22</mn> </msub> <mo>+</mo> <msub> <mi>r</mi> <mn>33</mn> </msub> </msqrt> <mn>2</mn> </mfrac> <mo>,</mo> <msub> <mi>q</mi> <mn>1</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>r</mi> <mn>12</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>21</mn> </msub> </mrow> <msub> <mrow> <mn>4</mn> <mi>q</mi> </mrow> <mn>4</mn> </msub> </mfrac> <mo>,</mo> <msub> <mi>q</mi> <mn>2</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>r</mi> <mn>31</mn> </msub> <mo>+</mo> <msub> <mi>r</mi> <mn>13</mn> </msub> </mrow> <msub> <mrow> <mn>4</mn> <mi>q</mi> </mrow> <mn>4</mn> </msub> </mfrac> <mo>,</mo> <msub> <mi>q</mi> <mn>3</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>r</mi> <mn>23</mn> </msub> <mo>-</mo> <msub> <mi>r</mi> <mn>32</mn> </msub> </mrow> <msub> <mrow> <mn>4</mn> <mi>q</mi> </mrow> <mn>4</mn> </msub> </mfrac> </mrow>

从公式(10)计算得到新的旋转矩阵:

<mrow> <mi>R</mi> <mo>=</mo> <mfenced open='[' close=']'> <mtable> <mtr> <mtd> <mo>-</mo> <mn>1</mn> <mo>+</mo> <msubsup> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>2</mn> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>1</mn> <mn>2</mn> </msubsup> </mtd> <mtd> <msub> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>2</mn> </msub> <msub> <mi>q</mi> <mn>3</mn> </msub> <mo>-</mo> <msub> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>1</mn> </msub> <msub> <mi>q</mi> <mn>4</mn> </msub> </mtd> <mtd> <msub> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>1</mn> </msub> <msub> <mi>q</mi> <mn>3</mn> </msub> <mo>+</mo> <msub> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>2</mn> </msub> <msub> <mi>q</mi> <mn>4</mn> </msub> </mtd> </mtr> <mtr> <mtd> <msub> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>2</mn> </msub> <msub> <mi>q</mi> <mn>3</mn> </msub> <mo>+</mo> <msub> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>4</mn> </msub> <msub> <mi>q</mi> <mn>1</mn> </msub> </mtd> <mtd> <mo>-</mo> <mn>1</mn> <mo>+</mo> <msubsup> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>1</mn> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>3</mn> <mn>2</mn> </msubsup> </mtd> <mtd> <mo>-</mo> <msub> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>1</mn> </msub> <msub> <mi>q</mi> <mn>2</mn> </msub> <mo>+</mo> <msub> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>3</mn> </msub> <msub> <mi>q</mi> <mn>4</mn> </msub> </mtd> </mtr> <mtr> <mtd> <mo>-</mo> <msub> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>1</mn> </msub> <msub> <mi>q</mi> <mn>3</mn> </msub> <mo>+</mo> <msub> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>4</mn> </msub> <msub> <mi>q</mi> <mn>2</mn> </msub> </mtd> <mtd> <msub> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>1</mn> </msub> <msub> <mi>q</mi> <mn>2</mn> </msub> <mo>+</mo> <msub> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>3</mn> </msub> <msub> <mi>q</mi> <mn>4</mn> </msub> </mtd> <mtd> <mo>-</mo> <mn>1</mn> <mo>+</mo> <msubsup> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>1</mn> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mrow> <mn>2</mn> <mi>q</mi> </mrow> <mn>4</mn> <mn>2</mn> </msubsup> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>11</mn> <mo>)</mo> </mrow> </mrow>

由于四元数q1,q2,q3,q4的平方和等于一,因此旋转矩阵的待优化参数简化为三 个;将已求得的左右两个CCD摄像机标定参数值作为优化初始值,利用优化算 法对摄像机主参数和畸变参数进行整体优化;对左右两个CCD摄像机主参数和 畸变参数的优化是基于最大似然估计准则,对于给定的n个标定板图片提供的 m×n个标定点坐标,左右两个CCD摄像机参数的最优化问题通过下式的最小化 问题来表达:

<mrow> <mi>S</mi> <mrow> <mo>(</mo> <mi>&theta;</mi> <mo>)</mo> </mrow> <mo>=</mo> <munderover> <mi>&Sigma;</mi> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>n</mi> </munderover> <munderover> <mi>&Sigma;</mi> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>m</mi> </munderover> <msubsup> <mi>&epsiv;</mi> <mi>i</mi> <mn>2</mn> </msubsup> <mo>=</mo> <munderover> <mi>&Sigma;</mi> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>n</mi> </munderover> <munderover> <mi>&Sigma;</mi> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>m</mi> </munderover> <msup> <mrow> <mo>[</mo> <msub> <mi>y</mi> <mi>i</mi> </msub> <mo>-</mo> <msub> <mi>m</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <msub> <mi>C</mi> <mi>j</mi> </msub> <mo>,</mo> <msub> <mi>p</mi> <mi>j</mi> </msub> <mo>,</mo> <msub> <mi>X</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> <mo>]</mo> </mrow> <mn>2</mn> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>12</mn> <mo>)</mo> </mrow> </mrow>

其中j是指参与计算的第j个摄像机,i是指第j个摄像机获取的第i个点,Xi是 输入的空间点坐标,yi是空间第i点的图像坐标,Cj是固定不变的摄像机参数向 量,其长度为n0,pj是需要调整的摄像机参数向量,其长度为n1,n0+n1即为摄 像机所有参数向量长度,mi(Cj,pj,Xi)是摄像机的成像方程;本发明采用基于 LM算法的光束平差法来解决上式的最小化优化问题;利用光束平差法的编写优 化程序时除了要提供左右两个CCD摄像机的标定参数初值外还要提供一定数量 的和形式的三维空间坐标点以及这些三维坐标点的相应图像坐标,因此为了精 确得到整个视场空间内高精密标定板的角点坐标,采用数控机床z轴进给带动 显微镜纵向移动的方法;并采用激光干涉仪(测量精度可达到0.1um)对z轴进 行测量得到z轴移动的实际距离,保证所采集角点坐标点的z向实际精度;在光 束平差法中,优化算法的迭代因子是由迭代参数的一阶导数组成的雅克比矩阵 构成的,因此需要在计算过程中引入待成像方程中调整参数的雅克比矩阵;

<mrow> <mi>J</mi> <mo>=</mo> <mo>[</mo> <mfrac> <mi>df</mi> <mrow> <mi>dp</mi> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </mrow> </mfrac> <mo>,</mo> <mfrac> <mi>df</mi> <mrow> <mi>dp</mi> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </mrow> </mfrac> <mo>,</mo> <mo>.</mo> <mo>.</mo> <mo>.</mo> <mo>,</mo> <mfrac> <mi>df</mi> <mrow> <mi>dp</mi> <mrow> <mo>(</mo> <msub> <mi>n</mi> <mn>1</mn> </msub> <mo>)</mo> </mrow> </mrow> </mfrac> <mo>]</mo> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>13</mn> <mo>)</mo> </mrow> </mrow>

光束平差法以其待调整参数的雅克比矩阵组成迭代因子,对模型参数进行调整, 得到使成像误差S(θ)达到最小的模型参数最优解。

(2)待测区域特征点的提取

本发明针对待测区域特征点的提取采用的是Harris角点检测算法;Harris 提出了用高斯函数代替方形区域计算梯度的算法,若目标点的像素坐标为x,y, 其在x方向和y方向的位移分别为u和v,则该目标点的灰度变化由如下公式表 示:

<mrow> <mi>E</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>)</mo> </mrow> <mo>=</mo> <munder> <mi>&Sigma;</mi> <mrow> <mi>u</mi> <mo>,</mo> <mi>v</mi> </mrow> </munder> <mi>&omega;</mi> <mrow> <mo>(</mo> <mi>u</mi> <mo>,</mo> <mi>v</mi> <mo>)</mo> </mrow> <msup> <mrow> <mo>[</mo> <mi>I</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>+</mo> <mi>u</mi> <mo>,</mo> <mi>y</mi> <mo>+</mo> <mi>v</mi> <mo>)</mo> </mrow> <mo>-</mo> <mi>I</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>)</mo> </mrow> <mo>]</mo> </mrow> <mn>2</mn> </msup> <mo></mo> </mrow>

<mrow> <mo>&ap;</mo> <munder> <mi>&Sigma;</mi> <mrow> <mi>u</mi> <mo>,</mo> <mi>v</mi> </mrow> </munder> <mi>&omega;</mi> <mrow> <mo>(</mo> <mi>u</mi> <mo>,</mo> <mi>v</mi> <mo>)</mo> </mrow> <mrow> <mo>(</mo> <mi>u</mi> <mo>,</mo> <mi>v</mi> <mo>)</mo> </mrow> <mfenced open='[' close=']'> <mtable> <mtr> <mtd> <msubsup> <mi>I</mi> <mi>x</mi> <mn>2</mn> </msubsup> </mtd> <mtd> <msub> <mi>I</mi> <mi>x</mi> </msub> <msub> <mi>I</mi> <mi>y</mi> </msub> </mtd> </mtr> <mtr> <mtd> <msub> <mi>I</mi> <mi>x</mi> </msub> <msub> <mi>I</mi> <mi>y</mi> </msub> </mtd> <mtd> <msubsup> <mi>I</mi> <mi>y</mi> <mn>2</mn> </msubsup> </mtd> </mtr> </mtable> </mfenced> <mfenced open='(' close=')'> <mtable> <mtr> <mtd> <mi>u</mi> </mtd> </mtr> <mtr> <mtd> <mi>v</mi> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>14</mn> <mo>)</mo> </mrow> </mrow>

<mrow> <mi>M</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>)</mo> </mrow> <mo>=</mo> <munder> <mi>&Sigma;</mi> <mrow> <mi>u</mi> <mo>,</mo> <mi>v</mi> </mrow> </munder> <mi>&omega;</mi> <mrow> <mo>(</mo> <mi>u</mi> <mo>,</mo> <mi>v</mi> <mo>)</mo> </mrow> <mfenced open='[' close=']'> <mtable> <mtr> <mtd> <msubsup> <mi>I</mi> <mi>x</mi> <mn>2</mn> </msubsup> </mtd> <mtd> <msub> <mi>I</mi> <mi>x</mi> </msub> <msub> <mi>I</mi> <mi>y</mi> </msub> </mtd> </mtr> <mtr> <mtd> <msub> <mi>I</mi> <mi>x</mi> </msub> <msub> <mi>I</mi> <mi>y</mi> </msub> </mtd> <mtd> <msubsup> <mi>I</mi> <mi>y</mi> <mn>2</mn> </msubsup> </mtd> </mtr> </mtable> </mfenced> </mrow> 则:

<mrow> <mi>E</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>)</mo> </mrow> <mo>&ap;</mo> <mrow> <mo>(</mo> <mi>u</mi> <mo>,</mo> <mi>v</mi> <mo>)</mo> </mrow> <mi>M</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>)</mo> </mrow> <mfenced open='(' close=')'> <mtable> <mtr> <mtd> <mi>u</mi> </mtd> </mtr> <mtr> <mtd> <mi>v</mi> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>15</mn> <mo>)</mo> </mrow> </mrow>

计算矩阵M的特征值λ1,λ2,如果两个特征值都比较大,说明目标点的图像灰度 自相关函数的两个正交方向上的值均较大,则该点即为特征点;通过Harris算 法对特征点进行初提取后;采用亚像素角点提取算法提取出精度更高的角点坐 标;对于理想角点,它附近的像素点的灰度梯度方向均垂直于该点与理想角点 的连线;这一特点可以用公式表达为:

<mrow> <mo>&dtri;</mo> <mover> <mi>H</mi> <mo>&RightArrow;</mo> </mover> <mrow> <mo>(</mo> <mover> <mi>&alpha;</mi> <mo>&RightArrow;</mo> </mover> <mo>-</mo> <mover> <mi>&beta;</mi> <mo>&RightArrow;</mo> </mover> <mo>)</mo> </mrow> <mo>=</mo> <mn>0</mn> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>16</mn> <mo>)</mo> </mrow> </mrow>

其中是理想角点的灰度梯度方向,向量为图像原点指向理想角点的坐标,
向量为图像原点指向理想角点附近任一边缘点的坐标;实际上,由于受到图像
噪声的影响,公式(10)通常都不为零,可以将其视为误差e,即

<mrow> <mi>e</mi> <mo>=</mo> <mo>&dtri;</mo> <mover> <mi>H</mi> <mo>&RightArrow;</mo> </mover> <mrow> <mo>(</mo> <mover> <mi>&alpha;</mi> <mo>&RightArrow;</mo> </mover> <mo>-</mo> <mover> <mi>&beta;</mi> <mo>&RightArrow;</mo> </mover> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>17</mn> <mo>)</mo> </mrow> </mrow>

在以角点为中心的邻域内,对所有点按上式计算,误差和为E,则

<mrow> <mi>E</mi> <mo>=</mo> <munder> <mi>&Sigma;</mi> <mi>i</mi> </munder> <mi>e</mi> <mo>=</mo> <munder> <mi>&Sigma;</mi> <mi>i</mi> </munder> <mo>&dtri;</mo> <mover> <mi>H</mi> <mo>&RightArrow;</mo> </mover> <mrow> <mo>(</mo> <mover> <mi>&alpha;</mi> <mo>&RightArrow;</mo> </mover> <mo>-</mo> <mover> <mi>&beta;</mi> <mo>&RightArrow;</mo> </mover> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>18</mn> <mo>)</mo> </mrow> </mrow>

其中,i是邻域内第i个点;这样求使误差和E最小的点即是亚像素级别的角点 坐标。

(3)特征点的立体匹配

采用归一化交叉算法(NCC)实现特征点的初匹配,以CCD摄像机2拍得 的图片作为基准图片,并以特征点pl为中心,构造一个N*N(代表图像像素)的 局部图像块作为模版图T,并令模版图T在待匹配图像S的包含有特征点的一定 范围内进行遍历,搜索模版所覆盖的子图记做Si,j,i,j为子图中心点在匹配图 像S中的像素坐标,通过使用归一化交叉算法(NCC)来计算模版图T和搜索 子图Si,j之间的归一化互相关系数,算法如下:

<mrow> <mi>R</mi> <mrow> <mo>(</mo> <mi>i</mi> <mo>,</mo> <mi>j</mi> <mo>)</mo> </mrow> <mo>=</mo> <mfrac> <mrow> <munder> <mi>&Sigma;</mi> <mrow> <mi>m</mi> <mo>,</mo> <mi>n</mi> </mrow> </munder> <mo>|</mo> <msup> <mi>S</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> </msup> <mrow> <mo>(</mo> <mi>m</mi> <mo>,</mo> <mi>n</mi> <mo>)</mo> </mrow> <mo>-</mo> <mi>E</mi> <mrow> <mo>(</mo> <msup> <mi>S</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> </msup> <mo>)</mo> </mrow> <mo>|</mo> <mo>|</mo> <mi>T</mi> <mrow> <mo>(</mo> <mi>m</mi> <mo>,</mo> <mi>n</mi> <mo>)</mo> </mrow> <mo>-</mo> <mi>E</mi> <mrow> <mo>(</mo> <mi>T</mi> <mo>)</mo> </mrow> <mo>|</mo> </mrow> <msqrt> <munder> <mi>&Sigma;</mi> <mrow> <mi>m</mi> <mo>,</mo> <mi>n</mi> </mrow> </munder> <msup> <mrow> <mo>(</mo> <msup> <mi>S</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> </msup> <mrow> <mo>(</mo> <mi>m</mi> <mo>,</mo> <mi>n</mi> <mo>)</mo> </mrow> <mo>-</mo> <mi>E</mi> <mrow> <mo>(</mo> <msup> <mi>S</mi> <mrow> <mi>i</mi> <mo>,</mo> <mi>j</mi> </mrow> </msup> <mo>)</mo> </mrow> <mo>)</mo> </mrow> <mn>2</mn> </msup> <munder> <mi>&Sigma;</mi> <mrow> <mi>m</mi> <mo>,</mo> <mi>n</mi> </mrow> </munder> <msup> <mrow> <mo>(</mo> <mi>T</mi> <mrow> <mo>(</mo> <mi>m</mi> <mo>,</mo> <mi>n</mi> <mo>)</mo> </mrow> <mo>-</mo> <mi>E</mi> <mrow> <mo>(</mo> <mi>T</mi> <mo>)</mo> </mrow> <mo>)</mo> </mrow> <mn>2</mn> </msup> </msqrt> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>19</mn> <mo>)</mo> </mrow> </mrow>

其中m,n代表各个像素的坐标,E(Si,j)和E(T)分别是搜索子图Si,j与模版图T 的灰度平均值;互相关系数R(i,j)的值越大,则搜索子图与模板图的匹配程度 越高;设定互相关系数R(i,j)的阈值后,选出特征点pl的候选匹配点;使用NCC 算法进行初匹配后还需要利用外极限约束和距离约束来对匹配的特征点对进行 校正;通过使用Longguet‑Higgins提出的归一化8点算法在摄像机标定时计算左 右两个CCD摄像机的基本矩阵F,根据外极线理论,对基准图像特征点pl,在 待匹配图像S中对应的外极线lr可以表示为:

lr=F·pl                       (20)

初匹配的候选特征点若在外极线lr附近则可进一步确定该候选特征点为匹配点; 最后统计每一候选特征点到其他特征点的距离和基准图像中特征点到其他特征 点点距离相等的个数num;num最大值时对应的候选特征点则为基准图像中特 征点的正确匹配点。

(4)三维坐标的求取

首先计算出双目图像中特征点经亚像素算法提取出的图像坐标值,并大致 估算出特征点在世界坐标系中的估计值;在优化程序中将左右两个CCD摄像机 参数对应雅克比矩阵设为零,同时在优化程序计算迭代因子部分中嵌入关于空 间坐标x、y、z的雅克比矩阵方程;将已经优化好的左右两个CCD摄像机参数 以及特征点点图像坐标值和世界坐标估计值输入光束平差法优化程序,就可以 对空间点的三维坐标进行优化,得到特征点的三维空间测量值,从而完成三维 坐标点的求取并最终实现特征点的定位。

本发明的有益效果是解决了双目显微视觉标定时透镜多、视场狭窄、景深 小、具有一定放大倍率并包含有较多的畸变因素的技术困难,实现了双目显微 立体视觉的高精度标定,从而完成了精密零件的精确定位。

图1所示为基于双目显微立体视觉的精密零件定位的装置模型图。其中, 1—为双目显微立体视觉系统提供照明的LED光源,2—左CCD摄像机,2’—右 CCD摄像机,3—立体显微镜,4—待测量的工件;5—可倾数控转台;6—待测 量区域;7—精密数控位移平台;8—床身;9—工业控制计算机。

图2所示为基于双目显微立体视觉的精密零件定位测量方法的流程图。

图3所示为基于光束平差算法的双目显微立体视觉的标定和三维坐标点的 求取的流程图。

以下结合技术方案和附图详细叙述本发明的具体实施方式。附图1为基于 双目显微立体视觉的精密零件定位的装置模型图。本装置通过左右两个CCD摄 像机2、2’采集待测工件中待测量区域的位置信息,通过已建立好的图像坐标系 和世界坐标系之间的关系,到待测量点的三维空间坐标即实现精密定位。其 装置的安装方式如下:机床床身8置于地面;精密数控位移平台7通过导轨与 机床床身安装在一起,其可以满足标定实验中对棋盘格标定板的高精度平移和 角度旋转需求;可倾数控转台5通过螺栓固定在铸铁平台的T形槽上;待测工 件4通过专用夹具安装在可倾数控转台上,其待测区域6位于双目显微立体视 觉系统的视场内;立体显微镜3通过专用的显微镜夹具安装在机床床身上;LED 光源1固定在夹具上,夹具通过螺栓固连在床身上;左右两个CCD摄像机2和 2’通过螺纹连接安装在立体显微镜上;左右两个CCD摄像机采集到的图像分别 通过1394图像采集卡传递到工业控制计算机9内进行图像数据处理。

附图1为本发明的一个实施例,采用左右两个CCD摄像机2、2’拍摄待测 物体的图像位置信息,左右两个CCD摄像机采用的是奥林巴斯DP26摄像机, 图像分辨率:2448*1920,采集速度:7fps,芯片尺寸:2/3英寸。立体显微镜3 采用的是奥林巴斯生产的SZX‑16研究级立体显微镜。变焦倍率:0.7‑11.5,工 作距离:60mm,采集图像最大视场:12.5mm×16.6mm,采集图像最小视场: 0.76mm×1.02mm。为实现双目显微立体视觉系统的标定,标定所需的标定板采 用的是深圳科创时代生产的CG‑050‑T‑0.5型玻璃基底棋盘格标定板。该型标定 板方格大小为0.5mm×0.5mm,棋盘格图案整体宽度为51mm×51mm,制造精度 为1μm,能够满足立体显微镜标定的需要。

附图2表示了基于双目显微立体视觉的精密零件定位测量方法的流程图, 精密定位的主要流程包括左右两个CCD摄像机2、2’的标定,双目图像特征点 的提取,双目图像特征点的匹配与校正,双目图像特征点的三维坐标求取。其 中,左右两个CCD摄像机的标定是利用高精密棋盘格标定板来实现对左右两 CCD摄像机2、2’的主参数和畸变参数求解的,主要包括采用张氏标定法对尺度 因子的标定、采用变倍率法实现主点坐标的标定和采用光束平差法对主参数、 畸变参数的优化;双目图像的特征点利用Harris角点检测算法和亚像素提取算 法进行提取;双目图像的标记点匹配是通过对双目图像的标记点进行初匹配及 特征点对的校正来实现的,最终通过光束平差法实现三维坐标点的求取。

附图3所示为基于光束平差算法的双目显微立体视觉的标定和三维坐标点 求取的流程图,通过拍摄得到的不同倍率下的图像使用变倍率法获得主点坐标, 通过30副各种角度的图像使用张氏标定法获得尺度因子,并通过已获得的内参 数求取左右两个CCD摄像机的外参数,将内参数和外参数作为光束平差法优化 的初始值,并使用激光干涉仪获取高精度的立体点阵作为优化的约束条件,以 待调整参数的雅克比矩阵组成迭代因子,通过光束平差法获得优化后的参数, 将优化后的参数带入到成像模型中得到最终成像模型,最后把成像模型、待测 点的图像坐标及三维坐标估计值输入到光束平差法的程序中求取待测点的空间 三维坐标。

测量方法的具体步骤如下:

(1)左右两个CCD摄像机的标定

本发明采用摄像机相对固定的方式,以精密加工的棋盘格标定板对摄像机 进行标定的方法。张氏摄像机标定法要求不同标定板图片有较大的角度关系, 通常最佳角度为45°;而我们限于立体显微镜的景深因素,当标定板与水平面角 度超过20°,视场内会出现许多模糊的角点,为了能够稳定的得到摄像机的尺度 因子,并对张氏标定的实验过程进行总结,发现标定板图片数量越多,相互平 行的标定板图片数量越少,就越能稳定的求得摄像机的尺度因子,因此对左右 两CCD摄像机2、2’的尺度因子的标定的具体步骤如下:将标定板旋转分为沿 水平轴旋转和竖直轴旋转两种方式,竖直轴旋转每次旋转40°,这样标定板在就 存在九个水平角位置。在每个水平角位置进行水平轴旋转,将标定板水平轴旋 转位置分为+10°,0°,10°三个位置,在每个位置上拍摄一张标定板图像。整个 对尺度因子的标定过程共拍摄54张标定板图片。针对主点坐标,考虑到主点实 际上是光轴与CCD成像面的交点,而在显微镜放大倍率发生改变时,光轴的位 置实际上是不变的也就是摄像机的主点坐标始终保持不变,在任何放大倍率下, 点(x1,y1,z1)的图像坐标都在同一直线上,且该直线一定经过摄像机主点(u0,v0); 这样,n个在图像平面上的投影不重合的空间点,通过在不同放大倍率下进行投 影,就能够形成n条相交于同一点的直线;这个交点的坐标就是摄像机主点坐 标。采用变倍率法实现主点坐标的标定具体操作步骤:确定12个在图像投影平 面上不重合的点,令这12个点分别在不同放大倍率下0.7、1、2、2.4、3、5、 5.8、8进行投影,利用最小二乘法对获得的12条相交于一点的直线,求取交点 坐标值。外参数矩阵在标定完内参数及单应性矩阵利用公式7进行求解。采用 光束平差法对内外标定参数的优化的过程:在优化程序中输入成像方程、空间 点图像坐标观察值、待优化的参数初值,程序的输出即为能使S(θ)输出最小的优 化后的标定参数值。

为了验证标定结果的精度,利用标定后的模型参数,对标定板提供的三维 空间点进行三维重建,得到三维空间点的测量值。通过将测量值与三维空间点 坐标的实际值进行比较分析,来评价参数的标定精度。标定结果显示,在1倍 率下,x方向、y方向、z方向的平均标定精度分别为:2.3um、1.4um、7.9um; 在2倍率下,x方向、y方向、z方向的平均标定精度分别为:1.5um、1.2um、 5.4um;在5倍率下,x方向、y方向、z方向的平均标定精度分别为:0.7um、 0.5um、3.2um。标定结果表明标定后的模型参数具有较高的横向重建精度和较 好的纵向重建精度。

(2)特征点的提取

对左右两个CCD摄像机2和2’采集到的图片通过使用采用Harris角点检测 算法进行特征点X1,X2……Xn初提取,Harris角点检测算法是通过计算特征值 λ1,λ2,如果两个特征值都比较大,说明目标点的图像灰度自相关函数的两个正 交方向上的值均较大,则该点即为特征点。以提取后的特征点X1,X2……Xn为 中心,分别计算每一个特征点的误差和E,求使误差和E最小的点即是亚像素 级别的角点位置。通过这种方式,就能够在Harris角点提取的基础上得到亚像 素角点坐标,提高特征点的提取精度。

(3)特征点的匹配

采用归一化交叉算法(NCC)实现特征点的初匹配,根据算法公式(19),以 CCD摄像机2拍得的图片作为基准图片,并以特征点pl为中心,构造一个N*N 的局部图像块作为模版图T,并令模版图T在待匹配图像S的包含有特征点的一 定范围内进行遍历,搜索模版所覆盖的子图记做Si,j,i,j为子图中心点在匹配 图像S中的像素坐标,通过使用归一化交叉算法(NCC)来计算模版图T和搜 索子图Si,j之间的归一化互相关系数。互相关系数R(i,j)的值越大,则搜索子图 与模板图的匹配程度越高。设定互相关系数R(i,j)的阈值后,选出特征点pl的候 选匹配点。使用NCC算法进行初匹配后还需要利用外极线约束和距离约束来对 匹配的特征点对进行校正;通过使用Longguet‑Higgins提出的归一化8点算法在 摄像机标定时计算左右两个CCD摄像机的基本矩阵F,根据外极线理论,对基 准图像特征点pl,使用公式(20)计算出基准图像中的特征点pl在待匹配图像中 的外极线位置,初匹配的候选特征点若在外极线附近则可进一步确定该候选特 征点为匹配点。最后是通过距离约束进一步校正匹配点对。采用归一化交叉算 法(NCC)进行特征点的初匹配,外极线约束和距离约束来校正匹配的特征点 对,左右两个图像匹配的特征点对分别为:Xl,Xr',X2l,X2r'……X3l,X3r'。

(4)三维坐标的求取

光束平差法除了能够对摄像机参数进行优化外,还能同时进行三维坐标点 的优化;此处三维坐标点的优化即是三维坐标的求取;求取过程如下:将优化 程序中摄像机参数对应雅克比矩阵设为零,同时在优化程序计算迭代因子部分 中嵌入关于空间坐标x、y、z的雅克比矩阵方程;将已经优化好的左右两个CCD 摄像机参数以及空间点的图像坐标和世界坐标估计值输入光束平差法优化程 序,就可以对空间点的三维坐标进行优化,得到三维空间点的测量值,从而完 成三维坐标点的求取。

本发明较好的解决了由于待测目标区域小、定位精度要求高、非接触等问 题所产生的测量难题。采用基于双目显微立体视觉的非接触式测量方法,很好 的完成了精密零件的精确定位。

本文发布于:2024-09-23 22:38:31,感谢您对本站的认可!

本文链接:https://www.17tex.com/tex/1/73655.html

版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系,我们将在24小时内删除。

留言与评论(共有 0 条评论)
   
验证码:
Copyright ©2019-2024 Comsenz Inc.Powered by © 易纺专利技术学习网 豫ICP备2022007602号 豫公网安备41160202000603 站长QQ:729038198 关于我们 投诉建议