Digital Image Fundamentals

Posted on 21/09/2018

数字图像处理基础

对一幅连续图像$f$，对坐标值进行数字化称为取样，对幅值数字化称为量化。
数字图像的表示：通过一个二维阵列$f(x,y)$表示图像，该阵列有M行和N列，其中$(x,y)$是离散坐标。以下是图像的矩阵表示： $f(x,y)= \begin{bmatrix} f(0,0) & f(0,1) & \cdots & f(0,N-1) \\ f(1,0) & f(1,1) & \cdots & f(1,N-1) \\ \vdots & \vdots & \ddots & \vdots \\ f(M-1,0) & f(M-1,1) & \cdots &f(M-1,N-1) \\ \end{bmatrix}$
灰度图像（Gray scale）是每个像素只有一个采样颜色的图像。这类图像通常显示为从最暗黑色到最亮的白色的灰度。在数字图像处理领域，灰度图像与黑白图像有区别。
1.1 像素间的基本关系

1.1.1 相邻像素表示
4邻域：位于坐标$(x,y)$的一个像素$p$有4个水平和垂直的相邻像素，这个像素集称为p的4邻域，用$N_4(p)$ 表示。
$p$的4个对角像素表示：$p$的4个对角的相邻像素用$N_D(p)$ 表示，
$p$的4邻域$N_4(p)$ 和4个对角像素$N_D(p)$ 一起称为p的8邻域，用$N_8(p)$ 表示。
1.1.2 3种邻接类型
令$V$是用于定义邻接性的灰度值集合，例如：灰度值范围为0到255的邻接像素中，集合V可能是这256个值的任何一个子集。
1. 4邻接：如果$q$在 $N_4(p)$ 集合中，则具有V中数值的两个像素$p$和$q$是4邻接的。
2. 8邻接: 如果$q$在 $N_8(p)$ 集合中，则具有V中数值的两个像素$p$和$q$是4邻接的。
3. m邻接（混合邻接）:如果（1）：$q$在$N_4(p)$ 中，或者（2）$q$在$N_D(p)$ 中，且集合$N_4(p)\cap N_4(q)$ 没有$V$值的像素，则具有$V$值的像素$p$和$q$是$m$邻接的。

2. 图像增强

图像增强的方法分为两大类：空间域方法和频域方法。空间域指图像平面本身，这类方法是以对图像像素的直接处理为基础。频域处理技术是以修改图像的傅里叶变换为基础。
2.1 灰度增强
按一定变换关系逐点改变原来图像中每个像素的灰度值：$s=T(r)$ ####2.1.1 空间域图像增强常用的三类函数
线性函数
对数函数
幂次函数
分段线性变换（任意合并多种变换，缺点在于需要更多的用户输入） ####2.1.2 直方图处理
灰度级为[0，L-1]范围的数字图像的直方图是离散函数：$g(r_k)=n_k$，其中$r_k$代表灰度值为k，$n_k$则代表灰度值为k的像素个数。
通常以灰度值为k的像素个数$nk$除以图像总像素个数n，以得到均一化的直方图：$p(r_k)=n_k/n$，其中$p(r_k)$ 代表灰度k的概率。
一个均一化的直方图的所有部分之和为1。
直方图是多种空间域图像处理技术的基础。
若一幅图像的像素占有全部可能的灰度级并且分布均匀，则这样的图像有高对比度和多变的灰度色调。
#####2.1.2.1 直方图均衡化（Histogram Equalization）
通过调节图像的直方图来增强对比度。
对于输入图像的每个具有r值的像素产生一个输出灰度值s：
$s=T(r)$
$T(r)$ 在$0 \leq r \leq L - 1$上为单调递增函数
当$0 \leq r \leq L - 1$时，$0 \leq T(r) \leq L - 1$
图像处理中，特别重要的变换函数为：$T(r)=(L-1)\int_0^rp_r(w)dw$
$p_s(s)=p_r(r) \frac{dr}{ds} \implies p_s(s)=\frac{1}{L-1}$

