diff --git a/CHANGELOG.md b/CHANGELOG.md
index 24c1a33..d6e6518 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -4,6 +4,13 @@ All notable changes to this project will be documented in this file.
The format is based on [Keep a Changelog](http://keepachangelog.com/)
and this project adheres to [Semantic Versioning](http://semver.org/).
+## [2.5.0] - yyyy-mm-dd
+
+### Added
+ - Implements natural neighbor interpolation. See type `NaturalNeighbor` for more on this!
+### Fix
+ - Fixes `PointProjection::reversed()` - the method would return a projection that gives sometimes incorrect results.
+
## [2.4.1] - 2023-12-01
### Fix
diff --git a/Cargo.toml b/Cargo.toml
index 863b1d0..9b4274f 100644
--- a/Cargo.toml
+++ b/Cargo.toml
@@ -38,8 +38,11 @@ rand = "0.8.3"
cgmath = "0.18.0"
svg = "0.14.0"
float_next_after = "1"
+image = "0.24.7"
+tiny-skia = "0.11.3"
criterion = { version = "0.5.1", features = ["html_reports"] }
-
+base64 = "0.21.5"
+anyhow = "1.0.75"
[[bench]]
name = "benchmarks"
diff --git a/README.md b/README.md
index 1146322..f39fc3d 100644
--- a/README.md
+++ b/README.md
@@ -6,14 +6,15 @@
Delaunay triangulations for the rust ecosystem.
-- 2D [Delaunay triangulation](https://en.wikipedia.org/wiki/Delaunay_triangulation), optionally backed by a hierarchy
- structure for improved nearest neighbor and insertion performance.
-- Allows both incremental and bulk loading creation of triangulations
-- Support for vertex removal
-- 2D constrained Delaunay triangulation (CDT)
-- [Delaunay refinement](https://en.wikipedia.org/wiki/Delaunay_refinement)
-- Uses precise geometric predicates to prevent incorrect geometries due to rounding issues
-- Supports extracting the Voronoi diagram
+ - 2D [Delaunay triangulation](https://en.wikipedia.org/wiki/Delaunay_triangulation), optionally backed
+ by a hierarchy structure for improved nearest neighbor and insertion performance.
+ - Allows both incremental and bulk loading creation of triangulations
+ - Support for vertex removal
+ - 2D constrained Delaunay triangulation (CDT)
+ - [Delaunay refinement](https://en.wikipedia.org/wiki/Delaunay_refinement)
+ - Uses precise geometric predicates to prevent incorrect geometries due to rounding issues
+ - Supports extracting the Voronoi diagram
+ - Natural neighbor interpolation
---------------------------
@@ -31,19 +32,10 @@ Project goals, in the order of their importance:
# Roadmap
-For Spade 2.x:
- - Add back the removed interpolation methods (natural neighbor interpolation, #67)
-
For Spade 3:
- - Possibly base `spade` on `nalgebra` as underlying vector and matrix library. Not much else planned yet!
-
-# Project state and contributing
-
-Looking for co-maintainers! Projects with just a single maintainer can be a little unreliable due to the single point of failure. I would love to see this burden being distributed on more shoulders! This is less about *implementing* things but rather about general maintenance tasks, e.g. package updates, minor bug fixes, reviewing PRs, etc...
-
-If you want to contribute, please consider opening an issue first. I'm happy to help out with any questions!
+ - Possibly API simplification by un-supporting non-f64 outputs.
-# Performance and comparison to other Delaunay crates
+# Performance and comparison to other crates
Refer to [the delaunay_compare readme](./delaunay_compare/README.md) for some benchmarks and a comparison with other crates.
diff --git a/benches/benchmarks.rs b/benches/benchmarks.rs
index 4b11b85..707cc13 100644
--- a/benches/benchmarks.rs
+++ b/benches/benchmarks.rs
@@ -5,15 +5,17 @@ use criterion::*;
mod benchmark_utilities;
mod bulk_load_benchmark;
mod bulk_load_vs_incremental;
+mod interpolation_benchmark;
mod locate_benchmark;
criterion_group! {
name = benches;
config = Criterion::default().sample_size(50).warm_up_time(Duration::from_secs(1));
targets =
- locate_benchmark::locate_benchmark,
bulk_load_benchmark::bulk_load_benchmark,
bulk_load_vs_incremental::bulk_load_vs_incremental_benchmark,
+ interpolation_benchmark::interpolation_benchmark,
+ locate_benchmark::locate_benchmark,
}
criterion_main!(benches);
diff --git a/benches/interpolation_benchmark.rs b/benches/interpolation_benchmark.rs
new file mode 100644
index 0000000..a352ab8
--- /dev/null
+++ b/benches/interpolation_benchmark.rs
@@ -0,0 +1,47 @@
+use criterion::*;
+use spade::Triangulation;
+use spade::{DelaunayTriangulation, FloatTriangulation, Point2};
+
+use crate::benchmark_utilities::*;
+
+pub fn interpolation_benchmark(c: &mut Criterion) {
+ let mut group = c.benchmark_group("interpolation benchmarks");
+
+ let points = uniform_distribution(*SEED, 1.0)
+ .take(50)
+ .collect::>();
+ let triangulation = DelaunayTriangulation::<_>::bulk_load(points).unwrap();
+
+ let query_point = Point2::new(0.5, -0.2);
+
+ group.bench_function("nearest neighbor interpolation", |b| {
+ b.iter(|| {
+ triangulation
+ .nearest_neighbor(query_point)
+ .unwrap()
+ .position()
+ });
+ });
+
+ let barycentric = triangulation.barycentric();
+ group.bench_function("barycentric interpolation", |b| {
+ b.iter(|| barycentric.interpolate(|v| v.position().x + v.position().y, query_point));
+ });
+
+ let natural_neighbor = &triangulation.natural_neighbor();
+ group.bench_function("natural neighbor interpolation (c0)", |b| {
+ b.iter(|| natural_neighbor.interpolate(|v| v.position().x + v.position().y, query_point));
+ });
+
+ group.bench_function("natural neighbor interpolation (c1)", |b| {
+ b.iter(|| {
+ natural_neighbor.interpolate_gradient(
+ |v| v.position().x + v.position().y,
+ |_| [0.0, 0.0],
+ 0.5,
+ query_point,
+ )
+ });
+ });
+ group.finish();
+}
diff --git a/examples/interpolation/main.rs b/examples/interpolation/main.rs
new file mode 100644
index 0000000..1fc9526
--- /dev/null
+++ b/examples/interpolation/main.rs
@@ -0,0 +1,214 @@
+use anyhow::*;
+use image::{ImageBuffer, Rgba};
+use std::io::{Cursor, Write};
+
+use base64::Engine;
+use tiny_skia::*;
+
+use spade::{DelaunayTriangulation, FloatTriangulation, HasPosition, Point2, Triangulation};
+
+#[derive(Debug, Copy, Clone)]
+struct PointWithHeight {
+ position: Point2,
+ height: f64,
+}
+
+impl PointWithHeight {
+ const fn new(x: f64, y: f64, height: f64) -> Self {
+ Self {
+ position: Point2::new(x, y),
+ height,
+ }
+ }
+}
+
+impl HasPosition for PointWithHeight {
+ type Scalar = f64;
+
+ fn position(&self) -> Point2 {
+ self.position
+ }
+}
+
+type TriangulationType = DelaunayTriangulation;
+
+const VERTICES: &[PointWithHeight] = &[
+ PointWithHeight::new(20.0, 20.0, 0.9),
+ PointWithHeight::new(20.0, -20.0, 0.9),
+ PointWithHeight::new(-20.0, 20.0, 0.9),
+ PointWithHeight::new(-20.0, -20.0, 0.9),
+ PointWithHeight::new(20.0, 10.0, 0.0),
+ PointWithHeight::new(10.0, 10.0, 0.5),
+ PointWithHeight::new(0.0, 23.0, 0.1),
+ PointWithHeight::new(0.0, -23.0, 0.1),
+ PointWithHeight::new(-10.0, 10.0, 0.0),
+ PointWithHeight::new(-10.0, -10.0, 0.1),
+ PointWithHeight::new(-15.0, 20.0, 0.5),
+ PointWithHeight::new(-20.0, 20.0, 0.0),
+ PointWithHeight::new(-5.0, 0.0, 0.25),
+ PointWithHeight::new(5.0, 7.0, 0.75),
+ PointWithHeight::new(12.0, -10.0, 0.4),
+ PointWithHeight::new(5.0, 10.0, 0.3),
+ PointWithHeight::new(5.0, 1.0, 0.2),
+];
+
+pub fn main() -> anyhow::Result<()> {
+ let mut t = TriangulationType::default();
+ for vertex in VERTICES {
+ t.insert(*vertex)?;
+ }
+
+ const OFFSET: f32 = 25.0;
+ const SCALAR: f32 = 50.0 / 512.0;
+
+ fn px_to_coords(x: u32, y: u32) -> Point2 {
+ Point2::new(
+ x as f64 * SCALAR as f64 - OFFSET as f64,
+ y as f64 * SCALAR as f64 - OFFSET as f64,
+ )
+ }
+
+ let dimensions = 512;
+ let mut nn_c0_pixmap =
+ Pixmap::new(dimensions, dimensions).context("Failed to allocate image")?;
+ let mut nn_c1_pixmap =
+ Pixmap::new(dimensions, dimensions).context("Failed to allocate image")?;
+ let mut barycentric_pixmap =
+ Pixmap::new(dimensions, dimensions).context("Failed to allocate image")?;
+ let mut nearest_neighbor_pixmap =
+ Pixmap::new(dimensions, dimensions).context("Failed to allocate image")?;
+
+ fn set_pixel(pixmap: &mut Pixmap, x: u32, y: u32, value: Option) {
+ let background_color = [255, 255, 255];
+ let [r, g, b] = value.map(float_to_color).unwrap_or(background_color);
+ let base = (y * pixmap.height() + x) as usize * 4;
+ pixmap.data_mut()[base] = r;
+ pixmap.data_mut()[base + 1] = g;
+ pixmap.data_mut()[base + 2] = b;
+ // Alpha
+ pixmap.data_mut()[base + 3] = 255;
+ }
+
+ let nn = t.natural_neighbor();
+ let barycentric = t.barycentric();
+
+ for y in 0..dimensions {
+ for x in 0..dimensions {
+ let coords = px_to_coords(x, y);
+ let value_nn_c0 = nn.interpolate(|v| v.data().height, coords);
+ set_pixel(&mut nn_c0_pixmap, x, y, value_nn_c0);
+
+ let value_nn_c1 =
+ nn.interpolate_gradient(|v| v.data().height, |_| [0.0, 0.0], 1.0, coords);
+ set_pixel(&mut nn_c1_pixmap, x, y, value_nn_c1);
+
+ let value_barycentric = barycentric.interpolate(|v| v.data().height, coords);
+ set_pixel(&mut barycentric_pixmap, x, y, value_barycentric);
+
+ let value_nearest_neighbor = t.nearest_neighbor(coords).map(|v| v.data().height);
+ set_pixel(&mut nearest_neighbor_pixmap, x, y, value_nearest_neighbor);
+ }
+ }
+
+ let path = {
+ let mut pb = PathBuilder::new();
+
+ for edge in t.undirected_edges() {
+ let [from, to] = edge.positions();
+ pb.move_to(from.x as f32, from.y as f32);
+ pb.line_to(to.x as f32, to.y as f32);
+ }
+
+ pb.finish().context("Could not build path")?
+ };
+
+ let stroke = Stroke {
+ width: 0.15,
+ ..Default::default()
+ };
+
+ let mut paint = Paint::default();
+ paint.set_color_rgba8(0, 0, 0, 70);
+ paint.anti_alias = true;
+
+ for pixmap in [
+ &mut nn_c0_pixmap,
+ &mut nn_c1_pixmap,
+ &mut barycentric_pixmap,
+ &mut nearest_neighbor_pixmap,
+ ] {
+ pixmap.stroke_path(
+ &path,
+ &paint,
+ &stroke,
+ Transform::from_translate(OFFSET, OFFSET).post_scale(1.0 / SCALAR, 1.0 / SCALAR),
+ None,
+ );
+ }
+
+ fn save_pixmap(pixmap: Pixmap, name: &str) -> anyhow::Result<()> {
+ // tiny_skia doesn't support jpg encoding which is required for small file size when embedding this into
+ // the documentation. We'll have to convert the data into ImageBuffer from the image crate and then do
+ // the jpeg encoding.
+ let (width, height) = (pixmap.width(), pixmap.height());
+ let data = pixmap.take();
+ let buffer = ImageBuffer::, _>::from_vec(width, height, data)
+ .context("Failed to convert to ImageBuffer")?;
+
+ let mut data_jpeg: Cursor> = Cursor::new(Vec::new());
+ buffer.write_to(&mut data_jpeg, image::ImageFormat::Jpeg)?;
+
+ std::fs::write(format!("images/{}.jpg", name), data_jpeg.get_ref())?;
+
+ // Encode image as ![]()
tag for inclusion into the documentation
+ let encoded = base64::engine::general_purpose::STANDARD_NO_PAD.encode(data_jpeg.get_ref());
+ let mut file = std::fs::File::create(format!("images/{}.img", name))?;
+ write!(file, r#"
"#, encoded)?;
+ Ok(())
+ }
+
+ save_pixmap(nn_c0_pixmap, "interpolation_nn_c0")?;
+ save_pixmap(nn_c1_pixmap, "interpolation_nn_c1")?;
+ save_pixmap(barycentric_pixmap, "interpolation_barycentric")?;
+ save_pixmap(nearest_neighbor_pixmap, "interpolation_nearest_neighbor")?;
+
+ Ok(())
+}
+
+fn float_to_color(value: f64) -> [u8; 3] {
+ // mostly AI generated..
+ // Converts a hue value in the range 0.0 ..= 1.0 from HLS to RGB
+ let value = value.clamp(0.0, 1.0);
+
+ const LIGHTNESS: f64 = 0.45;
+ const SATURATION: f64 = 0.55;
+
+ let c = (1.0 - (2.0 * LIGHTNESS - 1.0).abs()) * SATURATION;
+
+ let hue = value * 360.0;
+ let hs = hue / 60.0;
+
+ let x = c * (1.0 - (hs % 2.0 - 1.0).abs());
+ let m = LIGHTNESS - c * 0.5;
+
+ let (r, g, b) = if (0.0..60.0).contains(&hue) {
+ (c, x, 0.0)
+ } else if (60.0..120.0).contains(&hue) {
+ (x, c, 0.0)
+ } else if (120.0..180.0).contains(&hue) {
+ (0.0, c, x)
+ } else if (180.0..240.0).contains(&hue) {
+ (0.0, x, c)
+ } else if (240.0..300.0).contains(&hue) {
+ (x, 0.0, c)
+ } else {
+ (c, 0.0, x)
+ };
+
+ // Convert RGB to 8-bit values
+ let r = ((r + m) * 255.0) as u8;
+ let g = ((g + m) * 255.0) as u8;
+ let b = ((b + m) * 255.0) as u8;
+
+ [r, g, b]
+}
diff --git a/examples/svg_renderer/main.rs b/examples/svg_renderer/main.rs
index 2b8a8db..febb032 100644
--- a/examples/svg_renderer/main.rs
+++ b/examples/svg_renderer/main.rs
@@ -1,12 +1,26 @@
pub mod quicksketch;
mod scenario;
mod scenario_list;
-
-type Result = core::result::Result<(), Box>;
+use anyhow::Result;
/// Used for rendering SVGs for documentation. These are inlined (via #[doc = include_str!(...)])
/// into the doc comment of a few items. That makes sure they will be visible even for offline users.
-fn main() -> Result {
+fn main() -> Result<()> {
+ scenario_list::natural_neighbor_area_scenario(false)?.save_to_svg(
+ "natural_neighbor_insertion_cell",
+ "images/natural_neighbor_insertion_cell.svg",
+ )?;
+
+ scenario_list::natural_neighbor_area_scenario(true)?.save_to_svg(
+ "natural_neighbor_polygon",
+ "images/natural_neighbor_polygon.svg",
+ )?;
+
+ scenario_list::natural_neighbors_scenario().save_to_svg(
+ "natural_neighbor_scenario",
+ "images/natural_neighbor_scenario.svg",
+ )?;
+
scenario_list::refinement_maximum_area_scenario(None).save_to_svg(
"refinement_maximum_area_no_limit",
"images/refinement_maximum_area_no_limit.svg",
diff --git a/examples/svg_renderer/quicksketch/mod.rs b/examples/svg_renderer/quicksketch/mod.rs
index 1909ca4..dd86461 100644
--- a/examples/svg_renderer/quicksketch/mod.rs
+++ b/examples/svg_renderer/quicksketch/mod.rs
@@ -125,6 +125,7 @@ pub struct Style {
stroke_color: Option,
stroke_style: Option,
fill: Option,
+ opacity: Option,
marker_start: Option,
marker_end: Option,
}
@@ -144,6 +145,11 @@ impl Style {
.as_ref()
.map(|color| format!("stroke: {}", color));
+ let opacity = self
+ .opacity
+ .as_ref()
+ .map(|opacity| format!("opacity: {:.2}", opacity));
+
let fill = format!(
"fill: {}",
self.fill
@@ -175,6 +181,7 @@ impl Style {
stroke_color,
marker_start,
marker_end,
+ opacity,
stroke_dash_array,
])
.flatten()
@@ -356,6 +363,11 @@ impl SketchPath {
self
}
+ pub fn opacity(mut self, opacity: f64) -> Self {
+ self.style.opacity = Some(opacity);
+ self
+ }
+
pub fn fill(mut self, fill: SketchFill) -> Self {
self.style.fill = Some(fill);
self
diff --git a/examples/svg_renderer/scenario.rs b/examples/svg_renderer/scenario.rs
index c7e8b51..5791676 100644
--- a/examples/svg_renderer/scenario.rs
+++ b/examples/svg_renderer/scenario.rs
@@ -3,9 +3,10 @@ use spade::{ConstrainedDelaunayTriangulation, DelaunayTriangulation, HasPosition
use crate::convert_point;
-#[derive(Copy, Clone)]
+#[derive(Copy, Clone, PartialEq, Debug)]
pub struct VertexType {
position: Point,
+ pub color: Option,
pub radius: f64,
}
@@ -21,6 +22,7 @@ impl VertexType {
Self {
position: Point::new(x, y),
+ color: None,
radius: DEFAULT_CIRCLE_RADIUS,
}
}
@@ -34,6 +36,7 @@ impl HasPosition for VertexType {
}
}
+#[derive(Clone, Copy, Debug)]
pub struct UndirectedEdgeType {
pub color: SketchColor,
}
@@ -52,7 +55,7 @@ impl Default for UndirectedEdgeType {
}
}
-#[derive(Default)]
+#[derive(Clone, Copy, Debug, Default)]
pub struct DirectedEdgeType {}
#[derive(Clone, Copy, Debug)]
@@ -166,7 +169,9 @@ where
for vertex in triangulation.vertices() {
sketch.add_with_layer(
SketchElement::circle(convert_point(vertex.position()), vertex.data().radius)
- .fill(SketchFill::solid(options.vertex_color))
+ .fill(SketchFill::solid(
+ vertex.data().color.unwrap_or(options.vertex_color),
+ ))
.stroke_width(0.5)
.stroke_color(options.vertex_stroke_color),
crate::quicksketch::SketchLayer::VERTICES,
diff --git a/examples/svg_renderer/scenario_list.rs b/examples/svg_renderer/scenario_list.rs
index 042848f..55ea56b 100644
--- a/examples/svg_renderer/scenario_list.rs
+++ b/examples/svg_renderer/scenario_list.rs
@@ -3,13 +3,17 @@ use super::quicksketch::{
SketchLayer, StrokeStyle, Vector,
};
+use anyhow::{Context, Result};
+
use cgmath::{Angle, Bounded, Deg, EuclideanSpace, InnerSpace, Point2, Vector2};
+
use spade::{
handles::{
FixedDirectedEdgeHandle,
VoronoiVertex::{self, Inner, Outer},
},
- AngleLimit, FloatTriangulation as _, InsertionError, RefinementParameters, Triangulation as _,
+ AngleLimit, FloatTriangulation as _, HasPosition, InsertionError, RefinementParameters,
+ Triangulation as _,
};
use crate::{
@@ -974,6 +978,254 @@ pub fn shape_iterator_scenario(use_circle_metric: bool, iterate_vertices: bool)
result
}
+pub fn natural_neighbors_scenario() -> Sketch {
+ let triangulation = big_triangulation().unwrap();
+
+ let nn = triangulation.natural_neighbor();
+
+ let mut nns = Vec::new();
+ let query_point = spade::Point2::new(-1.0, 4.0);
+ nn.get_weights(query_point, &mut nns);
+
+ let mut result = convert_triangulation(&triangulation, &Default::default());
+
+ let offsets = [
+ Vector2::new(2.0, 5.0),
+ Vector2::new(-4.0, 6.0),
+ Vector2::new(0.0, 6.0),
+ Vector2::new(1.0, 6.0),
+ Vector2::new(0.0, 8.0),
+ Vector2::new(3.0, 4.0),
+ Vector2::new(5.0, 2.0),
+ ];
+
+ for (index, (neighbor, weight)) in nns.iter().enumerate() {
+ let neighbor = triangulation.vertex(*neighbor);
+ let position = convert_point(neighbor.position());
+ result.add(
+ SketchElement::circle(position, 2.0)
+ .fill(SketchFill::Solid(SketchColor::SALMON))
+ .stroke_width(0.5)
+ .stroke_color(SketchColor::BLACK),
+ );
+ result.add(
+ SketchElement::text(index.to_string())
+ .position(position)
+ .font_size(3.0)
+ .horizontal_alignment(HorizontalAlignment::Middle)
+ .dy(1.15),
+ );
+
+ result.add(
+ SketchElement::text(format!("{:.2}", weight))
+ .position(position + offsets[index])
+ .font_size(5.0),
+ );
+ }
+
+ result.add(
+ SketchElement::circle(convert_point(query_point), 1.5)
+ .fill(SketchFill::Solid(SketchColor::ROYAL_BLUE))
+ .stroke_width(0.2)
+ .stroke_color(SketchColor::BLACK),
+ );
+
+ result.set_view_box_min(Point2::new(-40.0, -30.0));
+ result.set_view_box_max(Point2::new(30.0, 45.0));
+ result.set_width(400);
+ result
+}
+
+/// Only used for internal documentation of natural neighbor area calculation
+pub fn natural_neighbor_area_scenario(include_faces: bool) -> Result {
+ let mut triangulation: Triangulation = Default::default();
+
+ triangulation.insert(VertexType::new(45.0, 30.0))?;
+ triangulation.insert(VertexType::new(7.5, 40.0))?;
+ triangulation.insert(VertexType::new(-45.0, 42.0))?;
+ triangulation.insert(VertexType::new(-55.0, 0.0))?;
+ triangulation.insert(VertexType::new(-32.0, -42.0))?;
+ triangulation.insert(VertexType::new(-2.0, -32.0))?;
+ triangulation.insert(VertexType::new(25.0, -32.0))?;
+
+ triangulation.insert(VertexType::new(70.0, 40.0))?;
+ triangulation.insert(VertexType::new(20.0, 65.0))?;
+ triangulation.insert(VertexType::new(-60.0, 50.0))?;
+ triangulation.insert(VertexType::new(-70.0, -15.5))?;
+ triangulation.insert(VertexType::new(-50.0, -60.0))?;
+ triangulation.insert(VertexType::new(-9.0, -70.0))?;
+ triangulation.insert(VertexType::new(42.0, -50.0))?;
+
+ for edge in triangulation.fixed_undirected_edges() {
+ triangulation.undirected_edge_data_mut(edge).color = SketchColor::LIGHT_GRAY;
+ }
+
+ for face in triangulation.fixed_inner_faces() {
+ triangulation.face_data_mut(face).fill = SketchFill::solid(SketchColor::ANTIQUE_WHITE);
+ }
+
+ let query_vertex = VertexType::new(-5.0, -5.0);
+ let mut nns = Vec::new();
+ triangulation
+ .natural_neighbor()
+ .get_weights(query_vertex.position(), &mut nns);
+
+ for (vertex, _) in &nns {
+ triangulation.vertex_data_mut(*vertex).color = Some(SketchColor::CRIMSON);
+ }
+
+ let mut result = convert_triangulation(&triangulation, &Default::default());
+
+ for (index, (vertex, _)) in nns.iter().enumerate() {
+ let vertex = triangulation.vertex(*vertex);
+ result.add(
+ SketchElement::text(format!("{index}"))
+ .position(convert_point(vertex.position()))
+ .font_size(2.5)
+ .horizontal_alignment(HorizontalAlignment::Middle)
+ .dy(0.9),
+ );
+ }
+
+ for edge in triangulation.undirected_voronoi_edges() {
+ let [v0, v1] = edge.vertices();
+
+ if let (Some(v0), Some(v1)) = (v0.position(), v1.position()) {
+ result.add(
+ SketchElement::line(convert_point(v0), convert_point(v1))
+ .stroke_color(SketchColor::SALMON)
+ .stroke_width(0.5),
+ );
+ }
+ }
+
+ let mut circumcenters = Vec::new();
+ for face in triangulation.inner_faces() {
+ let circumcenter = convert_point(face.circumcenter());
+ circumcenters.push(circumcenter);
+
+ result.add(
+ SketchElement::circle(circumcenter, 1.0)
+ .fill(SketchFill::Solid(SketchColor::ROYAL_BLUE)),
+ );
+ }
+
+ if !include_faces {
+ result.add(
+ SketchElement::circle(convert_point(query_vertex.position()), 1.0)
+ .fill(SketchFill::Solid(SketchColor::RED)),
+ );
+ }
+
+ let mut insertion_cell_points = Vec::new();
+ let inserted = triangulation.insert(query_vertex)?;
+ for edge in triangulation
+ .vertex(inserted)
+ .as_voronoi_face()
+ .adjacent_edges()
+ {
+ let context = "Edge was infinite - is insertion position correct?";
+ let from = edge.from().position().context(context)?;
+ let to = edge.to().position().context(context)?;
+ insertion_cell_points.push(convert_point(from));
+
+ result.add(
+ SketchElement::line(convert_point(from), convert_point(to))
+ .stroke_color(SketchColor::ORANGE_RED),
+ );
+ }
+ for cell_point in &insertion_cell_points {
+ result.add(
+ SketchElement::circle(*cell_point, 1.0)
+ .fill(SketchFill::solid(SketchColor::DARK_GREEN)),
+ );
+ }
+
+ if include_faces {
+ let nn3 = convert_point(triangulation.vertex(nns[3].0).position());
+ let nn4 = convert_point(triangulation.vertex(nns[4].0).position());
+ let nn5 = convert_point(triangulation.vertex(nns[5].0).position());
+
+ let last_edge = SketchElement::line(nn3, nn4)
+ .with_arrow_end(ArrowType::FilledArrow)
+ .stroke_color(SketchColor::ROYAL_BLUE)
+ .shift_from(-2.3)
+ .shift_to(-7.0);
+
+ result.add(
+ last_edge
+ .create_adjacent_text("last_edge")
+ .font_size(3.0)
+ .dy(3.5),
+ );
+ result.add(last_edge);
+ let stop_edge = SketchElement::line(nn4, nn5)
+ .with_arrow_end(ArrowType::FilledArrow)
+ .stroke_color(SketchColor::ROYAL_BLUE)
+ .shift_from(-2.3)
+ .shift_to(-7.0);
+ result.add(
+ stop_edge
+ .create_adjacent_text("stop_edge")
+ .font_size(3.0)
+ .dy(3.5),
+ );
+ result.add(stop_edge);
+
+ result.add(
+ SketchElement::text("first")
+ .position(insertion_cell_points[6] + Vector2::new(0.0, 0.0))
+ .dy(3.3)
+ .font_size(2.5),
+ );
+ result.add(
+ SketchElement::text("last")
+ .position(insertion_cell_points[0] + Vector2::new(-5.5, 0.0))
+ .dy(3.0)
+ .font_size(2.5),
+ );
+
+ result.add(
+ SketchElement::text("c2")
+ .position(circumcenters[4] + Vector2::new(-3.5, 0.0))
+ .dy(-1.0)
+ .font_size(2.5),
+ );
+
+ result.add(
+ SketchElement::text("c1")
+ .position(circumcenters[2])
+ .dy(-2.0)
+ .font_size(2.5),
+ );
+
+ result.add(
+ SketchElement::text("c0")
+ .position(circumcenters[3] + Vector2::new(1.5, 0.0))
+ .dy(0.5)
+ .font_size(2.5),
+ );
+
+ let path = SketchElement::path()
+ .move_to(insertion_cell_points[6])
+ .line_to(insertion_cell_points[0])
+ .line_to(circumcenters[4])
+ .line_to(circumcenters[2])
+ .line_to(circumcenters[3])
+ .close()
+ .fill(SketchFill::solid(SketchColor::ORANGE))
+ .opacity(0.75);
+
+ result.add_with_layer(path, SketchLayer::EDGES);
+ }
+
+ result
+ .set_view_box_min(Point2::new(-60.0, -45.0))
+ .set_view_box_max(Point2::new(45.0, 60.0));
+
+ Ok(result)
+}
+
pub fn refinement_scenario(do_refine: bool) -> Sketch {
let mut cdt = create_refinement_cdt();
diff --git a/images/grid_3d.svg b/images/grid_3d.svg
new file mode 100644
index 0000000..cfcfbe6
--- /dev/null
+++ b/images/grid_3d.svg
@@ -0,0 +1,12 @@
+
\ No newline at end of file
diff --git a/images/interpolation_barycentric.img b/images/interpolation_barycentric.img
new file mode 100644
index 0000000..91a7dab
--- /dev/null
+++ b/images/interpolation_barycentric.img
@@ -0,0 +1 @@
+
\ No newline at end of file
diff --git a/images/interpolation_barycentric.jpg b/images/interpolation_barycentric.jpg
new file mode 100644
index 0000000..ea36e85
Binary files /dev/null and b/images/interpolation_barycentric.jpg differ
diff --git a/images/interpolation_nearest_neighbor.img b/images/interpolation_nearest_neighbor.img
new file mode 100644
index 0000000..aa91cc0
--- /dev/null
+++ b/images/interpolation_nearest_neighbor.img
@@ -0,0 +1 @@
+
\ No newline at end of file
diff --git a/images/interpolation_nearest_neighbor.jpg b/images/interpolation_nearest_neighbor.jpg
new file mode 100644
index 0000000..8a935bf
Binary files /dev/null and b/images/interpolation_nearest_neighbor.jpg differ
diff --git a/images/interpolation_nn_c0.img b/images/interpolation_nn_c0.img
new file mode 100644
index 0000000..9412ea4
--- /dev/null
+++ b/images/interpolation_nn_c0.img
@@ -0,0 +1 @@
+
\ No newline at end of file
diff --git a/images/interpolation_nn_c0.jpg b/images/interpolation_nn_c0.jpg
new file mode 100644
index 0000000..1541960
Binary files /dev/null and b/images/interpolation_nn_c0.jpg differ
diff --git a/images/interpolation_nn_c1.img b/images/interpolation_nn_c1.img
new file mode 100644
index 0000000..0c16345
--- /dev/null
+++ b/images/interpolation_nn_c1.img
@@ -0,0 +1 @@
+
\ No newline at end of file
diff --git a/images/interpolation_nn_c1.jpg b/images/interpolation_nn_c1.jpg
new file mode 100644
index 0000000..bb196da
Binary files /dev/null and b/images/interpolation_nn_c1.jpg differ
diff --git a/images/natural_neighbor_insertion_cell.svg b/images/natural_neighbor_insertion_cell.svg
new file mode 100644
index 0000000..caf83c6
--- /dev/null
+++ b/images/natural_neighbor_insertion_cell.svg
@@ -0,0 +1,156 @@
+
\ No newline at end of file
diff --git a/images/natural_neighbor_polygon.svg b/images/natural_neighbor_polygon.svg
new file mode 100644
index 0000000..a31d799
--- /dev/null
+++ b/images/natural_neighbor_polygon.svg
@@ -0,0 +1,185 @@
+
\ No newline at end of file
diff --git a/images/natural_neighbor_scenario.svg b/images/natural_neighbor_scenario.svg
new file mode 100644
index 0000000..3aa4db6
--- /dev/null
+++ b/images/natural_neighbor_scenario.svg
@@ -0,0 +1,140 @@
+
\ No newline at end of file
diff --git a/images/preview.html b/images/preview.html
index 0097e83..a7da86c 100644
--- a/images/preview.html
+++ b/images/preview.html
@@ -23,7 +23,7 @@
}
-
+