CN109712195A - 利用球像的公共自极三角形进行单应性估计的方法 - Google Patents

Abstract

本发明涉及一种利用球像的公共自极三角形进行单应性估计的方法，以空间中的一个球作为靶标，在抛物折反射系统下,首先,从不同位置拍摄一个空间球的两幅图像，并且拟合出镜面轮廓投影方程、靶标投影方程及对拓球像方程；其次,分别计算两幅图像中球像与对拓球像构成的矩阵的特征值与特征向量；然后通过匹配特征值，求出三对线对应，这三对线对应即公共自极三角形的三对对应边；接着，连接公共自极三角形与二次曲线的交点，求得第四对线对应，这对线对应位于公共自极三角形内；最后，对这四对线对应采用SVD分解求出单应性矩阵。

Description

利用球像的公共自极三角形进行单应性估计的方法

技术领域

本发明属于计算机视觉领域，一个空间球在像平面上的投影是球像与对拓球像，利用球像与对拓球像的公共自极三角形进行单应性估计。

背景技术

计算机视觉是一门关于如何运用摄像机和电脑来获取被拍摄物体的数据与信息的科学，并且是人类观察世界和认知世界的重要手段，计算机视觉以视觉处理为中心，是人工智能范畴的一个新领域。它也是以图像处理、模式识别、计算机技术和生理学心理学为基础的信息处理科学中的一个重要分支。计算机视觉技术集数字图像处理、数字信号处理、光学、物理学、几何学、应用数学、模式识别及人工智能等知识于一体，主要应用于机器人导航、医学辅助诊断、监控、跟踪、制作三维环境模型等。由于这些应用对摄像机的可视范围有较高的要求，目前符合这一要求的摄像机是折反射摄像机。而折反射摄像机又根据是否有固定的单视点，分为中心和非中心折反射摄像机两种类型。文献“A theory of single-viewpoint catadioptric image formation”(Baker S.,Nayar K.,InternationalJournal of Computer Vision,35(2):175-196,1999)根据镜面参数ξ的值，将中心折反射镜面分为4种类型：ξ＝0，平面镜；0＜ξ＜1，椭圆或双曲镜面；ξ＝1，抛物镜面.从过去到目前为止中心折反射摄像机的标定方法总共可以分为两类：自标定和用标定物标定。

在中心折反射摄像机系统下，由于球像的轮廓可以被完全提取，这有利于提高标定的精确性，从而本文利用单个球在抛物折反射摄像机单位视球模型下的投影性质，根据球像与对拓球像的公共自极三角形进行单应性估计。

