【数字图像处理】（四）空间域图像增强

原创

二进制人工智能 2021-06-21 15:34:25 ©著作权

©著作权归作者所有：来自51CTO博客作者二进制人工智能的原创作品，请联系作者获取转载授权，否则将追究法律责任

文章目录

空间域图像增强

空间域图像增强

图像增强是数字图像处理相对简单却最具艺术性的领域之一，增强的目的是消除噪声，显现那些被模糊了的细节或简单地突出一幅图像中读者感兴趣的特征。一个简单的例子是增强图像的对比度，使其看起来更加一目了然。应记住，增强是图像处理中非常主观的领域，它以怎样构成好的增强效果这种人的主观偏好为基础，也正是这一点为其赋予了艺术性。这与图像复原技术刚好相反，图像复原也是改进图像外貌的一个处理领域，但它是客观的。

1 图像增强基础

1.1 为什么要进行图像增强

图像增强是指根据特定的需要突出一幅图像中的某些信息，同时，削弱或去除某些不需要的信息的处理方法。其主要目的是使处理后的图像对某种特定的应用来说，比原始图像更适用。因此，这类处理是为了某种应用目的而去改善图像质量的。处理的结果是使图像更适合于人的观察或机器的识别系统。
应该明确的是增强处理并不能增强原始图像的信息，其结果只能增强对某种信息的辨别能力，而同时这种处理有可能损失一些其他信息。正因如此很难找到一个评价图像增强效果优劣的客观标准，也就没有特别通用模式化的图像增强方法，需要根据具体期望的处理效果做出取舍。

1.2 图像增强的分类

图像增强技术基本上可分成两大类：一类是空间域增强，另一类是频率域增强。

空间域图像增强与频率域图像增强不是两种截然不同的图像增强技术，实际上在相当程度上说它们是在不同的领域做同样的事情，是殊途同归的，只是有些滤波更适合在空间域完成，而有些则更适合在频率域中完成。

空间域图像增强技术主要包括直方图修正、灰度变换增强、图像平滑化以及图像锐化等。在增强过程中可以采用单一方法处理，但更多实际情况是需要采用几种方法联合处理，才能达到预期的增强效果（永远不要指望某个单一的图像处理方法可以解决全部问题）。

通过灰度变换改善图像外观的方法，以及直方图灰度修正技术（即直方图均衡化和直方图规定化）都是图像增强的有效手段，这些方法的共同点是变换是直接针对像素灰度值的，与该像素所处的邻域无关，而空间域增强则是基于图像中每一个小范围（邻域）内的像素进行灰度变换运算，某个点变换之后的灰度由该点邻域之内的那些点的灰度值共同决定，因此空间域增强也称为邻域运算或邻域滤波。空间域变换可使用下式描述：
【数字图像处理】（四）空间域图像增强_数字图像处理

2 空间域滤波

滤波是信号处理中的一个概念，是将信号中特定波段频率滤除的操作，在数字信号处理中通常通过傅里叶变换及其逆变换实现。由于下面的内容实际上和通过傅里叶变换实现的频域下的滤波是等效的，故而也称为滤波。不同的是空间域滤波主要直接基于邻域（空间域）对图像中像素执行计算。

2.1 空间域滤波和邻域处理

对图像中的每一点 $(x, y)$ ，重复下面的操作。
（1）对预先定义的以 $（ x, y ）$ 为中心的邻域内的像素进行运算。
（2）将（1）中运算的结果作为 $（ x, y ）$ 点新的响应。

上述过程就称为邻域处理或空间域滤波。一幅数字图像可以看成一个二维函数 $f (x, y)$ ，而 $x - y$ 平面表明了空间位置信息，称为空间域，基于 $x - y$ 空间邻域的滤波操作叫作空间域滤波。如果对于邻域中的像素计算为线性运算，则又称为线性空间域滤波，否则称为非线性空间域滤波。

图5.1直观地展示了用一个3× 3的模板（又称为滤波器、模板、掩模、核或者窗口）进行空间滤波的过程，模板为w，用黑笔圈出的是其中心。
【数字图像处理】（四）空间域图像增强_数字图像处理_02

滤波过程就是在图像 $f (x, y)$ 中逐点地移动模板，使模板中心和点 $(x, y)$ 重合。在每一点 $（ x, y ）$ 处，滤波器在该点的响应是根据模板的具体内容并通过预先定义的关系来计算。一般来说模板中的非0元素指出了邻域处理的范围，只有那些当模板中心与点 $（ x, y ）$ 重合时，图像 $f$ 中和模板中非0像素重合的像素参与了决定点 $（ x, y ）$ 像素值的操作。在线性空间滤波中模板的系数则给出了一种加权模式，即 $（ x, y ）$ 处的响应由模板系数与模板下面区域的相应 $f$ 的像素值的乘积之和给出。例如，对于图5.1而言，此刻对于模板的响应R为：
$R = w (- 1, - 1) f (x - 1, y - 1) + w (- 1, 0) f (x - 1, y) + \dots + w (0, 0) f (x, y) + \dots + w (1, 0) f (x + 1, y) + w (1, 1) f (x + 1, y + 1)$

