Hello Mat

 找回密码
 立即注册
查看: 8453|回复: 4

SOA下的kmeans图像分割

[复制链接]

1323

主题

1551

帖子

0

金钱

管理员

Rank: 9Rank: 9Rank: 9

积分
22647
发表于 2019-10-21 23:58:27 | 显示全部楼层 |阅读模式
SOA下的kmeans图像分割
链接:https://pan.baidu.com/s/1hjGrKwBi56hhqT6l91w7ow 提取码:1tah

具体链接在halcom.cn论坛,联系人QQ:3283892722
该论坛是一个学习交流平台,我会逐一的和大家分享学习。
欢迎大家录制视频,你可在论坛进行打赏分享。
视频专用播放器:http://halcom.cn/forum.php?mod=viewthread&tid=258&extra=page%3D1\

数据获取(get_cluster_data函数):链接:https://pan.baidu.com/s/1r2AngQUwp8AP_8Qh89py_g 提取码:ejpf
聚类中心center_nearest:
  1. function [minDist,index] = center_nearest(center,ClusterData)
  2. % ClusterData:待聚类的数据
  3. % center:聚类的中心
  4. % 函数目的:求解每一个数据隶属于的中心
  5. minDist = zeros(size(ClusterData,1),1);
  6. index = zeros(size(ClusterData,1),1);
  7. dist = zeros(size(center,1),1);
  8. for i=1:length(ClusterData)
  9.     for j=1:size(center,1)
  10.         dist(j) = norm(ClusterData(i,:) - center(j,:), 2);
  11.     end
  12.     [ minDist(i), index(i) ]=  min(dist);
  13. end
复制代码
Kmeans_fun函数(适应度函数)如下:
游客,如果您要查看本帖隐藏内容请回复

kmeans分割函数如下:ClusterData是一个75088x2的矩阵;也就是图像数据;
  1. % K均值聚类算法
  2. clc,clear,close all;
  3. warning off;
  4. I_rgb=imread('2.jpg');
  5. I_rgb=im2double(I_rgb);   %双精度
  6. if size(I_rgb,3)==1
  7.     msgbox('请输入一张彩色图像!!!');
  8.     return;
  9. end
  10. R=I_rgb(:,:,1);
  11. G=I_rgb(:,:,2);
  12. B=I_rgb(:,:,3);
  13. [M,N]=size(R);   % 行、列
  14. ClusterData = get_cluster_data(R,G,B);
  15. % 分类数
  16. [len, depth]= size(ClusterData);
  17. Kc = 3;
  18. % 初始化中心
  19. index = ceil( unifrnd(1,len, 1, Kc) );   % 均匀分布
  20. center = ClusterData(index,:);
  21. center_new = center;   % 初始化新的中心
  22. % maim loop
  23. maxgen = 1e3;  % 最大迭代次数
  24. for i=1:maxgen
  25.     % 计算 每一个待分类的数据隶属于的中心编号
  26.     [center_minDist,center_index] = center_nearest(center,ClusterData);
  27.     % 更新center
  28.     for j=1:Kc
  29.         index = find(center_index==j);
  30.         center_new(j, :) = mean( ClusterData(index,:) );
  31.     end
  32.     error = norm( center_new-center, 2 );  % 相连两次中心点的距离值
  33.     center = center_new;   % 更新中心
  34.     if error<0.000001
  35.         fprintf('\n');
  36.         disp(['收敛迭代次数:', num2str(i)])
  37.         disp( '聚类中心为:' )
  38.         disp( num2str(center) )
  39.         break;
  40.     end
  41. end
  42. %% 分割显示
  43. % a = reshape(ClusterData(:,1),M,N);  % 第1列
  44. % b = reshape(ClusterData(:,2),M,N);  % 第2列
  45. A = zeros(M*N, 1);
  46. A(find(center_index==1), :) = 1;
  47. A1 = reshape(A,M,N);  % 第1类

  48. A = zeros(M*N, 1);
  49. A(find(center_index==2), :) = 1;
  50. A2 = reshape(A,M,N);  % 第2类

  51. A = zeros(M*N, 1);
  52. A(find(center_index==3), :) = 1;
  53. A3 = reshape(A,M,N);  % 第3类

  54. A4 = reshape(center_index,M,N);  % 第2列

  55. figure(1),
  56. subplot(221),imshow(A1,[]);xlabel('第1类')
  57. subplot(222),imshow(A2,[]);xlabel('第2类')
  58. subplot(223),imshow(A3,[]);xlabel('第3类')
  59. subplot(224),imshow(A4,[]);xlabel('整体显示')

  60. % 重构显示彩色图
  61. R1 = immultiply(R, logical(A1));
  62. G1 = immultiply(G, logical(A1));
  63. B1 = immultiply(B, logical(A1));
  64. RGB1 = cat(3,R1,G1,B1);

  65. R2 = immultiply(R, logical(A2));
  66. G2 = immultiply(G, logical(A2));
  67. B2 = immultiply(B, logical(A2));
  68. RGB2 = cat(3,R2,G2,B2);

  69. R3 = immultiply(R, logical(A3));
  70. G3 = immultiply(G, logical(A3));
  71. B3 = immultiply(B, logical(A3));
  72. RGB3 = cat(3,R3,G3,B3);
  73. figure(2),
  74. subplot(221),imshow(I_rgb,[]);xlabel('原始图像')
  75. subplot(222),imshow(RGB1,[]);xlabel('第1类')
  76. subplot(223),imshow(RGB2,[]);xlabel('第2类')
  77. subplot(224),imshow(RGB3,[]);xlabel('第3类')
