基于光流模型的视频卫星影像密集匹配方法及系统

摘要

本发明公开了一种基于光流模型的视频卫星影像密集匹配方法及系统，包括：步骤1，采用SGM匹配法对视频卫星影像进行粗匹配；步骤2，采用PMVS匹配法对视频卫星影像进行精匹配；步骤3，采用光流法精化步骤2的匹配结果本发明结合SGM和PMVS两种密集匹配法，选择SGM法粗匹配，然后进行PMVS算法精匹配的形式，尽可能发挥两者的优势，从精度和效率两方面提高影像数据的处理能力。

说明书

技术领域

本发明涉及影像密集匹配技术，特别是一种基于光流模型的视频卫星影像密集匹配方法及系统。

背景技术

经过30多年发展，我国航天技术取得了巨大进步，已形成资源、气象、海洋、环境、国防系列等构成的对地观测遥感卫星体系。特别是在“高分辨率对地观测系统”国家科技重大专项建设的推动下，通过在平台传感器研制、多星组网、地面数据处理等方面的创新,我国遥感卫星的空间分辨率、时间分辨率、数据质量大幅提升,为我国现代农业、防灾减灾、资源环境、公共安全等重要领域提供了信息服务和决策支持。随着遥感应用的深入，应用需求已从定期的静态普查向实时动态监测方向发展，利用卫星对全球热点区域和目标进行持续监测,获取动态信息已经成为迫切需求。由于视频卫星可获得一定时间范围内目标的时序影像，具备了对运动目标的持续监视能力，视频卫星成像技术已成为遥感卫星发展的一大热点。

影像密集匹配在机器人视觉、视频监控、导弹制导、无人机侦察等许多领域中都得到了广泛应用，是模式识别与计算机视觉中的研究热点。使用视频卫星影像进行密集匹配技术是一个对精度和实时性都要求较高的高科学技术，尽管人们对它的研究已经取得了不错的成绩，但是在很多性能指标上仍然有待提高。在自然场景不断变化时，从三维世界空间投影到二维图像空间会造成目标信息丢失，因此设计一个鲁棒的密集匹配方法是非常困难的。因此，如何提高密集匹配方法的适应性，满足对精度和实时性的要求，对保障视频卫星在动态观测领域的应用效果具有十分重要的意义。

发明内容

本发明的目的是提供一种基于光流模型的视频卫星影像密集匹配方法及系统，本发明可实现高精度和高可靠性的匹配。

本发明思路为：

结合SGM粗匹配策略和PMVS精匹配方法，利用视频卫星影像多帧影像间的几何关系，使用光流模型对匹配过程进行约束。

本发明的技术方案为：

一、基于光流模型的视频卫星影像密集匹配方法，包括：

步骤1，采用SGM匹配法对视频卫星影像进行粗匹配，本步骤进一步包括：

1.1对视频卫星影像逐像素视差计算得到视差影像；

1.2计算各像素的最小视差路径，即各像素的正确匹配，以此作为二维约束匹配；

1.3在子步骤1.2所得二维约束匹配下，对子步骤1.1所得视差影像进行SGM匹配，完成视频卫星影像左右影像同名点的匹配，获得匹配点；

步骤2，采用PMVS匹配法对视频卫星影像进行精匹配，本步骤进一步包括：

2.1对影像划分格网，在各格网的匹配点中选取兴趣值为局部极大值的点作为特征点；

2.2将视频卫星影像中各影像轮流作为参考影像，其他影像中选出主光轴与参考影像间夹角小于60°的影像，记为计算影像，基于特征点将参考影像分别与计算影像匹配，并生成种子patch集合；

2.3判断种子patch的邻域中是否存在与种子patch距离小于一个像素的patch或已存在平均相关系数大于阈值的patch，若有，则不向邻域扩散，对下一种子patch执行本步骤；否则，将种子patch向种子patch的邻域扩散；采用扩散后的新patch对参考影像和计算影像进行匹配；所述的阈值为经验值；

