Given a set X, a topology T on X is a collection of subsets of X satisfying (Ti) X, 0 ∈ T; (Tii) for U _α  ∈ T; and (Tiii) for U _i  ∈ T. X together with T is called a topological space. By convention, sets in T are referred to as open sets and their complements as closed sets.

A neighborhood N(x) of x is a set containing an open set U ∈ T such that x ∈ U. Points with neighborhood are called interior.

Two topological spaces X and Y are homeomorphic if there exists bijection which is continuous and has continuous inverse α ⁻¹. Since homeomorphisms copy topologies, homeomorphic spaces are topologically equivalent [1].

A function which is (Mi) positive-definite (d(x, y) > 0 for all x ≠ y and d(x, y) = 0 for x = y) and (Mii) subadditive (d(x, z) ≤ d(x, y) + d(y, z) for all x, y, z) is called a metric on X. The metric is an abstraction of the notion of distance between pairs of points on X. Property (Mii) is called triangle inequality and generalizes the known fact: the sum of the lengths of two edges of a triangle is greater or equal to the length of the third edge. The pair (X, d) is called a metric space.

A metric induces topology through the definition of open metric ball B _r (x) = { x ^′  ∈ X: d(x, x ^′ ) < r}. A neighborhood of x in a metric space is a set containing a metric ball B _r (x) [34].

Given two metric spaces (X, d) and (Y, δ), the set of pairs such that for every x ∈ X there exists at least one y ∈ Y such that (x, y) ∈ C, and similarly, for every y ∈ Y, there exists an x ∈ X such that (x, y) ∈ C is called a correspondence between X and Y. Note that a correspondence C is not necessarily a function. The correspondence is called bijective if every point in X has a unique corresponding point in Y and vice versa.

The discrepancy of the metrics d and δ between the corresponding points is called the distortion of the correspondence,

Metric spaces (X, d) and (Y, δ) are said to be ε–isometric if there exists a correspondence C with dis(C) ≤ ε. Such a C is called an ε–isometry.

A particular case of a 0-isometry is called an isometry. In this case, the correspondence is a bijection and X and Y are called isometric.

A Hausdorff space X which is locally homeomorphic to (i.e., for every x in X, there exists a neighborhood U and a homeomorphism ) is called an n–manifold or an n–dimensional manifold. The function α is called a chart. A collection of neighborhoods that cover X together with their charts is called an atlas on X. Given two charts α and β with overlapping domains U and V, the map is called a transition function. A manifold whose transition functions are all differentiable is called a differentiable manifold. More generally a –manifold has all transition maps k-times continuously differentiable. A -manifold is called smooth.

A manifold with boundary is not a manifold in the strict sense of the above definition. Its interior points are locally homeomorphic to , and every point on the boundary ∂ X is homeomorphic to .

Of particular interest for the discussion in this chapter are two-dimensional (n = 2) manifolds, which model boundaries of physical objects in the world surrounding us. Such manifolds are also called surfaces. In the following, when referring to shapes and objects, the terms manifold, surface, and shape will be used synonymously.

Locally, a manifold can be represented as a linear space, in the following way. Let be a chart on a neighborhood of x and be a differentiable curve passing through x = γ(0). The derivative of the curve is called a tangent vector at x. The set of all equivalence classes of tangent vectors at x forming an n-dimensional real vector space is called the tangent space T _x X at x.

A family of inner products depending smoothly on x is called Riemannian metric tensor. A manifold X with a Riemannian metric tensor is called a Riemannian manifold.

The Riemannian metric allows to define local length structures and differential calculus on the manifold. Given a differentiable scalar-valued function , the exterior derivative (differential) is a form on the tangent space TX. For a tangent vector , . ∇f is called the gradient of f at x and is a natural generalization of the notion of the gradient in vector spaces to manifolds. Similarly to the definition of Laplacian satisfying
for differentiable scalar-valued functions f and h, the operator is called the Laplace–Beltrami operator, a generalization of the Laplacian. Here, μ denotes the measure associated with the n-dimensional volume element (area element for n = 2). The Laplace–Beltrami operator is (Li) symmetric (), (Lii) of local action ( is independent of the value of f(x ^′ ) for x ^′ ≠ x), and (Liii) positive semi-definite () (in many references, the Laplace–Beltrami is defined as a negative semi-definite operator) and (Liv) coincides with the Laplacian on Euclidean domains, such that if f is a linear function and X is Euclidean.

Another important use of the Riemannian metric tensor is to measure the length of paths on the manifold. Given a continuously differentiable curve , its length is given by
For the set of all continuously differentiable curves between the points x, x ^′ ,
defines a metric on X referred to as length or geodesic metric. If the manifold is compact, for any pair of points x and x ^′ , there exists a curve called a minimum geodesic such that .

