GST  >> Vol. 7 No. 3 (July 2019)

    鱼眼镜头的畸变纠正分析
    The Analysis of Fish-Eye Lens Distortion Calibration

  • 全文下载: PDF(828KB) HTML   XML   PP.139-149   DOI: 10.12677/GST.2019.73020  
  • 下载量: 344  浏览量: 1,361  

作者:  

赵 栋:沈阳市勘察测绘研究院有限公司,辽宁 沈阳

关键词:
鱼眼镜头相机标定影像纠正Fish-Eye Lens Camera Calibration Images Calibration

摘要:

本文所指超广角镜头是一种焦距为16 mm或更短的并且视角接近或等于180˚的镜头,俗称“鱼眼镜头”。为使镜头达到最大的摄影视角,这种摄影镜头的前镜片直径很短且呈抛物状向镜头前部凸出,与鱼的眼睛颇为相似,“鱼眼镜头”因此而得名。鱼眼镜头的测量范围大,但是由于其变形严重,并没有广泛应用在摄影测量中,如果能有一种针对鱼眼镜头拍摄影像行之有效的纠正方法,则可以在一些测量中应用。目前对数码相机标定理论的研究已经相当成熟,国内外许多摄影测量和计算机视觉专家也提出了各种行之有效的标定方法,然而,针对鱼眼镜头的标定算法和软件仍然不成熟。本文将重点分析比较现有的鱼眼镜头标定模型,并通过编程实现鱼眼镜头的标定算法,最后以立体量测方式说明鱼眼镜头的标定精度。

The ultra-wide-angle lens’ focal length is 16 mm or shorter, and the angle of view is close to or equal to 180 degrees, commonly known as “fish-eye lens”. In order to make the lens reach the maximum viewing angle, the front lens of this lens is short in diameter and parabolic, projecting toward the front of the lens, which is similar to the fish’s eye, so “fish-eye lens” is named. Fish-eye lens has a large measurement range, but also large distortion, so it is not widely used in photogrammetry. If there is an effective correction method for the fish-eye lens image, it can be used in some measurements. Nowadays, the study of digital camera calibration is already quite mature, and the photogrammetry and computer vision experts both at home and abroad put forward many kinds of effective calibration method. Although there is a lot of camera calibrations method, to the fish-eye lens it is still immature. In this paper, I will focus on analyzing and comparing the existing fish-eye lens calibration models for fish-eye lens and achieve the calibration algorithm by programming, and finally use the stereoscopic measurement way to illustrate the accuracy of fish-eye lens calibration.

1. 引言

二十世纪八十年代,得益于Brown对相机标定的研究结果,尤其是光束法平差以及附加参数的引入,很多研究者进行了各种改进,也相应提出了各种方法。比如Ziemann提出采用多步法进行相机标定,同时提出了将铅垂线和光束法平差相结合的两步法,即在标定过程中将畸变差和主距、像主点分开求解,通过反复迭代完成标定。目前摄影测量和计算机视觉界中使用比较多、精度较高的相机标定方法还包括室内三维控制场标定,大量实验证明,基于高精度室内三维控制场的标定精度是最高的,而基于平面格网的标定方法又是最方便的。

随着普通数码相机的普及,研究各种方便实用、灵活和较高精度的相机标定系统已经成为了当今摄影测量界和计算机视觉界的研究热点之一。并且问题主要集中在使用什么标定参考对象和采用什么算法上。如果通过事先标定或者在测量过程中通过简单的计算就能够获得准确的相机内参数,那么就能够方便和快速地测量。其中,利用同名点和约束条件进行相机标定的自标定算法和利用高精度控制点的自检校光束法平差是两种非常有效的标定方法;前者更加灵活方便,可以用于较低精度的在线(实时)测量,而后者则可以用于特殊环境下的高精度测量。

2. 研究趋势