更一般的情况，对于一个大小为 $m \times n$ 的模板，其中 $m = 2 a + 1, n = 2 b + 1$ , $a, b$ 均为正整数，即模板长与宽均为基数，且可能的最小尺寸为3× 3（偶数尺寸的模板由于其不具有对称性因而很少被使用，而1× 1大小的模板的操作不考虑邻域信息，退化为图像点运算），可以将滤波操作形式化地表示为：
【数字图像处理】（四）空间域图像增强_数字图像处理_03
对于大小为 $M \times N$ 的图像 $f (0, \dots, M - 1, 0, \dots, N - 1)$ ，对 $x = 0, 1, 2, \dots, M - 1$ 和 $y = 0, 1, 2, \dots, N - 1$ 依次应用公式，从而完成了对于图像f所有像素的处理，得到新的图像 $g$ 。

2.2 边界处理

执行滤波操作要注意的一点是当模板位于图像边缘时，会产生模板的某些元素很可能会位于图像之外的情况，这时，对于在边缘附近执行滤波操作需要单独处理，以避免引用到本不属于图像的无意义的值。
以下3种策略都可以用来解决边界问题。
（1）收缩处理范围——处理时忽略位于图像f边界附近会引起问题的那些点，如对于图5.1中所使用的模板，处理时忽略图像 $f$ 四周一圈1个像素宽的边界，即只处理从 $x = 1, 2, 3, \dots, M - 2$ 和 $y = 1, 2, 3, \dots N - 2$ （在MATLAB中应为 $x = 2, 3, 4, \dots, M - 1$ 和 $y = 2, 3, 4, \dots N - 1$ ）范围内的点，从而确保了滤波过程中模板始终不会超出图像f的边界。
【数字图像处理】（四）空间域图像增强_数字图像处理_04

2）使用常数填充图像——根据模板形状为图像f虚拟出边界，虚拟边界像素值为指定的常数，如0，得到虚拟图像 $f^{'}$ ，保证模板在移动过程中始终不会超出 $f^{'}$ 的边界。
【数字图像处理】（四）空间域图像增强_数字图像处理_05

（3）使用复制像素的方法填充图像——和（2）基本相同，只是用来填充虚拟边界像素值的不是固定的常数，而是复制图像f本身边界的模式。
【数字图像处理】（四）空间域图像增强_数字图像处理_06

2.3 相关和卷积

相关和卷积与卷积的略有不同，表示如下。
【数字图像处理】（四）空间域图像增强_数字图像处理_07
尽管差别细微，但有本质不同，相关和卷积卷积时模板是相对其中心点做镜像后再对 $f$ 位于模板下的子图像做加权和的，或者说在做加权和之前模板先要以其中心点为原点旋转180°。如果忽略了这一细微差别将导致完全错误的结果，只有当模板本身是关于中心点对称时，相关和卷积的结果才会相同。

2.4 滤波操作的MATLAB实现

MATLAB中与滤波相关的函数主要有imfilter()和fspecial()。Imfilter()完成滤波操作，而fspecial()可以为使用者创建一些预定义的二维滤波器，直接供imfilter()函数使用。

（1）滤波函数imfilter()的原型如下。
【数字图像处理】（四）空间域图像增强_数字图像处理_08
参数说明：

f是要进行滤波操作的图像；
w是滤波操作所使用的模板，为一个二维数组；
option1, option2, …是可选项，具体可以包括以下内容。

– 边界选项：主要针边界处理问题：

【数字图像处理】（四）空间域图像增强_数字图像处理_09
采用第一种方式固定值填充虚拟边界的问题是在边缘附近会产生梯度，采用后面3种方式填充可让边缘显得平滑。

– 尺寸选项：由于滤波中填充了边界，有必要指定输出图像g的大小，如：
【数字图像处理】（四）空间域图像增强_数字图像处理_10

模式选项：指明滤波过程是相关还是卷积：
【数字图像处理】（四）空间域图像增强_数字图像处理_11

返回值：

g为滤波后的输出图像。

（2）可创建预定义的二维滤波器的fspecial()函数的常见调用格式如下。

【数字图像处理】（四）空间域图像增强_数字图像处理_12
参数说明：

参数type指定了滤波器的类型：
可选输入parameters是和所选定的滤波器类型type相关的配置参数，如尺寸和标准差等，如果不提供则函数使用该类型的默认参数配置。

返回值：

返回值h为特定的滤波器。

下面结合以下几种代表性的情况具体说明。

h=fspecial（'average', hsize），返回一个大小为hsize的平均模板滤波器h。参数hsize可以是一个含有两个分量的向量，指明h的行和列的数目；也可以仅为一个正整数，此时对应于模板为方阵的情况。hsize的默认值为[3 3]。
h=fspecial（'gaussian', hsize, sigma），返回一个大小为hsize，标准差σ=sigma的高斯低通滤波器。hsize的默认值为[3 3],sigma的默认值为0.5。
h=fspecial（'sobel'），返回一个加强水平边缘的竖直梯度算子。（如果需要检测竖直边缘，则使用h'。）

