Matlab图像处理基础(2)：区域处理，边沿检测

目录

0. 概要

2. 卷积和相关

2.1 卷积

2.2 相关

2.3 卷积与相关的关系 

2.4 Matlab函数 

2.5 2-D卷积/相关的分解

3. 高通滤波，edge detection

3.1 matlab edge()函数

3.2 各种算子简介

3.2.1 Robert算子

 3.2.2 Prewitt算子

3.2.3 Sobel算子 

3.4 Canny算子

3.5 Marr-Hildreth

3.6 Matlab代码示例

4. 低通滤波，平滑

---

0. 概要

         基于Matlab的图像处理基础介绍。Matlab提供了丰富的图像处理相关的工具箱，基于Matlab代码实验可以使得图像处理算法的理解变得更加简单直观。

        以下所介绍的一些matlab工具函数都限于篇幅只涉及其基本的用法，更丰富的参数设定等可以查阅matlab说明。

        本篇是本系列第二篇，介绍区域处理（Regional Processing），包括卷积、相关、高通滤波（边沿检测）、低通滤波（平滑）等。

        上一篇参见：Matlab图像处理基础（part 1）

        区域处理是指基于一个邻域内的像素值计算一个输出像素值。各种形状和大小的邻域可以使用，最简单的如基于周围8个相邻位置的像素的值计算一个输出像素值。计算方法则有比如说，对邻域内的像素值进行加权平均、从邻域内的像素值中根据某种规则选择某一个作为输出，等等。。。

       

2. 卷积和相关

2.1 卷积

        区域处理要达成的目的主要有两个：平滑或尖锐化（smoothing the data and sharpening it）。其中都涉及到卷积（convolution）处理。卷积处理的定义如下：

                

        其中，t[k,l]表示template（模板），也被称为卷积核（kernel）, d[i,j]表示输入数据.

        以上处理可以这样理解：将模板放置在输入图像上，并在图像上横向和纵向移动。在每个位置上，取卷积核与图像重叠部分的数据进行逐点相乘并累加得到输出值。

        当然，在图像边缘处（此时卷积和伸出图像之外了）的处理需要特别的注意，比如说要对输入图像采用padding处理等。

        滤波模板通常采用正方形的模板（m = n）。

2.2 相关

        与卷积相似的一个运算被称为相关（correlation），相关运算定义如下：

                

2.3 卷积与相关的关系 

        很显然，将卷积运算的模板旋转180度后做相关运算得到的结果与原卷积结果相同，反之亦然。而实际应用中使用的模板通常都是对称的（旋转180度后与原模板相同）或者反对称的（旋转180度后等于原模板乘以负一）。如果模板是对称的，则基于相同模板的卷积运算和相关运算的结果相等；如果模板是反对称的，则基于相同模板的卷积运算和相关运算的结果相反。 

        对于卷积和相关的效果的一个简单而直观的理解是，它给出了给定的核与图像之间的相似度的度量，尤其是对于特定的核，可以明显地看出这种效果。

2.4 Matlab函数 

        Matlab用于计算卷积或相关处理的通用函数为如下所示： 

B = imfilter(A,h)

B = imfilter(A,h,options,...)

        A: 输入图像

        B: 输出图像

        h: 即滤波器模板、滤波器核

        options: 控制滤波运算的选项，包括三个方面的选项：边界填充；输入大小；滤波类型，如下表所示：

例1：

clear; close all; clc
 
% Display a grayscale image 
A = imread('cameraman.tif'); % 'cameraman.tif'是matlab工具箱内置的样本图片
h = [1,0,-1;2,0,-2;1,0,-1];
B_same = imfilter(A,h);        % 缺省为corr
B_full = imfilter(A,h,'full'); % 缺省为corr
B_conv = imfilter(A,h,'conv');
figure;
subplot(1,3,1); imshow(A);
subplot(1,3,2); imshow(B_same);
subplot(1,3,3); imshow(B_conv);
 
