在航空影像拍摄过程中, 由于受到相机内部和外部因素的影响, 会造成影像之间色彩的不一致, 如影像的亮度、色调等。影像之间的色彩不一致直接影响镶嵌后正射影像的质量, 因此需要通过影像间匀光的方法来减弱色彩的差异^[14]。

Wallis滤波器是一种线性滤波函数, 其原理是通过局部影像变换, 使影像内部不同区域或不同影像之间的灰度均值和方差达到近似。Wallis变换可以表示为

f(x, y)=[g(x, y)-mg]csfcsg+(1-c)sf+bmf+(1-b)mg, (1)

式中: g(x, y)为原始影像的灰度值; f(x, y)为经过Wallis变换后影像的灰度值; m_g和s_g分别为原始影像的局部灰度均值与方差; m_f和s_f分别为目标影像的局部灰度均值和方差; c为影像的方差扩展常数; b为影像亮度系数, 且b和c大小均在[0, 1]。式(1)可以线性表示为

f(x, y)=r0g(x, y)+r1, (2)

r0=csfcsg+(1-c)sf, (3)

r1=bmf+(1-b-r0)mg, (4)

式中: 参数r₀和r₁为线性变换的乘性系数和加性系数。

本文基于基准色调的Wallis方法对影像间进行匀色处理(图1)。首先, 选择一幅标准模板影像, 统计其均值和标准差, 作为Wallis滤波处理时灰度均值和标准差变换的目标值; 然后, 统计待匀光影像的均值和标准差, 以目标值为基准, 来调整待匀光影像的灰度分布, 对所有影像进行Wallis匀光处理^[13]。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 影像间Wallis匀光流程Fig.1 Flowchart of Wallis dodging between images

中心投影的成像方式必然会造成高出地面的地物在影像上具有一定的投影变形, 且距离影像中心越远, 投影变形越大。无人机遥感影像有很高的重叠度, 在重叠区域合理地选择像素来源能有效地减小拼接结果的投影变形, 保证DOM精度。Voronoi图是计算几何中基于散点最邻近原则的空间划分方法。如果以每张影像像底点组成散点集生成Voronoi图对平面进行划分, 那么理论上拼接后DOM中每一点的投影变形最小。

传统镶嵌模式是先将测区影像纠正成正射影像, 然后通过边缘提取获得每张正射影像的成像区域, 最后获得成像区域中心位置生成测区Voronoi图网络。而本文从测区影像定向结果出发, 利用公式将影像的像主点投影到地面得到像底点(如图2所示), 然后以测区所有镶嵌点组成散点集, 利用生长法生成以测区所有影像像底点为剖分中心的Voronoi图镶嵌线网络^[16]。Voronoi图镶嵌线网络中以影像像底点为中心的多边形即为该影像对应的镶嵌后正射影像中的有效区域多边形, 公式为

X=Xs-(Zs-Z)a1x+a2y-a3fc1x+c2y-c3fY=Ys-(Zs-Z)b1x+b2y-b3fc1x+c2y-c3f, (5)

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 像底点的获取Fig.2 Photographic nadir

式中: (X, Y, Z)为地面点的物方坐标; (x, y, z)为对应地面点在像方坐标系中的坐标; (X_s, Y_s, Z_s)为影像外方位元素中的线元素; f为相机的焦距; a₁, a₂, a₃, b₁, b₂, b₃, c₁, c₂, c₃是根据影像外方位元素中3个角元素计算得到的旋转矩阵中的9个参数。

在CUDA架构下, CPU被视为主机, GPU被视为供主机调配的设备, 它们都有自己的存储地址空间并能互相通信协同, CPU负责逻辑运算和少量的串行计算, 而GPU则负责大量高度并行化的运算。CUDA将程序单次编译后, 将计算任务映射为大量可以并行执行的线程, 在GPU段Kernel函数是由线程格网单元进行执行。如图3所示, 每个线程格网(grid)由多个线程块(block)组成, 各个线程块并行执行, 且不能相互通信; 而每个线程块又由多个线程(thread)组成, 线程之间并行执行, 同属于相同线程块的线程之间可以相互通信。在CUDA架构下, GPU不仅具有良好的线性结构, 还提供了6种容量、特点和读写效率都不一样的存储器模型: 寄存器、局部存储器、共享存储器、全局存储器、常数存储器和纹理存储器, 使GPU具有了多层次性存储的特点。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 线程组织结构Fig.3 Thread organization pattern