【数字图像处理】（四）空间域图像增强_数字图像处理_14

例子：

f=imread('cameraman.tif');
w1=[1,1,1;
    1,1,1;
    1,1,1]/9;
w2=fspecial('average', 3);
g1=imfilter(f,w1,'corr','replicate');  % 均值滤波
g2=imfilter(f,w2,'corr','replicate');  % 均值滤波
subplot(1,3,1)
imshow(f)
title('原图像')
subplot(1,3,2)
imshow(g1)
title('w1')
subplot(1,3,3)
imshow(g2)
title('w2')

【数字图像处理】（四）空间域图像增强_数字图像处理_15

3 图像平滑

图像平滑是一种可以减少和抑制图像噪声的实用数字图像处理技术。在空间域中一般可以采用邻域平均来达到平滑的目的。

3.1 平均模板及其实现

【数字图像处理】（四）空间域图像增强_数字图像处理_16

从图可以看出滤波后的图 $g$ 有平滑或者说模糊的效果，这完全是模板 $w$ 作用的结果。2.4中的w提供了一种平均的加权模式，首先在以点 $（ x, y ）$ 为中心 $3 \times 3$ 邻域内的点都参与了决定在新图像 $g$ 中 $（ x, y ）$ 点像素值的运算；而且所有系数都为1表示它们在参与决定 $g （ x, y ）$ 值的过程中贡献（权重）都相同；最后前面的系数是要保证整个模板元素和为1，这里应为1/9，这样就能让新图像同原始图像保持在一个灰度范围中（如[0,255]）。这样的 $w$ 叫作平均模板，是用于图像平滑的模板中的一种，相当于一种局部平均。更一般的平均模板为：
【数字图像处理】（四）空间域图像增强_数字图像处理_17

3.1.1 工作原理

一般来说，图像具有局部连续性质，即相邻像素的数值相近，而噪声的存在使得在噪声点处产生灰度跳跃，但一般可以合理地假设偶尔出现的噪声影响并没有改变图像局部连续的性质，例如下面的局部图像f_sub，灰色底纹标识的为噪声点，在图像中表现为亮区中的2个暗点。
【数字图像处理】（四）空间域图像增强_数字图像处理_18
对f用3×3的平均模板进行平滑滤波后，得到的平滑后图像为g_sub。

【数字图像处理】（四）空间域图像增强_数字图像处理_19
【数字图像处理】（四）空间域图像增强_数字图像处理_20
显然，通过平滑滤波，原局部图像f_sub中噪声点的灰度值得到了有效修正，像这样将每一个点用周围点的平均替代，从而达到减少噪声影响的过程就称为平滑或模糊。

3.1.2 MATLAB实现

利用imfilter()和fspecial()，并以不同尺寸的平均模板实现平均平滑的MATLAB示例代码如下。

I=imread('coins.png');
I1=imnoise(I, 'gaussian', 0,0.02);     % 添加高斯噪声
h1=fspecial('average',3);  % 3x3平均模板
h2=fspecial('average',5);  % 5x5平均模板
h3=fspecial('average',7);  % 7x7平均模板
I2=imfilter(I1,h1,'corr','replicate');  % 相关滤波，重复填充边界
I3=imfilter(I1,h2,'corr','replicate');
I4=imfilter(I1,h3,'corr','replicate');
subplot(2,2,1)
imshow(I1)
title('噪声图像')
subplot(2,2,2)
imshow(I2)
title('3x3模板')
subplot(2,2,3)
imshow(I3)
title('5x5模板')
subplot(2,2,4)
imshow(I4)
title('7x7模板')

【数字图像处理】（四）空间域图像增强_数字图像处理_21

从图中可以看出随着模板的增大，滤波过程在平滑掉更多的噪声的同时也使得图像变得越来越模糊，这是由平均模板的工作机理所决定的。当模板增大到 $7 \times 7$ 时，图像中的某些细节，如相机结构已经难以辨识了。实际上，当图像细节与滤波器模板大小相近时就会受到比较大的影响，尤其当它们的灰度值又比较接近时，混合效应导致的图像模糊会更明显。随着模板的进一步增大，细节信息都会被当作噪声平滑掉。因此，在确定模板尺寸时应考虑好要滤除的噪声点的大小，有针对性地进行滤波。

3.2 高斯平滑及其实现

3.2.1 理论基础

平均平滑对于邻域内的像素一视同仁，为了减少平滑处理中的模糊，得到更自然的平滑效果，很自然地想到适当加大模板中心点的权重，随着距离中心点的距离增大，权重迅速减小，从而可以确保中心点看起来更接近于与它距离更近的点，基于这样的考虑得到的模板即为高斯模板。

