《增强现实：原理、算法与应用》读书笔记（3）特征与匹配

《增强现实：原理、算法与应用》读书笔记（3）特征与匹配

这一部分真的非常非常非常长，但是我个人感觉真的很重要，要不插个目录吧

第三章：实景的三维结构恢复和重建

特征

特征提取的目的是什么不用我多说了，一般图像中的点和边缘是我们最常用的特征，特征点可定义为在两个或两个以上方向颜色迅速变化的像素的集合（Forsyth et al., 2011）。由于常常出现在两条边的汇聚处，又称为角点（Corner point）。

特征最重要的性质是可重复性，指一张图像在进行旋转、尺度（分辨率）、光照、视角变化等操作后仍能在相同的位置重复提取出同一特征。由此引申出特征的尺度不变性、旋转不变形和仿射不变性等。

图像预处理与梯度提取

特征一般存在于图像梯度模较大的地方。由于特征提取过程中经常使用灰度图像，所以本节采用标量函数来表示。图像的梯度是两个正交方向的一阶导数，定义为：
$$▽I\left(x,y\right)={\left(\frac{\partial I\left(x,y\right)}{\partial x},\frac{\partial I\left(x,y\right)}{\partial y}\right)}^{\top }$$

因为图像是由间隔为1的离散值，所以水平和竖直方向上的梯度分量可以表示为：
$$\frac{\partial I\left(x,y\right)}{\partial x}\approx I\left(x+1,y\right)-I\left(x,y\right)$$$$\frac{\partial I\left(x,y\right)}{\partial y}\approx I\left(x,y+1\right)-I\left(x,y\right)$$

这样一来，就比较好理解Sobel算子的由来了。比如，我们计算某点的水平梯度时，应该使用如下算子：$$\left[-1,1\right]$$

但是，想想这个算子算的好像是中间的$0.5$像素处的梯度，不太方便，所以在中间加个0：$$\left[-1,0,1\right]$$

想想这样似乎不能体现出周围像素的影响（噪声、柔化），那就在上下再加一层：
$$\left[\begin{array}{ccc}-1& 0& 1\\ -1& 0& 1\\ -1& 0& 1\end{array}\right]$$

这样似乎太均匀化了，凸显不出中间像素的重要性，那就乘个2，然后我们就得到了Sobel算子：
$$\left[\begin{array}{ccc}-1& 0& 1\\ -2& 0& 2\\ -1& 0& 1\end{array}\right]$$

如果想要加强对中间层的强调，也可以进一步增大中间层的占比，就得到了衍化Sobel算子（大误）：
$$\left[\begin{array}{ccc}-3& 0& 3\\ -10& 0& 10\\ -3& 0& 3\end{array}\right]$$

图像的二阶导数原则上是矩阵，但一般采用拉普拉斯算子获得的标量来简单表示：$${▽}^{2}I\left(x,y\right)=\frac{{\partial }^{2}I\left(x,y\right)}{\partial x^2}+\frac{{\partial }^{2}I\left(x,y\right)}{\partial y^2}$$另一种表达方式是：$${▽}^{2}I\left(x,y\right)=I\left(x+1,y\right)+I\left(x-1,y\right)+I\left(x,y+1\right)+I\left(x,y-1\right)-4I\left(x,y\right)$$

从中可以看出拉普拉斯算子的由来：$$\left[\begin{array}{ccc}0& 1& 0\\ 1& -4& 1\\ 0& 1& 0\end{array}\right]$$

对图像求导过程实际是一个高通滤波过程，图像高频部分（颜色变化迅速）的响应值更大，由于图像中难免存在大量高频噪声，而导数计算器会放大噪声，因此往往需要先用一个低通滤波器平滑图像，以抑制噪声。一般采用圆对称滤波器，如高斯滤波器进行平滑，再取导数。高斯平滑与拉普拉斯求导可以合二为一，形成高斯拉普拉斯卷积核：
$$\left[\begin{array}{ccccc}0& 0& 1& 0& 0\\ 0& 1& 2& 1& 0\\ 1& 2& -16& 2& 1\\ 0& 1& 2& 1& 0\\ 0& 0& 1& 0& 0\end{array}\right]$$

Harris特征点提取