步骤3，采用光流法精化步骤2的匹配结果，本步骤进一步包括：

3.1构建相邻两帧影像间的光流方程：

其中，V_x和V_y分别表示x和y方向上的光流向量；q₁,q₂,…,q_n表示特征点(x,y)邻域内的像素点；I_x(q₁)、I_x(q₂)、……I_x(q_n)表示第t帧影像上像素点q₁,q₂,…,q_n的灰度值在x方向上的偏导数；I_y(q₁)、I_y(q₂)、……I_y(q_n)表示第t帧影像上像素点q₁,q₂,…,q_n的灰度值在y方向上的偏导数；I_T(q₁)、I_T(q₂)、……I_T(q_n)表示第t帧影像上像素点q₁,q₂,…,q_n的灰度值关于时间T的偏导数；

3.2对各影像，基于光流方程分别迭代计算其与相邻影像各匹配点间的光流向量；

3.3计算各光流向量的模长d(Vt)，统计影像上所有光流向量的模长的中误差，设阈值W为三倍中误差，将所有误差大于阈值W的匹配点作为粗差剔除。

步骤1.1具体为：

分别计算p′_i和其左临近点p′_i-1、右临近点p′_i+1的中点的灰度值得到p′_i邻域内的最小灰度值和最大灰度值则视差为max[0,I_L(p_i)-I_max,I_min-I_L(p_i)]；

其中，p_i和p′_i分别表示左右扫描线上的像素，I_L(p_i)表示左扫描线上点p_i处的灰度，I_R(p′_i)表示经过右扫描线上对采样点线性内插得到的p′_i的灰度值。

步骤1.2具体为：

对各像素p，在其周围设置若干方向的路径，分别计算各路径上的视差路径，若干路径的视差路径之和为各像素p的视差路径，以最小的视差路径作为各像素p的正确匹配。

步骤1.3具体为：

1.3a在半分辨率下进行SGM匹配，获得匹配后的视差影像；

1.3b提高分辨率，在当前分辨率下对当前视差影像进行SGM匹配；

1.3c不断提高分辨率，并重复步骤(2)直至分辨率达到原始分辨率。

步骤1还包括剔除误匹配的步骤，具体为：

利用右影像中同名点匹配左影像中同名点，如果两次匹配得到的视差不同，则该同名点视为无效匹配并剔除。

步骤2.2具体为：

2.2a对参考影像上各特征点f，寻找其在计算影像上的候选匹配点f′；

2.2b每一对(f,f′)前方交会获得模型点，计算各模型点与参考影像摄像机中心的距离，按距离由近到远对模型点排序，参考影像摄像机后文简记为参考相机；

2.2c按模型点排序顺序，取模型点轮流初始化匹配面元patch的中心c(p)，patch的一条边与参考相机的x轴平行，patch的单位法向量n(p)为由中心c(p)指向参考相机中心的单位向量；初始化后的patch应符合要求：n(p)和patch所在光线的夹角小于60度，且patch投影到计算影像和参考影像之间的相关系数大于预设的第一相关系数阈值；若不符合要求，该模型点失效，取下一个模型点重新执行本步骤；否则，执行步骤2.2d；第一相关系数阈值为经验值；

2.2d优化patch的中心坐标和单位法向量，使patch投影到各计算影像和参考影像之间的相关系数的平均值最大；若优化后，patch的平均相关系数大于第一相关系数阈值，则认为成功生成patch，将所有影像划分格网，将patch投影到影像，记录patch所在格网的坐标，同时去掉patch所在格网中所有特征点，考虑下一个模型点，重新执行步骤2.2c；否则放弃该patch，考虑下一个模型点重新执行步骤2.2c。

7、如权利要求1所述的基于光流模型的视频卫星影像密集匹配方法，其特征是：

步骤2.3中所述的采用扩散后的新patch对参考影像和计算影像进行匹配，具体为：

扩散后的新patch初始的单位法向量和种子patch的单位法向量相同，新patch的中心为通过种子patch邻域格网中心的光线与种子patch所在平面的交点；

