Maximum ball transform

Question

Consider a finite uniform grid $G$ in three dimensions with a function $f$ mapping integer grid positions $p$ to a boolean value $f(p)$ (i.e., a black/white volume image.)

A ball in $f$ is a set of grid points with center $c$ and radius $r$ where $f(p)$ is true for all $p \in G$ with euclidean distance $|c-p| \leq r$.

The maximum ball transform $MB$ finds for each grid point $p$ the radius $r$ of the largest ball in $f$ containing that point. (Note: This is probably not the actual, commonly used name for this. However I don't know what else to call it.)

Are there good algorithms to compute this transform?

---

I assume this is related to the (euclidean) distance transform. (Find for each grid point $p$ the radius of the largest disk in $f$ with center $p$.)

$DT(p) = min_{x \in f} |p-x|$

$MB(p) = max_{|x-p| < DT(x)} DT(x)$

In the literature I can find algorithms that compute the euclidean $DT$ (and its reverse) in linear time. There are also algorithms to compute the medial axis transform $MAT$ (the set of all center points of the maximum balls).

Is there a faster algorithm than computing the $MAT$, iterating over all found maximum balls and naively tracking the largest one for each grid point?

---

Edit: It seems this is sometimes known as Thickness Distribution or Euclidean Opening Function.

A recent paper computes an approximation: Fast and Accurate Approximation of the Euclidean Opening Function in Arbitrary Dimension (D. Coeurjolly, 2010), (presentation)

Answer 1

This answer has a problem; see comments below.

---

Following up on my comments, here is the rest of the full answer:

First, take a look at the 2008 paper Euclidean Skeletons of Digital Image and Volume Data in Linear Time by the Integer Medial Axis Transform by Hesselink and Roerdink.

From my reading of your question, Section 2 of the paper is answering your question fairly explicitly. In particular, the process they use is essentially a scan across grid points.

Additionally, I wouldn't worry about the runtime being excessive. The algorithms involved are all linear time, and in fact experimental results from the paper show that the $IMA$ transform can be computed in a matter of seconds on a modern computer, even on data sets in the tens of millions. An especially nice piece of this: It appears $MB$ is computed along the way (again, Section 2), so if you simply stop the algorithm at that point, you should have an even faster runtime than "a matter of seconds" in practice.