【计算机视觉】6. 图像处理的基础知识
图像的基本处理
在计算机视觉这一领域诞生的初期，一种普遍的研究范式，是将图像看做二维的数字信号，然后借用数字信号处理中的方法，这就是数字图像处理（digital image processing）。本节数字图像处理方面的算法进行介绍。
图像直方图与像素距离
传感器获取的图像是平面上的连续函数，将连续函数采样(sampled)为M行N列的矩阵，将每个连续样本量化(quantization)为一个整数值，即图像函数的连续范围被分成了K个区间，采样后得到的矩阵构成了离散图像，栅格中无限小的采样点对应于数字图像中的像元，即像素(pixel) 。
灰度图像中，最低值对应黑，最高值对应白；黑白之间的亮度值是灰度阶(gray-level) ，彩色图像则通过矢量函数（三阶张量）描述，可以将一幅彩色图像看做由R,G,B三种基础色进行堆叠形成，而这三种基础色又对应了三个大小相同的矩阵，矩阵的数值表征这一通道颜色的深浅。有时除了考虑RGB三种颜色外，还考虑像素的透明度a，称为RGBA描述。
色彩在人类视觉感知中极其重要，色彩与物体反射不同波长的电磁波的能力相关，一般将这三种颜色（三种不同波长的光）作为三原色：红(700nm)、绿(546.1nm)、蓝(438.5nm)，灰度图像的矩阵元素数值与彩色图像间满足Y=0.299R+0.587G+0.114B；RGB数字图像中，以(0,0,0)表示黑色，(255,255,255)表示白色；灰度图像中，以0表示黑色，以255表示白色；二值图像中，以0表示黑色，1表示白色。
图像每个位置[i,j]必定对应一个[0,255]的数值，统计每个数值所对应的像素点个数可以得到图像的亮度直方图，直方图（亮度直方图，brightness histogram）给出了图像中各个亮度值出现的概率，一幅k阶图像的直方图由具有k个元素的一维数组表示。
Def.(灰度级变换) 灰度级变换是这样一种（线性）变换，它将原来在范围[p0,pk]内的亮度p变换为一个新范围[q0,qk]内的亮度。
直方图均衡化（histogram equlization）的目标是创建一幅在整个亮度范围内具有相同亮度分布的图像，输入直方图H[p]，输入亮度范围为[p0,pk]，直方图均衡化的目标是找到一个单调的像素亮度变换q=T(p)，使输出直方图G[q]在整个输出亮度范围[q0,qk]内是均匀的；增强了靠近直方图极大值附近的亮度的对比度，减小了极小值附近的对比度。在对图像做进一步处理之前，直方图均衡化通常是对图像灰度值进行归一化的一个非常好的方法，并且可以增强图像的对比度，原先图像灰色区域的细节变得清晰。
直方图均衡化的步骤如下：
1. 对于k(256)个亮度级、大小为M × N的图像，创建长为k的数组H，初始化为0
2. 形成图像直方图H
3. 形成累计直方图Hc，Hc[p]=Hc[p-1]+H[p], Hc[0]=H[0]
4. 设置 T(p)=\text{int} \frac{(k-1)Hc[p]}{MN}
5. 重新扫描图像，根据查找表获得变换结果
（更加详细的解释参考：说说直方图均衡化）
坐标为(m,n)和(h,k)两点之间的距离有如下定义
- [欧氏距离（Euclid distance）] d_e=\sqrt{(m-h)^2+(n-k)^2 }
- [城区距离（曼哈顿距离，Manhattan distance）] d_4=|m-h|+|n-k|
- [棋盘距离（chess distance）] 因相当于国际象棋中王走一步能到达的位置而得名： d_8=max (|m-h|,|n-k|)
Def.[像素邻接性(adjacency)]
4－邻接：任意两像素之间的距离为D4＝1
8－邻接：任意两像素之间的距离为D8＝1
Def.[区域（region）] 由一些彼此邻接的像素组成的集合
Def.[连通性（continuous）] 一幅图像的两个像素之间存在一条路径；连通关系具有自反性、对称性和传递性；区域中没有孔，称为简单连通，有孔的区域称为复连通
距离变换（distance transform）给出图像中的每个像素与某个图像子集的距离，步骤如下：
1. 创建大小为M 乘 N的矩阵F，子集S初始化为0，其他位置初始化为INF
2. 按行遍历图像，从上到下，从左到右，对上方和左侧邻接像素（AL集合）计算 F(p)=\min_{q \in AL} [F(p),D(p,q)+F(q)]
3. 按行遍历图像，从下到上，从右到左，对下方和右侧邻接像素（BR集合）计算 F(p)=\min_{q \in BR} [F(p),D(p,q)+F(q)]
4. 返回数组F的结果
图像阈值化与二值图像的处理
Def.[图像二值化（阈值化）] 设输入图像为f，输出为g，则
$$g(m,n)=
\begin{cases}
1 & f(m,n)>t \ 0 & f(m,n) \ge t
\end{cases}
$$
其中t称为阈值（thresh）。一种最优阈值选择算法，通过不断迭代直至收敛的方法，寻找一个最优的阈值。
另一种选择方法，Otsu方法，选择使得类间方差最大的值作为阈值，设图像的归一化直方图为p，亮度阶为k（例如元素值取0~255，则k=256），则设
s_0=\sum_{i=0}^{t-1} p_i ,s_t=\sum_{i=t}^{k-1} p_i ,s=\sum_{i=0}^{k-1} p_i =1
令期望 \mu_0=\sum_{i=0}^{t-1} i \frac{p_i}{s_0},\mu_t=\sum_{i=t}^{k-1} i \frac{p_i}{s_t},\mu=\sum_{i=0}^{t-1} i \frac{p_i}{s}=\sum_{i=0}^{t-1}ip_i
则类间方差 \sigma_b^2=s_0 (\mu_0-\mu)^2+s_t (\mu_t-\mu)^2
数学形态学（morphology）又称图像代数（image algebra），是以形态为基础对二值图像进行分析的数学工具，基本思想是用具有一定形态的结构元素去度量和提取图像中的对应形状，以达到图像分析和识别的目的。
Def.[膨胀（dilation）] 将结构元素B的原点移至集合A的某一点，将结构元素B中的点与该点取“或”（只要两者中有一个值为1，结果就为1），得到对集合中该点的膨胀结果，对集合A中所有元素重复该过程，用$A \oplus B$表示。
Def.[腐蚀（erosion）] 将结构元素B的原点移至集合A的某一点，对两者进行“与”运算，用$A \ominus B$表示。
以下对比了膨胀和腐蚀操作的结果：
Def.[开运算（opening）] 开运算（opening）是指先腐蚀再膨胀的过程。开运算是结构元素B在集合A内完全匹配的并集，完全删除了不能包含结构元素的对象区域。开运算使图像的轮廓变得光滑，断开狭窄的间断，去掉细小的突出物。
Def.[闭运算（closing）] 闭运算（closing）是先膨胀再腐蚀的过程。结构元素B在集合A外侧平移，这些平移的并集就是闭运算的结果。闭运算使图像的轮廓变得光滑，但与开运算不同的是，它会将狭窄的缺口连接形成细长的弯口，并填充比结构元素小的洞。
Def.[形态学梯度运算] 求取膨胀与腐蚀结果间的差分，其效果是突出边缘轮廓。
Def.[顶帽（top hat）变换] 求取输入图像与其开运算间的差分。
Def.[黑帽（black hat）变换] 求取输入图像闭运算与其本身间的差分。
在灰度图像中，一幅图像的膨胀/腐蚀运算定义为：对每个像素赋值为其邻域内输入图像灰度级的最大值/最小值，二值形态学的结构元素代表一个邻域，灰度图像形态学的结构元素是一个二元函数，它规定了希望的局部灰度级性质。
连通区域检测是图像处理、模式识别中常用的一个基本方法，在目标分割，边缘检测，区域检测中有着广泛的应用。
首先，连通区域是对二值图像进行处理的，即，该图像只有黑（0）和白（255）两种颜色，这里，假设目标为白色，背景为黑色。标记算法首先对二值图像进行一次完整的扫描，标记所有目标像素点的同时，得到并记录等价标记对。等价标记对（以下简称等价对）的产生是由于扫描次序的不同，导致开始时认为是两个不同的连通区域，后来随着扫描的深入，又发现这两个区域是连通的。所以，需要记录等价对，以表明它们隶属于同一个连通区域，以便第一次扫描结束后进行修正。标记算法首先对二值图像的每一个像素进行8连通区域的标记，即:对任意一个像素的上、下、左、右、左上、右上、右下、左下，共8个相邻像素进行比较。
由于不是每个像素都有8个相邻像素，对于一些特殊位置的像素点需要特殊考虑，其中包括：
- 二值图像左上角的像素，由于是第一个要扫描的像素，无需进行8连通区域的检测，也无需考虑记录等价对的问题。
- 二值图像第一行的像素，只需要考虑左边相邻像素的连通性，无需考虑记录等价对。
- 二值图像第一列的像素，只需要考虑上和右上2个相邻像素的连通性。
- 二值图像最后一列的像素，只需要考虑左、左上、上3个相邻像素的连通性。
除了以上4种情况，其它像素，都需要考虑其8个相邻像素的连通性，如果出现不同连通标记的相邻像素，还需要考虑记录等价对的问题。步骤如下：
1. 标记图像左上角，即，第一行第一列的像素。如果其像素值为255，则标记该点的值为1，否则，开始扫描第一行第二列的像素。
2. 标记第一行的其它像素，此时，不会产生等价对的情况，不必考虑记录等价对。对该行的每一个像素，如果其值为255，检测左边像素是否为255，若是，则该点标记为左边像素点的标记；否则，该点的标记为前一个标记值加一；若该点的像素值为0，继续扫描下一个像素。
3. 对除了第一行以外的像素行进行标记，此时会出现等价对的情况，需要进行记录。
4. 首先对第一列进行处理，若该点像素值为0，则扫描该行下一个像素，否则，检测上、右上两个像素位置的像素值。若上被标记过，该点标记为上像素点的标记值。这时，再看右上是否被标记过，若也被标记过，比较上和右上的标记值是否相等，如果不相等，则记录上和右上为一个等价对，并将其记录在等价对记录表中。若上没有被标记，而右上被标记了，则该点标记为右上的标记值。如果上和右上都没有被标记，该点的标记值为上一个标记值加一。
5. 对中间列进行处理，若该像素的像素值为255，则检测左、左上、上、右上位置的像素值。若上述四个位置的像素值都为0，则该点的标记值为上一个标记值加一。如果上述四个位置中只有一个的像素值为255，则该点就标记为那个像素点的标记值。如果其中有m（m大于1，小于等于4）个像素点的像素值为255，则按照左、左上、上、右上的优先顺序来确定该点的标记值，然后对这m个像素位置的标记值进行等价对的分析，并进行相应的记录。
6. 对最后一列进行处理，步骤同上。
7. 依次扫描，直到所有像素值都被扫描。
8. 对等价记录表中的所有等价对进行处理，得到最终的连通区域标记。
图像卷积与滤波
一维信号的卷积（convolution） g(t)=\int_0^{\infty} f(t-\tau)h(\tau) d \tau
离散化得 g(m)=\sum\limits_{k=0}^{M_1-1} f(k)h(m-k)
定义二维函数f(x,y)与h(x,y)的卷积 f*h=\int_{-\infty}^{+\infty} \int_{-\infty}^{+\infty} f(x-a,y-b)h(a,b)dadb
在离散系统中，定义矩阵卷积
g(m,n)=\sum_{k=0}^{M_2-1} \sum_{l=0}^{N_2-1} f(m-k,n-l)h(k,l)=\sum_{k=0}^{M_1-1} \sum_{l=0}^{N_1-1} f(k,l)h(m-k,n-l)
其中宽（矩阵列数）M1=f.shape[1]，高（矩阵行数）N1=f.shape[0]，输出矩阵尺寸是M1+M2-1,N1+N2-1，称为full形式；same返回尺寸为M1,N1的矩阵；valid形式输出矩阵尺寸是M1-M2+1,N1-N2+1
利用傅里叶变换的卷积定理（频域相乘相当于时域卷积，时域相乘相当于频域卷积）可以实现卷积与 FFT 的互相转化。既然卷积与傅里叶变换可以互化，那么将图像与某些特定的卷积核（kernel）进行卷积（或将这些矩阵看做卷积掩膜， mask），就可以实现滤波的效果。
图像差分（图像梯度，gradient）显然可以由 I_x(x,y)=I(x+1,y)-I_x(x-1,y)
计算，如果写成卷积形式，即 I_x=I*h_x, I_y=I*h_y
就可以得到Roberts算子
h_x=\begin{pmatrix} -1 & 0 \\ 0 & 1 \end{pmatrix},h_y=\begin{pmatrix} 0 & -1 \\ 1 & 0 \end{pmatrix}
假如同时考虑用周围8邻域的数据来计算，就可以得到Prewitt算子
h_x=\begin{pmatrix} -1 & -1 & -1 \\ 0 & 0 & 0\\ 1 & 1 & 1 \end{pmatrix},h_y=\begin{pmatrix} -1 & 0 & 1 \\-1 & 0 & 1 \\ -1 & 0 & 1 \end{pmatrix}
假如对中间一行（列）进行加权，就得到Sobel算子
h_x=\begin{pmatrix} -1 & -2 & -1 \\ 0 & 0 & 0\\ 1 & 2 & 1 \end{pmatrix},h_y=\begin{pmatrix} -1 & 0 & 1 \\-2 & 0 & 2 \\ -1 & 0 & 1 \end{pmatrix}
如不声明，图像差分是指图像与Sobel算子的卷积，其幅值与相角
G(x,y)=\sqrt{I_x^2 (x,y)+I_y^2 (x,y)}, \theta (x,y)=\tan^{-1} \frac{I_y^2 (x,y)}{I_x^2 (x,y)}
图像积分
J(x,y)=\sum_{a<x,b<y} I(a,b)
假设每个像素上的噪声是一个均值为0，标准差为σ的独立随机变量，则可通过多次采集相同的静态景物来获得一幅平均图像，平均图像中的噪声仍是随机变量，均值为0，标准差为$ \sigma/\sqrt{k} $，若只能获得一幅带有噪声的图像，则通过图像的局部邻域平均完成滤波，这就是均值滤波，均值滤波属于线性滤波，如果噪声大小小于图像中感兴趣的最小尺寸，处理结果是可以接受的，但仍存在边缘模糊的问题，卷积核
k=\frac{1}{9}\begin{pmatrix} 1 & 1 & 1 \\ 1 & 1 & 1\\ 1 & 1 & 1 \end{pmatrix}
有时也采用
k=\frac{1}{10}\begin{pmatrix} 1 & 1 & 1 \\ 1 & 2 & 1\\ 1 & 1 & 1 \end{pmatrix}\quad \text{or} \quad k=\frac{1}{16}\begin{pmatrix} 1 & 2 & 1 \\ 2 & 4 & 2\\ 1 & 2 & 1 \end{pmatrix}
Def.(噪声) 实际图像常受一些随机误差的影响而退化，这些退化称为噪声(noise)，在图像获取、传输或处理中都可能出现噪声，噪声可能依赖于图像内容，也可能与其无关。
Note.
- [白噪声(white noise)] 具有常量的功率谱，强度不随频率的增加而衰减
- [高斯噪声(Gaussian noise)] 服从正态分布的随机变量
- [加性噪声(additive noise)] 噪声与出现的图像信号无关，g=f+noise
- [乘性噪声(additive noise)] 噪声与出现的图像信号有关，g=f(1+noise)
- [信噪比] 信号平方和与误差平方和之比 SNR=\frac{\sum_{(m,n)} f^2 (m,n)}{\sum_{(m,n)} noise^2 (m,n)}
椒盐噪声（salt-and-pepper noise）又称脉冲噪声，它随机改变一些像素值，在二值图像上表现为使一些像素点变白，一些像素点变黑；是由图像传感器，传输信道，解码处理等产生的黑白相间的亮暗点噪声，往往由图像切割引起。
椒盐噪声是指两种噪声，一种是盐噪声（salt noise），另一种是胡椒噪声（pepper noise）。盐为白色，胡椒为黑色，因此盐噪声又称白噪声，胡椒噪声又称黑噪声。前者是高灰度噪声，后者属于低灰度噪声。一般两种噪声同时出现，呈现在图像上就是黑白杂点。
Note. 其他均值滤波器
- [几何均值滤波器] f(x,y)=\left(\prod_{(s,t)\in S_{xy}} g(s,t)\right)^{1/mn} 与算术均值滤波器相比，它的优势在于在滤波过程中丢失的图像细节更少
- [谐波均值滤波器（调和均值滤波器）] f(x,y)={mn}{\prod_{(s,t)\in S_{xy}} \frac{1}{g(s,t)}} 处理高斯噪声和盐噪声效果好，但不适合胡椒噪声
- [逆谐波均值滤波器] f(x,y)=\frac{\prod_{(s,t)\in S_{xy}} g(s,t)^{Q+1}}{\prod_{(s,t)\in S_{xy}} g(s,t)^{Q}} 其中Q为滤波的阶数。它适合消除或减少椒盐噪声，当Q>0时，用于消除胡椒噪声，当Q<0时，用于消除盐噪声，因此不能同时消除这两种噪声
关于高斯滤波和拉普拉斯滤波，本文参考了高斯金字塔与拉普拉斯金字塔 和 拉普拉斯算子
高斯滤波基于高斯函数的形状形成卷积掩模，对于去除服从正态分布的噪声很有效；高斯滤波的特点：
- 高斯函数具有旋转对称性，高斯滤波在各个方向上的平滑程度是相同的，从而使后续的边缘检测不会偏向某一方向；
- 高斯函数是单值函数，高斯滤波器用像素邻域的加权均值来代替该点的像素值，而每一邻域像素点权值是随该点与中心点的距离单调增减的；
- 高斯函数具有特殊的形状，可通过调节高斯函数的标准差，在图像过于模糊（过平滑）与图像噪声/细纹理过多（欠平滑）之间进行折中；
- 由于高斯函数的可分离性，大高斯滤波器可以得以有效地实现。二维高斯函数卷积可以分两步来进行，首先将图像与一维高斯函数进行卷积，然后将卷积结果与方向垂直的相同一维高斯函数卷积．因此，二维高斯滤波的计算量随滤波模板宽度成线性增长而不是成平方增长；
- 高斯函数的DFT结果仍是高斯函数，因此，高斯函数的傅里叶变换频谱是单瓣的；这意味着，平滑图像不会被不需要的高频信号所污染，高斯滤波器是一种低通滤波器，图像中的噪声一般是高频信号，高斯滤波可以在一定程度上抑制较高频率的噪声，同时保留大部分所需的低频信号。
图像的二阶导数一般用拉普拉斯算子表示，即 \nabla^2 I=\frac{\partial^2 I}{\partial x^2 }+\frac{\partial^2 I}{\partial y^2 }
由一维情况下 \frac{\partial^2f}{\partial x^2} = f(x+1)-f(x)-f(x)+f(x-1)
\frac{\partial^2f}{\partial y^2} = f(y+1)-f(y)-f(y)+f(y-1)
那么得到上述二阶微分：
\nabla^2f = \frac{\partial^2f}{\partial x^2}+\frac{\partial^2f}{\partial y^2}\\=[f(x+1,y)+f(x-1,y)+f(x,y+1)+f(x,y-1)]-4f(x,y)
上述是数学表达形式的拉普拉斯算子，那么我们可以将其表达为卷积核的形式，它相当于I与如下算子进行卷积
L= \begin{pmatrix} 0 & 1 & 0\\ 1 & -4 & 1\\ 0 & 1 & 0 \end{pmatrix}
以上是基本的拉普拉斯算子，其变形有
L_1=\begin{pmatrix} 1 & 1 & 1\\ 1 & -8 & 1\\ 1 & 1 & 1 \end{pmatrix}\quad L_2= \begin{pmatrix} 0 & -1 & 0\\ -1 & 4 & -1\\ 0 & -1 & 0 \end{pmatrix}\quad L_3=\begin{pmatrix} -1 & -1 & -1\\ -1 & 8 & -1\\ -1 & -1 & -1 \end{pmatrix}
如果在图像中的一个比较暗的区域出现了一个亮点，那么经过拉普拉斯算子处理后，这个亮点将变得更亮，因为在一个很暗的区域内，很亮的点和其周围的点属于差异比较大的点，基于二阶微分的拉普拉斯算子就是求取这种像素值发生突然变换的点或线.此算子却可用二次微分正峰和负峰之间的过零点来确定，对孤立点或端点更为敏感，因此特别适用于以突出图像中的孤立点、孤立线或线端点为目的的场合。
同梯度算子一样，拉普拉斯算子也会增强图像中的噪声，有时用拉普拉斯算子进行边缘检测时，可将图像先进行平滑处理。但是在进行锐化的过程中，我们又不希望这个滤波器改变了图像中其他pixel的信息，所以保证了每个滤波器的元素数值和加起来为0
如果拉普拉斯掩模中心元素为负，那么拉普拉斯锐化 f(x,y)\leftarrow f(x,y)-\nabla^2 f(x,y) 如果拉普拉斯掩模中心元素为正，那么拉普拉斯锐化 f(x,y)\leftarrow f(x,y)+\nabla^2 f(x,y) 这种简单的锐化方法既可以产生拉普拉斯锐化处理的效果，同时又能保留背景信息，将原始图像叠加到拉普拉斯变换的处理结果中去，可以使图像中的各灰度值得到保留，使灰度突变处的对比度得到增强，最终结果是在保留图像背景的前提下，突现出图像中小的细节信息。
常用的平滑滤波器：高斯函数，满足
h(x,\sigma )=\frac{1}{\sqrt{2\pi} \sigma } \exp \left(-\frac{x^2}{2\sigma ^2} \right)
推广到二维
h(x,y,\sigma )=\frac{1}{2\pi \sigma^2 } \exp \left(-\frac{x^2+y^2}{2\sigma ^2} \right)
那么一阶微分滤波器
\begin{gather} \frac{\partial h}{\partial x}=-\frac{1}{\sqrt{2 \pi} \sigma }\exp \left(-\frac{y^2}{2\sigma^2 } \right) \cdot \frac{x}{\sqrt{2 \pi} \sigma ^3 } \exp \left(-\frac{x^2}{2\sigma^2 }\right)\\ \frac{\partial h}{\partial y}=-\frac{1}{\sqrt{2 \pi} \sigma } \exp \left(-\frac{x^2}{2\sigma ^2 } \right) \cdot \frac{x}{\sqrt{2 \pi} \sigma ^3 } \exp \left(-\frac{y^2}{2\sigma^2 }\right) \end{gather}
二阶微分滤波器（高斯-拉普拉斯算子，Laplacian of Gaussian, LoG）
\nabla^2 h=\frac{1}{2\pi \sigma ^4 } \left( \frac{x^2+y^2}{\sigma ^2} -2\right) \exp \left(-\frac{x^2+y^2}{2\sigma ^2} \right)
Note.
- 线性平滑的问题：图像的边缘被模糊了
- 非线性平滑方法：根据指定标准，将当前像素邻域分成两个子区域，仅使用邻域中与被处理像素有类似性质的点完成平滑
- 限制数据有效性下的平均仅使用满足某种标准的像素做平均，从而避免模糊；设定非法数据范围[min,max]，只有在[min,max]内的像素值才被其邻域的平均所取代，只有有效的数据才对邻域的平均有贡献
使用旋转掩模的平均通过搜索当前像素邻域的一致性部分来避免边缘模糊，像素的均值操作只在具有一致性的区域内进行；设区域R的像素数目是n，且输入图像是f，用亮度散布度量区域的一致性:
\sigma^2 = \frac{1}{n} \sum_{(i,j) \in R} [f(i,j)-\frac{1}{n} \sum_{(i,j) \in R} f(i,j)]^2
这样，比掩模小的细节信息将会损失，算法可以迭代使用，迭代过程会较快收敛到一个稳定状态。
中值滤波用像素点邻域灰度值的中值来代替该像素点的灰度值，不依赖于邻域内与典型值差别很大的灰度值，可以减少边缘的模糊，对去除椒盐噪声非常有效。
中值滤波的步骤：
- 确定掩模大小和掩模中心
- 在掩模内将像素点按亮度值大小排序
- 选择序列的中间值作为掩模中心的新像素值，如果像素点数为偶数，中值取排序像素中间两点的平均值。
不同于图像积分，箱式滤波的数组A中的每个元素的值是该像素邻域内的像素和（或像素平方和），在需要求某个矩形内像素和的时候，直接访问数组中对应的位置就可以了。因此可以看出它的复杂度是O(1)，步骤如下：
1. 给定一张图像，宽高为（M,N），确定待求矩形模板的宽高(m,n)
2. 开辟一段大小为M的数组，记为buff, 用来存储计算过程的中间变量；
3. 将矩形模板（紫色）从左上角（0，0）开始，逐像素向右滑动，到达行末时，矩形移动到下一行的开头（0，1），如此反复，每移动到一个新位置时，计算矩形内的像素和，保存在数组A中；
4. 计算sum[i] = sum[i-1] - buff[x-1] + buff[x+m-1]
5. 对buff进行更新。加上一个新进来的像素，再减去一个出去的像素，然后便开始新的一行的计算了。
Note.(ROF（Rudin-Osher-Fatemi） 去噪模型) 一幅（灰度）图像 I 的全变差（Total Variation， TV）定义为梯度范数之和。在连续表示的情况下，全变差表示为： J(I)=\int \nabla I(x,y)dxdy
在 Chambolle 提出的 ROF 模型里，目标函数为寻找降噪后的图像 U，使下式最小： \min_U ||I-U||^2+2\lambda IU
Note.(方向滤波器) 设(x',y')是(x,y)经仿射变换旋转θ之后的坐标，则 I(x,y)*f(\theta)=\frac{1}{c} \int_{-\infty}^{+\infty} h(t)I(x',y')dt \\\text{where }h(t)=\exp (-\frac{t^2}{2\sigma^2}),c=\int_{-\infty}^{+\infty} h(t)dt
Note.(自适应滤波器)
f(x,y)\leftarrow f(x,y)-\frac{\sigma_1^2}{\sigma_2^2}[f(x,y)-\mu]
其中 \sigma_1^2,\sigma_2^2 分别表示噪声方差和图象局部方差，μ为均值。
图像的分裂、归并、分割
在处理图像的过程中，由于图像中某个像素与相邻像素之间的有很强的相关性，即不管是从纹理还是从灰度级都很相似。所以如果物体的尺寸很小或者说对比度不高，通常则需要采用较高的分辨率来观察。如果物体的尺寸很大或者说对比度很强，那么就仅仅需要较低的分辨率就能够来传观了。如果现在物体的尺寸有大有小，对比度有强有弱，这些关系同时存在，这个时候，只能用多分辨率进行处理。图像金字塔就是为了以多分辨率来解释图像而诞生的一种简单有效的方法。
一幅图像的金字塔，是以一系列以金字塔形状排列的分辨率初步降低的图像的集合。金字塔的底部是待处理图像的高分辨率的表示，而顶部是低分辨率的表示；当金字塔向上层移动的时候，尺寸和分辨率都会降低。
假设基础级（也就是最底层）的尺寸为 2^J \times 2^J 即 N \times N ,J = \log_2 N ，那么金字塔中间任意一级j的尺寸大小为： 2^j \times 2^j, j \in [0,J] ，所以说，一个完整的图像金字塔可以由 J+1 个分辨率集所组成。
Def.[区域归并] 根据预先定义的标准，把像素或者子区域集合成较大区域，特殊的算法区别在于初始分割的定义和归并标准，区域归并的结果通常与归并的顺序有关。步骤如下：
- 定义某种初始化方法将图像分割很多具有一致性的小区域
- 为归并两个邻接区域定义一个标准
- 将满足归并标准的所有邻接区域归并起来，如果不再有两个区域归并后满足一致性要求，则停止
区域归并的基本方法：使用区域 2\times 2, 4\times 4, 8\times 8 像素的区域作为起始，使用区域灰度直方图描述区域一致性，将区域描述与相邻区域的描述进行比较，如果相互匹配，就将它们归并成一个较大的区域且计算新区域的描述；否则，区域标记为非匹配。对所有区域的相邻区域（包括新生成的区域）连续地进行归并，如果一个区域不能和其任何邻接的区域归并，则将其标记为“最终”，当所有区域都被标记为最终的，归并结束。
区域分裂是区域归并的反过程，首先将整个图像作为一个区域，然后连续地对每个存在的区域进行一致性测试，对不满足测试条件的区域进行分裂，直到所有区域均满足一致性条件。注意，即使使用相同的一致性测试条件，也不能保证区域归并和区域分裂的结果相同。
区域分裂与归并结合的方法可以集中两种方法的优点。分裂与归并方法使用金字塔形图像表示，区域是方形的，对应着金字塔表示的某层元素：
- 分裂 如果在任一层次上（不包括最低层）的任何区域不满足一致性测试条件，则将其分裂为四个子区域，它们是下一层较高分辨率元素，分裂过程对应于金字塔的降采样（pyramid-down）过程，即自上而下（top-down）的采样，可以使用差分高斯算子（DoG）作为掩膜来实现，称为高斯金字塔。
- 归并 如果在任一层次上存在四个具有接近相同的一致性度量的相邻区域，则将其归并为金字塔中上一层中的单个区域，归并过程对应于金字塔的升采样（pyramid-up）过程，即自下而上（bottom-up）的采样，可以使用差分拉普拉斯算子作为掩膜来实现，称为拉普拉斯金字塔。拉普拉斯金字塔是高斯金字塔的逆过程，作用是重建高斯金字塔。
- 降采样 对第i层图像金字塔进行高斯内核卷积，将所有偶数行和列去除，得到的图像即为第i+1层的图像，显而易见，结果图像只有原图的四分之一，向下取样会逐渐丢失图像的信息，缩小图像。
- 升采样 将图像在每个方向扩大为原来的两倍，新增的行和列以0填充，使用先前同样的内核(乘以4)与放大后的图像卷积，获得 “新增像素”的近似值，得到的图像即为放大后的图像，但是与原来的图像相比会发觉比较模糊，因为在缩放的过程中已经丢失了一些信息，如果想在缩小和放大整个过程中减少信息的丢失，就需要与原图进行比较。
可以采用分割四叉树的数据结构描述这一过程，叶子结点表示一致区域，区域分裂与归并对应于分割四叉树的删除或建立部分，分割结束后，树的叶子节点数对应于分割后的区域数。步骤如下：
1. 初始化：定义一个划分为区域的初始分割、一致性准则、金字塔数据结构
2. 在金字塔数据结构中，如果任意一个区域不是一致的，将其分裂为四个子区域；如果具有相同父节点的四个区域具有一致性，则归并它们；如果没有区域可以分裂或归并，继续执行下一步骤
3. 如果任意两个邻接区域具有一致性，（即使它们在金字塔的不同层或没有相同的父节点），则归并它们
4. 如果必须删除小尺寸区域，则将小区域与其最相似的邻接区域归并
Def.[图像分割（image segment）] 将图像划分为与其中含有的真实世界的物体或区域有强相关性的组成部分。
- [完全分割] 分割结果是一组唯一对应于输入图像中物体的互不相交的区域集合
- [部分分割] 分割后的区域并不直接对应于物体，图像被划分一些同性质区域
若将图像数据看做地形表面，梯度图像的灰度值表示高度，则区域边缘对应着高的分水岭（watershed），低梯度值的区域内部对应着集水盆地。显然，集水盆地区域具有如下性质：同一集水盆地的所有像素与该盆地的最小高程（灰度）区域有一条像素的简单路径相连，沿着该路径的灰度是单调递减的。 分水岭分割算法的目标，就是找到集水盆地和分水岭脊线。
分水岭图像分割方法之一
- 从图像的每个像素，开始寻找到图像表面高度局部最小值的下降路径；
- 定义集水盆地为满足如下条件的像素集合：下降路径终止于相同的高度最小点。
分水岭图像分割方法之二（部分文献中称为泛洪填充算法，即floodfill）
- 从底部填充集水盆地，从而标识下降路径
- 设每个局部最小值处存在一个孔，将地形表面浸泡在水中
- 水从极小值位置开始，填充所有的集水盆地，如果两个集水盆地交汇，需要在交汇处修建一条高达最高高度的水坝，水坝表示分水岭脊线
步骤如下：
1. 对输入梯度图像计算亮度直方图，创建一张具有灰度值的像素指针表
2. 设填充过程达到k层，则每个灰度小于等于k的像素被分配了唯一的集水盆地标号
3. 通过指针表获取灰度为k+1的元素，如果像素邻域至少有一个已经具有标号l，将其放入FIFO队列等待处理
4. 计算测定的集水盆地的测地学影响区域，即邻近集水盆地却未标注的图像像素所在地，进行标号
5. 队列中的元素按照顺序处理，没有分配已有标号的像素代表新发现的集水盆地，用新的标号标注。
泛洪填充算法又称洪水填充算法，是在很多图形绘制软件中常用的填充算法，如同windows画图中的油漆桶功能。算法的原理很简单，就是从一个点开始附近像素点，填充成新的颜色，直到封闭区域内的所有像素点都被填充新颜色为止。泛洪填充实现最常见有四邻域像素填充法，八邻域像素填充法，基于扫描线的像素填充方法。根据实现又可以分为递归与非递归（基于栈）。
