一种遥感视差图的快速、鲁棒生成方法与流程

1.本发明涉及遥感测绘技术领域，具体涉及一种遥感视差图的快速、鲁棒生成方法。


背景技术：

2.遥感测绘领域进行大范围数字表面模型dsm的计算是十分重要的工作，dsm计算的好坏直接影像遥感测绘的成败。其中视差图计算又是整个dsm计算的核心，依据视差图和外部参数就能直接计算dsm。如何既快速、又准确地计算视差图是当前遥感测绘领域的一个重要研究方向。传统方法在计算速度、鲁棒性上均有较大缺陷，无法达到目前计算要求。传统方法的缺陷根源在于其采用的盲搜索策略，在没有先验知识的支撑下，使用大搜索空间进行盲搜索以期能找到有效视差。这样做的后果是：(1)引入大量无效计算，严重影响计算速度；(2)无效数据占据大量内存，有效数据周转困难；(3)对存在相似纹理的图像，计算失败概率大幅增加；(4)没有迭代修正机制，结果还需要大量人工处理。


技术实现要素：

3.有鉴于此，本发明提供了一种遥感视差图的快速、鲁棒生成方法，能够有效加快视差图的计算速度，同时解决因相似纹理引起的视差计算异常问题，使遥感测绘人员快速、高质量地生成数字表面模型(dsm)。
4.本发明的遥感视差图的快速、鲁棒生成方法，包括：
5.步骤1，对左、右灰度立体遥感影像进行平行双目校正；
6.步骤2，对左、右图像进行降采样，采样图像按从上到下、分辨率从低到高的方式建立左、右图像金字塔，金字塔底层为原始图像；
7.步骤3，提取左、右金字塔顶层图像的特征点并匹配，构建特征点视差直方图，直方图主峰对应的视差范围作为初始视差搜索空间；
8.步骤4，从上至下依次提取金字塔各层图像，根据上一层确定的本层像素视差搜索空间，进行sgm计算，得到视差图，并将视差图分割为一般视差区域和异常区域；其中，顶层图像的像素视差搜索空间为初始视差搜索空间；
9.根据图像分辨率以及异常区域内异常像素的个数，将异常区域划分为小异常区域和大异常区域；
10.步骤5，对于一般视差区域内的像素，将该像素的视差值加上冗余a后扩充至下一层，作为下一层对应像素的视差搜索空间；对于小异常区域或大异常区域内的像素，首先对像素视差进行修正，将修改后视差加上冗余b或冗余c后，扩充至下一层，作为下一层对应像素的视差搜索空间；
11.步骤6，重复执行步骤4和步骤5，直到获得金字塔底层图像的区域视差图；
12.步骤7，对底层图像的小异常区域视差图进行像素视差修正，对底层图像的大异常区域视差图进行人工干预修正，最终得到遥感视差图。
13.进一步的，所述步骤2中，采用1/2降采样。
14.进一步的，所述特征点为sift、surf或orb特征点。
15.进一步的，使用ransac进行特征点匹配。
16.进一步的，所述步骤4中，采用区域跟踪算法对视差图进行分割，生成一般视差区域、小异常区域和大异常区域。
17.进一步的，所述步骤5中，对于一般视差区域内的像素，采用如下方法确定下一层的视差搜索空间：
18.对于有效视差点pi(x,y)，其视差是disp
pi
，灰度是gray
pi
；在pi(x,y)周围奇数窗口内搜索最小视差mindisp、最大视差maxdisp、最小灰度mingray和最大灰度maxgray；根据下式确定其在下一层对应的视差搜索空间：
[0019][0020]
其中，[reverse_limit[0],reverse_limit[1])即是有效点pi(x,y)的修正视差搜索空间；α是降采样系数，abs为取绝对值。
[0021]
进一步的，所述步骤5中，对于小异常区域或大异常区域内的像素，采用如下方法对像素视差进行修正：
[0022]
对小异常区域或大异常区域中某点o，在其周围8个方向寻找最近有效视差，设寻找到的最近有效视差disp＝{dispi},i∈[1,8]，对应的寻找距离dist＝{disti},i∈[1,8]，则有：
[0023][0024]
则o点的修正视差为
[0025]
进一步的，所述步骤5中，对于小异常区域内的像素，采用如下方法确定下一层的视差搜索空间：
[0026]
对于小异常区域中某点o，对{dispo，disp1，
…
，disp8}按从小到大重新排列，重新排列后dispo的前一个值是dispm，后一个值是dispn，则有：
[0027]
[0028]
[reverse_limit[0],reverse_limit[1])即是小异常区域某点o在下一层对应像素的修正视差搜索空间。
[0029]
进一步的，所述步骤5中，对于大异常区域内的某点像素o，由反距离加权计算的修正视差搜索空间和反相似度加权计算的修正视差搜索空间共同确定下一层的视差搜索空间：
[0030]
其中，反距离加权计算的修正视差搜索空间limd采用如下方式确定：
[0031]
对{dispo，disp1，
…
，disp8}按从小到大重新排列，重新排列后dispo的前一个值是dispm，后一个值是dispn，则有：
[0032][0033]
limd＝[limitd[0],limitd[1])即是大异常区域某点o在下一层对应像素的反距离加权计算的修正视差搜索空间；
[0034]
反相似度加权计算的修正视差搜索空间lims采用如下方式确定：
[0035]
计算当前异常区域外接矩形，以外接矩形大小将邻域划为3
×
3邻域；计算每个区域的hog特征和灰度直方图特征，将两特征向量连接成区域特征向量，定义中心区域与i区域相似度si＝norm_l2(f,fi)，norm_l2是2范数；则有
[0036][0037]
计算各邻域平均视差反相似度加权视差反相似度加权视差将从小到大重新排序，重现排列后disps的前一个值是后一个值是则
[0038][0039]
lims＝[limits[0],limits[1])是反相似度加权计算的修正视差搜索空间；
[0040]
则：
[0041]
(1)如果以[min(dispd,disps)-len,max(dispd,disps) len)为下
一层的视差搜索空间；
[0042]
(2)如果设左交点disp
l
,右交点dispr，以[disp
l-len,dispr len)为下一层的视差搜索空间；
[0043]
(3)如果取limd为下一层的视差搜索空间；反之，取lims为下一层的视差搜索空间；
[0044]
其中，len＝min(lens,lend)。
[0045]
进一步的，所述步骤7中，对于小异常区域内的像素，采用如下方法对像素视差进行修正：
[0046]
对小异常区域中某点o，在其周围8个方向寻找最近有效视差，设寻找到的最近有效视差disp＝{dispi},i∈[1,8]，对应的寻找距离dist＝{disti},i∈[1,8]，则有：
[0047][0048]
则o点的修正视差为
[0049]
有益效果：
[0050]
本发明针对传统方法存在的缺陷，提出了一种新的遥感视差图计算方法。利用不变特征、图像金字塔和修正准则层层约束视差搜索空间，使得每次计算都在一个有效的小范围内进行，最大程度减少无效计算和内存占用，极大提高了计算速度，降低了计算设备要求；同时对由相似纹理引起的视差计算异常问题也得到了有效解决，方法鲁棒性大幅提高。
[0051]
使用者不需要了解每个参数的含义和调整办法，只需设置其中几个参数即可获得较高质量视差图，降低使用门槛，提升用户体验。
附图说明
[0052]
图1为本发明流程图。
[0053]
图2为特征点视差直方图和初始视差范围确定准则。
[0054]
图3为图像金字塔与sgm迭代计算过程。
[0055]
图4为小异常区视差填充计算示意图。
[0056]
图5为大异常区修正视差搜索空间计算示意图。
具体实施方式
[0057]
下面结合附图并举实施例，对本发明进行详细描述。
[0058]
本发明提供了一种遥感视差图的快速、鲁棒生成方法，方法流程图如图1所示，包括：
[0059]
(1)输入左、右灰度立体遥感影像，图像对必须经过平行双目校正，以保证任意匹配点的行坐标相同。
[0060]
(2)使用1/2降采样方法对左、右图像进行降采样，采样图像按从上到下、分辨率从
低到高的方式建立左、右图像金字塔。金字塔底层就是原始图像。
[0061]
(3)对左、右顶层图像计算sift、surf或orb特征点，使用ransac进行特征点匹配。统计特征点视差直方图，将主峰对应的视差范围[low_lim,high_lim)作为每像素初始视差搜索空间。即首次sgm计算时每像素视差d都在low_lim≤d《high_lim间寻找最优值。
[0062]
(4)根据各像素视差搜索空间，进行首次sgm计算，得到视差图，并使用区域跟踪算法对视差图进行分割生成区域视差图。区域视差图中可能会出现异常区域，这里的异常指的是与周围存在明显差异。
[0063]
(5)判断当前sgm计算的是否是底层图像，如果是，对区域视差图中小异常区域进行视差填充，大异常区域需要人工以某种准则干预填充，算法结束，获得最终的视差图。如果不是，分3种情况处理，1)一般视差区域：根据视差加上冗余a后作为下一层修正视差搜索空间；2)小异常区域(根据图像分辨率确定，例如本实施例原始图像5800
×
5100，50以下认为小异常，大于50认为是大异常)：首先对像素视差进行修正，然后将修改后视差加上冗余b后作为下一层的视差搜索空间；3)大异常区域：首先对像素视差进行修正，然后将修改后视差加上冗余c后作为下一层的视差搜索空间。
[0064]
(6)根据视差扩充准则将(5)计算的修正视差搜索空间扩充到下一层，做为下一层sgm计算时的初始视差搜索空间。
[0065]
(7)迭代计算金字塔每一层，直到底层，算法结束。
[0066]
下面具体介绍每个模块的功能设计
[0067]
1、左、右立体遥感影像输入
[0068]
原始的左、右立体遥感影像是没有经过校正的，必须要使用相机标定参数进行平行双目校正，校正后的左、右图像才能作为输入图像。
[0069]
2、特征点约束初始视差搜索空间
[0070]
特征点约束初始视差搜索空间主要是依据特征点的视差直方图。如图2所示，横轴代表视差，纵轴代表特征点数量。本发明选择直方图的主峰为参考，设置一个阈值ε(实验中取0.8)，当水平直线下降到某一位置时此时对应左、右两侧视差就是初始视差搜索空间[low_lim,high_lim)。
[0071]
3、sgm算法
[0072]
使用经典的sgm算法计算视差图，但不同的是sgm要在图像金子塔上迭代计算。如图3所示，sgm会根据每像素自己的视差搜索空间进行计算。顶层图像每像素的视差搜索空间是一样的，之后下层图像像素的视差搜索空间要根据上层图像相关像素的视差计算得到，也就是除顶层外，其余各层每像素的视差搜索空间都是不一样的。
[0073]
sgm计算后生成视差图，使用区域跟踪算法进行视差图分割。设置分割阈值θ，则有v
(pi,pj∈n)dpi-d
pj
≤θ，其中pi、pj是区域n中的两点，d
pi
、d
pj
是pi、pj的视差。分割后的视差图可能会产生异常区域，异常指的是与周围区域视差存在明显差异。这类异常区域需要经过修正才能输出或进入下一层迭代。本发明将异常区域按区域内异常像素的个数，划分为小异常区域和大异常区域。
[0074]
4、小异常区视差填充修正
[0075]
如果当前层是金字塔底层，对小异常区域需进行视差填充以修正无效视差。方法
如图4所示，对异常区域中某点o，在8方向(正上、正下、正左、正右、左上、左下、右上、右下)寻找最近有效视差，设寻找到的最近有效视差disp＝{dispi},i∈[1,8]，对应的寻找距离dist＝{disti},i∈[1,8]，则有：
[0076][0077]
5、大异常区人工干预
[0078]
如果当前层是金字塔底层，对大异常区域需进行人工干预的方法进行视差填充。这是因为经过多层迭代后依然能不能确定有效视差的区域可能存在先天的不足(原始图像区域可能被云层遮挡或被积雪覆盖或其他原因)，这个时候需要使用先验知识(粗略底图)进行人工干预视差填充。
[0079]
6、一般视差区域，利用有效视差计算冗余a，获得修正视差搜索空间
[0080]
如果当前层不是金字塔底层，则需要对所有像素的视差搜索空间进行修正，将得到的修正视差搜索空间进行视差扩充后传入下一层，作为下一层sgm的初始视差搜索空间。
[0081]
对于有效视差点pi(x,y)，其视差是disp
pi
，使用如下方法计算其修正视差搜索空间：
[0082]
(1)在pi(x,y)周围7
×
7窗口内搜索最小视差mindisp和最大视差maxdisp。
[0083]
(2)在pi(x,y)周围7
×
7窗口内搜索最小灰度mingray和最大灰度maxgray。
[0084]
(3)若pi(x,y)的视差是disp
pi
，灰度是gray
pi
，则有：
[0085][0086]
[reverse_limit[0],reverse_limit[1])即是有效点pi(x,y)的修正视差搜索空间。α是修正系数，实验中取值0.5。
[0087]
7、小异常区计算，利用反距离加权和修正规则计算冗余b，获得修正视差搜索空间
[0088]
使用反距离加权修正法计算修正视差搜索空间。如图4所示，对异常区域中某点o，在8方向寻找最近有效视差，设寻找到的最近有效视差disp＝{dispi},i∈[1,8]，对应的寻找距离dist＝{disti},i∈[1,8]，则有：
[0089][0090]
o点预估视差
[0091]
对{dispo，disp1，
…
，disp8}进行从小到大排列，设排列后序列{
…
，dispm，dispo，dispn，
…
}，dispo在序列中的前一个值是dispm，后一个值是dispn，则有：
[0092][0093]
[reverse_limit[0],reverse_limit[1])即是小异常区域某点o的修正视差搜索空间。
[0094]
8、大异常区计算，利用反相似度加权、反距离加权和空间联合修正规则计算冗余c，获得修正视差搜索空间
[0095]
对于大异常区，要联合使用反相似度加权修正、反距离加权修正和空间联合修正规则。反距离加权修正和7是一样的，设由反距离加权计算的修正视差搜索空间是limd＝[limitd[0],limitd[1])，预估视差dispd，修正步长lend。
[0096]
反相似度加权修正计算过程如5所示，中心区域为异常区域。计算当前异常区域外接矩形，以外接矩形大小将邻域划为3
×
3邻域。计算每个区域的hog特征和灰度直方图特征(hog i是i区域hog特征，gray_hist i是i区域灰度直方图特征)，将两特征向量连接成区域特征向量，定义中心区域与i区域相似度si＝norm_l2(f,fi)，norm_l2是2范数。
[0097][0098]
计算各邻域平均视差反相似度加权视差反相似度加权视差将从小到大排序，设排列后序列从小到大排序，设排列后序列设disps在序列中的前一个值是后一个值是令
[0099][0100]
lims＝[limits[0],limits[1])是反相似度加权计算的修正视差搜索空间。令len＝min(lens,lend)，最终修正视差搜索空间按下面条件确定：
[0101]
(1)如果以[min(dispd,disps)-len,max(dispd,disps) len)为修正视差搜索空间。
[0102]
(2)如果设左交点disp
l
,右交点dispr，以[disp
l-len,dispr len)为修正视差搜索空间。
[0103]
(3)如果取limd为修正视差搜索空间；反之，取lims为修正视差搜索空间。
[0104]
9、修正视差搜索空间扩充
[0105]
预估准则：如果金子塔第i层(x,y)处的点的视差为d，则第i 1层(2x,2y),(2x 1,2y),(2x,2y 1),(2x 1,2y 1)处的视差可预估为2d。
[0106]
若第i层像素pi(x,y)的修正视差搜索空间是[limiti[0],limiti[1])，则i 1层像素p
i 1
(2x,2y)，p
i 1
(2x 1,2y)，p
i 1
(2x,2y 1)，p
i 1
(2x 1,2y 1)修正视差搜索空间可设为[2limiti[0]-σ,2limiti[1] σ)，其中σ是松弛量，实验取5。
[0107]
综上，本发明针对传统遥感视差图计算方法存在计算速度慢，相似纹理计算视差异常等情况，提出了不变特征约束下的sgm金字塔结合搜索空间修正策略的方法。首先对图像降采样生成图像金子塔，对金字塔顶层图像进行特征点检测和匹配，根据特征点视差直方图形成初始视差搜索空间；然后迭代进行sgm计算，除底层外，其他各层均以相应策略修正视差搜索空间并扩充，作为下一层sgm的视差搜索空间。最后计算底层时，对视差图小异常区和大异常区分别以不同策略处理，输出最终视差图。不同于传统方法的盲搜索策略，本方法在进行每层计算时视差搜索空间都被限制在一个有效的小范围内进行，大幅减少了无效计算和内存占用，极大提高了计算速度，降低了硬件要求；范围精确化会抑制周围相似纹理引起的计算异常，再通过迭代修正，使得计算鲁棒性和结果质量大幅提高，提升用户体验。
[0108]
1、对异常区的精细化处理方式
[0109]
对计算修正视差搜索空间进行分层化处理，使用三种方法分别处理有效视差、小异常区、大异常区的修正视差搜索空间计算问题。这种分层化处理机制得到了区别一般的更加精细化结果，大幅改善计算性能。
[0110]
2、对小异常区的处理方法
[0111]
使用反距离加权结合近中值排序联合约束的方法计算小异常区每像素的视差搜索范围，能将搜索约束在一个更加有效的小范围内进行。
[0112]
3、对大异常区的处理方法
[0113]
联合使用反相似度加权、反距离加权并集合运算确定最终修正视差搜索空间。反相似性度加权中相似度量化使用了异常区外接矩形一致性划分3
×
3邻域，使用hog-grayhist特征与norm_l2量化区域相似度，生成反相似度权值，计算反相似加权视差，结合近中值排序联合约束方法计算反相似度修正视差搜索空间。两者融合确定最终视差搜索空间。
[0114]
综上所述，以上仅为本发明的较佳实施例而已，并非用于限定本发明的保护范围。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。