常用的 $3 \times 3$ 的高斯模板如下所示：
【数字图像处理】（四）空间域图像增强_数字图像处理_22

高斯模板名字的由来是二维高斯函数，即读者熟悉的二维正态分布密度函数，回忆一下，一个均值为 $0$ 、方差为 $σ^2$ 的二维高斯函数如下。其3维示意图，如图5.4所示。
【数字图像处理】（四）空间域图像增强_数字图像处理_23
高斯模板正是将连续的二维高斯函数的离散化表示，因此任意大小的高斯模板都可以通过建立一个 $（ 2 k + 1 ） \times （ 2 k + 1 ）$ 的矩阵M得到，其 $（ i, j ）$ 位置的元素值可如下确定：
【数字图像处理】（四）空间域图像增强_数字图像处理_24

σ选择的小技巧
当标准差 $σ$ 取不同的值时，二维高斯函数的形状会有很大的变化，因而在实际应用中选择合适的σ值非常重要：如果 $σ$ 过小，偏离中心的所有像素权重将会非常小，相当于加权和响应基本不考虑邻域像素的作用，这样滤波操作退化为图像的点运算，无法起到平滑噪声的作用；相反如果σ过大，而邻域相对较小，这样在邻域内高斯模板将退化为平均模板；只有当σ取合适的值时才能得到一个像素值的较好估计。MATLAB中 $σ$ 的默认值为0.5，在实际应用中，通常对3× 3的模板取 $σ$ 为0.8左右，对于更大的模板可以适当增大 $σ$ 的值。

3.2.2 MATLAB实现

采用不同的 $σ$ 实现高斯平滑的MATLAB代码如下。

I=imread('cameraman.tif');
I1=imnoise(I, 'gaussian', 0,0.02);     % 添加高斯噪声
h3_5=fspecial('gaussian',3,0.5);  % sigma=0.5的3x3高斯模板
h3_8=fspecial('gaussian',5,0.8);   % sigma=0.8的3x3模板
h3_18=fspecial('gaussian',5,1.8);  % sigma=1.8的3x3模板，接近平均模板
h5_8=fspecial('gaussian',5,0.8);   % signma=0.8的5x5模板
h7_12=fspecial('gaussian',7,1.2);  % signma=1.2的7x7模板
I3_5=imfilter(I1,h3_5,'corr','replicate');  
I3_8=imfilter(I1,h3_8,'corr','replicate');
I3_18=imfilter(I1,h3_18,'corr','replicate');
I5_8=imfilter(I1,h5_8,'corr','replicate');
I7_12=imfilter(I1,h7_12,'corr','replicate');
subplot(2,3,1)
imshow(I1)
title('噪声图像')
subplot(2,3,2)
imshow(I3_5)
title('sigma=0.5的3x3高斯模板')
subplot(2,3,3)
imshow(I3_8)
title('sigma=0.8的3x3模板')
subplot(2,3,4)
imshow(I3_18)
title('sigma=1.8的3x3模板')
subplot(2,3,5)
imshow(I5_8)
title('signma=0.8的5x5模板')
subplot(2,3,6)
imshow(I7_12)
title('signma=1.2的7x7模板')

【数字图像处理】（四）空间域图像增强_数字图像处理_25

从结果图可以看到，当σ偏小时，平滑效果不明显；当σ增大至1.8时，高斯平滑效果类似于平均平滑效果。随着模板的增大，原图中的噪声得到了更好的抑制。还可以注意到模板大小越小，高斯滤波后的图像中图像细节被较好地保留。

上面介绍的平均平滑滤波器和高斯平滑滤波器都是线性平滑滤波器，在学习频率域滤波之后，还可以为它们赋予另外一个名字——低通滤波器。

3.3 自适应平滑滤波

利用平均模板的平滑在消除噪声的同时也使图像变得模糊，高斯平滑在一定程度上缓解了这些现象，但由平滑滤波机理可知这种模糊是不可避免的。这当然是使用者所不希望的，于是想到选择性的进行平滑，即只在噪声局部区域进行平滑，而在无噪声局部区域不进行平滑，将模糊的影响降到最少，这就是自适应滤波的思想。

那么怎样判断该局部区域是包含噪声的需要平滑的区域还是无明显噪声的不需平滑的区域呢？这要基于噪声的性质来考虑，噪声的存在会使得在噪声点处产生灰度跳跃，从而使噪声点局部区域灰度跨度较大。因此可以选择如下两个标准中的1个作为局部区域存在噪声的判据。

（1）局部区域最大值与最小值之差大于某一阈值 $T$ ，即： $m a x (R) - m i n (R) > T$ ，其中 $R$ 代表该局部区域。

（2）局部区域方差大于某一阈值 $T$ ，即： $D (R) > T, D (R)$ 表示区域 $R$ 中像素的方差。

自适应滤波算法的实现逻辑：

【数字图像处理】（四）空间域图像增强_数字图像处理_26

对于那些噪声位置具有随机性和局部性的图像，自适应的滤波具有非常好的效果。