计算机视觉是二十世纪六十年代中期迅速发展起来的一门新学科.主要任务是通过对采集的图片或视频进行处理以获得相应场景的三维信息。随着科技的发展和进步，摄像机标定在计算机视觉领域扮演着重要角色,是计算机视觉领域工作的基础，其标定结果与后续工作息息相关，并且标定的精确性直接影响着后续的工作。在一些应用中，由于球的外轮廓具有从任意位置观察都是可见和自身无遮挡的优势，所以利用球作为标定物使得标定结果比其他几何体更加精确。文献“Catadioptric projective geometry”(Geyer C.,Daniilidis K.，International Journal of Computer Vision,45(3):223-243,2001)和文献“A Unifying Theory for Central Panoramic Systems and PracticalApplications”(Geyer C.,Daniilidis K.,1843:445-461,2000)中研究了折反射摄像机的投影几何理论，证明了中心折反射摄像机的成像过程等价于经过单位视球的两步成像过程,并且证明在中心折反射摄像机下,直线在图像平面上的像为一条二次曲线。由于单位视球模型的提出为中心折反射摄像机的研究及应用提供了良好的数学依据,从而文献“Catadioptric camera calibration using geometric invariants”(Ying X.,Hu,Z.,IEEE Transactions on Pattern Analysis and Machine Intelligence,26(10):1260-1271,2004)首次提出利用球标定中心折反射摄像机,并且证明了空间直线和空间球在图像平面上的投影都是二次曲线，也证明了一条直线的投影提供三个不变量，一个球的投影提供两个不变量，由这些不变量，就可以获得中心折反射摄像机内参数的约束方程，并且得到一个重要结论：球的投影比线的投影鲁棒性好，精确度高.文献“Camera calibration withspheres:linear approaches”(Zhang H.,Zhang G.,Wong K.Y.,InternationalConference on Image Processing,IEEE,II-1150-53,2005)通过利用一幅图像上至少三个球的表面轮廓，提出了两种标定算法，即标量和正交方法。实验表明这两种方法的结果与以前提出的方法得到的结果相似甚至优于以前的结果。文献“Calibration of aparacatadioptric camera by projection imaging of a single sphere”(Li Y.,ZhaoY.,Applied Optics,56(8):2230,2017)根据一个球在单位视球上的投影是两个平行圆，从而利用三种不同的方法标定中心折反射摄像机的内参数，分别是:基于正交消失点的摄像机标定；基于平行圆的摄像机标定；基于消失线的摄像机标定。近年来，球已经广泛应用于摄像机标定。根据双接触理论“The seven circles theorem and other new theorems”(Evelyn C.,Money–Coutts G.，London:Stacey International Press,1974),文献“Geometric interpretations of the relation between the image of the absoluteconic and sphere images”(Ying X.,Zha H.，IEEE Transactions on Pattern Analysisand Machine Intelligence,28(12):2031–2036,2006)对球像与IAC的关系提出了两种几何解释。对于一个球像的几何解释是球像与IAC在两个双接触像点处相切，这两个切点是球的轮廓所在支撑平面的两个圆环点的像。对于三个球像的几何解释是三个球像和IAC满足双接触定理,IAC可直接由双接触定理确定。在文献“Geometric interpretations of therelation between the image of the absolute conic and sphere images”(Ying X.,Zha H.，IEEE Transactions on Pattern Analysis and Machine Intelligence,28(12):2031–2036,2006)的基础上,文献“Identical projective geometric properties ofcentral catadioptric line images and sphere images with applications tocalibration”(Ying X.,Zha H.,International Journal of Computer Vision,78(1):89-105,2008)发现在中心折反射摄像机系统下,绝对二次曲线的像与线像或球像没有双接触，但是存在一条虚构的二次曲线，定义为改进的绝对二次曲线的像(MIAC)，它在中心折反射摄像机系统下与线像或球像均双接触，从而得到一些线像和球像都具有的新的射影几何性质，由这些关于MIAC的几何性质，从而提出了两种用线像和球像线性标定中心折反射摄像机的方法。文献“A novel linear approach to camera calibration from sphereimages”(Ying X.,Zha H.,1(1):535-538,2006)基于球像矩阵分量的恒等约束提出了一种新的线性标定方法。文献“A calibration method for paracatadioptric camera fromsphere images”(Duan H.,Wu Y.,Pattern Recognition Letters,33(6):677-684,2012)和文献“Paracatadioptric camera calibration using sphere images”(Duan H.,WuY.,模式识别国家重点实验室,641-644.2011)中提出了一种基于球像的抛物折反射摄像机的标定方法，从而完善了基于球的中心折反射摄像机的标定方法。文献“Cameracalibration from images of spheres”(Zhang H.,Wong K.Y.K.,Zhang G.,IEEETransactions on Pattern Analysis&Machine Intelligence,29(3):499-503,2007)利用球像的对偶和绝对二次曲线的像的对偶之间的关系，表明两个球在图像平面上的投影的公共极点和极线也是绝对二次曲线的像的极点和极线。文献“Intrinsic parameterdetermination of a paracatadioptric camera by the intersection of two sphereprojections”(Zhao Y.,Wang Y.,Optical Society of America,32(11),2015)利用两个球在单位视球上的两组投影圆相交于四个点，且这四个点构成了一个矩形，从而确定了一组正交方向，然后根据正交消失点与绝对二次曲线的像的关系，至少用五幅图像求解摄像机内参数矩阵。近几年,公共自极三角形已经被用于摄像机标定。文献“Cameracalibration based on the common self-polar triangle of sphere images”(HuangH.,Zhang H.,Cheung Y.M.,In Computer Vision–ACCV 2014,pages 19–29.Springer,2014)提出同一平面上任意两个不相交的球像有且仅有一个公共自极三角形，从而利用三个两两不相交球像之间的公共自极三角形，通过用位于二次曲线外的顶点恢复图像平面的消失线，进而求解圆环点的像来标定摄像机。该方法从现有的射影几何理论出发恢复摄像机内参数，不需要计算球心的投影，从而比以前的方法更简单。文献“The common self-polar triangle of concentric circles and its application to cameracalibration”(Huang H.,Zhang H.,Cheung Y.M.,Computer Vision and PatternRecognition,IEEE,4065-4072,2015)研究并提出两个同心圆有无数多个公共自极三角形，通过研究这些公共自极三角形,发现它们具有两个重要性质：第一，所有公共自极三角形都有一个公共顶点和一条公共边，公共顶点是同心圆的圆心，公共边是同心圆所在平面的无穷远直线；第二，所有公共自极三角形都是直角三角形。基于这两个性质，圆心的像和消失线可以同时被恢复，从而由绝对二次曲线的像与圆环点的像之间的线性约束，求得摄像机的内参数矩阵。单应性矩阵是计算机视觉算法的一个基本工具，它已经成功地应用于许多不同的领域。比如场景重构，摄像机标定，视觉测量，位姿估计和对象识别等。文献“Homography estimation”(Dubrofsky E.,Master’s thesis,University of BritishColumbia(Vancouver),2009)中描述了单应性矩阵的定义，即单应性矩阵是从二维射影空间到二维射影空间的一个可逆映射，并且给出了求解单应性矩阵的方法。通常单应性矩阵估计用到的标定类型有：点，线和二次曲线。由于用二次曲线估计单应性矩阵的方法过于复杂。文献“Conics-based homography estimation frominvariant points and pole-polar relationships”(Conomis C.,International Symposium on 3D DataProcessing,Visualization,and Transmission,IEEE,908-915,2007)提出了一种基于两条无模型的共面二次曲线的单应性矩阵估计的方法。即给定两幅图像中一对非退化的共面二次曲线,它们之间与单应性矩阵有关。在每幅图像中,当由两条共面二次曲线的矩阵组成的矩阵的秩为3时,有三个特征向量,如果把特征向量看作公共自极三角形的极点，则这两幅图像提供三对点对应。寻找另一对点对应的主要想法是：求极点所对应的极线，极线与二次曲线的交点,即为第四对点对应，为了确保是实交点，选取的交点至少在一条二次曲线外。文献“Homography estimation from the common self-polar triangle of separateellipses”(Huang H.,Zhang H.,Cheung Y.M.,Computer Vision and PatternRecognition,IEEE,1737-1744,2016)提出用两个分离椭圆的公共自极三角形进行单应性估计。在每幅图像中,如果把由两条不相交的二次曲线的矩阵组成的矩阵的特征向量看作公共自极三角形的极线，则这两幅图像提供三对线对应，通过研究公共自极三角形的位置，第四对线对应由公共自极三角形与两条二次曲线的交点得到,这对线对应分别位于两条二次曲线的内部。一对线对应提供关于单应性矩阵的两个元素的两个独立方程，四对线对应提供八个独立方程，就可求得单应性矩阵。

