How to compare an image with another image?

I want to determine the degree of mixing of 2 different colored beads. Below image represents the ideal mixing case.

I want to use the above ideal image as a standard and use it to compare other images (like the one below) to determine how well the beads are mixed.

Please let me know how it can be done? Is there any factor/variable I can assign for the ideal case and compare it with the non-deal cases?

Thank You.

Answer 1

[Image Analyst]

First, put the beads on a black background rather than a white background. Next, use crossed polarizers to knock out the specular reflections. Then you can just convert to HSV color space and compute the area fraction of white, and of blue. See the Color Thresholder App on the Apps tab of the tool ribbon. Pretty easy - give it a try.

[Shashank]

Thank you Image Analyst. But isn't this code to identify the different colored beads and compute their area fraction.

I want to calculate the degree of mixing of the beads, that is, how well they are uniformly distributed.

For eg:

In the above image, area fraction of blue and white beads are same, that is, 50%. But their degree of mixing is zero, as they are in unmixed state.

If I want to calculate the degree of mixing, I thought we have to assign a variable (such as deg_mix) and equate it to 100 for ideal mixing case (image given in above thread) and to 0 for no mixing case. We could then compare randomly mixed image with ideal case and unmixed case and then assign a value to deg_mix such that 0<=deg_mix<=100. My question is, how should I compute the deg_mix variable?

Please let me know if this is the correct approach or is there any other better method to find the degree of mixing?

Thank you

[Image Analyst]

Well I'm not going to do the whole project for you. Don't you think that first identifying the blue and white pixels is the first step? And I did that.

It is not correct to say that in the above image where the two bead colors are separated into bins has a mixing of zero. To illustrate, what if you had 4 bins each with a single color? Would that still be zero? What if you have an 8-by-8 grid of bins, like a chessboard? Is that still zero? What if you have 600 bins so that each bead was in its own bin in a perfect alternating color chessboard? Is that still zero? What if the same bins and beads were there, but you shuffled the bins around so that the pattern was random? Is that still zero? No, of course not. So where exactly is the cutoff between mixed and not mixed? There is none - it has to be a continuum. And the only way you'll have no mixing is if there is a single color only in the image.

What you might want to look at are uniformity metrics. Like one way is to divide the image into tiles, like a chessboard. Then compute the mean area fraction in each tile. Then take the standard deviation of those. If the same area fraction is present in every tile (sub-image), then the standard deviation of the area fractions will be zero.

[Shashank]

Thank you sir. I have noted all the points mentioned above and have prepared a code which functions satisfactorily.

I would like you to have a look at it.

For the above image where the 2 bead colors are separated into bins, I get the standard deviation to be 0. What does it mean? Thanks.

[Image Analyst]

But you didn't note all the points mentioned below. You seemed to ignore my suggestions in that code. For example, I don't see anything in the code to deal with the white specular reflections in the blue beads. Also, you don't even display the masks so how can you tell if they're correct or not?

Plus the comments are nowhere near sufficient in terms of number of them or explaining what the code is doing.

I'd also advise to pick more descriptive variable names, like rgbImage instead of I, upperLeftQuadrant instead of A, etc. You don't want your code to look like an alphabet soup mess with no comments or descriptive variable names to help the reader understand what's going on.

[Shashank]

Sorry sir, I will make the necessary changes and follow proper etiquette from now.

Thank you.

Answer 2

[Image Analyst]

ALright, looks like you must be having trouble. So I made this code to get you started.

clc;    % Clear the command window.
close all;  % Close all figures (except those of imtool.)
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 16;

%===============================================================================
% Get the name of the image the user wants to use.
baseFileName = 'beads.jpg';
folder = pwd;
fullFileName = fullfile(folder, baseFileName);

% Check if file exists.
if ~exist(fullFileName, 'file')
	% The file doesn't exist -- didn't find it there in that folder.
	% Check the entire search path (other folders) for the file by stripping off the folder.
	fullFileNameOnSearchPath = baseFileName; % No path this time.
	if ~exist(fullFileNameOnSearchPath, 'file')
		% Still didn't find it.  Alert user.
		errorMessage = sprintf('Error: %s does not exist in the search path folders.', fullFileName);
		uiwait(warndlg(errorMessage));
		return;
	end
