Digital Geometry: Volumetric Analysis§

author:	David Coeurjolly

Distance Transformation/Voronoi Map§

Problem setting§

Distance transformation

Def.

At each point $p$ of a shape $X$ , we want the smallest distance to $q\in\bar{X}$

$\Rightarrow$ We need

A metric definition on $\mathbb{Z}^n$
Algorithms to compute the DT

Why the DT ?§

Very usefull

To extract metric information from digital objects
Implicit representation of surfaces
Differential estimators form DT and Voronoi Diagram
Representation as a union of balls

…for many applications

Shape description
Motion planning
Image synthesis
Material Sciences
…

Example: curvature from DT§

Idea

Shape $\rightarrow$ Boundary $\rightarrow$ Distance field $\rightarrow$ Implicit representation $f: (x,y)\rightarrow \mathbb{R}$

$f(x,y)&=0\\ \vec{g}(x,y)&=(I_x,I_y)^T\\ \vec{t}(x,y)&=\frac{(-I_x,I_y)^T}{\sqrt{I_x^2+I_y^2}}\\ k&=\frac{-t^THt}{\|\vec{g}\|}, \quad H= \left [ \begin{array}{cc} I_{xx} & I_{xy}\\ I_{yx} & I_{yy} \end{array}\right ]$

Application example 1: Snow Micro-structure Analysis§

Objectives

X-tomographic 3D images of snow samples
We need to extract geometrical information such as
- Density, specific surface
- Curvature map to guide metamorphosis
- Thickness, percolation path, … ( $\Rightarrow$ DT)

_images/Courb08iso_600_300_256p2_000.png

Application example 2: Metallic Foam Description§

Objectives

X-ray tomographic 3D images of metallic/polymer foam (used for catalytic processes)
We need to extract geometrical information such as
- Structural Complexity of the 1D structure
- Thickness, percolation path, thermodynamic diffusion paths, … ( $\Rightarrow$ DT)

Metric in Digital Space§

Def.

Metric $(d,E,F)$ is a map $d: E\times E \rightarrow F$ such that

$\forall p,q\in E\,, d(p,q)\geq 0$ (separation axiom)
$\forall p,q\in E\,, d(p,q) =0 \Rightarrow p=q$ (coincidence axiom)
$\forall p,q\in E\,, d(p,q) = d(q,p)$ (symmetry)
$\forall p,q,r\in E\,, d(p,q) \leq d(p,r) + d(r,q)$ (triangular inequality)

Def.

$(d,E,F)$ is a norm (metric induced by “normed” vector space) iff

$\forall p,q,r\in E\,, d(p,q) = d(p+r,q+r)\geq 0$ (translation invariant)
$\forall p,q,\lambda\in E\,, d(\lambda p, \lambda q) = |\lambda|d(p,q)\geq 0$ (homogeneity)

Examples§

E.g.

$d_1(p,q) = \| p-q\|_1 = \sum_{i=1}^n |p_i - q_i |$
$d_\infty(p,q) = \| p-q\|_\infty = \max_{i=1..n} \{|p_i - q_i |\}$
$d_2(p,q) = \| p-q\|_2 = \sqrt{\sum_{i=1}^n (p_i - q_i)^2}$
Weighted $l_p$ metric: $d_p(p,q) = \| p -q\|_p = \left ( \sum_{i=1}^n w_i|p_i-q_i |^p \right )^{\frac{1}{p}}$

l_p metrics Illustration§

Discrete Metrics§

Definition

We focus on triplets $(d,\mathbb{Z}^n,\mathbb{Z})$

Hence

$d_1$ , $d_\infty$ are Discrete Metrics
$d_2$ is not a discrete metric
$d_2^2$ is not a metric
$\lceil d_2 \rceil$ is a discrete metric
$\lfloor d_2 \rfloor$ is a not discrete metric
$[d_2]$ is not a metric

Hints for last two results use $p(2,3)\, q(-1,-1)\, r(0,0)$ and $p(1,1)\, q(-1,-1)\, r(0,0)$

Chamfer Mask§

Weigthed vector

$M = (\vec{v},\omega)$

Chamfer Mask

Set of weighted vector