#####2.1.2.2 直方图匹配（Histogram Matching）

直方图均衡化可以产生有均匀直方图的输出图像，但是对于某些应用，采用均匀直方图的基本增强并不是最好的方法，尤其是可以指定图像处理后具有的直方图形状。这种用于产生具有特定直方图形状的处理方法称为直方图匹配（Matching）或规定化（Specification）。
直方图均衡是直方图规定化的特例（均衡即目标直方图为均匀分布）。

####2.1.2.3局部直方图处理

以每个像素邻域的灰度分布为基础设计变换函数，变换函数有直方图均衡化、直方图匹配。

2.2 空间滤波

滤波一词借用于频率域处理，“滤波”是指接受（通过）或（拒绝）一定的频率成分。例如，通过低频率的滤波器称为“低通滤波器”。
空间滤波器 由（1）一个邻域（通常是一个较小的奇数像素值的矩形），（2）对该邻域所包围的图像像素进行预定义操作组成。滤波产生一个新像素，新像素的坐标等于邻域中心的坐标，像素的值是滤波操作的结果。滤波器的中心访问输入图像中的每个像素后，就生成类处理（滤波）后的图像。分为“线性空间滤波器”和“非线性滤波器”。
相关是滤波器模板移动过图像并计算每个位置乘积之和的操作。
卷积与相关的机理相似，但滤波器首先要旋转180度。卷积是线性系统的理论的基石。卷积的基本特性是函数与单位冲激卷积操作之后，得到的是该函数在冲激处的副本。
离散单位冲激 包含单个1而其余都是0的函数，一个函数与离散单位冲激相关，会在该冲激位置产生一个该函数的翻转版本。
空间域增强是指增强构成图像的像素,空间域处理可以由下式定义：g(x,y)=T[f(x,y)]。T操作最简单的形式是领域为1X1的尺度（即单个像素）。在这种情况下g仅仅依赖于f在（x,y)点的值，T操作称为 灰度级变换函数（又称强度映射），形式为：$s=T(r)$
变换函数的值通常存储在一个一维阵列中，通过查表可以得到从r到s的映射。

2.2.1 空间滤波（Spatial Filtering）基础

邻域处理工作是操作邻域图像的像素值，以及与邻域有相同维数的子图像的值。
掩模（mask）的大小(m*n)满足一下条件：$m=2a+1$,$n=2b+1$,其中a、b为非负数。
某个像素点的线性空间滤波处理的表达式为： $g(x,y)=\sum_{s=-a}^{a}\sum_{t=-b}^{b}w(s,t)f(x+s,y+t)$ 其中：g(x,y)为滤波处理后的灰度值，w(s,t)为掩模系数，f(x+s,y+s)为邻域像素的灰度值。
上述线性空间滤波（linear filtering）处理与频域（frequency domain)中的卷积（convolution）处理概念很相似，因此常常称为“掩模与图像的卷积”

2.2.2 平滑空间滤波器

平滑滤波器用于模糊处理和减小噪声。
平滑线性空间滤波器的输出是包含在滤波掩模邻域内像素的简单平均值。因此，这些滤波器也称为均值滤波器。
一幅$M\times N$ 的图像经过一个$m \times n$的加权均值滤波器滤波的过程如下： $g(x,y)= \frac {\displaystyle\sum_{s=-a}^{a}\sum_{t=-b}^{b}w(s,t)f(x+s,y+t)}{\displaystyle\sum_{s=-a}^{a}\sum_{t=-b}^{b}w(s,t)}$
一幅图像的完全滤波图像是对每个像素点执行上式得到。

2.3 频域滤波

2.4 中值滤波