单应性分两种情况，一种是空间平面与图像平面之间的单应性；如果空间平面在世界坐标系的坐标已知，则该单应性提供关于摄像机内参数的两个约束；如果空间平面在世界坐标系的坐标未知，则该单应性不能提供关于摄像机内参数的两个约束；另一种是两幅图像平面之间的单应性；两幅图像平面之间的单应性不能对摄像机内参数构成任何约束；本文求的是两幅图像平面之间的单应性，两幅图像平面之间的单应性可用于图像配准、图像矫正和图像拼接，也可以用于寻求图像间的稠密匹配点。

发明内容

本发明利用空间中的一个球作为靶标。在抛物折反射系统下,首先,从不同位置拍摄一个空间球的两幅图像，并且拟合出镜面轮廓投影方程、靶标投影方程及对拓球像方程；其次,分别计算两幅图像中球像与对拓球像构成的矩阵的特征值与特征向量；然后通过匹配特征值，求出三对线对应，这三对线对应即公共自极三角形的三对对应边；接着，连接公共自极三角形与二次曲线的交点，求得第四对线对应，这对线对应位于公共自极三角形内；最后，对这四对线对应采用SVD分解求出单应性矩阵。

本发明采用如下技术方案：

本发明以一个空间球作为靶标，用抛物折反射摄像机从不同位置拍摄这个空间球的两幅图像，通过Canny边缘检测算子分别对这两幅图像提取球像和镜面轮廓的边缘点，接着利用最小二乘法拟合获得每幅图像中球像和抛物镜面轮廓的投影方程，再根据抛物镜面轮廓的投影方程计算出摄像机内参数的初始值，最后，用优化方法求出优化后球像及对拓球像的方程。根据射影几何的知识，由于球像与对拓球像互不相交，从而它们有且仅有一个公共自极三角形，因此两幅图像提供两个公共自极三角形，通过计算球像与对拓球像构成的矩阵的特征值与特征向量，从而得到三对线对应,第四对线对应由公共自极三角形与二次曲线的交点的连线得到,并且这对线对应位于公共自极三角形内。具体的步骤包括：拟合出镜面轮廓投影方程、球像方程，根据像点与其对拓像点的关系,获得对拓像点，从而拟合出球像的对拓球像方程,再根据奇异值分解求得球像与对拓球像构成的矩阵的特征值与特征向量,以特征向量作为公共自极三角形的对应边，求出三对线对应，第四对线对应是连接公共自极三角形与二次曲线的交点得到的，这对线对应位于公共自极三角形内。然后由四对线对应的单应性关系,求得单应性矩阵，最后再通过公共自极三角形与二次曲线的交点的另一对连线验证本文求得的单应性矩阵是有效可行的。

