how to measure the diameter of a circle

im trying to write a script to measure the diameter of this photo and dont know how to go about it.

[Rik]

I think I would go for an iterative approach where I would increase the circle radius. For every iteration, check if moving the center would decrease the number of white pixels inside the circle. When the number of white pixels starts to increase dramatically, stop the process.

What would be a good cutoff for that increase probably depends on the seed location and how much noise you expect inside the circle.

Once you have the circle, calculating the diameter is of course trivial.

[Image Analyst]

Rik, you put this up in the comments, when I think it should be down in the Answers section.

Answer 1

[Image Analyst]

Try this code:

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

%=======================================================================================
% Have user browse for a file, from a specified "starting folder."
% For convenience in browsing, set a starting folder from which to browse.
startingFolder = pwd;  % or 'C:\wherever';
if ~exist(startingFolder, 'dir')
	% If that folder doesn't exist, just start in the current folder.
	startingFolder = pwd;
end
% Get the name of the file that the user wants to use.
defaultFileName = fullfile(startingFolder, 'n*.*');
[baseFileName, folder] = uigetfile(defaultFileName, 'Select a file');
if baseFileName == 0
	% User clicked the Cancel button.
	return;
end
fullFileName = fullfile(folder, baseFileName);
rgbImage = imread(fullFileName);

% rgbImage = decorrstretch(rgbImage,'tol',0.01);
% imwrite(rgbImage, 'test.png');

% Get the dimensions of the image.
[rows1, columns1, numberOfColorChannels1] = size(rgbImage)
% Display the original image.
subplot(2, 2, 1);
imshow(rgbImage, []);
axis('on', 'image');
caption = sprintf('Original Color Image "%s"', baseFileName);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;
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 Image Analyst', 'NumberTitle', 'Off')

% Do a color segmentation
[mask,maskedRGBImage] = createMask(rgbImage);
% Display the mask image.
subplot(2, 2, 2);
imshow(mask, []);
axis('on', 'image');
caption = sprintf('Color Segmentation');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;

% Delete blobs touching the border.
mask = imclearborder(mask);

% Fill holes, and take the largest blob.
subplot(2, 2, 3);
mask = imfill(mask, 'holes');
mask = bwareafilt(mask, 1);

% Display the image.
subplot(2, 2, 3);
imshow(mask, []);
axis('on', 'image');
title('Final Mask Image', 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;
hp = impixelinfo(); % Set up status line to see values when you mouse over the image.

% Measure the equivalent circular diameter.
props = regionprops(mask, 'EquivDiameter', 'Centroid');
% Plot these over the upper left image.
subplot(2, 2, 1);
hold on;
xCentroid = props.Centroid(1)
yCentroid = props.Centroid(2)
diameter = props.EquivDiameter
viscircles([xCentroid, yCentroid], props.EquivDiameter/2);
plot(xCentroid, yCentroid, 'r+', 'LineWidth', 2, 'MarkerSize', 100);

% Mask the image using bsxfun() function to multiply the mask by each channel individually.
maskedRgbImage = bsxfun(@times, rgbImage, cast(mask, 'like', rgbImage));
% Display the image.
subplot(2, 2, 4);
imshow(maskedRgbImage, []);
axis('on', 'image');
caption = sprintf('Masked Image.  Spot Diameter = %.2f', diameter);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;
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 Image Analyst', 'NumberTitle', 'Off')

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 19-Nov-2019
%------------------------------------------------------


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

% Define thresholds for channel 1 based on histogram settings
channel1Min = 0.000;
channel1Max = 0.175;

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

% Define thresholds for channel 3 based on histogram settings
channel3Min = 0.512;
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

Answer 2

[Walter Roberson]

bwareafilt to remove the small noise objects inside the circle. Invert the result so the open area becomes white. bwareafilt for the largest area. Now regionprops equivalent diameter.

[Victor Falola]

I should have mentioned that i am new to matlab so not sure how to use these functions, i would very much appreiciate it if you could write a simple code or send me a link that shows how to use that function. Thanks

[Rik]

Have you tried reading the documentation pages that Walter mentioned?

[Victor Falola]

where?

[Image Analyst]

Something like (untested):

binaryImage = bwareafilt(~binaryImage);
binaryImage = imfill(binaryImage, 'holes');
binaryImage = bwareafilt(binaryImage, 1);
props = regionprops(binaryImage, 'EquivDiameter');
diameter = props.EquivDiameter

There might be a better way. I highly doubt that the best first step was to take an edge image. Can you post the original gray scale image.

[Victor Falola]

Answer 3

[Image Analyst]

See attached demo.

[ren chongyu]

Its very good but how about this pic, the contrast ratio of the picture is low and I cant use your method to determine the circle center and diameter

[Image Analyst]

ren, you need to create a background image and divide by it. See attached.

clc;    % Clear the command window.
close all;  % Close all figures (except those of imtool.)
imtool close all;  % Close all imtool figures if you have the Image Processing Toolbox.
clear;  % Erase all existing variables. Or clearvars if you want.
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 16;
% Check that user has the Image Processing Toolbox installed.
hasIPT = license('test', 'image_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;
baseFileName = 'rsz_122.png';
% 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
grayImage = imread(fullFileName);
% Get the dimensions of the image.  numberOfColorBands should be = 3.
[rows, columns, numberOfColorBands] = size(grayImage);
if numberOfColorBands > 1
	grayImage = rgb2gray(grayImage);
end
% Display the original color image.
subplot(2, 3, 1);
imshow(grayImage);
axis on;
title('Original Image', 'FontSize', fontSize);
impixelinfo;
% Enlarge figure to full screen.
set(gcf, 'Units', 'Normalized', 'Outerposition', [0, 0.1, 1, 0.9]);

subplot(2, 3, 2);
imhist(grayImage);
grid on;
title('Histogram of Gray Scale Image', 'FontSize', fontSize);

% Illumination is not the same from top to bottom so we need to do a background correction:
verticalProfile = mean(grayImage(:, 1:100), 2);
% Smooth profile a bit with movmean.
verticalProfile = movmean(verticalProfile, 9);
subplot(2, 3, 3);
plot(verticalProfile, 'b-', 'LineWidth', 2);
grid on;
xlabel('Row', 'FontSize', fontSize);
ylabel('Mean Gray Level', 'FontSize', fontSize);
% Fit a cubic through it
coefficients = polyfit(1 : rows, verticalProfile, 3);
yFit = polyval(coefficients, 1 : rows);
hold on;
plot(yFit, 'r-', 'LineWidth', 2);
legend('Smoothed Vertical Profile', 'Fit', 'Location', 'southeast');
title('Vertical Profile Through Background', 'FontSize', fontSize);
hold off;
% Normalize 0-1 so we can divide by it.
yFit = yFit / max(yFit);
% Create background image.
backgroundImage = repmat(yFit', [1, columns]);
subplot(2, 3, 4);
imshow(backgroundImage);
axis on;
title('Background Image', 'FontSize', fontSize);
impixelinfo;

% Divide the original gray scale image by the background image
grayImage = uint8(double(grayImage) ./ backgroundImage);

% Update the histogram.
subplot(2, 3, 2);
imhist(grayImage);
grid on;
title('Histogram of Gray Scale Image', 'FontSize', fontSize);

% Get the binaryImage
binaryImage = grayImage >= 97; % imbinarize(grayImage);
% There is a white frame around it.  Get rid of that.
binaryImage = imclearborder(binaryImage);
% Display the image.
subplot(2, 3, 5);
imshow(binaryImage);
axis on;
title('Initial Binary Image', 'FontSize', fontSize);

subplot(2, 3, 6);
% It seems to be one connected giant blob.
% Use morphological opening to dissolve the tendrils that
% connect the large circular blob to all the other noise lines.
se = ones(5);
binaryImage = imopen(binaryImage, se);
binaryImage = bwareafilt(binaryImage, 1);
binaryImage = bwareafilt(binaryImage, 1);
imshow(binaryImage);
binaryImage = imfill(binaryImage, 'holes');
imshow(binaryImage);
binaryImage = bwareafilt(binaryImage, 1);
imshow(binaryImage);
axis('on', 'image');
% Now make the measurements.
props = regionprops(binaryImage, 'EquivDiameter', 'Centroid');
diameter = props.EquivDiameter
xCenter = props.Centroid(1);
yCenter = props.Centroid(2);
% Show the circle over the image
viscircles([xCenter, yCenter], diameter/2);
hold on;
plot(xCenter, yCenter, 'r+', 'LineWidth', 2, 'MarkerSize', 50);
caption = sprintf('Filled Binary Image.  The equivalent diameter is %f pixels.', diameter);
title(caption, 'FontSize', fontSize);
message = sprintf('The equivalent diameter is %f pixels.', diameter);
msgbox(message);

[ren chongyu]

Thank you very much sir, you saved my life : )

[ren chongyu]

you are such a genius

[ren chongyu]

Hello sir, i tried your method on this picture, but it doesnt work, no matter how i adjust the picture(increase contrast, reduce the noise first or even stright up draw a circle on it), do you know the reason and how can i fix it? thank you very much!

[Image Analyst]

If you don't have a background gradient (like that image), then just skip the background correction.