统计排序滤波器 是一种非线性空间滤波器，这种滤波器的响应以滤波器包围的图像区域中所包含的像素排序为基础，然后使用统计排序结果决定的值替代中心像素的值。
中值滤波器 :一个数值集合的中值 $\xi$是这样的数值：数值集合中一半小于等于 $\xi$，另一半大于大于 $\xi$。中值滤波器对随机噪声提供类优秀的去噪能力。
中值滤波器是使用最广泛的统计排序滤波器，除此之外还有 最大值滤波器 和 最小值滤波器。
2.5 锐化空间滤波器
锐化处理的主要目的是突出灰度的过渡部分。锐化处理可由空间微分来实现，图像微分会增强边缘和其他突变（如噪声），削弱灰度变化缓慢的区域。
图像边缘处的灰度值来说，通常有两种突变形式：
- 边缘两边图像灰度差异较大，这就形成了灰度台阶。在台阶处，一阶微分和二阶微分的值都不为0.
- 边缘两边图像灰度变化不如台阶那么剧烈，会形成一个缓慢变换的灰度斜坡。在斜坡的起点和终点一阶微分和二阶微分的值都不为0，但是沿着斜坡一阶微分的值不为0，而二阶微分的值为0。
对于图像的边缘来说，通常会形成一个斜坡过度。一阶微分在斜坡处的值不为0，那么用其得到的边缘较粗；而二阶微分在斜坡处的值为0，但在斜坡两端值不为0，且值得符号不一样，这样二阶微分得到的是一个由0分开的一个像素宽的双边缘。也就说，二阶微分在增强图像细节方面比一阶微分好得多，并且在计算上也要比一阶微分方便。
图像处理中的一阶微分通常使用梯度的幅值来实现。对于图像f(x,y),f在坐标(x,y)处的梯度是一个列向量： $\nabla f = grad(f) = \left[ \begin{array}{c} g_x \\ g_y \end{array} \right] = \left[\begin{array}{c} \frac{\partial f} {\partial x} \\ \frac{\partial f} {\partial y} \end{array} \right]$
要得到图像的梯度图，有以下步骤：
1. 图像在x方向的梯度$g_x$
2. 图像在y方向的梯度$g_y$
3. $M(x,y)=∣g_x∣+∣g_y∣$
Sobel和Laplace算子的特征：
- 一阶微分算子Sobel得到的边缘较粗
- 二阶微分算子Laplace得到的边缘则较细，并且边缘是双边缘Lpalace算子对噪声比较敏感，得到的边缘图像上噪声较明显

3 频率域滤波

3.1 背景知识

3.1.1 基本概念

复数$C=R+jI$是复平面指教坐标系统中的点 $(R,I)$。
极坐标下表示复数： $C= C (cos\theta + j*sin\theta)$
欧拉公式（Euler’s formula）：$e^{j\theta}=cos\theta+j*sin\theta$
结合前2个公式：$C= C e^{j\theta}$
冲激：线性系统和傅里叶变换研究的核心是冲激及其取样特性。
3.1.2 傅里叶分析
傅里叶级数(Fourier series):法国数学家傅里叶指出任何周期函数都可以表示为不同频率的正弦和/或余弦之和的形式，每个正弦项和/或余弦项都乘以不同的系数（傅里叶级数）,具有以下形式： $f(t)=\sum_{n=-\infty}^{\infty}c_ne^{j \frac{2\pi n}{T}t}$ 上式可展开为正弦与余弦之和这一事实来自欧拉公式。
傅里叶变换（Fourier transform）:非周期函数（该曲线下的面积是有限的）也可以用正弦和/或余弦乘以加权函数的积分表示，这种情况下的公式称为傅里叶变换。
傅里叶变换公式$F(\mu)$： $F(\mu)=\int_{-\infty}^{\infty}f(t)e^{-j2\pi\mu t}dt$
给定$F(\mu)$，通过傅里叶反变换可以得到$f(t)$: $f(t)=\int_{-\infty}^{\infty}F(\mu)e^{j2\pi\mu t}d\mu$
$F(\mu)$ 和$f(t)$ 称为 傅里叶变换对，表明一个函数可以由其变换来恢复。
使用欧拉公式，可以把傅里叶变换函数表示为： $F(\mu)=\int_{-\infty}^{\infty}f(t)[cos(2\pi\mu t)-jsin(2\pi\mu t)]dt$

