《数字图像处理》笔记
二、数字图像处理基础
2.1.数字图像采集
- 当前普遍采用的图像传感器：CCD、CMOS
- 共同点：利用感光二极管进行光电转换，获取图像并转换为数字数据
- 采样：空间坐标离散化
- 均匀采样：等间隔划分成$N$行$M$列的网格
- 量化：灰度值离散化
- 均匀量化：将灰度值划分为$L$（通常为$2^k$）个等间隔的区间
- 通常取$k=8,L=256$，灰度从小到大，像素从暗到亮
2.2.数字图像表示
- 通常采用$N\times M$大小的矩阵表示具有$N$行$M$列的数字图像
- 左上角坐标为$(0,0)$，每个坐标的函数值$f(x,y)$表示该点的灰度值
- 像素：每个坐标表示的单元格
2.3.分辨率与图像存储
- 空间分辨率：一般指数字图像的阵列大小$N\times M$
- 灰度分辨率：指数字图像的灰度级别$L$（通常为$2^k$）
- 灰度：图像中每个像素的亮度值，通常用0-255表示。灰度级别越多，图像越清晰。
- 图像存储：长×宽×位数×通道数/8
- 一幅空间分辨率为$N\times M$，灰度分辨率为$2^k$的RGB数字图像，需要$N\times M\times k\times 3 / 8 B$存储（$1B=8bit$）
2.4.像素间关系
- 相邻与邻域（以$p(x,y)$为例）
- 4邻域：上下左右四个像素，记为$N_4§$
- 对角邻域：上左、上右、下左、下右四个像素，记为$N_D§$
- 8邻域：九宫格，记为$N_8§$
- 像素的邻接性
- 4邻接：$p$与$N_4§$中的像素相邻
- 8邻接：$p$与$N_8§$中的像素相邻
- m邻接（混合邻接）：先进行4邻接，若存在四宫格中仅对角存在像素，则对角连接
- 像素的连通性
- 按照对应邻接性连接像素，形成连通块
- 像素的距离
- 欧几里得距离（欧氏距离，范数2）：$d(p_1,p_2)=\sqrt{(x_1-x_2)^2+(y_1-y_2)^2}$
- 曼哈顿距离（街区距离，范数1）：$d(p_1,p_2)=|x_1-x_2|+|y_1-y_2|$
- 棋盘距离（范数$\infty$）：$d(p_1,p_2)=\max(|x_1-x_2|,|y_1-y_2|)$
- $D_m$距离：m邻接后，$p_1$到$p_2$的最短路径长度
2.5. 习题
- 数字图像处理的五个经典应用领域：媒体通信、宇宙探索、气象预报、侦缉破案、考古
- 二值图像：具有两种灰度值的图像，但不一定是0和1。黑白图像一定是二值图像，但二值图像不一定是黑白图像。
- 视觉适应性：人眼对突然的亮度变化需要时间适应，亮适应性的时间比暗适应性短得多。
- 马赫带效应：人眼对两个颜色之间的边界，亮度感知受影响。较暗的色块靠近边界的一侧看起来更暗。
- 为什么在某些图像处理中，需要对图像的灰度进行对数运算？
- 大量实验表明人的视觉感知到的主观亮度与进入人眼的光强度成对数关系。如此处理可以达到更好的视觉效果。
三、数字图像的基本运算
3.1.矩阵统计量
1 | mean, stddev = cv2.meanStdDev(img) # 均值、标准差 |
3.2.灰度直方图
- $H§=[h(r_0),h(r_1),\cdots,h(r_{L-1})]$
- $h(r_i)=n_i$
- 统计图像中每个灰度级别的像素个数
- 用于分析图像的对比度、亮度等信息
1 | hist = cv2.calcHist([img], [0], None, [256], [0, 256]) # 计算直方图 |
- 对比度越小，直方图中像素分布越集中；对比度越大，直方图中像素分布越分散均匀
- 对比度：最亮与最暗像素之间的差异，可以表示为最大灰度和最小灰度的比值
- 亮度越低，直方图中像素分布越靠近左侧
灰度直方图的特征
- 所有函数值的总和是MN
- 只反映每种像素值出现的次数，不反映像素的空间分布
- 一幅图像的直方图是唯一的，但不同图像的直方图可以相同（由第二点决定）
归一化直方图
- $P(r_k)=\dfrac{h(r_k)}{MN}$
- 将直方图中的像素个数转换为频率
1 | hist = hist.ravel()/hist.sum() |
3.3.图像的几何运算
修改像素之间的空间关系，由两种基本操作组成：
- 空间变换：改变像素之间的空间关系
- 灰度级插值：确定空间变换后的像素灰度值
3.3.1.空间变换
3.3.1.1.平移变换
$$
\begin{cases}
x_1=x_0+\Delta x\\
y_1=y_0+\Delta y
\end{cases}
$$
3.3.1.2.旋转变换
记旋转前（角坐标系）：
$$
\begin{cases}
x_0=r\cos\theta\\
y_0=r\sin\theta
\end{cases}
$$
顺时针旋转$\alpha$角度：
$$
\begin{cases}
x_1=r\cos(\theta-\alpha)\\
y_1=r\sin(\theta-\alpha)
\end{cases}\\
=\begin{cases}
x_1=r\cos\theta\cos\alpha+r\sin\theta\sin\alpha\\
y_1=r\sin\theta\cos\alpha-r\cos\theta\sin\alpha
\end{cases}\\
=\begin{cases}
x_1=x_0\cos\alpha+y_0\sin\alpha\\
y_1=-x_0\sin\alpha+y_0\cos\alpha
\end{cases}
$$
3.3.1.3.镜像变换
- 水平镜像：$x_1=M-x_0,y_1=y_0$
- 垂直镜像：$x_1=x_0,y_1=N-y_0$
1 | img_flip = cv2.flip(img, 1) # 1水平镜像，0垂直镜像，-1水平+垂直镜像 |
3.3.1.4.缩放变换
$$
\begin{cases}
x_1=kx_0 \\
y_1=ky_0
\end{cases}
$$
3.3.1.5.转置变换
$$
\begin{cases}
x_1=y_0\\
y_1=x_0
\end{cases}
$$
3.3.2.灰度级插值
3.3.2.0.前向映射与反向映射
- 前向映射：原图像的一个像素映射到目标图像的多个像素（放大）
- 反向映射：原图像的多个像素映射到目标图像的一个像素（缩小）
3.3.2.1.最近邻插值/零阶内插
取前向映射变换所得位置，最近的的整数像素位置作为灰度值。
计算简单，容易产生锯齿。
3.3.2.2.双线性插值
记：原像素位置$(0,0),(0,1),(1,0),(1,1)$，反向映射变换所得位置$(x,y)$在它们中间。
根据像素相对四点的位置，求取算术平均值作为灰度值。
$$
f(x,0)=(1-x)f(0,0)+xf(1,0)\\
f(x,1)=(1-x)f(0,1)+xf(1,1)\\
f(x,y)=(1-y)f(x,0)+yf(x,1)
$$
四、空间域图像增强
图像增强技术可以分为两大类：空间域方法和频率域方法。
空间域方法以对图像的像素进行直接处理为基础，基于点运算的方法，在单个像素层面上进行。
统一表示：$g(x,y)=T[f(x,y)]$，其中算子T是一种点操作。
4.1.基于点运算的图像增强方法
4.1.1.对比度拉伸/分段线性变换
一种通过增强图像的对比度来改善图像质量的方法。
- 一种典型处理：压缩最暗和最亮部分，拉伸中间部分。
- 均匀变亮/均匀变暗：与原始灰度成线性关系。
4.1.2.窗切片
一种强制某区间的灰度值为给定值的方法。
- 区外灰度值不变
4.1.3.修剪变换
去除较暗和较亮部分，拉伸中间部分，特殊的对比度拉伸。
4.1.4.其他变换
- 对数变换：$s=c\log(1+r)$
- 幂次变换：$s=cr^\gamma$
- 离散化：$s=0\ if\ r<b\ else\ L-1$
4.2.基于直方图的图像增强方法
4.2.1.直方图均衡化
将一副已知图像变换成灰度具有均匀概率分布的图像，通过增大灰度值范围和对比度实现图像增强。
- 记原图像$f$的灰度级为$L$，对原图像$f$进行直方图均衡化得到的新图像为$g$。图像直方图均衡化的步骤如下：
- 统计原图像$f$中各个灰度级的像素个数，得到原图像的直方图$H(f)$。
- 给定 $16\times16$的图像$f$灰度级$L=8$
- 给定 $H(f)=[29,40,32,27,13,10,6,3]$
- 对原图像的直方图$H(f)$进行归一化处理，得到归一化直方图$P(f)$。
- 归一化直方图$P(f)=[0.18,0.25,0.20,0.17,0.08,0.06,0.04,0.02]$
- 计算归一化直方图$P(f)$的累积分布函数$C(f)$（前缀和）。
- 累积分布函数$C(f)=[0.18,0.43,0.63,0.80,0.88,0.94,0.98,1.00]$
- 计算新灰度级$t_i=round((L-1)c(i))$（四舍五入）
- $[1,3,4,6,6,7,7,7]$
- 将原图像$f$中的灰度级$i$替换为新灰度级$t_i$，得到直方图均衡化后的新图像$g$。
- 统计新图像$g$中各个灰度级的像素个数，得到新图像的直方图$H(g)$。
- $H(g)=[0,29,0,40,32,0,40,19]$
- 对新图像的直方图$H(g)$进行归一化处理，得到新图像的归一化直方图$P(g)$。
- $P(g)=[0.00,0.18,0.00,0.25,0.20,0.00,0.25,0.12]$
- 通过直方图均衡化，可以使得图像的灰度值分布更加均匀，从而增强图像的对比度。
1 | img_eq = cv2.equalizeHist(img) # 直方图均衡化 |
4.2.2.直方图规定化
直方图均衡化能显著的增强整个图像的对比度，但增强效果不易控制。通过直方图规定化，可以人为地改变直方图的形状。
- 图像直方图规定化的步骤如下：
- 统计原图像$f$中各个灰度级的像素个数，得到原图像的直方图$H(f)$。
- 给定$64\times64$的图像$f$灰度级$L=8$
- 给定 $H(f)=[790,1023,850,656,329,245,122,81]$
- 对原图像的直方图$H(f)$进行归一化处理，得到归一化直方图$P(f)$。
- $P(f)=[0.19,0.25,0.21,0.16,0.08,0.06,0.03,0.02]$
- 计算归一化直方图$P(f)$的累积分布函数$C(f)$（前缀和）。
- $C(f)=[0.19,0.44,0.65,0.81,0.89,0.95,0.98,1.00]$
- 对于给定的目标信息，求出目标累积分布函数$C_1(f_1)$。
- 假设给定$P_1(f_1)=[0,0,0,0.15,0.2,0.3,0.2,0.15]$
- $C_1(f_1)=[0,0,0,0.15,0.35,0.65,0.85,1.00]$
- 比较$C(f)$和$C(f_1)$，对于原图像$f$每个灰度级$i$，找到$|C(i)-C_1(j)|$最接近的$j$，得到映射关系。
- $C(f)=[0.19,0.44,0.65,0.81,0.89,0.95,0.98,1.00]$
- $C_1(f_1)=[0,0,0,0.15,0.35,0.65,0.85,1.00]$
- $T(i)=[3,4,5,6,6,7,7,7]$
- 将原图像$f$中的灰度级$i$替换为新灰度级$t_i$，得到直方图均衡化后的新图像$g$。
- 统计新图像$g$中各个灰度级的像素个数，得到新图像的直方图$H(g)$。
- $H(g)=[0,0,0,790,1023,850,985,448]$
- 对新图像的直方图$H(g)$进行归一化处理，得到新图像的归一化直方图$P(g)$。
- $P(g)=[0,0,0,0.19,0.25,0.21,0.24,0.11]$
1 | f = cv2.imread('lenna.jpg', cv2.IMREAD_GRAYSCALE) |
4.*.滤波
- 利用像素和其邻域像素的灰度值进行图像增强的方法。
- 空间滤波分为平滑滤波和锐化滤波。
- filter：滤波器，用于对图像进行滤波处理
- mask：滤波器的模板，用于计算滤波后的像素值
4.3.基于空间平滑滤波的图像增强方法
- 平滑滤波
- 目的：去除（细小）噪声、平滑图像
- 本质：对像素和其邻域像素的灰度值进行加权平均
4.3.1.线性滤波方法——邻域平均滤波
- mask：$h=\dfrac{1}{9}\begin{bmatrix}1&1&1 \\ 1&1&1 \\ 1&1&1\end{bmatrix}$
- 可以调整mask的大小和权重，实现不同的平滑效果。
- 参数越大，平滑效果越明显，图像越模糊。
- 一般中心像素权重最大，边缘像素权重较小。
- 可以根据Gaussian分布确定mask的权重，实现高斯滤波。
1 | ksize = (5, 5) # mask大小，奇数 |
4.3.2.非线性滤波方法——中值滤波
- 中值滤波窗口形状不定：矩形、圆形、十字形等
- 选取窗口内像素的中值作为中心像素的灰度值
- 优点：在减少边缘模糊方面比邻域平均效果好，对消除冲激噪声效果好
- 缺点：对消除高斯噪声效果差
1 | img_median = cv2.medianBlur(img, 5) # 中值滤波 |
4.4.基于空间锐化滤波的图像增强方法
- 锐化滤波
- 目的：增强图像的边缘、轮廓
- 要求：图像具有较高的信噪比
- 锐化处理通常和边缘检测结合使用
- 梯度
- 函数$z=f(x,y)$在点$(x,y)$处的梯度：$\nabla f(x,y)=\left[\dfrac{\partial f}{\partial x},\dfrac{\partial f}{\partial y}\right]^T$
- 梯度是一个向量，具有方向和模
- 梯度的方向$\theta$：$\theta=\arctan\left(\dfrac{\partial f}{\partial y}/\dfrac{\partial f}{\partial x}\right)$
- 梯度的模$G$：
- 按2-范数计算，对应欧几里得距离（欧氏距离）：$|\nabla f_{(2)}|=\sqrt{\left(\dfrac{\partial f}{\partial x}\right)^2+\left(\dfrac{\partial f}{\partial y}\right)^2}$
- 按1-范数计算，对应曼哈顿距离（街区距离）：$|\nabla f_{(1)}|=|\dfrac{\partial f}{\partial x}|+|\dfrac{\partial f}{\partial y}|$
- 按$\infty$-范数计算，对应棋盘距离：$|\nabla f_{(\infty)}|=\max(|\dfrac{\partial f}{\partial x}|,|\dfrac{\partial f}{\partial y}|)$
- 由于数字图像的离散性，偏导使用差分表示：
- 一阶水平：$\Delta_y f(x,y)=f(x,y+1)-f(x,y)$
- 一阶垂直：$\Delta_x f(x,y)=f(x+1,y)-f(x,y)$
- 二阶水平：$\Delta_{yy} f(x,y)=f(x,y+1)+f(x,y-1)-2f(x,y)$
- 二阶垂直：$\Delta_{xx} f(x,y)=f(x+1,y)+f(x-1,y)-2f(x,y)$
4.4.1.罗伯特算子（Roberts Operator）
- $G=|\Delta_x f(x,y)|+|\Delta_y f(x,y)|$
- $h=\begin{bmatrix}1&0\\0&-1\end{bmatrix}或\begin{bmatrix}0&1\\-1&0\end{bmatrix}$
- 优点：交叉差分法，简单，计算量小。
- 缺点：使用像素少，对噪声敏感。
1 | img_roberts = cv2.filter2D(img, -1, np.array([[1, 0], [0, -1]])) # 罗伯特算子 |
4.4.2.拉普拉斯算子（Laplacian Operator）
- $G=|\Delta_{xx} f(x,y)|+|\Delta_{yy} f(x,y)|$
- 常用mask：
- $h_{4,1}=\begin{bmatrix}0&1&0\\1&-4&1\\0&1&0\end{bmatrix}$
- $h_{4,2}=\dfrac{1}{2}\begin{bmatrix}1&0&1\\0&-4&0\\1&0&1\end{bmatrix}$
1 | img_laplacian = cv2.Laplacian(img, cv2.CV_64F, ksize=3) # 拉普拉斯算子，参数越大，线条越粗 |
4.4.3.普瑞维特算子（Prewitt Operator）
- 近似一阶导数
- 常用mask:
- $h_1 = \begin{bmatrix}1&1&1\\0&0&0\\-1&-1&-1\end{bmatrix}$
- 其余通过旋转得到
4.4.4.索贝尔算子（Sobel Operator）
- 近似一阶导数，往往用于估计水平和垂直方向的梯度
- 常用mask:
- $h_1 = \begin{bmatrix}-1&-2&-1\\0&0&0\\1&2&1\end{bmatrix}$
- $h_2 = \begin{bmatrix}-1&0&1\\-2&0&2\\-1&0&1\end{bmatrix}$
1 | img_sobelx = cv2.Sobel(img, cv2.CV_64F, 1, 0, ksize=3) # x方向梯度（竖直），ksize越大，成块的色块越大 |
五、频率域图像增强
- 三大类正交变换：正弦型变换、方波型变换、基于特征向量的变换
- 卷积：一种数学算子，施加于两个函数$f,g$，产生第三个函数。
- 定义：$(f*g)(t)=\int_{-\infty}^{+\infty}f(\tau)g(t-\tau)d\tau=\int_{-\infty}^{+\infty}f(t-\tau)g(\tau)d\tau$
- 频率域图像处理流程：
- $原图像f(x,y)\xrightarrow{傅里叶变换}F(u,v)\xrightarrow{H(u,v)滤波}G(u,v)\xrightarrow{傅里叶逆变换}结果图像g(x,y)$
5.1.离散傅里叶变换（Discrete Fourier Transform，DFT）
- 傅里叶变换将一个信号从时间域转换到频率域，反傅里叶变换将一个信号从频率域转换到时间域。
- 离散傅里叶变换是正弦型变换
5.1.1.一维连续傅里叶变换
- 一阶连续函数$f(x)$的傅里叶变换：$F(u)=\int_{-\infty}^{+\infty}f(x)e^{-j2\pi ux}dx$（j为虚数单位）
- 一阶连续函数$f(x)$的反傅里叶变换：$f(x)=\int_{-\infty}^{+\infty}F(u)e^{j2\pi ux}du$
- 欧拉公式：$e^{j\theta}=\cos\theta+j\sin\theta$
5.1.2.一维离散傅里叶变换
- 给定时间序列$x[n],n=0,1,2,\cdots,N-1$
- 离散傅里叶变换为：
- $X[k]=\sum\limits_{n=0}^{N-1}x[n]e^{-j2\pi kn/N}$
- 离散傅里叶逆变换为：
- $x[n]=\dfrac{1}{N}\sum\limits_{k=0}^{N-1}X[k]e^{j2\pi kn/N}$
- 离散傅里叶变换为：
- 以$N=4$为例，利用欧拉定理可知：
- $\begin{bmatrix}X_0\\X_1\\X_2\\X_3\end{bmatrix}=\begin{bmatrix}e^{-0\cdot i\pi/2}&e^{-0\cdot i\pi}&e^{-0\cdot i3\pi/2}&e^{-0\cdot i2\pi}\\e^{-0\cdot i\pi/2}&e^{-1\cdot i\pi}&e^{-2\cdot i3\pi/2}&e^{-3\cdot i2\pi}\\e^{-0\cdot i\pi/2}&e^{-2\cdot i\pi}&e^{-4\cdot i3\pi/2}&e^{-6\cdot i2\pi}\\e^{-0\cdot i\pi/2}&e^{-3\cdot i\pi}&e^{-6\cdot i3\pi/2}&e^{-9\cdot i2\pi}\end{bmatrix}\begin{bmatrix}x_0\\x_1\\x_2\\x_3\end{bmatrix}=\begin{bmatrix}1&1&1&1\\1&-i&-1&i\\1&-1&1&-1\\1&i&-1&-i\end{bmatrix}\begin{bmatrix}x_0\\x_1\\x_2\\x_3\end{bmatrix}$
5.1.3.二维离散傅里叶变换
- 对于给定的二维离散信号$f(x,y)$，$x\in[0,M-1],y\in[0,N-1]$
- 二维离散傅里叶变换：
- $F(u,v)=\sum\limits_{x=0}^{M-1}\sum\limits_{y=0}^{N-1}f(x,y)e^{-j2\pi(\dfrac{ux}{M}+\dfrac{vy}{N})}$
- 二维离散傅里叶逆变换：
- $f(x,y)=\dfrac{1}{\sqrt{MN}}\sum\limits_{u=0}^{M-1}\sum\limits_{v=0}^{N-1}F(u,v)exp({j2\pi(\dfrac{ux}{M}+\dfrac{vy}{N})})$
- 二维离散傅里叶变换：
- 可分离性：傅里叶变换中的指数项可分离成只含有x,u的一项和只含有y,v的一项的乘积
- $exp(-j2\pi\dfrac{ux+vy}{N})=exp(-j2\pi\dfrac{ux}{N})\cdot exp(-j2\pi\dfrac{vy}{N})$
- 即$F(u,v)=\sum\limits_{x=0}^{N-1}F(x,v)exp(-j2\pi\dfrac{ux}{N})$
- 结论：对每个x，对行进行一次一维DFT，对每个y，对列进行一次一维DFT，可以得到二维DFT
- 零频率分量与平均值
- $F(0,0)=\dfrac{1}{\sqrt{MN}}\sum\limits_{x=0}^{M-1}\sum\limits_{y=0}^{N-1}f(x,y)$
- 结论：$F(0,0)$与图像的平均值成正比，$k=\sqrt{MN}$
- 共轭对称性
- $F(u,v)=F^*(-u,-v)$
- 结论：对于实数图像，其频谱是共轭对称的
- 平移不变性
- 空间域平移：$f(x-x_0,y-y_0)\leftrightarrow F(u,v)exp(-j2\pi\dfrac{ux_0}{M})exp(-j2\pi\dfrac{vy_0}{N})$
- 频率域平移：$f(x,y)exp(j2\pi\dfrac{ux_0}{M})exp(j2\pi\dfrac{vy_0}{N})\leftrightarrow F(u-u_0,v-v_0)$
- 结论：在空间域中平移图像，频率域只发生相移，不发生幅度变化；在频率域中平移图像，空间域只发生相移，图像幅值不变。
- 周期性
- $|F(u,v)|=|F(u+M,v+N)|=|F(-u,-v)|=|F(M-u,N-v)|$
5.1.4.OpenCV
1 | img_dft = cv2.dft(np.float32(img), flags=cv2.DFT_COMPLEX_OUTPUT) # 二维离散傅里叶变换 |
5.2.图像的傅里叶频谱分析
5.2.1.频谱图像关于$(M/2,N/2)$对称
5.2.2.频谱中心化前后的频谱图像
- 中心化前：频谱图像的原点位于左上角。低频分量位于图像的四角，高频分量位于图像的中心。
- 中心化后：频谱图像的原点位于中心。低频分量位于图像的中心，高频分量位于图像的四角。
- 由于傅里叶变换的周期性，可以将中心化前的图像组成$2\times2$的周期图像，取中心部分即为中心化后的图像。
5.2.3.傅里叶变换的意义
- 简化计算：傅里叶变换将空间域图像转换到频率域，在空间域中处理图像时用到的复杂卷积运算，等同于在频率域中进行简单的乘法运算。
- 用频谱图表示的频率域图像中：
- 中心部位是能量集中的低频特征，反映图像的平滑部分
- 边缘部位是能量分散的高频特征，对应于细节、边缘、结构复杂区域、突变部位和噪声
- 因此，频率域中滤波的概念更直观、更易理解
- 某些特定的应用需求只能在频率域进行处理，如频率域图像特征提取、数据压缩、纹理分析等
5.3.基于频率域滤波的图像噪声消除——频率域低通滤波
- 图像的边缘和噪声对应于傅里叶频谱的高频成分。
- 低通滤波器：保留低频成分，抑制高频成分。
- 可以削弱噪声影响，但是会模糊边缘细节，降低图像清晰度。
- 与空间域平滑滤波效果类似。
5.3.1.理想低通滤波器
设置一个截止频率$D_0$，只保留频率小于$D_0$的成分。
- 理想低通滤波器：$H(u,v)=\begin{cases}1&D(u,v)\leq D_0\\0&D(u,v)>D_0\end{cases}$
- 其中，$D(u,v)$是频率域中点$(u,v)$到中心的距离。若变换被中心化，注意调整$D(u,v)$的计算方式。
- 理想低通滤波器无法通过电子器件实现。
5.3.2.巴特沃斯低通滤波器（Butterworth Lowpass Filter）
- n阶巴特沃斯低通滤波器：$H(u,v)=\dfrac{1}{1+[D(u,v)/D_0]^{2n}}$
- 巴特沃斯低通滤波器能够物理实现。
5.4.3.其他低通滤波器
- 梯形低通滤波器
- n指数型低通滤波器
5.4.基于频率域滤波的图像增强——频率域高通滤波
- 图像的边缘和细节对应于傅里叶频谱的高频成分。
- 高通滤波器：保留高频成分，抑制低频成分。
5.4.1.理想高通滤波器
与理想低通滤波器相反，只保留频率大于$D_0$的成分。
- 理想高通滤波器：$H(u,v)=\begin{cases}1&D(u,v)>D_0\\0&D(u,v)\leq D_0\end{cases}$
5.4.2.巴特沃斯高通滤波器（Butterworth Highpass Filter）
- n阶巴特沃斯高通滤波器：$H(u,v)=1-\dfrac{1}{1+[D_0/D(u,v)]^{2n}}$
5.4.3.高频增强滤波器
- 高通滤波会滤除低频成分，导致整体形象变暗。
- 高频增强滤波器：$H_e(u,v)=kH(u,v)+c$
- $k>1$：增强高频成分
- $c$：适当保留低频成分
5.5.带通滤波和带阻滤波
- 带通滤波：只保留频率在两个截止频率之间的成分。
- 带阻滤波：只保留频率在两个截止频率之外的成分。
5.6.习题
- 空间域图像和模板之间的运算是一种卷积运算。由卷积定理，空域的卷积可以通过频域中图像的傅立叶变换和模板的傅立叶变换间的乘积运算来实现。所以，对频域中的转移函数求傅立叶反变换即可得到用于空域卷积的模板。
七、图像压缩编码
7.1.离散余弦变换（Discrete Cosine Transform，DCT）
- DCT变换避免了复数运算。若图像矩阵是实数矩阵，则DCT变换后的矩阵也是实数矩阵。
- DCT是一种正交变换，变换矩阵是正交矩阵，变换核是可分离的。
- DCT有快速算法。
- DCT具有更强的信息集中能力。
7.1.1.一维离散余弦变换
- 给定长度为N的一维信号$x[n],n=0,1,2,\cdots,N-1$
- 一种常用的1D-DCT形式：$X[u]=a(u)\sum\limits_{m=0}^{N-1}x[m]cos(\dfrac{(2m+1)u\pi}{2N})$，其中$a(u)=\begin{cases}\sqrt{\dfrac{1}{N}}&u=0\\\sqrt{\dfrac{2}{N}}&u>0\end{cases}$
- 以$N=4$为例
- 1D-DCT：$\begin{bmatrix}X_0\\X_1\\X_2\\X_3\end{bmatrix}=\begin{bmatrix}0.5&0.5&0.5&0.5\\0.65&0.27&-0.27&-0.65\\0.5&-0.5&-0.5&0.5\\0.27&-0.65&0.65&-0.27\end{bmatrix}\begin{bmatrix}x_0\\x_1\\x_2\\x_3\end{bmatrix}$
- 1D-IDCT：$\begin{bmatrix}x_0\\x_1\\x_2\\x_3\end{bmatrix}=\begin{bmatrix}0.5&0.65&0.5&0.27\\0.5&0.27&-0.5&-0.65\\0.5&-0.27&-0.5&0.65\\0.5&-0.65&0.5&-0.27\end{bmatrix}\begin{bmatrix}X_0\\X_1\\X_2\\X_3\end{bmatrix}$
7.1.2.二维离散余弦变换
- 给定二维信号$x[m,n],m=0,1,2,\cdots,M-1,n=0,1,2,\cdots,N-1$
- 2D-DCT：
- $C(u,v)=a(u)a(v)\sum\limits_{m=0}^{M-1}\sum\limits_{n=0}^{N-1}f(x,y)cos(\dfrac{(2x+1)u\pi}{2M})cos(\dfrac{(2y+1)v\pi}{2N})$
- 2D-IDCT：
- $f(x,y)=\sum\limits_{u=0}^{M-1}\sum\limits_{v=0}^{N-1}a(u)a(v)C(u,v)cos(\dfrac{(2x+1)u\pi}{2M})cos(\dfrac{(2y+1)v\pi}{2N})$
- 其中，$a(u)=\begin{cases}\sqrt{\dfrac{1}{N}}&u=0\\\sqrt{\dfrac{2}{N}}&u>0\end{cases}$
- 2D-DCT：
7.1.3.OpenCV
1 | img_dct = cv2.dct(img.astype(np.float32)) # 二维离散余弦变换 |
7.2.数字图像压缩编码基础
7.2.1.信息冗余
- 图像压缩的原理
- 图像信号存在大量的冗余
- 人眼对图像的感知有限，可以以一定失真为代价换取数据量减少
- 压缩率：$压缩率=\dfrac{原始数据量}{压缩后数据量}$
- 信息冗余
- 空间冗余：相邻像素之间具有相关性
- 时间冗余：视频序列中相邻帧之间具有相关性
- 编码冗余：由编码方式导致的冗余
- 心理-视觉冗余：人眼对图像的感知有限，某些信息不重要
- 图像质量的主观评价：人的主观感受，如图像的清晰度、色彩、对比度等。
- 图像质量的客观评价：信噪比（SNR）、峰值信噪比（PSNR）、均方误差（MSE）等。
7.2.2.图像质量的客观评价
原始图像$f$，压缩后的图像$g$，图像大小为$M\times N$
7.2.2.1.均方误差（Mean Square Error，MSE）
- $MSE=\dfrac{1}{MN}\sum\limits_{x=0}^{M-1}\sum\limits_{y=0}^{N-1}[f(x,y)-g(x,y)]^2$
7.2.2.2.信噪比（Signal-to-Noise Ratio，SNR）
- $SNR=10\log_{10}\dfrac{\sigma^2}{MSE}$
- 其中$\sigma^2=\dfrac{1}{MN}\sum\limits_{x=0}^{M-1}\sum\limits_{y=0}^{N-1}f^2(x,y)$（原始图像的平方平均）
7.2.2.3.峰值信噪比（Peak Signal-to-Noise Ratio，PSNR）
- $PSNR=10\log_{10}\dfrac{max(I)^2}{MSE}$
- 通常$max(I)=L-1$
7.3.图像压缩编码
- 图像压缩编码处理流程：$f(x,y)\xrightarrow{信源编码器}\xrightarrow{信道编码器}信道传输\xrightarrow{信道解码器}\xrightarrow{信源解码器}g(x,y)$
- 信源编码：通过减少冗余来压缩数据的过程。
- 信源编码器模型：$f(x,y)\rightarrow映射变换器\rightarrow量化器\rightarrow符号编码器\rightarrow 信道$
- 信源解码器模型：$信道\rightarrow符号解码器\rightarrow反量化器\rightarrow反向映射变换器\rightarrow g(x,y)$
- 信道编码：也称差错控制编码，是一种在发送端给原数据添加与原数据相关的冗余信息，以便在接收端检测和纠正错误的编码方式。用于对抗传输过程中的噪声干扰。
- 汉明码、循环冗余校验码（CRC）
7.3.1.信息量与信息熵
- 事件提供的信息量：$I(x)=-\log_2P(x)$
- 信源：$X,x_i=p_i$，$x_i$发生的概率为$p_i$
- 信源熵：$H(X)=-\sum\limits_{i=1}^n p_i\log_2p_i$（比特/符号）
- 信源的平均信息量，即信源发出一个符号所携带的平均信息量
- 非负，在所有事件概率均等时，熵最大
7.3.2.编码
- 编码：由信源消息集到码字集的映射
- 根据长度是否相等分为：等长码和变长码
- 香农定理：无干扰下，平均码长的下限为信源熵
- 压缩比：$压缩比C=\dfrac{信源的平均比特率n}{编码后平均码长n_d}$
- 最大压缩比：$C_{max}=\dfrac{n}{H(X)}$
7.3.3.哈夫曼编码
- 带权路径长度（WPL）：$WPL=\sum\limits_{i=1}^n p_i l_i$
- 算法：
- 每次选取权重最小的两个节点，合并到一个新的节点，新节点的权重为两个节点的权重之和
- 直到所有节点合并到一起，最终得到一棵树
- 每个叶子节点代表一个源符号，从根节点到叶子节点的路径上的编码即为哈夫曼编码
7.3.4.算术编码
- 根据权重分配概率
- 符号序列映射为$[0,1)$的一个子区间
- 设$a:30%,b:70%$，则$a=[0,0.3),b=[0.3,1),aa=[0,0.09),ab=[0.09,0.3),ba=[0.3,0.51),bb=[0.51,1)$
7.4.位平面编码
- 位平面分解：一幅$m$比特表示灰度的图像，可以看作$m$个二值图像序列
- 在位平面中存在大量具有相同值的区域，可以利用连续的0和1进行高效的压缩编码和传输
7.4.1.格雷码分解编码
- 格雷码：相邻两个数的二进制码只有一位不同
- 设源图像中像素灰度值原$m$位为$x_{m-1}\cdots x_1x_0$，则有：$\begin{cases}g_i=x_i\oplus x_{i-1}\\g_{m-1}=x_{m-1}\end{cases}$
- 这种方式保证了相邻像素的格雷码值只有一位不同，在位平面上的连续0和1的序列更多，有利于压缩编码
7.5.游程编码/行程编码
- 适用于连续的0和1较多的序列
- 缺点：对于交替的0和1，效果不好甚至会增加编码长度
7.5.1.基本思想
- 对于连续的原符号，用一个符号和一个计数值来表示
- 如：$aaaaabbbcccc\rightarrow a5b3c4$
7.5.2.表0的游程编码
- 用一个$k$位二进制数表示相邻两个$1$之间$0$的个数
- 如：$0(14)10(9)110(20)10(30)110(11)1（括号内为连续0的数量）\rightarrow 1110(14),1001(9),0000(0),1111 0101(15+5=20),1111 1111(15+15=30),0000(0),1011(11),0000(0)（括号内为二进制表示数的含义）$
7.5.3.国际传真标准CCITT T.4(G3)
- 采用了霍夫曼编码和游程编码相结合的方法，对每种白长码字和黑长码字进行霍夫曼编码。
- 每一行总是以白长开始，允许为0；以唯一行尾（EOL）码字结束。
7.6.变换编码
九、图像分割
9.1.图像分割的基本概念
- 图像分割：依据图像的灰度、颜色、纹理、边缘等特征，把图像分成各自满足某种特征的连通区域的集合。
- 特征：灰度、颜色、纹理、形状、边缘轮廓等
9.2.基于边缘检测的图像分割——哈夫变换（Hough Transform）
- 边缘：具有不同灰度的区域的边界
- 常用算子：Roberts、Sobel、Prewitt、Laplacian、Canny
1 | edge = cv2.Canny(img, 80, 240) # Canny边缘检测，建议参数3是参数2的3倍 |
- 哈夫变换：用于在数字图像中检测直线、圆等形状特征
1 | lines = cv2.HoughLinesP(image, rho, theta, threshold, lines=None, minLineLength=None, maxLineGap=None) |
9.3.基于阈值的图像分割方法
提取物体与背景在灰度上的差异，将图像分为具有不同灰度级目标区域和背景区域
9.3.1.基于单一阈值的图像分割
- 找到一个介于目标和背景灰度之间的阈值$T$，将图像分为两部分
- $m(x,y)=\begin{cases}1&f(x,y)>T\\0&f(x,y)\leq T\end{cases}$
9.3.2. 一种阈值确定方法的步骤
- 选择一个初始估计值$T_0$
- 以$T_0$进行图像分割，计算两个像素集合各自的平均灰度值$M_1$和$M_2$
- 计算新的阈值$T_1=\dfrac{M_1+M_2}{2}$
9.3.3.自适应阈值
将图像分为$N\times N$个小区域，对每个小区域分别计算一个阈值
9.4.基于区域的图像分割方法
根据图像的灰度、颜色、纹理和图像像素统计特征的均匀性等图像的空间局部特征，将图像中的像素划归到各个物体或区域中。
9.4.1.区域生长法
- 根据事先定义的相似性准则，将具有相似性质的像素点合并成一个区域
- 步骤：
- 选择和确定一组能够正确代表所需区域的种子
- 原则：接近聚类中心的像素/红外图像目标检测中最亮的像素/按位置要求/根据经验
- 确定在生长过程中合并相邻像素的相似性原则
- 颜色（彩色图像）/灰度值差/和已知区域构成某种形状或尺寸
- 确定终止生长过程的条件或规则
- 一般：没有满足合并条件的像素
- 其他：区域大小、灰度差、形状等触发一定条件
- 从种子像素开始，逐步合并相邻像素，直到满足终止条件
- 选择和确定一组能够正确代表所需区域的种子
9.4.2.分裂合并法
- 从整幅图像开始，根据图像中区域间某种不一致性，逐步分裂成更小的区域；再根据相邻子区域的某种一致性准则，逐步合并。
- 图像四叉树表示法：
- 步骤：
- 对初始图像$R_0$进行一次四分裂得到四个子区域$R_1,R_2,R_3,R_4$
- 对于所有图像$R_i$进行检测，若满足分裂规则（$P(R_i)=TRUE$），则对该区域进行四分裂
- 对于所有图像$R_i$进行检测，若相邻区域满足合并规则（$M(R_i,R_j)=TRUE$），则合并这两个区域
- 重复步骤2和3，直到不再有区域满足分裂或合并规则
十二、形态学图像处理
12.1.二值形态学的基本运算
12.1.1.膨胀运算
$A\oplus B={x|(\hat{B})_y\cap A\neq\emptyset}$
计算过程：
创建一个与原图像$A$形状同的全0图像$G$。
对于原图像$A$中的每一个1所在的像素，将结构元素$B$的原点放在该位置，将$B$中的1对应的位置在$G$中置1。
- 图$G$中：白0表示原来是1新图像是0；灰1表示原来是1新图像还是1；黑2表示原来是0，新图像是1
1 | kernel = np.ones((5,1), np.uint8) # 5x1的结构元素 |
12.1.2.腐蚀运算
$A\ominus B={x|(\hat{B})_y\subseteq A}$
计算过程：
创建一个与原图像$A$形状同的全0图像$G$。
对于原图像$A$中的每一个1所在的像素，将结构元素$B$的原点放在该位置。
若$B$中的1对应的位置在$A$中都是1，则将$G$中该位置置1；
否则，$G$中该位置为0。
- 图$G$中：白0表示原来是1新图像是0；灰1表示原来是1新图像还是1
1 | kernel = np.ones((3,1), np.uint8) # 3x1的结构元素 |
12.1.3.开运算
- 先腐蚀后膨胀
- $A\circ B=(A\ominus B)\oplus B$
- 能够消除外部噪点，平滑轮廓
12.1.4.闭运算
- 先膨胀后腐蚀
- $A\bullet B=(A\oplus B)\ominus B$
- 能够消除内部噪点，保持物体的面积和形状
12.1.5.四种基本运算的性质
单调性、扩展性、交换性、结合性、平移不变性
12.2.二值图像的形态学处理
12.2.1.形态滤波
形态滤波器：先开后闭/先闭后开
12.2.2.边界提取
原图像-腐蚀图像=边界图像
12.2.3.区域填充
- 先对原图像$f$取反得到$f_c$
- 对$f_c$进行膨胀运算，得到$f_1$
- 对$f_1$和$f_c$求交集，得到$f_2$
- 对$f_2$和$f$求并集，得到目标图像$g$
12.2.4.骨架提取
通过迭代腐蚀和差分运算，得到图像的骨架
12.2.5.物体识别/击中击不中变换
利用具有特定形状的结构元素，对图像进行腐蚀运算，根据内部像素的变化判断物体的形状
零一、彩色图像与颜色模型
01.1.彩色图像的概述
- 颜色模型：也称作颜色空间，常见的包括：HSI、HSV、RGB、NTSC、YCbCr等
- 光谱分布和彩色感觉是多对一关系：为了得到某一种彩色感觉，可以用不同光谱分布的光以适当比例混合而成，参与混合的光的成分不唯一
- 三基色原理：合理地选取三种基本颜色，几乎自然界中所有颜色都可以用这三种颜色的适当混合来表示
- 彩色电视及图像处理系统的彩色显示器，都采用红绿蓝三色作为三基色
01.2.HSI
01.2.1.HSI
- HSI模型：色调（Hue）、饱和度（Saturation）、亮度/灰度（Intensity）
- 色调：反映颜色的类别，决定于彩色光的光谱成分，是“质”的特征
- 饱和度：反映颜色的纯度，即色彩的浓淡。纯色光谱色的含量越多，饱和度越高
- 亮度/灰度：决定了彩色光的强度
- 色调和饱和度合称“色度”
01.2.2.HSI模型
- HSI模型可以表示为一个纺锤体
- 色调$H$对应极坐标的$\rho$，饱和度$S$对应极坐标的$\theta$，亮度$I$对应垂直坐标$z$
- $\theta=0,120,240$对应红、绿、蓝三原色
- $I$越大，越靠近白色；$I$越小，越靠近黑色
01.2.3.HSI与RGB之间的转换
- 求色调$H$
- $\theta=\arccos(\dfrac{[(R-G)+(R-B)]}{2\sqrt{(R-G)^2+(R-B)(G-B)}})$
- $H=\begin{cases}\theta&G\geq B\\2\pi-\theta&G<B\end{cases}$
- 求饱和度$S$
- $S=1-\dfrac{3min(R,G,B)}{R+G+B}$
- 求亮度$I$
- $I=\dfrac{R+G+B}{3}$
01.3.HSV
01.3.1.HSV
- 与HSI类似，只是计算亮度的方式不同
- V：Value，表示颜色的亮度
01.3.2.HSV与RGB之间的转换
- $V=\dfrac{max(R,G,B)}{255}$
01.4.YUV
01.4.1.YUV
- 欧洲电视系统所采用的一种颜色编码
- Y：亮度信号；U、V：色差信号
01.4.2.YUV与RGB之间的转换
- $Y=0.299R+0.587G+0.114B$
- $U=-0.147R-0.289G+0.436B$
- $V=0.615R-0.515G-0.100B$
01.4.3.YUV的优点
- 采用YUV色彩空间，可以将亮度信号和色度信号分开处理，使得黑白电视机也能接收彩色信号
- 与RGB相比，最大的优点在于只需要占用极少的带宽，适合于传输和存储
01.5.YCbCr
01.5.1.YCbCr
- DVD、DV摄像机、数字电视等消费类视频产品中广泛使用的一种颜色编码
- Y：亮度信号；Cb：蓝色与参考值的差；Cr：红色与参考值的差
01.5.2.YCbCr与RGB之间的转换
$$
\begin{bmatrix} Y\\Cb\\Cr\end{bmatrix}=\begin{bmatrix}16\\128\\128\end{bmatrix}+\begin{bmatrix}0.257&0.504&0.098\\-0.148&-0.291&0.439\\0.439&-0.368&-0.071\end{bmatrix}\begin{bmatrix}R\\G\\B\end{bmatrix}
$$
01.5.3.YCbCr与YUV的对比
- YUV：适合于模拟电视系统
- YCbCr：适合于数字系统
01.6.彩色图像增强
- 单分量变化，亮度增强
- 单分量变化，饱和度增强
- 每个像素看作向量，以向量方式处理
零二、图像特征与理解
02.1.边界表示：链码
02.1.1.链码
- 用一组数字序列来表示边界的形状
- 以右向为0，顺时针方向为正方向，为每个方向编码。分为：4方向链码、8方向链码
- 从边界最左上角的像素开始，按照逆时针方向，将边界上的像素点用链码表示
- 链码是一个循环序列
- 上图的链码
- 4方向链码：$00333332322121110101$
- 8方向链码：$07666553321202$
02.1.2.链码归一化
对链码进行循环右移，直到链码的字典序最小
02.1.3.链码一阶差分
对链码进行模下差分，即$x_i=(x_i-x_{i-1}+n)\mod n$（n为链码的方向数）
02.2.二值图像的几何特征
02.2.1.位置
一般以物体面积的重心表示物体的位置
$$
\begin{cases}
x_c=\dfrac{1}{N}\sum\limits_{i=1}^N x_i\\\\
y_c=\dfrac{1}{N}\sum\limits_{i=1}^N y_i
\end{cases}
$$
02.2.2.方向
- 物体的方向：物体的主轴方向
- 确定主轴：（旋转）寻找最小外接矩形，主轴与矩形的长边平行
02.2.3.面积
物体的像素个数
02.2.4.周长
- 链码计算：$C=n_0+\sqrt{2}n_1$（$n_0$为偶数方向的个数，$n_1$为奇数方向的个数）
- 边界所占面积
- 边界的长度
02.2.5.形状
- 矩形度：$R=\dfrac{A}{A_r}$（$A$为物体面积，$A_r$为最小外接矩形的面积）
- 宽长比：$W=\dfrac{W_r}{L_r}$（$W_r$为最小外接矩形的宽度，$L_r$为最小外接矩形的长度）
- 圆形度
- 周长平方与面积的比值：$\dfrac{C^2}{A}$
- 面积与平均半径的比值：$\dfrac{A}{\bar{r}^2}$
- 偏心率
02.3.形状描述子
- 傅里叶描述子
- 边界链码
- 微分链码
02.4.矩描述
02.5.纹理分析：灰度共生矩阵
- 数组合数量
考试重点
- 3.2 灰度直方图
- 4.2 基于直方图的图像增强方法
- 直方图均衡化
- 4.3 基于空间平滑滤波的图像增强方法
- 4.4 基于空间锐化滤波的图像增强方法
- 5.1 2D-DFT
- 2D-DFT的性质
- 图像的傅里叶频谱图特性分析
- 7.1 DCT
- 7.3 图像压缩编码
- Huffman编码
- 7.4 位平面编码
- 9.3 基于阈值的图像分割
- 9.5 基于区域的图像分割
- 12.2 二值形态学基本运算
- 12.2 二值图像的形态学处理
cv2/np
1 | import cv2 |
plt
1 | from matplotlib import pyplot as plt |
评论