Matlab实现medfilt2函数功能

算法分析

1. 对原图像矩阵边界填充（填充的行列根据滤波器模板大小来设置），这里使用3*3的模板，所以向外扩充一圈，行和列增加2，这里使用的是复制边界的填充方式'replicate'(直接调用的padarray函数)，也可以直接通过矩阵赋值的方法实现边界填充

2. 遍历原图像每个通道的每个像素点，取扩充边界后图像中滤波器大小的像素点个数，求得中值，将中值赋给当前遍历的像素点

3. 输出图像的大小是新开辟的和原图像大小相等的矩阵，求得的中值应该赋值给新矩阵中每个元素

function outputimg = mymedfilt2(A)
% 利用中值滤波去噪
% 参数A是加噪的图像矩阵
% outputimg是输出中值滤波后的图像

[H,W,CH] = size(A);
N = 3;  %设置中值滤波器的大小为3*3
M = zeros(H+2,W+2,CH);

%对彩色图的每一个通道都进行边界填充
for i = 1 : CH   
    %直接调用padarray，在原图的每一维度的第一个元素和最后一个元素后，以复制边界的方式进行填充
    M(:,:,i) = padarray(A(:,:,i),[1,1],'replicate'); 
end

%矩阵赋值法，复制边界型填充
% M(2:height+1,2:width+1) = A(1:height,1:width);%原矩阵直接复制到新矩阵中心
% %使用直接复制最边界的层
% M(1,2:width+1) = A(1,:);  %将M的第一行
% M(height + 2,2:width+1) = A(height,:);  %M的最后一行
% M(2:height+1,1) = A(:,1);  %M的第一列
% M(2:height+1,width+2) = A(:,width);  %M的最后一列

outputimg = zeros(H,W,CH);
for k = 1 : CH
   for i = 1 : H  
       for j = 1 : W  
           c = M(i:i + (N-1),j:j + (N-1),k); %在C中从头取模板大小的块赋给c  
           e = c(1,:);      %e中存放是c矩阵的第一行  
           for u = 2:N  %  将c中的其他行元素取出来接在e后使e为一个行矩阵 
               e = [e,c(u,:)];          
           end  
           med = median(e);  %n*n这几个像素点的中值  
           outputimg(i,j,k) = med;   %将模板各元素的中值赋给模板中心位置的元素  
       end   
   end
end    

outputimg = uint8(outputimg);
end

实验结果

% 调用示例：
I = imread('cameraman.tif');
M = imnoise(I,'salt & pepper');
B = mymedfilt2(I);
subplot(1,2,1),imshow(M),title('加椒盐噪声');
subplot(1,2,2),imshow(B),title('中值滤波');

% 调用示例：
I = imread('cameraman.tif');
M = imnoise(I,'gaussian',0,0.02);
B = mymedfilt2(M);
subplot(1,2,1),imshow(M),title('加高斯噪声');
subplot(1,2,2),imshow(B),title('中值滤波');

实验分析

• 由实验结果可见，中值滤波是一个非线性滤波。它对椒盐噪声这类随机出现的噪点有比较好的平滑效果，但对于线性的噪声（如高斯噪声）效果不佳

• 如下图，使用Matlab自带的中值滤波函数medfilt2，生成的图像的四个边缘会有黑点,原因是直接调用，默认是边界0填充，当加入椒盐噪声时，会因求得的中值还是0，而得到黑点

  B = medfilt2(M);

• 解决办法，增加参数PADOPT控制边界填充的方式，如果PADOPT为' 0 '，则填充边界处都是0。如果PADOPT是‘symmetric’，则是对称地在边界上扩展的。如果PADOPT为'indexed'，如果A是双值，则用1填充;否则它被0填充

B = medfilt2(M,[3,3],'symmetric');