1.拟合镜面轮廓投影方程及靶标投影方程

在中心折反射摄像机下，由于镜面轮廓可以被完全提取，利用MATLAB中的Canny算子提取镜面轮廓投影边缘点,并对每幅靶标图像进行边缘检测，提取靶标图像边缘点的像素坐标，然后对获取的边缘点利用最小二乘法拟合得到镜面轮廓投影方程、球像方程。

2.拟合球像的对拓球像

一个空间球Q在抛物折反射摄像机下的投影过程等价于经过单位视球的两步投影过程：第一步，球Q以单位视球球心O为透视中心投影到单位视球球面上形成两个平行圆S_i+和S_i-(i＝1,2表示拍摄两幅图像)，并且S_i+和S_i-以单位视球球心O对称，其中S_i+是可见的，S_i-是不可见的。第二步，以单位视球上的点O_c为虚拟摄像机的光心将单位视球上的两个平行圆S_i+和S_i-分别投影到图像平面上形成两条二次曲线C_i+和C_i-，S_i+的像C_i+称为球像，S_i-的像C_i-称为对拓球像，它们是对拓关系，并且C_i+是可见的，C_i-是不可见的。利用MATLAB中的Canny算子提取第一幅图像中的镜面轮廓投影的边缘点和两幅图像平面上的球像与对拓球像的边缘点，通过最小二乘法拟合得到相应的二次曲线方程。这里用C₀表示第一幅图像镜面轮廓的系数矩阵。本文为了简化表述，用相同字母表示曲线和它的系数矩阵。通过C₀可得到折反射摄像机内参数矩阵K的一个初始值K₀，从而得到绝对二次曲线的像ω的初始值这里:ω＝K^-TK^-1，其中r是纵横比的初始值，f是有效焦距的初始值，s是倾斜因子的初始值，[u₀ v₀ 1]^T是折反射摄像机主点的初始齐次坐标矩阵形式，记为p＝[u₀ v₀ 1]^T。取C_i+(i＝1,2)上的一组点m_ij+(i＝1,2,j＝1,2,…,n,n≥5)，则与它相对应的一组对拓像点为m_ij－(i＝1,2,j＝1,2,…,n,n≥5)。因为本文只讨论抛物折反射摄像机的情况，从而关于对拓像点的关系式如下：这里m_ij+，m_ij－用齐次坐标表示。根据对拓像点的定义，点m_ij－在球像C_i+(i＝1,2)的对拓球像C_i-(i＝1,2)上，因此用最小二乘法拟合得到对拓球像C_i-的方程。在抛物折反射摄像机下，球像及其对拓球像满足下列表达式：其中这里C'_i+是小圆S_i+的系数矩阵，C'_i-是对拓小圆S_-的系数矩阵，C_i+是球像的系数矩阵，C_i-是对拓球像的系数矩阵,[n_x n_y n_z]^T是小圆所在基础平面的单位法向量。

3.公共自极三角形的恢复

设两条共面二次曲线C_i+,C_i-(i＝1，2)的公共极点极线是点m_i(i＝1,2)和线l_i(i＝1,2)，且满足如下关系式这里，μ_i是非零比例因子，通过化简得

其中E是单位矩阵。由上式,得到C_i+和C_i－的公共极线l_i是的特征向量。

