Detect circles with various radii in grayscale image via Hough Transform

Detect circular shapes in a grayscale image. Resolve their center positions and radii.

Description

Circular Hough transform based on the gradient field of an image.
1. Operates on grayscale images, NOT B/W bitmaps.
2. NO loops in the implementation of Circular Hough transform, which means faster operation but at the same time larger memory consumption.
3. Detection of radii works for concentric circles as well.

• Acknowledgements: This submission has inspired the following:; Coin Segmentation
• Required Products: Image Processing Toolbox; Signal Processing Toolbox
• MATLAB release: MATLAB 7.0.4 (R14SP2)

[e d]

 

[Gijs H]

Excellent work, well documented, easy to use.
Thanks!

[Iftach Nachman]

Excellent! Works amazingly well in detecting yeast cells in microscopy images.

[ali volkan]

Perfect

[Nurettin Senyer]

Thank you for EXCELLENT work. Surely, well document and easy to use.

[veli birol]

clearly and fast.

[Dany Simard]

Clear and efficient

[Sander L]

This code can be of great help in my project on pollen recognition in microscopic images. Thanks!

There is one problem however. For some(!) images I get the error message:
-----------------------------
??? Attempted to access SgnCv_LMPos(0); index must be a positive integer or logical.

Error in ==> CircularHough_Grd at 630
    cirrad(k) = SgnCv_LMPos(end);
-----------------------------
Can anyone help?

[Erinc D]

Thanks for sharing this file, it helps a lot!
However, I'm having the same strange problem too (see below). I understand that this is an issue about radius range selection, but I cannot predict which values would make it work.

[Tao Peng]

Hi, Sander and Erinc: The problem you have experienced is a bug in the code and it has been fixed (at least I think so :)). The revised version (March 7, 2007) will soon be posted by Matlab Central. Thank you for pointing out the problem.
Although I have requested that all user comments to be forwarded to my email, a couple of messages here seemed to have slipped. Therefore I would appreciate everyone who has questions about the code does send me emails. I can't guarantee a timely and satisfactory response but I'll do my best.

[Sander L]

I've solved it by changing the code around line 650 into this:

% Save the detected radii
if (length(SgnCv_LMPos) ~= 0)
    cirrad(k) = SgnCv_LMPos(end);
    for i_radii = (length(SgnCv_LMPos) - 1) : -1 : 1,
        circen = [ circen; circen(k,:) ];
        cirrad = [ cirrad; SgnCv_LMPos(i_radii) ];
    end
end

[Paul Premakumar]

Kudos to you!Your program is amazing. Your examples were great.
Thanks!

[hani ali]

good article

[Dr Samir badawy]

it straightforward. easy to follow and operate. not tested the programs on the MATLAB yet.

[Robert Keck]

Note that this also requires the Signal Processing Toolbox since it uses the filtfilt function.

[Martin Ohmes]

Excellent code with great documentation - this code saved my butt. It is a little memory intensive but with a few tweaks I got it working great for some detailed image analysis.

[Slobodan Vrhovac]

This code can be of great help in my project on disk recognition in experiments with granular materials.

[Nat Angarita]

Very good code, robust,...the only downside...it provides the values of the estimated radious as whole integers...could not see where I could change this in the code!

[Joan Wang]

Excellent performance!

[Vinh Nguyen]

 

[Ohil K]

I tried out this code, and I dont know if this problem can be solved. When I run the code on an image that has some straight lines too, all along the edge of the straight line, circles are detected. Is there a way to prevent this? I tried changing the threshold but to no avail.

[Lee Borland]

Hi all. This may seem like a silly question, but how does one run an add on file such as this? I am attempting to detect circles in images. So far, the only code I have in matlab is open/close code to enhance my edges, then blob analysis to detect my circles and finally convex hull code to fill in my circles. So I now have a BW image containing 4 white circles, very much like the images in the examples.

Can someone advise as to how i now implement this file for radius calculations? I attempted to copy and paste but got the error: Error: Function definitions are not permitted at the prompt or in scripts.

Obviously missing something silly here, can anyone help??

[Roxy]

I have been getting an error on one of your files when i try running it, there is no function name filtfilt.

[Artesis]

Great work, excellent performance and good results. The algorithm is perfectly customizable by its number of parameters. Thanks, this was the important part of a solution for a problem of weeks.

[Ken Jarman]

Excellent code and helpful comments on test cases. If you don't have IPToolbox or SPToolbox, you can still use this if you create or pick up a third-party version of bwlabel.m (labels blobs of dark/light pixels) and filtfilt.m (fwd/bkwd filter).

[Pavithra]

thanks for sharing

[Samuele Lilliu]

 

[Truc Phan]

Excellent and easy-to-use code with very useful documentation.
Works well in detecting small circles, but less robust for larger ones.
Thanks a lot for the great work.

[Fredrik]

 

[Fredrik]

Very well written code with extremely good examples and very well commented. Can't thank you enough!

[bass sawet]

Excellent code!! THX isn t enough

[Angelo]

Hi! why do i get this error message?

??? Maximum recursion limit of 500
reached. Use set(0,'RecursionLimit',N)
to change the limit. Be aware that
exceeding your available stack space
can
crash MATLAB and/or your computer.

Error in ==>
imbmpinfo>decodeCompression

[Brett Shoelson]

PLEASE NOTE: As written, this code requires the Signal Processing Toolbox in addition to the Image Processing Toolbox noted as a requirement above. (There's a single call to FILTFILT, which one may be able to replace with a couple calls to MATLAB's FILTER function.)

[Pascal]

this is great work.

However, circle detection has been investigated to perfusion. Ellipse detection is crucial, especially ellipses oriented at specific angles.....

[Yassine]

Good work! Thanks a lot!

[Nagaraju Annamdevula]

 

[White Apple]

same problem with Angelo. What does the error message means?

[scholar fathi]

excellent job!!
Thanks a lot!!

[Valentin Romanov]

same problem with Angelo too =((
please if someone knows how avid this - help!

[Karen ortega]

help error :S

[Julie]

I keep getting fewer circles when I increase the radius range
example
range =[1 30] # of circles=10
range =[1 50] # of circles=2

I am changing nothing else, is this a bug or am I doing something wrong?

[Joseph Corcoran]

For my application I am trying to identify circles that have different radii but similar (but not the same) centre points- They are sufficiently different that I would of thought that it would be possible to seperate them.
Anyone know of a way that I can differentiate the two circles?

[Kfir]

For some reason, I can't get this to work. I can get the 'accum' matrix and can visually find the local maxima, but the 'circen' and 'cirrad' remain empty.
Any suggestions as to what I should do or which parameters I need to experiment with?
(and I get an error message many times:
"??? Undefined function or method 'filtfilt' for input arguments of type 'double'.

Error in ==> CircularHough_Grd at 606
    SgnCv = filtfilt( fltr4SgnCv , [1] , SgnCv );"

[Paulo]

nice! thanks!

[Sarah]

 

[Ralvin Kumar]

For my application I am using this CircularHough_Grd function to detect COLORED circles.It works fine for red,green,blue,however, it does not get yellow circles or lighter circles.
Here is the image i'am working on:
http://i.imgur.com/OOrmZ.png

And this is what i get with my current codes:
http://i.imgur.com/tP1BY.jpg

My code is :

rgb = imread('capture4.bmp');
rawimg = rgb2gray(rgb);
 
% rawimg = rgb(:,:,1);
 
 fltr4img = [1 1 1 1 1; 1 2 2 2 1; 1 2 4 2 1; 1 2 2 2 1; 1 1 1 1 1];
 fltr4img = fltr4img / sum(fltr4img(:));
 imgfltrd = filter2( fltr4img , rawimg);
 tic;
 [accum, circen, cirrad] = CircularHough_Grd(imgfltrd, [2 20], 0, 3, 0.1);
 toc;
 figure(1); imagesc(accum); axis image;
 
 figure(2); imagesc(rgb); colormap('gray'); axis image;
 hold on;
 
 plot(circen(:,1), circen(:,2), 'r+');axis image;
 impixelinfo;
 for k = 1 : size(circen, 1),
     DrawCircle(circen(k,1), circen(k,2), cirrad(k), 32, 'r-');
 end
 hold off;
 title(['Raw Image with Circles Detected ', ...
     '(center positions and radii marked)']);
 

[Serena]

 

[Vince]

I'm fairly new to image processing, but I can follow your code halfway decent.

However, I am capturing an image from a webcam, and cant seem to get the function to work. What format should the image be in? I've tried to save it as several formats, but either get a lime green image as the result, or an error saying the image must be 2-dimensional. Any suggestions?

[Vilém Rechtacek]

Hello,

function CircularHough_Grd works very well for images, but I need to process video stream (for example from web camera). I want to detect eyes in video stream. I want to use CircularHough Tansform for detection of iris, becouse iris have circular shape. I assume, that the face fills the entire scanned image space, so I do not need to detect face.

Can anyone help me? When you send me the complete code, it would be best for me.

Thank you.
Vilem Rechtacek

[ketan patil]

I didn't understand following part in code
.
rr_4linaccum = double( prm_r_range ); %prm_r_range contain all the sort values of radius in ascending order
linaccum_dr = [ (-rr_4linaccum(2) + 0.5) : -rr_4linaccum(1) , ...
    (rr_4linaccum(1) + 0.5) : rr_4linaccum(2) ];

lin2accum_aJ = floor( ...
double(grdx(grdmasklin)./grdmag(grdmasklin)) * linaccum_dr + ...
repmat( double(grdmask_IdxJ)+0.5 , [1,length(linaccum_dr)] ) ...
);
lin2accum_aI = floor( ...
double(grdy(grdmasklin)./grdmag(grdmasklin)) * linaccum_dr + ...
repmat( double(grdmask_IdxI)+0.5 , [1,length(linaccum_dr)] ) ...
);
plese tell me

[Fer]

 

[Daniel Delgado]

Hello, I really need your help. I'm trying to detect the circle center of those red ball images but I can't find the correct parameters to make CircularHough_Grd to work. Thank you in advance, it's kinda urgent. It used to work perfectly on 160x120 images but now that those are 640x480 it fails.

Here is my code:

x = 90;

%% Read images
%images[];
for n=1:(numel(x)/3)
    images{n} = imread(sprintf('%d.jpg',n));
end

Temp = imfinfo('1.jpg');
width =Temp.Width;
height = Temp.Height;

%% Detect circles centers
tic;
for n=1:length(images)
    img = images{n};
    img2 = imadjust(img, [0 .7], [0 .5]);
    %img = imresize(img, 2);
    red = img2(:,:,1);
% red = img(:,:,1);
    %red = adapthisteq(red);
    fltr4img = [1 1 1 1 1; 1 2 2 2 1; 1 2 4 2 1; 1 2 2 2 1; 1 1 1 1 1];
% fltr4img = [1 1 1 1 1; 1 4 4 4 1; 1 4 12 4 1; 1 4 4 4 1; 1 1 1 1 1];
% fltr4img = [1 1 1 1 1; 1 4 8 4 1; 1 4 8 4 1; 1 4 8 4 1; 1 1 1 1 1];
% fltr4img = fltr4img / sum(fltr4img(:));
% imgfltrd = filter2( fltr4img , red );
% imgfltrd = filter2( fltr4img , img2(:,:,1) );
    imgfltrd = filter2( fltr4img , red );
% imgfltrd = red;

    %tic;
% [accum, circen, cirrad] = CircularHough_Grd(imgfltrd, [15 80], 8, 10);
% [accum, circen, cirrad] = CircularHough_Grd(imgfltrd, [0 26], 3, 4, 1);
    % [accum, circen, cirrad] = CircularHough_Grd(red, [11 50], 20, 13, 1);
    [accum, circen, cirrad] = CircularHough_Grd(imgfltrd, [90 120], 50, 5, 1);
    %toc

    if(length(circen)==2)
        %circenaccum {n} = circen(1,:)';
        circenaccum{n} = [circen(1,1);circen(1,2)];
        %circenaccum{n} = [(circen(1,1)-(width/2));(circen(1,2))-(height/2)];
    else
        circenaccum{n} = [0;0];
    end
   
    % if any(cirrad <= 0)
    % inds = find(cirrad>0);
    % cirrad = cirrad(inds);
    % circen = circen(inds,:);
    % end

    %imshow(img);
    figure(1); imagesc(accum); axis image;
    title('Accumulation Array from Circular Hough Transform');
% figure(2); imagesc(img2); colormap('gray'); axis image;
% figure(2); imagesc(img); colormap('gray'); axis image;
    figure(2); imagesc(imgfltrd); colormap('gray'); axis image;
% figure(9); surf(accum, 'EdgeColor', 'none'); axis ij;
% title('3-D View of the Accumulation Array');
    hold on;
    plot(circen(:,1), circen(:,2), 'r+');
    % for ii = 1 : size(circen, 1)
    % rectangle('Position',[circen(ii,1) - cirrad(ii), circen(ii,2) - cirrad(ii), 2*cirrad(ii), 2*cirrad(ii)],...
    % 'Curvature', [1,1], 'edgecolor', 'b', 'linewidth', 1.5);
    % end
    for k = 1 : size(circen, 1),
       DrawCircle(circen(k,1), circen(k,2), cirrad(k), 32, 'b-');
    end
    hold off;
    title(['Raw Image with Circles Detected ', ...
         '(center positions and radii marked)']);
% pause;
end
toc

%% Draw center points
figure
circenaccummat = cell2mat(circenaccum);
scatter(circenaccummat(1,:),circenaccummat(2,:))
% plot(circenaccummat(1,:), circenaccummat(2,:),'*b');
axis ij;

[Amalia adiba]

 

[Amalia adiba]

Dear Mr Taopeng

would you please give me a reason, why while i compile from any example on your program, i got this error

??? Maximum recursion limit of 500 reached. Use set(0,'RecursionLimit',N)
to change the limit. Be aware that exceeding your available stack space can
crash MATLAB and/or your computer.

Error in ==> imbmpinfo>decodeCompression

I'm using matlab 7.8 (2009a)
there's no problem while i try to call

rawimg = imread('TestImg_CHT_a3.bmp');

Thank you for yor reply

[Shubham]