size(A)
size(B_same)
size(B_full)

         由于该滤波核是反对称的，所以可以看出conv和corr的输出结果是有差异的，上图中和右两张图像恰好是“互补”的。

        A, B_same, B_full的大小分别为[256,256]，[256,256]，[258,258]。可以看到full模式让输出图像大小（H、W两个方向相同）增大了（L-1）,L表示滤波在各方向上的长度。 回想一维卷积（或者相关）的情况，两个输入长度分别为N1,N2，则卷积/相关输出的长度变为（N1+N2-1），此处的情况只是一维情况的自然推广。

例2： 

originalRGB = imread('peppers.png');
h = fspecial('motion', 50, 45); %使用 fspecial 函数创建一个运动模糊滤波器
filteredRGB = imfilter(originalRGB, h); % 将该滤波器应用于原始图像，以创建一个具有运动模糊的图像
figure;
subplot(1,2,1); imshow(originalRGB)
subplot(1,2,2); imshow(filteredRGB)

 

2.5 2-D卷积/相关的分解

        2-D的正方形的核有可能（并非总是可能，有一定的约束）分解为两个1-D的核，从而可以将2-D的卷积/相关分解为两个1-D的卷积/相关运算。这样做可以获得巨大的运算量的收益，因为2-D卷积/相关的运算复杂度是，而1-D卷积/相关的运算复杂度是。

        卷积神经网络中的“Depthwise separable convolution”的降低运算量的思路与此类似（可能根源就在这里），有兴趣者可以参考（A brief introduction to Depthwise Separable Convolution）。

3. 高通滤波，edge detection

        “edge”定义为图像强度的局部、剧烈的变化。最简单的edge是分割相邻的均一强度区域的分界线，但是实际情况永远不会这么简单。由于有噪声的影响、成像系统的非理想的光学特性、非理想的数字化处理等等都会使得图像出现模糊。因此即便在真正的edge出现的地方，观测到的也应该是局部的像素强度的相对“缓慢”的变化。

        从时域的观点来看，edge代表着梯度最大的地方，可以用差分的方法来进行检测。

        从频域的观点来看，edge对应于图像的2-D频谱中的高频成分，因此可以通过高通滤波的方式来检测图像中的edge。

3.1 matlab edge()函数

        matlab提供了edge()函数用于图像中的边沿检测。

BW = edge(I)

        % 返回二值图像 BW，其中的值 1 对应于灰度或二值图像 I 中函数找到边缘的位置，值 0 对应于其他位置。默认情况下，edge 使用 Sobel 边缘检测方法。
BW = edge(I,method)

        % 使用 method 指定的边缘检测算法检测图像 I 中的边缘。
BW = edge(I,method,threshold)

        % 返回强度高于 threshold 的所有边缘。
BW = edge(I,method,threshold,direction)

        % 指定要检测的边缘的方向。Sobel 和 Prewitt 方法可以检测垂直方向和/或水平方向的边缘。Roberts 方法可以检测与水平方向成 45 度角和/或 135 度角的边缘。仅当 method 是 'Sobel'、'Prewitt' 或 'Roberts' 时，此语法才有效。
BW = edge(___,'nothinning')

        % 跳过边缘细化阶段，这可以提高性能。仅当 method 是 'Sobel'、'Prewitt' 或 'Roberts' 时，此语法才有效。
BW = edge(I,method,threshold,sigma)

        % 指定 sigma，即滤波器的标准差。仅当 method 是 'log' 或 'Canny' 时，此语法才有效。
BW = edge(I,method,threshold,h)

        % 使用 'zerocross' 方法和您指定的滤波器 h 检测边缘。仅当 method 是 'zerocross' 时，此语法才有效。
[BW,threshOut] = edge(___) 

        % 除了图像结果BW外，还返回阈值。
[BW,threshOut,Gv,Gh] = edge(___) 

        % 除了图像结果BW和阈值外，还返回定向梯度幅值。对于 Sobel 和 Prewitt 方法，Gv 和 Gh 对应于垂直和水平梯度。对于 Roberts 方法，Gv 和 Gh 分别对应于与水平方向成 45° 和 135° 角的梯度。仅当 method 是 'Sobel'、'Prewitt' 或 'Roberts' 时，此语法才有效。

        边缘检测方法参数method，指定为下列方法之一。

