Using regionprops to detect circles in an image

I'm trying to detect circles in an image and find their centerpoints from that to find their deviations from each other. The code I have so far can loosely detect the circles that I want, however the outer circle is not being detected correctly even though in the image it's quite distinct. Right now its using regionprops, imfindcircles did not detect any circles at all. Where in my code can I improve things? I included the original picture for reference.

pic = imread("piercedc.jpg"); % Load image 
pic_gray = rgb2gray(pic); % Converts to grayscale 
pic_contrast = imadjust(pic_gray); % Improves contrast  
pic_filtered = medfilt2(pic_contrast); % Filter removes noise  
imshow(pic_filtered); 

%% BINARY IMAGE  
binarythreshold = graythresh(pic_filtered); % Calculate the global threshold using Otsu's method  
BW = imbinarize(pic_filtered, binarythreshold); % Binarize the image using calculated threshold  
inverted = imcomplement(BW); % invert image  
BW6 = bwareaopen(inverted, 50000); %Eliminate small black specks  
% Show Original Grayscale vs Binary  
figure; 
subplot(1,2,1); 
imshow(pic_filtered); 
title('Original Grayscale Image'); 
subplot(1,2,2); 
imshow(BW6); 
title('Binary Image'); 

%% EXTRACT ELEMENTS  
% Region properties 
stats = regionprops(BW6, 'Centroid', 'Area', 'Eccentricity', ... 
    'EquivDiameter', 'MajorAxisLength', 'MinorAxisLength'); 
% Filter by circularity 
circular_regions = stats([stats.Eccentricity] < 0.5); % Lower eccentricity = more circular 
% Sort circular regions by area (largest to smallest) 
[~, sortIdx] = sort([circular_regions.Area], 'descend'); 
sortedRegions = circular_regions(sortIdx); 
% Extract centroids 
centroids_sorted = cat(1, sortedRegions.Centroid); 
% Assign names based on area 
cartridge = centroids_sorted(1,:); % Largest 
piercing = centroids_sorted(2,:); % Second largest 
%% DISPLAY OVER ORIGINAL IMAGE WITH ELLIPSE OUTLINES 
figure; 
imshow(pic_histeq); hold on; 

Unrecognized function or variable 'pic_histeq'.

title('Detected Elliptical Regions Overlaid'); 
% Plot centroids and text 
plot(cartridge(1), cartridge(2), 'r*'); 
text(cartridge(1)+10, cartridge(2), 'cartridge', 'Color', 'b'); 
plot(piercing(1), piercing(2), 'r*'); 
text(piercing(1)+10, piercing(2), 'piercing', 'Color', 'b'); 
% Draw Elliptical Outlines 
theta = linspace(0, 2*pi, 100); 
% Draw outlines differently for outer (circle) and inner (ellipse) 
theta = linspace(0, 2*pi, 100); 
for i = 1:length(sortedRegions) 
    region = sortedRegions(i); 
    center = region.Centroid; 
    if i == 1 
        % Outer circle (cartridge) – assume perfect circle 
        radius = region.EquivDiameter / 2; 
        x = center(1) + radius * cos(theta); 
        y = center(2) + radius * sin(theta); 
        plot(x, y, 'g-', 'LineWidth', 2); 
    elseif i == 2 
        % Inner ellipse (piercing) 
a = region.MajorAxisLength / 2; 
b = region.MinorAxisLength / 2; 
x = center(1) + a * cos(theta); 
y = center(2) + b * sin(theta); 
plot(x, y, 'm-', 'LineWidth', 2); % Use different color for clarity 
end 
end 
%% CALCULATE CENTROID DISTANCES  
center_error = sqrt((piercing(1) - cartridge(1))^2 + (piercing(2) - cartridge(2))^2); % Uses Euclidean distance formula 
fprintf('Center error (distance between centroids): %.2f pixels\n', center_error); % Prints center Error  
% Find Relative Error  
outer_radius = sortedRegions(1).EquivDiameter / 2; 
relative_error = center_error / outer_radius; 
fprintf('Relative center error: %.2f%% of outer radius\n', relative_error * 100); 

[Image Analyst]

We don't see the outer circle and inner ellipse on your image. We can't run your code. Please attach pic_histeq so we can run your code and see the plots. Do you want the circle to run on the outside of the "gear" teeth, or at the bottom of the gear teeth? Or are you looking at the groove scored just inside the outer rim? Is the inner ellipse the region with the foil punched out?

[Rebecca]

This is what I want the image to look like, except I want the outer green line to be able to match up. Also pic_histeq wasn't supposed to be included, this is the code fixed. Thank you!

