《Mathlab数学图像处理与识别》的读书笔记

2017-02-02 fishedee 数学

1 概述

数学图像处理,是研究计算机如何处理图像的,在研读玩三大基础学科,高等数学,线性代数和概率论后,这本书就显得很简单的了。从中也能看到,其实计算机处理图像的方法其实还很初级。

另外,图像处理的很大方法与机器学习中也是相通的,PCA降维,神经网络与支持向量机等等。

2 点运算

点运算是在每个灰度点进行原地灰度变换的运算

2.1 灰度直方图

灰度直方图就是显示各灰度级数的图像

2.1.1 灰度直方图

假设一副图像的像素总和为N,灰度级数为L,灰度级数为g的像素总数为\(N_g\),那么灰度直方图的方程为

\[ P_g = N_g,g\in[0,L-1] \]

2.1.2 归一化灰度直方图

归一化的直方图就是让直方图的像素总数归一化为1,则归一化直方图的方程为:

\[ P_g = \frac {N_g}{N} = \frac {N_g}{\sum_{i=0}^{L-1} N_g},g \in [0,L-1] \]

2.2 灰度变换

灰度变换是将灰度图像增强对比度或降低对比度的算法

2.2.1 灰度的线性变换

灰度线性变换的方程为:

\[ D_B = f(D_A) = f_AD_A+f_B \]

参数\(D_A\)是输入图像的灰度,\(D_B\)是输出图像的灰度,\(f_A\)是线性函数的斜率,\(f_B\)是线性函数的截矩

可以看出,灰度线性变换是整体灰度图像的变换,当

  • \(f_A\)控制对比度和反色,\(\lvert f_A \rvert >1\)时,整体对比度增强,\(\lvert f_A \rvert <1\)时,整体对比度减少,\(\lvert f_A \rvert =1\)时,整体对比度不变。\(f_A<0\)时,反色处理,暗区域变亮,亮区域变暗
  • \(f_B\)控制整体的亮度,\(f_B>0\)时整体亮度上升,\(f_B<0\)时整体亮度降低。

特殊情况下,\(f_A=-1,f_B=255\)时,整体对比度不变,完全反色。

2.2.2 灰度的对数变换

灰度的对数变换为:

\[ D_B=c \log(1+D_A) \]

参数\(D_A\)是输入图像灰度,\(D_B\)是输出图像灰度,c是控制参数

可以看出,灰度对数变量是只对灰度较暗的区域进行亮度提升,而亮区域基本不变的变换。由于在灰度为255的区域,输出值仅仅为6多一点,所以一般有一个c值用来将灰度值映射到[0,255]的区域中。如果要改变灰度图像的弯曲程度,则需要改变对数的底,对数的底越小,灰度变换越弯,对数的底越大,灰度变换越平坦。

2.2.3 灰度的伽马变换

灰度的伽马变换为:

\[ D_B = (D_A+esp)^{\gamma} \]

其中\(D_A\)\(D_B\)都是[0,1]之间的灰度值,分别是输入像素与输出像素值,

可以看出,伽马变换是一种可以随意提升高灰度或低灰度对比度的方法

  • \(\gamma=1\)时,对比度不变,\(\gamma>1\)时,高灰度对比度增强,\(\gamma<1\)时,低灰度区域对比度增强
  • esp是整体亮度的偏移

2.2.4 灰度阀值变换

灰度阀值变换的方程为

\[ D_B=\begin{cases} 0,D_A<T\\ 255,D_A>T\\ \end{cases} \]

T为指定的阀值

我想,这是最简单的灰度变换了,没什么好说的。

2.2.5 分段线性变换

分段灰度变换是从线性变换与阀值变换组合而来的,线性变换只能整体对比度上升或下降,阀值能对部分进行操作。那么分段灰度变换的思路是,先将图像灰度按阀值划分到不同的灰度变换中,然后每个段用各自的灰度线性变换方法

