要实现一个核半径无关计算复杂度的排序滤波器其实非常简单,所利用的不过是像素值的离散性和有限性而已。简单地说,就是滤波窗口内的像素值取值范围就是0~255(对于8bit图像而言),取同一个值的这些元素之间是没有差别的,所以根本没有必要真的进行排序,使用一个直方图统计出每种灰度值出现的次数,最后从前往后累加到指定位置(对于中值滤波而言就是大于等于一半元素位置)的个数时停止,此时对应bin的灰度值就是指排序的值。而进行直方图统计又有一个小技巧,那就是利用上相邻窗口之间大量的重叠元素,统计时只考虑改变的部分即可。下面给出个人实现的中值滤波matlab程序:
function fimg=fastMedianFilter(img,r)
[m,n]=size(img);
rowHist1=zeros(n,256);
rowHist2=zeros(n,256);
fimg=zeros(m,n);
%构建第一行的邻域直方图
for i=1:r+1
for j=1:r+1
rowHist1(1,img(i,j)+1)=rowHist1(1,img(i,j)+1)+1;
end
end
for j=2:n
rowHist1(j,:)=rowHist1(j-1,:);
if j<=r+1
for i=1:r+1
rowHist1(j,img(i,j+r)+1)=rowHist1(j,img(i,j+r)+1)+1;
end
elseif j<=n-r
for i=1:r+1
rowHist1(j,img(i,j-r-1)+1)=rowHist1(j,img(i,j-r-1)+1)-1;
rowHist1(j,img(i,j+r)+1)=rowHist1(j,img(i,j+r)+1)+1;
end
else
for i=1:r+1
rowHist1(j,img(i,j-r-1)+1)=rowHist1(j,img(i,j-r-1)+1)-1;
end
end
end
%求第一行滤波结果
for j=1:n
sum=0;
jlow=max(1,j-r);
jhigh=min(n,j+r);
total=(jhigh-jlow+1)*(r+1);
for k=1:256
sum=sum+rowHist1(j,k);
if sum>=total/2
fimg(1,j)=k-1;
break;
end
end
end
%计算余下行的邻域直方图
rowHist2=rowHist1;
for i=2:m
rowHist1=rowHist2;
for j=1:n
jlow=max(1,j-r);
jhigh=min(n,j+r);
for k=jlow:jhigh
if i<=r+1
rowHist1(j,img(i+r,k)+1)=rowHist1(j,img(i+r,k)+1)+1;
elseif i<=m-r
rowHist1(j,img(i-r-1,k)+1)=rowHist1(j,img(i-r-1,k)+1)-1;
rowHist1(j,img(i+r,k)+1)=rowHist1(j,img(i+r,k)+1)+1;
else
rowHist1(j,img(i-r-1,k)+1)=rowHist1(j,img(i-r-1,k)+1)-1;
end
end
end
%求中值
for j=1:n
sum=0;
ilow=max(1,i-r);
ihigh=min(m,i+r);
jlow=max(1,j-r);
jhigh=min(n,j+r);
total=(jhigh-jlow+1)*(ihigh-ilow+1);
for k=1:256
sum=sum+rowHist1(j,k);
if sum>=total/2
fimg(i,j)=k-1;
break;
end
end
end
rowHist2=rowHist1;
end
fimg=uint8(fimg);
脚本文件:
clear all;close all;clc;
img=rgb2gray(imread('lena.tiff'));
fimg=fastMedianFilter(img,10);
figure,imshow(fimg);
运行结果(lena图原图就不上了):
滤波半径参数设为30 的结果如下: