%景曉軍文章中的二維直方圖
clear;
clc;
tic
a=imread（‘d:\image\xx.jpg‘‘jpg‘）;???????%Read?a?image?
%a=imread（‘d:\image\clock.tif‘‘tif‘）;
%a=imread（‘d:\image\zhuangjia.tiff‘‘tiff‘）;
%a=a（::1）;
%?%?a=M;
%?%load（‘source_N0.02.mat‘）;
%?%a=X;
%a=imread（‘syn1-g2.gif‘）;?
%load（‘裝甲車N_0.01.mat‘）;
%a=b;
%a=noise_h;
subplot（121）;
imshow（a）;??
[mn]=size（a）;
?
%b=imnoise（a‘salt?&?pepper‘0.003）;
%b=imnoise（b‘gaussian‘00.0015）;
%b?=?IMNOISE（a‘speckle‘0.09）;
%b=a;
a0=double（a）;
h=1;??????????????????????????
a1=zeros（mn）;
%?計算平均領域灰度的一維灰度直方圖
for?i=1:m
????for?j=1:n
????????for?k=-h:h
????????????for?w=-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（ij）=a0（pq）+a1（ij）;
?????????????end
????????end
????????a2（ij）=uint8（1/9*a1（ij））;
????end
end
fxy=zeros（256256）;
%?計算二維直方圖
for?i=1:m
????for?j=1:n
????????c=a0（ij）;
????????d=double（a2（ij））;
????????fxy（c+1d+1）=fxy（c+1d+1）+1;
?????end
??end
%??figure
%??mesh（fxy）;
%??title（‘二維灰度直方圖‘）;
Pxy=fxy/m/n;
P0=zeros（256256）;
Ui=zeros（256256）;
Uj=zeros（256256）;
P0（11）=Pxy（11）;
for?i=2:256
????P0（1i）=P0（1i-1）+Pxy（1i）;
end
for?i=2:256
????P0（i1）=P0（i-11）+Pxy（i1）;
end
for?i=2:256
????for?j=2:256
????????P0（ij）=P0（ij-1）+P0（i-1j）-P0（i-1j-1）+Pxy（ij）;
????end
end
P1=ones（256256）-P0;
Ui（11）=0;
for?i=2:256
????Ui（1i）=Ui（1i-1）+（1-1）*Pxy（1i）;
end
for?i=2:256
????Ui（i1）=Ui（i-11）+（i-1）*Pxy（i1）;
end
for?i=2:256
????for?j=2:256
????????Ui（ij）=Ui（ij-1）+Ui（i-1j）-Ui（i-1j-1）+（i-1）*Pxy（ij）;
????end
end
Uj（11）=0;
for?i=2:256
????Uj（1i）=Uj（1i-1）+（i-1）*Pxy（1i）;
end
for?i=2:256
????Uj（i1）=Uj（i-11）+（1-1）*Pxy（i1）;
end
for?i=2:256
????for?j=2:256
????????Uj（ij）=Uj（ij-1）+Uj（i-1j）-Uj（i-1j-1）+（j-1）*Pxy（ij）;
????end
end
uti=0;
utj=0;
for?i=1:256
????for?j=1:256
????????uti=uti+（i-1）*Pxy（ij）;?
????????utj=utj+（j-1）*Pxy（ij）;
????end
end
%?for?i=1:256
%?????for?j=1:256
%???????Ui1（ij）=Ui（256256）-Ui（ij）;
%???????Uj1（ij）=Uj（256256）-Uj（ij）;
%?????end
%?end
%計算類間方差
hmax=0;
for?i=1:256
????for?j=1:256
???????if?P0（ij）~=0&P1（ij）~=0
????????h（ij）=（（uti*P0（ij）-Ui（ij））^2+（utj*P0（ij）-Uj（ij））^2）/（P1（ij）*P0（ij））;
???????else
????????h（ij）=0;
???????end
????end
end
hmax=max（h（:））;
for?i=1:256
????for?j=1:256
????????if?h（ij）==hmax
?????????????s=i-1;
?????????????t=j-1;
?????????????continue;
?????????end
?????end
?end
hmax
s
t
z=ones（mn）;?
%c=zeros（mn）;??
for?i=1:m
????for?j=1:n
??????if?a（ij）<=s&a2（ij）<=t?
??????%if?double（b（ij））+double（a2（ij））<=s+t??
???????z（ij）=0;??????
??????end
??end
end
%imwrite（z‘D:\MATLAB7\work\feng\文章涉及到的\曲線otsu\曲線otsu實驗結果\baboon\wuprocess.jpg‘）;?
%imwrite（z‘D:\MATLAB7\work\feng\文章涉及到的\曲線otsu\曲線otsu實驗結果\baboon\you_process.jpg‘）;
%imwrite（uint8（b）

?屬性????????????大小?????日期????時間???名稱
-----------?---------??----------?-----??----

?????文件???????3208??2007-05-22?16:43??otsu_2d.m

-----------?---------??----------?-----??----

?????????????????3426????????????????????2