基于matlab的图像阈值分割算法---参数法（自动阈值选择）

1 迭代法

迭代法的的设计思想是，开始时选择一个阈值作为初始估计值，然后按某种策略不断的改进这一估计值，直到满足给定的准则为止。在迭代过程中，关键之处在于选择什么样的阈值改进策略。好的改进策略应该具备两个特征：一是能够快速收敛，二是在每一个迭代过程中，新产生的阈值优于上一次的阈值。

代码：

%迭代法

clc;
clear all;
f=imread('H:\数字图像处理\图片\标准图片集\baboon.bmp');
f=rgb2gray(f);
f=im2double(f);
T=0.5*(min(f(:))+max(f(:)));
done=false;
while ~done
g=f>=T;
Tn=0.5*(mean(f(g))+mean(f(~g)));
done=abs(T-Tn)<0.1;
T=Tn;
end
T
r=im2bw(f,T);
figure,imshow(f),title('原图');
figure,imshow(r),title('迭代法');

2 Otsu

该算法是日本人Otsu提出的一种动态阈值分割算法。它的主要思想是按照灰度特性将图像划分为背景和目标2部分，划分依据为选取门限值，使得背景和目标之间的方差最大。（背景和目标之间的类间方差越大，说明这两部分的差别越大，当部分目标被错划分为背景或部分背景错划分为目标都会导致这两部分差别变小。因此，使用类间方差最大的分割意味着错分概率最小。）这是该方法的主要思路。

代码：

%Ostu
clc;
clear all;
I=imread('H:\数字图像处理\图片\标准图片集\baboon.bmp');
I=rgb2gray(I);
I=double(I);
figure,imshow(uint8(I)),title('原图');
[m,n]=size(I);
Th=Otsu(I);
Th
for i=1:m
for j=1:n
if I(i,j)>=Th
I(i,j)=255;
else
I(i,j)=0;
end
end
end
figure,imshow(I),title('Otsu法');

function T=Otsu(I)
[m,n]=size(I);
I=double(I);
count=zeros(256,1);
pcount=zeros(256,1);
for i=1:m
for j=1:n
pixel=I(i,j);
count(pixel+1)=count(pixel+1)+1;
end
end
dw=0;
for i=0:255
pcount(i+1)=count(i+1)/(m*n);
dw=dw+i*pcount(i+1);
end
Th=0;
Thbest=0;
dfc=0;
dfcmax=0;
while (Th>=0 & Th<=255)
dp1=0;
dw1=0;
for i=0:Th
dp1=dp1+pcount(i+1);
dw1=dw1+i*pcount(i+1);
end
if dp1>0
dw1=dw1/dp1;
end
dp2=0;
dw2=0;
for i=Th+1:255
dp2=dp2+pcount(i+1);
dw2=dw2+i*pcount(i+1);
end
if dp2>0
dw2=dw2/dp2;
end
dfc=dp1*(dw1-dw)^2+dp2*(dw2-dw)^2;
if dfc>=dfcmax
dfcmax=dfc;
Thbest=Th;
end
Th=Th+1;
end
T=Thbest;

3 均匀性度量法

均匀性度量法的设计思想是，假设当图像被分为目标物和背景两个类别时，属于同一类别内的像素值分布应该具有均匀性。在这里采用方差来度量像素间的均匀性。

代码：

clc;
clear all;
I=imread('H:\数字图像处理\图片\标准图片集\baboon.bmp');
I=rgb2gray(I);
I=double(I);
[m,n]=size(I);
Smin=-1;
for T=0:255
sum1=0; num1=0;
sum2=0; num2=0;
for i=1:m
for j=1:n
if I(i,j)>=T
sum2=sum2+I(i,j);
num2=num2+1;
else
sum1=sum1+I(i,j);
num1=num1+1;
end
end
end
ave1=sum1/num1;
ave2=sum2/num2;
d1=-1;
d2=-1;
for i=1:m
for j=1:n
if I(i,j)>=T
d=(I(i,j)-ave2)^2;
if d2==-1
d2=d;
else
d2=d2+d;
end
else
d=(I(i,j)-ave1)^2;
if d1==-1
d1=d;
else
d1=d1+d;
end
end
end
end
p1=num1/(m*n);
p2=num2/(m*n);
S=p1*d1+p2*d2;
if(Smin==-1)
Smin=S;
else
if(S)
Smin=S;
Th=T;
end
end
end
Th
figure,imshow(uint8(I)),title('原图');
for i=1:m
for j=1:n
if I(i,j)>=Th
I(i,j)=255;
else
I(i,j)=0;
end
end
end
figure,imshow(I),title('均匀性度量法');

4 类间最大距离法