设Δefg是C₁₊,C_1-的公共自极三角形，e,f,g是顶点，Δe'f'g'是C₂₊,C_2-的公共自极三角形，e',f',g'是顶点。第一幅图像和第二幅图像与一个单应性H有关，从而C₁₊与C₂₊，C_1-与C_2-满足如下关系式：令从而有Q₂＝(H^T)^{- 1}Q₁H^T。故矩阵Q₁与Q₂是相似矩阵，根据相似矩阵具有相同的特征值，而相同特征值对应的特征向量不相同。

4.求四对线对应

Δefg是C₁₊,C_1-的公共自极三角形，Δe'f'g'是C₂₊,C_2-的公共自极三角形,第一幅图像和第二幅图像与一个单应性H有关,从而C₁₊,C_1-与C₂₊,C_2-的公共自极三角形提供三对线对应，即公共自极三角形的三对对应边。线段ef与C_1-交于点n，线段eg与C₁₊交于点d，线段e'f'与C_2-交于点n'，线段e'g'与C₂₊交于点d'，根据射影变换保持共线性，第一幅图像中的线段nd与第二幅图像中的线段n'd'构成第四对线对应。

5.单应性估计

设四对线对应的齐次线坐标为:l_j＝[x_j1 y_j1 z_j1]^T和l_j'＝[x_j2 y_j2 z_j2]^T(j＝1,2,3,4)它们与一个单应性有关,并且满足关系式:t_jl_j＝H^Tl_j'，写成矩阵的形式Au＝0，这里u＝[h₁₁ h₁₂ h₁₃ h₂₁ h₂₂h₂₃ h₃₁ h₃₂ h₃₃]^T。然后用SVD分解就可以求出H，即得到单应性估计。

本发明优点：

(1)该靶标制作简单，对该靶标的物理尺度没有要求，无需知道球心在世界坐标系下的坐标。

(2)该靶标的图像边界点几乎可以全部提取，这样可以提高曲线拟合的精确度。

(3)将单应性估计用于中心折反射摄像机单位视球模型。

附图说明

图1是靶标在抛物折反射摄像机单位视球模型的示意图。

图2是两幅图像之间的单应性估计。

具体实施方式

本发明提供了一种利用一个空间球在抛物折反射摄像机单位视球模型下的投影，进行单应性估计的方法，如图1。一个空间球进行单应性估计需要经过以下步骤：从图像中提取镜面轮廓投影边缘点、球像的边缘点，然后使用最小二乘法拟合获得镜面轮廓投影方程、球像方程，根据像点及对拓像点之间的关系获得对拓像点，从而求解出球像的对拓球像方程。设Δefg是C₁₊,C_1-的公共自极三角形，e,f,g是其顶点，Δe'f'g'是C₂₊,C_2-的公共自极三角形，e',f',g'是其顶点。用奇异值分解矩阵因为矩阵Q₁与Q₂是相似矩阵，根据相似矩阵具有相同的特征值，但相同特征值对应的特征向量不相同，所以通过匹配特征值，得到三对线对应，即公共自极三角形的三对对应边，第四对线对应nd,n'd'由公共自极三角形与二次曲线的交点的连线得到，并且这对线对应位于公共自极三角形内。因为单应性矩阵有八个自由度，一对线对应提供关于单应性矩阵元素的两个方程，四对线对应提供八个关于单应性矩阵元素的方程，就足以估计单应性矩阵，接着，再通过公共自极三角形与二次曲线的交点的另一对连线wq,w'q'，验证本文求得的单应性矩阵是有效可行的。具体步骤如下：

1.拟合镜面轮廓投影方程及靶标投影方程

在中心折反射摄像机下，由于镜面轮廓可以被完全提取，利用MATLAB中的Canny算子提取镜面轮廓投影边缘点,并对每幅靶标图像进行边缘检测，提取靶标图像边缘点的像素坐标，然后对获取的边缘点利用最小二乘法拟合得到镜面轮廓投影方程、球像方程。

2.拟合球像的对拓球像

