Here's something to get you started. Still needs some work to clean up some misidentified pixels and noise, but it's a very good start for you. It background corrects the image with CLAHE (imadapthisteq function), then converts to LAB color space and uses kmeans to classify the image into 3 colors.
clc; % Clear the command window.
close all; % Close all figures (except those of imtool.)
clear; % Erase all existing variables. Or clearvars if you want.
workspace; % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 20;
% Check that user has the Image Processing Toolbox installed.
hasIPT = license('test', 'image_toolbox'); % license('test','Statistics_toolbox'), license('test','Signal_toolbox')
if ~hasIPT
% User does not have the toolbox installed.
message = sprintf('Sorry, but you do not seem to have the Image Processing Toolbox.\nDo you want to try to continue anyway?');
reply = questdlg(message, 'Toolbox missing', 'Yes', 'No', 'Yes');
if strcmpi(reply, 'No')
% User said No, so exit.
return;
end
end
%===============================================================================
% Read in a standard MATLAB color demo image.
folder = pwd; % Determine where demo folder is (works with all versions).
baseFileName = 'a.jpg';
% Get the full filename, with path prepended.
fullFileName = fullfile(folder, baseFileName);
if ~exist(fullFileName, 'file')
% Didn't find it there. Check the search path for it.
fullFileName = baseFileName; % No path this time.
if ~exist(fullFileName, 'file')
% Still didn't find it. Alert user.
errorMessage = sprintf('Error: %s does not exist.', fullFileName);
uiwait(warndlg(errorMessage));
return;
end
end
rgbImage = imread(fullFileName);
% Get the dimensions of the image. numberOfColorBands should be = 3.
[rows, columns, numberOfColorChannels] = size(rgbImage);
% Display the original color image.
subplot(3, 3, 1);
imshow(rgbImage);
axis on;
title('Original Color Image', 'FontSize', fontSize, 'Interpreter', 'None');
% Enlarge figure to full screen.
set(gcf, 'Units', 'Normalized', 'Outerposition', [0, 0, 1, 1]);
drawnow;
% Extract the individual red, green, and blue color channels.
redChannel = rgbImage(:, :, 1);
greenChannel = rgbImage(:, :, 2);
blueChannel = rgbImage(:, :, 3);
redChannel = adapthisteq(redChannel, 'clipLimit',0.02,'Distribution','rayleigh');
greenChannel = adapthisteq(greenChannel, 'clipLimit',0.02,'Distribution','rayleigh');
blueChannel = adapthisteq(blueChannel, 'clipLimit',0.02,'Distribution','rayleigh');
% Recombine separate color channels into a single, true color RGB image.
rgbImage = cat(3, redChannel, greenChannel, blueChannel);
subplot(3, 3, 2);
imshow(rgbImage);
axis on;
title('Background corrected Color Image', 'FontSize', fontSize, 'Interpreter', 'None');
% Convert to grayscale.
grayImage = rgb2gray(rgbImage);
% It's still a bit noisy from the CLAHE operation so do a median filter
grayImage = medfilt2(grayImage, [15, 15]);
subplot(3, 3, 3);
imshow(grayImage);
axis on;
title('Gray Scale Image', 'FontSize', fontSize, 'Interpreter', 'None');
% Take the histogram
subplot(3, 3, 4);
histogram(grayImage);
xlabel('Gray Level', 'FontSize', fontSize);
ylabel('Count', 'FontSize', fontSize);
grid on;
ax = gca
ax.XTick = [0:16:255];
ax.GridAlpha = 0.75;
drawnow;
% Binarize the image
binaryImage = grayImage < 112;
% Just keep the larger blobs larger than 1400
binaryImage = bwareafilt(binaryImage, [1400, inf]);
% Fill Holes
binaryImage = imfill(binaryImage, 'holes');
% Get rid of any blobs touching the border of the image.
binaryImage = imclearborder(binaryImage);
subplot(3, 3, 5);
imshow(binaryImage);
axis on;
title('Binary Image', 'FontSize', fontSize, 'Interpreter', 'None');
% Label
[labeledImage, numBlobs] = bwlabel(binaryImage);
% Let's assign each blob a different color to visually show the user the distinct blobs.
coloredLabels = label2rgb (labeledImage, 'hsv', 'k', 'shuffle'); % pseudo random color labels
% coloredLabels is an RGB image. We could have applied a colormap instead (but only with R2014b and later)
subplot(3, 3, 6);
imshow(coloredLabels);
axis image; % Make sure image is not artificially stretched because of screen's aspect ratio.
caption = sprintf('Pseudo colored labels, from label2rgb().\nBlobs are numbered from top to bottom, then from left to right.');
title(caption, 'FontSize', fontSize);
% Get all the blob properties. Can only pass in originalImage in version R2008a and later.
blobMeasurements = regionprops(labeledImage, 'Area');
allAreas = sort([blobMeasurements.Area], 'descend')
% Dilate some and mask out the original image.
binaryImage = imdilate(binaryImage, true(9));
% Find the mean gray level of the background.
backgroundMask = ~binaryImage;
meanGL = mean(grayImage(backgroundMask))
% Make color channels the mean where the background is so we can use kmeans() or rgb2ind().
redChannel(backgroundMask) = meanGL;
greenChannel(backgroundMask) = meanGL;
blueChannel(backgroundMask) = meanGL;
% Recombine separate color channels into a single, true color RGB image.
rgbImage = cat(3, redChannel, greenChannel, blueChannel);
subplot(3, 3, 7);
imshow(rgbImage);
title('Masked RGB Image', 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;
% Now do kmeans color classification
% http://www.mathworks.com/help/images/examples/color-based-segmentation-using-k-means-clustering.html?prodcode=IP&language=en&requestedDomain=www.mathworks.com
labImage = rgb2lab(rgbImage);
[rows, columns, numberOfColorChannels] = size(labImage);
labList = reshape(labImage, (rows * columns), 3);
nColors = 3;
% repeat the clustering 3 times to avoid local minima
[cluster_idx, cluster_center] = kmeans(labList, nColors,...
'distance','sqEuclidean', 'Replicates',3);
% Show the results.
subplot(3, 3, 8:9);
classifiedImage = reshape(cluster_idx, rows, columns);
imshow(classifiedImage, []);
title('Image labeled by cluster index');
