摘  要
 
 中值滤波技术是一种在去除噪声的同时能较好保护图像边缘细节的非线性技术，在图像增强和恢复等领域中得到了广泛的应用。文章阐述了中值滤波的原理和特点，并使用软件工具MATLAB实现了图像的中值滤波。结果表明，利用中值滤波法消除图像中的随机噪声是图像噪声处理的最佳方法。
 
 关键词
 
 中值滤波，图像处理，MATLAB
 
 
 
 1. 引言
 
 数字图像是建立在人眼的感光特性和大量变换、处理、传输等应用技术基础之上的，在其形成、传输、接收和处理过程中，不可避免地存在着一定程度的噪声干扰，而图像的噪声是影响人眼观察的最直接、最有害和最关键的问题，如何使得图像在终端画面上即保持显示图像的各类要素信息、并得以加强，又适合人眼观察，还能合理地消除大量的噪声是图像处理研究领域进行的主要研究之一[] [] 。去噪即滤波，通常可分为线性滤波和非线性滤波两类。对受到噪声污染的退化图像的复原可以采用线性滤波方法来处理，在许多情况下是很有效的，但是多数线性滤波具有低通特性，在去除噪声的同时也使图像的边缘变模糊，不能获得较好的复原效果[] 。而中值滤波是一种去除噪声的非线性处理方法，在某些条件下既可去除噪声又能保护图像边缘，加之，中值滤波对图像信号的衰减不影响人眼对图像的感觉，因此通过中值滤波能获得较好的图像复原效果。自Turky在1971年提出中值滤波技术以来，该技术就得到了广泛的应用，并被应用到语音处理、图像恢复等领域中[] -[] 。近年，江景涛教授课题组[] 从人眼的结构和视觉特性出发深入研究了图像中噪声产生的机理和消除的方法，提出了利用中值滤波法消除图像噪声的实用方法，通过对影响图像信号的不同噪声采用了不同的滤波方法，并对信号的频谱和滤波的效果进行了研究和分析。
 
 2. 中值滤波技术
 
 中值滤波是基于排序统计理论的一种能有效抑制噪声的非线性信号处理技术。这种滤波器的优点是运算简单而且速度快，在滤除叠加白噪声和长尾叠加噪声方面显示了极好的性能。而且中值滤波器在滤除噪声的同时能很好地保护图像边缘，使图像较好地复原。另外，中值滤波器很容易自适应化，从而可以进一步提高其滤波性能。因此，它就非常适应于一些线性滤波器无法胜任的数字图像处理应用场合。
 
 2.1. 中值滤波的基本原理[] []
 
 中值滤波其基本原理是把数字图像或数字序列中一点的值用该点的一个领域中各点值的中值代替。中值的定义如下：一数组，把n个数按值的大小顺序排列于下：，
 
 
 
 y称为序列的中值。把一个点的特定长度或形状的领域称为窗口。在一维情况下，中值滤波器是一个含有奇数个像素的滑动窗口，窗口正中间那个像素的值用窗口内各像素值的中值代替。
 
 设输入序列为，I为自然数集合或子集，窗口长度为n，则滤波器输出为
 
 
 
 式中，；u = (n − 1)/2。
 
 很容易将中值滤波的概念推广到二维，此时可以利用某种形式的二维窗口。设表示数字图像各点的灰度值，滤波窗口为A的二维中值滤波可定义为二维中值滤波可以取方形，也可以取近似圆形或十字形。
 
 简言之，中值滤波的主要原理是：首先确定一个以某个像素为中心点的领域，一般为方形领域；然后将领域中的各个像素的灰度值进行排序，取其中间值作为中心点像素灰度的新值，这里的领域通常被称为窗口；当窗口在图像中上下左右进行移动后，利用中值滤波算法可以很好地对图像进行平滑处理。具体步骤如下：
 
 (1) 将模板在图像中漫游，并将模板中心与图像中心某个像素的位置重合；
 
 (2) 读取模板下各对应像素的灰度值；
 
 (3) 将这些灰度值从小到大排列成一列；
 
 (4) 找出排在中间的一个值；
 
 (5) 将这个中间值赋给对应模板中心位置的像素。
 
 中值滤波的输出像素是由领域图像的中间值决定的，因而中值滤波对极限像素值(与周围像素灰度值差别较大的像素)远不如平均值那么感，从而可以消除孤立的噪声点，又可以使图像产生较少的模糊。
 
 2.2. 中值滤波的主要特性[]-[]
 
 2.2.1. 对某些输入信号中值滤波的不变性
 
 对某些特定的输入信号，如在窗口2n + 1内单调增加或单调减少的序列，中值滤波输出信号仍保持输入信号不变，利用这个特点，可以使中值滤波既能去除图像中的噪声，又能保持图像中一些物体的边缘。
 
 2.2.2. 中值滤波去噪声性能
 
 中值滤波是非线性运算，因此对于随机性质的噪声输入，数学分析是相当复杂的。由大量实验可得，对于零均值正态分布的噪声输入，中值滤波输出与输入噪声的密度分布有关，输出噪声方差与输入噪声密度函数的平方成反比。
 
 对随机噪声的抑制能力，中值滤波性能要比平均值滤波差些。但对于脉冲干扰来讲，特别是脉冲宽度小于m/2，相距较远的窄脉冲，中值滤波是很有效的。
 
 2.2.3. 中值滤波的频谱特性
 
 由于中值滤波是非线性运算，输入和输出之间在频率上不存在一一对应关系。故不能用一般线性滤波器频率特性的研究方法。为了能够直观地定性地看出中值滤波输入和输出频谱变化情况，我们采用总体试验观察方法。设G为输入信号频谱，F为输出信号频谱，定义：H = |G/F|实验表明，H是与G有关的，呈不规则波动不大的曲线。
 
 因此，中值滤波频谱特性H起伏不大，其均值比较平坦，可以认为信号经中值滤波后，频谱基本不变，这点认识对从事设计和使用中植滤波器的工作是很有意义的。
 
 2.3. 视觉与中值滤波的关系[]
 
 图像最终需要由人来观察、辨别、理解。图像与视角之间有着重要的关系。图像质量的好坏和人眼的感光特性是密切相关的。由于人眼对图像的感觉是以亮度感觉为主的，而亮度感觉主要是由人眼的杆状细胞来完成的，杆状细胞对显示图像的感光是人眼多个细胞的平均光。因此，人眼感受的是若干个杆状细胞，对平均光的刺激，而中值滤波所取的图像信号的中值正好适合于人眼杆状细胞的感光。对图像进行中值滤波后，可以很好地消除分布在图像上的各类随机噪声，且不影响图像给人的感觉。在图像传输过程中，外部干扰和系统内部干扰会给图像带来很多的随机噪声干扰，利用中值滤波就能很好地消除这些干扰，而中值滤波对图像信号的衰减又不影响人眼对图像的感觉。因此，利用中值滤波法来消除图像的随机干扰噪声是图像随机噪声消除的最佳方法。
 
 3. MATLAB软件介绍
 
 MATLAB是Maths Works公司的商业软件，是一种直译式语言，其最重要的功能就是进行矩阵的数值运算，它的数值分析、模拟与运算功能非常强大，而且程序结构完整，具有很强的平行移植性，不仅如此，MATLAB语言编程与C语言十分类似，语法简单、易学易用[] -[] 。因此，在图像处理、自动控制、语音处理、信号分析等工业领域，MATLAB是研究、开发和分析首选的计算工具。MATLAB提供了强大的矩阵运算和绘图功能，其中，图像处理工具箱(Image Processing Toolbox)提供了将近200种最基本的图像处理函数，利用这些图像处理工具箱， 结合其强大的数据处理能力，我们可不必关心图像文件的格式、读写、显示等细节，而把精力集中在算法研究上，大大提高了工作效率。同时，在测试这些算法时既可方便地得到统计数据，又可得到直观图示。
 
 利用MATLAB语言及其所提供的函数进行图像处理，编程简单、操作方便、可靠性强，能够达到预期的效果。MATLAB工具箱功能非常强大，简化了繁琐的数学计算工作，有效的促进了数字图像处理相关的研究开发工作。
 
 4. 中值滤波在图像处理中的MATLAB实现
 
 源程序如下：
 
 I=imread('D:\myimages\chost.jpg')；%从D盘名为myimages的文件夹中读取。格式为jpg的图像文件chost J=imnoise(I,'salt & pepper',0.02)；%给图像加入均值为0，方差为0.02的淑盐噪声
 
 subplot(2,4,1);
 
 imshow(I);
 
 title('原始图像');
 
 subplot(2,4,2);
 
 imshow(J);
 
 title('加入椒盐噪声之后的图像');
 
 %h=ones(3,3)/9; %产生3 × 3的全1数组
 
 %B=conv2(J,h); %卷积运算
 
 %采用MATLAB中的函数对噪声干扰的图像进行滤波
 
 Q=wiener2(J,[3 3]); %对加噪图像进行二维自适应维纳滤波
 
 P=filter2(fspecial('average',3),J)/255; %均值滤波模板尺寸为3 K1=medfilt2(J,[3 3]); %进行3 × 3模板的中值滤波
 
 K2= medfilt2(J,[5 5]); %进行5 × 5模板的中值滤波
 
 K3= medfilt2(J,[7 7]); %进行7 × 7模板的中值滤波
 
 K4= medfilt2(J,[9 9]); %进行9 × 9模板的中值滤波
 
 %显示滤波后的图像及标题
 
 subplot(2,4,3);
 
 imshow(Q);
 
 title('3 × 3模板维纳滤波后的图像');
 
 subplot(2,4,4);
 
 imshow(P);
 
 title('3 × 3模板均值滤波后的图像');
 
 subplot(2,4,5);
 
 imshow(K1);
 
 title('3 × 3模板的中值滤波的图像');
 
 subplot(2,4,6);
 
 imshow(K2);
 
 title('5 × 5模板的中值滤波的图像');
 
 subplot(2,4, 7);
 
 imshow(K3);
 
 title('7 × 7模板的中值滤波的图像');
 
 subplot(2,4,8);
 
 imshow(K4);
 
 title('9 × 9模板的中值滤波的图像');
 
 得到图像结果如图1所示。
 
 
 
 图1. 滤波效果对比图
 
 
 ating5. 结论
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 (1) 采用不同的滤波方法对图像进行滤波，通过进行比较，可以看到：均值滤波和维纳滤波在衰减噪声的同时也会使图像的边界等细节模糊；而采用中值滤波能较好保护图像边缘细节，消除图像的随机噪声效果很好，同时，滤波后的图像又能符合人眼的视觉感觉。这也是中值滤波器受欢迎的主要原因。
 
 (2) 由上图可见，同一种滤波方法并不是窗口越大滤波效果越好。通过比较，3 × 3模板的中值滤波的效果是最好的，随着所使用的中值滤波器尺寸的增大，消除噪声的效果得到增强；但是，图像的细节锐化程度相应降低，图像变得模糊起来。从而造成图像细节信息的丢失，也破坏了图像的精度。因此选择适当大小的中值滤波窗口可以在最大限度的保持图像精度的基础上去除图像噪声。
 
 综上所述，中值滤波在去除噪声的同时能较好保护图像边缘细节信息，解决了多数线性滤波在去噪的同时模糊图像这一缺点，获得了较好的图像复原效果。因此，利用中值滤波法消除图像中的随机噪声是图像噪声处理的最佳方法。

