相机检校是摄影测量及计算机视觉至关重要的一个问题，相机检校是否精确，会影响到后续测图及三维重建的精度。影响相机检校结果的因素很多，包括相机自身因素和外部因素等，如相机制造工艺、相机焦距、调节光圈、更换镜头、相机晃动、更换电池、光照条件、检校靶标数量、靶标点提取精度、安装角、测量误差、相机模型选择等^[1-4]，因此在进行摄影测量时，在获得相机的初始参数后，还需进行相机检校，减少相机参数对测量精度的影响^[5-6]。摄影测量领域的传统方法通常是在专用检校场(鉴定场)进行相机检校^[7-8]，但检校场建造维护代价很高，且要综合考虑气候条件、地形地貌等的影响。为了满足相机鉴定的迫切需求，使相机鉴定更加方便灵活，文献[9-12]使用靶标进行相机鉴定。基于靶标的相机鉴定结果精度很高，稳健性强，但需要制作特定靶标，而且需要精密仪器测量世界坐标。

最好、最适宜的相机检校方法取决于实际应用类型^[13]，立体相机根据两个相机之间固定相对关系建立约束进行标定^[14]，存在“空间平行线”等建筑物场景时，则可利用灭点进行相机检校^[15]。为了真实精确反映航空摄影测量时所获取图像的相机参数和畸变系数，提高解算精度，进行精确三维建模。本文提出基于航空影像的相机自检校方法，联合所有航空影像同名点基于共面约束进行相机检校。对于高分辨率大尺寸航空影像，中心及边缘的畸变差异较大，为了提高畸变参数解算精度，提出分区域网格划分检校法。

1 相机检校

航空影像是现代信息化战场的重要信息来源，为了提高数据处理精度，保证数据的可信性，关键一步就是对相机进行严格的检校，求出相机的内方位元素及各项畸变系数。在数据处理中消除影像的畸变差，提高同名像点的空间解算精度，实现高精度三维重建。本文提出基于航空影像的相机自检校方法，具体步骤如下：首先进行航空影像同名点匹配，构建立体像对；然后基于共面约束使用直接解法和迭代优化进行相对定向，解算相机位置与姿态；最后对航空影像进行区域划分，使用最小二乘优化方法解算相机内方位元素和畸变系数。由于该方法联合所有影像同名点进行解算，为了避免大量同名点对存储和计算速度的影响，优化了矩阵解算策略，使方程解算与同名点数量无关。该方法操作简单，对图像无特殊要求，不需要制作精良靶标，也不需要使用特定检校场。使用大比例尺航空影像进行解算能真实精确反映摄影测量时所获取图像的相机参数和畸变系数，避免检校环境与使用环境不同解算得到的相机畸变参数不能真实反映所获取影像的畸变问题；使用所有同名点解算，避免由于选择不同特征点或控制点对检校精度的影响^[16]，同时，焦距、像主点和相机畸变系数可以独立解算，降低了参数之间的耦合度^[12]。

在进行相机检校前首先必须根据实际应用选择合理的相机畸变模型，模型的选择要比参数估计更重要，同时还需要选择合适的评价指标衡量选择的模型^[17]。本文研究了小孔模型、经典的畸变模型及在x、y方向使用不同的径向畸变系数的畸变模型，同时提出区域网格划分法进行畸变解算，并对检校结果进行了比较。同名点误差及共面约束误差对三维重建精度影响很大，因此本文选择同名点误差及共面误差衡量各模型的性能。

1.1 相机模型

被广泛使用的相机模型为小孔成像模型，如式所示

式中，s为空间点深度或尺度因子；(u, v)为以像素为单位的图像坐标；(X_w, Y_w, Z_w)为世界坐标系中的坐标；(u₀, v₀)为像主点坐标；α是u轴和v轴不垂直因子；(f_u, f_v)分别是相机以像素为单位的横坐标轴和纵坐标轴上的焦距；R为旋转矩阵；t为平移向量。

本文设u轴和v轴具有相同的尺度因子，即f_u=f_v，并假设成像平面的横、纵轴垂直，则最终的相机投影变换模型为

式中，f为等效焦距；K为内方位元素矩阵或投影变换矩阵；M为外参数矩阵，包含旋转矩阵R和平移向量t。

由于制造工艺等多种因素，相机通常都会存在畸变，如式(3)所示

式中，(x, y)为中心化的像点坐标，x=x₀-c_x, y=y₀-c_y；(x₀, y₀)为实际像点坐标；c_x、c_y为光心坐标；r为像点到光心的距离，；k₁、k₂为径向畸变系数；p₁、p₂为偏心畸变系数；s₁、s₂为薄棱镜畸变系数。

对于高分辨率大尺寸航空影像，x、y方向的径向畸变系数可能并不相同，因此在进行解算时，设x、y方向具有不同的径向畸变系数k_x1、k_x2及k_y1、k_y2，如式(4)所示