end

%=======================================================================================
% Read in demo image.
rgbImage = imread(fullFileName);
% Get the dimensions of the image.
[rows, columns, numberOfColorChannels] = size(rgbImage)

% Display image.
subplot(2, 2, 1);
imshow(rgbImage, []);
impixelinfo;
axis on;
caption = sprintf('Original Color Image\n%s', baseFileName);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
hp = impixelinfo(); % Set up status line to see values when you mouse over the image.

% Set up figure properties:
% Enlarge figure to full screen.
set(gcf, 'Units', 'Normalized', 'OuterPosition', [0 0.05 1 0.95]);
% Get rid of tool bar and pulldown menus that are along top of figure.
% set(gcf, 'Toolbar', 'none', 'Menu', 'none');
% Give a name to the title bar.
set(gcf, 'Name', 'Demo by ImageAnalyst', 'NumberTitle', 'Off')
drawnow;

% Get the blue mask.
[blueMask, maskedRGBImage] = createMask(rgbImage);

% Display the image.
subplot(2, 2, 2);
imshow(blueMask);
caption = sprintf('Initial Blue Mask Image');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
impixelinfo;
axis('on', 'image');

% Need to get rid of "holes" caused by not using crossed polarizing filters.
binaryImage = imclearborder(~blueMask);
% Display the image.
subplot(2, 2, 3);
imshow(binaryImage);
caption = sprintf('Blue Mask Image');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
impixelinfo;
axis('on', 'image');

% % Compute the sizes of the holes.
% props = regionprops(binaryImage, 'Area');
% allAreas = sort([props.Area], 'descend')
% subplot(2, 2, 4);
% histogram(allAreas);
% grid on;
% title('Histogram of "holes"');

% Get small blobs smaller than 1000 pixels.
binaryImage = bwareafilt(binaryImage, [0, 1000]);
% Use that to fill in "holes" in the blue mask.
blueMask = blueMask | binaryImage;
% Compute the area fraction of blue.
blueAreaFraction = sum(blueMask(:)) / numel(blueMask)

% Display the image.
subplot(2, 2, 3:4);
imshow(blueMask);
caption = sprintf('Final Blue Mask Image, Blue Area Fraction = %.3f', blueAreaFraction);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
impixelinfo;
axis('on', 'image');

% Tell user via a popup message.
uiwait(helpdlg(caption));

function [BW,maskedRGBImage] = createMask(RGB)
%createMask  Threshold RGB image using auto-generated code from colorThresholder app.
%  [BW,MASKEDRGBIMAGE] = createMask(RGB) thresholds image RGB using
%  auto-generated code from the colorThresholder app. The colorspace and
%  range for each channel of the colorspace were set within the app. The
%  segmentation mask is returned in BW, and a composite of the mask and
%  original RGB images is returned in maskedRGBImage.

% Auto-generated by colorThresholder app on 02-Dec-2018
%------------------------------------------------------


% Convert RGB image to chosen color space
I = rgb2hsv(RGB);

% Define thresholds for channel 1 based on histogram settings
channel1Min = 0.529;
channel1Max = 0.652;

% Define thresholds for channel 2 based on histogram settings
channel2Min = 0.559;
channel2Max = 1.000;

% Define thresholds for channel 3 based on histogram settings
channel3Min = 0.000;
channel3Max = 1.000;

% Create mask based on chosen histogram thresholds
sliderBW = (I(:,:,1) >= channel1Min ) & (I(:,:,1) <= channel1Max) & ...
	(I(:,:,2) >= channel2Min ) & (I(:,:,2) <= channel2Max) & ...
	(I(:,:,3) >= channel3Min ) & (I(:,:,3) <= channel3Max);
BW = sliderBW;

% Initialize output masked image based on input image.
maskedRGBImage = RGB;

% Set background pixels where BW is false to zero.
maskedRGBImage(repmat(~BW,[1 1 3])) = 0;

end