\[ D_B = \begin{cases} \frac {y_1} {x_1} D_A,D_A <x_1\\ \frac {y_2-y_1}{x_2-x_1}(D_A-x_1)+y_1,x_1<=D_A<X_2\\ \frac {255-y_2}{255-x_2}(D_A-x_2)+y_2,D_A>x_2\\ \end{cases} \]

其中\(D_A,D_B\)分别是输入和输出的灰度值,而\([x_1,y_1]\)\([x_2,y_2]\)是两个灰度控制参数,按照需要,你可以设置三个或以上的灰度控制参数。

2.3 灰度自适应变换

2.3.1 直方图均衡化

直方图均衡化就是将集中在某个区域的灰度图像进行自适应映射到全区域的灰度上。

其变换方程为:

\[ D_B = f(D_A) = D_{max}\sum\limits_{i=0} ^{D_A} \frac {N_i}{N} \]

其中\(D_A,D_B\)分别是输入和输出的灰度值,\(D_{max}\)是灰度的最大值,\(N_i\)是灰度级数为i的像素个数,\(N\)是总的像素个数

证明:

灰度归一化直方图为

\[ P_g = \frac {N_g} {N} \]

每个灰度等级为g的频率为\(P_g\)看成是随机变量X的概率密度值域,则直方图均衡化的问题就是找出一个随机变量函数\(Y=f(x)\),使得随机变量Y为均衡分布。

为了简化问题,我们将随机变量X和随机变量Y都看成是[0,1]上的概率分布,设X的概率密度为g(x),Y的概率密度为h(y),x与y之间概率分布函数为Y=f(x)。

则按题意可知,Y的概率密度函数为\(h(y)=1\),按照概率密度函数变换公式,有

\[ h(y)=g(x)\frac {dx} {d(f(x))}\\ 1 = g(x)\frac {dx} {df(x)}\\ \frac {df(x)} {dx} = g(x)\\ \int_0^{x} \frac {df(x)} {dx}dx = \int_{0}^{x} g(x) dx\\ f(x) = \int_0^{x}g(x)dx \]

因此,我们找出了随机变量Y的概率密度函数公式,而实际情况中,X的概率密度公式就是\(P_g\),而且灰度分布图为离散值,则

\[ f(x) = \sum\limits_{g=0}^{x} P_g \\ = \sum\limits_{g=0}^{x} \frac {N_g} {N}\\ \]

然后,我们将灰度值映射回\([0,D_{max}]\)的灰度范围中,得到

\[ D_B = f(D_A) = D_{max}\sum\limits_{i=0} ^{D_A} \frac {N_i}{N} \]

所以得证

这就是直方图均衡化的全部了,里面包含了很重要的将归一化灰度图看成是概率分布,以及将离散型看成是连续型,然后再将连续型看回成离散型的思想,非常典型。

2.3.2 直方图规定化

直方图规定化,就是将指定图像的灰度特征迁移到另外一副图像中。其方法就是利用均衡化将图像的灰度特征与均衡化一一对应,然后做反变换得到的。

假设原图像的灰度为r,均衡化后灰度为s,带匹配图像的灰度灰度为z,均衡化后灰度为v。则我们要求的是r到将原图像的颜色r,变化匹配图像的灰度z的对应关系。

因为

\[ s = f(r) = \int_{0}^{r} p_r(\mu)d\mu\\ v = g(z) = \int_{0}^{z} p_z(\lambda)d\lambda \]

假设均衡化后图像颜色一致,即\(s=v\),那么

\[ \int_{0}^{r} p_r(\mu)d\mu =\int_{0}^{z} p_z(\lambda)d\lambda\\ f(r) = g(z)\\ z = g^{-1}[f(r)] \]

因此,我们的得到了颜色r对应到颜色z的关系。这个算法看起来比较神奇,可以将任意一张色系的图片套入到另外一张图片上。

3 几何变换

