[starlove1006]

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%聚类算法%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%不清楚malab查找表的函数
%%%%%%%%%%%%改用两个存储矩阵,IsSame存储同类像素点坐标,NotSame存储不同类的像素点坐标
%%%%%%%%%%%但不知道为什么会出现死循环%%%%%%%%%%%%%%%%%%%

clc
a2=imread('D:/1/frame46.bmp');
figure;imshow(a2);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
p=im2double(a2);
r=p(:,:,1);
g=p(:,:,2);
b=p(:,:,3);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%RGB转HSI
num=0.5*((r-g)+(r-b));
den=sqrt((r-g).^2+(r-b).*(g-b));
theta=acos(num./(den+eps));
H=theta;
H(b>g)=2*pi-H(b>g);
num=min(min(r,g),b);
den=r+g+b;
den(den==0)=eps;
S=1-3.*num./den;
H(S==0)=0;
I=(r+g+b)/3;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%求H值的直方图
% H=H/(2*pi);
H=H/(2*pi);
H=round(H*360);
%hist(H,360);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[WW,HH]=size(P); %取图像大小
[x,y]=find(H==90); % 暂时先确定为90 取峰值对应h值（即第一个聚类中心）的任意一个坐标（x1,y1）
x1=x(1);
y1=y(1);
seed=90; %将生长起始点的h值存入seed seed=firstpeak;
threshold=5; %阈值 待求

MarkP=zeros(WW,HH); % 作一个与原图像等大的图像矩阵MarkP，作为输出图像矩阵，同一类的赋予相同的标号
lab=1; % lab是标号
MarkP(x1,y1)=lab; % 对MarkP中与所取点的相对位置的点标号

count=0; % 记录 每判断一点周围八点符合条件的新点的数目
samp=zeros(100,1);%存放最终显示的H值，（聚类中心的H值）
a=1;
IsSame=zeros(WW*HH,2); %IsSame存放聚类中属同类的点的坐标，IsSame(a,1)放横坐标，IsSame(a,2)放纵坐

标
b=1;
NotSame=zeros(WW*HH,2); %NotSame存放聚类中属同类的点的坐标，

while b>0
b=0;
while a>=0
for u=-1:1
for v=-1:1 %u、v为偏移量
if (x1+u-1)>0 & (x1+u+1)<(WW+1) & (y1+v-1)>0 & (y1+v+1)<(HH+1) %判断此点是否为图像

边界上的点
d= abs(H(x1+u,y1+v)-seed);
if MarkP(x1+u,y1+v)==0 && d<=threshold %判断是否为未标记，且符合阈值条件的点
a=a+1; %统计同类的点的个数
IsSame(a,1)=x1+u; %存放点的坐标，为区域生长做准备
IsSame(a,2)=y1+v;
count=count+1;
MarkP(x1+u,y1+v)=lab; %标记为一类
seed=(seed+H(x1+u,y1+v))/2; %取平均值为新的seed
else if MarkP(x1+u,y1+v)==0 %点未标记，但不属于该类
b=b+1;
NotSame(b,1)=x1+u;
NotSame(b,2)=y1+v;
end
end
end
end
end
if count>0 & a>0 %区域生长
x1= IsSame(a,1);
y1= IsSame(a,2);
IsSame(a,1)=0; %取出点后清空
IsSame(a,2)=0;
a=a-1;
end
end %while 1次聚类结束
if b>0 %2次聚类 聚集可能被1次聚类忽略的点
for n=1:b %遍历NotSame中的点，属同类的在MakP中标记（但没进行再一次的生长）并在NotSame清空

，否则取其中一点作为另一次的聚类中心 前面第一个while后的 'b=0'相当于把NotSame清空
if abs(H(NotSame(n,1),NotSame(n,2))-seed)<=threshold
MarkP(NotSame(n,1),NotSame(n,2))=lab; %标记为一类
seed=(seed+H(NotSame(n,1),NotSame(n,2)))/2 %取平均值为新的seed
NotSame(n,1)=0; %取出点后清空
NotSame(n,2)=0;
else x1=NotSame(n,1); %为另一次聚类准备
y1=NotSame(n,2);
end
end
end
samp(lab)=seed;
%为另一次聚类准备
lab=lab+1;
seed=H(x1,y1); %x1,y1的取值在上一个if结构中 这里直接这样似有点不妥
end %while 全部聚类结束