class Delaunator::Triangulator
Constants
- EDGE_STACK
- EPSILON
Attributes
coords[R]
halfedges[R]
hull[R]
triangles[R]
Public Class Methods
new(coords)
click to toggle source
# File lib/delaunator/triangulator.rb, line 9 def initialize(coords) @coords = coords @n = @coords.length >> 1 @max_triangles = [2 * @n - 5, 0].max @triangles_len = 0 @_triangles = Array.new(@max_triangles * 3, 0) @_halfedges = Array.new(@max_triangles * 3, 0) @hash_size = Math.sqrt(@n).ceil @hull_prev = Array.new(@n, 0) @hull_next = Array.new(@n, 0) @hull_tri = Array.new(@n, 0) @hull_hash = Array.new(@hash_size, -1) end
Public Instance Methods
add_triangle(i0, i1, i2, a, b, c)
click to toggle source
# File lib/delaunator/triangulator.rb, line 321 def add_triangle(i0, i1, i2, a, b, c) i = @triangles_len @_triangles[i] = i0 @_triangles[i + 1] = i1 @_triangles[i + 2] = i2 link(i, a) link(i + 1, b) link(i + 2, c) @triangles_len += 3 i end
circumcenter(ax, ay, bx, by, cx, cy)
click to toggle source
# File lib/delaunator/triangulator.rb, line 262 def circumcenter(ax, ay, bx, by, cx, cy) dx = bx - ax dy = by - ay ex = cx - ax ey = cy - ay bl = dx * dx + dy * dy cl = ex * ex + ey * ey d = 0.5 / ( dx * ey - dy * ex).to_f x = ax + (ey * bl - dy * cl) * d y = ay + (dx * cl - ex * bl) * d [x, y] end
circumradius(ax, ay, bx, by, cx, cy)
click to toggle source
# File lib/delaunator/triangulator.rb, line 277 def circumradius(ax, ay, bx, by, cx, cy) dx = bx - ax dy = by - ay ex = cx - ax ey = cy - ay bl = dx * dx + dy * dy cl = ex * ex + ey * ey d = 0.5 / (dx * ey - dy * ex).to_f x = (ey * bl - dy * cl) * d y = (dx * cl - ex * bl) * d r = x * x + y * y return 0 if [bl, cl, d, r].include? 0 r end
dist(ax, ay, bx, by)
click to toggle source
# File lib/delaunator/triangulator.rb, line 244 def dist(ax, ay, bx, by) dx = ax - bx dy = ay - by dx * dx + dy * dy end
hash_key(x, y)
click to toggle source
# File lib/delaunator/triangulator.rb, line 333 def hash_key(x, y) angle = (pseudo_angle(x - @cx, y - @cy) * @hash_size).floor() angle % @hash_size end
in_circle(ax, ay, bx, by, cx, cy, px, py)
click to toggle source
# File lib/delaunator/triangulator.rb, line 297 def in_circle(ax, ay, bx, by, cx, cy, px, py) dx = ax - px dy = ay - py ex = bx - px ey = by - py fx = cx - px fy = cy - py ap = dx * dx + dy * dy bp = ex * ex + ey * ey cp = fx * fx + fy * fy ( dx * (ey * cp - bp * fy) - dy * (ex * cp - bp * fx) + ap * (ex * fy - ey * fx) ) < 0 end
legalize(a)
click to toggle source
# File lib/delaunator/triangulator.rb, line 343 def legalize(a) # if the pair of triangles doesn't satisfy the Delaunay condition # (p1 is inside the circumcircle of [p0, pl, pr]), flip them, # then do the same check/flip recursively for the new pair of triangles # # pl pl # /||\ / \ # al/ || \bl al/ \a # / || \ / \ # / a||b \ flip /___ar___\ # p0\ || /p1 => p0\---bl---/p1 # \ || / \ / # ar\ || /br b\ /br # \||/ \ / # pr pr i = 0 ar = 0 loop do b = @_halfedges[a] a0 = a - a % 3 ar = a0 + (a + 2) % 3 # convex hull edge if b == -1 break if i.zero? i -= 1 a = EDGE_STACK[i] next end b0 = b - b % 3 al = a0 + (a + 1) % 3 bl = b0 + (b + 2) % 3 p0 = @_triangles[ar] pr = @_triangles[a] pl = @_triangles[al] p1 = @_triangles[bl] illegal = in_circle(@coords[2 * p0], @coords[2 * p0 + 1], @coords[2 * pr], @coords[2 * pr + 1], @coords[2 * pl], @coords[2 * pl + 1], @coords[2 * p1], @coords[2 * p1 + 1]) if illegal @_triangles[a] = p1 @_triangles[b] = p0 hbl = @_halfedges[bl] # edge swapped on the other side of the hull (rare) # fix the halfedge reference if hbl == -1 e = @hull_start loop do if @hull_tri[e] == bl @hull_tri[e] = a break end e = @hull_prev[e] break if e != @hull_start end end link(a, hbl) link(b, @_halfedges[ar]) link(ar, bl) br = b0 + (b + 1) % 3 if i < EDGE_STACK.length EDGE_STACK[i] = br i += 1 end else break if i.zero? i -= 1 a = EDGE_STACK[i] end end ar end
link(a, b)
click to toggle source
# File lib/delaunator/triangulator.rb, line 338 def link(a, b) @_halfedges[a] = b @_halfedges[b] = a if b != -1 end
orient(rx, ry, qx, qy, px, py)
click to toggle source
a more robust orientation test that's stable in a given triangle (to fix robustness issues)
# File lib/delaunator/triangulator.rb, line 258 def orient(rx, ry, qx, qy, px, py) (orient_if_sure(px, py, rx, ry, qx, qy) || orient_if_sure(rx, ry, qx, qy, px, py) || orient_if_sure(qx, qy, px, py, rx, ry)) < 0 end
orient_if_sure(px, py, rx, ry, qx, qy)
click to toggle source
return 2d orientation sign if we're confident in it through J. Shewchuk's error bound check
# File lib/delaunator/triangulator.rb, line 251 def orient_if_sure(px, py, rx, ry, qx, qy) l = (ry - py) * (qx - px) r = (rx - px) * (qy - py) ((l - r).abs >= (3.3306690738754716e-16 * (l + r).abs)) ? l - r : 0 end
pseudo_angle(dx, dy)
click to toggle source
# File lib/delaunator/triangulator.rb, line 316 def pseudo_angle(dx, dy) p = dx / (dx.abs + dy.abs).to_f (dy > 0 ? 3 - p : 1 + p) / 4.to_f end
triangulate()
click to toggle source
# File lib/delaunator/triangulator.rb, line 24 def triangulate dists = Array.new(@n, 0.0) ids = Array.new(@n, 0) min_x = Float::INFINITY min_y = Float::INFINITY max_x = -Float::INFINITY max_y = -Float::INFINITY # compute bounds (0..@n - 1).each do |i| x = @coords[2 * i] y = @coords[2 * i + 1] min_x = x if x < min_x min_y = y if y < min_y max_x = x if x > max_x max_y = y if y > max_y end cx = (min_x + max_x).to_f / 2 cy = (min_y + max_y).to_f / 2 min_dist = Float::INFINITY i0 = 0 i1 = 0 i2 = 0 # pick a seed point close to midpoint i0 (0..@n - 1).each do |i| d = dist(cx, cy, @coords[2 * i], @coords[2 * i + 1]) if d < min_dist i0 = i min_dist = d end end i0x = @coords[2 * i0] i0y = @coords[2 * i0 + 1] min_dist = Float::INFINITY (0..@n - 1).each do |i| next if i == i0 d = dist(i0x, i0y, @coords[2 * i], @coords[2 * i + 1]) if (d < min_dist) && (d > 0) i1 = i min_dist = d end end i1x = @coords[2 * i1] i1y = @coords[2 * i1 + 1] min_radius = Float::INFINITY # find the third point which forms the smallest circumcircle with the first two (0..@n - 1).each do |i| next if (i == i0) || (i == i1) r = circumradius(i0x, i0y, i1x, i1y, @coords[2 * i], @coords[2 * i + 1]) if r < min_radius i2 = i min_radius = r end end i2x = @coords[2 * i2] i2y = @coords[2 * i2 + 1] if min_radius == Float::INFINITY # order collinear points by dx (or dy if all x are identical) # and return the list as a hull (0..n - 1).each do |i| dists[i] = (@coords[2 * i] - @coords[0]) || (@coords[2 * i + 1] - @coords[1]) end ids = dists.map.with_index.sort.map(&:last) @hull = Array.new(@n, 0) j = 0 d0 = -Float::INFINITY (0..@n - 1).each do |i| id = ids[i] next unless dists[id] > d0 @hull[j] = id j += 1 d0 = dists[id] end @hull = hull.slice(0, j) @triangles = [] @halfedges = [] return end # swap the order of the seed points for counter-clockwise orientation if orient(i0x, i0y, i1x, i1y, i2x, i2y) i = i1 x = i1x y = i1y i1 = i2 i1x = i2x i1y = i2y i2 = i i2x = x i2y = y end center = circumcenter(i0x, i0y, i1x, i1y, i2x, i2y) @cx, @cy = center (0..@n - 1).each do |k| dists[k] = dist(@coords[2 * k], @coords[2 * k + 1], center[0], center[1]) end # sort the points by distance from the seed triangle circumcenter ids = dists.map.with_index.sort.map(&:last) @hull_start = i0 @hull_size = 3 @hull_next[i0] = @hull_prev[i2] = i1 @hull_next[i1] = @hull_prev[i0] = i2 @hull_next[i2] = @hull_prev[i1] = i0 @hull_tri[i0] = 0 @hull_tri[i1] = 1 @hull_tri[i2] = 2 @hull_hash.fill(-1) @hull_hash[hash_key(i0x, i0y)] = i0 @hull_hash[hash_key(i1x, i1y)] = i1 @hull_hash[hash_key(i2x, i2y)] = i2 @triangles_len = 0 add_triangle(i0, i1, i2, -1, -1, -1) xp = 0 yp = 0 (0..ids.length - 1).each do |k| i = ids[k] x = @coords[2 * i] y = @coords[2 * i + 1] if (k > 0 && (x - xp).abs <= EPSILON && (y - yp).abs <= EPSILON) next end xp = x yp = y #skip seed triangle points if i == i0 || i == i1 || i == i2 next end # find a visible edge on the convex hull using edge hash start = 0 key = hash_key(x, y) (0..@hash_size - 1).each do |m| start = @hull_hash[(key + m) % @hash_size] break if (start != -1) && (start != @hull_next[start]) end start = @hull_prev[start] e = start q = nil loop do q = @hull_next[e] break if orient(x, y, @coords[2 * e], @coords[2 * e + 1], @coords[2 * q], @coords[2 * q + 1]) e = q if e == start e = -1 break end end if e == -1 next end t = add_triangle(e, i, @hull_next[e], -1, -1, @hull_tri[e]) @hull_tri[i] = legalize(t + 2) @hull_tri[e] = t @hull_size += 1 n = @hull_next[e] loop do q = @hull_next[n] break if !orient(x, y, @coords[2 * n], @coords[2 * n + 1], @coords[2 * q], @coords[2 * q + 1]) t = add_triangle(n, i, q, @hull_tri[i], -1, @hull_tri[n]) @hull_tri[i] = legalize(t + 2) @hull_next[n] = n @hull_size -= 1 n = q end if e == start loop do q = @hull_prev[e] break if !orient(x, y, @coords[2 * q], @coords[2 * q + 1], @coords[2 * e], @coords[2 * e + 1]) t = add_triangle(q, i, e, -1, @hull_tri[e], @hull_tri[q]) legalize(t + 2) @hull_tri[q] = t @hull_next[e] = e @hull_size -= 1 e = q end end @hull_start = e @hull_prev[i] = e @hull_next[e] = i @hull_prev[n] = i @hull_next[i] = n # save the two new edges in the hash table @hull_hash[hash_key(x, y)] = i @hull_hash[hash_key(@coords[2 * e], @coords[2 * e + 1])] = e end @hull = Array.new(@hull_size, 0) e = @hull_start (0..@hull_size - 1).each do |k| @hull[k] = e e = @hull_next[e] end # trim typed triangle mesh arrays @halfedges = @_halfedges.slice(0, @triangles_len) @triangles = @_triangles.slice(0, @triangles_len) @triangles end