$\mathcal{M} = \{ M_i\in \mathbb{Z}^n\times \mathbb{N}^*\}_{1\leq i \leq m}$

with central symmetry
which contains at least a basis of $\mathbb{Z}^n$

Usually, chamfer masks are G-symmetric, i.e. restricted to

$\mathcal{M} = \{ M_i\in \mathcal{G}\times \mathbb{N}^*\}_{1\leq i \leq m}$

with

$\mathcal{G} = \{ (x_1,\ldots,x_n)\in\mathbb{Z}^n\,|\, x_n\geq \ldots \geq x_1\geq 0 \}$

Chamfer Distances§

Chamfer path

$k$ -Path based on vectors from a chamfer mask

$\mathcal{P} =\{ \alpha_1\vec{v}_{i_1}, \ldots, \alpha_k \vec{v}_{i_k} \}$

Length of a chamfer path

$d_\mathcal{M}(\mathcal{P}) = \sum_k \alpha_k\omega_{i_k}$

Chamfer distance

Minimal length of chamfer path between $p$ and $q$

All chamfer distances induced distances, not necessarily norm

Simple examples§

Path based distance

Matrix representation for masks

$\mathcal{M}_{ab}=\{ (a,(0,1)^T) , (b,(1,1)^T) \}$ , $\mathcal{M}_{abc}=\{ (a,(0,1)^T) , (b,(1,1)^T), (c,(2,1)^T)\}$

For example:

$\mathcal{M}_{3,4} = \{ (3,(0,1)^T), (4,(1,1)^T) \}$

(distances must be divided by 3 at the end)

Chamfer balls§

$\Rightarrow$ We need constraints on $\{\omega_i\}$ to induce norms

e.g.

$0 < a \leq b \leq 2a$

$0 < 2a \leq c \leq a+b\quad\text{and}\quad 3b\leq 2c$

Mask Construction§

We construct the mask to approximate the Euclidean Metric

We first fix a set of vectors (usually, Farey fraction vectors in $\mathcal{G}$ )
We find optimal weights to minimize the error (uniformly, average error, …) with respect to on specific configuration
- We minimize the error on the column $x=N$
- We minimize the error on the circle or radius $N$
We approximate optimal weights for integer numbers (plus scaling integer)

Drawbacks

Just an approximation of the Euclidean metric
Isotropic error distribution: error is maximized on specific orientations
Increasing the mask size reduces the errors but increases the computational cost

Distance Transformation algorithm with Chamfer Masks§

Propagation using Dijkstra’s algorithm

Implicitly construct a regular graph from $\mathbb{Z}^n$ vertices and edges taken from $\mathcal{M}$
Use Dijkstra’s like algorithm to propagate distances from background points

$\Rightarrow$ Computation cost in $O(mn\log n)$ for $n$ grid points and $|\mathcal{M}|=m$

Raster Scan Algorithm§

Split the mask into two sub-masks and perform forward/backward scans with “min” operations.

Init: $DT(p) = 0 \quad \text{if} \quad p\not\in X\\ DT(p) = +\infty \quad \text{if} \quad p\in X$
Then: $DT(p) = min( DT(p), min_{(\omega_i,\vec{v}_i) \text{ in sub-mask}} ( DT(p+\vec{v}_i) + \omega_i ))$

$\Rightarrow$ Computational cost in $O(nm)$

Other path-based distances§

Neighborhood sequence

We consider a sequence of Chamfer masks $\mathcal{M}_i$
At each step in the path construction, we consider another mask $\mathcal{M}_i\rightarrow \mathcal{M}_{i+1}$
Idea each mask has directional error, changing masks reduce the error propagation

Example

“Octogonal” distance with infinite sequence $\{ d_1, d_\infty, d_1,\ldots, \}$

Sometimes, explicit forms exist

$d_{oct}(p,q)=\max \left \{ \left \lfloor\frac{2}{3} d_1(p,q) + 1\right\rfloor, d_\infty(p,q) \right\}$

Bit more difficult to ensure that a sequence/weighted sequence induce a metric

Euclidean metric§

Idea

Still consider $(d,\mathbb{Z}^n, \mathbb{R})$ distances but with integer based representations and algorithmic