如图1，一个空间球Q在抛物折反射摄像机下的投影过程等价于经过单位视球的两步投影过程：第一步，球Q以单位视球球心O为透视中心的世界坐标系O-x_wy_wz_w投影到单位视球球面上形成两个平行圆S_i+和S_i-(i＝1,2表示拍摄两幅图像)，并且S_i+和S_i-以单位视球球心O对称，其中S_i+是可见的，S_i-是不可见的，图1上用i＝1为例。第二步，在透视坐标系O_c-x_cy_cz_c，以单位视球球面上的点O_c为虚拟摄像机的光心的将单位视球球面上的两个平行圆S_i+和S_i-分别投影到图像平面Π上，其中图像平面Π与光轴O_cO垂直交于主点p，即z_w,z_c轴与光轴O_cO重合，x_w,x_c轴和y_w,y_c轴与图像平面的u,v轴平行。所形成两条二次曲线C_i+和C_i-(图1，以i＝1为例)，S_i+的像C_i+称为球像，S_i-的像C_i-称为对拓球像，它们是对拓关系，并且C_i+是可见的，C_i-是不可见的，即下标+表示可见，-表示不可见。利用MATLAB中的Canny算子提取第一幅图像中的镜面轮廓投影的边缘点和两幅图像平面上的球像与对拓球像的边缘点，通过最小二乘法拟合得到相应的二次曲线方程。这里用C₀表示第1幅图像镜面轮廓的系数矩阵。本文为了简化表述，用相同字母表示曲线和它的系数矩阵。通过C₀可得到折反射摄像机内参数矩阵K的一个初始值K₀，其中

其中,r是纵横比的初始值，f是有效焦距的初始值，s是倾斜因子的初始值，[u₀ v₀1]^T是折反射摄像机主点的初始齐次坐标矩阵形式，记为p＝[u₀ v₀ 1]^T。从而得到绝对二次曲线的像ω的初始值这里:

取C_i+(i＝1,2)上的一组点m_ij+(i＝1,2,j＝1,2,…,n,n≥5)，则与它相对应的一组对拓像点为m_ij－(i＝1,2,j＝1,2,…,n,n≥5)。因为本文只讨论抛物折反射摄像机的情况，从而关于对拓像点的关系式如下：

这里m_ij+，m_ij－用齐次坐标表示。根据对拓像点的定义，点m_ij－在球像C_i+(i＝1,2)的对拓球像C_i-(i＝1,2)上，因此用最小二乘法拟合得到对拓球像C_i-的方程。

在抛物折反射摄像机下，球像及其对拓球像满足下列表达式：

其中

这里C'_i+是小圆S_i+的系数矩阵，C'_i-是对拓小圆S_-的系数矩阵，C_i+是球像的系数矩阵，C_i-是对拓球像的系数矩阵,[n_x n_y n_z]^T是小圆所在基础平面的单位法向量，“±”表示+和-的组合，即可见与不可见。

3.公共自极三角形的恢复

如图2，一个空间球Q在图像平面上的投影成两条分离的二次曲线，从不同位置拍摄球Q的两幅图像，设这两幅图像上的二次曲线C_i+,C_i-(i＝1,2)的公共极点极线是点m_i(i＝1,2)和线l_i(i＝1,2)，且满足如下关系式

这里，μ_i是非零比例因子，上式相减得

两边同时乘C_i+得

其中E是单位矩阵。由上式,我们得到C_i+和C_i－的公共极线l_i是的特征向量。

图2所示，设Δefg是C₁₊,C_1-的公共自极三角形，其中e,f,g是公共自极三角形顶点；Δe'f'g'是C₂₊,C_2-的公共自极三角形，e',f',g'是公共自极三角形顶点。第一幅图像和第二幅图像与一个单应性H有关,从而C₁₊与C₂₊，C_1-与C_2-满足如下关系式：

令通过整理(9)式,得到

由(14)式知，矩阵Q₁与Q₂是相似矩阵，根据相似矩阵具有相同的特征值，而相同特征值对应的特征向量不相同。设u₁是矩阵的特征向量，v₁是与特征向量u₁相对应的特征值矩阵，u₂是矩阵的特征向量，v₂是特征向量u₂对应的特征值矩阵，则有

Q₁u₁＝u₁v₁, (15)

这里u₁是可逆矩阵(实际上,在特征空间的维数小于3时,u₁将是奇异的)，从而有

由(10)—(13)式得

Q₂＝(H^T)^-1Q₁H^T， (17)

将(16)式代入(17)式有

由(18)式和(19)式得，矩阵Q₁与矩阵Q₂的特征向量也与该平面射影变换H有关。

4.求四对线对应