3.2 各种算子简介

3.2.1 Robert算子

        Robert算子:

                 

   

Advantages:

1. Detection of edges and orientation are very easy
2. Diagonal direction points are preserved

Limitations:

1. Very sensitive to noise
2. Not very accurate in edge detection

 3.2.2 Prewitt算子

         Prewitt算子：

                

Advantages:

1. Good performance on detecting vertical and horizontal edges
2. Best operator to detect the orientation of an image

Limitations:

1. The magnitude of coefficient is fixed and cannot be changed
2. Diagonal direction points are not preserved always

3.2.3 Sobel算子 

        Sobel算子：

                

        This operator provides greater resilience to noise and is the best estimator of edge orientation and strength of all the “small” kernels.  

Advantages:

1. Simple and time efficient computation
2. Very easy at searching for smooth edges

Limitations:

1. Diagonal direction points are not preserved always
2. Sensitive to noise
3. Not very accurate in edge detection
4. Detect with thick and rough edges does not give appropriate results

3.4 Canny算子

        Canny针对边沿检测采取了基于信息论的方法（an information theoretic approach），认为一个边沿检测器应该满足以下要求：

        1. Detect an edge

        2. Should give a response in the correct location

        3. Have a single response to an edge

        Canny基于高斯噪声的假设定义了一个匹配滤波器，这个匹配滤波器可以用高斯差分近似。该算法不仅在高斯噪声条件下给出了最优性能，在其它噪声条件下，也同样能给出充分的检测性能。Canny Edge detection包含以下四个主要步骤：

• Reduce Noise using Gaussian Smoothing.
• Compute image gradient using Sobel filter.
• Apply Non-Max Suppression or NMS to just jeep the local maxima
• Finally, apply Hysteresis thresholding with 2 threshold values T_upper and T_lower.

 

3.5 Marr-Hildreth

        Marr-Hildreth算子别名“LoG”（Laplacian of Gaussian），matlab.edge()中边沿检测方法参数method设为“log”即指这种算法。

        基于差分的边沿检测方法的缺陷在于难以精确地定义边沿的位置。这个问题可以用双重差分来解决，Marr-Hildreth算子表达式如下所示：

        Edges are at the 'zero crossings' of the LoG, which is where there is a change in gradient. 

 

3.6 Matlab代码示例

close all; clear; clc
%%  读取图像并转换成灰度图
% I = imread( "plane.jpg" ) ;
% T_B = rgb2gray( I ) ;
% I_B = im2double ( T_B ) ;
I = imread( "pout.tif" ) ;
I_B = im2double ( I ) ;
 
%% sharp（图像锐化）
%% H 为常用的拉普拉斯算子
H = [ 0 ,1 , 0 ; 0 , -4 , 0 ; 0 , 1 ,0 ] ;
T_S = conv2 ( I_B ,H ,"same" ) ;
K = T_S - I_B ; 
%% prewitt 算子计算图像边缘
[ I_E , thresh ] =edge (K , "prewitt" , [] , "both" ) ;
%% sharp（图像再次锐化）、保存
T_E = conv2 ( I_E ,H ,"same" ) ;
K_1 = T_E - I_B ;
% imwrite ( I_E , "plane_edge.jpg" ) ;
imwrite ( I_E , "pout_edge.jpg" ) ;
%
 
figure;
subplot(1,2,1); imshow(I);
subplot(1,2,2); imshow(I_E);
 
I = imread('circuit.tif');
figure;
subplot(2,2,1); imshow(I); 
BW1 = edge(I,'Canny'); % 使用 Canny 方法查找边缘。
BW2 = edge(I,'Sobel'); % 使用 Sobel 方法查找边缘。
BW3 = edge(I,'LoG');   % 使用 LoG 方法查找边缘。
subplot(2,2,2); imshow(BW1)
subplot(2,2,3); imshow(BW2)
subplot(2,2,4); imshow(BW3)

 

 

4. 低通滤波，平滑

        coming soon