Fitting circle inside image. Exit condition for loop

Question

I am trying to write some code to check if my circle fits inside a binary mask.

So, I have a circle with a given center and a radius. I also have a binary image and I can plot it and see that the circle fits inside the mask. The plot is attached:

So what I am trying to do is find the best circle that fits in this binary mask. So, what I have so far is as follows:

import numpy as np
# Start with a random radius
radius = 100
# Keep track of the best radius
best_radius = 0
# Mask is the binary mask
x, y = np.ogrid[:mask.shape[0], :mask.shape[1]]
# distance computation
distance = np.sqrt((x-r)**2 + (y-c)**2)

while True:
    # check if the distance < radius
    m = distance < radius
    # if this is true everything in the circle is in the image
    if np.all[mask[m] > 0]:
        # Update best radius
        best_radius = radius
        # Increase radius for next try
        radius += radius/2
    else:
        # decrease radius
        radius -= radius/2
        if radius <= best_radius:
            break

So, the first issue I am having is the break condition. While this works as in it always gets me a circle, which is inside the mask, it is not always optimum. So, running it on the example from above, I get the following best circle:

As you can see, the radius can be increased still. So, I have a feeling that the way I increment and decrement radius is not right.

The other big question is that this is a very hacky way to do this in a brute barbarian like approach. If anyone has more elegant solution to share, I would be very grateful. Perhaps, some numerical optimization routine could be employed here?

Answer 1

• Compute the distance from each non-masked point to (r,c).
• The maximal radius is the minimum of the distances.

---

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.patches as patches

h, w = 600, 700
r, c = 250, 350

def make_mask():
    x, y = np.ogrid[:h, :w]
    y -= c
    x -= r
    theta = np.arctan2(y, x)
    rad = 200 + 100*(np.sin(theta+np.pi/2)**2) 
    mask = (x**2 + y**2) < rad**2
    return mask

mask = make_mask()
x, y = np.ogrid[:h, :w]
distance_sq = (x-r)**2+(y-c)**2
radius = np.sqrt((distance_sq[~mask]).min())

fig, ax = plt.subplots()
ax.imshow(mask)
ax.add_patch(patches.Circle((c,r), radius=radius))
plt.show()


yields