4 中值滤波

中值滤波本质上是一种统计排序滤波器。对于原图像中某点 $(i, j)$ ，中值滤波以该点为中心的邻域内的所有像素的统计排序中值作为 $（ i, j ）$ 点的响应。
中值不同于均值，是指排序队列中位于中间位置的元素的值，例如：采用 $3 \times 3$ 中值滤波器，某点 $（ i, j ）$ 的8个邻域的一系列像素值为： $12 、 18 、 18 、 11 、 23 、 22 、 13 、 25 、 118$ 。统计排序结果为： $11 、 12 、 13 、 18 、 18 、 22 、 23 、 25 、 118$ 。排在中间位置（第5位）的18即作为 $（ i, j ）$ 点中值滤波的响应 $g (i, j)$ 。显然，中值滤波为非线性滤波器。

4.1 性能比较

中值滤波对于某些类型的随机噪声具有非常理想的降噪能力，对于线性平滑滤波而言，在处理的像素邻域之内包含噪声点时，噪声的存在总会或多或少的影响该点的像素值的计算，（对于高斯平滑，影响程度同噪声点到中心点的距离成正比），但在中值滤波中，噪声点则常常是直接被忽略掉的；而且同线性平滑滤波器相比，中值滤波在降噪的同时引起的模糊效应较低。中值滤波的一种典型应用是消除椒盐噪声。

4.1.1 图片加噪声

MATLAB中为图片加噪声的语句代码如下：
【数字图像处理】（四）空间域图像增强_数字图像处理_27
参数说明：

$I$ 为原图像；
可选参数type指定了噪声类型

返回值：
$J$ 为添加了噪声后的图像。

提示：
使用imnoise（'gaussian', m, v）添加高斯噪声时，相当于对原图像中每一个像素叠加一个从均值为m、方差为v的高斯分布中产生的随机样本值。当m=0时，较小的方差v通常保证了高斯分布在0附近的随机样本（高斯分布密度函数f(x)在x=0附近具有最大值）有一个较大的概率产生值，从而大部分的像素位置对原图像影响较小。

4.1.2 中值滤波的MATLAB实现

MATLAB提供了medfilt2()函数实现中值滤波，原型如下：
【数字图像处理】（四）空间域图像增强_数字图像处理_29
参数说明：

I1是原图矩阵；
m和n是中值滤波处理的模板大小，默认3×3。

返回值：
输出I2是中值滤波后的图像矩阵。

对比：

I=imread('coins.png');
J=imnoise(I, 'salt & pepper');     % 添加高斯噪声
w1=[1,2,1;
    2,4,2;
    1,2,1]/16;  % 高斯滤波3x3模板
w2=[1,1,1;
    1,1,1;
    1,1,1]/9;   % 平均平滑3x3滤波模板
I1=imfilter(J,w1,'corr','replicate');  
I2=imfilter(J,w2,'corr','replicate');
I3=medfilt2(J,[3,3]);   % 中值滤波
subplot(2,2,1)
imshow(J)
title('噪声图像')
subplot(2,2,2)
imshow(I1)
title('高斯滤波3x3模板')
subplot(2,2,3)
imshow(I2)
title('平均平滑3x3滤波模板')
subplot(2,2,4)
imshow(I3)
title('中值滤波')

【数字图像处理】（四）空间域图像增强_数字图像处理_30
程序的运行结果如图所示，从中可见线性平滑滤波在降噪的同时不可避免地造成了模糊，而中值滤波在有效抑制椒盐噪声的同时模糊效应明显低得多，因而对于椒盐噪声污染的图像，中值滤波要远远优于线性平滑滤波。

4.2 一种改进的中值滤波策略

中值滤波效果依赖于滤波窗口的大小，太大会使边缘模糊，太小则去噪效果不好。因为噪声点和边缘点同样是灰度变化较为剧烈的像素，普通中值滤波在改变噪声点灰度值的时候，会一定程度地改变边缘像素的灰度值。但是噪声点几乎都是邻域像素的极值，而边缘往往不是，因此可以利用这个特性来限制中值滤波。

具体的改进方法如下：逐行扫描图像，当处理每一个像素时，判断该像素是否是滤波窗口所覆盖下邻域像素的极大或者极小值，如果是，则采用正常的中值滤波处理该像素，如果不是，则不予处理。在实践中这种方法能够非常有效地去除突发噪声点，尤其是椒盐噪声，而几乎不影响边缘。

由于算法可以根据局部邻域的具体情况而自行选择执行不同的操作，因此改进的中值滤波也称为自适应中值滤波。

4.3 中值滤波的工作原理

与线性平滑滤波考虑邻域中每个像素的作用不同，中值滤波在每个 $n \times n$ 邻域内都会忽略掉那些相对于邻域内大部分其余像素更亮或更暗，并且所占区域小于像素总数一半 $n^2/2）$ 的像素的影响，而实际上满足这样条件被忽略掉的像素往往就是噪声。

5 图像锐化

