MATLAB realization of Hough Circle detection

Source:

M file 1:

function Newpic=hough_circle_main (oldpic)
Close all
Clear
Clc


[Filename,pathname]=uigetfile ({' *.bmp;*.jpg;*.tif;*.gif '}, ' newpic ');
If IsEqual (filename,0)
Disp (' User selected Cancel ')
Else
Disp ([' User selected: ', FullFile (pathname,filename)])
End
Root=fullfile (Pathname,filename);
Oldpic=imread (root);
I=oldpic;


% Close All
% clear
% CLC
%
% i=imread (' circle.bmp ');
[M,n,k] = size (I);
If k>1
i = Rgb2gray (i);
End
BW = Edge (I, ' Sobel ');
Step_r = 10;
Step_angle = 0.1;% 5.7296 degrees
MINR = 100; %need to modify by fact
MAXR = 200; %need to modify by fact
Thresh = 0.5; % 76.8000
[Hough_space,hough_circle,para] = hough_circle (Bw,step_r,step_angle,minr,maxr,thresh);

Subplot (221), Imshow (I), title (' original ')
Subplot (222), Imshow (BW), title (' Edge ')
Subplot (223), Imshow (hough_circle), title (' Test result ')

M File 2:

function [Hough_space,hough_circle,para] = hough_circle (bw,step_r,step_angle,r_min,r_max,p)

% %%%%%%%%%%%%%%%%%%%%%%%%%%
% input
% BW: two value image;
% Step_r: Round radius Step for detection
% Step_angle: Angle step in radians
% R_min: Minimum circle radius
% R_max: Maximum circle radius
% P: The maximum value of p*hough_space is the threshold value, p takes 0, 1 number


% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% output
% hough_space: Parameter space, H (a,b,r) indicates the number of points on a circle with a radius of R (A, B)
% HOUGH_CIRCL: Two value image, detected circle
% para: The center and radius of the detected circle

[M,n] = size (BW);
Size_r = Round ((r_max-r_min)/step_r) +1;
Size_angle = round (2*pi/step_angle);
Hough_space = zeros (M,n,size_r);
[Rows,cols] = find (BW);
Ecount = size (rows);

% Hough Transform
% The image space (x, y) corresponds to the parameter space (A,B,R)
% A = X-r*cos (angle)
% B = y-r*sin (angle)
For I=1:ecount
For R=1:size_r
For K=1:size_angle
A = round (rows (i)-(r_min+ (r-1) *step_r) *cos (K*step_angle));
B = Round (cols (i)-(r_min+ (r-1) *step_r) *sin (K*step_angle));
if (a>0&a<=m&b>0&b<=n)
Hough_space (a,b,r) = Hough_space (a,b,r) +1;
End
End
End
End

% search for aggregation points exceeding threshold values
Max_para = Max (max (hough_space));
index = find (HOUGH_SPACE&GT;=MAX_PARA*P);
Length = size (index);
Hough_circle=zeros (M,n);


For I=1:ecount
For K=1:length
PAR3 = Floor (index (k)/(M*n)) +1;
PAR2 = Floor ((index (k)-(par3-1) * (m*n))/m) +1;
PAR1 = index (k)-(par3-1) * (m*n)-(par2-1) *m;
if ((Rows (i)-par1) ^2+ (cols (i)-par2) ^2< (r_min+ (par3-1) *step_r) ^2+5& ....
(Rows (i)-par1) ^2+ (cols (i)-par2) ^2> (r_min+ (par3-1) *step_r) ^2-5)
Hough_circle (Rows (i), cols (i)) = 1;
End
End
End

% Print Results
For K=1:length
PAR3 = Floor (index (k)/(M*n)) +1;
PAR2 = Floor ((index (k)-(par3-1) * (m*n))/m) +1;
PAR1 = index (k)-(par3-1) * (m*n)-(par2-1) *m;
PAR3 = r_min+ (par3-1) *step_r;
fprintf (1, ' k=%d Center (%d,%d) Radius%d\n ', K,PAR2,PAR1,PAR3);
Para (:, k) = [PAR1,PAR2,PAR3];
End
End