注：目前开通个人网站朝思录，之后的博文将在上面更新，CSDN博客会滞后一点

---

前言

中值滤波算法是一种非线性的滤波算法，其中心思想是采用模板内所有像素的排序中值作为目标像素的值，进行滤波。一般情况下，若模板大小为m×m，可有效滤除面积小于m2/2的脉冲像素团。比如3×3中值滤波模板可有效滤除面积为1的椒盐噪声。
最近在写中值滤波算法时，因为一些bug偶然发现了一些有趣的现象。中值滤波算法的大致算法流程为： % 算法1 % 读取图像 im = imread('某个地址'); % 中值滤波 im_backup = im; %备份图像 for i = size(im_backup,1) for j = size(im_backup,2) value = GetMedian(im_backup,i,j,3); %获取3*3模板的中值 im(i,j) = value; end end 注：GetMedian为自定义的获取模板中值的函数。
其中比较关键的部分是备份图像，得到im和im_backup，并采用im_backup逐像素得出中值修改im。备份图像的目的是防止已修改的像素对下一像素的模板中值产生影响。（在线性滤波中这种备份很重要，否则滤波就不再具有线性性）。
但是，刚开始写程序时，笔者忘记备份图像，此时笔者算法流程如下所示： % 算法2 % 读取图像 im = imread('某个地址'); % 中值滤波 for i = size(im,1) for j = size(im,2) value = GetMedian(im,i,j,3); %获取3*3模板的中值 im(i,j) = value; end end 但是采用该“错误”的算法取得的结果出乎意料的好，椒盐噪声滤除非常干净，甚至比原来的中值滤波算法好。因为算法2中，除了第一次迭代，之后所有迭代中模板里总包含被更新的像素。为了方便区分，故称算法2为更新式中值滤波算法，算法1为传统中值滤波算法。

对比

1、椒盐噪声滤除效果的对比

原始图像采用冈萨雷斯第三版中受椒盐噪声污染的PCB图。

采用的模板大小为3×3，对比结果如表所示： 3×3传统中值滤波算法 3×3更新式中值滤波算法 可以发现在同样的滤波条件下，更新式中值滤波算法在牺牲轻微锐度（如QFP下的总线边缘略不清晰）的基础上滤除了更多噪声。而牺牲的锐度可以通过锐化来补偿。
当然，在模板大小为5×5时，传统中值滤波算法也会滤除掉大部分噪声，但此时图像引入了非常多的blurring，导致图像变得更加模糊。对比如下图所示。 5×5传统中值滤波 3×3更新中值滤波 所以在某些应用中，可能采用更新式中值滤波算法效果较好。除了椒盐噪声的去除效果更好之外，其不需要备份图像更节省了存储空间。

2、彩 {MOD}图像滤波效果对比

待滤波图像如下图所示。该图一度作为笔者的壁纸使用。出处不详，往知情者告知。

这里都采用7×7模板进行滤波。滤波结果图如下所示。

传统中值滤波算法

更新式中值滤波算法

可以看到更新式中值滤波算法滤除了更多的高频信息，如图片左侧中部的梯田：传统中值滤波算法导致了不自然的高频横纹；而更新式滤波算法并没有相应的纹路。而且更新式中值滤波算法得到的图像 {MOD}块偏浓重，略有野兽派的画风；传统中值滤波算法在近景出还有一定的写实特点。

附代码

脚本

1、Median_Filter_Script.m

clc clear close all %% 定义 % 滤波器模板大小 s = 5; % 取值位置 % 注（中值滤波为0.5，最大值滤波为1，最小值滤波为0） p = 0.5; % 更新式滤波 % 1：非更新式滤波，0：更新式滤波 method = 1; %% 图像输入 im = imread('circuit_board.tif'); im = im(:,:,1); %% 图像扩展 % 根据滤波器模板大小定义扩展图像 % Define augmented figure according to the filter mask scale [H,W] = size(im); [im_augmented, index_x, index_y] = Image_Augment(im, s); %% 图像滤波 % 注意，如果采用不更新方式进行滤波，出来的效果貌似更好 im_filtered_augmented = Median_Filter(im_augmented,s,p,method); im_filtered = im_filtered_augmented(index_x,index_y); %% 图像显示 figure subplot(1,2,1 ) imshow(im); subplot(1,2,2) imshow(uint8(im_filtered));

2、Color_Image_Median_Filter_Script.m

