霍夫圆检测的matlab实现

close all

clc

clear

img=imread('cir1.bmp');
% img=imread('imagefff.bmp');
img=imcrop(img);

if ndims(img) == 3
    img = rgb2gray(img);
end
% figure,imhist(img)
figure,imshow(img),title('ROI')

I = im2uint8(img(:));
num_bins = 256;
counts = imhist(I,num_bins);
p = counts / sum(counts);
omega = cumsum(p);
mu = cumsum(p .* (1:num_bins)');
mu_t = mu(end);

sigma_b_squared = (mu_t * omega - mu).^2 ./ (omega .* (1 - omega));
maxval = max(sigma_b_squared);
isfinite_maxval = isfinite(maxval);

if isfinite_maxval
    idx = mean(find(sigma_b_squared == maxval));
    threshold = round(idx - 1);
    level = (idx - 1) / (num_bins - 1);    
else
    level = 0.0;
end
th=[threshold,level]
    
bw1 = im2bw(img, level);
bw2=edge(bw1,'log');
figure,imshow(bw2); 
title('ostu阈值分割 & sobel edge');
bw2 = bwareaopen(bw2,50);
figure,imshow(bw2);
title('bwareaopen=50');
hold on
[m n]=size(img);
str= sprintf('threshold=%2.2d',threshold);
text(20,20,str,'Color','red');

% [centers1, r1] = imfindcircles(img,[30 150],'ObjectPolarity','bright')
% str1 = sprintf('r1=%2.3f',r1);
% text(centers1(1)-20,centers1(2)-30,str1,'Color','White');

[centers2, r2] = imfindcircles(img,[30 150],'ObjectPolarity','dark','Method','twostage')
str2 = sprintf('r2=%2.3f',r2);
text(centers2(1)-20,centers2(2)-30,str2,'Color','White');
hold off