一、获取代码方式
 
获取代码方式1：
 通过订阅紫极神光博客付费专栏，凭支付凭证，私信博主，可获得此代码。
 
获取代码方式2：
 完整代码已上传我的资源：【图像去噪】基于matlab中值+均值+高斯+Laplacian+Sobel+Prewitt图像去噪【含Matlab源码 025期】
 
备注：
 订阅紫极神光博客付费专栏，可免费获得1份代码（有效期为订阅日起，三天内有效）；
 
二、简介
 
1 中值滤波
 （1）概念：
 
 （2）原理解释：
 
 2 均值滤波
 均值滤波是指任意一点的像素值，都是周围 N \times M 个像素值的均值。例如下图中，红色点的像素值是其周围蓝色背景区域像素值之和除25，25=5\times5 是蓝色区域的大小。
 
 均值滤波详细的计算方法如下图所示：
 
 其中5\times5的矩阵称为核，针对原始图像内的像素点，采用核进行处理，得到结果图像，如下图所示：
 
 
 提取 1/25 可以将核转换为如下形式：
 
 3 高斯滤波
 
 （1）核大小为 33
 
 （2）核大小为 55
 
 高斯滤波让临近的像素具有更高的重要度，对周围像素计算加权平均值，较近的像素具有较大的权重值。如下图所示，中心位置权重最高为0.4。
 
 
 