复制代码
SOA优化下的Kmeans图像分割:
ClusterData是一个75088x2的矩阵;也就是图像数据
  1. clc,clear,close all;
  2. warning off;
  3. format longG;
  4. %% 获取图像数据
  5. I_rgb=imread('2.jpg');
  6. I_rgb=im2double(I_rgb);   %双精度
  7. if size(I_rgb,3)==1
  8.     msgbox('请输入一张彩色图像!!!');
  9.     return;
  10. end
  11. R=I_rgb(:,:,1);
  12. G=I_rgb(:,:,2);
  13. B=I_rgb(:,:,3);
  14. [M,N]=size(R);   % 行、列
  15. ClusterData = get_cluster_data(R,G,B);
  16. % 分类数
  17. [len, depth]= size(ClusterData);
  18. Kc = 3;
  19. %% SOA 参数
  20. maxiter = 20;  % 迭代次数
  21. sizepop = 20;  % 种群数量
  22. Umax=0.9500;   % 最大隶属度值
  23. Umin=0.0111;   % 最小隶属度值
  24. popmin = min(min(ClusterData));  popmax = max(max(ClusterData)); % x

  25. fai1 = 0.5;
  26. fai2 = 0.5;
  27. w1 = 0.5;
  28. % w2 = 0.5;
  29. w2max = 0.7;
  30. w2min = 0.2;
  31. nvar = 2*Kc;      % 变量的个数
  32. %% 初始化种群
  33. for i=1:sizepop
  34.     index = ceil( unifrnd(1,len, 1, Kc) );   % 均匀分布
  35.     center = ClusterData(index,:);
  36.     pop(i,:) = center(:);
  37.     [fitness(i), center_new]= Kmeans_fun( pop(i,:), ClusterData, Kc );
  38.     pop(i,:) = center_new(:);
  39. end
  40. %% 记录一组最优值
  41. [bestfitness,bestindex]=min(fitness); % 索引最小值
  42. zbest=pop(bestindex,:);   % 全局最佳
  43. gbest=pop;                % 个体最佳
  44. fitnessgbest=fitness;     % 个体最佳适应度值
  45. fitnesszbest=bestfitness; % 全局最佳适应度值
  46. Frontpop = pop;
  47. %% 迭代寻优
  48. for i=1:maxiter
  49.     disp(['当前迭代次数: ', num2str(i)])
  50.     for j=1:sizepop
  51.        dego = gbest(j,:)-pop(j,:);
  52.        dalt = zbest-pop(j,:);
  53.        dpro = Frontpop(j,:) - pop(j,:);
  54.        Frontpop = pop;  % 记录上一时刻的最优种群
  55.        dij = sign( w1*dpro+fai1*dego+fai2*dalt );   % 确定方向
  56.        % 确定步长
  57.        [maxfitness,index] = sort(fitness,'descend');
  58.        u = Umax - index(j)/sizepop *(Umax-Umin);
  59.        u = u.*rand(1,nvar);       % nvar个变量
  60.        T = sqrt(-log(u) );
  61.        xmin = pop(index(1),:);    % 最差的个体
  62.        xmax = pop(index(end),:);  % 最优的个体
  63.       
  64.        w2 = w2max - j/sizepop*(w2max-w2min);
  65.        delta = w2*abs( xmin-xmax );
  66.       
  67.        alpha = delta.*T;
  68.       
  69.        % 步长进行限制
  70.        for k=1:length(alpha)
  71.            if alpha(k)>0.6
  72.                alpha(k)=0.6;
  73.            elseif alpha(k)<-0.6
  74.                alpha(k) = -0.6;
  75.            end
  76.        end
  77.       
  78.        % 位置更新
  79.        pop(j,:) = pop(j,:) + dij.*alpha;        

  80.         % x 越界限制
  81.         for k =1:nvar
  82.             if pop(j,k)>popmax
  83.                 pop(j,k)=popmax;
  84.             end
  85.             if pop(j,k)<popmin
  86.                 pop(j,k)=popmin;
  87.             end
  88.         end

  89.         % 适应度更新
  90.         [fitness(j), newpop]= Kmeans_fun(pop(j,:), ClusterData, Kc );
  91.         pop(j,:) = newpop(:);
  92.         
  93.         % 比较  个体间比较
  94.         if fitness(j)<fitnessgbest(j)
  95.             fitnessgbest(j) = fitness(j);
  96.             gbest(j,:) = pop(j,:);
  97.         end
  98.         if fitness(j)<bestfitness
  99.             bestfitness = fitness(j);
  100.             zbest =  pop(j,:);
  101.         end
  102.         
  103.     end
  104.     fitness_iter(i) = bestfitness;
  105. end
  106. %% 结果显示
  107. fprintf('\n')
  108. disp(['最优适应度值: ', num2str(bestfitness)])
  109. fprintf('\n')
  110. zbest_center = reshape(zbest, Kc, 2 );
  111. disp( '聚类中心为:' )
  112. disp( num2str(zbest_center) );
  113.   
  114. figure('color',[1,1,1])
  115. plot(fitness_iter,'ro-','linewidth',2);
  116. xlabel('迭代次数'); ylabel('适应度值');
  117. axis tight;grid on;

  118. %% 分割显示
  119. % 计算 每一个待分类的数据隶属于的中心编号
  120. [center_minDist,center_index] = center_nearest(zbest_center,ClusterData);
  121. % a = reshape(ClusterData(:,1),M,N);  % 第1列
  122. % b = reshape(ClusterData(:,2),M,N);  % 第2列
  123. A = zeros(M*N, 1);
  124. A(find(center_index==1), :) = 1;
  125. A1 = reshape(A,M,N);  % 第1类

  126. A = zeros(M*N, 1);
  127. A(find(center_index==2), :) = 1;
  128. A2 = reshape(A,M,N);  % 第2类

  129. A = zeros(M*N, 1);
  130. A(find(center_index==3), :) = 1;
  131. A3 = reshape(A,M,N);  % 第3类

  132. A4 = reshape(center_index,M,N);  % 第2列

  133. figure(1),
  134. subplot(221),imshow(A1,[]);xlabel('第1类')
  135. subplot(222),imshow(A2,[]);xlabel('第2类')
  136. subplot(223),imshow(A3,[]);xlabel('第3类')
  137. subplot(224),imshow(A4,[]);xlabel('整体显示')

  138. % 重构显示彩色图
  139. R1 = immultiply(R, logical(A1));
  140. G1 = immultiply(G, logical(A1));
  141. B1 = immultiply(B, logical(A1));
  142. RGB1 = cat(3,R1,G1,B1);

  143. R2 = immultiply(R, logical(A2));
  144. G2 = immultiply(G, logical(A2));
  145. B2 = immultiply(B, logical(A2));
  146. RGB2 = cat(3,R2,G2,B2);

  147. R3 = immultiply(R, logical(A3));
  148. G3 = immultiply(G, logical(A3));
  149. B3 = immultiply(B, logical(A3));
  150. RGB3 = cat(3,R3,G3,B3);
  151. figure(2),
  152. subplot(221),imshow(I_rgb,[]);xlabel('原始图像')
  153. subplot(222),imshow(RGB1,[]);xlabel('第1类')
  154. subplot(223),imshow(RGB2,[]);xlabel('第2类')
  155. subplot(224),imshow(RGB3,[]);xlabel('第3类')
