高分三号顺轨影像快速匹配并行计算方法

耿 勇 胜，洪 中 华*，童 小 华，刘 世 杰，冯 永 玖，张 云，周 汝 雁，韩 彦 岭，潘 海 燕，王 静，杨 树 瑚，徐 利 军

(1.上海海洋大学信息学院，上海 201306；2.同济大学测绘与地理信息学院，上海 200092)

0 引言

星载SAR影像广泛应用于图像融合与拼接、三维重建、变化检测等方面[1-4]，而SAR影像匹配是这类应用的重要基础。星载SAR影像匹配主要有基于强度的方法和基于特征的方法[5,6]，前者需要计算大量的相关性矩阵，时间复杂度高，后者包含特征提取、特征描述子构建和相似性准则匹配等过程，其速度更快、精度更高。SIFT(Scale Invariant Feature Transform)方法是具有代表性的特征匹配方法之一[7]，具有良好的尺度不变性和抗噪性[8]，已广泛应用于各种遥感影像匹配中[9,10]。然而，SAR影像不同于光学影像，存在相干斑噪声，导致SIFT算法难以提取高质量的特征点。针对这一难题，Dellinger等[11]提出了SAR-SIFT算法，通过图像强度指数加权算子(ROWEA)和梯度影像对数计算以降低相干斑噪声的影响，但由于SAR-SIFT采用CPU串行计算处理SAR卫星影像效率较低，无法满足地震等自然灾害发生后快速响应场景之需。

针对遥感影像匹配效率不高的难题，当前处理方法主要有4种：分块匹配方法、特征尺度降维法、CPU并行计算、GPU并行计算。分块匹配策略的核心是减少影像匹配中非重叠区域的冗余计算，如杨佳宾等[12]通过无人机影像的POS信息提取重叠区域以减少Dense SIFT快速拼接的冗余计算，但该类方法仅采用分块和CPU串行的方式效率仍不高；魏天华[13]采用特征尺度降维和多核DSP方法提高了SAR-SIFT匹配效率，但该方法是通过降低金字塔维数损失匹配精度以提高效率；利用多核CPU并行计算可加倍加快处理效率，在影像匹配中应用更为广泛，如程栋等[14]提出利用OpenMP基于多核CPU并行计算对SAR影像匹配中的互相关模板匹配进行加速(加速比为8)，宋伟等[15]利用OpenMP基于多核CPU并行计算对无人机影像匹配方法进行加速(加速比优于2)，但CPU核数有限导致该类方法效率提升并不明显；GPU核数远大于CPU核数，更适用于影像匹配效率提升，如许川佩等[16]利用OpenCL基于GPU并行加速计算提高SIFT算法图像匹配效率(加速比优于10)，张杰等[17]提出利用CUDA基于GPU对特征提取中的图像金字塔和梯度计算进行加速(加速比为40)，肖汉等[18]利用基于CUDA的GPU并行计算对SIFT算法加速(加速比约为30)。

综合上述研究可以发现，目前针对SAR影像SIFT匹配加速的研究较多，且取得了明显进展，但专门针对SAR-SIFT的并行计算方法未见报道。因此，本文提出一种基于CPU-GPU协同并行计算的SAR-SIFT快速匹配方法，以期提高SAR影像匹配速度。

1 SAR影像快速匹配并行计算方法

本文提出的SAR-SIFT快速匹配方法主要步骤为：首先对SAR影像进行分块预处理并提取主副影像间的重叠区域，然后采用Fork/Join模式的CPU-GPU协同并行计算以提升SAR-SIFT匹配效率。

1.1 SAR影像分块预处理与重叠区域提取

针对星载SAR影像文件大导致匹配效率低问题，本文先对影像进行预处理，求出重叠区域，再进行后续匹配。1)分别利用主副影像有理函数模型参数(RPCs)计算影像4个角点的经纬度坐标，然后通过求交判断该影像是否与其他影像相交；2)将每幅影像依次设为主影像并进行分块，求出与相交副影像的重叠区域，以减少影像分块时不重叠区域的冗余计算；3)将主影像分成N×N块，并求得每一子块4个角点的经纬度坐标，进而求得主副影像块之间的重叠区域；4)采用RPCs反算出主副影像块中重叠区域的行列坐标，而后对主副影像块中的重叠区域进行匹配。

