近年来，鱼眼相机由于其超大的视场范围，在全景视觉、机器人导航、虚拟现实及视觉监控等领域中得到广泛应用。然而大视场同时也带来严重的图像畸变，影响人眼直观的视觉感受以及图像信息的利用。为了对图像畸变进行校正，需要进行鱼眼相机标定^[1]。标定选用的成像模型方面，大多根据径向畸变和切向畸变建模，也有多项式逼近模型、抛物面投影模型、球面投影模型等等^[2-5]。标定方法方面，主要分为传统标定法、主动视觉法和自标定法。其中传统标定法标定过程相对复杂，但精度最高^[6-10]。鱼眼相机标定对精度要求较高，因而采用传统标定法。

用传统标定法对鱼眼相机进行标定时，由于鱼眼相机的大视场，单幅图像中标定板难以覆盖整幅图像从而导致标定误差增大。为弥补这一缺陷，一般采用多图标定法，将平面标定板放置于相机前不同的位置采集多张图像，以获得分布范围较大的标定原始数据。该方法需多次摆放棋盘并采集图像，对于某些需要快速安装并完成标定的场合并不适用，比如在生产线上批量装配车载相机时。针对这一问题，本文提出一种立体标定法，利用三块平面黑白棋盘搭建立体标定板，采集一幅图像即可实现快速标定。

1 鱼眼相机标定原理 1.1 成像模型

相机标定首先要建立成像模型。本文选用Scaramuzza提出的折反射模型^[5]作为鱼眼相机成像模型。

折反射模型坐标系设定如图 1。其中，(X, Y, Z)为世界坐标，(x, y, z)为相机坐标，(u", v")为感光面坐标，(u′, v′)为输出图像坐标。由世界坐标系至相机坐标系的映射模型为外部模型，与相机的安装位置有关，涉及参数称为外参；由相机坐标系至图像坐标系的映射模型为内部模型，与相机的内部构造有关，涉及参数称为内参。

世界坐标系中的一点P经相机成像落于感光面的P", 由相机坐标系原点O至P点的向量记为p, 折反射模型对于如何由P"逆投影至P做以下假设：

$\mathit{\boldsymbol{p=}}\left[ \begin{matrix} x \\ y \\ z \\ \end{matrix} \right]=\lambda \left[ \begin{matrix} {{u}''} \\ {{v}''} \\ f\left( \rho \right) \\ \end{matrix} \right]$

(1)

式中：λ为尺度变换因子，f(ρ)为包含相机畸变参数的非线性投影函数，用泰勒级数展开多项式表示为

$f\left( \rho \right)={{a}_{0}}+{{a}_{1}}\rho +{{a}_{2}}{{\rho }^{2}}+{{a}_{n}}{{\rho }^{n}}$

(2)

其中，$\rho =\sqrt{{{{{u}''}}^{2}}+{{{{v}''}}^{2}}}$表示投影点至成像面中心的距离；a₀，a₁，…，a_n为多项式系数，表示相机畸变参数。

成像面坐标到图像坐标的映射需经过仿射变换和中心偏移，关系式为

$\left[ \begin{matrix} {{u}'} \\ {{v}'} \\ \end{matrix} \right]=\left[ \begin{matrix} c&d \\ e&1 \\ \end{matrix} \right]\left[ \begin{matrix} {{u}''} \\ {{v}''} \\ \end{matrix} \right]+\left[ \begin{matrix} {{u}_{c}} \\ {{v}_{c}} \\ \end{matrix} \right]$

(3)

式中：c、d、e为仿射变换系数，u_c、v_c为图像中心偏移量。

以上参数为相机内参。

相机坐标(x, y, z)与世界坐标(X, Y, Z)间存在旋转和平移，其映射关系为

$\left[ \begin{matrix} x \\ y \\ z \\ \end{matrix} \right]=\left[ \begin{matrix} \mathit{\boldsymbol{R}}&\mathit{\boldsymbol{T}} \\ \end{matrix} \right]\left[ \begin{matrix} x \\ y \\ z \\ 1 \\ \end{matrix} \right]$

(4)

其中，旋转矩阵：

$\mathit{\boldsymbol{R}}=\left[ \begin{matrix} {{\mathit{\boldsymbol{r}}}_{1}}&{{\mathit{\boldsymbol{r}}}_{2}}&{{\mathit{\boldsymbol{r}}}_{3}} \\ \end{matrix} \right]=\left[ \begin{matrix} {{r}_{11}}&{{r}_{12}}&{{r}_{13}} \\ {{r}_{21}}&{{r}_{22}}&{{r}_{23}} \\ {{r}_{31}}&{{r}_{32}}&{{r}_{33}} \\ \end{matrix} \right]$

平移向量:

$\mathit{\boldsymbol{T}}={{\left[ \begin{matrix} {{t}_{1}}&{{t}_{2}}&{{t}_{3}} \\ \end{matrix} \right]}^{\rm{T}}}$

二者构成相机外参。

1.2 模型参数求取

要进行模型参数求取，首先需获取多个样点的世界坐标(X, Y, Z)和图像坐标(u′, v′)。对标定板拍摄图像，其上样点世界坐标已知，通过图像检测可获得样点图像坐标。

文献[4]中详细介绍了平面多图标定法的模型参数求解方法。参考文献中的方法，首先忽略仿射变换及图像中心偏移量，以(u, v)表示图像坐标，同时也可表示感光面坐标。以下标j代表样点序号。将式(1)两边叉乘并展开后得到如下公式：

$\begin{align} &\ \ \ \ \ \ {{v}_{j}}\left( {{r}_{31}}{{X}_{j}}+{{r}_{32}}{{Y}_{j}}+{{r}_{33}}{{Z}_{j}}+{{t}_{3}} \right)- \\ &f\left( {{\rho }_{j}} \right)\left( {{r}_{21}}{{X}_{j}}+{{r}_{22}}{{Y}_{j}}+{{r}_{23}}{{Z}_{j}}+{{t}_{2}} \right)=0 \\ \end{align}$

(5)

$\begin{align} &f\left( {{\rho }_{j}} \right)\left( {{r}_{11}}{{X}_{j}}+{{r}_{12}}{{Y}_{j}}+{{r}_{13}}{{Z}_{j}}+{{t}_{1}} \right)- \\ &\ \ \ \ \ {{u}_{j}}\left( {{r}_{31}}{{X}_{j}}+{{r}_{32}}{{Y}_{j}}+{{r}_{33}}{{Z}_{j}}+{{t}_{3}} \right)=0 \\ \end{align}$

(6)

$\begin{align} &\ \ {{u}_{j}}\left( {{r}_{21}}{{X}_{j}}+{{r}_{22}}{{Y}_{j}}+{{r}_{23}}{{Z}_{j}}+{{t}_{2}} \right)- \\ &{{v}_{j}}\left( {{r}_{11}}{{X}_{j}}+{{r}_{12}}{{Y}_{j}}+{{r}_{13}}{{Z}_{j}}+{{t}_{1}} \right)=0 \\ \end{align}$

(7)

其中式(7)为线性方程，以此为突破口进行参数求解。由于旋转矩阵R内部存在约束关系，不能直接简单地对式(7)求解，需改进求解算法。

首先从标定板的三个平面中选取样点较多的一个平面，利用该平面内的样点坐标求外参初值。因为无论选定哪个平面，都只涉及r₁、r₂、r₃中的某两个向量，不存在内部约束，可简单求解。然后令所有样点的三维坐标都参与计算，对外参进行优化。

以下按选定X-O-Y面为例进行说明。在该平面上，所有点的坐标Z=0，式(7)简化为

${{u}_{j}}\left( {{r}_{21}}{{X}_{j}}+{{r}_{22}}{{Y}_{j}}+{{t}_{2}} \right)-{{v}_{j}}\left( {{r}_{11}}{{X}_{j}}+{{r}_{12}}{{Y}_{j}}+{{t}_{1}} \right)=0$

(8)

该式中包含r₁₁、r₁₂、r₂₁、r₂₂、t₁、t₂六个参数，可构建超定线性方程组并求解^[11]。然后根据旋转矩阵R内部存在的约束关系，即r₁、r₂、r₃两两相互正交，且|| r₁||=|| r₂||=1，求得除t₃外的其他外参。

然后用所有样点数据根据式(7)对外参进行最小二乘迭代优化^[11-13]。迭代步长为

${{\delta }_{k}}={{\left( {{\mathit{\boldsymbol{J}}}_{k}}^{\rm{T}}{{\mathit{\boldsymbol{J}}}_{k}} \right)}^{-1}}{{\mathit{\boldsymbol{J}}}_{k}}^{\rm{T}}{{\varepsilon }_{k}}$

(9)

式中：ε表示上一次的估计误差，J表示Jacob矩阵，下标k代表第k次迭代。以函数g表示式(7)左部，Jacob矩阵如式(10)所示，式中n代表样点数量：