%% LOAD AND EDIT IMAGE
pic = imread("piercedc.jpg");                                                  % Load image
pic_gray = rgb2gray(pic);                                                   % Converts to grayscale
pic_contrast = imadjust(pic_gray);                                          % Improves contrast 
pic_filtered = medfilt2(pic_contrast);                                      % Filter removes noise 
imshow(pic_filtered);
%% BINARY IMAGE 
binarythreshold = graythresh(pic_filtered); % Calculate the global threshold using Otsu's method  
BW = imbinarize(pic_filtered, binarythreshold); % Binarize the image using calculated threshold 
inverted = imcomplement(BW); % invert image  
BW6 = bwareaopen(inverted, 50000); %Eliminate small black specks  
% Show Original Grayscale vs Binary 
figure;
subplot(1,2,1);
imshow(pic_filtered);
title('Original Grayscale Image');

subplot(1,2,2);
imshow(BW6);
title('Binary Image');
%% EXTRACT ELEMENTS 
% Region properties
stats = regionprops(BW6, 'Centroid', 'Area', 'Eccentricity', ...
    'EquivDiameter', 'MajorAxisLength', 'MinorAxisLength');
% Filter by circularity
circular_regions = stats([stats.Eccentricity] < 0.5);  % Lower eccentricity = more circular
% Sort circular regions by area (largest to smallest)
[~, sortIdx] = sort([circular_regions.Area], 'descend');
sortedRegions = circular_regions(sortIdx);
% Extract centroids
centroids_sorted = cat(1, sortedRegions.Centroid);
% Assign names based on area
cartridge = centroids_sorted(1,:);  % Largest
piercing  = centroids_sorted(2,:);  % Second largest
%% DISPLAY OVER ORIGINAL IMAGE WITH ELLIPSE OUTLINES
figure;
imshow(pic_filtered); hold on;
title('Detected Elliptical Regions Overlaid');
% Plot centroids and text
plot(cartridge(1), cartridge(2), 'r*');
text(cartridge(1)+10, cartridge(2), 'cartridge', 'Color', 'w');
plot(piercing(1), piercing(2), 'r*');
text(piercing(1)+10, piercing(2), 'piercing', 'Color', 'w');
% Draw Elliptical Outlines
theta = linspace(0, 2*pi, 100);
% Draw outlines differently for outer (circle) and inner (ellipse)
theta = linspace(0, 2*pi, 100);
for i = 1:length(sortedRegions)
    region = sortedRegions(i);
    center = region.Centroid;
    if i == 1
        % Outer circle (cartridge) – assume perfect circle
        radius = region.EquivDiameter / 2;
        x = center(1) + radius * cos(theta);
        y = center(2) + radius * sin(theta);
        plot(x, y, 'g-', 'LineWidth', 2);
    elseif i == 2
        % Inner ellipse (piercing)
        a = region.MajorAxisLength / 2;
        b = region.MinorAxisLength / 2;
        x = center(1) + a * cos(theta);
        y = center(2) + b * sin(theta);
        plot(x, y, 'm-', 'LineWidth', 2); % Use different color for clarity
    end
end
%% CALCULATE CENTROID DISTANCES 
center_error = sqrt((piercing(1) - cartridge(1))^2 + (piercing(2) - cartridge(2))^2); % Uses Euclidean distance formula
fprintf('Center error (distance between centroids): %.2f pixels\n', center_error); % Prints center Error 
% Find Relative Error 
outer_radius = sortedRegions(1).EquivDiameter / 2;
relative_error = center_error / outer_radius;
fprintf('Relative center error: %.2f%% of outer radius\n', relative_error * 100);

[Matt J]

@Rebecca you have that now, don't you? If the outer circle was the missing piece, then my answer gives it to you.

Answer 1

[Matt J]

For this, you will need to download,

https://www.mathworks.com/matlabcentral/fileexchange/87584-object-oriented-tools-to-fit-plot-conics-and-quadrics

pic = imread("piercedc.jpg"); % Load image 
pic_gray = rgb2gray(pic); % Converts to grayscale 
pic_contrast = imadjust(pic_gray); % Improves contrast  
pic_filtered = medfilt2(pic_contrast); % Filter removes noise  

B=bwboundaries( bwconvhull( imbinarize(pic_filtered)) );

circ=circularFit(fliplr(B{1})')

circ = 
  circularFit with properties:

    center: [1.0869e+03 1.0379e+03]
    radius: 821.3029

imshow(pic_filtered); hold on
showfit(circ); hold off

●