1.2 基于SAR-SIFT算法的连接点提取

SAR-SIFT算法在特征尺度空间和特征描述方面对SIFT进行了改进，提高了SAR影像的匹配精度，但随着数据规模增大，计算效率逐渐降低。SAR-SIFT生成尺度空间响应函数R(x,y,σ)计算如下[11]：

(1)

R(x,y,σ)=Det(CSH(x,y,σ))-t·Trace(CSH(x,y,σ))

(2)

首先，通过尺度空间响应函数R(x,y,σ)在每个尺度空间上寻找3×3范围内的极大值进行特征点提取，其中响应函数值大于0.8。在确定特征点坐标后，通过特征点邻域的点经过高斯加权求得特征点主方向，然后采用GLOH[20]方法构建特征描述子。利用特征描述子初匹配后，采用RANSAC筛选出正确的匹配点对，即连接点。在连接点提取过程中，通过Fork/Join模式对金字塔影像进行分层并行处理，该并行模式在CPU端和GPU端分别采用OpenMP和CUDA技术对SAR-SIFT匹配方法进行加速。

1.3 基于OpenMP和CUDA的SAR-SIFT协同并行计算方法

本文提出的SAR-SIFT并行计算方法主要技术流程如图1所示。

图1 SAR-SIFT协同并行计算流程Fig.1 Framework of collaborative parallel computing for SAR-SIFT

(1)在获取输入影像后，采用GPU构建特征金字塔，以充分发挥GPU并行处理和浮点运算力的优势；同时利用CPU对SAR-SIFT协同并行计算，其遵循Fork/Join并发模式。在构建8层特征金字塔时，由OpenMP编译指令分配8个独立线程处理对应子任务中的响应函数。在每个子任务中，影像从内存传入显存中，由GPU完成SAR-HRRIS响应函数的计算。该计算首先对影像做高斯模糊操作并计算其梯度幅值和梯度方向，然后计算像素点协方差矩阵，最后计算当前层的SAR-HARRIS响应函数。主线程等待8个工作线程处理完毕后一起将该层的SAR-HRRIS响应函数返回CPU内存中，该过程中仅需将影像传入显存，由GPU统一计算特征尺度，以减少显存与内存之间数据传输耗时问题。

(2)完成特征金字塔构建后，利用OpenMP在CPU端并行提取特征点和构建特征描述子。在提取特征点时，由OpenMP编译指令分配8个工作线程用于提取每一层的特征点；在构建特征描述子时，由OpenMP分配N个工作线程用于特征描述子构建；匹配阶段应用knnmatch方法(基于欧氏距离匹配)[21]，采用RANSAC[22]三倍中误差筛除错误的匹配点。这种CPU-GPU协同并行计算方法充分利用了计算机硬件性能，可显著提高SAR影像匹配效率。

2 实验与分析

2.1 实验数据与环境

本文实验数据为2021年5月26日采集的覆盖云南省大理市漾濞县的5景高分三号L1A级全极化条带影像，该影像工作频段为C频段，方位向空间分辨率为8 m，距离向空间分辨率为6～9 m，幅宽为30 km。漾濞县地势北高南低，海拔最高为4 122 m，最低为1 174 m，总面积为1 954 km2，其中山区面积占98.4%。

实验硬件环境：CPU I-9-10900X，显卡 RTX2080，内存 64 GB，显存 11 GB；软件环境：Windows 10 64位操作系统、Visual Studio 2017、Opencv 3.4.11、GDAL 2.3.1、Eigen 3、CUDA Toolkit 10.1。所有实验在同一环境下运行，采用C++编写代码，在x64debug模式下记录程序运行时间。

算法参数设置如下：在构建金字塔时，金字塔层数设为8；在提取特征点时，SAR-HARRIS响应函数的阈值设为0.8；在构建特征描述子时，在提取的所有特征点中保留响应值前4 000个特征点用于构建特征描述子匹配，从而提高匹配速度；在匹配阶段，欧氏距离匹配阈值设为0.8；最后保存每块影像上最优的10个连接点。

