图像处理的一些基础内容

数字图像的表示

数字图像通常是利用矩阵来表示。矩阵中的每个数字代表了图像中该点的灰度值，灰度值的范围通常用 0-255 或者 0-65536 的整数表示，数值越大表示该像素点越亮。

对于灰度图像，只需要一个矩阵就能表示，对于一个彩色图像则需要 3 个矩阵来描述，每个矩阵分别代表 RGB 三种颜色的强度值。下图为 RGB 彩色立方体。

灰度变换

灰度变换是图像处理中最简单也是最常用的处理方法。灰度变换可表示为 ，其中 为输入像素， 为输出像素。其中，主要包含：反转变换、对数变换、指数变换、恒等变换等。下图显示了一些常见灰度变换的函数图像。

图像反转

反转变换表达式如下，通常用于展现黑色背景中白色或者灰色物体的细节。


上图中是一张原始 X 光图片，右图为经过反转变换后的图片。

对数变换

对数变换可以将低亮度像素区域拉伸、高亮度区域压缩，提高图像中暗处的对比度，使低亮度区域展现得更加清晰。变换函数为：


上图中左图为某图片的傅里叶频谱图片，右图为经过对数变换后的傅里叶频谱图片。

指数变换

指数变换既可以用于拉伸低亮度区域也可以拉伸高亮度区域，其变换函数如下。


下面我们看几个使用指数变换进行对比度增强的例子。
压缩高亮度区域，将低亮度区域展现得更加清晰。

上图中从左到右，从上到下依次为：人体脊椎的核磁共振图像，经过指数变换 之后的图像，经过指数变换 之后的图像，经过指数变换 之后的图像。从中可以看出，低亮度区域展现得越来越明显。
看另一个对比度增强的例子：压缩低亮度区域，将高亮度区域展现得更加清晰。

上图中从左到右，从上到下依次为：原图，经过指数变换 之后的图像，经过指数变换 之后的图像，经过指数变换 之后的图像。从中可以看出，高亮度区域展现得越来越明显。

分段线性变换

分段线性变换就是利用自定义的分段线性函数作为灰度变换。一般可以用于对比度拉伸：压缩低亮度与高亮度区域，拉伸中等亮度区域。

上图左上显示了一个分段线性函数，右上是原始图像，左下是经过灰度变换所得到的图像，右下是利用阈值分割得到的图像。
分段线性函数也可以用于灰度切片。下图显示了两个灰度变换函数。


上图中左图为主动脉血管造影图像，中间图片为使用第一个分段线性函数变换后的图像，右侧图像为使用第二个分段线性变换后的图像。

直方图处理

直方图是众多空间域处理技术的基础，它也可以用于图像增强。一幅灰度图的直方图表示了所有像素值的统计信息。对于较暗的图片，像素值分布在左侧；对于较亮的图片，像素值分布在直方图的右侧；对于对比度较高的图像，像素的分布则比较均匀。通常来说，我们想要改变一张图片的直方图，使直方图更加均匀，从而使图片的对比度更高，可视效果更好。

直方图均衡

我们先假设一张图片中所有像素点的灰度值都是连续型随机变量，记 分别为随机变量 的概率密度函数。我们令：

根据积分上限函数的微分有：

为了确保变换前后，两个分布在坐标轴上所围成的面积不变，我们有：

有：

这表明，任意的分布 经过公式 (1) 变换得到的分布 ，分布 是一个均匀分布。我们这就完成了直方图均衡的工作。

对于离散的情况，有


式中 为总像素个数， 为灰度值 出现的次数。
下图给出了不同亮度的图经过直方图均衡后的效果，可以明显的看出，图像的对比度更大了，视觉效果更好了。

直方图匹配

直方图均衡中，经过变换的像素值的分布为均匀分布，我们也可以让经过变换后的像素值的分布为一个特定的分布。这个过程称之为直方图匹配。
下面，我们令原始像素为 ，经过变换后的像素为 。我们先对原始像素进行均衡处理，有

其次，我们定义变换 ，它为对像素 进行均衡变换，有：

那么，有：

