Maximum Inscribed Circle using Distance Transform

Approximately computes the largest inner circle of a contour/region using distance transform

Description

Maximum Inscribed Circle
 
Or in other words, "largest inner circle" , "maximum empty circle" etc.

This is a very common problem in computational geometry, and it is not simple to solve efficiently.

Addressing 2D image/contour processing, I couldn't find a good implementation on the web. Generally, the reasonable way of solving this problem is to make use of Voronoi Diagrams, which are generally O(nlogn).

After analyzing the problem a bit, I noticed that the problem can easily and approximately be solved using well-known distance transform.

Here is how:
 
The computational aim can be written as:

(x, y) maximizes r = min_{i} r_{i}
where r_i = ||(x_i, y_i) − (x, y)|| and d_i = r_i − r
     
(x_i, y_i): Pairs data points
(x, y), r : Pair, scalar circle centre and radius
     
In non-mathematical terms:
     
  1. The center of the maximum inscribed circle will lie inside the polygon
  2. The center of such a circle will be furthest from any point on the edges of the polygon.
 
So we seek for the point that lies inside the polygon and has maximal distance to the closest edge. This is exactly the maximum value of the pixel (of DT) that lies inside the contour.
 
Notice that this approach is completely in-precise. It is only pixel-precise and never subpixel accurate. However, unlike optimization approaches, it does guarantee a global convergence. In the case of ambiguity, any of the solutions will be valid.
 
To detect the points inside the region, inpolygon remains very slow. So, I make use of the great code of Darren Engwirda, here. As well as being contained in this package, it can also be downloaded from:
http://www.mathworks.com/matlabcentral/fileexchange/10391-fast-points-in-polygon-test
 
Here are other implemnatations, which are more accurate, but much slower than my approach (only slower in Matalb of course!)
     
http://www.mathworks.com/matlabcentral/fileexchange/2794-cvoronoi
 
Using "Largest pixel that doesn't cross any point" approach:
     http://www.mathworks.com/matlabcentral/newsreader/view_thread/283296

--------

Here is a sample call:

I=imread('hand_contour.png');
[R cx cy]=max_inscribed_circle(I)

The max_inscribed_circle function, finds the boundaries of the image, traces them and retrieves a boundary, where neighboring pixels follow each other. It uses the points on the boundary to compute the maximum inscribed circle

Cheers,

Acknowledgements

The author wishes to acknowledge the following in the creation of this submission:
cvoronoi

• Required Products: Image Processing Toolbox
• MATLAB release: MATLAB 7.9 (2009b)

[osunmbs unduo]

problem solved . tnks, however i want to ask , is the use of bwtraceboundary only for unit8 images cos i guess that was my problem. After converting my logical im back to unit8, d code ran perfectly.
Also, how do i crop out the content of the circle with a link to d initial binary image or the original image.
Thanks for ur time

[Tolga Birdal]

There is nothing wrong with the code. The contour start parameter in bwtraceboundary is hardcoded.

The start point should be a white pixel on the boundary of the image for bwtraceboundary to work. If you are using [Y(1), X(1)] there is very little chance for you to get a result. This is because the point [Y(1), X(1)] probably contains a black pixel instead of a white one.

One common approach to find this pixel automatically is to loop through all the points until you hit a white pixel. As soon as you hit that pixel (x,y) you break out of the loop and input (x,y) to bwtraceboundary.

In the example I didn't bother with this. I left the implementation to the user as the purpose of this algorithm is circle fitting not contour segmentation.

And yes, it is of course orientation invariant (given that you obtain a white pixel for bwtraceboundary).

Regards,

[osunmbs unduo]

@ Tolga;the code is worhing perfectly but pls i'm still expecting ur reply to my latest question tanks.

[Tolga Birdal]

bwtraceboundary works for binary images. Typically a logical or appropriate uint8 image will do (see that in the code I'm using logicals). For further info, you could read matlab's manual.

--

For your second question, for a very dirty solution:

You can make a mask having value 1 for the points in the circle (eg. use inpoly), and 0 outside. Then just multiply (or mask) your image with that mask.

[osunmbs unduo]

@Tolga , can u be more explicit on creating the mask. Hoe do i use inpoly for dis. Though i tried using the createMask function but realized that d input to match with represent the total circle corodinates, thus not working.
Kindly shed more light
Thanks S M

[Tolga Birdal]

Check every pixel in the image for whether it is contained in the circle or not. You can use inpoly for this. The polygon should be composed of all the points on the circle. For the pixels contained in the circle make the value 1, if not contained make it 0. There is your mask.

[osunmbs unduo]

@ Tolga, can u pls check the last two lines if correct in creatind the mask to crop out the portion in the circle. Although, the last line is giving error
<??? Attempt to reference field of non-structure array.

Error in ==> max_inscribed_circle at 106
binaryImage =Cin.createMask();>

figure,imshow(BW,[]);
    hold on
    plot(X,Y,'r','LineWidth',2);
   theta = [linspace(0,2*pi) 0];
  hold on
AA=cos(theta)*R+cx; BB=sin(theta)*R+cy;
  plot(AA,BB,'color','g','LineWidth', 2);

Cin =inpoly([X,Y],[AA(:),BB(:)]);
binaryImage = Cin.createMask();

thanks

[Tolga Birdal]

If your circle points are correct, finding the polygon seems ok. Then you can do sth like

ind= sub2ind(size(I) , Y(Cin), X(Cin));
Mask=zeros(size(I));
Mask(ind)=1;

I have no idea about the createMask function.

[osunmbs unduo]

ok, tnks,i'll try dis out

[Sebastian W]