A particular realization of a Riemannian manifold called embedded manifold (or embedded surface for n = 2) is a smooth submanifold of (m > n). In this case, the tangent space is an n-dimensional hyperplane in , and the Riemannian metric is defined as the restriction of the Euclidean inner product to the tangent space, .

The length of a curve on an embedded manifold is expressed through the Euclidean metric,
and the geodesic metric d _X defined according to (1) is said to be induced by . (Repeating the process, one obtains that the metric induced by d _X is equal to d _X . For this reason, d _X is referred to as intrinsic metric [34].)

Though apparently embedded manifolds are a particular case of a more general notion of Riemannian manifolds, it appears that any Riemannian manifold can be realized as an embedded manifold. This is a consequence of the Nash embedding theorem [84], showing that a Riemannian manifold can be isometrically embedded in a Euclidean space of dimension m = n ² + 5n + 3. In other words, any smooth Riemannian manifold can be defined as a metric space which is isometric to a smooth submanifold of a Euclidean space with the induced metric.

Riemannian manifolds do not have a unique realization as embedded manifolds. One obvious degree of freedom is the set of all Euclidean isometries: two congruent embedded manifolds are isometric and thus are realizations of the same Riemannian manifold. However, a Riemannian manifold may have two realizations which are isometric but incongruent. Such manifolds are called nonrigid. If, on the other hand, a manifold’s only isometries are congruences, it is called rigid.

A function is called a diffusion kernel if it satisfies the following properties: (Ki) nonnegativity: k(x, x) ≥ 0; (Kii) symmetry: k(x, y) = k(y, x); (Kiii) positive-semidefiniteness: for every bounded f,
(Kiv) square integrability: ; and (Kv) conservation: (y) = 1. The value of k(x, y) can be interpreted as a transition probability from x to y by one step of a random walk on X.

Diffusion kernel defines a linear operator
which is known to be self-adjoint. Because of (Kiv), has a finite Hilbert norm and therefore is compact. As the result, it admits a discrete eigendecomposition with some eigenfunctions and eigenvalues . α _i  ≥ 0 by virtue of property (Kiii), and α _i  ≤ 1 by virtue of (Kv) and consequence of the Perron–Frobenis theorem.

By the spectral theorem, the diffusion kernel can be presented as . Since form an orthonormal basis of L ²(X),
a fact sometimes referred to as Parseval’s theorem. Using these results, properties (Kiii–Kv) can be rewritten in the spectral form as 0 ≤ α _i  ≤ 1 and .

An important property of diffusion operators is the fact that for every t ≥ 0, the operator is also a diffusion operator with the eigenbasis of and corresponding eigenvalues . The kernel of expresses the transition probability by random walk of t steps. This allows to define a scale space of kernels, , with the scale parameter t.

There exist a large variety of possibilities to define a diffusion kernel and the related diffusion operator. Here, we restrict our attention to operators describing heat diffusion. Heat diffusion on surfaces is governed by the heat equation
where u(x, t) is the distribution of heat on the surface at point x in time t, u ₀ is the initial heat distribution, and is the positive-semidefinite Laplace–Beltrami operator, a generalization of the second-order Laplacian differential operator to non-Euclidean domains. (If X has a boundary, boundary conditions should be added.)

On Euclidean domains (), the classical approach to the solution of the heat equation is by representing the solution as a product of temporal and spatial components. The spatial component is expressed in the Fourier domain, based on the observation that the Fourier basis is the eigenbasis of the Laplacian , and the corresponding eigenvalues are the frequencies of the Fourier harmonics. A particular solution for a point initial heat distribution u ₀(x) = δ(x − y) is called the heat kernel , which is shift invariant in the Euclidean case. A general solution for any initial condition u ₀ is given by convolution , where is referred to as heat operator.

In the non-Euclidean case, the eigenfunctions of the Laplace–Beltrami operator can be regarded as a “Fourier basis” and the eigenvalues given the “frequency” interpretation. The heat kernel is not shift invariant but can be expressed as .

It can be shown that the heat operator is related to the Laplace–Beltrami operator as , and as a result, it has the same eigenfunctions ϕ _i and corresponding eigenvalues . It can be thus seen as a particular instance of a more general family of diffusion operators diagonalized by the eigenbasis of the Laplace–Beltrami operator, namely, ’s as defined in the previous section but restricted to have the eigenfunctions ψ _i  = ϕ _i . The corresponding diffusion kernels can be expressed as
where K(λ) is some function such that α _i  = K(λ _i ) (in the case of , K(λ) = e ^{−t λ} ). Since the Laplace–Beltrami eigenvalues can be interpreted as frequency, K(λ) can be thought of as the transfer function of a low-pass filter. Using this signal processing analogy, the kernel k(x, y) can be interpreted as the point spread function at a point y, and the action of the diffusion operator on a function f on X can be thought of as the application of the point spread function by means of a non shift-invariant version of convolution. The transfer function of the diffusion operator is K ^t (λ), which can be interpreted as multiple applications of the filter K(λ). Such multiple applications decrease the effective bandwidth of the filter and, consequently, increase its effective support in space. Because of this duality, both k(x, y) and K(λ) are often referred to as diffusion kernels.