由于径向畸变的影响远大于其他畸变^[18-19]，因此本文只考虑径向畸变。

1.2 相机检校

1.2.1 同名点匹配

同名点匹配是进行相机检校的基础。输入图像为无序的遥感影像集，因此必须首先提取特征点，进行同名点匹配，进而确定各影像间的邻接关系，建立立体模型。从无序影像集构建立体模型的策略如下：首先提取所有影像的特征点；对于任意两幅影像，进行特征匹配，因为从不同视点拍摄的两幅影像存在缩放、旋转等仿射变换，基于RANSAC算法求解它们之间的仿射变换关系，以判断匹配的有效性：若两幅图像无重叠或误点太多，则无法得到一个有效的仿射变换，该像对被剔除；对于有重叠并匹配成功的两幅影像，通过RANSAC算法可以剔除少量误匹配点，最终获得一定数量的可靠同名像点，对应的两幅影像即构成了一个立体像对^[20]。在立体像对模型内进行稠密同名点匹配，获取大量的观测值，可以为后续的最小二乘优化提供大量的观测方程。

1.2.2 相对定向

根据解析法相对定向，同名射线在空间对对相交，射线op、o′p′及摄影基线B共面，如图 1所示，根据三矢量共面，它们的混合积为0，可得约束公式

式中，op、o′p′为同名射线矢量；B为基线矢量。

首先采用直接解法相对定向获取各相机位置与姿态初值，然后使用最小二乘进行迭代优化，依据共面约束，建立优化函数如式(6)所示

式中，F=B·op×o′p′。

具体过程如下：①选择系列图像中中间位置的图像作为主图，相邻的一幅图作为辅图，邻接图为相互有重叠的图，从其他未定位图中找一幅邻接点最多的邻接图作为候选图。通过提取主图与邻接图间的同名点，由相对定向过程，即可确定邻接图对应相机的相对位置和姿态；重复上述过程，即可确定全部相机的相对位置和姿态。②以位置与姿态参数为优化目标，基于共面约束建立优化函数，解算各相机位置与姿态参数。对于多幅航空影像，为了快速解算位置与姿态，迭代解算能够收敛，避免参数过多相互干扰，采用如下策略：对主、辅图定位后，主、辅图只绕Z轴旋转κ₁、κ₂，不进行平移变换，其他图像进行平移调整(t_xi, t_yi, t_zi)，同时绕Z轴旋转κ_i，经过迭代调整后使同名像点重合。此次位置与姿态解算不加入相机内方位元素。主、辅图各对应1个位置与姿态参数，其他图像对应4个位置与姿态参数，对于M幅图像N对同名点，则可生成包含(4M-6)个未知数N个方程的线性约束方程组，如式(7)所示

式中

经过直接解算与迭代优化后可以解算出相机相对位置与姿态。

1.2.3 相机内方位元素解算

在相对定向解算的位置与姿态基础上加入相机内方位元素(焦距、像主点、畸变系数)进行联合解算，进一步优化位置与姿态参数，解算相机内方位元素。约束方程和优化目标不变，以相机相对于公共坐标系的旋转角参数、平移向量及相机内方位元素为未知量，作为最小二乘参数迭代优化，如式(8)所示

式中

在调整过程中，主图绕X、Y轴旋转ω₁、φ₁，辅图绕3个轴旋转ω₂、φ₂、κ₂，其他图像进行平移调整(t_xi, t_yi, t_zi)，同时绕3个轴旋转ω_i、φ_i、κ_i。主、辅图对应5个位置与姿态参数，其他图像对应6个位置与姿态参数，共有7个相机参数(f, c_x, c_y, k_x1, k_x2, k_y1, k_y2)，对于M幅图像N对同名点，则可生成包含6M个未知数、N个方程的线性约束方程组。

经过最小二乘迭代优化，可以解算出相机的内方位元素和畸变系数及更精确的位置与姿态参数。

1.2.4 分区域畸变校正

对于高分辨率大尺寸航空影像，图像中心及边缘的畸变差异较大，为了准确地解算相机畸变系数，进行畸变校正，对图像进行了区域网格划分，分区域进行畸变校正。将航空影像划分为m×n矩形网格，m、n分别为x、y方向的网格数。联合所有影像同名点使用双线性插值和最小二乘迭代优化算法解算网格交点畸变量，然后根据网格交点处畸变量校正该网格区域内像点。通过联合计算所有航空影像并通过合理划分矩形网格，能够最真实地反映摄影时相机的内方位元素及畸变系数，对像点进行畸变校正。

在图 2中，设网格(i, j)交点p_ij的畸变量为δ_xp_{i, j}、δ_yp_{i, j}，则根据双线性插值，同名像点Q_k处的畸变量为

式中，u、v分别为Q_k到直线P_{i+1, j}P_{i+1, j+1}和P_{i, j+1}P_{i+1, j+1}的归一化距离。然后进行最小二乘优化迭代，在位置与姿态参数、相机内方位元素的基础上加入所有网格交点处畸变量进行迭代优化。x、y方向的网格数分别为m、n，则网格交点畸变量参数个数为2×m×n，对于N对同名点，构成的线性约束方程组如式(10)所示

式中

2.2.5 优化策略

联合所有航空影像大规模同名点进行解算，对存储和计算要求较高。本文采用最小二乘方法进行迭代优化，方程AX=B与A^TAX=A^TB同解，为了节省存储空间，不保存系数矩阵A而直接构造A^TA，方程右边直接构造A^TB，对于式(8)最终所需的存储单元仅为6M×(6M+1)。常规方法需要的存储单元为2N×(6M+1)，这样的处理使此方法的方程解算规模与引入的同名点数量无关，从而极大地减少了存储量，提高了计算速度。多次迭代直至||δX₂||＜ε(ε为某一小值)，进而得到未知量的最优解。

2 试验 2.1 试验数据

试验数据1采用威海天福山试验场拍摄的航空影像^[21]，试验场面积为2 km²，场区内均匀布设了183个精确量测的地面控制点，每个控制点都制作了明显的地面标志，控制点坐标精度优于0.01 m。用有人驾驶飞机对场区进行航拍，得到了5条航线、共115幅航空影像，原始影像数据量约为20 GB；图像分辨率为8956×6708像素。图像的航向重叠率约为60%，旁向重叠率大于20%，影像空间拍摄相机为哈苏H4D-60，镜头焦距为50 mm，像元大小为6 μm。无POS记录。图 3(a)为试验用的部分影像。

试验数据2采用使用佳能5DⅡ相机对青岛王家村进行航空摄影拍摄的航空影像，全航区共拍摄了78(13×6航线)幅影像，镜头焦距为35 mm，像元大小为6.4 μm，图像大小为5616×3744像素。部分试验数据如图 3(b)所示。

2.2 试验结果

从数据1威海天福山试验场的22幅航空影像64个同名点文件中匹配14幅图像32个同名点文件进行了相机内方位元素及畸变联合解算，同名点对数量约为640 000对。参数初始值(区域网格划分时网格交点畸变量初值为0)及试验结果如表 1所示。畸变量(δ_x，δ_y)为位于图片角上的像点在x、y方向的径向畸变，用于评估畸变大小。共面误差(MSE)为左像同名点和基线构成平面与右像同名点的距离的均方误差。方法1和方法2为文献[18]基于室外检校场联立全部影像和4幅影像的检校结果；方法*-1(*=1, 2, 3)为本文方法单独解算焦距(f)和像主点(c_x, c_y)；方法*-2(*=1, 2, 3)为使用传统畸变模型解算相机内方位元素和畸变系数；方法*-3为x、y方向使用不同的径向畸变参数解算相机内方位元素和畸变系数；方法*-4、*-5和*-6为采用区域网格划分进行畸变解算，网格数分别为16×12，32×24，64×48。

从表 1试验结果可以看出，本文方法可以有效求出相机内方位元素，且精度较高。像主点解算结果非常稳定，与基于室外检校场的解算精度相当。基于区域网格划分的畸变校正方法共面误差明显小于其他畸变校正方法。相机焦距的检校结果与初始值及文献[18]相差较大，经分析主要原因如下：一是由于相机镜头加工、装配过程中引入的误差，使实际值与标称值存在误差；二是在不同的环境下拍摄，摄影测量时外部条件变化等引入的误差^{[1, 22-27]}；如不同的摄影距离下主距和像主点也会发生变化，不同的检校场、图像集检校结果也会有差异^[27-28]。

表 2为本文方法焦距和结果分析，可以看出该方法对优化算法初值不敏感，解算的结果非常稳定。表 3和图 4为各图定向后同名点误差，使用区域网格划分的畸变解算算法精度明显高于其他方法，同名点误差显著减少，解算结果对初值不敏感，所有最小二乘迭代优化均小于7次收敛。

表 4-表 6和图 5为畸变校正前后共面误差，采用本文方法进行畸变校正后共面误差显著减小，使用分区域畸变校正后，共面误差在0.4像素以内。

为了进一步验证本文方法对普通数码相机的标定精度，使用数据2中10幅图像进行了试验，试验结果如表 7所示。从试验结果可以看出，该方法对普通数码相机标定精度也很高。

3 结论

本文提出了基于大比例尺航空影像利用共面约束的相机自检校方法。与传统方法相比，该方法操作简单，精度较高，能真实反映航空摄影测量时相机参数。通过对高分辨率大尺寸航空影像进行区域划分解算相机畸变量，进一步提高了解算精度。通过实际图像运算结果可以看出，畸变校正后同名点匹配精度有了显著提高，共面误差在0.4像素以内，该方法已应用于大规模航空影像拼接，进行三维重建。