几何变换是将图像的坐标移动到另外一个坐标的变换,一般情况下,我们采用的方法不是从原坐标中逐个复制到新坐标去。而是遍历新坐标,找到对应的原坐标,然后取出原坐标的灰度复制到新坐标去。

3.1 线性变换

3.1.1 缩放

\((x_0,y_0)\)为原图像上的一点,图像水平缩放量为\(s_x\),垂直缩放量为\(s_y\),则缩放后的点坐标\((x_1,y_1)\)为:

\[ \begin{bmatrix} x_0 & y_0 & 1 \end{bmatrix}=\begin{bmatrix} x_1 & y_1 & 1 \end{bmatrix}\begin{bmatrix} \frac {1}{s_x} & 0 & 0\\ 0 & \frac {1} {s_y} & 0\\ 0 & 0 & 1\\ \end{bmatrix} \]

证明:

从题意可以得到,

\[ \begin{cases} x_1 = x_0 \cdot s_x\\ y_1 = y_0 \cdot s_y \\ \end{cases} \]

用矩阵表示为:

\[ \begin{bmatrix} x_1 & y_1 & 1 \end{bmatrix}=\begin{bmatrix} x_0 & y_0 & 1 \end{bmatrix}\begin{bmatrix} s_x & 0 & 0\\ 0 & s_y & 0\\ 0 & 0 & 1\\ \end{bmatrix} \]

将矩阵求逆矩阵为:

\[ \begin{bmatrix} x_0 & y_0 & 1 \end{bmatrix}=\begin{bmatrix} x_1 & y_1 & 1 \end{bmatrix}\begin{bmatrix} \frac {1}{s_x} & 0 & 0\\ 0 & \frac {1} {s_y} & 0\\ 0 & 0 & 1\\ \end{bmatrix} \]

所以得证

3.1.2 最近插值

\((x,y)\)为几何变换中求得的原坐标值,这个坐标值可能为小数,故不存在对应的原坐标值\((x',y')\),最近插值的算法就是简单地对这个小数四舍五入得到的,即

\[ f(x,y)=g(round(x),round(y)) \]

其中f为所求的坐标的灰度值,g为原来图像坐标的灰度值。简单暴力,也是最容易想到的。

3.1.3 双线性插值

双线性插值就是根据该点周边的四个点来线性插值这个点的灰度,其中\((x,y)\)为坐标中的小数部分,f(0,0),f(0,1),f(1,0),f(1,1)分别是这个坐标周边的四个点,则其方程为:

\[ f(x,y)=[f(1,0)-f(0,0)]x+[f(0,1)-f(0,0)]y+[f(1,1)+f(0,0)-f(0,1)-f(1,0)]xy+f(0,0) \]

证明:

当点落在\((0,0)与(1,0)\)之间的直线时,

\[ f(x,0)=f(0,0)+x[f(1,0)-f(0,0)] \]

当点落在\((0,1)与(1,1)\)之间的直线时,

\[ f(x,1)=f(0,1)+x[f(1,1)-f(0,1)] \]

那么,我们在点\((x,y)\)的点,相当于落在\((x,0)\)\((x,1)\)的直线上,那我们可以再用一次一元插值的方法

\[ f(x,y) = f(x,0)+y[f(x,1)-f(x,0)] \]

然后将\(f(x,1)\)\(f(x,0)\)代入就可以得到方程了,简单。

这个方法本质就是二元函数上的线性插值算法,它将一元上进行推广得到的,推导的方法也很巧妙,固定一方不动,改变另外一边就可以用一元插值了。

3.1.4 三次插值

高阶插值就是用根据该点邻近的4x4邻域内的像素取样得到的,设要计算的坐标为

\[ f(x,y)=f(i+u,j+v) \]

其中,i,j是坐标的整数部分,u,v是坐标的小数部分,则