那么直方图匹配过程也就比较简单了：

1. 对原始像素 做直方图均衡得到 .
2. 根据目标分布的均衡过程得到 .
3. 将第一步中的 代入 ，得到 .

由此，我们可以看出：直方图均衡是直方图匹配的桥梁。

空间滤波

相关操作与卷积操作

相关操作与卷积操作是最常用的两个空间滤波操作，其原理可用下图以及以下公式概括：滤波器模板中的元素与原始图像中的元素对应相乘后相加得到一个值，然后移动滤波器继续相同操作，直到图像中所有的像素点都进行过这个操作。


卷积操作与相关操作略有不同，两者公式分别为：


其中 为原始图像在 处的灰度值， 表示相关或卷积 (滤波器) 模板。为了便于理解，我们给出对一维信号的相关与卷积过程。

对于图像的相关与卷积操作如下图所示：

平滑空间滤波器

使用平滑空间滤波器可以对一幅图像进行模糊操作也能去除图像中的一些噪声。一个 的平滑滤波器可以如下图所示：

由于这个滤波器是对称的，因此经过卷积操作和相关操作的结果是一样的。下图是经过不同大小的平滑滤波器滤波后的图片。

中值滤波器

中值滤波器是一种非线性滤波器，滤波操作为：取模板范围内所有像素灰度值的中间值为输出值，通常用于去除图像中的噪声。下图左图为一幅电路板的 X 光图像，中间图为使用平滑滤波器滤波的效果，右图为使用中值滤波器滤波的效果。可以看出，中值滤波器对于去噪有很良好的表现。

锐化空间滤波器

通常，我们使用二阶微分算子—拉普拉斯算子进行图像锐化操作。离散情况下的拉普拉斯算子可以表示为：

具体的，锐化操作的公式可表达为：

常见的锐化滤波器 ()，为：

下图左上为原始图片，上中图为经过拉普拉斯算子得到的图片，上右为上中图乘以某个系数得到的图片。下图分别为两张锐化后的图片与原图相加所得到的锐化后的图片。

除了二阶微分算子，也有使用一阶微分算子进行锐化图片的操作。常用的一阶微分算子为：

下图显示了一阶微分算子进行锐化操作的结果，用于提取物体的边缘。

边缘检测

我们知道锐化操作是可以得到物体的边缘的，当前还有许多进行边缘检测的算子，这些算子基本都是微分算子。我们现在来介绍一下。

常见边缘检测方法

对于一个二维的离散的函数 ，对两个方向分别求梯度 (差分)，有

根据差分公式我们就能构建不同的边缘检测算子，下面介绍三种边缘检测算子。
Roberts 算子

我们定义滤波器每个值为：

那么 Roberts 算子的计算公式为：

Prewitt 算子 (横竖方向)

该滤波器计算公式为：

Sobel 算子 (横竖方向)

该滤波器计算公式为：

当我们有了边缘检测的算子，滤波后的像素值就是两个方向梯度的幅值：

下面我们看一个使用 Sobel 算子的例子：

上图中从上到下，从左到右依次为：原图，经过 滤波后的图像， 经过 滤波后的图像， 图像。可以明显的看到， 能得到横线， 能得到竖线， 既能得到横线也能得到竖线。
此外，Sobel 算子和 Prewitt 算子都有对角线检测版本，如下图所示，分别是：Sobel 算子 (对角线方向)，Prewitt 算子 (对角线方向)。

一般来说，在边缘检测里，我们会先对图像做平滑操作，然后利用边缘检测算子，最后利用阈值分割得到最终的边缘图像。下图中我们提前对原图使用了 的均值滤波器进行滤波，然后再利用 Sobel 算子检测边缘。可以看到，墙上的瓷砖已经被抹去了，它们本来就不应该是边缘。

在配上阈值分割，可以得到下图，屋顶上的瓦片也被抹去了。

Canny边缘检测

Canny 边缘检测的步骤大致分为4步：

1. 使用高斯滤波器对原图进行平滑操作
2. 利用边缘检测算子 (梯度算子) 计算原图的梯度的幅值以及相位
3. 对梯度幅值图像进行非极大值抑制，让边缘更窄
4. 利用双阈值进行边缘连通性分析和连接边缘