根据目前的研究现状,找到能适应不同投影类型(包括中心投影),标定简单、可靠性高、精度高、算法快速的软件是主要的研究趋势。传统的标定模型的畸变项已经不适合用来矫正鱼眼影像的畸变,目前的讨论集中在用新的模型和畸变项取代传统的标定模型。其中,投影模型的确定是核心,除了使用不同的投影模型,一些文献还提到用物体几何特征,立体标定来确定相机内参数。用投影模型的文章例如:英向华等提出的一种基于球面透视投影约束的鱼眼镜头校正方法,此方法通过拟合一个与选取的离散点距离最近的一个圆曲面,辅以多项式的切向和径向畸变项来确定相机内参数。但其却缺陷是不同相机的投影方式各异,投影面可能并不是一个圆曲面,而是一个不规则曲面 [1] 。基于几何特征的标定方法通过提取特征物体并矫正到实际的几何特征来标定内参数,如邱志强等提出的用射影不变性纠正鱼眼镜头畸变,利用直线物体在纠正影像中也应该是直线这样一种约束来标定相机内参,这对特征的选取和提取精度要求较高 [2] 。还有通过立体标定来确定内参,这种方法对像点匹配精度要求很高。在目前成本低廉的非专业测量用相机大量普及的情况下,用其部分替代成本高昂的传统测量用相机可以大大降低测量成本,研究适用于非专业测量相机的标定的高精度和简便实用的算法可使普通相机在测量领域的应用更加广泛。

3. 标定模型

3.1. 经典标定模型

经典的模型以中心投影为基础,加上线性和非线性的畸变项对图像边缘的径向和切向畸变进行纠正。经典模型的求解首先通过二维直接线性变换(DLT)或者求同形矩阵给出外方位元素的初值,然后用光束平差进行参数的优化,最后得到可靠的相机内参数 [3] [4] [5] ,(1)式描述了经典模型所依赖的共线方程。

x + Δ x x 0 = C i { m 11 ( X X 0 ) + m 12 ( Y Y 0 ) + m 13 ( Z Z 0 ) } { m 31 ( X X 0 ) + m 32 ( Y Y 0 ) + m 33 ( Z Z 0 ) } y + Δ y y 0 = C j { m 11 ( X X 0 ) + m 12 ( Y Y 0 ) + m 13 ( Z Z 0 ) } { m 31 ( X X 0 ) + m 32 ( Y Y 0 ) + m 33 ( Z Z 0 ) } (1)

通过平面格网控制点将其转化为8个参数相关的二维DLT模型,如下:

x = L 1 X + L 2 Y + L 3 L 7 X + L 8 Y + 1 y = L 4 X + L 5 Y + L 6 L 7 X + L 8 Y + 1 (2)

式(2)中: L 1 , , L 8 为二维DLT模型的8个参数,由于8个参数相关,故得到的参数精度不够高,故只能作为迭代解算的初始值。由于切向畸变和径向畸变的影响对高精度测量不可忽略,如公式(3)所示Δx,Δy为像点坐标的畸变改正数,其中: k 1 , k 2 , k 3 , 为径向畸变系数,p1,p2为切向畸变系数。

Δ x = ( x x 0 ) ( k 1 r 2 + k 2 r 4 + ) + p 1 ( r 2 + 2 ( x x 0 ) 2 ) + 2 p 2 ( x x 0 ) ( y y 0 ) Δ y = ( y y 0 ) ( k 1 r 2 + k 2 r 4 + ) + p 2 ( r 2 + 2 ( y y 0 ) 2 ) + 2 p 1 ( x x 0 ) ( y y 0 ) (3)

3.2. 鱼眼相机成像模型

鱼眼相机由于其超大视场角的原因,不遵守透视投影,故物方点、投影中心和像方点不在一条直线上(除了与主光轴平行的那条光线),基于共线方程的经典模型不适用于鱼眼相机的成像模型,故有必要讨论鱼眼相机的投影模型的性质。

投影模型

人工选取的区域,并由计算机自动识别的圆心

i) 中心投影: r = f × tan ( θ ) (4)

ii) 立体投影: r = 2 × f × tan ( θ / 2 ) (5)

iii) 等距投影: r = f × θ (6)

iv) 等角投影: r = 2 × f × sin ( θ / 2 ) (7)