复制代码
输出:
  1. 收敛迭代次数:19
  2. 聚类中心为:
  3. 21.515      32.6575
  4. 3.61474       4.6504
  5. 15.0719      16.9821
复制代码

















本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有帐号?立即注册

x
算法QQ  3283892722
群智能算法链接http://halcom.cn/forum.php?mod=forumdisplay&fid=73
回复

使用道具 举报

0

主题

4

帖子

18

金钱

新手上路

Rank: 1

积分
22
发表于 2020-1-25 18:02:58 | 显示全部楼层
看一下大佬的代码
回复 支持 反对

使用道具 举报

0

主题

59

帖子

1

金钱

注册会员

Rank: 2

积分
72
发表于 2020-6-20 21:55:51 | 显示全部楼层
非常好的代码
回复 支持 反对

使用道具 举报

0

主题

54

帖子

1

金钱

注册会员

Rank: 2

积分
55
发表于 2021-1-10 11:40:44 | 显示全部楼层
111111111111
回复 支持 反对

使用道具 举报

0

主题

27

帖子

0

金钱

注册会员

Rank: 2

积分
54
发表于 2021-12-19 01:55:10 | 显示全部楼层
666666666666666666
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

Python|Opencv|MATLAB|Halcom.cn ( 蜀ICP备16027072号 )

GMT+8, 2024-11-26 07:45 , Processed in 0.211283 second(s), 22 queries .

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表