\[ f(x,y)=f(i+u,j+v)=ABC\\ A = \begin{bmatrix} S(1+v)\\ S(v)\\ S(1-v)\\ S(2-v)\\ \end{bmatrix}^T\\ C = \begin{bmatrix} S(1+u)\\ S(u)\\ S(1-u)\\ S(2-u)\\ \end{bmatrix}\\ B=\begin{bmatrix} f(i-1,j-1)&f(i-1,j)&f(i-1,j+1)&f(i-1,j+2)\\ f(i,j-1)&f(i,j)&f(i,j+1)&f(i,j+2)\\ f(i+1,j-1)&f(i+1,j)&f(i+1,j+1)&f(i+1,j+2)\\ f(i+2,j-1)&f(i+2,j)&f(i+2,j+1)&f(i+2,j+2)\\ \end{bmatrix}\\ S(x)=\begin{cases} 1-2\lvert x \rvert^2+\lvert x \rvert^3,0<=\lvert x \rvert < 1\\ 4 - 8 \lvert x \rvert+5\lvert x \rvert^2-\lvert x \rvert^3,1<=\lvert x \rvert < 2\\ 0,\lvert x \rvert >= 2\\ \end{cases} \]

证明略

这是理论上比双线性插值更好的差值算法,只是计算量太大了。

3.2 仿射变换

3.2.1 平移

\((x_0,y_0)\)为原图像上的一点,图像水平偏移量为\(t_x\),垂直偏移量为\(t_y\),则平移后的点坐标\((x_1,y_1)\)为:

\[ \begin{bmatrix} x_0 & y_0 & 1 \end{bmatrix}=\begin{bmatrix} x_1 & y_1 & 1 \end{bmatrix}\begin{bmatrix} 1 & 0 & 0\\ 0 & 1 & 0\\ -t_x & -t_y & 1\\ \end{bmatrix} \]

证明:

根据题意得

\[ \begin{cases} x_1 = x_0 + t_x\\ y_1 = y_0 + t_y \end{cases} \]

用矩阵表示为:

\[ \begin{bmatrix} x_1 & y_1 & 1 \end{bmatrix}=\begin{bmatrix} x_0 & y_0 & 1 \end{bmatrix}\begin{bmatrix} 1 & 0 & 0\\ 0 & 1 & 0\\ t_x & t_y & 1\\ \end{bmatrix} \]

那么逆矩阵为:

\[ \begin{bmatrix} x_0 & y_0 & 1 \end{bmatrix}=\begin{bmatrix} x_1 & y_1 & 1 \end{bmatrix}\begin{bmatrix} 1 & 0 & 0\\ 0 & 1 & 0\\ -t_x & -t_y & 1\\ \end{bmatrix} \]

因此我们得到了图像平移公式,我们遍历一下新图像,将新图像的每一点\((x_1,y_1)\)代入公式得到原图像坐标\((x_0,y_0)\),然后将原图像坐标的颜色套入新图像上就可以了

3.2.2 镜像

\((x_0,y_0)\)为原图像上的一点,图像的宽度为Width,高度为Height,则水平镜像后的点坐标\((x_1,y_1)\)为:

\[ \begin{bmatrix} x_0 & y_0 & 1 \end{bmatrix}=\begin{bmatrix} x_1 & y_1 & 1 \end{bmatrix}\begin{bmatrix} -1 & 0 & 0\\ 0 & 1 & 0\\ Width & 0 & 1\\ \end{bmatrix} \]

垂直镜像后的点坐标\((x_1,y_1)\)为:

\[ \begin{bmatrix} x_0 & y_0 & 1 \end{bmatrix}=\begin{bmatrix} x_1 & y_1 & 1 \end{bmatrix}\begin{bmatrix} 1 & 0 & 0\\ 0 & -1 & 0\\ 0 & Height & 1\\ \end{bmatrix} \]

证明:

显然,水平镜像时,有:

\[ \begin{cases} x_1 = width - x_0\\ y_1 = y_0\\ \end{cases} \]