E.g.

Represent $d_2(p,q)\in\mathbb{R}$ by $d_2^2(p,q)\in \mathbb{Z}$
Represent $d_2(p,q)\in\mathbb{R}$ by vector $(p,q)^T\in \mathbb{Z}^2$
Similarly, $d_p(p,q)\in\mathbb{R}$ by $d_p^p(p,q)\in \mathbb{Z}$ or even $(p,q)^T\in \mathbb{Z}^2$ for $l_p$ metrics

Nice but are there fast algorithms for such exact metrics ?

Separable Approach For Squared Euclidean Distance Transform§

We want to compute (for all $p\in X$ )

$DT_2(p) = \min_{q\in\bar{X}} \{ d_2(p,q)\} =\sqrt{ \min_{q\in\bar{X}} \{ (p_1 - q_1)^2 + (p_2 - q_2)^2\}}$

$DT_2(p) = \sqrt { \min_{q\in\bar{X}} SEDT(p) }$

Separable approach with intermediate map

$g( i,j) = \min_{x} \{ (x-i)^2\}$

$SEDT( p(i,j) ) = \min_{y} \{ (y-j)^2 + g(i,y)\}$

in dimension 3, we would have

$g(i,j,k) = \min_x \{(x-i)^2 \}, h(i,j,k) = \min_y \{(y-j)^2 + g(i,y,k) \}\\ SEDT( p(i,j,k) ) = \min_{z} \{ (z-k)^2 + h(i,j,z)\}$

First Step§

Simple two-scan propagation

$\Rightarrow$ $O(N^2)$ in 2D for NxN image

$\Rightarrow$ $O(N^d)$ in d-D for N^d image

Second Step§

$g( i,j) = \min_{x} \{ (x-i)^2\}$ and $SEDT( p(i,j) ) = \min_{y} \{ (y-j)^2 + g(i,y)\}$

Key-point Lower envelope computation of a set of parabolas

Lower Envelope Computation§

Consider the set of parabolas $\{ (x-k)^2 + g_k \}_{k=1\ldots N}$

Any two parabolas have single point intersection
For lower envelope computation, the intersection point acts as a pivot
- If $P_k$ and $P_{k'}$ are two parabolas with intersection $q$ and $k<k'$
- $P_k$ cannot appear in the lower envelope for abscissa greater than $q$

$\Rightarrow$ Lower envelope computation in $O(N)$ using stack based approach ;)

Overall SEDT Algorithm§

Given a $N^d$ image

Algorithm

First Step: perform two-scan propagation
For all $(d-1)$ remaining dimensions: compute independent lower envelope which are in $O(N)$

$\Rightarrow\quad O(d\cdot N^d)$ algorithm for error free Euclidean metric DT

Generalizations§

Thanks to separability

We have independent 1D problem to solve
Synchronization steps only occurs when we go to the next dimension

Optimal multi-thread implementation

Generalization to toric domains

Useful to characterize periodic structures in arbitrary dimensions

Generalization to other metrics§

Principle

Separable decomposition works for any metric satisfying the monotonicity property:

Def.

We consider $p(x,y)$ , $q(x',y')$ with $x<x'$

$r( x'',O)$ be a point on the x-axis such that $d(p,r) = d(q,r)$

Let $s(u,0)$ be another point on the x-axis A metric $d$ is monotonic if

$u < x'' \implies d(p,s) \leq d(q,s)$

$u > x'' \implies d(p,s) \geq d(q,s)$

Result

All $l_p$ metrics are monotonic
All Chamfer masks induced by norms are monotonic
All path based distances inducing norms with axis-symmetric unit ball are monotonic

$\Rightarrow$ Let’s use the separable approach for other metrics !

Voronoi Diagram§

Definition

Given a set of sites $S=\{ s_i\in \mathbb{R}^d\}$ , the Voronoi Diagram is a decomposition of the space into closed cells ${c_i}$ such that: $Voro_{S}(s_i) = \{ x\in\mathbb{R}^d,\, d(x,s_i) \leq d(x,s_j),\, \forall s_j\in S\}$