矩阵Q₁与Q₂都是非奇异的，计算出矩阵Q₁的特征值为α_i(i＝1,2,3)，l_i(i＝1,2,3)是对应于α_i(i＝1,2,3)的特征向量，矩阵Q₂的特征值为β_i(i＝1,2,3)，l'_i(i＝1,2,3)是对应于β_i(i＝1,2,3)的特征向量。根据相似矩阵的定义，由(14)式可得，矩阵Q₁与Q₂相似，再由相似矩阵的性质，矩阵Q₁与Q₂具有相同的特征值，但相同特征值对应的特征向量不相同，从而有α_i＝β_i(i＝1,2,3)且l'_i＝H^-Tl_i(i＝1,2,3)。因为矩阵Q₁的特征向量l_i(i＝1,2,3)是第一幅图像中公共自极三角形的公共极线，矩阵Q₂的特征向量l'_i(i＝1,2,3)是第二幅图像中公共自极三角形的公共极线，由图2知,C₁₊,C_1-与C₂₊,C_2-的公共自极三角形提供三对线对应，即公共自极三角形的三对对应边。

因为单应性矩阵有八个自由度，所以至少需要四对线对应来估计单应性矩阵。现在给出第四对线对应。如图2，C₁₊,C_1-的公共自极三角形是Δefg，C₂₊,C_2-的公共自极三角形是Δe'f'g'，线段ef与C_1-交于点n，线段eg与C₁₊交于点d，线段e'f'与C_2-交于点n'，线段e'g'与C₂₊交于点d'，根据射影变换保持共线性，第一幅图像中的线段nd与第二幅图像中的线段n'd'构成第四对线对应。

5.单应性估计

设四对线对应的齐次线坐标为:l_j＝[x_j1 y_j1 z_j1]^T和l_j'＝[x_j2 y_j2 z_j2]^T(j＝1,2,3,4)它们与一个单应性H有关,并且满足关系式:

t_jl_j＝H^Tl'_j， (21)

其中t_j是非零比例因子。令

从而由(21)式得

(23)用非齐次坐标(x'_j1＝x_j1/z_j1,y'_j1＝y_j1/z_j1)表示得

写成矩阵的形式

Au＝0, (25)

其中

u＝[h₁₁ h₁₂ h₁₃ h₂₁ h₂₂ h₂₃ h₃₁ h₃₂ h₃₃]^T。

然后用SVD分解就可以求出H，即单应估计矩阵。最后再通过另一对连线wq,w'q'验证本文求得的单应性矩阵是有效可行的。

实施例

本发明提出了一种利用一个空间球作为标定物，用两幅图像的球像与对拓球像的公共自极三角形进行单应性估计。本发明采用的实验模板结构示意图如图1所示。下面用一实例对本发明的实施方案做出更为详细的描述。

如图1所示，本发明采用的实验模板是一个空间球，记为Q。利用本发明的方法进行单应性估计的具体步骤如下：

1.拟合镜面轮廓投影方程、球像方程

用抛物折反射摄像机从不同位置拍摄一个空间球的两幅图像，读入图像，利用MATLAB中的Canny算子提取镜面轮廓投影边缘点,并对每幅靶标图像进行边缘检测，提取靶标图像边缘点的像素坐标，然后对获取的边缘点利用最小二乘法拟合得到镜面轮廓投影方程、球像方程。设镜面轮廓投影方程的系数矩阵为C₀，两幅图像中球像方程的系数矩阵分别为C_i+(i＝1,2)，其结果如下：

2.拟合球像的对拓球像

将(26)代入(1)和(2)可得的系数矩阵，结果如下：

先分别在球像C₁₊和C₂₊上取至少五个互异的点，再将所取的点和(29)分别带入(3)得对拓球像上点的坐标，通过最小二乘法拟合获得对拓球像C_1-和C_2-的估计，系数矩阵结果分别如下：

3.公共自极三角形的恢复

令由(27)、(28)和(30)、(31)得

4.求四对线对应

由图2知,C₁₊,C_1-与C₂₊,C_2-的公共自极三角形提供三对线对应，即公共自极三角形的三对对应边，从而分别对Q₁和Q₂奇异值分解得到公共自极三角形Δefg和Δe'f'g'的三对线对应：

ef＝[823.3631448892227 -334.9154785909649 1.000000000000000]^T， (34)

fg＝[-64.7370672847409 -284.3456442063169 1.000000000000000]^T， (35)

eg＝[253.1800257577513 401.9887420262962 1.000000000000000]^T； (36)

