立体匹配算法

立体匹配算法#

简介#

在立体匹配中，匹配问题可以看成是寻找两组数据相关程度的过程。立体匹配算法由多种分类

根据算法运行时约束的作用范围：分为局部（local）匹配算法和全局（Global）匹配算法
基于生成的视差图：可分为稠密（Dense）匹配和稀疏（Sparse）匹配

1. 全局匹配算法#

全局/半全局立体匹配算法主要是采用了全局的优化理论方法估计视差，建立一个全局能量函数，其包含一个数据项和平滑项，通过最小化全局能量函数得到最优的视差值。其中，求解能量最小化的优化算法包括：

图割(Graph cuts， GC）
置信传播（Belief Propagation，BP）
动态规划（Dynamic Programming,DP）
粒子群算法（Particle Swarm Optimization,PSO）
遗传算法（Genetic Algorithm，GA）

全局匹配算法一般定义如下能量函数

E (d) = E_{d a t a} (d) + E_{s m o o t h} (d) = \sum_{p \in R} C (p, d) + \sum_{q, p \in R} P (d_{q} - d_{p})

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

考虑到能量优化问题在一维空间的复杂度是多项式级的，因此一些研究试图做一些近似来降低算法的复杂度。例如，半全局算法（SGM）就利用了这一特性将二维问题简化为8到16个一维问题，以实现一种较好的近似。其在各个方向上计算累积代价后，将各方向代价相加得到总代价，这样就模拟了二维的优化问题。SGM是立体匹配逐渐取代激光雷达生成视差图的技术关键，同时也是商业软件中应用最多的立体匹配算法。

2. 局部匹配算法#

局部立体匹配算法又称基于窗口的方法或基于支持区域的方法，算法对参考图像中的每个像素计算一个合适大小、形状和权重的窗口，然后对这个窗口内的视差值进行加权平均。理想的支持窗口应该完全覆盖弱纹理区域，并在窗口内深度连续。与全局立体匹配算法相似，通过优化一个代价函数的方法计算最佳视差。但是，在局部立体匹配算法的能量函数中，只有基于局部区域的约束数据项，没有平滑项。局部匹配算法仅利用某一点邻域的灰度、颜色、梯度等信息进行计算匹配代价，计算复杂度较低，大多实时的立体匹配算法都属于局部立体匹配的范畴，但局部立体匹配算法对低纹理区域、重复纹理区域、视差不连续和遮挡区域匹配效果不理想。

基于区域的局部立体算法是最早开始研究，算法成熟、计算简单、速度快，匹配精度较高; 基本原理：在参考图像中选择一个点，选择该点邻域内一个支持窗口，然后依据一定的相似性判断准则，在待匹配图像中寻找与支持窗口最相似的子窗口，该子窗口所对应的像素点即为对应的匹配点。

固定窗口代价聚合使用固定大小和形状的窗口作为代价聚合的基元，通常是一个矩形，并假设支持窗口内的其它像素点与待匹配点具有相同的视差。固定窗口法精度不高，但易实现、耗时短，在一些对实时性要求极高的场合得到了应用。

基于双边滤波的代价聚合算法仍然使用固定大小和形状的窗口，但窗口内的元素权重不同，权重由目标图像在该窗口内像素与窗口中心的灰度差和距离计算。基于双边滤波的代价聚合算法精度高，但计算复杂，实时性差，算法性能随窗口尺寸指数增加。

基于分割的代价聚合算法的主要思想是：预先将作为参考图像的左图进行分割，对于支持窗口内与窗口中心处于同一分割的像素，对应的权值取1，否则为一个远小于1的正数。但是图像分割是一个非常耗时的操作，同样无法在实时性要求较高的场合使用。

基于十字的代价聚合算法（Cross-based Cost Aggregation,CBCA）的支持窗口形状并不确定，会根据匹配点邻域的灰度值而改变，该方法可以使用GPU并行计算，具有较好的实时性，现广泛应用于各种算法的代价聚合步骤。

3.立体匹配步骤(半全局)#

3.1.匹配代价计算#

计算匹配代价，即计算图像上每个像素点 $L_{p}$ ，以所有视差可能性去匹配目标图像上对应点 $R_{p}$ 的代价值，因此计算得到的代价值可以存储在一个 $h * w * d_{m a x}$ 的三维数组中，通常称这个三维数组为视差空间图（Disparity Space Image，DSI）。匹配代价是立体匹配的基础，设计抗噪声干扰、对光照变化不敏感的匹配代价，能提高立体匹配的精度。匹配代价的设计在全局算法和局部算法中都是研究的重点。

主要方法：

AD(absolute diff)
AD算法是基于单个像素点计算的匹配代价，受光照不均、图像噪声影响较大，但对纹理丰富区域有较好的匹配效果。

C_{A D} (p, q) = | I_{L} (p) - I_{R} (q) |

SAD(sum of absolute diff)
Np、Nq分别表示p、q周围的像素点,m,n 对应patch in L&R

C_{SAD} (p, q) = \sum_{m \in N_{p}, n \in N_{q}} | I_{L} (m) - I_{R} (n) |

Cencus
基于Census变换的匹配代价计算方法是计算左右影像对应的两个像素的Census变换值的汉明（Hamming）距离，即匹配代价为:

C_{Cencus} (p, q) = Haming (T (p), T (q))

其中Cencus变换公式定义如下:

T (p) = \underset{q \in N_{p}}{\otimes} ξ (I (p), I (q))

其中 $p$ 是窗口中心像素， $q$ 是窗口中心像素以外的其他像素， $N p$ 表示中心像素 $p$ 的邻域。 $I (*)$ 表示像素点*处的灰度值， $\otimes$ 为bit的逐位连接运算。 $ξ ()$ 运算则由下面公式定义：

ξ (I (p), I (q)) = {\begin{cases} 0 I (q) \leq I (p) \\ 1 I (q) \geq I (p) \end{cases}

AD和SAD算法都对光照较为敏感，这里的Cencus变换则对图片的明暗变化并不敏感，因为Cencus算法是比较的相对灰度关系，所以即使左右影像亮度不一致，也能得到较好的匹配效果。但是Cencus变换对重复区域的匹配效果不好。

AD-Cencus
前文已经介绍了AD算法和Cencus变换，显而易见，AD-Cencus是将AD和Census结合,这样就能对两种方式起到一个互补作用。Cencus算法对重复纹理的效果不好，而AD算法是基于单像素的，可以在一定程度上缓解Cencus算法对重复纹理处理棘手的问题。但是将两种算法结合存在算法结果尺度不一致的问题，需要进行归一化处理。，AD的结果是亮度差，范围是[0,255]，而Census是比特串对应位值不相同的个数，范围为[0,N]（N等于比特串的位数,跟设置的窗口大小有关系）。因此，需要通过归一化，将两者的结果归一化到相同的范围区间，AD-Census所采用的方法是一个值区间在[0,1]的自然指数函数：

ρ (c, λ) = 1 - e^{- \frac{c}{λ}}

其中c是代价值，λ是控制参数，当c和λ都为正值时，这个函数的值区间在[0,1]的。并且c即代价值越大，函数值越大。因此可以通过该函数将任意代价值归一化到[0,1]的范围。
最终，AD-Census的代价计算公式为：

C (p, q) = ρ (C_{Cencus} (p, q), λ_{cencus}) + ρ (C_{A D} (p, q), λ_{A D})

3.2.代价聚合#

通常全局算法不需要代价聚合，而局部算法需要通过求和、求均值或其他方法对一个支持窗口内的匹配代价进行聚合而得到参考图像上一点 $p$ 在视差d处的累积代价 $C_{A_{(p, d)}}$ ，这一过程称为代价聚合。通过匹配代价聚合，可以降低异常点的影响，提高信噪比，进而提高匹配精度。代价聚合策略通常是局部匹配算法的核心，策略的好坏直接关系到最终视差图的质量。
由于噪声等因素，基于像素的代价计算通常是不明确的，并错误匹配比正确匹配的代价更低。因此，增加了一个附加的约束，通过惩罚相邻差的变化来支持平滑性。

像素代价和平滑项约束可以用能量函数 $E (D)$ 来表示，能量函数 $E (D)$ 依赖于深度图像 $D$ ：

\begin{aligned} E (D) = \sum_{p} (C (p, D_{p}) + \sum_{q \in N_{p}} P_{1} T [| D_{p} - D_{q} | = 1] & + \sum_{q \in N_{p}} P_{2} T [| D_{p} - D_{q} | > 1]) . \end{aligned}

代价聚合公式:
像素p沿着某条路径r的路径代价计算公式

L_{r} (p, d) = C (p, d) + min {\begin{cases} L_{r} (p - r, d) \\ L_{r} (p - r, d - 1) + P_{1} \\ L_{r} (p - r, d + 1) + P_{1} \\ min_{i} L_{r} (p - r, i) + P_{2} \end{cases}} - min_{i} L_{r} (p - r, i)

其中，第一项为匹配代价值 $C$ 属于数据项；第二项是平滑项，表示的是和式2相同的含义，累加到路径代价上的值取不做惩罚、做 $p_{1}$ 惩罚和做 $p_{2}$ 惩罚三种情况中代价最小的值；第三项是为了保证新的路径代价值不超过一定数值上限.

关于公式中第二项中 $m i n_{i} L_{r} (p - r, i)$ 的i值，理解是它应该满足 $| i - d | > 1$ ,即相关像素的视差不是第二项中前三种情况，也即像素的视差>1.

聚合的具体怎么进行的...
我们现在需要计算像素p点的聚合代价值，在计算之前肯定是知道p的各个视差下的初始匹配代价，所有对于p来说，要求p的路径聚合代价，那么根据公式 $L_{r} (p, d)$ ，其第一部分C(p,d)已知。那么对于图像中的任意一个像素，都可以在初始代价矩阵C中找到各个视差下的初始匹配代价，再来看公式4的第二项，不难看出，这一项主要是要寻找前一个像素的路径聚合代价, 上一个像素的路径聚合初始代价=初始的代价。

Copy Highlighter-prismjs prismjs-lines-highlighted language-c
memcpy(cost_aggr_row, cost_init_row, disp_range * sizeof(uint8_t));
memcpy(&cost_last_path[1], cost_aggr_row, disp_range * sizeof(uint8_t));

最后总路径代价

S (p, d) = \sum_{r} L_{r} (p, d)

3.3.视差计算#

局部立体匹配算法的思想，在滑动窗口内聚合完匹配代价后，获取视差的过程就比较简单。通常采用'赢者通吃'策略（WTA，Winner Take All）,即在视差搜索范围内选择累积代价最优的点作为对应匹配点，与之对应的视差即为所求的视差。即 $p$ 点的视差为 $d = \underset{d}{a r g m i n} C_{A_{(p, d)}}$ 。

3.4.后处理#

视差图还存在遮挡点视差不准确、噪声点、误匹配点等问题存在，因此还需要对视差图进行优化，进一步执行后处理步骤对视差图进行修正。常用的方法有插值（Interpolation）、亚像素增强（Subpixel Enhancement）、精细化（Refinement）、图像滤波（Image Filtering）等操作。

优化手段：