v) 正交投影: r = f × sin ( θ ) (8)

式中:f为焦距, θ 为入射角,r为投影半径。

五种投影模型中投影半径r与入射角 θ 的函数关系如图1所示。

Figure 1. Functional relationship between projection radius r and incident angle theta in the above five projection models

图1. 以上五种投影模型中投影半径r与入射角 θ 的函数关系

图2可以看出,P为物方点,经过投影中心,弯曲投射到了p,而不是沿直线到达 p

目前鱼眼镜头的成像模型大都可以用以上其中一种公式近似表达,等距投影是理想的球面投影,但是鱼眼镜头的等效投影面是一个近似球面的不规则曲面而非理想曲面,考虑到拟合度,可以选择更一般的多项式作为投影模型即:

r = k 1 × θ + k 2 × θ 2 + k 3 × θ 3 + k 4 × θ 4 + k 5 × θ 5 + + k n × θ n (9)

将以上5个典型的投影模型展开为泰勒多项式:

中心投影: r = f × ( θ + 1 3 × θ 3 + 2 15 × θ 5 + 17 315 × θ 7 + + k n × θ 2 n 1 ) (10)

Figure 2. Diagram of the relationship between projection radius and incident angle in projection model of fish-eye camera

图2. 鱼眼相机投影模型中投影半径与入射角的关系图示

立体投影: r = f × ( θ + 1 12 × θ 3 + 1 120 × θ 5 + 34 40320 × θ 7 + + k n × θ 2 n 1 ) (11)

等距投影: r = f × θ (12)

等角投影: r = f × ( θ 1 24 × θ 3 + 1 1920 × θ 5 2 645120 × θ 7 + + k n × θ 2 n 1 ) (13)

正交投影: r = f × ( θ 1 6 × θ 3 + 1 120 × θ 5 1 5040 × θ 7 + + k n × θ 2 n 1 ) (14)

由以上式子可以看出理想的鱼眼成像模型的泰勒展开式都是奇次项,在 θ 接近0时,三次项以上的项可以忽略不计。随着 θ 增大,高次项的作用越来越明显。以上投影模型的泰勒展开式呈现出高度的一致性,故可以选择一个一般的奇次多项式开拟合各种成像模型的参数,定为 r ( θ ) ,k1、k2、k3等代表成像模型系数。选择多少次的模型取决于控制点数量,精度和特征点提取的精度。对于多项式拟合,次数越高,理论上也就越精确,但是需要一定数量和足够高精度的观测值。否则,高次项没有意义。

鱼眼镜头也存在非线性的对称和非对称畸变,对称畸变可以由成像模型来修复,非对称畸变可以通过增加适当的参数来加以改正。

这里自然而然的从形式上给出不考虑非线性畸变的不同于经典模型的鱼眼相机成像模型:

u = r ( θ ) × cos ( φ ) × m u + u 0 (15)

v = r ( θ ) + sin ( φ ) × m v + v 0 (16)

以上两式中: ( u , r ) 为像点坐标, θ 为像点对应的光线入射角, φ 为像点与主点的连线与像坐标系中u方向轴的夹角, ( u 0 , v 0 ) 为主点坐标, m u 为影像水平方向上单位长度所包含的像素个数, m v 为影像垂直方向上单位长度所包含的像素个数 [6] [7] [8] [9] 。

4. 实验过程

4.1. 特征点提取

特征点的提取在模式识别中有很多方法例如:二值分割,边缘算子提取等,这里我们选用计算机视觉开源库里OPENCV里的Canny算子进行圆形物边缘的提取以便计算圆的中心。Canny算子是John F. Canny 1986年开发出来的一个多级边缘检测算法,实践证明,这种边缘检测算法是众多边缘检测算子中表现较突出的一种方法,下面我们详细阐述特征点提取的步骤:

1) 对原图进行中值滤波,然后用Canny算子提取边缘,如图3所示。

Figure 3. Edge binary images obtained by canny operator

图3. 使用canny算子得到的边缘二值图像

2) 人工选择区域并且由计算机将识别圆心,如图4所示。