$\mathit{\boldsymbol{J}}=\left[ \begin{matrix} \frac{\partial {{g}_{1}}}{\partial {{r}_{11}}}&\frac{\partial {{g}_{1}}}{\partial {{r}_{12}}}&\frac{\partial {{g}_{1}}}{\partial {{r}_{21}}}&\frac{\partial {{g}_{1}}}{\partial {{r}_{22}}}&\frac{\partial {{g}_{1}}}{\partial {{t}_{1}}}&\frac{\partial {{g}_{1}}}{\partial {{t}_{2}}} \\ \frac{\partial {{g}_{2}}}{\partial {{r}_{11}}}&\frac{\partial {{g}_{2}}}{\partial {{r}_{12}}}&\frac{\partial {{g}_{2}}}{\partial {{r}_{21}}}&\frac{\partial {{g}_{2}}}{\partial {{r}_{22}}}&\frac{\partial {{g}_{2}}}{\partial {{t}_{1}}}&\frac{\partial {{g}_{2}}}{\partial {{t}_{2}}} \\ \vdots &\vdots &\vdots &\vdots &\vdots &\vdots \\ \frac{\partial {{g}_{n}}}{\partial {{r}_{11}}}&\frac{\partial {{g}_{n}}}{\partial {{r}_{12}}}&\frac{\partial {{g}_{n}}}{\partial {{r}_{21}}}&\frac{\partial {{g}_{n}}}{\partial {{r}_{22}}}&\frac{\partial {{g}_{n}}}{\partial {{t}_{1}}}&\frac{\partial {{g}_{n}}}{\partial {{t}_{2}}} \\ \end{matrix} \right]$

(10)

每次迭代过程中，r₁₁、r₁₂、r₂₁、r₂₂、t₁、t₂由式(9)计算迭代步长后更新，而r₁₃、r₂₃则由更新后的r₁₁、r₁₂、r₂₁、r₂₂计算得到。

外参优化完成后，利用式(5)、(6)，采用最小二乘法^[12]求解畸变参数a₀，a₁，…，a_n以及t₃。

最后，采用Levenberg-Marquardt迭代优化算法^[14-15]对所有参数进行优化。其中，c、d、e初值为1、0、0，u_c、v_c初值分别为图像高宽的一半。至此，所有参数求解完成。

2 鱼眼相机快速标定方法 2.1 标定板设计

立体标定板由两两互相垂直的三块黑白棋盘模板构成，如图 2所示。此标定板有如下优点：棋盘格成像后角点(即黑白方格的顶点，作为参数求取所需的样点)邻域内灰度梯度较大，易检测；立体标定板可包围相机，使棋盘角点遍布整幅图像；平面交界处绘有红色边界线，方便角点检测前进行所属平面划分，确定三维世界坐标。

2.2 标定过程

1)图像采集

将标定板置于待标定相机前，调整标定板位置及角度，使图像中遍布棋盘角点，采集标定图像一张，如图 3。

2)图像分割与角点检测

标定图像中没有复杂的干扰物，只有黑白方格与红色曲线，因而RGB三通道中，边界线处像素点的R通道灰度值明显高于G、B通道，其他区域三通道灰度值接近。求R通道与G或B通道的灰度差，采用自适应阈值进行二值化处理，可将边界线检测出来^[16]，结果如图 4。

对三条边界线分别进行二次多项式拟合，依据边界线将原始图像分割为三个区域，形成三幅子图像，如图 5。

对三幅子图像分别进行角点检测并排序^{[3, 6]}。由于角点在子图像中的坐标与其在原始图像的坐标相同，因而可以得到原始图像中所有角点的图像坐标，如图 6所示。

3)坐标映射建立

各区域的角点都有一维世界坐标为0，对应图 5中三幅子图像分别为Y=0，X=0，Z=0。根据棋盘方格边长以及与世界坐标原点间隔的方格数，可得到所有角点的世界坐标。从而建立起二维图像坐标与三维世界坐标的一一映射，用于模型参数的求解。

4)参数求解

依据标定原理，求解鱼眼相机参数。

以上标定过程在Visual Studio 2010平台下，用C/C++语言结合OpenCV编写标定程序，输入标定图像，得到鱼眼相机参数。

3 实验验证

图 3的标定结果如下

$\begin{align} &\ \ \ \ \ \ \ \ \ \ {{a}_{0}}=-1.5362\times {{10}^{2}},{{a}_{1}}=0, \\ &{{a}_{2}}=1.5294\times {{10}^{-3}},{{a}_{3}}=-5.1205\times {{10}^{-6}} \\ &\ \ \ \ \ \ \ {{a}_{4}}=1.8467\times {{10}^{-8}},c=0.9992, \\ &d=-1.3012\times {{10}^{-5}},e=2.7501\times {{10}^{-3}} \\ &\ \ \ \ \ \ \ \ {{u}_{c}}=243.5637,{{v}_{c}}=323.4531 \\ \end{align}$