三、部分源代码
 
clc
clear all
close all
img = imread('testPic.tif');
img=img(1:128*16,1:128*16);
imwrite(img,'Source0.tif')
%% 均值滤波部分
figure(1);
subplot(3,4,5)
subimage(img);
title('原始图像');
Byave = filter2(fspecial('average',3),img)/255; %进行3*3均值滤波
%空域锐化
p1 =fspecial('sobel');
subplot(3,4,2)
add1 =imfilter(Byave,p1);
out1=add1+Byave;
subimage(out1);
title('sobel锐化结果(均值滤波）');
[SNR,MSE,PSNR]=getSNR(img,out1);
fprintf('均值滤波+sobel算子锐化:\nSNR:%f\nMSE:%f\nPSNR:%f\n',SNR,MSE,PSNR)
subplot(3,4,3)
p2 =fspecial('prewitt');
add2 =imfilter(Byave,p2);
out2=add2+Byave;
subimage(out2);
title('prewitt锐化结果(均值滤波）');
[SNR,MSE,PSNR]=getSNR(img,out2);
fprintf('均值滤波+prewitt算子锐化:\nSNR:%f\nMSE:%f\nPSNR:%f\n',SNR,MSE,PSNR)
subplot(3,4,4)
p3 =fspecial('laplacian');
add3 =imfilter(Byave,p3);
out3=add3+Byave;
imwrite(out3,'均值滤波+laplacian算子锐化.tif');
subimage(out3);
title('laplacian锐化结果(均值滤波））');
[SNR,MSE,PSNR]=getSNR(img,out2);
fprintf('均值滤波+laplacian算子锐化:\nSNR:%f\nMSE:%f\nPSNR:%f\n',SNR,MSE,PSNR)
%% 中值滤波部分
Bymed=medfilt2(img,[3,3]);
%空域锐化
p1 =fspecial('sobel');
subplot(3,4,6)
add1 =imfilter(Bymed,p1);
out1=add1+Bymed;
subimage(out1);
title('sobel锐化结果(中值滤波）');
[SNR,MSE,PSNR]=getSNR(img,out1);
fprintf('中值滤波+sobel算子锐化:\nSNR:%f\nMSE:%f\nPSNR:%f\n',SNR,MSE,PSNR)
subplot(3,4,7)
p2 =fspecial('prewitt');
add2 =imfilter(Bymed,p2);
out2=add2+Bymed;
subimage(out2);
title('prewitt锐化结果(中值滤波）');
[SNR,MSE,PSNR]=getSNR(img,out2);
fprintf('中值滤波+prewitt算子锐化:\nSNR:%f\nMSE:%f\nPSNR:%f\n',SNR,MSE,PSNR)
subplot(3,4,8)
p3 =fspecial('laplacian');
add3 =imfilter(Bymed,p3);
out3=add3+Bymed;
imwrite(out3,'中值滤波+laplacian算子锐化.tif');
subimage(out3);
title('laplacian锐化结果(中值滤波）');
[SNR,MSE,PSNR]=getSNR(img,out3);
fprintf('中值滤波+laplacian算子锐化:\nSNR:%f\nMSE:%f\nPSNR:%f\n',SNR,MSE,PSNR)
%% 高斯滤波部分
W = fspecial('gaussian',[3,3],1); 
Bygaussian = imfilter(img, W, 'replicate');
 p1 =fspecial('sobel');
