一种用于遥感影像几何校正的在线自动匹配方法
【技术领域】
[0001] 本发明涉及一种实用的遥感影像自动匹配方法，能够利用网络航空、航天影像地 图对遥感影像进行在线自动匹配。可应用于遥感、摄影测量、测绘、图像处理等领域。
【背景技术】
[0002] 影像自动匹配和几何校正是摄影测量与遥感任务中关键的环节，它们是影像融 合、镶嵌、变化检测、地图更新等高级应用的基础。尽管在过去几十年中对影像自动匹配的 研究已经非常多，但是遥感影像的自动匹配仍然十分具有挑战性。实用的自动匹配方法应 该在效率、稳健性、精度等方面均具有好的表现，但由于遥感影像数据量大、场景大、获取条 件多变、几何变形复杂等特点，现有的方法很难兼顾这几方面。此外，参考影像的准备也是 遥感影像自动匹配和几何校正的一个难点，尤其是高分辨率参考影像的获取往往需要很高 的成本。
[0003] 鉴于现有自动匹配方法的局限性，一种用于遥感影像几何校正的在线自动匹配方 法具有重要的实用价值。
【发明内容】
[0004] 本发明的目的在于解决现有技术的不足，提出一种针对遥感影像几何校正的快 速、稳健、精确的在线自动匹配方法，该方法能够处理任意大小的遥感影像，利用网络影像 地图在较短的时间内自动匹配到一定数量的、分布良好的高精度控制点，直接用于遥感影 像的几何精校正。该方法的优势主要在于充分利用遥感影像的先验几何信息，借助免费或 低成本的网络影像地图资源，高效可靠地获取均匀分布的控制点。
[0005] 为解决上述问题，本发明提供了一种用于遥感影像几何校正的在线自动匹配方 法，该方法包括步骤：
[0006] SI.根据需要采集的控制点的个数将待匹配遥感影像均匀划分为若干个区域；
[0007] S2.如果所有影像区域均已处理，则完成整个影像匹配过程；否则，开始处理下一 个影像区域；
[0008] S3.将当前影像区域按256像素X 256像素的大小划分为若干个影像单元；
[0009] S4.如果当前影像区域中的所有影像单元均已处理，则标记该影像区域已完成处 理，并转入步骤S2 ;否则，开始处理下一个影像单元；
[0010] S5:根据当前待匹配影像单元的范围和初始成像模型计算网络影像地图中参考影 像块的大致范围，然后通过网络下载对应的参考影像块并将其重采样为与待匹配影像单元 相近的分辨率；
[0011] S6 :利用SIFT匹配算子对待匹配影像单元和参考影像块进行匹配，并剔除粗差 点，如果该步骤得到4个以上的匹配点，则进行步骤S7,否则转入步骤S4 ;
[0012] S7.从步骤S6得到的匹配点中挑选最佳匹配点对并对其进行最小二乘匹配，修正 SIFT特征点的点位坐标，将该对匹配点添加到匹配结果中。
[0013] 其中，步骤S5进一步包括：
[0014] S5. 1根据待匹配影像单元的分辨率计算网络影像地图的最临近缩放级别；
[0015] S5. 2计算对应参考影像块的宽度和高度；
[0016] S5. 3计算对应参考影像块的中心点的经炜度坐标；
[0017] S5. 4发送静态地图服务请求并下载对应参考影像块；
[0018] S5. 5将下载的参考影像块重采样为与待匹配影像单元相近的分辨率。
[0019] 其中，步骤S6进一步包括：
[0020] S6. 1对待匹配影像单元和参考影像块进行SIFT匹配；
[0021] S6. 2通过尺度约束剔除粗差点；
[0022] S6. 3利用旋转角度约束剔除粗差；
[0023] S6. 4利用RANSAC估计相似变换约束剔除粗差；
[0024] S6. 5利用仿射变换约束剔除粗差点。
【附图说明】
[0025] 图1为依照本发明一种实施方式的在线自动匹配方法流程图；
[0026] 图2为本发明方法中影像区域和影像单元示意图；
[0027] 图3为依照本发明一种实施方式的在线自动匹配方法中获取网络影像地图参考 影像块流程图；
[0028] 图4为依照本发明一种实施方式的在线自动匹配方法中对待匹配影像单元和参 考影像块匹配流程图；
【具体实施方式】
[0029] 本发明提出的在线自动匹配方法，结合附图详细说明如下。
[0030] 如图1所示，依照本发明一种实施方式的自动匹配方法包括步骤：
[0031] SI.根据需要采集的控制点的个数将待匹配遥感影像均匀划分为若干个区域，并 将各影像区域标记为未处理状态，影像区域示意图见图2 ;
[0032] S2.如果所有影像区域均已处理，则完成并结束整个影像匹配过程；否则，开始处 理下一个影像区域；
[0033] S3.将当前影像区域按256像素X256像素的大小划分为若干个影像单元，并将各 影像单元记为未处理状态，影像单元示意图见图2 ;
[0034] S4.如果当前影像区域中的所有影像单元均已处理，则标记该影像区域已完成处 理，并转入步骤S2 ;否则，开始处理下一个影像单元；
[0035] S5.根据当前待匹配影像单元的范围和初始成像模型计算网络影像地图（目前可 用的网络影像地图包括Google卫星影像图、Bing航空影像地图、MapQuest卫星影像地图和 Mapbox卫星影像地图）中参考影像块的大致范围，然后通过网络下载对应的参考影像块并 将其重采样为与待匹配影像单元相近的分辨率；
[0036] S6.利用SIFT匹配算子对待匹配影像单元和参考影像块进行匹配，并剔除粗差 点，如果该步骤得到4个以上的匹配点，则进行步骤S7,否则转入步骤S4 ;
[0037] S7.从步骤S6得到的匹配点中，选取局部对比度最大的一对，并将SIFT匹配得到 的局部几何变换模型作为初值对该对匹配点进行最小二乘匹配，精化SIFT特征点的点位 坐标，最后将该对匹配点添加到匹配结果中，最小二乘匹配对某一点的条件方程用式（1) 表不，
[0038] kjs(xs, ys)+k2-Ir(xr, yr) = 0 (I)
[0039] 其中xs，ysS待匹配影像单元像素坐标，x y ι为参考影像块像素坐标，
[0040] xr= a o+a^.+a^, yr= b o+biX.+b^,
[0041] a。，B1, a2, b。，bp b2为6个几何变换参数，k p 1^2为2个辐射变换参数，
[0042] Is (xs，ys)和Ir (xr，yr)分别为待匹配影像单元和参考影像块的灰度值，
[0043] 对待匹配点周围11像素Xll像素的区域中各点依式（1)建立误差方程，然后用 Levenberg-Marquardt算法进行最优化求解，得到最优几何变换参数，从而确定精化后的匹 配点坐标。
[0044] 将待匹配影像划分为影像单元后，需要根据当前待匹配影像单元的范围和初始成 像模型计算网络影像地图中参考影像块的大致范围。具体地，如图3所示，步骤S5进一步 包括：
[0045] S5. 1根据待匹配影像单元的分辨率计算网络影像地图的最临近缩放级别，缩放级 别根据待匹配影像单元所在处的经炜度坐标和分辨率利用式（2)计算，影像单元所在处的 经度λ和炜度φ通过待匹配影像的初始成像模型计算得到，
(2)
[0047] 式中R^th米为地球半径，其近似值为6378137米，
[0048] GSD为影像单元所在处的分辨率，
[0049] η为缩放级别，
[0050] [·]表示取最邻近整数的运算；
[0051] S5. 2计算对应参考影像块的宽度width和高度height :给定缩放级别η，网络影 像地图影像坐标X，y与经度λ、炜度φ的换算方法如式（3)和式（4)所示，
[0054] 分别利用待匹配影像单元4个顶点的经炜度坐标通过式（3)计算相应的4个图像 点坐标，然后再求这4个点的最小外接矩形，该矩形的宽度和高度即为参考影像块的宽度 width 和高度 height ;
[0055] S5. 3计算对应参考影像块的中心点的经炜度坐标λ rc，φ ，根据S5. 2中得到的 最小外接矩形计算其中心点图像坐标，然后代入式（4)计算对应的经炜度坐标λκ，（K。；
[0056] S5. 4发送静态地图服务请求并下载对应参考影像块，静态地图服务请求以统一资 源定位符（URL)的方式发出，常用网络影像地图的静态地图服务请求格式如下：
[0057] (