1. 使用高斯滤波器对原图进行平滑操作
高斯滤波器的模板，可用以下公式获得：

利用滤波器模板对原图进行卷积操作

2. 利用边缘检测算子 (梯度算子) 计算原图的梯度的幅值以及相位
利用上面说的边缘检测算子 (Sobel等)，计算 的梯度幅值以及相位图片：

3. 对梯度幅值图像进行非极大值抑制，让边缘更窄
我们先规定边缘的四个方向：水平方向 ()、竖直方向 ()、+45°方向 ()、-45°方向 ()。假设有一个边缘它的方向角度为 ，取距离 最近的方向为该边缘的方向。这四个方向的范围如下图所示。

那么非极大值抑制过程如下：

1. 找到四个方向 中与 最近的方向
2. 如果 小于之一 ，其中是 沿着 方向相邻的两个点的梯度幅值。那么令 (抑制)，否则 。这么做的目的就是在一个区域沿着某个方向，只挑选梯度幅值最大的那一个点作为边缘，其他的点归零，最终目的是使边缘变窄。

4. 利用双阈值进行边缘连通性分析和连接边缘
使用阈值分割梯度幅值图像时，阈值选得太大就会漏掉很多边缘，选得太小又会增加一些不必要的边缘。为了解决这个问题，可以使用双阈值方法。我们先定义两个阈值 ，利用这两个阈值分别对 进行阈值分割，得到 ，由于 中包含了 中所有的点，我们令

那么 中不再包含 中的点，此时 表示强边缘点 (这个点一定是边缘)， 表示弱边缘点 (这个点可能是边缘)。

1. 定位到下一个中没有被访问过的边缘点 ，并将 标记为边缘点。
2. 找到在 中所有与 连通 (8连通) 的非零点 (弱边缘点)，将找到的这些点标记为边缘点。
3. 如果 中所有非零点都被访问过了，跳到步骤 d，否则回到步骤 a。
4. 将 中没有被标记为边缘点的点置为 0。

在上述过程中，所有标记为边缘点的点就是最终得到的边缘。

下面我们看一个例子，下图中从左到右从上到下依次为：原图，经过高斯滤波后再利用梯度算子操作再经过阈值分割得到的图像，使用 Marr-Hildreth 算法 (本节没讲) 得到的边缘，使用 Canny 算法所得到的边缘。可以看到 Canny 边缘检测算法效果还是很不错的。

一些重要的特征点

特征点由关键点（Key-point）和描述子（Descriptor）两部分组成。比方说，当我们谈论 SIFT 特征时，是指“提取 SIFT 关键点，并计算 SIFT 描述子”两件事情。

Harris 角点检测

角点是图像中的重要特征，角点这么重要的原因在于：对于边缘来说灰度值沿着边缘的方向是不会变化的，而角点附近沿着任意方向的灰度值都会急剧的变化。因此，角点是一个非常重要的特征，它具有旋转不变性。

那么，我们应该如何检测角点呢？上面说了，在角点附近灰度值都会急剧的变化，那么我们只要找到一个点它附近的灰度值急剧变化就行了。我们先用一个 的窗口，如果这个窗口沿着某个方向移动，例如， 方向，无论从哪个方向移动，整个窗口中 9 个位置的灰度值都发生了很大的变化，说明这个窗口中有角点。

利用二维的泰勒展开公式：

对 线性展开，有

有：

当无论这个窗口沿着什么方向移动时灰度值的变化都很大，那么此时 的特征值一定很大。所以我们只需要判断 的特征值大小就行了。下面是如何根据特征值判断该窗口中是否存在角点的方法。

FAST 角点检测

FAST (Features from Accelerated Segment Test) 检测也用于检测角点，它的优势在于计算复杂度很低，一般用于实时检测中。FAST 主要检测局部像素灰度变化明显的地方。它的思想是：如果一个像素与它邻域的像素差别较大（过亮或过暗）, 那它更可能是角点。检测过程如下：