3.1.3 卷积（Convolution）

卷积指函数$f(t)$和$g(t)$通过运算得到第三个函数的一种数学算子，表征函数$f$与经过翻转和平移的函数$g$的乘积函数所围成曲边梯形的面积。术语卷积表示卷积的运算过程和结果函数。
连续变量$t$的的两个连续函数$f(t)$和$g(t)$卷积的积分形式： $f(t)\star g(t)=\int_{-\infty}^{\infty}f(\tau)g(t-\tau)d\tau \\ \qquad\qquad=\int_{-\infty}^{+\infty}g(\tau)f(t-\tau)d\tau$
卷积图示解释：
$t$从$-\infty$移动到$+\infty$过程中，计算$f$与$g$相交时的积分值。计算结果即为在位置$t$处，函数$f(\tau)$的加权平均值，加权函数为$g(-\tau)$。
如果$f(t)=\delta(t)$是一个单位冲激，上式的结果就是$g(t)$，称为冲激响应。
卷积定理（Convolution Theorem）
1. 空间域中两个函数的卷积的傅里叶变换，等于两个函数的傅里叶变换在频率域中的乘积；反过来，如果有两个变换的乘积，可以通过傅里叶反变换得到空间域的卷积。这一结果是卷积定理的一半，可以表示为（其中“ $\star$”表示为求两个函数的卷积） $f(t)\star g(t) \Leftrightarrow G(\mu)F(\mu)$
2. 卷积定理的另一半为： $f(t)g(t) \Leftrightarrow G(\mu)\star F(\mu)$ 其中$F(\mu)$和$G(\mu)$ 分别是$f(t)$、$g(t)$的傅里叶变换，“$\Leftrightarrow$”表示从左往右为傅里叶变换，从右往左为傅里叶逆变换。
卷积定理是频率域滤波的基础。
3.1.4 二维离散傅里叶变换及其反变换
二维离散傅里叶变换(DFT)： $F(u,v)=\sum_{x=0}^{M-1}\sum_{y=0}^{N-1}f(x,y)e^{-j2\pi(ux/M+vy/N)}$ 式中，$f(x,y)$是大小为$M \times N$的数字图像。
傅里叶反变换（IDFT): $f(x,y)=\frac{1}{MN}\sum_{u=0}^{M-1}\sum_{v=0}^{N-1}F(u,v)e^{j2\pi(ux/M+vy/N)}$ 以上两式构成二维离散傅里叶变换对。
3.1.4.1 二维离散傅里叶变换的性质
平移和旋转：对图像进行平移，不影响其傅里叶变换频谱。若图像f(x,y)旋转 $\theta$角度，则F(u,v)也旋转相同的角度，反之，若F(u,v)旋转一个角度，f(x,y)也旋转相同的角度。
3.1.4.2 二维卷积定理
二维循环卷积表达式： $f(x,y)\star g(x,y)=\sum_{m=0}^{M-1}\sum_{n=0}^{N-1}f(m,n)g(x-m,y-n)$
二维卷积定理： $f(x,y)\star g(x,y)\Leftrightarrow F(u,v)G(u,v)$ $f(x,y)g(x,y)\Leftrightarrow F(u,v) \star G(u,v)$

3.2 使用频率域滤波器平滑图像