Each cell can be further decomposed into sub-dimensional i-facets taking into account cases where $d(x,s_i)= d(x,s_j)$

Voronoi Diagram $\equiv$ Distance Transformation

$DT(p) = d(p,q)\, \text{ with } q\in\bar{X}\text{ such that }p\in Voro_{\bar{X}}(q)$

$\Rightarrow$ Getting the distance value is equivalent to localizing a point in a Voronoi diagram

Separable Voronoi Map§

Input set: $X\subset\mathbb{Z}^2$ , we construct $Voro_{\bar{X}}\cap\mathbb{Z}^2$

Separable Voronoi Map§

Generic Algorithm§

Main Result

For any monotonic metric and an image $[1\ldots n]^d\rightarrow \{0,1\}$ , the Voronoi Map (and the distance transformation) can be obtained by the separable algorithm in $O( d\cdot n^d\cdot (C + H) )$

C: Closest(u,v, p), decide whether u or v is closest to p
H: HiddenBy(u,v,w, 1D-line), decide if Voronoi cells of u and w hide the Voronoi cell of v on the 1D-line

Metric	C	H	Total
$l_2$	$O(1)$	$O(1)$	$\Theta(d\cdot n^d)$
$l_1$	$O(1)$	$O(1)$	$\Theta(d\cdot n^d)$
$l_\infty$	$O(1)$	$O(1)$	$\Theta(d\cdot n^d)$
Exact $l_p$	$O(log(p))$	$O(log(p).log(n))$	$O(d\cdot n^d\log(p)\cdot\log(n))$
Chamfer Norms	$O(log(m))$	$O(log^2(m))$	$O(d\cdot n^d\cdot\log^2(m))$
Neigh. Seq. Norms	O(1)	$O(\log(n))$	$O(d\cdot n^d\cdot\log(n))$

Examples§

$l_2$ $l_6$

Path based Metrics§

Better expected bounds for path based norms

Metric	C	H	Total
Chamfer with adapter	$O(m)$	$O(m\cdot log(m))$	$O(d\cdot m\cdot n^d\cdot\log(n))$
Chamfer Norms	$O(\log(m))$	$O(\log^2(m))$	$O(d\cdot n^d\cdot\log^2(m))$

Similar expected results for neighborhood sequences

Examples§

$l_2$ $l_6$

Subquadratic Algorithm for path based distances§

Path Based Metric and Rational balls§

Notations

Chamfer masks: $\mathcal{M} = \{ (\vec{v}_i,\omega_i) \in \mathbb{Z}^n\times \mathbb{N}^*\}_{1\leq i \leq m}$ (we consider only chamfer masks inducing norms)
Rational ball: $\mathcal{B}_{\mathcal{M}} = Conv\left ( \{ \frac{\vec{v}_i}{\omega_i} \} \right )$ [Normand, Strand,…]
Rational ball faces have normal vector $\mathcal{F}_i$

$\mathcal{M}_{7,8,11,14}$

Distance Evaluation§

[Normand et al.]

$d(O,p) = \max_i \{ \mathcal{F}_i\cdot \vec{Op}\}$

$\Rightarrow$ $O(m)$

Can be generalized to other path based distances to get similar expression

$d(O,p) = \max_i \{ f_i(\mathcal{F}_i\cdot \vec{Op}) \}$

for some function $f_i$ (based Lambek-Moser inverse sequences)

Optimized Distance Evaluation§

Computational Geometry setting

The facet inducing the $\max$ is given by the facet pierced by the straight line $(Op)$
$\Rightarrow$ Ray shooting problem in convex polytopes

Fast Distance computation

Following [Matousek and Schwarzkpof]
$\Rightarrow$ $O(m^{\lfloor d/2\rfloor})$ space/pre-processing and $O(\log m)$ per query

Separable Predicates for Chamfer Masks§

Goal

If we could have predicates in $O(h)$ then we have exact Voronoi Map/DT in $O(d\cdot h\cdot n^d)$
Raster scan is in $O(m\cdot n^d)$

Main Result

$h = O(\log^2m)$

==> First sub-quadratic DT algorithm for Chamfer metrics

Dimension 2§

Key point