1. 在图像中选取像素点 ，其灰度值为
2. 设置一个阈值
3. 以像素 为中心，选取半径为 3 的圆上的 16 个像素点。
4. 假如选取的圆上，有连续的 N 个点的亮度大于 或小于 ，那么像素 可以被认为是特征点，N 通常取 12，即为 FAST-12。其它常用的 N 取值为 9 和 11，他们分别被称为 FAST-9，FAST-11。
5. 循环以上四步，对每个像素执行相同操作。

由于 FAST 角点检测通常会在同一个区域检测到多个点，所以通常需要非极大抑制将重复的点删除。
由于 FAST 角点并没有旋转不变性和尺度不变性，因此学者提出了 Oriented FAST， 它利用图像金字塔以及灰度质心法来实现旋转不变性和尺度不变性。利用图像金字塔就很好理解了，检测角点是时候检测不同尺度图像的角点就行，匹配的时候也同时跟多个尺度的特征点进行匹配就行了。

对于旋转不变性，假设我们已经在图像中找到了一个关键点 ，我们定义在 附近的大小为 的邻域为 ，首先计算 的矩：

表示这块区域的灰度值的和 (总质量)， 表示这块区域沿 轴方向的矩 (力矩)， 表示这块区域沿 轴方向的矩 ，那么这块区域的重心点为：

由于这块邻域的中心点为 ，那么向量 的夹角为

有了这个夹角，对角点进行匹配的时候，提前把两个待匹配的区域旋转成一个角度就行了。这就满足了旋转不变性的要求。有了旋转不变性和尺度不变性的 FAST 算法称之为 Oriented FAST。

BRIEF 特征描述子

假设现在我们已经通过某种方法找到了图像中所有的感兴趣的点，例如角点。我们如何描述这些点，使它能够快速的与另一张图中的感兴趣的点进行匹配呢？
特征点描述子
假设我们已经在图像中找到了一个关键点 ，我们定义在 附近的大小为 的邻域为 。在 内随机选 个点对 ， 可以为 128、256、512。随机选择的方法可以用多种，下面介绍 5 中随机在 中选择这些 pairs 的方法。

1. 都服从二维均匀分布
2. 都服从二维高斯分布
3. 服从高斯分布 ， 服从高斯分布 ，先采样 后采样
4. 在极坐标下服从均匀分布
5. 固定不动， 在极坐标下满足均匀分布

下图展示了 5 种随机采样方法的结果：点-点对。

再 中找到了 个点对 后，依次比较这些点对的灰度值，得到 ：

将这 个 0 或 1，依次放进一个二进制数组中作为 BRIEF 描述子。

可以直接用一个十进制的数去存这个数组，计算方式为：

匹配的时候，利用汉明距离对不同的描述子进行匹配就行了。
注意：利用 Oriented FAST 提取提取特征点，并利用 BRIEF 作为特征描述子的特征点方法称为 ORB 方法，它在视觉 SLAM 中应用非常广泛，因为 ORB (Oriented FAST and Rotated BRIEF) 算法能做到实时检测与匹配。

SIFT 特征

下面我们一次介绍 SIFT 特征的处理过程。

高斯模糊
对图像金字塔中的每一层图片都进行不同程度的高斯模糊。


例如，金字塔中的每一层图片都利用 进行高斯模糊买得到几张不同程度的模糊图片。降采样时，高斯金字塔上一组图像的初始图像(底层图像)是由前一组图像的倒数第三张图像隔点采样得到的。

在边缘检测中，我们用了拉普拉斯算子对图像进行卷积操作，但一般在边缘检测之前要用高斯滤波器先对原始图像进行一次卷积，去除一些噪声。我们知道，卷积操作具有结合律，因此我们先将高斯平滑滤波器与拉普拉斯滤波器进行卷积，然后利用得到的混合滤波器去对图片进行卷积以得到所需的结果。这样就可以减少很多计算量。这个混合滤波器就是：高斯-拉普拉斯算子。
以 0 为均值， 为标准差的二维 LoG 函数表达式为：

