两个matlab实现最大熵法图像分割程序.doc
《两个matlab实现最大熵法图像分割程序.doc》由会员分享，可在线阅读，更多相关《两个matlab实现最大熵法图像分割程序.doc（5页珍藏版）》请在冰点文库上搜索。
%两个程序，亲测可用
clearall
a=imread（'moon.tif'）;
figure,imshow（a）
count=imhist（a）;
[m,n]=size（a）;
N=m*n;
L=256;
count=count/N;%%每一个像素的分布概率
count
fori=1:
L
ifcount（i）~=0
st=i-1;
break;
end
end
st
fori=L:
-1:
1
ifcount（i）~=0
nd=i-1;
break;
end
end
nd
f=count（st+1:
nd+1）;%f是每个灰度出现的概率
size（f）
E=[];
forTh=st:
nd-1%%%设定初始分割阈值为Th
av1=0;
av2=0;
Pth=sum（count（1:
Th+1））;
%%%第一类的平均相对熵为
fori=0:
Th
av1=av1-count（i+1）/Pth*log（count（i+1）/Pth+0.00001）;
end
%%%第二类的平均相对熵为
fori=Th+1:
L-1
av2=av2-count（i+1）/（1-Pth）*log（count（i+1）/（1-Pth）+0.00001）;
end
E（Th-st+1）=av1+av2;
end
position=find（E==（max（E）））;
th=st+position-1
fori=1:
m
forj=1:
n
ifa（i,j）>th
a（i,j）=255;
else
a（i,j）=0;
end
end
end
figure,imshow（a）;
%%%%%%%%%%%%%%%%%%%%%2-d最大熵法（递推方法）%%%%%%%%%%%
clearall;
clc;
tic
a=imread（'trial2_2.tiff'）;
figure,imshow（a）;
a0=double（a）;
[m,n]=size（a）;
h=1;
a1=zeros（m,n）;
%计算平均领域灰度的一维灰度直方图
fori=1:
m
forj=1:
n
fork=-h:
h
forw=-h:
h;
p=i+k;
q=j+w;
if（p<=0）|（p>m）
p=i;
end
if（q<=0）|（q>n）
q=j;
end
a1（i,j）=a0（p,q）+a1（i,j）;
end
end
a2（i,j）=uint8（1/9*a1（i,j））;
end
end
fxy=zeros（256,256）;
%计算二维直方图
fori=1:
m
forj=1:
n
c1=a0（i,j）;
d=double（a2（i,j））;
fxy（c1+1,d+1）=fxy（c1+1,d+1）+1;
end
end
Pxy=fxy/m/n;
%figure,
%mesh（Pxy）;
%title（'二维灰度直方图'）;
%计算Hl
Hl=0;
fori=1:
256
forj=1:
256
ifPxy（i,j）>0.00001
Hl=Hl-Pxy（i,j）*log（Pxy（i,j））;
else
Hl=Hl;
end
end
end
%计算PA,HA
PA=zeros（256,256）;
HA=zeros（256,256）;
PB=zeros（256,256）;
PA（1,1）=Pxy（1,1）;
ifPA（1,1）<1e-4
HA（1,1）=0;
else
HA（1,1）=-PA（1,1）*log（PA（1,1））;
end
fori=2:
256
PA（i,1）=PA（i-1,1）+Pxy（i,1）;
ifPxy（i,1）>0.00001
HA（i,1）=HA（i-1,1）-Pxy（i,1）*log（Pxy（i,1））;
else
HA（i,1）=HA（i-1,1）;
end
end
forj=2:
256
PA（1,j）=PA（1,j-1）+Pxy（1,j）;
ifPxy（1,j）>0.00001
HA（1,j）=HA（1,j-1）-Pxy（1,j）*log（Pxy（1,j））;
else
HA（1,j）=HA（1,j-1）;
end
end
fori=2:
256
forj=2:
256
PA（i,j）=PA（i-1,j）+PA（i,j-1）-PA（i-1,j-1）+Pxy（i,j）;
ifPxy（i,j）>0.00001
HA（i,j）=HA（i-1,j）+HA（i,j-1）-HA（i-1,j-1）-Pxy（i,j）*log（Pxy（i,j））;
else
HA（i,j）=HA（i-1,j）+HA（i,j-1）-HA（i-1,j-1）;
end
end
end
%计算最大熵
PB=1-PA;
h=zeros（256,256）;
hmax=0;
fori=1:
256
forj=1:
256
ifabs（PA（i,j））>0.00001&abs（PB（i,j））>0.00001
h（i,j）=log（PA（i,j）*PB（i,j））+HA（i,j）/PA（i,j）+（Hl-HA（i,j））/PB（i,j）;
else
h（i,j）=0;
end
ifh（i,j）>hmax
hmax=h（i,j）;
s=i-1;
t=j-1;
end
end
end
z=ones（m,n）;
fori=1:
m
forj=1:
n
ifa0（i,j）<=s&a2（i,j）<=t
%ifdouble（a（i,j））+double（a2（i,j））+a3（i,j）<=s+t+q
z（i,j）=0;
else
z（i,j）=255;
end
end
end
hmax
s
t
figure,imshow（z）;
toc
升级会员 