Given to points and a straight line, detect the position of the Voronoi edge on the line
we are looking for point $r\in l$ such that $d_\mathcal{M}(p,r) = d_\mathcal{M}(q,r)$

Warm up: Localizing a point§

Question Find the cone at $p$ containing a point $r$

==> Dichotomic/Binary search (thanks to convexity of the metric)

==> $O(\log m)$

Algorithm Overview§

Idea

If we have localized the Voronoi edge point, we are done (find the exact position given by linear system with one unknown)

Step 1: Shrinking $\mathcal{M}_p$ §

ShrinkMp( Mp, Mq )

   if |Mp| == 1
     return the cone in Mp
   else
     Split cones Mp -> { Mp, cone, M-} with |M+|~|M-|
     {v1,v2} = cone
     dp1 = distance d_M(p, v1 intersection l)   //O(1)
     dp2 = distance d_M(p, v1 intersection l)   //O(1)
     dq1 = localize and get the distance of d_M(q, v1 intersection l) //O(log(m))
     dq2 = localize and get the distance of d_M(q, v2 intersection l) //O(log(m))

     c1  = closest point between p and q at v1
     c2  = closest point between p and q at v2

     if (c1 == c2 == GREEN)
      return ShrinkMp(M+)

     if (c1 == c2 == BLUE)
      return ShrinkMp(M-)

      return cone

Correctness

Chamfer norm implies 2 connected sets (blue/green), maybe overlapping
When evaluating distance on the line, the function is convex (-> orientation)

Step 2: Shrinking $\mathcal{M}_q$ and final computation§

Shrinking $\mathcal{M}_q$

Similar algorithm
Can be speed-up using $\mathcal{M}_p$ cone
$\Rightarrow O(\log^2 m)$

Final step

Small computation in each cone
$\Rightarrow O(1)$

Fast computations in higher dimensions§

Basic Idea for $\mathcal{M}$ in $\mathbb{R}^d$

Each ${p,l}$ ${q,l}$ defines a plane ( $P_{pl},\, P_{ql}$ )
$\mathcal{B}_\mathcal{M}\cap P_{pl}$ induces a 2-dimensional polytope
$\Rightarrow$ 2D problem with $O(\log^2 m)$ computational cost

Conclusion

Closest() and HiddenBy() predicates can be implemented in $O(\log^2m)$
Exact Voronoi map/Distance transformation of Chamfer norms using separable approach in

$O(d\cdot \log^2 m\cdot n^d)$

Reverse Distance Transformation§

Reverse Transformation§

Problem setting

Def.

Given a metric $(d,E,G)$ and a set of balls $\mathcal{B}=\{ B_i=(p_i,r_i)\in E\times G\}_{i=1\ldots N}$ , reconstruct the binary shape $X$

$X = \bigcup_{i=1\ldots N} B_i$

Why?

Reverse operation of the Distance Transformation

$\text{If } \mathcal{B}=\left \{(p,DT(p))\, \forall p\in X\right \}, \text{ then } RDT( \mathcal{B} ) = X$
To reconstruct the shape $X$ if we characterize it as a union of balls (e.g. via medial axis)

Bruteforce approach

For $n\times m$ image

$O(Nnm)$

Separable Approach for $l_2$ §

W.l.o.g. we consider $d=2$

Let us denote $p_k=(x_k,y_k)$ for $k=1\ldots N$ , then

$X = \left \{ (i,j)\,|\, \exists k\in[1,N] \, (i - x_k)^2 + (j-y_k)^2 \leq r_k^2\right \}$

Which can be rewritten

$X =\left \{ (i,j)\,|\, \max_{k=1\ldots N}\{ r_k^2 -(i - x_k)^2 - (j-y_k)^2\} >0\right \}$

$\Rightarrow$ Separable decomposition

Start from a map $f: \mathbb{Z}^2\rightarrow \mathbb{Z}$ with $f(x,y) = r_k^2$ if $((x,y),r_k)\in\mathcal{B}$ ( $f(x,y) = 0$ otherwise)

$g(i,j) = \max_{x} \{ f(x,j) - (x-i)^2\}$

$REDT(i,j) = \max_{y} \{g(i,y) - (y-j)^2 \}$

Illustration§