Figure 4. Artificially selected area and automatically recognized by computer (green ring)

图4. 人工选取的区域,并由计算机自动识别的圆心(绿色环状物)

4.2. 内参数初始化

对于主点 ( u 0 , v 0 ) 的初始值,可以由像片的长度和高度来决定,即

u 0 = w i d t h 2 , v 0 = h e i g h t 2 , (17)

(17)式中: w i d h t 为影像以像素为单位的宽度, h e i g h t 为影像以像素为单位的高度,对于 m u , m v ,可以由像片大小,视场角和焦距算出,式中符号含义与前文相同

m u = w i d t h 2 / ( f × θ ) m v = h e i g h t 2 / ( f × θ ) (18)

对于5个投影参数k1、k2、k3、k4、k5,我们采用 r = k 1 × θ + k 2 × θ 3 进行拟合。 θ 为0˚~75˚的等间隔采样点(单位为弧度),投影半径 r = f × θ ,并且暂时忽略k3、k4、k5。得到至少2个方程,解出k1、k2作为整体迭代求解的初始值。

4.3. 外参数初始化

由于鱼眼影像的变形越靠近边缘越大,我们只将离主点较近的9个圆心点作为像方点用二维DLT进行外参数的初始化。

x = h 1 X + h 2 Y + h 3 h 7 X + h 8 Y + 1 y = h 4 X + h 5 Y + h 6 h 7 X + h 8 Y + 1 (19)

(X,Y)为控制点坐标,(x,y)为相应的像点坐标,当控制点大于等于4个时,将式(18)线性化,通过线性方程求解二维dlt的8个参数,由于8个参数是内外方为元素的函数值,由这8个参数解出外方位元素。当平面控制点的z坐标为0时,经典的共线方程可表示为式(19)。

x x 0 = f a 1 ( X X s ) + b 1 ( Y Y S ) + c 1 ( Z s ) a 3 ( X X s ) + b 3 ( Y Y S ) + c 3 ( Z s ) y y 0 = f a 2 ( X X s ) + b 2 ( Y Y S ) + c 2 ( Z s ) a 3 ( X X s ) + b 3 ( Y Y S ) + c 3 ( Z s ) (20)

将(20)转化成(21)

x = ( f a 1 λ a 3 λ x 0 ) X + ( f b 1 λ b 3 λ x 0 ) Y + ( x 0 f λ ( a 1 X s + b 1 Y s + c 1 Z s ) ) a 3 λ X b 3 λ Y + 1 y = ( f a 2 λ a 3 λ y 0 ) X + ( f b 2 λ b 3 λ y 0 ) Y + ( y 0 f λ ( a 2 X s + b 2 Y s + c 2 Z s ) ) a 3 λ X b 3 λ Y + 1 (21)

式中: λ = ( a 3 × X s + b 3 × Y s + c 3 × Z s )

比较式(19)与(21),可得出以下八式:

h 1 = f a 1 λ a 3 λ x 0 h 2 = f b 1 λ b 3 λ x 0 h 3 = x 0 f λ ( a 1 X S + b 1 Y S + c 1 Z s ) h 4 = f a 1 λ a 3 λ y 0 h 5 = f b 2 λ b 3 λ y 0 h 6 = y 0 f λ ( a 2 X S + b 2 Y S + c 2 Z s ) h 7 = a 3 λ h 8 = b 3 λ } (22)

由上式可得

h 1 h 7 x 0 f = a 1 λ h 2 h 8 x 0 f = b 1 λ h 4 h 7 y 0 f = a 2 λ h 5 h 8 y 0 f = b 2 λ h 7 = a 3 λ h 8 = b 3 λ } (23)

有上式最终可得到:

a 1 a 3 = h 1 h 7 x 0 f h 7 b 1 b 3 = h 2 h 8 x 0 f h 8 a 2 a 3 = h 4 h 7 y 0 f h 7 b 2 b 3 = h 5 h 8 y 0 f h 8 } (24)

主点 ( x 0 , y 0 ) 与4.2节提到的主点 ( u 0 , v 0 ) 等同,f在前文4.2节已有论述。