类间最大距离法德设计思想是，在某个适当的阈值下，图像分割后的前景目标与背景两个类之间的差异最大为最佳阈值。在这里两个类别（目标与背景）的差异，用两个类别中心与阈值之间的距离差来度量。

代码：

clc;
clear all;
I=imread('H:\数字图像处理\图片\标准图片集\baboon.bmp');
I=rgb2gray(I);
I=double(I);
[m,n]=size(I);
Smax=0;
for T=0:255
sum1=0; num1=0;
sum2=0; num2=0;
for i=1:m
for j=1:n
if I(i,j)>=T
sum2=sum2+I(i,j);
num2=num2+1;
else
sum1=sum1+I(i,j);
num1=num1+1;
end
end
end
ave1=sum1/num1;
ave2=sum2/num2;
S=((ave2-T)*(T-ave1))/(ave2-ave1)^2;
if(S>Smax)
Smax=S;
Th=T;
end
end
Th
figure,imshow(uint8(I)),title('原图');
for i=1:m
for j=1:n
if I(i,j)>=Th
I(i,j)=255;
else
I(i,j)=0;
end
end
end
figure,imshow(I),title('类间最大距离法');

5 最大熵法

熵是信息论中对不确定性的度量，是对数据中所包含信息量大小的度量。熵取最大值时，就表明获得的信息量为最大。最大熵法的设计思想是，选择适当的阈值将图像分为两类，两类的平均熵之和为最大时，可以从图像中获得最大信息量，以此来确定最佳阈值。

代码：

clc;
clear all;
I=imread('H:\数字图像处理\图片\标准图片集\baboon.bmp');
I=rgb2gray(I);
figure,imshow(I);title('原始图像');
h=imhist(I);
h1=h;
len=length(h);
[m,n]=size(I);
h1=h1/(m*n);
for i=1:(len-1)
if h(i)~=0
P1=sum(h1(1:i));
P2=sum(h1((i+1):len));
else
continue;
end
H1(i)=-(sum(P1.*log(P1)));
H2(i)=-(sum(P2.*log(P2)));
H(i)=H1(i)+H2(i);
end
m1=max(H);
Th=find(H==m1);
Th
for i=1:m
for j=1:n
if I(i,j)>=Th
I(i,j)=255;
else
I(i,j)=0;
end
end
end
figure,imshow(I);title('最大熵法');

6 最大类间、类内方差比法

从统计意义上讲，方差是表征数据分布不均衡性的统计量，要通过阈值对两类问题进行分割，显然，适当的阈值使得两类数据间的方差越大越好，表明该阈值的确将两类不同的区域分开了，同时希望属于同一类的方差越小越好，表明同一类区域有一定的相似性。因此可以采用类内，类间方差比作为选择阈值的评价参数。

代码：

clc;
clear all;
f=imread('H:\数字图像处理\图片\标准图片集\baboon.bmp');
imshow(uint8(f)),title('彩图');
I=rgb2gray(f);
I=double(I);
[m,n]=size(I);
Smax=-1;
for T=0:255
sum1=0; num1=0;
sum2=0; num2=0;
for i=1:m
for j=1:n
if I(i,j)>=T
sum2=sum2+I(i,j);
num2=num2+1;
else
sum1=sum1+I(i,j);
num1=num1+1;
end
end
end
ave1=sum1/num1;
ave2=sum2/num2;
ave=(sum1+sum2)/(m*n);
d1=-1;
d2=-1;
for i=1:m
for j=1:n
if I(i,j)>=T
d=(I(i,j)-ave2)^2;
if d2==-1
d2=d;
else
d2=d2+d;
end
else
d=(I(i,j)-ave1)^2;
if d1==-1
d1=d;
else
d1=d1+d;
end
end
end
end
p1=num1/(m*n);
p2=num2/(m*n);
S1=p1*(ave1-ave)^2+p2*(ave2-ave)^2;
S2=p1*d1+p2*d2;
S=S1/S2;
if S>Smax
Smax=S;
Th=T;
end
end
Th
figure,imshow(uint8(I)),title('原图');
for i=1:m
for j=1:n
if I(i,j)>=Th
I(i,j)=255;
else
I(i,j)=0;
end
end
end
figure,imshow(I),title('最大类内类间方差比法');

7 局部阈值法

代码：

clc;
clear all;
I=imread('H:\数字图像处理\图片\标准图片集\baboon.bmp');
I=rgb2gray(I);
I=im2double(I);
figure,imshow(I),title('原图');
se=strel('disk',10);
ft=imtophat(I,se);
Th=graythresh(ft);
Th
G=im2bw(ft,Th);
figure,imshow(G),title('局部阈值');

至于那种方法最好，这也没有一个统一的评判标准，大多是取决于实际需求和人的主观评判。