POS数据支持下的面阵航测影像GPU几何纠正方法

马远征，方 勇

1.长安大学，陕西 西安，710061;

2.西安测绘研究所，陕西 西安，710054

1 引 言

相比于模拟相机，面阵航测相机能够直接获取数字影像而无需数字化，摄影效率高、图像质量好[1]，可作为应急救灾数据源。当重大灾害发生时，其成为决策者快速了解灾情、及时制定救灾策略的重要依据。而数字航测影像往往数据规模大，迫切需要开发更为高效的影像数据处理算法，以满足应急测绘需要。随着图形处理器(graphics processing unit，GPU)计算性能的不断提升[2]以及定位定向系统(position orientation system，POS)数据在影像直接地理定位上的广泛应用，使得摄影数据快速处理成为可能。就此，本文对POS数据支持下的面阵航测影像快速几何纠正技术展开研究。

传统的影像几何纠正方法，大多运行在中央处理器(central processing unit，CPU)上，其效率提升主要依靠多CPU并行技术，但由于CPU时钟周期长期局限于4GHz，且多CPU散热问题限制了多核并行技术的进一步发展。为提升影像数据处理效率，许多学者开始尝试利用CPU-GPU协同技术处理遥感影像。文献[3]将CPU-GPU协同处理技术用于卫星影像正射校正;文献[4]分析了几种不同GPU架构下的影像正射校正效率;文献[5]利用CPU-GPU协同处理技术，使得影像正射校正效率提升了10.91倍。数字航测相机与POS系统集成，可实现数字航测影像快速地理定位。文献[6]指出经后处理的POS数据可直接作为影像的外方位元素;文献[7]设计基于POS数据的高光谱影像CPU快速几何纠正算法;文献[8]开展了基于POS数据的机载航空影像几何纠正试验;文献[9]基于POS数据使用ERDAS软件进行航空遥感影像正射校正;文献[10]提出结合POS数据、地形数据库和GDAL开源图像库的航空影像快速几何纠正方案。相机畸变也是数字航测影像几何纠正的一个重要方面。文献[11]分析了面阵航测相机的几何畸变误差传递模型;文献[12]利用POS计算影像外方位元素，对推扫式成像方式的OMIS影像进行几何畸变纠正;文献[13]认为多面阵组合相机的光学畸变差会导致像点相对于其理想位置发生径向和切向的偏移，并推导了数码相机畸变校正公式。

上述研究从不同角度论证基于GPU编程框架和POS数据对面阵航测影像几何纠正的可行性，并开展了许多有益的工作。但是，这些文献大多仅侧重于某一技术层面。如CPU-GPU协同处理方面的文献往往侧重于分析影像几何纠正算法的实现效率，而没有涉及对校正精度的详细分析;基于POS数据进行影像几何纠正的文献关注的重点是影像外方位元素求解精度及影像几何校正方法的精度，而较少考虑纠正效率。而在应急测绘保障需求条件下，纠正效率和纠正精度同等重要。因此，本文基于POS数据建立无控制点面阵影像对地定位关系，设计面阵航测影像GPU几何纠正算法;并以某新型国产面阵数字航测相机为例，验证所提方法的效率和精度。

2 原理与方法

2.1 GPU几何纠正原理

1)无控严密对地定位原理

由于影像倾斜和地形起伏等因素，使得获取的原始中心投影影像存在一定的变形，需要按照一定的数学模型实现原始影像与纠正后影像之间的几何变换。航测相机获取的数字影像，一般都满足中心投影成像特性，可以采用共线方程模型建立严密对地定位关系，公式表述为：

其中，XS、YS、ZS为摄站点坐标;X、Y、Z为地面点坐标;为3 个外方位元素组成的9个方向余弦，外方位元素由集成在相机上的POS设备直接测定。

2)CUDA并行程序设计原则