转化为矩阵形式,再求逆就好了,不难。

3.2.3 旋转

\((x_0,y_0)\)为原图像上的一点,绕原点逆时针\(\theta\)到点\((x_1,y_1)\),则

\[ \begin{bmatrix} x_0 & y_0 & 1 \end{bmatrix}=\begin{bmatrix} x_1 & y_1 & 1 \end{bmatrix}\begin{bmatrix} cos \theta & -sin\theta & 0\\ sin \theta & cos \theta & 0\\ 0 & 0 & 1\\ \end{bmatrix} \]

顺时针旋转\(\theta\)到点\((x_1,y_1)\),则

\[ \begin{bmatrix} x_0 & y_0 & 1 \end{bmatrix}=\begin{bmatrix} x_1 & y_1 & 1 \end{bmatrix}\begin{bmatrix} cos \theta & sin\theta & 0\\ -sin \theta & cos \theta & 0\\ 0 & 0 & 1\\ \end{bmatrix} \]

证明:

\(L=\lvert OP \rvert = \sqrt{x_0^2+y_0^2}\),则有

\[ sin \alpha= y_0 / L\\ cos \alpha = x_0 / L \\ \]

到达\(P_1\)点,有

\[ sin (\alpha +\theta) = y_1 / L \\ cos (\alpha +\theta) = x_1 / L \\ \]

联立上面方程,将\(\alpha\)消去,得到\(y_1,x_1\)\(x_0,y_0,\theta\)的关系为

\[ x_1 = cos \theta x_0 - sin \theta y_0\\ y_1 = cos \theta y_0 - sin \theta x_0\\ \]

将上面方程变为矩阵形式,然后进行矩阵求逆就可以求得

\[ \begin{bmatrix} x_0 & y_0 & 1 \end{bmatrix}=\begin{bmatrix} x_1 & y_1 & 1 \end{bmatrix}\begin{bmatrix} cos \theta & -sin\theta & 0\\ sin \theta & cos \theta & 0\\ 0 & 0 & 1\\ \end{bmatrix} \]

所以得证。

这个定理得出了绕原点旋转的公式,对于绕任意一点的方程,则只需要简单地,先平移,后旋转,再平移就能解决了,就不再啰嗦了。

4 空间域增强

空间域增强就是直接在图像的时域上进行滤波操作的行为。图像上从左到右,从上到下遍历每一个点,相当于随着时间的变化出现的每一个信号,所以我们也称图像为时域信号。在后面的频率域增强中,我们先将图像转化为频率域信号,然后在频率域进行操作的行为,我们称为频率域信号。

空间域增强与频率域增强气其实都是同一种操作在不同角度进行操作的结果而已。

滤波的操作,就是遍历原图像的每一个点,然后取该点的邻域NxN个点,将几个点的灰度和模板做加权叠加,输出的值作为新图像的灰度值。

滤波中的相关操作为:

\[ g(x,y)=\sum\limits_{s=-a}^{a}\sum\limits_{t=-b}^{b}w(s,t)f(x+s,y+t) \]

其中g是新图像的灰度值,f是原图像的灰度值,a和b是滤波模板的长宽,w是滤波模板各点的权值。

滤波中的卷积操作为:

\[ g(x,y)=\sum\limits_{s=-a}^{a}\sum\limits_{t=-b}^{b}w(-s,-t)f(x+s,y+t) \]

跟相关操作差不多,注意取滤波模板的值不同。频率域的公式大多都是和卷积有关,而不是相关,所以要注意区别。

4.1 图像平滑

4.1.1 平均平滑

设模板长度为\(2k+1\),则平均模板为:

\[ w=\frac {1}{(2k+1)^2}\begin{bmatrix} 1&1&\cdots&1\\ 1&1&\cdots&1\\ \cdots&\cdots&&\cdots\\ 1&1&\cdots&1\\ \end{bmatrix}_{(2k+1)x(2k+1)} \]

