I = imread('C:\1.bmp');

figure,imshow(I)

[x,y]=size(I);                % 求出图象大小

b=double(I);

z0=max(max(max(I)));                   % 求出图象中最大的灰度

z1=min(min(min(I)));                   % 最小的灰度

T=(z0+z1)/2;                      % T赋初值，为最大值和最小值的平均值

TT=0;                             %TT赋初值

S0=0.0; n0=0.0;                       %为计算灰度大于阈值的元素的灰度总值、个数赋值

S1=0.0; n1=0.0;                       %为计算灰度小于阈值的元素的灰度总值、个数赋值

allow=50;                       % 新旧阈值的允许接近程度

d=abs(T-TT);

count=0;                         % 记录几次循环

while(d>=allow)                 % 迭代最佳阈值分割算法

count=count+1;

for i=1:x,

for j=1:y,

if I(i,j)>=T,

S0=S0+b(i,j);

n0=n0+1;

end

if I(i,j)

S1=S1+b(i,j);

n1=n1+1;

end

end

end

T0=S0/n0;

T1=S1/n1;

TT=(T0+T1)/2;            %在阈值T下，迭代阈值的计算过程

d=abs(T-TT);

T=TT;

end

S0

i1=im2bw(I,T/255);          % 图像在最佳阈值下二值化

figure,imshow(i1)

i2=bwperim(i1);         % 图像细化过程

figure,imshow(i2)

J=1-i2

% subplot(1,2,2);imshow(J);title('反相');

for m=1:50

J1=zeros(285,339);

for i=1:285

for j=1:339

k=0;

s=0;

if J(i,j)==1&J(i-1,j)*J(i+1,j)*J(i,j-1)*J(i,j+1)*J(i-1,j-1)*J(i-1,j+1)*J(i+1,j+1)*J(i+1,j-1)==0;

%统计非零相邻像素个数*************

k=J(i-1,j)+J(i+1,j)+J(i,j-1)+J(i,j+1)+J(i-1,j-1)+J(i-1,j+1)+J(i+1,j+1)+ J(i+1,j-1);

%统计相邻像素由0到1的次数****************

if J(i-1,j)==0 &J(i-1,j+1)==1

s=s+1;

end

if J(i-1,j+1)==0 &J(i,j+1)==1

s=s+1;

end

if J(i,j+1)==0 &J(i+1,j+1)==1

s=s+1;

end

if J(i+1,j+1)==0 &J(i+1,j)==1

s=s+1;

end

if J(i+1,j)==0 &J(i+1,j-1)==1

s=s+1;

end

if J(i+1,j-1)==0 &J(i,j-1)==1

s=s+1;

end

if J(i,j-1)==0 &J(i-1,j-1)==1

s=s+1;

end

if J(i-1,j-1)==0 &J(i-1,j)==1

s=s+1;

end

if (J(i,j+1)*J(i+1,j)==0|(J(i-1,j)==0&J(i,j-1)==0))&k<=6&k>=2&s==1

J1(i,j)=1;

end

end

end

end

J=J-J1;

J2=zeros(285,339);

for i=1:285

for j=1:339

k=0;

s=0;

if J(i,j)==1&J(i-1,j)*J(i+1,j)*J(i,j-1)*J(i,j+1)*J(i-1,j-1)*J(i-1,j+1)*J(i+1,j+1)*J(i+1,j-1)==0;

%统计非零相邻像素个数*************

k=J(i-1,j)+J(i+1,j)+J(i,j-1)+J(i,j+1)+J(i-1,j-1)+J(i-1,j+1)+J(i+1,j+1)+ J(i+1,j-1);

%统计相邻像素由0到1的次数****************

if J(i-1,j)==0 &J(i-1,j+1)==1

s=s+1;

end

if J(i-1,j+1)==0 &J(i,j+1)==1

s=s+1;

end

if J(i,j+1)==0 &J(i+1,j+1)==1

s=s+1;

end

if J(i+1,j+1)==0 &J(i+1,j)==1

s=s+1;

end

if J(i+1,j)==0 &J(i+1,j-1)==1

s=s+1;

end

if J(i+1,j-1)==0 &J(i,j-1)==1

s=s+1;

end

if J(i,j-1)==0 &J(i-1,j-1)==1

s=s+1;

end

if J(i-1,j-1)==0 &J(i-1,j)==1

s=s+1;

end

if (J(i-1,j)*J(i,j-1)==0|(J(i,j+1)==0&J(i+1,j)==0))&k<=6&k>=2&s==1

J2(i,j)=1;

end

end

end

end

J3=J-J2;

if J1==0&J2==0;

break

else

J=J3;

end

end

figure;imshow(J);title('细化后');

K=1-J;

figure;imshow(K);title('最终图')