python-窗口代价计算视差

一、窗口代价计算视差的的基本原理

二、实现NCC视差匹配方法

三、图像集

四、总结

 

一、窗口代价计算视差的的基本原理

• 立体匹配算法的原理：就是找出两张图像的对应关系，根据三角测量原理，得到视差图；在获得了视差信息后，根据投影模型很容易地可以得到原始图像的深度信息和三维信息。
• 立体匹配算法的分类：全局匹配算法、局部匹配算法
• 局部立体匹配算法：又称基于窗口的方法或基于支持区域的方法，对图像集中每张图像的每个像素计算一个适合大小、形状和权重的窗口。然后对这个窗口内的视差值进行加权平均。理想的支持窗口应该完全覆盖弱纹理区域，并在窗口内深度连续。与全局立体匹配算法相似，通过优化一个代价函数的方法计算最佳视差。但是，在局部立体匹配算法的能量函数中，只有基于局部区域的约束数据项，没有平滑项。局部匹配算法仅利用某一点邻域的灰度、颜色、梯度等信息进行计算匹配代价，计算复杂度较低，大多实时的立体匹配算法都属于局部立体匹配的范畴，但局部立体匹配算法对低纹理区域、重复纹理区域、视差不连续和遮挡区域匹配效果不理想。

• 全局立体匹配算法： 主要是采用了全局的优化理论方法估计视差，建立一个全局能量函数，其包含一个数据项和平滑项，通过最小化全局能量函数得到最优的视差值。其中，图割(Graph cuts， GC）、置信传播（Belief Propagation，BP）、动态规划（Dynamic Programming,DP）,粒子群算法（Particle Swarm Optimization,PSO）、遗传算法（Genetic Algorithm，GA）等优化算法都是常用的求解能量最小化的方法。

         全局立体匹配算法一般用到的函数：

       其中数据项描述了匹配程度，平滑项体现了定义场景的约束，C是匹配代价，P是不同两像素p和q视差的函数，一般称之为惩罚项（penalty），当p点和q点视差不相等时，P>0，且与两者差值越大，P值越大。当p和q视差相等时，P=0。由于全局匹配算法在数学上是一个能量函数的优化问题，因此可以找到最优解。这个问题被证明在二维空间是NP困难的。因此，即使全局算法具有准确性较高的优点，其计算速度确非常慢，在实时性要求高的场合不适合使用全局立体匹配算法。

      比较两种算法之间的差异：局部算法运算速度快，但是弱纹理、遮挡和视差不连续区域错误匹配率高。全局算法精度高，但计算复杂度也相应的较高，耗时长、难以达到实时性。

      立体匹配的步骤：

①匹配代价计算（Cost Computation）：计算参考图像上每个像素点IR(P)，以所有视差可能性去匹配目标图像上对应点IT(pd)的代价值，因此计算得到的代价值可以存储在一个h*w*d(MAX)的三维数组中，通常称这个三维数组为视差空间图

②代价聚合（Cost Aggregation）：  通常全局算法不需要代价聚合，而局部算法需要通过求和、求均值或其他方法对一个支持窗口内的匹配代价进行聚合而得到参考图像上一点p在视差d处的累积代价CA(p,d)，这一过程称为代价聚合。通过匹配代价聚合，可以降低异常点的影响，提高信噪比（SNR,Signal Noise Ratio）进而提高匹配精度。

③视差计算（Disparity Computation）： 局部立体匹配算法的思想，在支持窗口内聚合完匹配代价后，获取视差的过程就比较简单。通常采用‘胜者为王’策略（WTA，Winner Take All）,即在视差搜索范围内选择累积代价最优的点作为对应匹配点，与之对应的视差即为所求的视差。即P点的视差为。

④后处理（Post Process）：  一般的，分别以左右两图为参考图像，完成上述三个步骤后可以得到左右两幅视差图像。但所得的视差图还存在一些问题，如遮挡点视差不准确、噪声点、误匹配点等存在，因此还需要对视差图进行优化，采用进一步执行后处理步骤对视差图进行修正。

二、实现NCC视差匹配方法

2、1原理：对于原始的图像内任意一个像素点(px,py) (p_x,p_y)(p x​ ,p y)构建一个n×n n imes nn×n的邻域作为匹配窗口。然后对于目标相素位置(px+d,py) (p_x+d, p_y)(p x+d,p y)同样构建一个n×n n imes nn×n大小的匹配窗口，对两个窗口进行相似度度量，注意这里的d dd有一个取值范围。对于两幅图像来说，在进行NCC NCCNCC计算之前要对图像处理，也就是将两帧图像校正到水平位置，即光心处于同一水平线上，此时极线是水平的，否则匹配过程只能在倾斜的极线方向上完成，这将消耗更多的计算资源。

NCC计算公式：

 其中NCC(P, d)得到的值得范围将在[-1, 1]之间。

Wp为之前提到的匹配窗口。

I(2, y)为原始图像的像素值。

Ii (PI, Py)为原始窗口内像素的均值。

I2(x + d, y)为原始图像在目标图像上对应点位置在方向上偏移d后的像素值。

2(P:+ d, Py)为目标图像匹配窗口像素均值。

.若NCC= -1,则表示两个匹配窗口完全不相关，相反，若NCC= 1时，表示两个匹配窗口相关程度非常高。匹配流程 

2、2匹配流程

• 采集图像
• 级线矫正
• 特征匹配
• 深度恢复

 

2、3代码实现：

# -*- coding: utf-8 -*-
from PIL import Image
from pylab import *
import cv2
from numpy import *
from numpy.ma import array
from scipy.ndimage import filters
def plane_sweep_ncc(im_l,im_r,start,steps,wid):
    """ 使用归一化的互相关计算视差图像 """
    m,n = im_l.shape
    # 保存不同求和值的数组
    mean_l = zeros((m,n))
    mean_r = zeros((m,n))
    s = zeros((m,n))
    s_l = zeros((m,n))
    s_r = zeros((m,n))
    # 保存深度平面的数组
    dmaps = zeros((m,n,steps))
    # 计算图像块的平均值
    filters.uniform_filter(im_l,wid,mean_l)
    filters.uniform_filter(im_r,wid,mean_r)
    # 归一化图像
    norm_l = im_l - mean_l
    norm_r = im_r - mean_r
    # 尝试不同的视差
    for displ in range(steps):
        # 将左边图像移动到右边，计算加和
        filters.uniform_filter(np.roll(norm_l, -displ - start) * norm_r, wid, s) # 和归一化
        filters.uniform_filter(np.roll(norm_l, -displ - start) * np.roll(norm_l, -displ - start), wid, s_l)
        filters.uniform_filter(norm_r*norm_r,wid,s_r) # 和反归一化
        # 保存 ncc 的分数
        dmaps[:,:,displ] = s / sqrt(s_l * s_r)
        # 为每个像素选取最佳深度
    return np.argmax(dmaps, axis=2)

def plane_sweep_gauss(im_l,im_r,start,steps,wid):
 """ 使用带有高斯加权周边的归一化互相关计算视差图像 """
 m,n = im_l.shape
 # 保存不同加和的数组
 mean_l = zeros((m,n))
 mean_r = zeros((m,n))
 s = zeros((m,n))
 s_l = zeros((m,n))
 s_r = zeros((m,n))
 # 保存深度平面的数组
 dmaps = zeros((m,n,steps))
 # 计算平均值
 filters.gaussian_filter(im_l,wid,0,mean_l)
 filters.gaussian_filter(im_r,wid,0,mean_r)
 # 归一化图像
 norm_l = im_l - mean_l
 norm_r = im_r - mean_r
 # 尝试不同的视差
 for displ in range(steps):
     # 将左边图像移动到右边，计算加和
     filters.gaussian_filter(np.roll(norm_l, -displ - start) * norm_r, wid, 0, s) # 和归一化
     filters.gaussian_filter(np.roll(norm_l, -displ - start) * np.roll(norm_l, -displ - start), wid, 0, s_l)
     filters.gaussian_filter(norm_r*norm_r,wid,0,s_r) # 和反归一化
     # 保存 ncc 的分数
     dmaps[:,:,displ] = s / np.sqrt(s_l * s_r)
 # 为每个像素选取最佳深度
 return np.argmax(dmaps, axis=2)

im_l = array(Image.open(r'jia/1.jpg').convert('L'), 'f')
im_r = array(Image.open(r'jia/2.jpg').convert('L'),'f')
# 开始偏移，并设置步长
steps = 12
start = 4
# ncc 的宽度
wid = 9
res = plane_sweep_ncc(im_l,im_r,start,steps,wid)
import scipy.misc
scipy.misc.imsave('depth.png',res)
show()

2、4结果截图

当wid=2时

运行结果：

当wid=9时

运行结果：

当wid=15时

运行结果：

小结：两张图片是从网上搜集的，图像集下面附有网址，1.jpg为模板，2.jpg为待测图像。我分别把wid设为2、9、15，运行出来的结果如上，当wid越小的时候，结果图的点就越密集，我所观察到的是有轮廓的地方，颜色就越深，而其他地方的颜色就比较浅，有的地方甚至呈现白色。我认为是窗口越小，像素点越密集，图像的清晰度就越高。

三、图像集

                   left                      right

图片链接：http://vision.middlebury.edu/stereo/data/scenes2003/

四、总结

这次试验我理解不是很精，只是跟着我的理解和参考一些博客来写的，图像是我网上收集的，打开我发的链接就可以搜到，代码是网上搜索到的，我运行出来的结果如上，我学到的知识是窗口代价计算视差的原理，也不是很精，就是大概知道一些原理和实现的步骤，从运行结果来看，我认为当窗口值越小时，匹配出的图像的点就越密集，像素点就越密集。有轮廓的地方，颜色就越黑，应该是因为有轮廓的区域或者很黑的的区域与与原模板相同位置的反差越小，反之，颜色越白的区域白表示周围跟模板相同位置反差越大。

参考博客：https://blog.csdn.net/qq_38898129/article/details/93074387