对图像进行中值滤波的几种写法以及基于排序阈值的开关中值滤波算法

自己编写的中值滤波算法和库函数调用的中值滤波算法 psnr值一样，因此不必多虑，放心调用.

以下文件：

1:自己编写的中值滤波

2:库函数中值滤

3：极值中值算法：中心像素等于窗口最大值或者最小值，即该点为噪声点的可能性很大，中心点则使用中值代替

4：一种基于排序阈值的开关中值滤波方法（秦鹏等）：先将窗口内所有像素分成三类：1） 噪声点 2） 平坦区域 3）边缘细节

5：极大中值滤波：水平、垂直、45°、135°方向各自取中值，再对4个值取最大值

6：med(水平中值，90°中值，中心点) ，中心点 ，med（垂直，135°，中心点），再对3个取中值

%commonfilt2_1.m
%自己编写的中值滤波
function [y]=commonfilt2_1(x,N)
y=x;
k=floor(N/2);%中心点旁边对称点的数目
k0=floor( (N*N)/2 )+1;
[m,n]=size(x);

for i=k+1:m-k
for j=k+1:n-k
z=x(i-k:i+k,j-k:j+k);
z1=z(1:end);%矩阵z转化为行向量z1
%a*a个数冒泡排序法
for i1=1:1:(N^2-1)
for j1=1:1:(N^2-i1)
if z1(1,j1) > z1(1,j1+1)
z0=z1(1,j1);
z1(1,j1)=z1(1,j1+1);
z1(1,j1+1)=z0;
end
end
end
y(i,j)=z1(1,k0);%z1中的中值
end
end

%commonfilt2_2.m
%库函数中值滤波
function [y]=commonfilt2_2(x,N)
y=x;
k=floor(N/2);%中心点旁边对称点的数目
[m,n]=size(x);
for i=k+1:m-k
for j=k+1:n-k
z=x(i-k:i+k,j-k:j+k);
y(i,j)=median(z(:));
end
end

%commonfilt2_3.m
%极值中值算法：中心像素等于窗口最大值或者最小值，即该点为噪声点的可能性很大，中心点则使用中值代替
function [y]=commonfilt2_3(x,N)
y=x;
k=floor(N/2);%中心点旁边对称点的数目
[m,n]=size(x);
for i=k+1:m-k
for j=k+1:n-k
z=x(i-k:i+k,j-k:j+k);
zmin=min(z(:));
zmax=max(z(:));
zmed=median(z(:));
if  x(i,j) == zmin | x(i,j) == zmax
y(i,j) = zmed;
end
%中心灰度值介于最大值最小值之间，则不变
end
end

%commonfilt2_4.m
%基于排序阈值的开关中值滤波算法：先将窗口内所有像素分成三类：1 噪声点 2 平坦区域 3边缘细节
%第一和第三要中值取代中心点灰度值
function [y]=commonfilt2_4(x,N)
[m,n]=size(x);
y=x;
k=floor(N/2);%中心点旁边对称点的数目
%按照文献取参数
a=5;b=7;
TA=20;TB=4;

for i=k+1:m-k
for j=k+1:n-k
z=x(i-k:i+k,j-k:j+k);
z1=z(1:end);%将矩阵z转化为列向量
z1=sort(z1);%对列向量z1排序，从小到大,1:N^2
zmin=min(z(:));
zmax=max(z(:));
zmed=median(z(:));
%若中心点处于最小值 最大值之间，则灰度值不变
%若中心点等于最小值，即疑似噪声点
if x(i,j) == zmin
if z1(1,1+a)-x(i,j) >= TA
y(i,j) = zmed;
else
if abs( z1(1,1+b)-z1(1,N^2-b) )<= TB %平坦区
y(i,j) = zmed;
%否则的话则为边缘细节区，灰度值不变
end
end
end
%若中心点等于最大值，即疑似噪声点
if x(i,j) == zmax
if x(i,j)-z1(1,N^2-a) >= TA
y(i,j) = zmed;
else
if abs( z1(1,1+b)-z1(1,N^2-b) )<= TB %平坦区
y(i,j) = zmed;
%否则的话则为边缘细节区，灰度值不变
end
end
end
end
end

