kmeansclusters.m

利用matlab实现彩色图像的分割。算法主要是利用聚类算法。

function [clusters, mus, gs, mapping] = kmeansClusters(image, numClusters)

% KMEANS k-Means clustering - used to initialize clusters
%
%  [CLUSTERS,MUS,ITERATIONS] = KMEANS(DATA,NUMCLUSTERS) 
%
%  Input:
%       IMAGE - MxNx3 image, where each pixel is a point in a 
%               perceptually unif color space (e.g. LAB)
%       NUMCLUSTERS - the number of clusters to use
%
%  Output: 
%       CLUSTERS -  MxNxC matrix, where C is the number of clusters (not necessarily=numClusters)
%                   If a pixel belongs to a cluster c, clusters(y,x,c)=1
%       MUS - Cx3 matrix containing the final set of mean values for each cluster 
%       GS - MxN matrix, where each element is the group number that this pixel
%            belongs to
%       MAPPING - Cx1 vector that tracks the group number of each cluster
%
% Jeff Walters & Angi Chau
% Feb 2003

% define stopping condition
EPSILON = 0.005;

% reshape the image into a vector for k-means clustering (scan down columns 
% like in matlab)
height = size(image,1);
width = size(image,2);
imageVect = [];
for x=1:width
    for y=1:height
        imageVect(end+1,:) = [x y image(y,x,1) image(y,x,2) image(y,x,3)];
    end
end
disp('image reshaped into vector for clustering');       
%imageVect = reshape(image, [], 3);

% uniformly pick out pixels to use as the initial mean points
%initIndices = 1:floor(height*width/numClusters):height*width;
%initIndices = initIndices(1:numClusters);

% to make sure we don't get pixels too close together, let's
% pick from the sorted vector. so unless, numClusters close to MxN
% we will hopefully have avoided this problem.
%sortedVect = sort(imageVect);
%currentmus = sortedVect(initIndices,:);

% pick out uniformly distributed random points in LAB/XY space
Lrand = rand(numClusters,1)*100;
Arand = rand(numClusters,1)*200-100;
Brand = rand(numClusters,1)*200-100;
xrand = rand(numClusters,1)*width;
yrand = rand(numClusters,1)*height;
currentmus = [xrand yrand Lrand Arand Brand];

% initialize variable to hold the old means
oldmus = zeros(numClusters, 3);

notdone = 1;
iterations = 0;

while notdone
  disp(sprintf('Iteration %d',iterations));
  
    % classify the points based on the current means
    groups = kclassify(imageVect, currentmus);
    
    % check out the groupings
    %colormap(hsv(numClusters));
    %imagesc(reshape(groups, height, width));
    
    % remember old means before resetting new means
    oldmus = currentmus;
    
    % calculate new means
    for j=1:numClusters
        ind = find(groups == j);   % get everyone in this cluster
        
        if (isempty(ind) ~= 1)
            currentmus(j,:) = mean(imageVect(ind, :),1);
        else
            % if no one's in this group, it's a dead class. let's just
            % leave the value alone and see if anyone gets assigned to
            % it on the next iteration.
            %disp('dead class alert!!');
            currentmus(j,:) = oldmus(j,:);
%             currentmus(j,:) = [ rand*width
%                                 rand*height
%                                 rand*100
%                                 rand*200-100
%                                 rand*200-100]';
        end    
    end
    
    % if the difference is small enough, we're done!
    delta = norm(currentmus - oldmus) / numClusters;
    notdone = (EPSILON < delta);
    disp(sprintf('Delta = %2.4f.\n', delta));

    iterations = iterations + 1;
    %pause
end

% shape groupings back into image size
groupsShaped = reshape(groups, height, width);
gs=groupsShaped;

% create masks for all non-empty groups
clusters = zeros(height,width,0);
mus = zeros(0,size(currentmus,2));
mapping = zeros(0,1);
for j=1:numClusters
    mask=(groupsShaped==j);
    if (sum(mask(1:end)) ~= 0)
        clusters(:,:,end+1) = mask;
        mus(end+1,:) = currentmus(j,:);  
        mapping(end+1) = j;
    end
end

% debuggy things
figure
colormap(bone(numClusters));
imagesc(groupsShaped);
colorbar;
title('kmeans clusters')