Extracting "endpoints" from a skeleton image to enable a circle to be defined.

Hi, I am trying to find the centre and radius of the circle (whilst ignoring the cross) so I can locate the centre and display a circle on the image:

After reading another matlab central thread, IA has suggested skeletonising and finding the endpoints, so I am trying this approach

So after binarising my raw tiff image and filling in the holes, I have created the skeleton.

But Im not sure how to extract say the 8 most extreme endpoints and use this to define a circle (or fit a circle to)

figure
subplot(1,4,1)
myImshow(IM,2)         %IM is my raw image  
title(['Raw Image'])

subplot(1,4,2)
imshow(BI,[0 1])       %BI is my binarised image
title(['Binary: Thresh=',num2str(level,'%.1f')])

%Fill holes
subplot(1,4,3)
bw = imfill(BI, 'holes');
imshow(bw,[0 1]) 
title(['Fill Holes'])

bw3 = bwmorph(bw2,'skel',Inf);
imshow(bw3,[0 1]); hold on
title(['Skeletonise'])

ep=bwmorph(bw2,'endpoints');

What to do with ep?

I have included the fill holes image incase its useful.

Thanks

Jason

Answer 1

[Image Analyst]

I don't remember seeing that image and I don't think I would have recommended that, anything with skeletonization and endpoints unless I totally misunderstood. What I'd for this image is to call bwconvhull() on it and then call regionprops to get the circle. Something like (untested):

mask = bwconvhull(mask, 'union');
props = regionprops(mask, 'Centroid', 'EquivDiameter');
viscircles(props.Centroid, props.EquivDiameter/2);

There are other ways but with only two or three lines of code I think this is the simplest. See if that works. Note that since the edge is fuzzy/blurry and there is apparently no ground truth answer to check accuracy against, you just have to see if whatever center and radius you get from whatever method suits your needs. It could be that consistency is more important than accuracy.

[Jason]

Both this and Rik's answer are much better than the skeleton idea. Im really sorry but I can't choose between I.A & Rik to accept the answer

[Rik]

Feel free to accept either and give a vote to the other if you think both would suit your needs.

fprintf('fminsearch: C=[%.1f %.1f], r=%.1f\n',[x,y,r])
fprintf('bwconvhull: C=[%.1f %.1f], r=%.1f\n',[props.Centroid props.EquivDiameter/2])
%returns:
fminsearch: C=[79.5 88.3], r=74.3
bwconvhull: C=[79.5 90.1], r=78.4

As you can see the centroids are very close and the radius is slightly smaller for my method (which is expected)

Answer 2

[Rik]

This can probably be improved, but this will find an approximate circle. It will slightly underestimate the radius, as it tries to maximize the number of correctly labeled pixels. You may find the perfomance improves if you manage to fill in that crosshair.

Watch out for the differences in convention for the orientation between image space and plotting space: y is flipped. The code below works also if you extend your masked image in one direction, so I'm reasonably confident I didn't make a mistake here.

%convert image back to binary
mask=imread('im154.png');mask=mask(:,:,1)>128;
[x,y,r]=fit_circle(mask);

t=linspace(0,2*pi,200);
x_c=r*cos(t)+x;
y_c=r*sin(t)+y;

figure(1),clf(1)
imshow(mask),hold on
plot(x_c,y_c,'r-*')
plot(x,y,'ro')

function [x,y,r]=fit_circle(mask)
[Y,X]=ndgrid(1:size(mask,1),1:size(mask,2));
count_of_correct_pixels=@(y,x,r) sum(sum( (sqrt((Y-y).^2+(X-x).^2)<=r) == mask));

%initialize as centered
y=size(mask,1)/2;x=size(mask,2)/2;
r=mean([y x]);

p_fitted=fminsearch(@(p) -count_of_correct_pixels(p(1),p(2),p(3)),[y,x,r]);

y=p_fitted(1);
x=p_fitted(2);
r=p_fitted(3);
end