is a planar triangulation. Tripple point: si,j,k Ci ⇥ Cj ⇥ Ck. Double point: wi,j = argmin x Ci ⇥Cj d(x, xi ). Geometric realization of DS : xj xk xi si,j,k wi,j Union of two geodesics starting from wi,j . Special case: boundary curves if xi and xj . Fast computation in O(n log(n)) with Fast Marching.
, xn } M. Parametric domain sampling: : [0, 1]2 ⇥ M, xi,j = (i/ ⇤ n, j/ ⇤ n). -covering: i B (xi ) = M, where B (x) def. = {y \ dM (x, y) }. Theorem: The sampling {x1, . . . , xn } is an -covering that is -separated for = max i=1,...,n min j=1,...,n dM (xi, xj ). 1. Initialization: x1 random, d(x) dM (x0, x), set i = 1. 2. Select point: xi+1 = argmax x d(x), = d(xi+1 ). 3. Local update of the distance: d(x) min(d(x), dM (xi+1, x)). 4. Stop: If i < n or > 0 , set i i + 1 and go back to 2. Farthest point -separated: min(dM (xi, xj )) .
of functions. 8 Jonathan Richard Shewchuk Figure 2: A visual illustration of how large angles, but not small angles, can cause the error ∇f − ∇g to explode. In each triangulation, 200 triangles are used to render a paraboloid. 35 50 65 40 20 20 40 40 40 20 Figure 3: As the large angle of the triangle approaches 180◦, or the sliver tetrahedron becomes arbitrarily ﬂat, the magnitude of the vertical component of ∇g becomes arbitrarily large. sphere, then perturbing one of the vertices just off the equator so that the sliver has some (but not much) volume. Because of this sensitivity, mesh generators usually choose the shapes of elements to control ∇f − ∇g ∞, and not f − g ∞, which can be reduced simply by using smaller elements. Section 6.1 presents quality measures that judge the shape of elements based on their ﬁtness for interpolation. Table 2 gives two upper bounds on ∇f − ∇g ∞ over a triangle. The ﬁrst upper bound is almost tight, to within a factor of two. The “weaker but simpler upper bound” of 3ctrcirc is not as good an indicator as the stronger upper bound, but it has the advantages of being smooth almost everywhere (and therefore more Isotropic Anisotropic, Anisotropic, bad shape good shape f(x + h) = f(x) + ⇤xf, h⇥ + Hx (f)h, h⇥ + O(||h||2) f(x + h) (f(x) + ⇥⌅xf, h⇤) = µ1 ⇥h, e1 ⇤2 + µ2 ⇥h, e2 ⇤2 + O(h2). = Local optimal shape: Hx = µ1e1e1 T + µ2e2e2 T with |µ1 | > |µ2 |. width/length = |µ2 |/|µ1 | x e1 ⇥ µ2 ⇥ µ1
image f, that is C2 outside C2 edges, one can however built an adapted triangulation that enhances signiﬁcantly the approximation performance of a wavelet basis. The following theorem sketch the construction of such a triangulation. f ∗ h f Figure 9.2: Approximation with ﬁnite elements on a triangulation for a function without and with additional blurring. Theorem 2. If f is C2 outside a set of C2 contours, then it exists a constant C such that for any M, it exists an adapted triangulation (V, T) over which the piecewise-linear interpolation fM satisﬁes ||f fM ||2 C M 2. (9.2) Proof. We give here only a sketch of the proof. [todo : show an image of the tubes] Near an edges, the lengths of the triangles should be of order M 1 and their widths should be of order M 2. We deﬁne a thin band of width M 2 around all the edge curves, see ﬁgure ??. Since the edge curve are C2, this band can be sub-divided in elongated tubes of length proportional to M 1, each of witch is decomposed in two elongated triangles. The number of such triangle is 2LM where L is the total length of the edges in the image. A special care should be taken near edge corner or edge crossings but we ignore these technicalities here. Over these tubes, the function is bounded and thus the approximation error ||f e f||2 L2( ) in the tubes is of the order of area( )||f||2 ⇥ = LM 2||f||2 ⇥ . In the complementary of the tubes c, one packs M large equilateral triangles of edge lengths approximatively ⌃ ⇥ M 1/2. Since f is C2 outside the edge curve, the approximation error of a linear interpolation e f over such The piecewise linear interpolation ⇥ fM of f on (V, T) is such that ⇥ fM (xi ) = f(xi ) and ⇥ fM is linear on each triangular face tj . This is a two-dimensional extension of the spline approximations studied in Section ??. This piecewise linear approximation depends on the position of the vertices V and the connectivity T of the triangulation. In order to e⇥ciently approximate a given function f with M triangles, one needs to ﬁnd the optimal shape of the triangles. For a geometrically regular image f, that is C2 outside C2 edges, one can however built an adapted triangulation that enhances signiﬁcantly the approximation performance of a wavelet basis. The following theorem sketch the construction of such a triangulation. f ∗ h f Figure 9.2: Approximation with ﬁnite elements on a triangulation for a function without and with additional blurring. Theorem 2. If f is C2 outside a set of C2 contours, then it exists a constant C such that for any M, it exists an adapted triangulation (V, T) over which the piecewise-linear interpolation ⇥ fM satisﬁes ||f ⇥ fM ||2 C M 2. (9.2) Proof. We give here only a sketch of the proof. [todo : show an image of the tubes] Near an edges, 9.4. Geometric Image Approximation As already seen in the proof of theorem ??, one can design triangles of width and legnth M the tube so that the error in this area is also of the order of M 2, which concludes the p The proof of this theorem shows that when an image is smoothed by an unknow width ⇤, the triangulation should depend on ⇤ in order to get a fast decay of the app error. To reach an error decay of O(M 2), in the neighborhood of a contour smoothed of width ⇤, the triangle should have a length of order ⇤1/4M 1/2 and a width of order ⇤ as shown on ﬁgure ??. The scale ⇤ is most of the time unknown and one thus needs an algorithm to devise the size of the triangles. s1/ 4M1/ 2 triangles s1/ 4M−1/ 2 s3/ 4M−1/ 2 Figure 9.4: Aspect ratio of triangles for the approximation of a blurred contou This shows that locally, an optimal triangle should be aligned with the direction µ function is the most regular. The local behavior of a smooth C2 image is well descr second order by a quadratic approximation f(x + h) = f(x) + ⌅⌥xf, h⇧ + ⌅Hx (f)h, h⇧ + O(||h||2) 2⇥2 ||f ⇥ fM ||2 C M 2. (9.2) . We give here only a sketch of the proof. [todo : show an image of the tubes] Near an edges, ngths of the triangles should be of order M 1 and their widths should be of order M 2. We e a thin band of width M 2 around all the edge curves, see ﬁgure ??. Since the edge curve are his band can be sub-divided in elongated tubes of length proportional to M 1, each of witch is mposed in two elongated triangles. The number of such triangle is 2LM where L is the total length e edges in the image. A special care should be taken near edge corner or edge crossings but we e these technicalities here. Over these tubes, the function is bounded and thus the approximation ||f e f||2 L2( ) in the tubes is of the order of area( )||f||2 ⇤ = LM 2||f||2 ⇤ . In the complementary e tubes c, one packs M large equilateral triangles of edge lengths approximatively ⇤ M 1/2. f is C2 outside the edge curve, the approximation error of a linear interpolation e f over such gles is ||f e f||⇤( c ||f||C 2 and thus the error satisﬁes ||f e f||2 = O(LM 2||f||2 ⇤ + ||f||2 C M 2). M−1 M−2 M−1 M−α Figure 9.3: Finite elements for the approximation around a singularity curve. theorem shows that an adapted triangulation should balance the approximation error d outside the singularity curves. The triangles located far away from these singularities e large and isotropic, while triangles that cover the edges should be stretched along the y curves. This construction can be generalized by replacing triangles by higher order c primitives whose boundaries are polynomial curves of degree , as shown on the right gure ??. The adapted approximation using polynomials deﬁned on M such higher order Sharp edges Blurred edges Sharp edges Blurred edges
Figure 16: Mandrill original, a 200×200 raster image. Figure 17: Mandrill approximated with Gouraud shaded triangles created by subsampling on a uniform 20×20 grid (400 vertices). Figure 16: Mandrill original, a 200×200 raster image. Figure 17: Mandrill approximated with Gouraud shaded triangles created by subsampling on a uniform 20×20 grid (400 vertices). Figure 18: Mandrill approximated with Gouraud Figure 19: Mesh for the image to the left. Original image Isotropic mesh (400 vert) Anisotropic mesh (400 vert)
A. Iske (a) (b) (c) (d) Fig. 8. Reflex. (a) Original image of size 128 × 128. Compression at 0.251 bpp and reconstruction by (b) SPIHT with PSNR 30.42 db, (d) AT∗ 2 with PSNR 41.73 db. 16 L. Demaret, N. Dyn, M.S. Floater, A. Iske (a) (b) (c) (d) Fig. 8. Reflex. (a) Original image of size 128 × 128. Compression at 0.251 bpp and Original image Wavelets @.25bpp Triangulation @.25bpp
computational geometry extend 105 to the setting of a Riemannian metric. 106 Starting from a set of points S = {xi }m i=1 , one can deﬁne graphs and trian- 107 gulations that reﬂect the geometry of the Riemannian manifold. These points 108 and the corresponding graphs are the basic building blocks of the algorithms for 109 perceptual grouping and planar domain meshing. 110 In the following, the boundary ⇧ of is assumed to be a set of closed smooth curves. At least one point of S is located on each curve, and these boundary points segment ⇧ as a set of sub-curves ⇧ = ⇥ i,j ⇥i,j with ⇥i,j ⇤ P(xi, xj ) ⇥i,j ↵ S = {xi, xj }. (one can have xi = xj if there is only one point on a curve). Ω Ωc x1 x2 x3 x4 111 2.1 Delaunay and Voronoi Graphs 112 The segmentation of the domain in Riemannian Voronoi cells is 113 = C0 ⇥ xi S Ci where Ci = {x ⇤ \ ⇧ j ⌅= i, d(xi, x) d(xj, x)} . (7) The outer Voronoi cell is deﬁned as C0 = Closure( c). 114 The Delaunay graph DS of S is a graph where two points are connected if their respective Voronoi cells are adjacent (xi, xj ) ⇤ DS ⇥ Ci ↵ Cj ⌅= ⌃. To each Delaunay edge (xi, xj ) ⇤ DS corresponds a double point wi,j , which is the closest point to xi and xj on the common Voronoi cell boundary wi,j = argmin x Ci ⌅Cj d(x, xi ). ECCV-08 submission ID1057 11 e to the boundary. A boundary sub-curve i,j is 245 S if it exists a triple point wi,k,0 i,j , see ﬁgure 246 nnot be part of the Delaunay triangulation, and 247 orithm by inserting a mid point. Similarly, triple 248 croach any boundary sub-curve (the sub-curve is 249 250 he Delaunay graph DS of S is not necessarily a 251 ng S is not dense enough, see [7]. This is because 252 connected to only one other point of S in DS . 253 dd points on the Voronoi cell boundary of such a 254 255 x2 x1 x2 x3 aches the boundary curve ⇥1,2. Right: the vertex x3 e (x1, x2) is a Delaunay edge. ECCV if a point xk is located too close to the boundary. 245 said to be encroached by xk S if it exists a triple 246 7. Such an encroached edge cannot be part of the 247 is automatically split by the algorithm by inserting 248 points are not added if they encroach any boundar 249 subdivided instead). 250 Another di⇥culty is that the Delaunay graph 251 valid triangulation if the sampling S is not dense en 252 of some isolated point, that is connected to only 253 The algorithm automatically add points on the Vo 254 point. 255 x1 x2 x3 x1 Fig. 7. Left: the vertex x3 encroaches the boundary does not encroach anymore because (x1, x2) is a Delau 2D Riemannian manifold: [0, 1]2 equipped with H(x) ⇥ R2 2. Boundaries: ⇥ = ⇥ i,j i,j with i,j P(xi, xj ) i,j ⇥ S = {xi, xj }. Requirement: DS . Good situation Bad: x3 encroaches 1,2
seeding strategy of [?] w 242 metrics and domains with arbitrary boundaries. 243 Our anisotropic meshing algorithm proceeds by iteratively i point si,j,k ⇥ S to an already computed set of points S. In or an anisotropic mesh with triangles of high quality with resp metric, one inserts si,j,k for a Delaunay triangle (xi, xj, xk ) w circumradius to shortest edge ratio ⇥(si,j,k ) = d(si,j,k, xi ) min(d(xi, xj ), d(xj, xk ), d(xk, xi )), which is a quantity computed for each triple point in parallel to 244 ing propagation. In the Euclidean domain, a triangle (xi, xj, xk ) 245 of ⇥(si,j,k ) is badly shaped since its smallest angle is close to 0 . 246 [?], this property extends to an anisotropic metric H(x) if angl 247 using the inner product deﬁned by H(x). 248 = circumradius shortest edge Fig. 7. Left: the vertex x3 encroaches the boundary curve ⇥1,2. Right: the vertex x3 does not encroach anymore because (x1, x2) is a Delaunay edge. Table 2 details this algorithm. A bound ⇥⇥ on ⇥ enforces the reﬁnement to reach some quality criterion, while a bound U⇥ enforces a uniform reﬁnement to match some desired triangle density. 1. Initialization: set S with at least one point on each curve of , compute US with a Fast Marching. 2. Boundary enforcement: while it exists ⇥i,j ⌃ encroached by some xk ⇤ S, subdivide: S ⇥ S ⌃ argmax x⇤⇤i,j US (x). Update US with a local Fast Marching. 3. Triangulation enforcement: while it exists (xi, xj) ⇤ DS with xi or xj isolated, insert w⇥ = argmax w⇤Ci⇧Cj d(xi, w). 4. Select point: s⇧ ⇥ argmin s⇤ S ⌃⇥ ⇤(s). – If in S ⌃ {s⇧}, s⇧ encroaches some ⇥i,j ⌃ , subdivide: S ⇥ S ⌃ argmax x⇤⇤i,j US (x). – Otherwise, add it: S ⇥ S ⌃ s⇧. Update US with a local Fast Marching. 5. Stop: while ⇤(s⇧) > ⇤⇧ or US (s⇧) > U⇧, go back to 2. x1 x3 s1,2,3 s1,2,3 x3 x2 x2 x1
metric tensor field generate anisotropic mesh. nisotropic mesh. If metric tensor derivatives, no triangle has angle < 20° as measured by any point in the triangle. Main Result M is smooth with bounded etric dual of an anisotropic Voronoi di- angulation. We describe conditions in guaranteed to be entirely visible from ction 5. For the special case of two di- ons also guarantee that the planar dual with no inverted triangles. ible an algorithm that generates high- isotropic meshes by reﬁ ning an aniso- nforce the conditions that guarantee that nd to remove any poor-quality elements
g(A) = argmin x A ||x y||2dy. A M, geodesic center of gravity: gM (A) = argmin x A dM (x, y)2dy Sampling {xi }m i=1 , Voronoi tessellation Voronoi M ({xi }i ) = {Vi }m i=1 . Centroidal tessellation i, xi = gM (Vi ), local minimizer of E({xi }, {Vi }i ) = m i=1 ⇥ Vi dM (xi, y)2dy. 1. Initialize: {xi }i random. 2. Update regions: ⇥ i, {Vi }i Voronoi M ({xi }i ). 3. Update centers: ⇥ i, xi gM (Vi ). 4. Stop: while not converged, go back to 2. Lloyd-Max x y ⇥ nx (y) EA (x) = A dM (x, y)2dy =⇤ ⌅x EA = A dM (x, y) ⇥ nx (y)dy