对单通道灰度图(32位)的自定义初始聚类中心灰度值的Kmeans聚类分割
- 最近学习需要使用Kmeans聚类分割一幅灰度图像,但是发现opencv自带的kmeans随机分配初始中心的方式KMEANS_PP_CENTERS并不好用,而且封装的函数看不到每次迭代的聚类中心变化过程,使用KMEANS_USE_INITIAL_LABELS自定义初始中心,发现相关的说明文档不多,且效果不好,必须手动先分配一波labels所有分类,很麻烦,效果也不好
所有就自己写了一个特别好用的基于opencv的Kmeans聚类分割灰度图像函数,可以用户自定义初始化聚类中的每各初始灰度值。
opencv中代码如下
void main()
{
Mat C=imread("灰度图(单通道)");
Mat Img = C.clone();
float c1[10] = { 0 }, c2[10] = { 0 }, c3[10] = { 0 }, c4[10] = { 0 }; //定义迭代十次,c存放每次迭代的聚类中心值
c1[1] = 5; //初始化第一类聚类中心灰度值
c2[1] = 50;
c3[1] = 100;
c4[1] = 200;
int kk = 4, times = 0; //times记录迭代次数,,kk代表分4类
Mat cla = C.clone(); //结果图
Mat result = C.clone();
int num1 = 0, num2 = 0, num3 = 0, num4 = 0; //存放每类的像素个数和元素总和
int sum1 = 0, sum2 = 0, sum3 = 0, sum4 = 0;
int d1 = 2, d2 = 2, d3 = 2, d4 = 2; //先任意给定可以使循环开始的大于1的聚类中心间距值
int dis1 = 0, dis2 = 0, dis3 = 0, dis4 = 0;
while ((d1 > 0.01 || d2 > 0.01 || d3 > 0.01 || d4 > 0.01) && (times < 10)) //迭代核心语句
{
times = times + 1;
for (int i = 0; i < Img.rows; i++)
{
uchar *data = Img.ptr<uchar>(i);
uchar *datacla = cla.ptr<uchar>(i);
for (int j = 0; j < Img.cols; j++)
{
int value = data[j];
int value2 = datacla[j];
dis1 = (value - c1[times])*(value - c1[times]); //欧式距离
dis2 = (value - c2[times])*(value - c2[times]);
dis3 = (value - c3[times])*(value - c3[times]);
dis4 = (value - c4[times])*(value - c4[times]);
int mindis = min_dis(dis1, dis2, dis3, dis4);
if (mindis == dis1) //分成4类直接分别赋值0,60,150,240灰度值
{
datacla[j] = 0;
num1 = num1 + 1;
sum1 = sum1 + data[j];
}
else if (mindis == dis2)
{
datacla[j] = 60;
num2 = num2 + 1;
sum2 = sum2 + data[j];
}
else if (mindis == dis3)
{
datacla[j] = 150;
num3 = num3 + 1;
sum3 = sum3 + data[j];
}
else if (mindis == dis4)
{
datacla[j] = 240;
num4 = num4 + 1;
sum4 = sum4 + data[j];
}
}
}
//调整新的聚类中心
c1[times + 1] = sum1 / num1;
c2[times + 1] = sum2 / num2;
c3[times + 1] = sum3 / num3;
c4[times + 1] = sum4 / num4;
//两次之间的聚类中心变化大小 以此做为终止迭代的标准
d1 = abs(c1[times] - c1[times + 1]);
d2 = abs(c2[times] - c2[times + 1]);
d3 = abs(c3[times] - c3[times + 1]);
d4 = abs(c4[times] - c4[times + 1]);
//在下一次迭代前将每类的加和置0
sum1 = 0;
sum2 = 0;
sum3 = 0;
sum4 = 0;
num1 = 0;
num2 = 0;
num3 = 0;
num4 = 0;
}
namedWindow("kk", 0);
imshow("kk", cla);
waitKey(0);
}
附上Matlab代码
%% 对功率灰度图Img进行K-means聚类 分四类
Img=imread('D:/image/C.bmp');
Img=Img(1:end,1:end);
imshow(Img);
Img=double(Img);
%% 各个参数
c1=zeros(1,15);c2=zeros(1,15);c3=zeros(1,15);c4=zeros(1,15); %c存放每次迭代的聚类中心值
% 初始聚类中心
c1(1)=5; %类1
c2(1)=50; %类2
c3(1)=100; %类3
c4(1)=200; %类4
%分为四类
k=4;
[m,n]=size(Img);
%记录迭代次数
times=0;
%存储分类标识 取值为1 2 3 4
class=zeros(m,n);
%存储聚类结果
result=zeros(m,n);
%存放每类的像素个数和元素总和
num1=0;num2=0;num3=0;num4=0;
sum1=0;sum2=0;sum3=0;sum4=0;
%先任意给定可以使循环开始的大于1的聚类中心间距值
d1=2;d2=2;d3=2;d4=2;
%% 迭代部分
tic
while((d1>0.01||d2>0.01||d3>0.01||d4>0.01)×<15) %限制迭代次数
times=times+1
%%%%%%%%%%%%%%%%%归类%%%%%%%%%%%%%%%%%%%%%%%%%
for i=1:m
for j=1:n
dis1=(Img(i,j)-c1(times))^2; %欧式距离
dis2=(Img(i,j)-c2(times))^2;
dis3=(Img(i,j)-c3(times))^2;
dis4=(Img(i,j)-c4(times))^2;
min_dis=min([dis1,dis2,dis3,dis4]);
switch(min_dis)
case dis1
class(i,j)=1;
num1=num1+1;
sum1=sum1+Img(i,j);
case dis2
class(i,j)=2;
num2=num2+1;
sum2=sum2+Img(i,j);
case dis3
class(i,j)=3;
num3=num3+1;
sum3=sum3+Img(i,j);
case dis4
class(i,j)=4;
num4=num4+1;
sum4=sum4+Img(i,j);
end
end
end
e_new=sqrt(min_dis);
%调整新的聚类中心
c1(times+1)=sum1/num1;
c2(times+1)=sum2/num2;
c3(times+1)=sum3/num3;
c4(times+1)=sum4/num4;
%两次之间的聚类中心变化大小 以此做为终止迭代的标准
d1=abs(c1(times)-c1(times+1));
d2=abs(c1(times)-c1(times+1));
d3=abs(c1(times)-c1(times+1));
d4=abs(c1(times)-c1(times+1));
%在下一次迭代前将每类的加和置0
sum1=0;
sum2=0;
sum3=0;
sum4=0;
num1=0;
num2=0;
num3=0;
num4=0;
end
toc
%%%%%%%%%%%%将聚类后的功率灰度图赋值%%%%%%%%%%%%%%%%%%%
%利用矩阵
result(class==1)=0;
result(class==2)=60;
result(class==3)=150;
result(class==4)=240;
%显示图像
result=uint8(result); %谨记图像的两种显示格式 uint8和[0,1]
%写入图像
imwrite(mat2gray(result),'D:/image/C-4.bmp');
imshow(result);title('聚类后的结果');
print('Done');