Similar algorithm

Lower Envelope computation of parabolas $\rightarrow$ Upper envelope computation of a set of parabolas
Generalization of arbitrary dimension

$\Rightarrow$ $O(d\cdot n^d)$ separable algorithm for REDT

Associated Structure from Computational Geometry§

Voronoi map –> Power map

Kind of Voronoi diagram with additive power metric. For example the power of a point x w.r.t. ball $(s,r)\in\mathbb{R}^d\times\mathbb{R}$

$\pi(x,(s,r)) = d^2(x,s) - r^2$

Def.

Given a set of weighted sites $S=\{ (s_i,w_i)\in \mathbb{R}^d\times\mathbb{R}\}$ , the Power Diagram is a decomposition of the space into closed cells ${c_i}$ such that: $Power_{S}((s_i,w_i)) = \{ x\in\mathbb{R}^d,\, \pi(x,(s_i,w_i)) \leq \pi(x,(s_j,w_j)),\, \forall (s_j,w_j)\in S\}$

Each cell can be further decomposed into sub-dimensional i-facets taking into account cases where $d(x,s_i)= d(x,s_j)$

_images/Radical_axis_intersecting_circles.png

Power Map§

Idea

REDT of a set of balls $\mathcal{B}$ is equivalent to the construction of $\left(Power_{\mathcal{B}}\cap\mathbb{Z}^d\right )$

Results

The separable algorithm for monotonic metrics can be extended to generate $\left(Power_{\mathcal{B}}\cap\mathbb{Z}^d\right )$
$O( d\cdot n^d\cdot (C + H) )$ computational cost for a large class of metrics

Medial Axis Extraction§

Problem Description§

Alternative Definitions

Self intersection of wavefronts in prairie model ( $\rightarrow$ PDE fromulation)
Center of maximal balls contained in the shape

$B$ is maximal in X if $\not\exists B'\subset X,\, B\subset B'$
Set of balls touching $\partial X$ at least twice
1D topological equivalent of the contour
…

Contact Points based Geometrical Definition§

Voronoi Based Approximation

Shape $\Rightarrow$ Point set approximation $\Rightarrow$ Voronoi Diagram $\Rightarrow$ Medial Axis approximation

Convergence results exists for various classes of Voronoi based medial axis

Maximal Ball based Definition§

Def.

A maximal ball is a ball contained in the shape not entirely covered by another ball contained in the shape

Def.

The medial axis of a shape is the set of maximal ball centers contained in the shape.

Digital Setting

$X\subset\mathbb{Z}^d$
We consider digital balls $B\cap\mathbb{Z}^d$
$\Rightarrow$ Finite set of digital balls contained in $X$
$\Rightarrow$ Medial Axis Extraction $\equiv$ Combinatorial Covering problem

Reversible Encoding of X

$X = \bigcup_{B_i\in MA(X)} B_i$

DT and Digital Medial Axis§

DT as preliminary step

Given $p\in X$ and $r\in\mathbb{R}$ such that Euclidean ball with $B(p,r)\cap X\in X$ , we have

$B(p,r) \subseteq B(p, DT(p))$

(defined for $l_2$ but trivial generalizations to other metrics)

$\Rightarrow$ Set of candidate balls $O(|X|)$

Digital Ball vs. Euclidean Balls§

Covering Test

Let us consider a IsCoveredBy(B,B’) a predicate returning true if $B\subseteq B'$

If B and B’ are Euclidean balls $\Rightarrow$ The predicate is in $O(1)$
If B and B’ are Digital balls $\Rightarrow$ The predicate is in $O(\max(|B|,|B'|))$
$B\subset B'$ $\Rightarrow$ $(B\cap\mathbb{Z}^2)\subset (B'\cap\mathbb{Z}^2)$

but

$B\subset B'$ $\not\Leftarrow$ $(B\cap\mathbb{Z}^2)\subset (B'\cap\mathbb{Z}^2)$

$\Rightarrow$ Bruteforce Digital Medial Axis Extraction $O(|X|^2r^2_{max})$ (with $r_{max}$ the maximal DT value)

Implementing IsCoveredBy()§

Goal