GPU-CPU协同处理就是将CPU看作控制中心负责数据的准备、传输、GPU的调度和逻辑性强的串行计算, 而GPU则负责接收CPU分配的耗时巨大、高度并行的大规模运算, 并将结果反馈给CPU。

航空影像的正射纠正是利用影像的内外方位元素和数字高程模型(digital elevation model, DEM), 基于构像方程将中心投影的原始影像纠正为正射影像的过程。传统正射纠正是将全幅影像按照反解法数字微分纠正基于CPU串行逐像素进行纠正。基于CPU串行的全幅纠正不仅会造成巨大的数据冗余和时间开销, 还会造成内存资源的浪费。为此, 本文根据已生成的镶嵌线网络, 确定每张影像在镶嵌后影像中的有效部分, 只对该部分进行正射纠正, 并将纠正后的结果直接赋给结果影像, 同时将正射纠正计算过程在GPU上运行进一步提高效率。GPU-CPU协同影像有效部分的纠正过程如图4所示。

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 GPU-CPU协同正射纠正流程Fig.4 Flowchart of GPU-CPU cooperation orthographic rectification

1)由Voronoi图镶嵌线网络得到当前影像的有效区域范围(该影像对应的Voronoi多边形), 设Rec为该多边形的最小外接矩形。

2)根据Rec计算该区域覆盖的DEM范围, 再根据公式(6)将Rec投影到原始影像上得到其对应的纹理范围, 将原始影像局部纹理、局部DEM和外参数等从内存拷贝到显存, 即

x=-fa1(X-Xs)+b1(Y-Ys)+c1(Z-Zs)a3(X-Xs)+b3(Y-Ys)+c3(Z-Zs)y=-fa2(X-Xs)+b2(Y-Ys)+c2(Z-Zs)a3(X-Xs)+b3(Y-Ys)+c3(Z-Zs)。 (6)

3)根据Rec的长宽与DOM地面空间分辨率的关系设置线程块的大小和数目, 每一个线程对应一个像素。

4)根据线程索引, 每条线程计算出对应像素的地面坐标(X, Y), 判断该像素是否在有效多边形内部, 若在则内插高程Z, 根据公式(6)计算其在原始影像的位置, 内插并进行灰度赋值; 若不在, 则该像素不做任何处理。

5)将GPU计算的结果从显存传到内存, 并写到硬盘上。

将传统的CPU程序, 按照CUDA的语法转换成可在GPU上执行的程序, 并不能很好的发挥GPU的高性能处理能力, 还需要根据GPU硬件和CUDA的优化策略对GPU程序进行性能优化。CUDA的基本思想就是最大程度的线程并行化, 程序在执行时线程块被映射到流多处理器中调度执行。麦克斯韦架构中, 满足以下条件时, 流多处理器的调度性能最优: ①线程块中线程数的大小不大于1 024; ②每个流多处理器中的线程数最大为2 048; ③一个流多处理器, 最多可以同时调用16个线程块。综合以上条件, 当线程块大小为1 024时, 流多处理器的调度性能最优。

CUDA模型中提供的每种存储器的大小和访问速度不同, 合理有效地利用GPU层次性存储的特点将极大提高程序执行效率。本文根据数据的大小和用途作以下优化:

1)对于数据量不大且每个像素都要做只读访问的影像定向参数、DEM文件参数, 将其用常数存储器存储可以减少大量访问造成的时间开销。

2)对于频繁用于灰度内插的原始影像纹理信息和用于高程内插的DEM栅格数据, 将其存储于纹理存储器之上, 可以有效利用纹理存储器的缓存加速和局部数据快速访问的特点。

3)对于同一线程块内所有线程都要频繁访问的数据, 可以存储在该线程块速度极快的共享存储器上。

4)寄存器虽然有很高的访问速度, 但是每个线程的寄存器数量非常有限, 如果线程中定义的变量过多, 过多的变量将会存储于访问效率低下的全局存储器, 因此在每个线程中尽可能少的定义变量也有助于提高程序的执行效率。

