errordiff.m

在二值化前,将像素与阈值的误差扩散到周围像素上,提高图像的二值化质量

%******************************************************
% filename:ErrorDiff.m
% function:在二值化前,将像素与阈值的误差扩散到周围像素上
% date:2006.03.23
% programmer:wt
%******************************************************
function ErrorDiff()
clc,clear

I=imread('E:\work\picture\L0.jpg');
if isrgb(I)
    I=rgb2gray(I);
end
I=double(I);
% 选取原始图像的均值作误差扩散的阈值
t=mean(mean(I));
[m n]=size(I);
% 误差扩散到右,左下,下,右下四个像素上,比例如下
%    * 7
%  3 5 1
t1=16;
I1=Diff_16(I,m,n,t,t1);

% 误差扩散比例
%     * 7 5
% 3 5 7 5 3
% 1 3 5 3 1
t2=48;
I2=Diff_48(I,m,n,t,t2);

% 误差扩散比例
%     * 8 4
% 2 4 8 4 2
% 1 2 4 2 1
t3=42;
I3=Diff_42(I,m,n,t,t3);
% 误差扩散比例
%          **    32
% 12    26    30    16
%    12    26    12
%  5    12    12     5
t4=200;
I4=Diff_200(I,m,n,t,t4);
% 显示图像
%imshow(uint8(I),[]);
figure,imshow(uint8(I1),[]),title('1975:Floyd&Steinberg');
imwrite(logical(I1),'E:\work\picture\L0_16.bmp');
figure,imshow(uint8(I2),[]),title('1976:Jarvis,Judice&Ninke');
imwrite(logical(I2),'E:\work\picture\L0_48.bmp');
figure,imshow(uint8(I3),[]);title('1981:Stucki');
imwrite(logical(I3),'E:\work\picture\L0_42.bmp');
figure,imshow(uint8(I4),[]);title('1985:Stevenson&Arce');
imwrite(logical(I4),'E:\work\picture\L0_200.bmp');
return
%****************************************
%sub1 图像二值化
%****************************************
function bw=ImageToBw(I,t)
Pixelhigh=find(I>=t);
Pixellow=find(I<t);
I(Pixelhigh)=1;
I(Pixellow)=0;
bw=I;
return
%****************************************
%sub2 误差扩散系数为16
%****************************************
function  A=Diff_16(I,m,n,t,t1)
col=zeros(m,1);                     %图像矩阵的左右各加一列,下加一行
row=zeros(1,n+2);
I=[col,I,col];
I=[I;row];
[m1 n1]=size(I);
for i=1:m1-1
    for j=2:n1-1
        I(i,j)=cut(I(i,j));
      if I(i,j)>=t                      %分为两种情况,大于阈值时
        e=(I(i,j)-255)/t1;
        I(i,j+1)=I(i,j+1)+7*e;
        I(i+1,j-1)=I(i+1,j-1)+3*e;
        I(i+1,j)=I(i+1,j)+5*e;
        I(i+1,j+1)=I(i+1,j+1)+e;
      else
        e=I(i,j)/t1;                       %小于阈值时
        I(i,j+1)=I(i,j+1)+7*e;
        I(i+1,j-1)=I(i+1,j-1)+3*e;
        I(i+1,j)=I(i+1,j)+5*e;
        I(i+1,j+1)=I(i+1,j+1)+e;
    end
end
end
I=I(1:m1-1,2:n1-1);    % 舍去增加的行和列
A=ImageToBw(I,t);
return
%****************************************
%sub3 误差扩散系数为48
%****************************************
function B=Diff_48(I,m,n,t,t2)
col=zeros(m,2);                     %图像矩阵的左右各加两列,下加两行
row=zeros(2,n+4);
I=[col,I,col];
I=[I;row];
[m1 n1]=size(I);
for i=1:m1-2
    for j=3:n1-2
        I(i,j)=cut(I(i,j));
      if I(i,j)>=t                      %分为两种情况,大于阈值时
        e=(I(i,j)-255)/t2;
        I(i,j+1)=I(i,j+1)+7*e;
        I(i,j+2)=I(i,j+2)+5*e;
        I(i+1,j-2)=I(i+1,j-2)+3*e;
        I(i+1,j-1)=I(i+1,j-1)+5*e;
        I(i+1,j)=I(i+1,j)+7*e;
        I(i+1,j+1)=I(i+1,j+1)+5*e;
        I(i+1,j+2)=I(i+1,j+2)+3*e;
        I(i+2,j-2)=I(i+2,j-2)+e;
        I(i+2,j-1)=I(i+2,j-1)+3*e;
        I(i+2,j)=I(i+2,j)+5*e;
        I(i+2,j+1)=I(i+2,j+1)+3*e;
        I(i+2,j+2)=I(i+2,j+2)+e;
      else
        e=I(i,j)/t2;                       %小于阈值时
        I(i,j+1)=I(i,j+1)+7*e;
        I(i,j+2)=I(i,j+2)+5*e;
        I(i+1,j-2)=I(i+1,j-2)+3*e;
        I(i+1,j-1)=I(i+1,j-1)+5*e;
        I(i+1,j)=I(i+1,j)+7*e;
        I(i+1,j+1)=I(i+1,j+1)+5*e;
        I(i+1,j+2)=I(i+1,j+2)+3*e;
        I(i+2,j-2)=I(i+2,j-2)+e;
        I(i+2,j-1)=I(i+2,j-1)+3*e;
        I(i+2,j)=I(i+2,j)+5*e;
        I(i+2,j+1)=I(i+2,j+1)+3*e;
        I(i+2,j+2)=I(i+2,j+2)+e;
    end