优化新patch的中心坐标和单位法向量，使新patch投影到各计算影像和参考影像之间的相关系数的平均值最大；采用优化后的新patch对计算影像和参考影像进行匹配，若匹配的计算影像数大于预设的数量阈值，则认为成功扩散了一个patch，否则失败，跳过该失效patch，对下一个种子patch执行步骤2.3。

步骤2还包括对匹配结果进行滤波，具体为：

(1)过滤掉同一格网中平均相关系数差异大于0.6的patch；

(2)找出格网中深度最小的patch，过滤掉格网中与该patch的距离大于平均距离、且与该patch的法向量夹角大于90°的patch，这里平均距离指该patch和格网中所有patch的距离的平均值；

(3)将patch和其邻域的patch拟合二次曲面，过滤掉邻域patch残差之和大于残差阈值的patch，残差阈值为经验值；

(4)根据patch在物方的相互距离对格网中patch聚类分组，去掉patch数小于预设数量阈值的patch分组，数量阈值为经验值。

步骤3.2具体为：

对相邻两帧影像，分别在x、y、T方向上求灰度梯度；设置原始光流场为零，采用光流方程循环迭代计算各匹配点的光流向量；当两次迭代中匹配点的光流向量误差小于预设的误差阈值，结束迭代。

二、基于光流模型的视频卫星影像密集匹配系统，包括：

(1)粗匹配模块，用来采用SGM匹配法对视频卫星影像进行粗匹配；

所述的粗匹配模块进一步包括：

视差计算模块，用来对视频卫星影像逐像素视差计算得到视差影像；

二维约束匹配模块，用来计算各像素的最小视差路径，即各像素的正确匹配，以此作为二维约束匹配；

SGM匹配模块，用来在二维约束匹配下，对视差影像进行SGM匹配，完成视频卫星影像左右影像同名点的匹配，获得匹配点；

(2)精匹配模块，用来采用PMVS匹配法对视频卫星影像进行精匹配；

所述的精匹配模块进一步包括：

特征点检测模块，用来对影像划分格网，在各格网的匹配点中选取兴趣值为局部极大值的点作为特征点；

种子patch集合生成模块，用来将视频卫星影像中各影像轮流作为参考影像，其他影像中选出主光轴与参考影像间夹角小于60°的影像，记为计算影像，基于特征点将参考影像分别与计算影像匹配，并生成种子patch集合；

PMVS匹配模块，用来判断种子patch的邻域中是否存在与种子patch距离小于一个像素的patch或已存在平均相关系数大于阈值的patch，若有，则不向邻域扩散，对下一种子patch执行本步骤；否则，将种子patch向种子patch的邻域扩散；采用扩散后的新patch对参考影像和计算影像进行匹配；所述的阈值为经验值；

(3)精化匹配模块，用来采用光流法精化步骤2的匹配结果；

所述的精化匹配模块进一步包括：

光流方程构建模块，用来构建相邻两帧影像间的光流方程：

其中，V_x和V_y分别表示x和y方向上的光流向量；q₁,q₂,…,q_n表示特征点(x,y)邻域内的像素点；I_x(q₁)、I_x(q₂)、……I_x(q_n)表示第t帧影像上像素点q₁,q₂,…,q_n的灰度值在x方向上的偏导数；I_y(q₁)、I_y(q₂)、……I_y(q_n)表示第t帧影像上像素点q₁,q₂,…,q_n的灰度值在y方向上的偏导数；I_T(q₁)、I_T(q₂)、……I_T(q_n)表示第t帧影像上像素点q₁,q₂,…,q_n的灰度值关于时间T的偏导数；

光流向量计算模块，用来对各影像，基于光流方程分别迭代计算其与相邻影像各匹配点间的光流向量；

粗差剔除模块，用来计算各光流向量的模长d(Vt)，统计影像上所有光流向量的模长的中误差，设阈值W为三倍中误差，将所有误差大于阈值W的匹配点作为粗差剔除。