这说明了平均模板的长宽都是相等,并且必须都是奇数。

4.1.2 高斯平滑

设模板长度为\(2k+1\),则高斯模板的矩阵为:

\[ M(i,j)=\frac {1} {2\pi \sigma^2}exp(-\frac {(i-k-1)^2+(j-k-1)^2}{2\sigma^2})),i,j\in [-k,k] \]

特别地,一个3x3的高斯模板为:

\[ w = \frac {1} {16} \begin{bmatrix} 1 & 2 & 1\\ 2 & 4 &2 \\ 1 & 2 & 1\\ \end{bmatrix} \]

显然,高斯模板就是一个套用二维高斯分布函数作为权值分布的公式。当\(\sigma\)较大时,权值非常分散,倾向于变为平均模板,当\(\sigma\)较小时,权值非常集中,倾向于变为点运算。经验估计,一般情况,\(\sigma\)为0.5,3x3模板时\(\sigma\)为0.8。

4.1.3 自适应平滑

自适应平滑就是在融合了阀值和普通算法平滑的算法。它的思想是,先对该店的邻域进行平滑判断,如果认为该点有噪音,就进行平滑,否则就不进行平滑。

判断噪音的两个方法为:

  • 该邻域的灰度最大值与最小值之差大于某个阀值T,即\(max(R)-min(R)>T\)
  • 该邻域的灰度方差大于某个阀值T,即\(D(R)>T\)

执行平滑的两个算法

  • 平均平滑
  • 高斯平滑

对于噪音是局部和随机的图像而言,自适应平滑的效果非常好

4.2 中值滤波

4.2.1 中值滤波

中值滤波是统计排序的滤波器,它的算法就是先将邻域的灰度取出来,然后所有灰度值从小到大排序,然后取中位数作为滤波响应。

可以看出,中值滤波对于周期性的强烈噪音有很好的效果。

4.2.2 改进的中值滤波

改进的中值滤波加入了噪音的判断标准而已,当该点的灰度值是邻域内的极值时,就执行中值滤波,否则就不执行。

相当于自适应平滑中的噪音判断算法而已

4.3 图像锐化

4.3.1 水平垂直差分算子

获取图像中水平边缘的算子为:

\[ w_1 = \begin{bmatrix} -1 & 0 \\ 1 & 0 \\ \end{bmatrix} \]

垂直边缘的算子为:

\[ w_2 = \begin{bmatrix} -1 & 1 \\ 0 & 0 \\ \end{bmatrix} \]

水平和垂直方向边缘的算子为:

\[ w = \lvert w_1 \rvert + \lvert w_2 \rvert \]

证明:

将图像看成是连续的二元函数图,求图像的水平边缘,就是相当于求二元函数中各点在y方向的偏导数,这个偏导数越大,代表这个点为水平方向的边缘的可能性越大。而,在连续函数中,y方向的偏导数为:

\[ \frac {\partial f}{\partial y} = \lim\limits_{\Delta y \to 0} \frac {f(x,y+\Delta y)-f(x,y)}{\Delta y} \]

对应图像中的离散形式为:

\[ \frac {\partial f}{\partial y} = \frac {f(x,y+1)-f(x,y)} {1} = f(x,y+1)-f(x,y) \]

改为矩阵形式为:

\[ w_1 = \begin{bmatrix} -1 & 0 \\ 1 & 0 \\ \end{bmatrix} \]

所以得证。并且求得图像中,每个点可能为正可能为负,所以我们一般将图像的每个点进行绝对值的处理,只关注大小不关注方向,这样,求得图像中,水平边缘越明显的地方越白,越不明显的地方越黑。

同理,检测垂直方向的算子为:

\[ w_2 = \begin{bmatrix} -1 & 1 \\ 0 & 0 \\ \end{bmatrix} \]

要求出两个方向都变化较大的点,则只需要将两个求得的图像相加就可以了。