Design a IsCoveredBy() predicate with cost as a function of $m$

Elementary Chamfer Masks $\mathcal{M}\in\{d_1, d_\infty\}$

$(p,DT(p)) \in MA \Leftrightarrow DT(p+\vec{v}) < DT(p) + \omega,\, \forall (\vec{v},\omega)\in\mathcal{M}$

Also true for $\mathcal{M}_{3,4}$ with the following rewriting rules of the DT map:

$3 \rightarrow 1$
$6 \rightarrow 5$

Other path-based distances: Look-up table approach

We pre-compute

$Lut(\vec{v},r) =\min \{ r'\,|\, B(O,r)\subseteq B(O+\vec{v},r'\}$
Then,

$(p,DT(p)) \in MA \Leftrightarrow DT(p+\vec{v}) < Lut(\vec{v},DT(p)),\, \forall \vec{v}\in\mathcal{V}$

$\mathcal{V}$ is the neighborhood test.

Bottlenecks Efficient computation of Lut, bounds on $|\mathcal{V}|$ , bounds on $|Lut|$ , …

Global approach using Power Map§

Idea

Get the Medial Axis as a by-product of the Power map

Lemma

Let $S\subset \mathbb{R}^d\times\mathbb{R}$ and $X=\bigcup_{B_i\in S} B_i$

$B\subset B' \implies Power_{S}(B) \cap X = \emptyset$

Non-empty power map cells are related to maximal balls

[Skipping details…]

$\implies$ Separable algorithm to extract the medial axis

$\implies$ $O( d\cdot n^d\cdot (C + H) )$ computational cost for a large class of metrics

One algorithm to rule them all§

Toward Minimal Medial Axis§

Question

Is the set of maximal balls a minimal representation of X as union of balls ?

Answer No