end
end
I=I(1:m1-2,3:n1-2);    % 舍去增加的行和列
B=ImageToBw(I,t);
return
%****************************************
%sub4 误差扩散系数为42
%****************************************
function C=Diff_42(I,m,n,t,t3)
col=zeros(m,2);                     %图像矩阵的左右各加两列,下加两行
row=zeros(2,n+4);
I=[col,I,col];
I=[I;row];
[m1 n1]=size(I);
for i=1:m1-2
    for j=3:n1-2
        I(i,j)=cut(I(i,j));
      if I(i,j)>=t                      %分为两种情况,大于阈值时
        e=(I(i,j)-255)/t3;
        I(i,j+1)=I(i,j+1)+8*e;
        I(i,j+2)=I(i,j+2)+4*e;
        I(i+1,j-2)=I(i+1,j-2)+2*e;
        I(i+1,j-1)=I(i+1,j-1)+4*e;
        I(i+1,j)=I(i+1,j)+8*e;
        I(i+1,j+1)=I(i+1,j+1)+4*e;
        I(i+1,j+2)=I(i+1,j+2)+2*e;
        I(i+2,j-2)=I(i+2,j-2)+e;
        I(i+2,j-1)=I(i+2,j-1)+2*e;
        I(i+2,j)=I(i+2,j)+4*e;
        I(i+2,j+1)=I(i+2,j+1)+2*e;
        I(i+2,j+2)=I(i+2,j+2)+e;
      else
        e=I(i,j)/t3;                       %小于阈值时
        I(i,j+1)=I(i,j+1)+8*e;
        I(i,j+2)=I(i,j+2)+4*e;
        I(i+1,j-2)=I(i+1,j-2)+2*e;
        I(i+1,j-1)=I(i+1,j-1)+4*e;
        I(i+1,j)=I(i+1,j)+8*e;
        I(i+1,j+1)=I(i+1,j+1)+4*e;
        I(i+1,j+2)=I(i+1,j+2)+2*e;
        I(i+2,j-2)=I(i+2,j-2)+e;
        I(i+2,j-1)=I(i+2,j-1)+2*e;
        I(i+2,j)=I(i+2,j)+4*e;
        I(i+2,j+1)=I(i+2,j+1)+2*e;
        I(i+2,j+2)=I(i+2,j+2)+e;
    end
end
end
I=I(1:m1-2,3:n1-2);    % 舍去增加的行和列
C=ImageToBw(I,t);
return
%****************************************
%sub5 误差扩散系数为200
%****************************************
function D=Diff_200(I,m,n,t,t4)
col=zeros(m,3);                     %图像矩阵的左右各加三列,上下各加三行
row=zeros(3,n+6);
I=[col,I,col];
I=[I;row];
[m1 n1]=size(I);
for i=1:m1-3
    for j=4:n1-3
      I(i,j)=cut(I(i,j));
      if I(i,j)>=t                      %分为两种情况,大于阈值时
        e=(I(i,j)-255)/t4;
        I(i,j+2)=I(i,j+2)+32*e;
        I(i+1,j-3)=I(i+1,j-3)+12*e;
        I(i+1,j-1)=I(i+1,j-1)+26*e;
        I(i+1,j+1)=I(i+1,j+1)+30*e;
        I(i+1,j+3)=I(i+1,j+3)+16*e;
        I(i+2,j-2)=I(i+2,j-2)+12*e;
        I(i+2,j)=I(i+2,j)+26*e;
        I(i+2,j+2)=I(i+2,j+2)+12*e;
        I(i+3,j-3)=I(i+3,j-3)+5*e;
        I(i+3,j-1)=I(i+3,j-1)+12*e;
        I(i+3,j+1)=I(i+3,j+1)+12*e;
        I(i+3,j+3)=I(i+3,j+3)+5*e;
      else
        e=I(i,j)/t4;                       %小于阈值时
        I(i,j+2)=I(i,j+2)+32*e;
        I(i+1,j-3)=I(i+1,j-3)+12*e;
        I(i+1,j-1)=I(i+1,j-1)+26*e;
        I(i+1,j+1)=I(i+1,j+1)+30*e;
        I(i+1,j+3)=I(i+1,j+3)+16*e;
        I(i+2,j-2)=I(i+2,j-2)+12*e;
        I(i+2,j)=I(i+2,j)+26*e;
        I(i+2,j+2)=I(i+2,j+2)+12*e;
        I(i+3,j-3)=I(i+3,j-3)+5*e;
        I(i+3,j-1)=I(i+3,j-1)+12*e;
        I(i+3,j+1)=I(i+3,j+1)+12*e;
        I(i+3,j+3)=I(i+3,j+3)+5*e;
    end
end
end
I=I(1:m1-3,4:n1-3);    % 舍去增加的行和列
D=ImageToBw(I,t);
return
%***********************************
%sub6 限幅
%***********************************
function E=cut(ee)
if ee>255
    ee=255; 
 elseif ee<0
    ee=0;
end
 E=ee;
return