计算统一设备架构(compute unified device architecture，CUDA)是NIVDIA公司推出的一款内嵌在GPU上的通用计算平台，支持使用标准C语言开发GPU代码，这些代码即可运行在CPU上，也可运行在GPU上。CPU负责派生出在GPU运行的多线程任务，即CUDA的内核程序;GPU则通过内部调度器将内核程序映射到GPU的流处理器簇(stream multiprocessors，SM)上，SM数目越多，程序运行速度提升比例就会越高。每个SM调度多个线程块，在每一个线程块内部又包含着许多线程。这些线程以批处理方式运行，称为线程束，是SM中线程执行的最小单元。基于CUDA层次化调度机制，进行GPU并行程序设计时需遵循以下三个原则：(1)要解决的任务可分解为多个可独立执行的子任务;(2)执行每个子任务的时间消耗不大于解决母任务的时间消耗;(3)每个子任务的解决可独立执行，且不需要复杂的流程控制。而面阵影像几何纠正中的逐像素几何变换及灰度重采样操作计算密集度高，无需进行流程控制，符合上述原则，故可将其映射到CUDA内核函数中。

因此，面阵影像GPU几何纠正原理为将CPU视为程序的流程控制端和串行代码执行端，负责面阵影像数据读入、旋转矩阵计算及纠正影像地理范围概算等任务;而将地面点坐标计算、像点坐标计算和灰度重采样三个环节映射到CUDA内核程序中，发挥GPU并行计算优势，高效完成面阵影像快速几何纠正任务。

图1 GPU几何纠正算法流程图

2.2 GPU几何纠正方法

使用共线方程建立无控严密对地定位关系的方法有两种，分别为正解法和反解法。正解法纠正是由原始影像的像素点位置利用共线方程直接解算出纠正后影像像素位置，并将原始影像上该像点的灰度值填充到纠正后影像的相应位置。反解法则是由纠正后的影像像素位置出发利用共线方程解算其在原始影像中的像素位置，并将原始影像的像素位置的灰度值填充到纠正后影像的相应位置。

在CPU平台下，程序是以串行方式执行，体现在影像几何纠正上就是像素的处理是顺序执行的，不同的采样方法之间的差异仅仅是计算量的问题。在GPU平台下，程序是以线程并发的方式执行，并行性能发挥的关键在于应尽量保证程序对每个像素的处理是独立进行的，且流程控制部分不宜过多。正解法纠正对每个像素的处理都需要迭代判断其高程值是否满足限差，计算量大，流程控制复杂，存储器寻址操作复杂度高，GPU并行计算优势无法充分发挥。而反解法纠正有效地避免了线程访问冲突问题及数据不一致性错误，求解过程明确，计算复杂度低，GPU并行易于实现。故本文设计的GPU几何纠正算法使用反解法，校正流程如图1所示。

详细步骤如下：

1)POS数据读入与旋转矩阵计算

在CPU端将待纠正影像和对象的POS数据读入到GPU的全局内存中。根据读入的外方位角元素计算影像旋转矩阵，计算公式为：

式中，φ、ω、κ分别为航向倾角、旁向倾角和旋偏角。

2)计算地面点坐标

根据影像尺寸及旋转矩阵，利用共线方程计算影像地理坐标的四周范围。假定得到的影像左下角地面点坐标为(Xmin，Ymin)，纠正后影像分辨率为m，则影像上位置为(i，j)的像素点对应的地面坐标为：

式中，blockDim、blockIdx、thread Idx均为 CUDA内建变量，数据类型均为dim3类型，分别标识线程块的维度、线程块索引及线程索引。其中，x方向和y方向分别对应面阵影像的宽度和高度。

根据设定的线程索引方式及规定的CUDA并行程序设计原则，可知面阵影像几何纠正算法需要的线程块和线程网格配置维度均为二维，即由二维线程块组成线程网格。在影像处理领域，往往以2的幂次配置每个线程块中的线程数目，其组织方式可采用以下几种方案：2 ×2、4 ×4、8 ×8、16×16、32×32。由于NIVDIA规定线程束的大小为32，前两种方案中线程束均小于该值，造成计算资源的大量浪费。8×8的配置方案，每个SM内的线程数为512，小于每个SM可处理768个线程的约束;32×32的配置方案调度每个线程块中的线程

3)CUDA线程块及线程网格配置

待纠正面阵影像仅存在x方向和y方向两个索引维度，则其在CUDA中的线程索引可以表示为：数目为1024，超出一个线程块只能包含512个线程的限制。而16×16的线程块组织方案，SM每次调度3个线程块，可充分利用SM的处理性能。故本文选择16×16的线程块调度方案，若纠正后的影像宽度为W，高度为H，则其对应的线程网格配置方案为：