subplot(3,4,10)
add1 =imfilter(Bygaussian,p1);
out1=add1+Bygaussian;
subimage(out1);

 
四、运行结果
 
 
五、matlab版本及参考文献
 
1 matlab版本
 2014a
 
2 参考文献
 [1] 蔡利梅.MATLAB图像处理——理论、算法与实例分析[M].清华大学出版社，2020.
 [2]杨丹,赵海滨,龙哲.MATLAB图像处理实例详解[M].清华大学出版社，2013.
 [3]周品.MATLAB图像处理与图形用户界面设计[M].清华大学出版社，2013.
 [4]刘成龙.精通MATLAB图像处理[M].清华大学出版社，2015.
 [5]陈浩,方勇,朱大洲,王成,陈子龙.基于蚁群算法的玉米植株热红外图像边缘检测[J].农机化研究. 2015,37(06)

图像的平滑、锐化都是利用掩模操作来完成的。通过掩模操作实现一种邻域运算，待处理像素点的结果由邻域的图像像素以及相应的与邻域有相同维数的子图像得到。这些子图像被称为滤波器、掩模、核、模板或窗口；
 掩模运算的数学含义是卷积（或互相关）运算；
 掩模子图像中的值是系数值，而不是灰度值；
 卷积示例图：
 
 一般来说，在MN的图像f(x,y)上，用mn大小的滤波器掩模进行线性滤波由下式给出：
 
 模板为1*5的中值滤波和均值滤波的对比：
 
 
均值滤波
 
简单来说就是对某个区域内的像素值取平均值代替原像素值
 常用的3*3的滤波器掩模为：
 
 
 一幅M×N的图像经过m×n的加权均值滤波器滤波的过程可由下式给出：
 

 一般选取n*n的模板，便于运算，下面给出示例代码：
 
img = imread('');

[M , N] = size(img);%图片尺寸
img_result = zeros(M, N);%预生成，提高速度

muban_size = 3;%模板尺寸
expand_size = floor(muban_size / 2);%扩展尺寸
muban = 1 / (muban_size * muban_size) .* ones(muban_size, muban_size);


expand_img = double(wextend('2D','zpd', img, expand_size));%扩展0，转double为了矩阵运算

for i=1:M
    for j=1:N
        ave = sum( sum( expand_img(i:i+muban_size-1,j:j+muban_size-1) .* muban)); %取出扩展元素与模板相乘,并求矩阵元素之和
        img_result(i,j) = ave;
    end
end

img_result = uint8(img_result);%转int8，图像
subplot(1 ,2, 1);
title('原图像')
imshow(img)
subplot(1 ,2, 2);
imshow(img_result)
da = ['模板大小为' num2str(muban_size) '，变化后的图像'];
title(da)
 
结果示例：
 
 可见均值滤波对于噪声有一定的抑制作用，但是会出现部分的涂抹感。
 
中值滤波
 
中值滤波和均值滤波不同的地方是，中值滤波是对图像的像素值进行排序，取中间的像素值赋给新的图像。
 主要功能：使拥有不同灰度的点看起来更接近于它的邻近值。
 主要用途：去除“椒盐”噪声
 
示例代码：
 
img = imread('');

[M , N] = size(img);%图片尺寸
img_result = zeros(M, N);%预生成，提高速度

muban_size = 3;%模板尺寸
expand_size = floor(muban_size / 2);%扩展尺寸
muban = ones(muban_size, muban_size);


expand_img = double(wextend('2D','zpd', img, expand_size));%扩展0，转double为了矩阵运算

for i=1:M
    for j=1:N
        mat = expand_img(i:i+muban_size-1,j:j+muban_size-1) .* muban; %取出x1中从(i,j)开始的n行n列元素与模板相乘
        mat = mat(:);%转数组
        mat = sort(mat);%排序
        if mod(muban_size, 2)==1
            img_result(i,j) = mat(floor(muban_size*muban_size/2)+1);%取中间
        else
            img_result(i,j) = (mat(muban_size*muban_size/2) + mat(muban_size*muban_size/2+1))/2;
        end
    end
end

img_result = uint8(img_result);%转int8，图像
subplot(1 ,2, 1);
title('原图像')
imshow(img)
subplot(1 ,2, 2);
imshow(img_result)
da = ['模板大小为' num2str(muban_size) '，变化后的图像'];
title(da)
 
结果示例：
 
 效果好像特别好。
 
一些其他的中值滤波器：
 

 另：
 中值滤波的窗口形状和尺寸对滤波效果影响较大，不同的图像内容和不同的应用要求，往往采用不同的窗口形状和尺寸。常用的二维中值滤波窗口有线状、方形、圆形、十字形以及圆环形等。
 窗口尺寸一般先用3X3，再取5X 5逐渐增大，直到滤波效果满意为止。就经验来讲，对于有缓变的较长轮廓线物体的图像，采用方形或圆形窗口为宜。对于包含有尖顶物体的图像，用十字形窗口，而窗口大小则以不超过图像中最小有效物体的尺寸为宜。如果图像中点、线、尖角细节较多，则不宜采用中值滤波
 
最后，点个赞？

 
图像的中值滤波、均值滤波及边缘检测
 