%commonfilt2_5.m
%极大中值滤波：水平、垂直、45°、135°方向各自取中值，再对4个值取最大值
function [y]=commonfilt2_5(x,N)
y=x;
k=floor(N/2);%中心点旁边对称点的数目
[m,n]=size(x);
for i=k+1:m-k
for j=k+1:n-k
%z=x(i-k:i+k,j-k:j+k);
%y(i,j)=median(z(:));
z1=zeros(1,2*k+1);
z2=zeros(1,2*k+1);
z3=zeros(1,2*k+1);
z4=zeros(1,2*k+1);
k1=1;k2=1;k3=1;k4=1;
%垂直方向上
for i1=i;
for j1=j-k:j+k
z1(1,k1) = x(i1,j1);
k1=k1+1;
end
end
%水平方向上
for j1=j;
for i1=i-k:i+k
z2(1,k2) = x(i1,j1);
k2=k2+1;
end
end
%45°方向
for i1=-k:k
z3(1,k3) = x(i+i1,j+i1);
k3=k3+1;
end
%135°方向
for i1=-k:k
z3(1,k4) = x(i+i1,j-i1);
k4=k4+1;
end
%求各个矩阵的中值
zmed1=median(z1);%垂直方向
zmed2=median(z2);%水平方向
zmed3=median(z3);%45°方向
zmed4=median(z4);%135°方向
y(i,j)=max([zmed1,zmed2,zmed3,zmed4]);
end
end

%commonfilt2_6.m
%多重中值滤波
%med(水平中值，90°中值，中心点) ，中心点 ，med（垂直，135°，中心点），再对3个取中值
function [y]=commonfilt2_6(x,N)
y=x;
k=floor(N/2);%中心点旁边对称点的数目
[m,n]=size(x);
for i=k+1:m-k
for j=k+1:n-k
%z=x(i-k:i+k,j-k:j+k);
%y(i,j)=median(z(:));
z1=zeros(1,2*k+1);
z2=zeros(1,2*k+1);
z3=zeros(1,2*k+1);
z4=zeros(1,2*k+1);
k1=1;k2=1;k3=1;k4=1;
%垂直方向上
for i1=i;
for j1=j-k:j+k
z1(1,k1) = x(i1,j1);
k1=k1+1;
end
end
%水平方向上
for j1=j;
for i1=i-k:i+k
z2(1,k2) = x(i1,j1);
k2=k2+1;
end
end
%45°方向
for i1=-k:k
z3(1,k3) = x(i+i1,j+i1);
k3=k3+1;
end
%135°方向
for i1=-k:k
z3(1,k4) = x(i+i1,j-i1);
k4=k4+1;
end
%求各个矩阵的中值
zmed1=median(z1);%垂直方向
zmed2=median(z2);%水平方向
zmed3=median(z3);%45°方向
zmed4=median(z4);%135°方向
y(i,j)=median( [median( [x(i,j),zmed2,zmed1 ]),median( [x(i,j),zmed3,zmed4 ]),x(i,j) ]);
end
end

I=imread('couple.png');
%I=rgb2gray(I);
J=imnoise(I,'salt & pepper',0.2);

k1=commonfilt2_1(J,5);%自己编写的中值滤波方法
k2=commonfilt2_2(J,5);%库函数中值滤波方法
k3=commonfilt2_3(J,5);%极值中值滤波方法
k4=commonfilt2_4(J,5);%极值中值滤波方法
%使用第四种方法进行迭代，处理密度较大的噪声很有效
%k4=commonfilt2_4(k4,5);%极值中值滤波方法
%k4=commonfilt2_4(k4,5);%极值中值滤波方法
%k4=commonfilt2_4(k4,5);%极值中值滤波方法
k5=commonfilt2_5(J,5);%极大中值滤波器
k6=commonfilt2_6(J,5);%多重中值滤波器

subplot(231),imshow( uint8(k1) ),title('自己编写的中值滤波');
subplot(232),imshow( uint8(k2)),title('库函数中值滤波');
subplot(233),imshow( uint8(k3)),title('极值中值滤波方法');
subplot(234),imshow( uint8(k4)),title('基于排序阈值的开关中值滤波');
subplot(235),imshow( uint8(k5)),title('极大中值滤波器');
subplot(236),imshow( uint8(k6)),title('多重中值滤波器');

% 计算三种算法的峰值信噪比
[e f]=size(I);
B=8;                %编码一个像素用多少二进制位
MAX=2^B-1;          %图像有多少灰度级

I=double(I);
k1=double(k1);
k2=double(k2);
k3=double(k3);
k4=double(k4);
k5=double(k5);
k6=double(k6);

%%%%%% psnr1=
MES1=sum(sum((I-k1).^2))/(e*f);
PSNR1=20*log10(MAX/sqrt(MES1));

%%%%%% psnr2=
MES2=sum(sum((I-k2).^2))/(e*f);
PSNR2=20*log10(MAX/sqrt(MES2));

%%%%%% psnr3=
MES3=sum(sum((I-k3).^2))/(e*f);
PSNR3=20*log10(MAX/sqrt(MES3));

%%%%%% psnr4=
MES4=sum(sum((I-k4).^2))/(e*f);
PSNR4=20*log10(MAX/sqrt(MES4));

%%%%%% psnr5=
MES5=sum(sum((I-k5).^2))/(e*f);
PSNR5=20*log10(MAX/sqrt(MES5));

%%%%%% psnr6=
MES6=sum(sum((I-k6).^2))/(e*f);
PSNR6=20*log10(MAX/sqrt(MES6));