4)计算像点坐标

当测量参数无误差时，依据共线方程能够得到无几何失真的地面影像。而在实际应用中，由于相机畸变的存在，使得纠正后的影像仍存在残余误差。因此，建立无控制严密对地定位关系时，首先校正相机畸变。

航测相机畸变主要由相机镜头径向畸变(k1，k2，k3)、切向畸变(p1，p2)、面阵 CCD非正方形比例因子b1、CCD阵列排列非正交性畸变系数b2等因素引起。像点坐标畸变差计算公式为：

式中，δx、δy为航测相机畸变导致的像点坐标偏移量;dx0、dy0、df为像主点及相机焦距改正量，由相机文件提供;dx、dy为像主点坐标相对于面阵影像中心点(0，0) 点的偏移量;Δx=x-x0、Δy=y-y0，分别为量测像点坐标与像主点坐标在x方向和y方向的坐标差值;。

将公式(4)代入到公式(1)，可推导无几何失真的无控严密对地定位关系表达式：

其中：

由于面阵影像的屏幕显示是以行列号表示像素位置，坐标系定义为以左上角为原点，x方向取向右为正，y轴方向取向下为正。因此，还需要建立像平面坐标到面阵影像坐标系之间的坐标变换关系，像平面坐标(x，y)到屏幕坐标(I，J)的转换公式为：

其中，Ic、Jc为像主点所在的行、列号;μ为像元尺寸。

5)灰度重采样

采用反解法计算得到的像点坐标不可能为整数，因此，存在灰度重采样问题。在CUDA程序执行过程中，每个线程块内所有线程共同拥有一块共享存储器，所有线程均可访问全局存储器。共享存储器访存速度优于全局存储器。在灰度重采样时，每个线程都需要从GPU的存储器中读取这组参数，如果使用全局存储器，则会存在严重的访存延迟。因此，本文将坐标变换参数从全局存储器读取到访存速度相对较高的共享存储器中，提升影像灰度重采样效率。

3 实验及结果分析

采用Microsoft Visual Studio 2010开发环境分别开发了在CPU平台和GPU平台上运行的C++程序，其中CPU程序设备环境为Intel®CoreTMi7-5950HQ CPU 2.9GHZ，内存为 32GB，64 位操作系统;GPU设备环境为GeForce GTX 980M，CUDA版本为CUDA7.0。

实验数据为某新型大面阵数字航测相机在中国西部某地拍摄的符合中心投影特性的面阵影像，测区处于平原与山地交界地带，地形起伏较大。单幅影像原始大小为6024×8008，影像位数为8位，数据格式为TIFF;相机主距45mm，像素大小为12μm，航高1500m，影像分辨率为0.19m。相机上集成有POSAV510系统，记录相机摄影时的位置和姿态参数;经POSPAC后处理软件处理后，可为面阵航测影像提供精确的摄站点坐标和相机姿态参数。相机参数文件包含相机畸变系数、像中心点的像素坐标和相机焦距等信息。在每幅影像上均匀分布有32个地面控制点或空三加密点，坐标系为2000国家大地坐标系。

3.1 POS数据定位精度验证

机载POS系统直接对地定位精度可由严密空三解算的地面控制点进行检验，检验公式为：

式中，ica、jca为本文推导公式计算控制点行列号;ico、jco为经严密空三解算得到控制点行列号。

使用本文方法及公式(7)，选取5幅实验影像，利用已有控制数据，验证经检校后的POS数据直接对地定位精度，结果见表1。

表1 1752号影像像点坐标计算精度(单位：像素)

由表1可知，经严格检校的POS数据无控定位精度优于0.25个像素。

3.2 无控制面阵影像几何纠正试验

1)纠正精度分析

在应急测绘应用场景下，无法在短时间内获取大量高精度的控制点。因此，以1752号影像(图2)为例，取测区平均高程开展无控面阵影像几何纠正试验。采用双线性插值法进行灰度重采样，几何纠正结果如图3所示。

图2 原始影像

图3 CPU/GPU纠正后影像

对图3中的32个控制点进行量测，其误差分布情况见表2。

由表2可知，32个检查点的误差最大值为7.79m，最小值为0.03m。可得出结论：使用平均高程进行无控几何纠正精度可达到米级。

表3 相机畸变对纠正精度的影响(单位：m)

表3对比了相机畸变对面阵影像几何纠正精度的影响。结果显示：加入畸变纠正后几何纠正精度分别提高了5.23倍(X方向)和3.92倍(Y方向)。说明航测相机畸变是面阵影像几何纠正的关键步骤，不可忽略。

对比表2和表3中校正相机畸变后的结果可知，所有检查点的误差值均在3倍中误差以内。

2)纠正效率分析

使用加速比Sr对本文提出快速几何纠正算法效率进行评价，定义如下：

式中，Tcpu为CPU串行算法进行面阵影像几何纠正所需的时间;Tgpu为GPU并行算法进行面阵影像几何纠正的所需的时间。

利用最邻近像元法和双线性插值法进行灰度重采样，随机选取5幅影像，其在CPU环境和GPU环境下的运行时间及加速比，见表4和表5。

表4 最邻近像元法所需时间(单位：m s)

表5 双线性插值法所需时间(单位：m s)

从上表可以看出：由于采用了GPU并行算法，面阵影像几何纠正时间大幅缩短。较最邻近像元法，双线性插值方法计算量更大，所需时间更长，但两种插值方法所需时间均小于400ms。较CPU串行程序而言，本文设计的GPU算法加速比在20倍以上。

4 结 论

本文从数字航测相机数据快速处理需求出发，提出了POS数据支持下无控面阵航测影像几何纠正算法，并进行了GPU编程实现。该算法使用经检校后的POS数据实现面阵影像的直接对地定位，定位精度优于0.25个像素;应急条件下，使用平均高程进行无控几何纠正精度可达到米级;较CPU串行纠正方法，GPU并行几何纠正效率提升20倍以上。本文方法同样适用于其他满足中心投影构像条件的遥感影像，具有一定的通用性。

本文研究主要实现了单景航摄影像的快速纠正，还有以下任务有待深入研究：其一，本文方法未使用DEM数据，对于如何开发加入DEM数据的GPU几何纠正算法，进一步提高影像几何纠正精度还需继续研究;其二，在本文研究的基础上，开发处理大数据量的影像快速拼接算法，提升获取大范围区域性影像的时效性，也是未来亟待解决的问题之一。

[1]张祖勋.航空数码相机及其有关问题[J].测绘工程，2004，13(4)：1-5.

[2]苏统华，马培军，刘曙等.CUDA专家手册[M].北京：机械工业出版社，2014.

[3]方留杨，王密，李德仁.CPU和GPU协同处理的光学卫星遥感影像正射校正方法[J].测绘学报，2013，42(5)：668-675.

[4]侯毅，沈彦男，王睿索等.基于GPU的数字影像的正射纠正技术的研究[J].现代测绘，2009，32(3)：10-11.

[5]杨靖宇，张永生，李正国等.遥感影像正射纠正的GPU-CPU协同处理研究[J].武汉大学学报·信息科学版，2011，36(9)：1043-1046.

[6]徐秋辉.无控制点的无人机遥感影像几何校正与拼接方法研究[D].南京：南京大学，2013.

[7]张浩，张兵，卫征.基于POS的高光谱影像快速几何纠正[J].测绘通报，2009(1)：14-17.

[8]程红，王志强，张耀宇.航空影像几何校正方法的研究[J].东北师大学报(自然科学)，2009，41(3)：50-54.

[9]袁国体.航空遥感影像正射校正关键技术研究[D].焦作：河南理工大学，2011.

[10]李德龙，邙晓宇，刘杨等.基于地形数据库与机载POS的影像正射纠正技术探讨[J].地理信息世界，2012，10(6)：58-61.

[11]周前飞，刘晶红，李刚.面阵CCD航空相机斜视图像几何畸变校正误差分析[J].仪器仪表学报，2014，35(6)：1-8.

[12]王伟，童小华，谢锋等.基于POS数据的OMIS影像几何校正[J].遥感信息，2010(1)：89-93.

[13]周国香.UAV载多面阵数码相机拼接技术的研究[D].青岛：山东科技大学，2009.