在以Y轴为主轴的转角系统下, tan ( κ ) = b 1 b 2 ,根据式(22)与式(24), tan ( κ ) = h 2 h 8 x 0 h 5 h 8 y 0 ,再由

b 1 2 + b 2 2 + b 3 2 = 1

b 3 2 = 1 1 + ( h 2 h 8 x 0 ) 2 f 2 h 8 2 + ( h 5 h 8 y 0 ) 2 f 2 h 8 2 (25)

通过 sin ( ω ) = b 3 可求得 ω ,再由旋转矩阵的正交性

( c 1 c 2 c 3 ) = ( a 1 a 2 a 3 ) × ( b 1 b 2 b 3 ) = ( a 2 b 3 a 3 b 2 a 3 b 1 a 1 b 3 a 1 b 2 a 2 b 1 ) (26)

tan ( ϕ ) = a 3 c 3 = 1 a 1 b 3 b 2 a 2 a 3 b 1 b 1 b 2 在上一步已求出, a 1 b 3 a 2 a 3 可由式(24)确定。

由式(22),可得到

h 3 = x 0 f λ ( a 1 X s + b 1 Y s + c 1 Z s ) h 6 = y 0 f λ ( a 2 X s + b 2 Y s + c 2 Z s ) λ = ( a 3 X s + b 3 Y s + c 3 Z s ) (27)

X s Y s Z s 的初值可通过解(27)方程组得到。

4.4. 整体求解

观测方程如下:

u = r ( θ ) × cos ( φ ) × m u + u 0 (28)

v = r ( θ ) × sin ( φ ) × m v + v 0 (29)

由4.2节求出的3个旋转角可得到 a 1 , a 2 , a 3 , b 1 , b 2 , b 3 ,加上3个平移参数 X s Y s Z s ,根据坐标齐次变换:

x = a 1 X + b 1 Y + X s y = a 2 X + b 2 Y + Y s z = a 3 X + b 3 Y + Z s (30)

(30)式(X,Y)为物方坐标点,(x,y,z)为近似像空间坐标点:

θ = acos ( z s q r t ( x 2 + y 2 + z 2 ) ) φ = asin ( y s q r t ( x 2 + y 2 ) ) (31)

式(19),(20)建立了物方点和像方点之间的关系,由上述关系,对于每一对物方点和像方点可得间接平差的误差方程:

v = A x L (32)

式(32)中:A为观测值对每个待求参数的偏导,v为残差,L为观测值减去近似值。对每一点对依次求出像点坐标对内参、外参、畸变系数的偏导数。但是这里对内外参求偏导数的解析式较为复杂。我们应用数值微分,通过8点的函数值求得近似偏导数,由数值微分理论,求偏导数的精度可达o(x4)。

最后我们用非线性邻域里较为可靠的优化方法L_M进行非线性最小二乘求解,此方法亦将“残差平方和最小”作为准则进行优化。

4.5. 影像纠正

4.5.1. 直接纠正

通过原影像坐标的反解,得到入射角 θ 和极角 φ ,再代入中心投影方程得到纠正后影像的坐标,由于重投影后,纠正影像的高度和宽度比原图像大,会有空白点,可以对其进行格网插值,选择双三次卷积进行插值。

具体步骤如下:

1) 依次由原图像的像点 ( u , v ) 得到球面坐标 ( θ , φ ) x = u u 0 m u y = v v 0 m v r = s q r t ( x 2 + y 2 ) 求解

多项式 r = k 1 × θ + k 2 × θ 3 + + k 5 × θ 9 的根,选择合理的值 φ = acos ( x / r ) ,注意 r = 0 的情形;