2.2 匹配结果分析

采用本文提出的SAR-SIFT快速匹配方法在云南漾濞震区5景高分三号顺轨影像的连接点提取结果如图2所示，其中重叠区域中的紫色十字符号代表相邻两景影像的连接点坐标，在左右两边对应的高分三号匹配影像中，紫线代表连接点连线。由图2可知，影像重叠区域约占每景影像的20%，提取的连接点分布较均匀且全部落入重叠区域内。本文采用RMSE[23]作为匹配精度指标，得到A、B、C、D 4个区域的RMSE分别为0.3076、0.2673、0.2478、0.2569，均为亚像元精度，满足图像拼接需求。

图2 高分三号遥感影像匹配结果Fig.2 Matching results of GF-3 remote sensing images

2.3 匹配效率分析

为验证本文提出的基于重叠区域分块匹配策略和CPU-GPU协同并行计算在SAR-SIFT匹配方法中的有效性，设计了两组对比实验：1)在分块匹配策略中，设计了未分块以及从4×4分块到10×10分块8个对比方案；在CPU-GPU协同并行计算方案中，设计了CPU单线程、CPU多线程、GPU并行、CPU-GPU协同并行4种对比方案，其中GPU并行是指在CPU串行执行每一层金字塔构建时采用GPU并行计算。本文重点关注星载SAR影像匹配效率，采用加速比作为量化参数[24]。实验结果如表1所示，其中耗时包括影像预处理、通讯消耗和特征点提取、特征描述子构建、特征匹配所需时间，加速比1为未分块影像匹配耗时除以不同分块影像匹配最低耗时，加速比2为单线程影像匹配耗时除以多线程影像匹配耗时。2)将本文方法与OpenMP加速[15]、分块匹配[25]、减层降维[13]不同分块策略进行对比，结果如表2所示。

由表1可知：1)影像匹配耗时随着CPU线程数增加而递减，当线程数大于5时耗时减少并不明显，主要原因是5个线程只能同时处理8层金字塔影像中的5层，导致后3层金字塔影像处理时出现2个线程等待，当线程数增至8时，加速比为5.74；GPU并行方案可大幅提升处理效率，加速比达12.67；CPU八线程+GPU协同并行计算在所有计算方案中耗时最少，加速比达14.44。2)同一线程数的大部分情况下，6×6分块方案的加速比最高。3)在所有计算方案中，分块方案比未分块方案快2倍以上。4)在分块方案中，匹配次数随着影像分块数增加而增加，进而增加耗时，而未分块影像只需匹配一次，耗时较少；此外，GPU并行计算时分块影像需多次传入GPU显存，导致通讯时间变长、加速比降低，但分块可减少冗余计算，使加速比下降不明显。5)“6×6分块+CPU八线程+GPU并行计算”为最优方案，其相对于未分块CPU串行计算方案的加速比达45.15。由表2中的加速比可知，本文方法在5景高分三号顺轨影像匹配中效率最高，加速比为8.90，验证了本文方法的有效性。

表1 影像不同分块和不同并行计算方案的实验结果Table 1 Experimental results of different image blocks by using various parallel computing schemes

表2 不同分块、并行策略耗时对比Table 2 Time-consuming comparison of different image blocks and parallel strategies

3 结语

本文提出一种基于重叠区域分块匹配策略和CPU-GPU协同计算的SAR-SIFT特征匹配方法，采用分块匹配策略以减少非重叠区域的冗余计算，并基于Fork/Join模式CPU-GPU协同并行计算以提高SAR-SIFT特征匹配效率。为验证本文方法的有效性，在云南漾濞县5景高分三号影像上开展影像匹配实验，并与当前主流的分块匹配策略、OpenMP多线程、减层降维等方法进行对比，结果表明本文方法效率最高，加速比为8.90，为解决SAR影像匹配效率低问题提供了一种解决方案，为进一步提升大区域SAR影像匹配性能提供了有效策略和方法。