图像中的边缘和其他急剧灰度变化（如噪声）主要影响其傅里叶变换的高频部分。因此，在频率域平滑（模糊）可通过衰减高频来达到，这就是低通滤波。
3.2.1 理想低通滤波器（ILPF)
以圆点为圆心、以$D_0$为半径的园内，无衰减的通过所有频率，而在该圆外阻断所有频率的二维低通滤波器，称为 理想低通滤波器(ILPF)。由以下函数确定： $G(u,v)= \begin{cases} 1, &D(u,v) \leq D_0 \\ 0, &D(u,v) > D_0 \end{cases}$ 式中，$D_0$是一个正常数，表示截止频率。$D(u,v)$是频率域中点(u,v)与频率矩形中心的距离，即： $D(u,v)=[(u-P/2)^2+(v-Q/2)^2]^{1/2}$
ILPF函数的的透视图、图像显示和径向剖面显示图：
使用理想低通滤波器可以产生明显的振铃状波纹, 因此应用中很少采用理想低通滤波器。
3.2.2 布特沃斯(Butterworth)低通滤波器(BLPF)
是一种平滑理想低通滤波器。布特沃斯滤波器有一个参数，称为滤波器的“阶数”。当阶数较高时，布特沃斯滤波器接近理想滤波器；对于较低的阶数，布特沃斯滤波器更像高斯滤波器。因此，布特沃斯滤波器可以视为两种极端滤波器的过渡。
n阶布特沃斯低通滤波器传递函数定义： $G(u,v)=\frac{1}{1+(\frac{D(u,v)}{D_0})^{2n}}$
BLPF函数的透视图、图像显示和径向剖面显示图：
3.2.3 高斯(Gauss)低通滤波器(GLPF)
高斯低通滤波器的定义为： $G(u,v)=e^{-D^2(u,v)/2\sigma^2}$ 式中 $\sigma$是关于中心的扩展度的度量，通过令 $\sigma=D_0$，可以表示为： $G(u,v)=e^{-D^2(u,v)/2D_0^2}$
GLPF函数的透视图、图像显示和径向剖面显示图：
高斯滤波后图像被平滑的程度取决于标准差。它的输出是邻域像素的加权平均，同时离中心越近的像素权重越高。因此，相对于均值滤波（mean filter）它的平滑效果更柔和，而且边缘保留的也更好。
OpenCV中的高斯滤波函数：
1. 复合形式

void cvSmooth(const CvArr* src, CvArr* dst,intsmoothtype=CV_GAUSSIAN,int param1=3, int param2=0,double param3=0, doubleparam4=0 );

单独的高斯滤波函数

void GaussianBlur(InputArray src, OutputArray dst, Size ksize, double sigmaX, double sigmaY=0, int borderType=BORDER_DEFAULT ) ;

3.3 使用频率域滤波锐化图像

在低通滤波中我们说明了通过衰减图像傅里叶变换的高频信号可以平滑图像。因为边缘和其他灰度急剧变化的区域与高频分量有关, 所以图像的锐化可以通过在频率域的高通滤波实现。
3.2.1 理想低通滤波器（IHPF)

$G(u,v)= \begin{cases} 0, &D(u,v) \leq D_0 \\ 1, &D(u,v) > D_0 \end{cases}$

3.2.2 布特沃斯(Butterworth)低通滤波器(BHPF)

$G(u,v)=\frac{1}{1+(\frac{D_0}{D(u,v)})^{2n}}$

3.2.3 高斯(Gauss)低通滤波器(GHPF)

$G(u,v)=1-e^{-D^2(u,v)/2D_0^2}$

3.4 选择性滤波器

3.4.1 带通/带阻滤波器

带通允许一定频率范围内的信号通过而阻止其他频率范围内的信号通过。
时域：描述数学函数或物理信号对时间的关系。例如一个信号的时域波形可以表达信号随着时间的变化。
频率域：自变量是频率,即横轴是频率,纵轴是该频率信号的幅度,也就是通常说的频谱图。频谱图描述了信号的频率结构及频率与该频率信号幅度的关系。
空间域：又称图像空间(image space)。由图像像元组成的空间。在图像空间中以长度(距离)为自变量直接对像元值进行处理称为空间域处理。时域
空间频率是指单位长度亮度（灰度）做周期性变化的次数。是图像中灰度变化剧烈程度的指标，可以理解为灰度在平面空间上的梯度。

The spectrum of a signal refers to the plot of the magnitudes and phases of different frequency components