与现有技术相比，本发明具有如下优点和有益效果：

(1)基于光流模型，估计卫星视频影像多帧之间几何关系，对影像匹配结果进行约束，可以极大提高匹配可靠性。

(2)结合SGM和PMVS两种密集匹配算法，选择SGM算法粗匹配、然后进行PMVS算法精匹配的形式，尽可能发挥两者的优势，从精度和效率两方面提高影像数据的处理能力。

附图说明

图1为SGM匹配路径选择示意图，其中，图(a)表示所有最小路径的最终累加结果，图(b)表示每个像素上最小路径选择方向示意图；

图2为PMVS匹配流程图；

图3为光流法多帧约束算法流程图。

具体实施方式

下面结合附图进一步说明本发明的具体实施方式，步骤如下：

步骤1，采用SGM匹配法对视频卫星影像进行粗匹配。

SGM(semi-global matching)是一种产生在计算机视觉领域的影像匹配法，利用互信息描述同名像点间复杂的对应关系，以实现影像匹配。SGM影像匹配的基本思想是：首先，基于互信息执行逐像素视差计算；然后，利用多方向上的一维约束近似二维约束。

本步骤进一步包括：

步骤1.1，对视频卫星影像逐像素视差计算得到视差影像。

视差函数d(p_i,p′_i,I_L,I_R)定义如下：

视差函数表示两像素为同名像点的可能性，其中，p_i和p′_i分别表示左右扫描线上的像素。视差计算时，p′取值范围在p′_i左右1/2个像素。I表示对应像素处的灰度值，左扫描线上点p_i处的灰度定义为I_L；I_R则是经过右扫描线上对采样点线性内插得到的p′_i的灰度值。

因为分段线性函数的极值点一定是断点，因此视差d的计算非常容易。实际视差d的计算过程中，根据其极值点的性质优化如下：

首先，计算p′_i和左临近点p′_i-1的中点的灰度值

式(2)中，I_R(p′_i)和分别表示p′_i和p′_i-1的灰度值。

同样的，计算p′_i和右临近点p′_i+1的中点的灰度值

式(3)中，I_R(p′_i)和分别表示p′_i和p′_i+1的灰度值。

得到y_i左右1/2个像素邻域内的最小灰度值和最大灰度值于是视差d(p_i,p′_i,I_L,I_R)＝max[0,I_L(p_i)-I_max,I_min-I_L(p_i)]。

步骤1.2，用一维约束近似二维约束。

对于一种匹配可能，最小视差路径函数L_r(p,d)定义为：

L_r(p,d)＝C(p,d)+min(L_r(p-r,d),L_r(p-r,d±1)+P₁,minL_r(p_i-r,i)+P₂)(4)

像素p的周围，以45°为间隔设置16个路径。通过16个路径分别计算各路径上的最小视差路径L_r(p,d)，d表示子步骤1.1中计算的视差d(p_i,p′_i,I_L,I_R)，r表示路径方向。

式(4)中：

第一项C(p,d)表示对像素p，视差为d时的不相似性代价；

第二项min(L_r(p-r,d_p),L_r(p-r,d±1)+P₁,minL_r(p_i-r,i)+P₂)，如图1所示，视差路径函数在16个方向进行累加。每条累积路径上，像素p当前视差的累积视差路径与像素p视差以及路径上前一个像素的累积视差路径相关。因此从前一个像素的视差和像素p的视差做差，取为0、1或大于1三种情况中挑出最小者作为像素p的累积量：

(1)如果结果为0，像素p的累积量为前一像素的视差路径L_r(p-r,d)；

(2)如果结果为1，像素p的累积量为前一像素的视差路径L_r(p-r,d±1)和惩罚因子P₁之和；

(3)如果结果大于1，像素p的累积量为前一像素所有视差路径的最小者min L_r(p_i-r,i)和惩罚因子P₂之和，p_i表示像素p的前一像素，i表示前一像素的路径方向。

对于每一像素，最后的视差路径为16个方向的视差路径和。视差路径和最小的路径被认为当前像素的正确匹配。以此来近似二维约束匹配计算。

步骤1.3，视差迭代计算。

本发明采用一种多级匹配策略，即首先从半分辨率情况下开始SGM匹配，将得到的视差影像作为初始值；然后，递归地利用上一次的视差影像再次进行SGM匹配，这一过程中不断提高分辨率直至原始分辨率。从而将每一级匹配范围限制在一定的限度内，通过逐级传递，自适应确定每个像素的匹配视差范围。可在减少内存的需求的同时，加快匹配的速度，并减少匹配的粗差。

步骤1.4，误匹配的剔除。

左右影像同名点匹配完成后，利用右影像中同名点匹配左影像中同名点。如果两次匹配得到的视差不相同，则该同名点视为无效匹配并剔除。

步骤2，采用PMVS匹配法对视频卫星影像进行精匹配。

PMVS是一个多视立体匹配方法，由一组照片以及相机参数，重建出照片中物体或场景的三维结构。PMVS只重建刚性结构，它会自动忽略非刚性结构，如建筑物前的行人。软件输出一组带方向的三维点，估计了每个点的三维坐标和法向量。

PMVS的匹配思路是：先找到稳健的、能可靠匹配的特征点，然后再按照区域增长的思路进行匹配传播。整个流程如附图2，具体过程如下：

步骤2.1，特征点提取。

每张影像上画格网，格网大小为32×32pixels。从步骤1中通过SGM粗匹配得到的匹配点中在每个格网中选η＝4个兴趣值为局部极大值的点，作为特征点。

步骤2.2，特征点匹配，生成种子patch集合。

每张影像轮流作为参考影像R(P)，在其他影像中选出主光轴与R(P)间夹角小于60°的影像I(P)，记为计算影像，将参考影像和计算影像匹配，具体为：

(1)对参考影像上每个特征点f，寻找其在计算影像上的候选匹配点f′，f′组成集合F。要求f′在离同名核线两个像素的范围之内，同名核线即同一核面与左右影像相交形成的两条核线，其中核面指物方点与摄影基线所确定的平面。

(2)每一对(f,f′)前方交会获得模型点，计算各模型点与参考影像摄像机(后文简记为“参考相机”)中心的距离，按距离由近到远对模型点排序。

(3)按模型点排序顺序，取模型点轮流初始化匹配面元patch。匹配面元patch是近似于物体表面的局部正切平面，包括中心c(p)、法向量n(p)、参考影像R(p)。patch的一条边与参考相机的x轴平行。patch上有μ×μ的格网，μ一般为5或7。

其中，patch的中心c(p)为模型点，单位法向量n(p)为由中心c(p)指向参考相机中心的单位向量。要求patch的法向量和patch所在光线的夹角要小于60度，且patch投影到计算影像和参考影像之间的相关系数要大于第一相关系数阈值。如果不满足要求，该模型点失败，跳过该点考虑下一个模型点。第一相关系数阈值为经验值，本实施例中，第一相关系数阈值取为0.4。

(4)优化patch的中心坐标和法向量，具体为：细微调整中心坐标和法向量，使平均相关系数最大化，即最大化patch的参考影像与计算影像之间的相关系数的平均值。该优化过程需要把patch的中心固定在参考影像的光线上，优化自由度为3。

如果优化后，patch的平均相关系数大于第一相关系数阈值，则认为成功生成了patch。在所有影像上打上大小为2×2pixels的格网cells Ci(x,y)，把patch投影到影像上，记录patch所在的格网坐标，格网坐标即格网在像片坐标系中坐标。每个格网有两个存储patch的集合Qi(x,y)和Qi*(x,y)，分别存储投影到空间参考坐标系V(p)和像方坐标系V*(p)的patch。同时去掉patch所在格网中所有特征点，接下来考虑下一个模型点。否则放弃该patch，考虑下一个模型点。