LoG 对于特征点检测是一个非常好的算子，但是使用一下 DoG 也能近似计算 LoG，DoG 的速度更快。

DOG (Difference of Gaussian)
现在我们假设图像金字塔中一共有 4 张图片，每张图片又被 4 个不同的高斯核进行模糊处理，那么一共有 张图片。现在我们定义 DOG 为

就是在一层金字塔中的 5 张不同模糊程度的图片中，相邻两张图片的差，如下图所示。

找到关键点
现在，我们要判断一个点是否为 DOG 空间中的局部的极值点。在一个点附近的 邻域中有 8 个像素点，它的上一层和下一层还有对应的 9 个像素点，因此一共要对比 26 个像素点。如果一个点在其立方体的邻域中是局部极值点，那么它就是潜在的特征点。

下图展示了上述三步的过程，注意：降采样时，高斯金字塔上一组图像的初始图像(底层图像)是由前一组图像的倒数第三张图像隔点采样得到的。


删除不好的极值点**
通过找到关键点步骤中所得到的所有关键点，有些可能位于边上，有些极值点可能不是那么强 (极值点与邻域的灰度值相差太小)，这些点是需要被去除的点。此外离散空间里找的极值点可能不是连续空间中真正的极值点，可能在真正的极值点的旁边。

为了提高关键点的稳定性，需要对当前尺度空间（这个尺度空间有图像金字塔和高斯核共同确定）的 DOG 函数进行曲线拟合。对 DOG 在极值点附近进行泰勒展开。注意到下式中 都为 0，因为在极值点展开的。

在原始论文里， 的点被剃除。此外还要删除，在边缘上的点。
我们知道：曲面在某点的主曲率与该点 Hessian 矩阵的特征值成正比：


其中 为矩阵的两个特征值。设 ，表示大特征值与小特征值的比值，则

当特征值的比值很大时，说明这个点是边 (参考 Harris 角点检测)，此时这个点要去除。因此，当满足下式时，特征值比值很大，这个点要删除：

原始论文里 。

求特征点主方向
通过以上步骤，我们已经得到了在不同尺度下所有图片的特征点，实现了尺度不变性，为了实现旋转不变性，需要定义每个特征点的方向。利用特征点邻域像素的梯度分布特性来确定其方向参数，再利用图像的梯度直方图求取关键点局部结构的稳定方向。
计算以特征点为中心、以 为半径的区域图像的每个像素点的梯度的幅值 与相位 ，计算公式如下：

计算得到梯度方向后，就要使用直方图统计特征点邻域内像素对应的梯度方向和幅值。梯度方向的直方图的横轴是梯度方向的角度，纵轴是梯度方向对应梯度幅值的累加，在直方图的峰值就是特征点的主方向。

直方图中，值大于 80% 就可以当做当前这个特征点的主方向。某个点也可能有多个主方向。有了特征点的主方向，匹配的时候先旋转再匹配，就可以满足旋转不变性了。

特征描述子
找到了特征点，就可以定义特征描述子了，得到特征描述子一般具有三个步骤：

1. 矫正旋转主方向，确保旋转不变性。
2. 得到描述子向量。
3. 归一化描述子向量，确保不受光照影响。

为了保证特征矢量的旋转不变性，要以特征点为中心，在附近邻域内将坐标轴旋转一定角度，即将坐标轴旋转为特征点的主方向。旋转后以主方向为中心取 的窗口。每 的窗口构成一个子窗口，子窗口中包含 8 个像素点，我们计算了每个像素点的梯度的幅值以及相位，每个像素点中的箭头方向代表该像素的梯度方向，长度代表梯度幅值，然后利用高斯窗口对其进行加权运算。最后在每个 的小块上绘制 8 个方向的梯度直方图，计算每个梯度方向的累加值，即可形成一个种子点，如右图所示。每个特征点由 16 个种子点组成，每个种子点有 8 个方向的向量信息。这种邻域方向性信息联合增强了算法的抗噪声能力，同时对于含有定位误差的特征匹配也提供了比较理性的容错性。