Maximal balls are defined by binary predicates $IsCoveredBy(B,B')$
A maximal ball B can be covered by union of balls, which makes B unnecessary in the representation of X

Toward minimal MA

Thm.

If we allow k-ary predicates $IsCoveredBy(B, S_k)$ with $|S_k|=k$ the minimal medial axis problem becomes NP-complete

Topological Skeleton§

Introduction§

Digital Medial Axis is defined as a set of balls without any topological information

We are thus looking for

a minimal subset of X (skeleton) with the same topological structure
as far as possible, we would like the skeleton to be centered

$\Rightarrow$ Iterative thinning via Simple Point Removal

Simple Point§

Def.

A point $p\in X$ is simple for $X$ if $X$ and $X\setminus\{p\}$ are in the same homotopy equivalence class

In dimension 2, same number of connected components of the object and of its complementary
In dimension 3, we must preserve holes/tunnels
…

Topological Transformations§

From Simple Point Definition

Let $\phi$ by any sequence of insertions/removals of simple points, then $X$ and $\phi(X)$ are in the same homotopy equivalence class

How to characterize simple points ?

$\alpha$ -simple points§

Definition

Def.

A point $p\in X$ is $(\kappa,\lambda)$ -simple for $X$ if

$X$ and $X\setminus\{p\}$ have the same number of $\kappa$ -components
$\overline{X}$ and $\overline{X\setminus\{p\}}$ have the same number of $\lambda$ -components

Example

(which are resp. (0,1)- and (1,0)-simple ?)

Local characterization§

Main Results

Thm.

In dimension 2 and 3, $(\kappa,\lambda)$ -simplicity of $p\in X$ can be decided locally at $p$

( $3\times 3$ neighborhood in 2D, $3\times 3\times 3$ )

E.g. 2D

$N^*_8(p)$ : 8-neighborhood around $p$ (without p)
$C_\kappa^p(X))$ : set of $\kappa$ -connected components in $X$ adjacent to $p$
$T_\kappa(p,X) = | C_\kappa^p(N^*_8(p)\cap X|$
$\Rightarrow$ $p\in X$ is $(\kappa,\lambda)$ -simple for $X$ $\Leftrightarrow$ $T_\kappa(p,X) =T_\lambda(p,\bar{X})=1$

In dimension 3, $T_\kappa(p,X)$ definition is a bit more complex but still local

Illustration§

All configurations in 2D


(0,1)	(1,0)

Homotopic thinning§

Idea

Iterate until stability over sequential simple points removal $\Rightarrow$ ultimate homotopic thinning

Algorithm§

P = { p in X | p is simple for X }
while ( P != empty )
   Q = emptyset
   for all points p in P
     if (p is simple for X)
       X = X \ {p}
       for all q in N(p)
          Q = Q + {q}

   P = emptyset
   for all points p in Q
     if (p is simple for X)
       P = P+ {p}

Homotopic thinning with anchor points§

Idea

Based on an Oracle, we decide to block some simple points during the thinning

Generic algorithm

Breadth first thinning if P is implemented as a queue

P = { p in X | p is simple for X }
while ( P != empty )
   Q = emptyset
   for all points p in P
     if (p is simple for X) and (p is not anchor point)
       X = X \ {p}
       for all q in N(p)
          Q = Q + {q}

   P = emptyset
   for all points p in Q
     if (p is simple for X)
       P = P+ {p}

E.g.

p is anchor point if it has only one neighbor in X

Curve or Surface based Skeleton§

Idea

Anchor points can be specified to generate surface based skeleton

Misc.§

Guided Thinning

Instead of using a queue for P, we consider a priority list with distance transformation values

$\Rightarrow$ Better geometry (central axis) of the skeleton

Parallel thinning

Sequential algorithm needs to test the simplicity twice
Idea of parallel thinning: we mark some simple points and remove them in parallel
- We add constraints on simple points to allow parallel removal
- Orientation based process (N,S,W,E)

$\Rightarrow$ usually, parallel thinning algorithms are more efficient and provide centered skeletons

Active works

Extensions to grayscale images
Optimize the simplicity test (binary decision tree, …)
Definition of simple point sets to allow more efficient parallel removal
Definition of simple structures in topological cellular spaces
…

Digital Geometry: Volumetric Analysis§

Distance Transformation/Voronoi Map§

Problem setting§

Why the DT ?§

Example: curvature from DT§

Application example 1: Snow Micro-structure Analysis§

Application example 2: Metallic Foam Description§

Metric in Digital Space§

Examples§

l_p metrics Illustration§

Discrete Metrics§

Chamfer Mask§

Chamfer Distances§

Simple examples§

Chamfer balls§

Mask Construction§

Distance Transformation algorithm with Chamfer Masks§

Raster Scan Algorithm§

Other path-based distances§

Euclidean metric§

Separable Approach For Squared Euclidean Distance Transform§

First Step§

Second Step§

Lower Envelope Computation§

Overall SEDT Algorithm§

Generalizations§

Generalization to other metrics§

Voronoi Diagram§

Separable Voronoi Map§

Separable Voronoi Map§

Generic Algorithm§

Examples§

Path based Metrics§

Examples§

Subquadratic Algorithm for path based distances§

Path Based Metric and Rational balls§

Distance Evaluation§

Optimized Distance Evaluation§

Separable Predicates for Chamfer Masks§

Dimension 2§

Warm up: Localizing a point§

Algorithm Overview§

Step 1: Shrinking §

Step 2: Shrinking and final computation§

Fast computations in higher dimensions§

Reverse Distance Transformation§

Reverse Transformation§

Separable Approach for §

Illustration§

Associated Structure from Computational Geometry§

Power Map§

Medial Axis Extraction§

Problem Description§

Contact Points based Geometrical Definition§

Maximal Ball based Definition§

DT and Digital Medial Axis§

Digital Ball vs. Euclidean Balls§

Implementing IsCoveredBy()§

Global approach using Power Map§

One algorithm to rule them all§

Toward Minimal Medial Axis§

Heuristics§

Topological Skeleton§

Introduction§

Simple Point§

Topological Transformations§

-simple points§

Local characterization§

Illustration§

Homotopic thinning§

In Dimension 3§

Algorithm§

Homotopic thinning with anchor points§

Illustration§

Curve or Surface based Skeleton§

Misc.§

Step 1: Shrinking $\mathcal{M}_p$ §

Step 2: Shrinking $\mathcal{M}_q$ and final computation§

Separable Approach for $l_2$ §

$\alpha$ -simple points§