The bandwidth of a signal is the spread of the frequency components with significant energy existing in a signal

Fourier series and Fourier transform are ways to find spectrums for periodic and aperiodic signals, respectively

高斯滤波（Gauss Filter）

将图像频域处理和时域处理相联系，作为低通滤波器使用，可以将低频能量（比如噪声）滤去，起到图像平滑作用。高斯模糊是低通滤波的一种，也就是滤波函数是低通高斯函数，但是高斯滤波是指用高斯函数作为滤波函数，至于是不是模糊，要看是高斯低通还是高斯高通，低通就是模糊，高通就是锐化。

对于 $(2k+1)\times (2k+1)$ 的模板，模板中各元素值的计算公式如下： $H_{i,j}=\frac{1}{2\pi\sigma^2}e^{-\frac{(i-k-1)^2+(j-k-1)^2}{2\sigma^2}}$

漫水填充（FloodFill）

漫水填充法是一种用特定的颜色填充联通区域，通过设置可连通像素的上下限以及连通方式来达到不同的填充效果的方法。漫水填充经常被用来标记或分离图像的一部分以便对其进行进一步处理或分析，也可以用来从输入图像获取掩码区域，掩码会加速处理过程，或只处理掩码指定的像素点，操作的结果总是某个连续的区域。
OpenCV提供了2个实现版本：

int floodFill(InputOutputArray image, Point seedPoint, Scalar newVal, Rect* rect=0, Scalar loDiff=Scalar(), Scalar upDiff=Scalar(), int flags=4 )

int floodFill(InputOutputArray image, InputOutputArray mask, Point seedPoint,Scalar newVal, Rect* rect=0, Scalar loDiff=Scalar(), Scalar upDiff=Scalar(), int flags=4 )

图像采样

下采样(subsampling)：对于一幅尺寸为$MN$的图像，对其进行s倍下采样，即得到$(M/s) * (N/s)$尺寸的得分辨率图像（s是M和N的公约数），如果是矩阵形式的图像，就是把原始图像$ss$窗口内的图像变成一个像素，这个像素点的值就是窗口内所有像素的均值。
上采样(upsampling)：图像放大几乎都是采用内插值方法，即在原有图像像素的基础上在像素点之间采用合适的插值算法插入新的元素。插值算法分为传统插值、基于边缘的插值和基于区域的插值等。

并查集

在计算机科学中，并查集是一种树型的数据结构，用于处理一些不交集（Disjoint Sets）的合并及查询问题。有一个联合-查找算法（union-find algorithm）定义了两个用于此数据结构的操作：
- Find：确定元素属于哪一个子集。它可以被用来确定两个元素是否属于同一子集。
- Union：将两个子集合并成同一个集合。
由于支持这两种操作，一个不相交集也常被称为联合-查找数据结构（union-find data structure）或合并-查找集合（merge-find set）。其他的重要方法，MakeSet，用于创建单元素集合。有了这些方法，许多经典的划分问题可以被解决。

线性移不变系统（Linear Spatially Invariant,LSI）

线性（Linear System） $T[\alpha_1x_1(n_1,n_2)+\alpha_2x_2(n_1,n_2)]=\alpha_1T[x_1(n_1,n_2)]+\alpha_2T[x_2(n_1,n_2)]$ 如果$\alpha_1=\alpha_2$，$x_1(n_1,n_2)=x_2(n_1,n_2)$： $T[\alpha_1x_1(n_1,n_2)+\alpha_2x_2(n_1,n_2)]=T[0]=0\label{linear1}$ 公式$\ref{linear1}$是线性系统的必要条件，可以用于判断一个系统是否为线性。

只有当输入为0时，线性系统的输出才为0？
移不变(Spatially Invariant System) $T[x(n_1,n_2)]=y(n_1,n_2)\\$

如果 $T[x_1(n_1-k_1,n_2-k_2)]=y(n_1-k_1,n_2-k_2)$ 则$T[]$为spatially invariant。