$\begin{align} &\ \ \ \ \ \mathit{\boldsymbol{R=}}\left[ \begin{matrix} 0.3657&0.3410&-0.8660 \\ -0.7142&0.6994&-0.0262 \\ 0.5968&0.6281&0.4780 \\ \end{matrix} \right] \\ &{{\mathit{\boldsymbol{T}}}^{\rm{T}}}=\left[ \begin{matrix} 30.5170&-3.8548&-125.0822 \\ \end{matrix} \right] \\ \end{align}$

上述相机参数中，只有平移向量T可以简单测得，旋转矩阵R由于旋转角较难测量而无法得到真实值，内参由于是抽象模型参数也无法获取真实值。因而难以通过所得参数本身来评价标定精度。

本文设计了四组对比实验对标定结果进行评价。

3.1 实验设计 3.1.1 标定实验设计

设计四组标定实验如下：

A、B两组采用立体标定法，标定板可在保证棋盘角点遍布整幅图像的前提下随意摆放；C、D两组采用平面多图标定法，各自用平面棋盘拍摄八幅图像进行标定，其中C组严格控制棋盘摆放位置，使棋盘绕相机光轴等角度旋转一周，D组棋盘摆放较为随意, 如图 7所示。

3.1.2 验证实验设计

由于四组标定所获图像的场景各不相同，要进行结果对比，需设计一个统一的验证实验。在杂物较少的水平地面上摆放一张棋盘，棋盘长约1 m，宽约0.6 m，相机距地面高度1.2 m。采集鱼眼图像作为测试图像，如图 8。

结合标定实验所得四组相机内参，对测试图像标定外参。此过程只需通过求解超定线性方程组即可得到外参^[15]。由此获得同一场景的四组相机内外参数，用于对比分析。

3.2 实验结果及分析 3.2.1 鱼眼图像畸变校正

首先利用四组参数对测试图像进行畸变校正，从视觉效果上进行主观评价。为了便于分析校正效果，在校正基础上进行俯视投影。结果如图 9所示。

在校正效果图中，A、B、C组的整体校正效果较好，将棋盘及地面的方格较好的还原；D组还原效果稍差，图像左下角和右下角可看出弯曲变形未完全校正。

3.2.2 世界坐标估计

利用标定所得相机参数，用角点的图像坐标逆投影得到其世界坐标的估计值，计算估计坐标与真实坐标的距离作为估计误差，如图 10。世界坐标估计误差分析见表 1。

由以上数据可看出四组实验中世界坐标的估计误差相近。粗略计算世界坐标估计值的相对误差，棋盘中心与相机距离取1 300 mm，估计误差取3 mm，则相对误差为0.23%，可认为误差较小，说明四组标定都有较高的标定精度。

其中，C组标定精度最高，A组与B组精度极为接近，D组精度相对略低。分析原因，A组与B组采用立体标定板，虽然标定板位置随意摆放，但都保证了棋盘角点遍布整幅图像，数据分布较好，因而求得的参数精度较高；C组标定板摆放位置经过精心设计以保证角点分布较均匀且遍布图像，D组为随意放置，角点分布散乱，无法保证获取的角点数据的分布情况，因而出现标定精度上的差异。

在标定过程的复杂程度上，A、B组采用的立体标定法，只需采集一幅图像，且因其标定板的立体设计，使得摆放棋盘时很容易将角点覆盖整幅图像，极易操作；C、D组采用的平面多图标定法，需多次采集图像，若要实现高精度标定，需要兼顾每幅图像中的棋盘位置，且若无辅助设备则人工操作较为繁琐。

4 结束语

本文提出了一种利用立体标定板对鱼眼相机进行快速标定的方法，根据鱼眼相机成像特点设计专用的立体标定板，并提出适用于立体标定的参数求解方法。经实验验证，与传统的平面标定法相比，该方法具有较高的标定精度和鲁棒性，同时简化了标定流程，操作简易，耗时较短。对于车载鱼眼相机安装等需要进行批量标定的场合，该方法具有实际应用意义。

无论是本文所述方法还是传统标定法，由于鱼眼图像边缘区域畸变较大，与图像中心区域相比样点较密集且棋盘方格变形较大，样点检测精度降低，对标定参数的精度存在一定影响。本文认为可以尝试样点非均匀分布的标定板，改善边缘样点检测精度以提高相机参数标定精度。