原理概述
图像滤波
中值滤波
均值滤波
边缘检测
  
代码
中值滤波
均值滤波
Laplace算子
prewitt算子
roberts算子
sobel算子
  
结果
总结
 

 
原理概述
 
图像滤波
 
图像滤波：即在尽量保留图像细节特征的条件下对目标图像的噪声进行抑制，是图像预处理中不可缺少的操作，其处理效果的好坏将直接影响到后续图像处理和分析的有效性和可靠性。消除图像中的噪声成分叫作图像的平滑化或滤波操作。信号或图像的能量大部分集中在幅度谱的低频和中频段是很常见的，而在较高频段，感兴趣的信息经常被噪声淹没。因此一个能降低高频成分幅度的滤波器就能够减弱噪声的影响。
 图像滤波的目的有两个:一是抽出对象的特征作为图像识别的特征模式;另一个是为适应图像处理的要求，消除图像数字化时所混入的噪声。
 而对滤波处理的要求也有两条:一是不能损坏图像的轮廓及边缘等重要信息;二是使图像清晰视觉效果好。
 关于滤波器，一种形象的比喻法是：我们可以把滤波器想象成一个包含加权系数的窗口，当使用这个滤波器平滑处理图像时，就把这个窗口放到图像之上，透过这个窗口来看我们得到的图像。举一个滤波在我们生活中的应用：美颜的磨皮功能。如果将我们脸上坑坑洼洼比作是噪声的话，那么滤波算法就是来取出这些噪声，使我们自拍的皮肤看起来很光滑。
 
中值滤波
 
中值滤波是一种典型的非线性滤波，是基于排序统计理论的一种能够有效抑制噪声的非线性信号处理技术，基本思想是用像素点邻域灰度值的中值来代替该像素点的灰度值，让周围的像素值接近真实的值从而消除孤立的噪声点。该方法在取出脉冲噪声、椒盐噪声的同时能保留图像的边缘细节。这些优良特性是线性滤波所不具备的。中值滤波首先也得生成一个滤波模板，将该模板内的各像素值进行排序，生成单调上升或单调下降的二维数据序列，二维中值滤波输出为
 g(x, y)=medf{f(x-k, y-1),(k, l∈w)}
 其中f(x,y)和g(x,y)分别是原图像和处理后图像, w为输入的二维模板，能够在整幅图像上滑动，通常尺寸为33或55区域，也可以是不同的形状如线状、圆形、十字形、圆环形等。通过从图像中的二维模板取出奇数个数据进行排序，用排序后的中值取代要处理的数据即可。
 下图为中值滤波示意图：
 四邻域八邻域
 
均值滤波
 
均值的原理是，一个规定的邻域内，所有像素的平局值作为最终计算的结果，每个像素的权值相同，为总像素的倒数。均值平滑都是线性的，一旦参数给定，模板就确定下来，不会因为位置和像素分布不同而改变，而线性模板的基本运算是卷积。
 
 
边缘检测
 
边缘检测是图像处理和计算机视觉中的基本问题，边缘检测的目的是标识数字图像中亮度变化明显的点。图像属性中的显著变化通常反映了属性的重要事件和变化。这些包括：深度上的不连续、表面方向不连续、物质属性变化和场景照明变化。 边缘检测是图像处理和计算机视觉中，尤其是特征提取中的一个研究领域。图像边缘检测大幅度地减少了数据量，并且剔除了可以认为不相关的信息，保留了图像重要的结构属性。有许多方法用于边缘检测，它们的绝大部分可以划分为两类：基于查找一类和基于零穿越的一类。基于查找的方法通过寻找图像一阶导数中的最大和最小值来检测边界，通常是将边界定位在梯度最大的方向。基于零穿越的方法通过寻找图像二阶导数零穿越来寻找边界，通常是Laplacian过零点或者非线性差分表示的过零点。
 人类视觉系统认识目标的过程分为两步:首先，把图像边缘与背景分离出来;然后，才能知觉到图像的细节，辨认出图像的轮廓。计算机视觉正是模仿人类视觉的这个过程。因此在检测物体边缘时，先对其轮廓点进行粗略检测，然后通过链接规则把原来检测到的轮廓点连接起来，同时也检测和连接遗漏的边界点及去除虚假的边界点。图像的边缘是图像的重要特征，是计算机视觉、模式识别等的基础，因此边缘检测是图象处理中一个重要的环节。然而，边缘检测又是图象处理中的一个难题，由于实际景物图像的边缘往往是各种类型的边缘及它们模糊化后结果的组合，且实际图像信号存在着噪声。噪声和边缘都属于高频信号，很难用频带做取舍。这就需要边缘检测来进行解决的问题了。边缘检测的基本方法有很多：一阶的有Laplace算子，Roberts算子，Prewitt算子，Sobel算子，Canny算子等。
 
