Understanding the "threshold" parameter of RadonTransformCircles

[anntzer]

I am trying to understand how to set the "threshold" parameter of RadonTransformCircles. My understanding is that this parameter sets a threshold on topographic prominence (via WatershedMaxima) in the Radon transform, but some questions remain. In particular, does that imply that the parameter should scale wrt. image intensity? (e.g., does it make sense to pick it twice higher for an image with twice higher intensity -- in which case the default value of "1" is arbitrary?) How does this scaling interact e.g. with the "normalize" option, which if I understand correctly will rescale the values in the Radon transform?

[crisluengo]

The Radon transform convolves the image with a circle whose integral is the length of the perimeter. With the normalize option its integral is 1 (it is the convolution kernel that is normalized, not the image). Either way the peak intensities in the Radon space scale with image intensities.

However, compared to the Hough transform, the Radon space is very smooth, and there are many fewer spurious local maxima. That is why the default threshold value is rather low. It really only comes into play when the circles are not really circular.

[crisluengo]

I think you're right! The maxima detection with WatershedMaxima works a lot better when connectivity is 2 rather than 1. We're getting artifacts here. With this change, your image yields exactly one local maximum.

With a much smaller threshold you do get spurious local maxima, but they don't form suspicious patterns like the one you show here.

[crisluengo]

Updated in 0ff10d8.

[anntzer]

Thank you for the very quick fix.

I guess I'll continue asking my circle-detection related issues on this thread (but I can also open another one if you prefer). Consider the following example:

import diplib as dip
import matplotlib as mpl
from matplotlib import pyplot as plt
import numpy as np


ii, jj = np.mgrid[-60:61, -60:61]
disks = 256 * (
    np.exp(-((ii+15)**2+(jj+15)**2)/(2*7**2))
    + np.exp(-((ii-15)**2+(jj-15)**2)/(2*7**2))
)
circles = dip.GradientMagnitude(disks)
params, rcps = dip.RadonTransformCircles(
    circles, slice(3, 10), options={"normalize", "correct", "filled"})

axs = plt.figure().subplots(1, 2)
axs[0].imshow(circles)
axs[1].imshow(np.asarray(params).max(0))
for ax in axs:
    for rcp in rcps:
        ax.add_patch(mpl.patches.Circle(rcp.origin, rcp.radius, ec="tab:red", fc="none"))
plt.show()

I don't really mind the two spurious "small" circles (or rather, I can understand that a small extra maximum is found on each circle's edge); I am more annoyed by the spurious detection between the two circles, that essentially tries to draw a circle that's tangent to the other two. A better choice of threshold would get rid of it, but this would not work if there are other circles in the image with wildly different intensities. Do you know a good way to avoid these spurious "tangent circles"? I guess one could check the gradient direction on the circle's edge (I am really trying to detect disks, not circles, so I could check whether the gradient is approximately radial), but perhaps there's something better available?

[crisluengo]

The typical solution is to require a minimum intensity for the peak in the Radon space. But indeed, that assumes the circles are of known and similar intensity.

Another solution could be to verify each detection by reading pixel values along the circle and seeing if they are uniform or vary strongly.

But yes, this is not a trivial problem.

[anntzer]

Thanks. I was mostly trying to make sure I didn't miss anything obvious here.