待所有影像的patch计算过程完成后，将这些patch放入种子patch集合中。

步骤2.3，匹配传播，由种子patch向种子patch所在格网的邻域扩散，若邻域中已存在与种子patch距离小于一个像素的patch或邻域中已存在一个平均相关系数大于第二相关系数阈值的patch，则不向该邻域扩散。第二相关系数阈值为经验值，本实施例中，第二相关系数阈值取为0.6。

扩散后的新patch初始的单位法向量和种子patch的单位法向量相同，新patch的中心为通过邻域格网中心的光线与种子patch所在平面的交点。优化新patch的中心坐标和单位法向量，并采用优化后的新patch对计算影像和参考影像进行匹配，若匹配的计算影像数大于预设的数量阈值，则认为成功扩散了一个patch，否则失败，跳过该失效patch对下一种子patch进行邻域扩散。

重复步骤2.3，直到无法扩散。

步骤2.4，匹配结果滤波。

匹配结果滤波包括以下四种情况：

(1)过滤掉同一格网中平均相关系数差异大于0.6的patch；

(2)找出格网中深度最小的patch，过滤掉格网中与该patch的距离大于平均距离、且与该patch的法向量夹角大于90°的patch，这里平均距离指该patch和格网中所有patch的距离的平均值。深度指物点与摄像机中心间的距离投影到主光轴方向上的长度。

(3)将patch和其八邻域中的neighbors patches拟合二次曲面，过滤掉neighborspatches残差之和大于残差阈值的patch，残差阈值为经验值。

(4)过滤掉数量较小的patch分组，具体为：根据patch在物方的相互距离对格网中patch聚类分组，去掉patch数小于20的patch分组。

步骤3，基于步骤2得到的匹配点，按照光流模型分别计算匹配点的光流向量，重复迭代然后得到所有匹配点的光流向量并统计其均值与中误差，根据结果将超出阈值的匹配点作为误匹配点剔除，得到最后的匹配点。

参考附图3，具体步骤如下：

步骤3.1，使用稀疏光流法估计步骤2获得的特征点的运动状态估计参数。光流法是一种估计图像序列中像素点在连续帧中运动情况的方法，稀疏光流法只需要处理图像中的某些像素，并不计算图像所有像素的光流估计，从而计算量要相对小。

假设相邻两帧影像之间仅存在微小的变化量，从而建立如下光流方程：

式(5)中：

V_x和V_y分别表示x和y方向上的光流向量，x方向为影像像素列方向，y方向为影像像素行方向；

q₁,q₂,…,q_n表示特征点(x,y)邻域内的像素点；

I_x(q₁)、I_x(q₂)、……I_x(q_n)表示第t帧影像上像素点q₁,q₂,…,q_n的灰度值在x方向上的偏导数；

I_y(q₁)、I_y(q₂)、……I_y(q_n)表示第t帧影像上像素点q₁,q₂,…,q_n的灰度值在y方向上的偏导数；

I_T(q₁)、I_T(q₂)、……I_T(q_n)表示第t帧影像上像素点q₁,q₂,…,q_n的灰度值关于时间T的偏导数。

令v＝(v_x,v_y)，v_x和v_y表示V_x和V_y中的元素，可通过最小二乘法获得近似解：

式(6)中，w_i是点q_i的权值，权值的初始值设为1。

步骤3.2，对各影像，分别迭代计算其与相邻影像各匹配点间的光流向量。

本步骤具体为：

读取多帧相邻原始影像，分别在x、y、T方向上求灰度梯度；设置原始光流场为零，采用光流方程循环迭代计算各匹配点的光流向量；当两次迭代中匹配点的光流向量误差小于预设的误差阈值，结束迭代。

步骤3.3，计算各光流向量的模长d(Vt)，并统计影像上所有光流向量的模长的均值和中误差，设阈值W为三倍中误差，将所有误差大于阈值W的匹配点作为粗差剔除。

最终所得的匹配点即为本方法获得的结果。