判断一个像素是否是特征所对应的位置点，可以对该点周围一个窗口$W$中的像素构成的自相关矩阵来确定（Szeliski et al., 2010）：
$$M=\sum \omega \left(x,y\right)\ast \left(\begin{array}{cc}\frac{\partial I\left(x,y\right)}{\partial x}\frac{\partial I\left(x,y\right)}{\partial x}& \frac{\partial I\left(x,y\right)}{\partial x}\frac{\partial I\left(x,y\right)}{\partial y}\\ \frac{\partial I\left(x,y\right)}{\partial y}\frac{\partial I\left(x,y\right)}{\partial x}& \frac{\partial I\left(x,y\right)}{\partial y}\frac{\partial I\left(x,y\right)}{\partial y}\end{array}\right)$$

式中，$\omega \left(x,y\right)$为定义在窗口$W$中的函数，一般为圆对称函数，如高斯函数。

图像中的大部分信息存在于梯度变化剧烈的地方，通常将自相关矩阵的最小特征值作为衡量在该点处提取特征的一个标准。

从自相关矩阵中引出三个结论：
（1）特征点的自相关矩阵的两个特征值都较大；
（2）边缘点的自相关矩阵的一个特征值较大，另一个接近0；
（3）平坦区域的点的自相关矩阵特征值都接近0。

说实话，很多综述里看到这儿的时候我都是懵逼的。？？？为啥就特征值了？？？第三条还好理解，毕竟$I_x,I_y$都接近0，所以下面我举一个$3\times 3$的窗口作为例子，
$$\left[\begin{array}{ccccc}1& 1& 1& 2& 2\\ 1& 1& 1& 2& 2\\ 1& 1& 1& 2& 2\\ 2& 2& 2& 2& 2\\ 2& 2& 2& 2& 2\end{array}\right]$$

显而易见，这是一个含有一个角点的$5\times 5$区域，取中间一块作为窗口，计算每个点的梯度，得：
$$\left[\begin{array}{ccc}(0,0)& (1,0)& (1,0)\\ (0,-1)& (1,-1)& (1,0)\\ (0,-1)& (0,-1)& (0,0)\end{array}\right]$$

然后，简单取$\omega =1$，计算$M$，得：
$$\sum _{3\times 3}\left(\begin{array}{cc}{I}_x^2& I_x{I}_y\\ I_x{I}_y& I_y^2\end{array}\right)=\left(\begin{array}{cc}4& -1\\ -1& 4\end{array}\right)$$

计算其特征值：$\left(4-\lambda \right)\left(4-\lambda \right)-1=0$
得：${\lambda }_0=15,{\lambda }_1=1$

再举一个边缘的例子：
$$\left[\begin{array}{ccccc}1& 1& 1& 2& 2\\ 1& 1& 1& 2& 2\\ 1& 1& 1& 2& 2\\ 1& 1& 1& 2& 2\\ 1& 1& 1& 2& 2\end{array}\right]$$取中间一块作为窗口，计算每个点的梯度，得：
$$\left[\begin{array}{ccc}(0,0)& (1,0)& (1,0)\\ (0,0)& (1,0)& (1,0)\\ (0,0)& (1,0)& (1,0)\end{array}\right]$$
简单取$\omega =1$，计算$M$，得：
$$\sum _{3\times 3}\left(\begin{array}{cc}{I}_x^2& I_x{I}_y\\ I_x{I}_y& I_y^2\end{array}\right)=\left(\begin{array}{cc}6& 0\\ 0& 0\end{array}\right)$$计算其特征值：$\left(6-\lambda \right)\left(-\lambda \right)=0$
得：${\lambda }_0=6,{\lambda }_1=0$

这样一来是不是直观很多。

但是，Harris角点检测方法采用的是另一个响应函数来替代最小特征值方法：
$$\xi =det\left(M\right)-\alpha \cdot trace{\left(M\right)}^2={\lambda }_0{\lambda }_1-\alpha \left({\lambda }_0+{\lambda }_1\right)$$

式中${\lambda }_0,{\lambda }_1$是$M$的两个特征值，$\alpha$取值范围通常设为$0.04\sim 0.06$。

