Robust checks for orientation

Cargo.toml

@@ -11,6 +11,9 @@ categories = ["algorithms", "data-structures"]
 keywords = ["delaunay", "triangulation", "tessellation", "spatial", "geometry"]
 authors = ["Vladimir Agafonkin <agafonkin@gmail.com>"]
 
+[dependencies]
+robust = "0.2.3"
+
 [dev-dependencies]
 criterion = "0.3.4"
 rand = "0.8.3"
@@ -26,3 +29,11 @@ bench = false
 [[bench]]
 name = "bench"
 harness = false
+
+[[example]]
+name = "triangulate"
+path = "examples/triangulate.rs"
+
+[[example]]
+name = "svg"
+path = "examples/svg.rs"

README.md

@@ -20,7 +20,7 @@ let points = vec![
     Point { x: 0., y: 1. },
 ];
 
-let result = triangulate(&points).expect("No triangulation exists.");
+let result = triangulate(&points);
 
 println!("{:?}", result.triangles); // [0, 2, 1, 0, 3, 2]
 ```

examples/svg.rs

@@ -0,0 +1,81 @@
+use std::{env, fs::File, io::Write};
+use delaunator::{EMPTY, Point, Triangulation, next_halfedge};
+const CANVAS_SIZE: f64 = 800.;
+const POINT_SIZE: usize = 4;
+const LINE_WIDTH: usize = 1;
+const HULL_COLOR: &str = "green";
+const LINE_COLOR: &str = "blue";
+const POINT_COLOR: &str = "black";
+const HULL_POINT_COLOR: &str = "red";
+
+/// Takes the first argument and use as path to load points data. If no argument provided, loads one of the test fixtures data file
+/// Example: cargo run --example svg -- tests/fixtures/issue24.json
+fn main() -> std::io::Result<()> {
+    // load points from file
+    let default_path = "tests/fixtures/robust4.json".to_string();
+    let args = env::args().collect::<Vec<String>>();
+    let path = args.get(1).unwrap_or(&default_path);
+    let points: Vec<Point> = serde_json::from_reader::<_, Vec<(f64, f64)>>(File::open(path)?)?
+        .iter().map(|p| Point { x: p.0, y: p.1 }).collect();
+
+    // triangulate and scale points for display
+    let triangulation = delaunator::triangulate(&points);
+    println!("{:#?}", triangulation);
+    let points = center_and_scale(&points, &triangulation);
+
+    // generate SVG
+    let contents = format!(
+        r#"
+<svg viewBox="0 0 {width} {height}" xmlns="http://www.w3.org/2000/svg">
+<rect width="100%" height="100%" fill="white" />
+    {circles}
+    {lines}
+</svg>"#,
+        width = CANVAS_SIZE,
+        height = CANVAS_SIZE,
+        circles = render_point(&points, &triangulation),
+        lines = (0..triangulation.triangles.len()).fold(String::new(), |acc, e| {
+            if e > triangulation.halfedges[e] || triangulation.halfedges[e] == EMPTY {
+                let start = &points[triangulation.triangles[e]];
+                let end = &points[triangulation.triangles[next_halfedge(e)]];
+                let color = if triangulation.halfedges[e] == EMPTY { HULL_COLOR } else { LINE_COLOR };
+                acc + &format!(r#"<line x1="{x0}" y1="{y0}" x2="{x1}" y2="{y1}" style="stroke:{color};stroke-width:{width}" />"#, x0 = start.x, y0 = start.y, x1=end.x, y1=end.y, width = LINE_WIDTH, color = color)
+            } else {
+                acc
+            }
+        })
+    );
+    File::create("triangulation.svg")?
+        .write_all(contents.as_bytes())
+}
+
+
+/// Finds the center point and farthest point from it, then generates a new vector of
+/// scaled and offset points such that they fit between [0..SIZE]
+fn center_and_scale(points: &Vec<Point>, t: &Triangulation) -> Vec<Point> {
+    let center = &points[*t.triangles.get(0).unwrap_or(&0)];
+    let farthest_distance = points.iter().map(|p| {
+        let (x, y) = (center.x - p.x, center.y - p.y);
+        x*x + y*y
+    }).reduce(f64::max).unwrap().sqrt();
+    let scale = CANVAS_SIZE / (farthest_distance * 2.0);
+    let offset = ((CANVAS_SIZE / 2.0) - (scale * center.x), (CANVAS_SIZE / 2.0) - (scale * center.y));
+    points.iter().map(|p| Point { x: scale * p.x + offset.0, y: scale * p.y + offset.1 }).collect()
+}
+
+
+fn render_point(points: &[Point], triangulation: &Triangulation) -> String {
+    let mut circles = points.iter().enumerate().fold(String::new(), |acc, (i, p)| {
+        let color = if triangulation.hull.contains(&i) { HULL_POINT_COLOR } else { POINT_COLOR };
+        acc + &format!(r#"<circle cx="{x}" cy="{y}" r="{size}" fill="{color}"/>"#, x = p.x, y = p.y, size = POINT_SIZE, color = color)
+    });
+
+    // show ids for points if input is relatively small
+    if points.len() < 100 {
+        circles = points.iter().enumerate().fold(circles, |acc, (i, p)| {
+            acc + &format!(r#"<text x="{x}" y="{y}" font-size="20" fill="black">{i}</text>"#, i = i, x = p.x + 10., y = p.y - 5.)
+        })
+    }
+
+    circles
+}

examples/triangulate.rs

@@ -1,6 +1,3 @@
-extern crate delaunator;
-extern crate rand;
-
 use std::iter::repeat_with;
 
 const N: usize = 1_000_000;

src/lib.rs

@@ -27,6 +27,7 @@ extern crate std;
 #[macro_use]
 extern crate alloc;
 
+use robust::orient2d;
 use core::{f64, fmt};
 use alloc::vec::Vec;
 
@@ -47,15 +48,24 @@ impl fmt::Debug for Point {
     }
 }
 
+impl Into<robust::Coord<f64>> for &Point {
+    fn into(self) -> robust::Coord<f64> {
+        robust::Coord::<f64> {
+            x: self.x,
+            y: self.y,
+        }
+    }
+}
+
 impl Point {
     fn dist2(&self, p: &Self) -> f64 {
         let dx = self.x - p.x;
         let dy = self.y - p.y;
         dx * dx + dy * dy
     }
 
-    fn orient(&self, q: &Self, r: &Self) -> bool {
-        (q.y - self.y) * (r.x - q.x) - (q.x - self.x) * (r.y - q.y) < 0.0
+    fn orient(&self, q: &Self, r: &Self) -> f64 {
+        orient2d(self.into(), q.into(), r.into())
     }
 
     fn circumdelta(&self, b: &Self, c: &Self) -> (f64, f64) {
@@ -130,6 +140,7 @@ pub fn prev_halfedge(i: usize) -> usize {
 }
 
 /// Result of the Delaunay triangulation.
+#[derive(Debug, Clone)]
 pub struct Triangulation {
     /// A vector of point indices where each triple represents a Delaunay triangle.
     /// All triangles are directed counter-clockwise.
@@ -347,7 +358,7 @@ impl Hull {
         start = self.prev[start];
         let mut e = start;
 
-        while !p.orient(&points[e], &points[self.next[e]]) {
+        while p.orient(&points[e], &points[self.next[e]]) <= 0. {
             e = self.next[e];
             if e == start {
                 return (EMPTY, false);
@@ -419,7 +430,7 @@ fn find_seed_triangle(points: &[Point]) -> Option<(usize, usize, usize)> {
         None
     } else {
         // swap the order of the seed points for counter-clockwise orientation
-        Some(if p0.orient(p1, &points[i2]) {
+        Some(if p0.orient(p1, &points[i2]) > 0. {
             (i0, i2, i1)
         } else {
             (i0, i1, i2)
@@ -516,7 +527,7 @@ pub fn triangulate(points: &[Point]) -> Triangulation {
         let mut n = hull.next[e];
         loop {
             let q = hull.next[n];
-            if !p.orient(&points[n], &points[q]) {
+            if p.orient(&points[n], &points[q]) < 0. {
                 break;
             }
             let t = triangulation.add_triangle(n, i, q, hull.tri[i], EMPTY, hull.tri[n]);
@@ -529,7 +540,7 @@ pub fn triangulate(points: &[Point]) -> Triangulation {
         if walk_back {
             loop {
                 let q = hull.prev[e];
-                if !p.orient(&points[q], &points[e]) {
+                if p.orient(&points[q], &points[e]) < 0. {
                     break;
                 }
                 let t = triangulation.add_triangle(q, i, e, EMPTY, hull.tri[e], hull.tri[q]);
@@ -610,7 +621,7 @@ fn f64_sqrt(f: f64) -> f64 {
 
     let sc = f64_sqrt(f / 4.0) * 2.0;
     let lc = sc + 1.0;
-    
+
     if lc * lc > f {
         sc
     } else {

tests/fixtures/grid.json

@@ -0,0 +1,34 @@
+[
+    [0.5071359338391455, -0.3236832077937956],
+    [0.3830151816637668, 0.5627457940628474],
+    [1.5196837554166323, -0.2974359857979777],
+    [-0.6149866437812772, 1.0450423584955024],
+    [-0.21916717030359378, -2.702309985665537],
+    [-0.3128297966002196, 0.32627037936130365],
+    [-0.2585335697093173, -3.2042352606594826],
+    [-1.05847858399307686, 2.3754785777859544],
+    [1.0, 1.0],
+    [-1.0, 1.0],
+    [-1.0, -1.0],
+    [1.0, -1.0],
+    [4.666589448684871, -1.160460175259872],
+    [4.95915032227072, 0.7037936855165575],
+    [-2.2472073777517436, 3.363632890638188],
+    [-0.30704845634464, -4.6397107833345075],
+    [-0.5217818413507629, 1.0],
+    [1.0, 0.5817628040530143],
+    [1.0, -0.3109072062236856],
+    [1.0, 0.09583827903291398],
+    [-0.5217818413507627, 1.0],
+    [-0.59605176954064, 1.0],
+    [-0.2718132304410775, -1.0],
+    [0.3006256102638203, -1.0],
+    [-0.7824622340091141, -1.0],
+    [-1.0, -0.08479064614066445],
+    [1.0, 0.09583827903291409],
+    [1.0, 0.5817628040530143],
+    [-1.0, -0.08479064614066623],
+    [-0.7824622340091145, -1.0],
+    [-0.2718132304410775, -1.0],
+    [-0.5960517695406401, 1.0]
+]

tests/fixtures/issue10.json

@@ -0,0 +1,22 @@
+[
+    [-1.0849334460551594, -1.0849334460551594],
+    [-1.5265847189170323, -1.5265847189170323],
+    [-2.095815826893781, -1.274850699584578],
+    [-2.770865690566727, -0.9763199041819535],
+    [-3.6843898126113546, -0.5723274031100705],
+    [-4.403658177176129, -0.25424163117441645],
+    [-5.02567003534954, 0.020833892521026076],
+    [-5.701084620214019, 0.3195259804640507],
+    [-6.561463270870451, 0.7000156846325452],
+    [-7.31105511135829, 1.0315115642859167],
+    [-8.0, 1.336187228503853],
+    [-7.339371017948894, 1.1048861305058357],
+    [-6.616986689211032, 0.8519630935920505],
+    [-5.816767071935726, 0.5717881676590966],
+    [-5.121128245254447, 0.3282293338915143],
+    [-4.512948676142796, 0.11529195763725286],
+    [-3.850960067633096, -0.11648517623155441],
+    [-3.1594534122386113, -0.3585972436874558],
+    [-2.288934450277422, -0.663385554827794],
+    [-1.6751076145244035, -0.8783001664294665]
+]

tests/fixtures/issue43js.json

@@ -0,0 +1,7 @@
+[
+    [-537.7739674441619, -122.26130468750004],
+    [-495.533967444162, -183.39195703125006],
+    [-453.29396744416204, -244.5226093750001],
+    [-411.0539674441621, -305.6532617187501],
+    [-164, -122]
+  ]

tests/fixtures/robust4.json

@@ -0,0 +1 @@
+[[-3.0381276552207055, 10.481881920449052], [-2.1931270567413446, 11.016647278872279], [-1.3481264582619854, 11.551412637295508], [-0.5031258597826245, 12.086177995718735], [0.3418747386967347, 12.620943354141964], [1.1868753371760938, 13.155708712565193], [-2.5033622967974765, 9.63688132196969], [-1.6583616983181173, 10.171646680392918], [-0.8133610998387582, 10.706412038816147], [0.03163949864060278, 11.241177397239376], [0.8766400971199619, 11.775942755662605], [1.721640695599322, 12.310708114085832], [-1.9685969383742474, 8.791880723490332], [-1.1235963398948883, 9.326646081913559], [-0.2785957414155291, 9.861411440336788], [0.5664048570638318, 10.396176798760017], [1.411405455543191, 10.930942157183246], [2.25640605402255, 11.465707515606473], [-1.4338315799510184, 7.9468801250109715], [-0.5888309814716592, 8.4816454834342], [0.2561696170076999, 9.016410841857429], [1.10117021548706, 9.551176200280658], [1.94617081396642, 10.085941558703885], [2.791171412445779, 10.620706917127112], [-0.8990662215277911, 7.1018795265316115], [-0.05406562304843021, 7.6366448849548405], [0.7909349754309281, 8.17141024337807], [1.635935573910288, 8.706175601801297], [2.4809361723896473, 9.240940960224526], [3.3259367708690073, 9.775706318647753], [-0.3643008631045621, 6.256878928052252], [0.48069973537479704, 6.7916442864754805], [1.3257003338541562, 7.326409644898709], [2.1707009323335162, 7.861175003321938], [3.0157015308128763, 8.395940361745165], [3.8607021292922354, 8.930705720168394]]

tests/fixtures/robust5.json

@@ -0,0 +1 @@
+[[1.2000108248775223,0.08939813876119695],[3.197878575122478,0.08939813876119695],[3.9347512248775223,0.08939813876119695],[5.932618975122478,0.08939813876119695],[6.669491624877522,0.08939813876119695],[1.2000108248775223,1.9106018612388114],[3.197878575122478,1.9106018612388114],[3.9347512248775223,1.9106018612388114],[5.932618975122478,1.9106018612388114],[6.669491624877522,1.9106018612388114],[2.0267626140791766,-0.21046296655057145],[2.371126785920824,-0.21046296655057145],[4.761503014079176,-0.21046296655057145],[5.105867185920824,-0.21046296655057145],[7.496243414079176,-0.21046296655057145],[2.0267626140791766,2.2104629665505797],[2.371126785920824,2.2104629665505797],[4.761503014079176,2.2104629665505797],[5.105867185920824,2.2104629665505797],[7.496243414079176,2.2104629665505797],[1.5802128246274927,0.22503184151033584],[2.8176765753725075,0.22503184151033584],[4.314953224627493,0.22503184151033584],[5.552416975372507,0.22503184151033584],[7.049693624627492,0.22503184151033584],[1.5802128246274927,1.7749681584896724],[2.8176765753725075,1.7749681584896724],[4.314953224627493,1.7749681584896724],[5.552416975372507,1.7749681584896724],[7.049693624627492,1.7749681584896724],[0.8315745000000003,1.000000000000004],[3.5663149,1.000000000000004],[3.5663149,1.000000000000004],[6.3010553,1.000000000000004],[6.3010553,1.000000000000004],[0.8315745000000003,1.0000000000000044],[3.5663149,1.0000000000000044],[3.5663149,1.0000000000000044],[6.3010553,1.0000000000000044],[6.3010553,1.0000000000000044],[0.8315745000000003,-0.488685099999997],[3.5663149,-0.488685099999997],[3.5663149,-0.488685099999997],[6.3010553,-0.488685099999997],[6.3010553,-0.488685099999997],[0.8315745000000003,2.4886851000000054],[3.5663149,2.4886851000000054],[3.5663149,2.4886851000000054],[6.3010553,2.4886851000000054],[6.3010553,2.4886851000000054],[2.1989447,-0.61],[2.1989447,-0.61],[4.9336851,-0.61],[4.9336851,-0.61],[7.6684255,-0.61],[2.1989447,2.6100000000000083],[2.1989447,2.6100000000000083],[4.9336851,2.6100000000000083],[4.9336851,2.6100000000000083],[7.6684255,2.6100000000000083],[2.1989447,0.878685100000001],[2.1989447,0.878685100000001],[4.9336851,0.878685100000001],[4.9336851,0.878685100000001],[7.6684255,0.878685100000001],[2.1989447,1.1213149000000073],[2.1989447,1.1213149000000073],[4.9336851,1.1213149000000073],[4.9336851,1.1213149000000073],[7.6684255,1.1213149000000073]]

tests/fixtures/robust6.json

@@ -0,0 +1,7 @@
+[
+    [0, 0],
+    [0.0000000000000002220446049250313, 0.0000000000000002220446049250313],
+    [0.0000000000000004440892098500626, 0.0000000000000002220446049250313],
+    [-0.0000000000000004440892098500626, 0.0000000000000002220446049250313],
+    [-0.0000000000000005551115123125783, 0.0000000000000004440892098500626]
+]