4.3.2 一阶微分的Robert算子

Robert算子,算\(+45^o\)边缘的算子为:

\[ w_1 = \begin{bmatrix} -1 & 0 \\ 0 & 1 \\ \end{bmatrix} \]

\(-45^o\)边缘的算子为:

\[ w_2 = \begin{bmatrix} 0 & -1 \\ 1 & 0 \\ \end{bmatrix} \]

证明方式跟水平垂直差分算子一样,就不啰嗦了,只是角度变为45度而已。

4.3.3 一阶微分的Prewitt算子

Prewitt算子,算垂直边缘的算子为:

\[ w_1 = \begin{bmatrix} -1 & -1 & -1 \\ 0 & 0 & 0 \\ 1 & 1 & 1\\ \end{bmatrix} \]

水平方向的算子为:

\[ w_2 = \begin{bmatrix} -1 & 0 & 1 \\ -1 & 0 & 1 \\ -1 & 0 & 1\\ \end{bmatrix} \]

显然是跟水平差分算子一样的算法,只是将边缘数值从中间方向推广到周边共三个方向的边缘叠加而已。同理,这个算子也可以推广到45度的边缘检测上。

\(+45^0\)边缘的算子为:

\[ w_3 = \begin{bmatrix} -1 & -1 & 0 \\ -1 & 0 & 1 \\ 0 & 1 & 1\\ \end{bmatrix} \]

\(-45^o\)边缘的算子为:

\[ w_4 = \begin{bmatrix} 0 & -1 & -1 \\ 1 & 0 & -1 \\ 1 & 1 & 0\\ \end{bmatrix} \]

4.3.4 一阶微分的Sobel算子

Sobel算子,算垂直边缘的算子为:

\[ w_1 = \begin{bmatrix} -1 & -2 & -1 \\ 0 & 0 & 0 \\ 1 & 2 & 1\\ \end{bmatrix} \]

水平方向的算子为:

\[ w_1 = \begin{bmatrix} -1 & 0 & 1 \\ -2 & 0 & 2 \\ -1 & 0 & 1\\ \end{bmatrix} \]

显然是跟Prewitt算子一样,只是让中间方向的权值提高到2,而其他保持不变。同理的,Sobel算子还有一个优化,就是将中间方向的权值变为\(\sqrt{2}\)

和Prewitt算子一样,Sobel算子也可以用来检测45度方向的边缘。

4.3.5 二阶微分的拉普拉斯算子

在拉普拉斯算子中,同时检测水平方向和垂直方向边缘的算子为:

\[ w_1 = \begin{bmatrix} 0 & 1 & 0\\ 1 & -4 & 1 \\ 0 & 1 & 0\\ \end{bmatrix} \]

证明:

在水平方向上,连续函数的二阶y偏导数为:

\[ \frac {\partial^2 f}{\partial y^2}= \lim\limits_{\Delta y \to 0} \frac {\frac {f(x,y+\Delta y)-f(x,y)}{\Delta y}-\frac {f(x,y)-f(x,y-\Delta y)}{\Delta y}}{\Delta y } \]

对应的离散形式为:

\[ \frac {\partial^2 f}{\partial y^2}= \frac {\frac {f(x,y+1)-f(x,y)}{1}-\frac {f(x,y)-f(x,y-1)}{1}}{1 }=f(x,y+1)+f(x,y-1)-2f(x,y) \]

所以,检测水平边缘的矩阵形式为:

\[ w_1 = \begin{bmatrix} 0 & 0 & 0\\ 1 & -2 & 1 \\ 0 & 0 & 0\\ \end{bmatrix} \]

同理,检测垂直边缘的矩阵形式为:

\[ w_1 = \begin{bmatrix} 0 & 1 & 0\\ 0 & -2 & 0 \\ 0 & 1 & 0\\ \end{bmatrix} \]

由于两个算子没有关联操作,故可以组合一起为:

\[ w_1 = \begin{bmatrix} 0 & 1 & 0\\ 1 & -4 & 1 \\ 0 & 1 & 0\\ \end{bmatrix} \]

所以得证,从算子可以看出,拉普拉斯算子即使在图像旋转90度后依然得出一样的结果,所以拉普拉斯算子具有各向同性的性质。

同理,我们可以推广出同时检测水平垂直和45度的拉普拉斯算子为:

\[ w_1 = \begin{bmatrix} 1 & 1 & 1\\ 1 & -8 & 1 \\ 1 & 1 & 1\\ \end{bmatrix} \]

更进一步的,我们根据高斯分布的性质优化拉普拉斯算子为:

\[ w_1 = \begin{bmatrix} 1 & 4 & 1\\ 4 & -20 & 4 \\ 1 & 4 & 1\\ \end{bmatrix} \]

因此,拉普拉斯的各向同性使得它在图像边缘检测中有更广泛的意义。

4.3.6 高斯-拉普拉斯(LoG)算子

高斯拉普拉斯算子为:

\[ h(r)=-[\frac {r^2-2\sigma^2}{\sigma^4}]e^{-\frac {r^2}{2\sigma^2}} \]

其中

\[ r^2 = x^2 + y^2 \]

证明:

拉普拉斯算子作为一个检测边缘不错的算子,其最大的问题是,检测边缘后会将噪音也放大了。因此更好的办法应该是,先将图像进行高斯平滑处理,然后再进行拉普拉斯算子锐化。

而高斯平滑的本质就是将高斯函数与原图像进行卷积处理,即

\[ f(x,y)\ast G(x,y) \]

而继续求拉普拉斯算子,相当于将得到的图像进行二次梯度辐值的结果,即

\[ \Delta^2(f(x,y)\ast G(x,y)) \]

根据数学推导,可知

\[ \Delta^2(f(x,y)\ast G(x,y))\\ = f(x,y) \ast \Delta^2[G(x,y)] \]

即我只需要对高斯分布进行求二次梯度辐值,然后跟原函数图像做卷积就可以了,即

\[ G(x,y)=-e^{-\frac {x^2+y^2}{2\sigma^2}}\\ \frac {\partial G(x,y)}{\partial x} =-\frac {x}{\sigma^2}e^{-\frac {x^2+y^2}{2\sigma^2}}\\ \frac {\partial^2 G(x,y)}{\partial x^2} =\frac {x^2-\sigma^2}{\sigma^4}e^{-\frac {x^2+y^2}{2\sigma^2}}\\ \frac {\partial^2 G(x,y)}{\partial y^2} =\frac {y^2-\sigma^2}{\sigma^4}e^{-\frac {x^2+y^2}{2\sigma^2}}\\ \Delta^2 G(x,y) =\frac {x^2+y^2-2\sigma^2}{\sigma^4}e^{-\frac {x^2+y^2}{2\sigma^2}}\\ \]

所以得证。离散形式的5x5的高斯-拉普拉斯算子为:

4.3.7 高提升滤波

以上的算法都是介绍了图像的边缘检测,进行图像锐化的过程,就是将图像锐化边缘明显的部分,而不改变边缘不明显的部分。则图像锐化的过程为:

  • 图像A经过以上算子获得边缘检测图像B,注意图像B不要进行绝对化操作
  • 将图像A与图像B根据比例混合起来
  • 混合后执行灰度自适应拉伸处理。

混合函数为

\[ g(i,j)=\begin{cases} Af(i,j)+h(i,j),h(i,j)的锐化算子中心系数>0\\ Af(i,j)-h(i,j),h(i,j)的锐化算子中心系数<0 \end{cases} \]

其中g(i,j)为混合后的图 像灰度,f(i,j)为原图像,h(i,j)为图像边缘检测后的结果,A为混合比例系数,A越大原图像的比例越高,反之越低。

相关文章