Laplace算子是最简单的各向同性微分算子。一个二维图像函数f(x,y) f(x, y)f(x,y) 定义为：
 
 
 Roberts算子是一种利用局部差分算子寻找边缘的算子，它有下式给出：
 
 其中f(x,y)、f(x+1,y+1)、f(x,y+1)、f(x+1,y) 分别为4领域的坐标，且是具有整数像素坐标的输人图像;其中的平方根运算使得该处理类似于人类视觉系统中发生的过程。Roberts算子是2X2算子模板。图1所示的2个卷积核形成了Roberts算子。图象中的每一个点都用这2个核做卷积。
 
 Sobel算子是一种一阶微分算子，它利用像素邻近区域的梯度值来计算1个像素的梯度，然后根据一定的绝对值来取舍。它由下式给出:
 
 Sobel算子是33算子模板。图2所示的2个卷积核dx 、 dy形成Sobel算子。一个核对通常的垂直边缘响应最大，而另一个核对水平边缘响应最大。2个卷积的最大值作为该点的输出值。运算结果是一幅边缘幅度图像。
 
 Prewitt算子由下式给出：
 
 Prewitt算子是33算子模板。图3所示的2个卷积核dx ,不要形成了Prewitt算子。与Sobel算子的方法一样，图像中的每个点都用这2个核进行卷积，取最大值作为输出值。Prewitt算子也产生一幅边缘幅度图像。
 
 Canny算子是是一阶算子。其方法的实质是用1个准高斯函数作平滑运算fs=f(x,y)G(x,y),然后以带方向的一阶微分算子定位导数最大值。
 平滑后fs(x，y)的梯度可以使用22的一阶有限差分近似式:
 P[i,j]≈(fs[i,j+1]-fs[i,j]+fs[i+1,j+1]-fs[i+1,j])/2
 Q[i,j]≈(fs[i,j]-fs[i+1,j]+fs[i,j+1]-fs[i+1,j+1])/2
 在这个2x2正方形内求有限差分的均值，便于在图像中的同一点计算二和y的偏导数梯度。幅值和方向角可用直角坐标到极坐标的坐标转化来计算:
 
 M[i,j]反映了图象的边缘强度;反映了边缘的方向。使得M}i,j}取得局部最大值的方向角，就反映了边缘的方向。
 Canny算子也可用高斯函数的梯度来近似，在理论上很接近4个指数函数的线性组合形成的最佳边缘算子。在实际工作应用中编程较为复杂且运算较慢。
 
代码
 
中值滤波
 
% 灰度图像中值滤波平滑处理
function result = median_filtering(image, width, leng)
result = im2uint8(zeros(size(image)));
[rows,cols] = size(image); % 图像行列数
temp = im2uint8(zeros([rows+2,cols+2]));
temp(2:rows+1,2:cols+1) = image;

for i = 1:(rows + 2 - leng + 1)
    for j = 1:(cols + 2 - width + 1)
        data = temp(i:i + leng - 1, j:j + width - 1); % 取3*3矩阵
        data = sort(data(:));                          % 返回排序结果
        result(i, j) = data((length(data)+1)/2);
    end                                                % 赋值
end

% image = imread('noise.jpg');
% width=3;leng=3;
% result = median_filtering(image, width, leng);
% imshow(result)
 
均值滤波
 
% 灰度图像均值滤波平滑处理
function result = average_filtering(image, width, leng)
result = im2uint8(zeros(size(image)));
[rows,cols] = size(image); % 图像行列数
temp = im2uint8(zeros([rows+2,cols+2]));
temp(2:rows+1,2:cols+1) = image;

for i = 1:(rows + 2 - leng + 1)
    for j = 1:(cols + 2 - width + 1)
        data = temp(i:i + leng - 1, j:j + width - 1); % 取3*3矩阵
        average = mean(data(:));                          % 返回排序结果
        result(i, j) = average; % 赋值
    end                                                
end

% image = imread('noise.jpg');
% width=3;leng=3;
% result = average_filtering(image, width, leng);
% imshow(result)
 
Laplace算子
 
% 灰度图像Laplace锐化处理
function laplaceion = laplace(image,threshold)

image = double(image);
[rows,cols,n] = size(image); % 图像行列数
temp1 = int16(zeros([rows+2,cols+2]));
temp2 = int16(zeros([rows+2,cols+2]));
temp3 = int16(zeros([rows+2,cols+2]));
temp1(2:rows+1,2:cols+1) = image(:, :, 1);
temp2(2:rows+1,2:cols+1) = image(:, :, 2);
temp3(2:rows+1,2:cols+1) = image(:, :, 3);

laplaceion = uint8(zeros(size(image)));

for i = 2:(rows + 1)
    for j = 2:(cols + 1)
        laplaceion(i-1,j-1,1) = temp1(i-1,j+1) + temp1(i+1,j+1) + temp1(i-1,j-1) + temp1(i+1,j-1) + temp1(i+1,j) + temp1(i-1,j) + temp1(i,j+1) + temp1(i,j-1) - 8 * temp1(i,j);
        if laplaceion(i-1,j-1,1) <= threshold
            laplaceion(i-1,j-1,1) = 0;
        else
            laplaceion(i-1,j-1,1) = 255;
        end
    end                                                
end

for i = 2:(rows + 1)
    for j = 2:(cols + 1)
        laplaceion(i-1,j-1,2) = temp2(i-1,j+1) + temp2(i+1,j+1) + temp2(i-1,j-1) + temp2(i+1,j-1) + temp2(i+1,j) + temp2(i-1,j) + temp2(i,j+1) + temp2(i,j-1) - 8 * temp2(i,j);
        if laplaceion(i-1,j-1,2) <= threshold  
            laplaceion(i-1,j-1,2) = 0;
        else
            laplaceion(i-1,j-1,2) = 255;
        end
    end                                                
end