为解决传统正射影像镶嵌方法自动化程度低、耗时长、数据冗余大和硬件内存开销大的问题, 本文提出了一种GPU-CPU协同无人机影像快速镶嵌的方法, 该方法以测区所有影像像底点为散点集, 自动生成测区Voronoi图镶嵌线网络; 从匀光后影像出发将正射纠正的过程嵌入影像镶嵌中, 只对每张影像有效贡献区域进行正射纠正直接赋值给镶嵌后的DOM; 同时将GPU-CPU协同处理技术运用于影像正射纠正和镶嵌过程, 快速获取测区全景DOM影像图。GPU-CPU协同无人机影像快速镶嵌流程如图5所示。

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 GPU-CPU协同无人机影像快速镶嵌流程Fig.5 Flowchart of GPU-CPU cooperation rapidly mosaicking of UAV images

图6给出了影像间匀光的结果, 可以看出在匀光前目标影像与模板影像之间存在着明显的色彩差异, 模板影像色彩均匀, 整体色调偏浅绿, 而目标影像整体色调偏暗, 经过Wallis匀光后, 结果影像的整体色调和模板影像基本趋于一致。

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 Wallis匀光处理结果Fig.6 Results of Wallis dodging

图7给出了用本文方法自动生成的镶嵌线网络和测区镶嵌后DOM全景图。其中图7(a)是将测区所有像主点投影到地面得到的像底点分布; 图7(b)则是以测区所有像底点为散点集自动生成的Voronoi图镶嵌线网络, 网络中每一个多边形对应一张影像的有效区域; 图7(c)和7(d)是镶嵌之后测区全景DOM图, 镶嵌后DOM大小为76 104像元× 66 850像元。从图7看出本文方法快速获得的DOM全景图没有拼接漏洞, 实现了无缝镶嵌。

图7

Fig.7

	Figure Option ViewDownloadNew Window
	图7 快速镶嵌结果Fig.7 Results of rapid mosaicking

图8则给出了测区DOM全景图的镶嵌细节。

图8

Fig.8

	Figure Option ViewDownloadNew Window
	图8 DOM镶嵌细节Fig.8 Mosaicking details of DOM

测区内地形和地物丰富, 从图8中可以看出, 基于DEM进行纠正的DOM不仅在测区平坦区域可以实现无缝镶嵌, 在起伏较大的山地区域也没有纹理错位; 由于DEM不包含建筑物等信息, 因此当镶嵌线穿过建筑物时会有几何错位, 本文使用Voronoi图镶嵌线网络每张影像的有效部分的投影变形都在理论上最小, 因此在镶嵌线穿过非高层建筑时, 也能取得较好的镶嵌效果。

从GPU-CPU协同对影像有效区域正射纠正和整个镶嵌过程对本文方法进行效率分析。由于每张影像的有效区域大小不一, 因此耗时也不同。对一张影像的有效区域正射纠正包括硬盘读取、数据传输、纠正计算和写回硬盘4个部分, 对硬盘的读写时间占据了整个过程的很大比重。为了更好地分析GPU-CPU协同计算的效率, 本研究将整个过程分为数据传输及计算、硬盘读写2个部分进行效率分析。

图9(a)为经过选择配置和存储器优化之后GPU-CPU协同和基于CPU对有效区域正射纠正耗时对比。基于CPU平均耗时为3 610 ms, 基于GPU-CPU协同经过选择配置优化后平均耗时为167 ms, 再进行存储器优化之后平均耗时进一步减少到149.8 ms, 较基于CPU串行计算效率平均提升25.8倍。影像的内部读写有3种格式, 按照行存储的BIL格式、按照像素存储的BIP格式和按照波段存储的BSQ存储, 大多影像默认为BIL格式。图9(b)为3种方式存储格式下, 对相邻15张影像有效区域读取和结果写入的耗时情况。图9(c)为基于GPU-CPU协同且运用BSQ存储策略下与传统基于CUP串行且BIL存储策略下对有效区域纠正并存储整个过程的耗时比较。传统方法对15张影像有效区域正射纠正平均耗时为14 389 ms, 而本文方法平均耗时仅为4 958 ms, 平均加速比为3.2倍。

图9

Fig.9

	Figure Option ViewDownloadNew Window
	图9 影像有效区域正射纠正效率分析Fig.9 Efficiency of ortho-rectification of the valid area

利用本文GPU-CPU协同无人机影像处理方法对230张无人机影像进行快速正射纠正处理直接输出测区DOM全景图只需耗时28.4 min。如果采用文献[7]和文献[14]等的方法, 先对测区所有影像进行全幅正射纠正, 再进行镶嵌线网络选择进而进行影像镶嵌的传统方法, 在全区影像全幅正射纠正上不仅耗费巨大的时间, 而且还会造成巨大的数据冗余, 不利于缩短应急测绘的响应时间。