每个关键点，可以用 128 维向量表示： 16 子窗口 (右上图)，每个子窗口由 8 个值描述 (8 个方向的统计值)。匹配的时候直接对比这 128 维的向量就行了。**

SURF 特征

SURF (Speeded-Up Robust Features) 是一种更快的特征点，我觉得比 SIFT 更好。下面我们顺序介绍 SURF 特征的流程。
积分图像
一幅图像 对应的积分图像 ，定义如下：

Hessian 矩阵行列式图片
对于经过不同的高斯核滤波的图片 ，它的 Hessian 矩阵定义为：

对于 的计算，可以使用高斯-拉普拉斯滤波器进行滤波，但是这个计算量有点大，在 SIFT 里使用了 DoG 代替这个步骤，而在 SURF 里则采用了更快的滤波器：盒子滤波器。盒子滤波器就是简化版本的 LoG 滤波器，而且非常简单粗暴，如下图所示：

上图，图 1 为使用 LoG 滤波器的 ，对应的盒子滤波器为图 3；图2中为使用 LoG 滤波器的 ，对应的盒子滤波器为图 4。因为使用了更加简化的 LoG 滤波器，所以 SURF 更快了。盒子滤波器我们记为：

那么 Hessian 矩阵的行列式图像为：

其中， 可以取 0.9。使用积分图像完成 Hessian 矩阵的行列式的计算会非常的快。下图展现了行列式图像的获得过程：


SURF 的尺度空间**
与 SIFT 不同的是，SIFT 在获取图像金字塔时，只能先得到最后一层，对原始图片进行不同程度的滤波 (改变 的值)，然后将倒数第三张滤波后的图片进行降采样，作为第二层的最后一张图，以此内推，最后得到图像金字塔。这个过程是串行的。
SURF算法先从 的盒子滤波器开始，对盒子滤波器的尺寸进行扩展。在 SURF 中，保持图像不变，仅仅改变盒子滤波器窗口的大小以及采样间隔来获得不同尺度的图像，即构成了尺度空间。不同组、不同层的图像都是可以并行获得的。所以，可以得到多个不同尺度的行列式图像。一般滤波器大小变化为：


这样也就构成了图像金字塔，与 SIFT 一样图像金字塔有多组，每组有多个层。要知道有了积分图像，滤波是相当快的。

找特征点
SURF 找特征点的方式与 SIFT 一模一样，同时在空间 邻域中采用了非极大值抑制，去除不必要的特征点。

此外，还设定了 Hessian 行列式的阈值，如果该点小于某个阈值，则该点不能选为特征点。

特征描述
为了使特征具备较好的旋转不变性，需分配给每个特征点一个主方向。在 SURF 中，主方向是对以特征点为中心的 6 倍特征尺度为半径的圆形区域内的Haar小波响应做统计运算得到的。
即以特征点为中心，计算半径为6s (s为特征点所在的尺度值) 的邻域内，统计 60 度扇形内所有点在水平和垂直方向的Haar小波响应总和 (Haar小波边长取4s)，并给这些响应值赋高斯权重系数，使得靠近特征点的响应贡献大，而远离特征点的响应贡献小，然后60度范围内的响应相加以形成新的矢量，遍历整个圆形区域，选择最长矢量的方向为该特征点的主方向。

有了主方向，就可以构造特征描述子了。在 SURF 中，我们在关键点周围选取一个正方形框，方向为关键点的主方向，边长为 20s。将其划分为16个区域（边长为5s），每个区域统计 25 个像素的水平方向和垂直方向的Haar小波特性（均相对于正方形框的主方向确定的），我们记 分别为某点沿与主方向水平的方向的相应、沿与主方向竖直方向的相应以及这两个相应的绝对值。统计他们的和得到

那么这个 的区域的统计信息可以由 4 个值来描述，一共有 16 个这样的区域，因此 SURF 可以用一个 128 维的向量表示。

Reference

• https://medium.com/data-breach/introduction-to-sift-scale-invariant-feature-transform-65d7f3a72d40
• https://www.cnblogs.com/wangguchangqing/p/4853263.html
• https://blog.csdn.net/molihong28/article/details/89296677