for i = 2:(rows + 1)
    for j = 2:(cols + 1)
        laplaceion(i-1,j-1,3) = temp3(i-1,j+1) + temp3(i+1,j+1) + temp3(i-1,j-1) + temp3(i+1,j-1) + temp3(i+1,j) + temp3(i-1,j) + temp3(i,j+1) + temp3(i,j-1) - 8 * temp3 (i,j);
        if laplaceion(i-1,j-1,3) <= threshold  
            laplaceion(i-1,j-1,3) = 0;
        else
            laplaceion(i-1,j-1,3) = 255;
        end
    end                                                
end

imwrite(laplaceion,'laplace.jpg');
% 与库函数的laplace算子处理效果进行比较
% H = [1 1 1; 1 -8 1; 1 1 1];
% image = uint8(image);
% res = imfilter(image,H);
% cmp=(res==laplaceion);
% length(cmp(cmp==0))
% 测试代码
% image = imread('test.jpg');
% laplaceion = laplace(image);
% result = image - laplaceion;
% imshow(result);
% imwrite(result,'laplace_result.jpg');

 
prewitt算子
 
% 灰度图像prewitt锐化处理
function result = prewitt(image,threshold)

image = double(image);
[rows,cols,n] = size(image); % 图像行列数
temp1 = int16(zeros([rows+2,cols+2]));
temp2 = int16(zeros([rows+2,cols+2]));
temp3 = int16(zeros([rows+2,cols+2]));
temp1(2:rows+1,2:cols+1) = image(:, :, 1);
temp2(2:rows+1,2:cols+1) = image(:, :, 2);
temp3(2:rows+1,2:cols+1) = image(:, :, 3);

result = uint8(zeros(size(image)));

for i = 2:(rows + 1)
    for j = 2:(cols + 1)
        result(i-1,j-1,1) = abs(temp1(i-1,j-1)+temp1(i,j-1)+temp1(i+1,j-1)-temp1(i-1,j+1)-temp1(i,j+1)-temp1(i+1,j+1))+abs(temp1(i+1,j-1)+temp1(i+1,j)+temp1(i+1,j+1)-temp1(i-1,j-1)-temp1(i-1,j)-temp1(i-1,j+1));
        if result(i-1,j-1,1) <= threshold
            result(i-1,j-1,1) = 0;
        else
            result(i-1,j-1,1) = 255;
        end
    end                                                
end

for i = 2:(rows + 1)
    for j = 2:(cols + 1)
        result(i-1,j-1,2) = abs(temp2(i-1,j-1)+temp2(i,j-1)+temp2(i+1,j-1)-temp2(i-1,j+1)-temp2(i,j+1)-temp2(i+1,j+1))+abs(temp2(i+1,j-1)+temp2(i+1,j)+temp2(i+1,j+1)-temp2(i-1,j-1)-temp2(i-1,j)-temp2(i-1,j+1));
        if result(i-1,j-1,2) <= threshold  
            result(i-1,j-1,2) = 0;
        else
            result(i-1,j-1,2) = 255;
        end
    end                                                
end

for i = 2:(rows + 1)
    for j = 2:(cols + 1)
        result(i-1,j-1,3) = abs(temp3(i-1,j-1)+temp3(i,j-1)+temp3(i+1,j-1)-temp3(i-1,j+1)-temp3(i,j+1)-temp3(i+1,j+1))+abs(temp3(i+1,j-1)+temp3(i+1,j)+temp3(i+1,j+1)-temp3(i-1,j-1)-temp3(i-1,j)-temp3(i-1,j+1));
        if result(i-1,j-1,3) <= threshold  
            result(i-1,j-1,3) = 0;
        else
            result(i-1,j-1,3) = 255;
        end
    end                                                
end

imwrite(result,'prewitt_result.jpg');

% 测试代码
% image = imread('test.jpg');
% result = prewitt(image,240);
% imshow(result)

 
roberts算子
 
% 灰度图像Laplace锐化处理
function result = roberts(image,threshold)

image = double(image);
[rows,cols,n] = size(image); % 图像行列数
temp1 = int16(zeros([rows+2,cols+2]));
temp2 = int16(zeros([rows+2,cols+2]));
temp3 = int16(zeros([rows+2,cols+2]));
temp1(2:rows+1,2:cols+1) = image(:, :, 1);
temp2(2:rows+1,2:cols+1) = image(:, :, 2);
temp3(2:rows+1,2:cols+1) = image(:, :, 3);

result = uint8(zeros(size(image)));

for i = 2:(rows + 1)
    for j = 2:(cols + 1)
        result(i-1,j-1,1) = abs(temp1(i+1,j+1)-temp1(i,j))+abs(temp1(i,j+1)-temp1(i+1,j));
        if result(i-1,j-1,1) <= threshold
            result(i-1,j-1,1) = 0;
        else
            result(i-1,j-1,1) = 255;
        end
    end                                                
end

for i = 2:(rows + 1)
    for j = 2:(cols + 1)
        result(i-1,j-1,2) = abs(temp2(i+1,j+1)-temp2(i,j))+abs(temp2(i,j+1)-temp2(i+1,j));
        if result(i-1,j-1,2) <= threshold  
            result(i-1,j-1,2) = 0;
        else
            result(i-1,j-1,2) = 255;
        end
    end                                                
end