2) 得到的球面坐标 ( θ , φ ) 应用中心投影模型 r = f × tan ( θ ) { X = r × cos ( φ ) + X 0 Y = r × sin ( φ ) + Y 0 进行重投影,灰度

值的确定可采用双三次卷积。

4.5.2. 间接纠正

直接纠正会导致纠正图像上存有空白点,间接纠正是从纠正影像上的坐标点通过中心投影反算得到 ( θ , φ ) ,再应用模型正算映射到原影像,通过双三次卷积插值恢复其灰度。具体介绍其步骤:

由相机视场角估算出纠正影像的大小,并算出主点坐标 ( u 0 , v 0 ) ,依次求出纠正影像的像素坐标所对

应的 ( θ , φ ) ,其中 x = u u 0 m u y = v v 0 m v r = s q r t ( x 2 + y 2 ) θ = atan ( r f ) φ = acos ( x / r ) 。注意 r = 0 的情形应用解出的参数确定的模型: r = k 1 × θ + k 2 × θ 3 + + k 5 × θ 9 { X = r × cos ( φ ) + X 0 Y = r × sin ( φ ) + Y 0 ,确定纠正影像

的坐标与原影像的坐标映射关系,灰度插值采用双3次卷积插值。

5. 实验和结果

我们基于MFC框架设计了一个从数据输入到结果输出这样一个操作过程流程化的工具,对算法进行了实现,我们采用了University of Oulu的Juho Kannala和Sami S. Brandt提供的影像数据,其相机为Watec 221S CCD color camera,鱼眼镜头的焦距为1.178 mm,纠正结果如图5图6所示。

Figure 5. Pictures taken with fish-eye lenses

图5. 鱼眼镜头拍摄的影像

Figure 6. Corrected image of one of the channels

图6. 其中一个通道进行纠正后的影像

对于纠正影像的精度,我们放弃了用观测值残差来评价,而直接用特殊物体几何特征。我们选取纠正影像上共线的一组圆心坐标,进行直线拟合,得到的标准差为0.4157个像素,纠正精度达到了摄影测量的要求。

文章引用:
赵栋. 鱼眼镜头的畸变纠正分析[J]. 测绘科学技术, 2019, 7(3): 139-149. https://doi.org/10.12677/GST.2019.73020

参考文献

[1] 英向华, 胡占义. 一种基于球面透视投影约束的鱼眼镜头校正方法[J]. 计算机学报, 2003, 26(12): 1702-1708.
[2] 邱志强, 等. 用神经网络变易有效焦距的摄像机标定法[J]. 国防科技大学学报, 2002, 24(5): 16-19.
[3] Zhang, Z. (2000) A Flexible New Technique for Camera Calibration. IEEE Transactions on Pattern Analysis and Machine Intelligence, 22, 1330-1334.
https://doi.org/10.1109/34.888718
[4] 詹总谦. 基于纯平液晶显示器的相机标定方法与应用研究[D]: [博士学位论文]. 武汉: 武汉大学, 2006.
[5] 张永军, 张祖勋, 张剑清. 利用二维DLT及光束法平差进行数字摄像机标定[J]. 武汉大学学报(信息科学版), 2002, 27(6): 566-571.
[6] Kannala, J. and Brandt, S.S. (2006) A Generic Camera Model and Calibration Method for Conventional, Wide-Angle, and Fish-Eye Lenses. IEEE Transactions on Pattern Analysis and Machine Intelligence, 28, 1335-1340.
https://doi.org/10.1109/TPAMI.2006.153
[7] Thirthala, S. and Pollefeys, M. (2005) Multi-View Geometry of 1D Radial Cameras and Its Application to Omnidirectional Camera Calibration. Tenth IEEE International Conference on Computer Vision, Beijing, 17-21 October 2005, 1539-1546.
https://doi.org/10.1109/ICCV.2005.158
[8] Barreto, J.P. and Daniilidis, K. (2005) Fundamental Matrix for Cameras with Radial Distortion. Tenth IEEE International Conference on Computer Vision, Beijing, 17-21 October 2005, 625-632.
https://doi.org/10.1109/ICCV.2005.103
[9] Kannala, J. and Brandt, S. (2004) A Generic Camera Calibration Method for Fish-Eye Lenses. Proceedings of the 17th International Conference on Pattern Recognition, Cambridge, 26-26 August 2004, 10-13.
https://doi.org/10.1109/ICPR.2004.1333993