图像锐化的目的是使模糊的图像变得更加清晰起来。

5.1 理论基础

图像锐化主要用于增强图像的灰度跳变部分，这一点与图像平滑对灰度跳变的抑制正好相反。事实上从平滑与锐化的两种运算算子上也能说明这一点，线性平滑都是基于对图像邻域的加权求和或者说积分运算的，而锐化则通过其逆运算导数（梯度）或者说有限差分来实现。

在平滑中要平滑的是噪声，希望平滑处理不要涉及边缘，而在锐化中要锐化的对象是边缘，希望处理不要涉及噪声。

5.2 基于一阶导数的图像增强——梯度算子

对于连续的二维函数 $f (x, y)$ ，其在点 $(x, y)$ 处的梯度是下面的二维列向量：
【数字图像处理】（四）空间域图像增强_数字图像处理_31
其中,
【数字图像处理】（四）空间域图像增强_数字图像处理_32

为在点（x, y）处f对x的偏导；
【数字图像处理】（四）空间域图像增强_数字图像处理_33
为在点 $（ x, y ）$ 处f对y的偏导。

梯度的方向就是函数 $f (x, y)$ 最大变化率的方向。梯度的幅值作为变化率大小的度量，其值为
【数字图像处理】（四）空间域图像增强_数字图像处理_34

对于离散的二维离散函数 $f （ i, j ）$ ，可以用有限差分作为梯度幅值的一个近似，如下式所示。
【数字图像处理】（四）空间域图像增强_数字图像处理_35
尽管梯度幅值和梯度两者之间有着本质的区别，但在数字图像处理中提到梯度时，往往不加区分，即将上式的梯度幅值称为梯度。

上式中包括平方和开方，不方便计算，因此可近似为绝对值的形式：
【数字图像处理】（四）空间域图像增强_数字图像处理_36

而在实际使用中，经常被采用的是另外一种近似梯度——Robert交叉梯度：
【数字图像处理】（四）空间域图像增强_数字图像处理_37

Robert交叉梯度
Robert交叉梯度对应的模板如下
【数字图像处理】（四）空间域图像增强_数字图像处理_38
其中，w1对接近45°边缘有较强响应；w2对接近-45°边缘有较强响应。

【例】基于Robert交叉梯度的图像锐化示例。
有了前面学习的滤波的知识，只要分别以 $w 1$ 和 $w 2$ 为模板，对原图像进行滤波就可得到 $G 1$ 和 $G 2$ ，最终的Robert交叉梯度图像为： $G = ∣ G 1 ∣ + ∣ G 2 ∣$ 。

在进行锐化滤波之前读者要将图像类型从uint8转换为double，这是因为锐化模板的负系数常常使得输出产生负值，如果采用无符号的uint8型，则负值会被截断。

在调用函数imfilter()时，还要注意不要使用默认的填充方式，因为MATLAB默认会在滤波时进行“0”填充，这会导致图像在边界处产生一个人为的灰度跳变，从而在梯度图像中产生高响应，而这些人为高响应值的存在将导致对图像中真正的边缘和其他使用者关心的细节的响应在输出梯度图像中被压缩在一个很窄的灰度范围，同时也影响显示的效果。本章这里采用了“replicate”的重复填充方式，也可采用“symmetric”的对称填充方式。

程序编程实现如下：

I=imread('coins.png');
I=double(I);      % 转换成double型，这样可以保存负值，否则uint8会把负值截掉。
w1=[-1,0;
    0,1];
w2=[0,-1;
    1,0];
G1=imfilter(I,w1,'corr','replicate');
G2=imfilter(I,w2,'corr','replicate');
G=abs(G1)+abs(G2);   % 计算Robert梯度
subplot(2,2,1)
imshow(I,[])
title('原图像')
subplot(2,2,2)
imshow(G,[])
title('Robert交叉梯度图像')
subplot(2,2,3)
imshow(G1,[])
title('w1滤波后取绝对值并重新标定')
subplot(2,2,4)
imshow(G2,[])
title('w2滤波后取绝对值并重新标定')

运行结果：
【数字图像处理】（四）空间域图像增强_数字图像处理_39

由于G1和G2中都可能有负值，（c）和图（d）分别是对G1和G2取绝对值后的图像，（c）中接近45°边缘较明显，而（d）中则突显出接近-45°方向的边缘，这与直接分析w1和w2模板结构得出的结论是一致的。

提示：
为便于观察效果，(a)、（b）、（c）、（d）都做了显示时的重新标定，即将图像的灰度范围线性变换到0～255之内并使得图像的最小灰度值为0，最大灰度值为255。在MATLAB中只需在用imshow()函数显示图像时加一个参数[ ]即可。

Sobel梯度

计算Sobel梯度的Sobel模板（奇数）：
【数字图像处理】（四）空间域图像增强_数字图像处理_40

w1对水平边缘有较大响应的竖直梯度，w2对竖直边缘有较大响应的水平梯度。

下面的MATLAB程序计算了一幅图像的竖直和水平梯度，它们的和可以作为完整的Sobel梯度：