for i = 2:(rows + 1)
    for j = 2:(cols + 1)
        result(i-1,j-1,3) = abs(temp3(i+1,j+1)-temp3(i,j))+abs(temp3(i,j+1)-temp3(i+1,j));
        if result(i-1,j-1,3) <= threshold  
            result(i-1,j-1,3) = 0;
        else
            result(i-1,j-1,3) = 255;
        end
    end                                                
end

imwrite(result,'roberts_result.jpg');

% 测试代码
% image = imread('test.jpg');
% result = roberts(image,240);
% imshow(result);

 
sobel算子
 
% 灰度图像sobel锐化处理
function result = sobel(image,threshold)

image = double(image);
[rows,cols,n] = size(image); % 图像行列数
temp1 = int16(zeros([rows+2,cols+2]));
temp2 = int16(zeros([rows+2,cols+2]));
temp3 = int16(zeros([rows+2,cols+2]));
temp1(2:rows+1,2:cols+1) = image(:, :, 1);
temp2(2:rows+1,2:cols+1) = image(:, :, 2);
temp3(2:rows+1,2:cols+1) = image(:, :, 3);

result = uint8(zeros(size(image)));

for i = 2:(rows + 1)
    for j = 2:(cols + 1)
        result(i-1,j-1,1) = abs(temp1(i-1,j-1)+2*temp1(i,j-1)+temp1(i+1,j-1)-temp1(i-1,j+1)-2*temp1(i,j+1)-temp1(i+1,j+1))+abs(temp1(i+1,j-1)+2*temp1(i+1,j)+temp1(i+1,j+1)-temp1(i-1,j-1)-2*temp1(i-1,j)-temp1(i-1,j+1));
        if result(i-1,j-1,1) <= threshold
            result(i-1,j-1,1) = 0;
        else
            result(i-1,j-1,1) = 255;
        end
    end                                                
end

for i = 2:(rows + 1)
    for j = 2:(cols + 1)
        result(i-1,j-1,2) = abs(temp2(i-1,j-1)+2*temp2(i,j-1)+temp2(i+1,j-1)-temp2(i-1,j+1)-2*temp2(i,j+1)-temp2(i+1,j+1))+abs(temp2(i+1,j-1)+2*temp2(i+1,j)+temp2(i+1,j+1)-temp2(i-1,j-1)-2*temp2(i-1,j)-temp2(i-1,j+1));
        if result(i-1,j-1,2) <= threshold  
            result(i-1,j-1,2) = 0;
        else
            result(i-1,j-1,2) = 255;
        end
    end                                                
end

for i = 2:(rows + 1)
    for j = 2:(cols + 1)
        result(i-1,j-1,3) = abs(temp3(i-1,j-1)+2*temp3(i,j-1)+temp3(i+1,j-1)-temp3(i-1,j+1)-2*temp3(i,j+1)-temp3(i+1,j+1))+abs(temp3(i+1,j-1)+2*temp3(i+1,j)+temp3(i+1,j+1)-temp3(i-1,j-1)-2*temp3(i-1,j)-temp3(i-1,j+1));
        if result(i-1,j-1,3) <= threshold  
            result(i-1,j-1,3) = 0;
        else
            result(i-1,j-1,3) = 255;
        end
    end                                                
end

imwrite(result,'sobel_result.jpg');

% 测试代码
% image = imread('test.jpg');
% result = sobel(image,240);
% imshow(result)

 
结果
 
原图（含椒盐噪声）
 
 中值滤波结果、均值滤波结果
  
 laplace算子边缘检测结果
 
 prewitt算子边缘检测结果
 
 roberts算子边缘检测结果
 
 sobel算子边缘检测结果
 
 
总结
 
均值滤波和和中值滤波都可以起到平滑图像，虑去噪声的功能。均值滤波采用线性的方法，平均整个窗口范围内的像素值，均值滤波本身存在着固有的缺陷，即它不能很好地保护图像细节，在图像去噪的同时也破坏了图像的细节部分，从而使图像变得模糊，不能很好地去除噪声点。均值滤波对高斯噪声表现较好，对椒盐噪声表现较差。中值滤波采用非线性的方法，它在平滑脉冲噪声方面非常有效,同时它可以保护图像尖锐的边缘，选择适当的点来替代污染点的值，所以处理效果好，对椒盐噪声表现较好，对高斯噪声表现较差。
 laplace算子是与方向无光的各向同性边缘检测算子,若只关心边缘点的位置而不顾其周围的实际灰度差时,一般选择该算子进行检测.特点:各向同性,线性和位移是不变的,对线性和孤立点检测效果好, 但边缘方向信息丢失,常产生双像素的边缘,对噪声有双倍的加强作用。Roberts算子边缘定位准，但是对噪声敏感。适用于边缘明显且噪声较少的图像分割。Roberts边缘检测算子是一种利用局部差分算子寻找边缘的算子,Robert算子图像处理后结果边缘不是很平滑。经分析，由于Robert算子通常会在图像边缘附近的区域内产生较宽的响应，故采用上述算子检测的边缘图像常需做细化处理，边缘定位的精度不是很高。Prewitt算子对噪声有抑制作用，抑制噪声的原理是通过像素平均，但是像素平均相当于对图像的低通滤波，所以Prewitt算子对边缘的定位不如Roberts算子。Sobel算子和Prewitt算子都是加权平均，但是Sobel算子认为，邻域的像素对当前像素产生的影响不是等价的，所以距离不同的像素具有不同的权值，对算子结果产生的影响也不同。一般来说，距离越远，产生的影响越小。