Farthest point sampling (FPS) is a greedy procedure constructing a sequence of samplings . A sampling is constructed from by adding the farthest point
The sequence {r _N } of the sampling radii associated with is nonincreasing, and, furthermore, each is also r _N -separated. The starting point x ₁ is usually picked up at random, and the stopping condition can be either the sampling size or radius.

Though FPS does not strictly minimize any of the error criteria defined in the previous section, in terms of ε _∞ , it is no more than twice inferior to the optimal sampling of the same size [59]. In other words, for produced using FPS,
This result is remarkable, as finding the optimal sampling is known to be an NP-hard problem.

In order for a sampling to be ε ₂-optimal, each sample x _i has to minimize
A point minimizing the latter quantity is referred to as the centroid of V _i . Therefore, an ε ₂-optimal sampling induces a so-called centroidal Voronoi tessellation (CVT), in which the centroid of each Voronoi region coincides with the sample point inducing it [45, 90]. Such a tessellation and the corresponding centroidal Voronoi sampling are generally not unique.

A numerical procedure for the computation of a CVT of a shape is known as the Lloyd–Max algorithm [68, 73]. Given some initial sampling of size N (produced, e.g., using FPS), the Voronoi tessellation induced by it is computed. The centroids of each Voronoi region are computed, yielding a new sampling of size N. The process is repeated iteratively until the change of becomes insignificant. While producing high-quality samplings in practice, the Lloyd–Max procedure is guaranteed to converge only to a local minimum of ε ₂. For computational aspects of CVTs on meshes, the reader is referred to [90].

The simplest representation of a shape is obtained by considering the points of the sampling as points in the ambient Euclidean space. Such a representation is usually referred to as a point cloud. Points in the cloud are called vertices and denoted by . The notion of a point cloud can be generalized using the formalism of simplical complexes. For our purpose, an abstract k–simplex is a set of cardinality k + 1. A subset of a simplex is called a face. A set K of simplices is said to be an abstract simplical complex if any face of σ ∈ K is also in K, and the intersection of any two simplices σ ₁, σ ₂ ∈ K is a face of both σ ₁ and σ ₂. A simplicial k-complex is a simplicial complex in which the largest dimension of any simplex is k. A simplicial k-complex is said to be homogeneous if every simplex of dimension less than k is the face of some k-simplex. A topological realization of a simplicial complex K maps K to a simplicial complex in , in which vertices are identified with the canonical basis of , and each simplex in K is represented as the convex hull of the corresponding points . A geometric realization is a map of the simplicial complex to defined by associating the standard basis vectors with the vertex positions .

In this terminology, a point cloud is a simplicial 0-complex having a discrete topology. Introducing the notion of neighborhood, we can define a subset of pairs of vertices that are adjacent. Pairs of adjacent vertices are called edges, and the simplicial 1-complex has a graph topology, i.e., the set of vertices forms an undirected graph with the set of edges E. A simplicial 2-complex consisting of vertices, edges, and triangular faces built upon triples of vertices and edges is called a triangular mesh. The mesh is called topologically valid if it is homeomorphic to the underlying continuous surface X. This usually implies that the mesh has to be a two manifold. A mesh is called geometrically valid if it does not contain self-intersecting triangles, which happens if and only if the geometric realization of the mesh is bijective. Consequently, any point x on a geometrically valid mesh can be uniquely represented as . The vector is called the barycentric coordinates of x and has at most three nonzero elements. If the point coincides with a vertex, is a canonical basis vector; if the point lies on an edge, has two nonzero elements; otherwise, has three nonzero elements and x lies on a triangular face.

A particular way of constructing a triangular mesh stems from the Voronoi tessellation induced by the sampling. We define the simplicial 3-complex as
in which a pair of vertices spans an edge and a triple of vertices spans a face if the corresponding Voronoi regions are adjacent. A mesh defined in this way is called a Delaunay mesh. (Unlike in the Euclidean case where every sampling induces a valid Delaunay triangulation, an invalid Voronoi tessellation results in a topologically invalid Delaunay mesh. In [66], Leibon and Letscher showed that under the same conditions sufficient for the existence of a valid Voronoi tessellation, the Delaunay mesh is also topologically valid.)

