Computational Geometry csci3250 Laura Toma Bowdoin - - PowerPoint PPT Presentation

Computational Geometry csci3250

• Laura Toma
• Bowdoin College

Today

• Triangulations
• Delaunay triangulation
• Properties
• How to construct it
• Applications

Triangulation

Input: P = {p1, p2,…,pn} set of points in the plane

• A triangulation T(P) is a maximal planar subdivision whose vertices are P.
• maximal planar subdivision: no edge can be added without distorting planarity

planar : edges can only intersect at endpoints

• Every bounded face must be a triangle.

The unbounded face : bounded by convex hull.

• Many triangulations possible

P = {p1, p2,…,pn} set of points in the plane

• In practice, want triangulations without small angles
• small angles cause numerical instability with geometric predicates, etc

Size of a triangulation

• Let k = number of points in P that lie on CH(P).
• Any triangulation of P has 2n - 2 -k triangles and 3n - 3 - k edges.
• Proof: Use Euler’s formula: e - v + f = 2;

P = {p1, p2,…,pn} set of points in the plane

Voronoi Diagram Vor(P)

• Vor(pi): all points in the plane that are closer to pi than to any other site

P = {p1, p2,…,pn} set of points in the plane

v p1 p2 p3

Voronoi Diagram Vor(P)

• The edges of Vor(P) are (segments of) perpendicular bisectors

P = {p1, p2,…,pn} set of points in the plane

v p1 p2 p3

Voronoi Diagram Vor(P)

• Every Voronoi vertex is the intersection of 3 Voronoi regions

(if no 4-points co-circular)

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

v p1 p2 p3

Voronoi Diagram Vor(P)

• Degeneracies
• more than 3 points lie on same circle
• collinear points

P = {p1, p2,…,pn} set of points in the plane

Voronoi Diagram Vor(P)

• Every Voronoi vertex is the intersection of 3 Voronoi regions
• v is equidistant from p1, p2, p3

v p1 p2 p3

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

Voronoi Diagram Vor(P)

• Every Voronoi vertex is the intersection of 3 Voronoi regions
• v is equidistant from p1, p2, p3
• C(v): circle of p1,p2,p3 has center v

v p1 p2 p3

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

Voronoi Diagram Vor(P)

• Every Voronoi vertex is the intersection of 3 Voronoi regions
• v is equidistant from p1, p2, p3
• C(v): circle of p1,p2,p3 has center at v
• No other site can be inside C(v) <———— empty circle property

v p1 p2 p3

• Every Voronoi vertex is the

center of a circle that has 3 sites on its boundary and no other sites inside

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

The dual of Vor(P)

• For every pair of Voronoi regions Vor(u) and Vor(v) that share and edge, draw

an edge between u and v ==> Each Voronoi vertex v defines a triangle v p1 p2 p3

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

Delaunay triangulation DT

Theorem [Delaunay 1934]: The straight-line dual of Vor(P) is planar, and is a triangulation. Proof: not trivial (see p.197 in 4M book) v p1 p2 p3

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

Delaunay triangulation DT

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

Delaunay triangulation DT

By definition, DT(P) is the dual of Vor(P). Would be nice to have a direct way to characterize DT(P) (without Vor(P)).

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

v p1 p2 p3 p1 p2 p3

Delaunay triangulation DT

What we know from Vor(P)

• p1p2p3 is a triangle in DT(P) <==> circumcircle of p1p2p3 is empty
• Theorem:

A triangulation T(P) is the DT(P) iff all triangles of T(P) have the empty-circle property. Proof: “==>” easy “<==“ requires a bit more work

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

Delaunay triangulation DT

• How to build a triangulation such that all triangles have empty-circle property?

p1 p2 p3

circumcircle of p1p2p3 is empty

• How to build a triangulation?
• Input: a set of points
• Output: a partition of CH(P) into triangles
• P = {p1, p2,…,pn} set of points in the plane

no 4 points co-circular

• Many triangulations possible

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

• Come up with an algorithm to compute a triangulation of P
• Possible ideas: incremental, plane sweep, divide-and-conquer, …

Class work:

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

Flipping edges

edge flips

How to get a DT?