clc clear close all %% 定义 % 滤波器模板大小 s = 7; % 取值位置 % 注（中值滤波为0.5，最大值滤波为1，最小值滤波为0） p = 0.5; % 更新式滤波 % 1：非更新式滤波，0：更新式滤波 method = 0; %% 图像输入 im = imread('test.jpg'); %% 图像增广 [im_R_augmented,index_R_x,index_R_y] = Image_Augment(im(:,:,1),s); [im_G_augmented,index_G_x,index_G_y] = Image_Augment(im(:,:,2),s); [im_B_augmented,index_B_x,index_B_y] = Image_Augment(im(:,:,3),s); %% 滤波 im_filtered_R_augmented = Median_Filter(im_R_augmented,s,p,method); im_filtered_G_augmented = Median_Filter(im_G_augmented,s,p,method); im_filtered_B_augmented = Median_Filter(im_B_augmented,s,p,method); im_filtered = cat(3,... im_filtered_R_augmented(index_R_x,index_R_y),... im_filtered_G_augmented(index_G_x,index_G_y),... im_filtered_B_augmented(index_B_x,index_B_y)); %% 显示图像 imshow(uint8(im_filtered));

函数

3、Image_Augment.m

function [im_augmented,index_h,index_w] = Image_Augment(im,s) % Image_Augment - 增广图像，以便于图像模板匹配 % 输入参数 % im - 待增广图像 % s - 模板边长 % 输出参数 % im_augmented - 增广后图像 % index_h,index_w - 图像在高和宽的嵌入位置索引 [H,W] = size(im); index_h = s:(s+H-1); index_w = s:(s+W-1); im_augmented = zeros(H + 2*(s-1),... W + 2*(s-1)); im_augmented(index_h,index_w) = double(im); im_augmented(index_h,1:(s-1)) = repmat(... im_augmented(index_h,s),1,s-1); im_augmented(index_h,(s+W):(W+2*(s-1))) = repmat(... im_augmented(index_h,s+W-1),1,s-1); im_augmented(1:(s-1),index_w) = repmat(... im_augmented(s,index_w),s-1,1); im_augmented((s+H):(H+2*(s-1)),index_w) = repmat(... im_augmented(s+H-1,index_w),s-1,1); for i = 1:(s-1) for j = 1:(s-1) if i > j im_augmented(i,j) = im_augmented(s,j); im_augmented(H+2*(s-1)-i+1,j) = ... im_augmented(s+H-1,j); im_augmented(i,W+2*(s-1)-j+1) = ... im_augmented(s,W+2*(s-1)-j+1); im_augmented(H+2*(s-1)-i+1,W+2*(s-1)-j+1) = ... im_augmented(s+H-1,W+2*(s-1)-j+1); elseif i < j im_augmented(i,j) = im_augmented(i,s); im_augmented(H+2*(s-1)-i+1,j) = ... im_augmented(H+2*(s-1)-i+1,s); im_augmented(i,W+2*(s-1)-j+1) = ... im_augmented(i,s+W-1); im_augmented(H+2*(s-1)-i+1,W+2*(s-1)-j+1) = ... im_augmented(H+2*(s-1)-i+1,s+W-1); else im_augmented(i,j) = im_augmented(s,s); im_augmented(H+2*(s-1)-i+1,j) = ... im_augmented(s+H-1,s); im_augmented(i,W+2*(s-1)-j+1) = ... im_augmented(s,s+W-1); im_augmented(H+2*(s-1)-i+1,W+2*(s-1)-j+1) = ... im_augmented(s+H-1,s+W-1); end end end end

4、Median_Filter.m

function im_filtered = Median_Filter(im,s,ratio,method) % Median_Filter - 图像中值滤波 % 输入参数 % im - 待滤波图像 % s - 滤波模板大小 % ratio - 滤波阈值 % method - 滤波方法 % 1：非更新式滤波，0：更新式滤波 % 注：非更新式滤波最好只采用0.5滤波阈值，采用其他阈值会出现不理想的效果 % 输出参数 % im_filtered - 滤波后图像 %% 基本参数求解 [H,W] = size(im); %% 图像增广 [im_augmented,index_x,index_y] = Image_Augment(im,s); im_augmented_temp = im_augmented; for i = (s-1):(H+2*(s-1)-(s-2)) for j = (s-1):(W+2*(s-1)-(s-2)) if method temp = im_augmented_temp(... (i-(s-2)):(i+(s-2)),(j-(s-2)):(j+(s-2))); else temp = im_augmented(... (i-(s-2)):(i+(s-2)),(j-(s-2)):(j+(s-2))); end temp = sort(temp(:)); index = ceil(length(temp)*ratio); if index == 0 index = 1; end im_augmented(i,j) = temp(index); end end im_filtered = im_augmented(index_x,index_y); end