e'f'＝[811.1691178418838 -245.0685946933112 1.000000000000000]^T， (37)

f'g'＝[-8.533739513317000 -200.3926614508949 1.000000000000000]^T, (38)

e'g'＝[214.0920089477229 350.3431516307489 1.000000000000000]^T。 (39)

然后根据公共自极三角形与二次曲线的交点，得到第四对线对应，并且这对线对应位于公共自极三角形内：

nd＝[878.1409386246340 -2129.294804927160 1.000000000000000]^T， (40)

n'd'＝[903.1925050159500 -1407.426148819028 1.000000000000000]^T。 (41)

5.单应性估计

对(25)用SVD分解就可求得单应性矩阵

如图2，根据公共自极三角形Δefg的边ef与eg与二次曲线C₁₊,C_1-的交点为w,q，同理公共自极三角形Δe'f'g'的边e'f'与e'g'与二次曲线C₂₊,C_2-的交点为w',q'，这一对线对应位于公共自极三角形外，其齐次坐标矩阵：

wq＝[-504.9725768929506 577.0724065709180 1.000000000000000]^T, (43)

w'q'＝[-523.3735233141710 481.2032270143448 1.000000000000000]^T。 (44)

求出单应性矩阵以后，通过匹配线对应的齐次线坐标矩阵，得到

H^Tw'q'＝[-504.9725768927912 577.0724065704330 1.000000000000000]^T。 (45)

通过公式(43)和(45)可得，利用单应性矩阵矫正及求解匹配的线w'q'。

Claims (1)

1.一种利用球像的公共自极三角形进行单应性估计的方法，其特征在于由空间中的单个球作为靶标；所述方法的具体步骤包括：首先,从不同位置拍摄一个空间球的两幅图像，并且拟合出镜面轮廓投影方程、靶标投影方程及对拓球像方程；其次,分别计算两幅图像中球像与对拓球像构成的矩阵的特征值与特征向量；然后通过匹配特征值，求出三对线对应，这三对线对应即公共自极三角形的三对对应边；接着，连接公共自极三角形与二次曲线的交点，求得第四对线对应，这对线对应位于公共自极三角形内；最后，对这四对线对应采用SVD分解求出单应性矩阵；

(1)公共自极三角形的恢复

设两条共面二次曲线C_i+,C_i-，其中i＝1，2，的公共极点极线是点m_i(i＝1,2)和线l_i，其中i＝1,2，且满足如下关系式这里，μ_i是非零比例因子，通过化简得其中E是单位矩阵；由上式,得到C_i+和C_i－的公共极线l_i是的特征向量；

设Δefg是C₁₊,C_1-的公共自极三角形，e,f,g是顶点，Δe'f'g'是C₂₊,C_2-的公共自极三角形，e',f',g'是顶点；第一幅图像和第二幅图像与一个单应性H有关，从而C₁₊与C₂₊，C_1-与C_2-满足如下关系式：令从而有Q₂＝(H^T)^{- 1}Q₁H^T；故矩阵Q₁与Q₂是相似矩阵，根据相似矩阵具有相同的特征值，而相同特征值对应的特征向量不相同；

(2)求四对线对应

Δefg是C₁₊,C_1-的公共自极三角形，Δe'f'g'是C₂₊,C_2-的公共自极三角形,第一幅图像和第二幅图像与一个单应性H有关,从而C₁₊,C_1-与C₂₊,C_2-的公共自极三角形提供三对线对应，即公共自极三角形的三对对应边；线段ef与C_1-交于点n，线段eg与C₁₊交于点d，线段e'f'与C_2-交于点n'，线段e'g'与C₂₊交于点d'，根据射影变换保持共线性，第一幅图像中的线段nd与第二幅图像中的线段n'd'构成第四对线对应；

(3)单应性估计

设四对线对应的齐次线坐标为:l_j＝[x_j1 y_j1 z_j1]^T和l_j'＝[x_j2 y_j2 z_j2]^T(j＝1,2,3,4)它们与一个单应性有关,并且满足关系式:t_jl_j＝H^Tl_j'，写成矩阵的形式Au＝0，这里

u＝[h₁₁ h₁₂ h₁₃ h₂₁ h₂₂ h₂₃ h₃₁ h₃₂ h₃₃]^T；然后用SVD分解就以求出H，即单应性矩阵，又称单应性估计。