• Compute an arbitrary triangulation T(P)
• T ————> DT

edge flips

a c b d a c b d a c b d

edge flips

Flipping edges and empty circles

C(abc) contains d b outside C(acd) flip ab to cd

WHY?

a c b d

Lemma: d inside C(abc) <==> b outside C(abd)

• Why?
• basic geometry of angles (Thales theorem)

a c b a > b > c

a c b d

Note condition is symmetrical: d inside C(abc) <==> c inside C(abs)

• Why?
• basic geometry of angles (Thales theorem)

a c b a > b > c

Notation:

• edge ab is illegal iff d is inside C(abc)

(a legal edge also called a locally Delaunay edge)

• edge ab is legal iff not illegal (d on or outside C(abc))

a c b d a c b d a c b d

edge flips

Flipping edges and empty circles

ab illegal cd legal flip

Theorem, revisited: A triangulation where all edges are locally Delaunay (legal) is the DT. (locally Delaunay ==> globally Delaunay) a c b d a c b d

Flipping edges and empty circles

ab illegal cd legal flip ab to cd

DT via edge flipping

• construct an arbitrary triangulation T
• mark all edges in T and put them in S
• while S not empty
• pop edge ab from S and un-mark it
• if ab not legal:
• flip edge and update T
• for each new edge ac, ad, bc, bc: if

not marked, push it on S and mark it a c b d

Delaunay triangulation

Does this ever terminate?.. If it terminates, T is DT (by theorem)

a1 a2 a3 a4 a5 a6

d a b c

Edge flipping and angles

angle vector before A(T) = [a1, a2, a3, a4, a5, a6] after A(T’) = [b1, b2, b3, b4, b5, b6]

• Claim: For each new angle, there is an old angle that’s <= to it.

Proof: b1 = a1+a4 > a1 b2 > a5 (b outside C(acd)) b3 > … similar b4 = a3 + a5 > a3 b5 > … similar b6 > a3 (b outside C(acd))

b1

c b d

b2 b3 b4 b5 b6

a

Edge flipping and angles

Flipping an illegal edge to a legal edge increases minimum angle. Proof: Follows because each angle in T’ is larger than an angle in T.

• DT maximizes minimum angle.

Proof: DT has only legal edges. Any other flip would cause the min angel to decrease.

a1 a2 a3 a4 a5 a6

d a b c

b1

c b d

b2 b3 b4 b5 b6

a

More generally..

A triangulation T. consider all angles in each of its m triangles => 3m angles sort angles increasingly: a1 ≦ a2 ≦ a3 ≦…… ≦ a3m A(T) = [a1, a2, a3, ……,a3m]

• We say A(T) > A(T’) if A(T) > A(T’) lexicographically.

[2,5,6] > [1,2,3] [1,2,3] > [1,1,9]

• Of all possible triangulations, DT maximizes the angle vector.

That is, for all T’, A(DT) > A(T’). Proof: An edge flip increases the angle vector.

a1 a2 a3 a4 a5 a6

d a b c

b1

c b d

b2 b3 b4 b5 b6

a

DT Applications

• DT used because it is angle-optimal.
• Small angles
• Numerical instability, Interpolation errors
• Issues with geometric predicates (leftOf, tangents, derivatives)

Computing DT: RIC

RIC (Randomized Incremental Construction)

• Start with a large triangle T that contains P.
• Compute a random permutation of the points in P
• For each point p in P do:

//insert p in T

• find triangle abc of T that contains p
• splitting into 3 triangles: abp, bcp, cap and update T
• LegalizeEdge(p, ab, T)
• LegalizeEdge(p, bc, T)
• LegalizeEdge(p, ca, T)
• discard the initial triangle T and its edges
• return T

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

a p c b

Computing DT: RIC

LegalizeEdge(p, uv, T) //the point being inserted is p, and edge uv may need to be flipped

• if uv is illegal
• let uvq be the triangle adjacent to uv, on the other side of p
• flip uv and update T
• LegalizeEdge(p, uq, T)
• LegalizeEdge(p, qv, T)

P = {p1, p2,…,pn} set of points in the plane no 4 points co-circular

a p c b

a p c b

a p c b

a p c b q

a p c b q

a p c b q