Shapes homeomorphic to a disk can be parametrized using a single global chart, e.g., on the unit square, . (Manifolds with more complex topology can still be parametrized in this way by introducing cuts that open the shape into a topological disk.) Such surfaces are called parametric and can be sampled directly in the parametrization domain. For example, if the parametrization domain is sampled on a regular Cartesian grid, the shape can be represented as three N × N arrays of x, y, and z values. Such a completely regular structure is called a geometry image [56, 69, 99] and can be thought indeed as a three-channel image that can undergo standard image processing such as compression. Geometry images are ideally suitable for processing by vector and parallel hardware.

Another way of representing a shape is by considering the isosurfaces of some function defined on a region of . Such a representation is called implicit, and it often arises in medical imaging applications, where shapes are two-dimensional boundaries created by discontinuities in volumetric data. Implicit representation can be naturally processed using level-set-based algorithms, and it easily handles arbitrary topology. A disadvantage is the bigger amount of storage commonly required for such representations.

A famous algorithm for the solution of the recursion (21) was proposed by Dijkstra. Dijkstra’s algorithm measures the distance map from the source vertex x _i to all the vertices in the graph.

Every vertex in Dijkstra’s algorithm is processed exactly once; hence, Nn outer iterations are performed. Extraction of vertex with smallest d is straightforward with complexity and can be reduced to using efficient data structures such as Fibonacci heap. In the inner loop, updating adjacent vertices in our case, since the graph is sparsely connected, is . The resulting overall complexity is .

Unfortunately, the graph distance d _L is an inconsistent approximation of d _X , in the sense that d _L usually does not converge to d _X when the sampling becomes infinitely dense. This phenomenon is called metrication error, and the reason is that the graph induces a metric inconsistent with d _X (Fig. 3). While metrication errors make in general the use of d _L an approximation of d _X disadvantageous, the Bernstein–de Silva–Langford–Tenenbaum theorem [14] states that under certain conditions, the graph metric d _L can be made as close as desired to the geodesic metric d _X . The theorem is formulated as a bound of the form
where λ ₁, λ ₂ depend on shape properties, sampling quality, and graph connectivity. In order for d _L to represent d _X accurately, the sampling must be sufficiently dense, length of edges in the graph bounded, and sufficiently close vertices must be connected, usually in a non-regular manner.

An alternative to computation of a discrete metric on a discretized surface is the discretization of the metric itself. The distance map d(x) = d _X (x ₀, x) (Fig. 4) on the manifold can be associated with the time of arrival of a propagating front traveling with unit speed (illustratively, imagine a fire starting at point x ₀ at time t = 0 and propagating from the source). Such a propagation obeys the Fermat principle of least action (the propagating front chooses the quickest path to travel, which coincides with the definition of the geodesic distance) and is governed by the eikonal equation
where ∇ _X is the intrinsic gradient on the surface X. Eikonal equation is a hyperbolic PDE with boundary conditions d(x ₀) = 0; minimum geodesics are its characteristics. Propagation direction is the direction of the steepest increase of d and is perpendicular to geodesics.

Since the distance map is not everywhere differentiable (in particular, at the source point), no solution to the eikonal equation exists in the classical sense, while there exist many non functions satisfying the equation and the boundary conditions. Among such functions, the largest d satisfying the boundary conditions and the inequality
at every point where ∇ _X d exists is called the viscosity solution [42]. The viscosity solution of the eikonal equation always exists and is unique, and its value at a point x coincides with d _X (x, x ₀). It is known to be monotonous, i.e., not having local maxima.

A family of algorithms for finding the viscosity solution of the discretized eikonal equation by simulated wavefront propagation is called fast marching methods [63, 100, 112]. Fast marching algorithms can be thought of as continuous variants of Dijkstra’s algorithm, with the notable difference that they consistently approximate the geodesic metric d _X on the surface.

The general structure of fast marching closely resembles that of Dijkstra’s algorithm with the main difference lying in the update step. Unlike the graph case where shortest paths are restricted to pass through the graph edges, the continuous approximation allows paths passing anywhere in the simplicial complex. For that reason, the value of d(x) has to be computed from the values of the distance map at two other vertices forming a triangle with x. In order to guarantee consistency of the solution, all such triangles must have an acute angle at x. Obtuse triangles are split at a preprocessing stage by adding virtual connections to nonadjacent vertices.

Given a triangle (x, x ₁, x ₂) with known values of d(x ₁) and d(x ₂), the goal of the update step is to compute d(x). The majority of fast marching algorithms do so by simulating the propagation of a planar wavefront in the triangle. The wavefront arrival time to x ₁ and x ₂ is set to d(x ₁) and d(x ₂), from which the parameters of the wave source are estimated. Generally, there exist two solutions for d(x) consistent with the input, the smallest corresponding to the wavefront first arriving to x and then to x ₁ and x ₂

Stay updated, free articles. Join our Telegram channel

Join