clc;clear
I=imread('rice.png');
I=double(I);      % 转换成double型，这样可以保存负值，否则uint8会把负值截掉。
w1=[1,2,1;
    0,0,0;
    -1,-2,-1];
w2=w1';
G1=imfilter(I,w1); % 水平Sobel梯度
G2=imfilter(I,w2); % 竖直Sobel梯度
G=abs(G1)+abs(G2);   % 计算Sobel梯度
subplot(2,2,1)
imshow(I,[])
title('(a)原图像')
subplot(2,2,2)
imshow(G,[])
title('(b)Soble梯度图像')
subplot(2,2,3)
imshow(G1,[])
title('(c)w1滤波后取绝对值并重新标定')
subplot(2,2,4)
imshow(G2,[])
title('(d)w2滤波后取绝对值并重新标定')

【数字图像处理】（四）空间域图像增强_数字图像处理_41

（c）中接近水平方向的边缘较明显，（d）中接近竖直方向的边缘较明显

还可以直接利用MATLAB梯度函数gradient()计算Sobel梯度：

I=imread('rice.png');
I=double(I);      % 转换成double型，这样可以保存负值，否则uint8会把负值截掉。
[Gx,Gy]=gradient(I) ; % 计算x,y方向梯度
G=abs(Gx)+abs(Gy);   % 计算整体梯度

5.3 基于二阶微分的图像增强——拉普拉斯算子

下面介绍一种对于图像锐化而言应用更为广泛的基于二阶微分的拉普拉斯(Laplacian)算子。

理论基础

二维函数 $f (x, y)$ 的二阶微分（拉普拉斯算子）定义为：
【数字图像处理】（四）空间域图像增强_数字图像处理_42
对于离散的二维图像f(i, j)，可以用下式作为对二阶偏微分的近似：
【数字图像处理】（四）空间域图像增强_数字图像处理_43
将上面两式相加就得到用于图像锐化的拉普拉斯算子：
【数字图像处理】（四）空间域图像增强_数字图像处理_44
对应的滤波模板如下：
【数字图像处理】（四）空间域图像增强_数字图像处理_45
因为在锐化增强中，绝对值相同的正值和负值实际上表示相同的响应，故也等同于使用如下的模板W2：
【数字图像处理】（四）空间域图像增强_数字图像处理_46
分析拉普拉斯模板的结构，可知这种模板对于90°的旋转是各向同性的。所谓对于某角度各向同性是指把原图像旋转该角度后再进行滤波与先对原图像滤波再旋转该角度的结果相同。这说明拉普拉斯算子对于接近水平和接近竖直方向的边缘都有很好的增强，从而也就避免读者在使用梯度算子时要进行两次滤波的麻烦。更进一步，读者还可以得到如下对于45°旋转各向同性的滤波器：
【数字图像处理】（四）空间域图像增强_数字图像处理_47

沿用高斯平滑模板的思想，根据到中心点的距离给模板周边的点赋予不同的权重，还可得到如下的模板W5
【数字图像处理】（四）空间域图像增强_数字图像处理_48
MATLAB实现

clc;clear
I=imread('riceblurred.png');
% h1=fspecial('average',7);  % 3x3平均模板
% I=imfilter(I,h1,'corr','replicate');  % 相关滤波，重复填充边界
I=double(I);      
w1=[0,-1,0;
    -1,4,1;
    0,-1,0];
w4=[-1,-1,-1;
    -1,8,-1;
    -1,-1,-1];
w5=[1,4,1;
    4,-20,4;
    1,4,1];
I1=imfilter(I,w1,'corr','replicate');
I2=imfilter(I,w4,'corr','replicate');
I3=imfilter(I,w5,'corr','replicate');
subplot(2,2,1)
imshow(I,[])
title('(a)原图像')
subplot(2,2,2)
imshow(abs(I1),[])
title('(b)使用W1模板拉普拉斯锐化')
subplot(2,2,3)
imshow(abs(I2),[])
title('(c)使用W4模板拉普拉斯锐化')
subplot(2,2,4)
imshow(abs(I3),[])
title('(d)使用W5模板拉普拉斯锐化')

【数字图像处理】（四）空间域图像增强_数字图像处理_49

拉普拉斯锐化效果与之前的Robert与Sobel梯度锐化明显不同的一点是输出图像中的双边缘。拉普拉斯锐化对一些离散点有较强的响应，当然由于噪声也是离散点，因此这个性质有时是读者所不希望的。

5.4 高提升滤波

高提升滤波的原理
无论是基于一阶微分的Robert、Sobel模板还是基于二阶微分的拉普拉斯模板，其中各系数和均为0。这说明算子在灰度恒定区域的响应为0，即在锐化处理后的图像中，原图像的平滑区域近乎于黑色，而原图中所有的边缘、细节和灰度跳变点都作为黑背景中的高灰度部分突出显示。在基于锐化的图像增强中常常希望在增强边缘和细节的同时仍然保留原图像中的信息，而不是将平滑区域的灰度信息丢失。因此可以把原图像加上锐化后的图像得到比较理想的结果。