如果响应$\xi$超过给定阈值并且是最大值，则该点为特征点（以刚才的角点区域为例，$\xi =14.2$）；
如果$\xi$的绝对值很小，则为平坦区域；
如果$\xi$的绝对值较大而且是负值，则为边缘点（以刚才的角点区域为例，$\xi =-0.3$）。

Harris特征对旋转、平移、光照变化和噪声都不太敏感，很长一段时间里都是非常重要的角点检测算子。

Fast特征点提取

Fast的特征是一种快速检测算法，其核心思想是，如果像素与周围邻域内足够多的像素的灰度差异较大，则该像素可能是角点。

Fast速度快的原因在于它会先选较少的点进行比较，再依次增多，如以点$p$为圆心，半径为3的圆周上选16个像素点为候选比较点，表示为$p_1\sim p_{16}$。

先计算中心点$p$与$p_1$和$p_9$的像素差，如果绝对值都小于阈值$t$，则该中心点非特征点；否则$p_1$和$p_9$作为候选点，并扩大比较范围：先比较中心点与$p_1,p_5,p_9,p_{13}$的差异，如果3个或以上的点满足阈值$t$，再比较中心点与所有点。

Fast特征的最大问题是可能在纹理丰富处大量聚集，有必要对所有计算得到的特征进行非极大值抑制，如要求特征的响应（如Fast特征中中心点与16个临近点的像素值差异之和）必须大于其周围半径$r$区域内的特征并超过一定阈值。

旋转不变性

旋转不变性指的是，物体发生旋转时，相关的特征仍然能够准确提取。主流的做法是对检测到的特征点周围的图像块计算一个主方向（dominant orientation）。

常见的主方向的计算方法有四种：

（1）计算特征点周围区域的平均梯度，以此梯度方向作为主方向。
（2）以梯度值最大的点的方向作为主方向；
（3）如SIFT特征，先对特征点周围梯度做一次高斯滤波，随后统计周围所有像素点对应的梯度，建立一个梯度直方图，将梯度直方图的峰值对应的梯度方向作为主方向；
（4）如ORB特征，通过求解灰度质心方向作为主方向。

尺度不变性

尺度不变性指的是图像的模糊程度对特征没有影响。实际拍摄时，可能使用不同的相机从不同的距离对图像进行拍摄，导致图像的分辨率存在很大的差异，需要排除分辨率对特征提取的影响。

实现尺度不变性的方法一般是基于原图像通过不同尺寸的高斯核卷积形成一系列不同模糊尺度的新图像，对这些图像统一进行特征提取，这样当对象出现在不同分辨率的图像中时，都能找到最佳的特征与其匹配。

图像的边缘检测常常采用高斯拉普拉斯算子操作(Laplacian of Gaussian,LoG)，LoG滤波器可以表示为：
$$LoG\left(x,y,\sigma \right)=\frac{{\partial }^2{G}_{\sigma }\left(x,y\right)}{\partial x^2}+\frac{{\partial }^2{G}_{\sigma }\left(x,y\right)}{\partial y^2}$$

但是LoG的计算比较复杂，一般采用高斯差分(difference of Gaussian,DoG)来近似LoG。高斯差分是用两个不同方差的高斯核对图像进行卷积后，再进行差分，其实就是不同尺度图像间的差分：
$$DoG\left(I,{\sigma }_1,{\sigma }_2\right)=G_{{\sigma }_1}\ast I-G_{{\sigma }_2}\ast I=\left(G_{{\sigma }_1}-G_{{\sigma }_2}\right)\ast I$$

想要保持尺度不变性，需要找到尺度空间下的极值点。在尺度空间，对图像的LoG滤波等价于DoG滤波。实际操作中，很难构建连续的尺度空间进行计算，SIFT特征通过构建高斯金字塔来近似计算极值点。

高斯金字塔的每个塔之间是降采样关系，塔内每一层都是原图像在不同$\sigma$下进行高斯卷积操作，塔内相邻层相减就形成了高斯差分金字塔，从而得到了若干个离散的高斯差分尺度空间。

此外，也可以保持图像分辨率不变，改变滤波器窗口的大小，来构造尺度空间金字塔。这种做法允许尺度空间内多层图像并行处理，并且省去了图像的降采样操作，从而提高了算法性能。