需要注意具有正的中心系数和具有负的中心系数的模板之间的区别，对于中心系数为负的模板(如w1,w3, w5)，要达到上述的增强效果，显然应当让原图像f(i, j)减去锐化算子直接处理后的图像，即：
【数字图像处理】（四）空间域图像增强_数字图像处理_50

其中：Sharpen(.)表示通用的锐化算子。

由于锐化后边缘和细节处的高灰度值的存在，经灰度伸缩后（归一化在[0,255]），原图灰度被压缩在一个很窄的范围内，整体上显得较暗。为了改善这种情况，对上面介绍的方法进行推广，具体的说就是在复合f(i, j)和Sharpen(f(i, j))时适当地提高f(i, j)的比重，形式化地描述如下。
【数字图像处理】（四）空间域图像增强_数字图像处理_51
MATLAB实现

clc;clear
I=imread('riceblurred.png');
% h1=fspecial('average',7);  % 3x3平均模板
% I=imfilter(I,h1,'corr','replicate');  % 相关滤波，重复填充边界
I=double(I);      
w1=[0,-1,0;
    -1,4,1;
    0,-1,0];
w4=[-1,-1,-1;
    -1,8,-1;
    -1,-1,-1];
w5=[1,4,1;
    4,-20,4;
    1,4,1];
I1=imfilter(I,w1,'corr','replicate');
I2=imfilter(I,w4,'corr','replicate');
I3=imfilter(I,w5,'corr','replicate');
subplot(2,2,1)
imshow(I,[])
title('(a)原图像')
subplot(2,2,2)
I1=I+abs(I1);
% 归一化到[0,255]
I1=255*((I1-min(I1(:)))/(max(I1(:))-min(I1(:))));
imshow(I1,[])
title('(b)使用W1模板拉普拉斯锐化')
subplot(2,2,3)
I2=6*I+abs(I2);
% 归一化到[0,255]
I2=255*((I2-min(I2(:)))/(max(I2(:))-min(I2(:))));
imshow(I2,[])
title('(c)使用W4模板拉普拉斯锐化')
subplot(2,2,4)
I3=6*I+abs(I3);
% 归一化到[0,255]
I3=255*((I3-min(I3(:)))/(max(I3(:))-min(I3(:))));
imshow(I3,[])
title('(d)使用W5模板拉普拉斯锐化')

【数字图像处理】（四）空间域图像增强_数字图像处理_52

一般来说权重系数A应为一个大于等于1的实数，A越大原图像所占比重越大，锐化效果越不明显。

6 高斯-拉普拉斯变换（Laplacianof a Gaussian, LoG）

6.1 原理

锐化在增强边缘和细节的同时往往也“增强”了噪声，因此如何区分开噪声和边缘是锐化中要解决的一个核心问题。
基于二阶微分的拉普拉斯算子对于细节（细线和孤立点）能产生更强的响应，并且各向同性。然而，它对于噪声点的响应也更强，经拉普拉斯锐化后噪声更明显。
考虑高斯型函数：
【数字图像处理】（四）空间域图像增强_数字图像处理_53
其中： $r^2=x^2+y^2$ , $σ$ 为标准差。
图像经该函数滤波将产生平滑效应，且平滑的程度由 $σ$ 决定。进一步计算h的拉普拉斯算子（ $h$ 关于 $r$ 求二阶导数），从而得到著名的高斯-拉普拉斯算子（Laplacianof a Gaussian, LoG）：
【数字图像处理】（四）空间域图像增强_数字图像处理_54

LoG函数的三维形状：
【数字图像处理】（四）空间域图像增强_数字图像处理_55

6.2 MATLAB实现

I=imread('rice.png');
I=imnoise(I,'salt & pepper',0.01);
Id=double(I);
% 拉普拉斯算子
h_lap=[-1,-1,-1;
    -1,8,-1;
    -1,-1,-1];
I_lap=imfilter(Id,h_lap,'corr','replicate'); %Laplacian锐化
I_lap=uint8(abs(I_lap));% 取绝对值并将255以上的响应截断
h_log1=fspecial('log',5,0.5); % 大小为5，sigma=0.5的LoG算子
I_log1=imfilter(Id,h_log1,'corr','replicate');
I_log1=uint8(abs(I_log1));
h_log2=fspecial('log',5,2);% 大小为5，sigma=2的LoG算子
I_log2=imfilter(Id,h_log2,'corr','replicate');
I_log2=uint8(abs(I_log2));
subplot(2,2,1)
imshow(I)
title('原图像')
subplot(2,2,2)
imshow(I_lap,[])
title('Laplacian锐化')
subplot(2,2,3)
imshow(I_log1,[])
title('大小为5，sigma=0.5的LoG算子')
subplot(2,2,4)
imshow(I_log2,[])
title('大小为5，sigma=2的LoG算子')