WIP

node-graph/gvector-nodes/Cargo.toml

@@ -0,0 +1,38 @@
+[package]
+name = "graphene-vector-nodes"
+version = "0.1.0"
+edition = "2024"
+description = "graphene vector nodes"
+authors = ["Graphite Authors <contact@graphite.rs>"]
+license = "MIT OR Apache-2.0"
+
+[features]
+default = ["serde"]
+serde = [
+	"dep:serde",
+	"bezier-rs/serde",
+]
+
+[dependencies]
+# Local dependencies
+dyn-any = { workspace = true }
+bezier-rs = { workspace = true }
+graphene-core = { workspace = true }
+graphene-vector = { workspace = true }
+node-macro = { workspace = true }
+
+# Workspace dependencies
+kurbo = { workspace = true }
+glam = { workspace = true }
+specta = { workspace = true }
+log = { workspace = true }
+rustc-hash = { workspace = true }
+petgraph = { workspace = true }
+rand = { workspace = true }
+
+# Optional workspace dependencies
+serde = { workspace = true, optional = true, features = ["derive"] }
+
+[dev-dependencies]
+# Workspace dependencies
+tokio = { workspace = true }

node-graph/gvector-nodes/src/algorithms/bezpath_algorithms.rs

@@ -0,0 +1,316 @@
+use super::poisson_disk::poisson_disk_sample;
+use crate::misc::{PointSpacingType, dvec2_to_point};
+use glam::DVec2;
+use kurbo::{BezPath, DEFAULT_ACCURACY, Line, ParamCurve, ParamCurveDeriv, PathEl, PathSeg, Point, Rect, Shape};
+
+/// Splits the [`BezPath`] at `t` value which lie in the range of [0, 1].
+/// Returns [`None`] if the given [`BezPath`] has no segments or `t` is within f64::EPSILON of 0 or 1.
+pub fn split_bezpath(bezpath: &BezPath, t: f64, euclidian: bool) -> Option<(BezPath, BezPath)> {
+	if t <= f64::EPSILON || (1. - t) <= f64::EPSILON || bezpath.segments().count() == 0 {
+		return None;
+	}
+
+	// Get the segment which lies at the split.
+	let (segment_index, t) = t_value_to_parametric(bezpath, t, euclidian, None);
+	let segment = bezpath.get_seg(segment_index + 1).unwrap();
+
+	// Divide the segment.
+	let first_segment = segment.subsegment(0.0..t);
+	let second_segment = segment.subsegment(t..1.);
+
+	let mut first_bezpath = BezPath::new();
+	let mut second_bezpath = BezPath::new();
+
+	// Append the segments up to the subdividing segment from original bezpath to first bezpath.
+	for segment in bezpath.segments().take(segment_index) {
+		if first_bezpath.elements().is_empty() {
+			first_bezpath.move_to(segment.start());
+		}
+		first_bezpath.push(segment.as_path_el());
+	}
+
+	// Append the first segment of the subdivided segment.
+	if first_bezpath.elements().is_empty() {
+		first_bezpath.move_to(first_segment.start());
+	}
+	first_bezpath.push(first_segment.as_path_el());
+
+	// Append the second segment of the subdivided segment in the second bezpath.
+	if second_bezpath.elements().is_empty() {
+		second_bezpath.move_to(second_segment.start());
+	}
+	second_bezpath.push(second_segment.as_path_el());
+
+	// Append the segments after the subdividing segment from original bezpath to second bezpath.
+	for segment in bezpath.segments().skip(segment_index + 1) {
+		if second_bezpath.elements().is_empty() {
+			second_bezpath.move_to(segment.start());
+		}
+		second_bezpath.push(segment.as_path_el());
+	}
+
+	Some((first_bezpath, second_bezpath))
+}
+
+pub fn position_on_bezpath(bezpath: &BezPath, t: f64, euclidian: bool, segments_length: Option<&[f64]>) -> Point {
+	let (segment_index, t) = t_value_to_parametric(bezpath, t, euclidian, segments_length);
+	bezpath.get_seg(segment_index + 1).unwrap().eval(t)
+}
+
+pub fn tangent_on_bezpath(bezpath: &BezPath, t: f64, euclidian: bool, segments_length: Option<&[f64]>) -> Point {
+	let (segment_index, t) = t_value_to_parametric(bezpath, t, euclidian, segments_length);
+	let segment = bezpath.get_seg(segment_index + 1).unwrap();
+	match segment {
+		PathSeg::Line(line) => line.deriv().eval(t),
+		PathSeg::Quad(quad_bez) => quad_bez.deriv().eval(t),
+		PathSeg::Cubic(cubic_bez) => cubic_bez.deriv().eval(t),
+	}
+}
+
+pub fn sample_polyline_on_bezpath(
+	bezpath: BezPath,
+	point_spacing_type: PointSpacingType,
+	amount: f64,
+	start_offset: f64,
+	stop_offset: f64,
+	adaptive_spacing: bool,
+	segments_length: &[f64],
+) -> Option<BezPath> {
+	let mut sample_bezpath = BezPath::new();
+
+	let was_closed = matches!(bezpath.elements().last(), Some(PathEl::ClosePath));
+
+	// Calculate the total length of the collected segments.
+	let total_length: f64 = segments_length.iter().sum();
+
+	// Adjust the usable length by subtracting start and stop offsets.
+	let mut used_length = total_length - start_offset - stop_offset;
+
+	// Sanity check that the usable length is positive.
+	if used_length <= 0. {
+		return None;
+	}
+
+	const SAFETY_MAX_COUNT: f64 = 10_000. - 1.;
+
+	// Determine the number of points to generate along the path.
+	let sample_count = match point_spacing_type {
+		PointSpacingType::Separation => {
+			let spacing = amount.min(used_length - f64::EPSILON);
+
+			if adaptive_spacing {
+				// Calculate point count to evenly distribute points while covering the entire path.
+				// With adaptive spacing, we widen or narrow the points as necessary to ensure the last point is always at the end of the path.
+				(used_length / spacing).round().min(SAFETY_MAX_COUNT)
+			} else {
+				// Calculate point count based on exact spacing, which may not cover the entire path.
+				// Without adaptive spacing, we just evenly space the points at the exact specified spacing, usually falling short before the end of the path.
+				let count = (used_length / spacing + f64::EPSILON).floor().min(SAFETY_MAX_COUNT);
+				if count != SAFETY_MAX_COUNT {
+					used_length -= used_length % spacing;
+				}
+				count
+			}
+		}
+		PointSpacingType::Quantity => (amount - 1.).floor().clamp(1., SAFETY_MAX_COUNT),
+	};
+
+	// Skip if there are no points to generate.
+	if sample_count < 1. {
+		return None;
+	}
+
+	// Decide how many loop-iterations: if closed, skip the last duplicate point
+	let sample_count_usize = sample_count as usize;
+	let max_i = if was_closed { sample_count_usize } else { sample_count_usize + 1 };
+
+	// Generate points along the path based on calculated intervals.
+	let mut length_up_to_previous_segment = 0.;
+	let mut next_segment_index = 0;
+
+	for count in 0..max_i {
+		let fraction = count as f64 / sample_count;
+		let length_up_to_next_sample_point = fraction * used_length + start_offset;
+		let mut next_length = length_up_to_next_sample_point - length_up_to_previous_segment;
+		let mut next_segment_length = segments_length[next_segment_index];
+
+		// Keep moving to the next segment while the length up to the next sample point is greater than the length up to the current segment.
+		while next_length > next_segment_length {
+			if next_segment_index == segments_length.len() - 1 {
+				break;
+			}
+			length_up_to_previous_segment += next_segment_length;
+			next_length = length_up_to_next_sample_point - length_up_to_previous_segment;
+			next_segment_index += 1;
+			next_segment_length = segments_length[next_segment_index];
+		}
+
+		let t = (next_length / next_segment_length).clamp(0., 1.);
+
+		let segment = bezpath.get_seg(next_segment_index + 1).unwrap();
+		let t = eval_pathseg_euclidean(segment, t, DEFAULT_ACCURACY);
+		let point = segment.eval(t);
+
+		if sample_bezpath.elements().is_empty() {
+			sample_bezpath.move_to(point)
+		} else {
+			sample_bezpath.line_to(point)
+		}
+	}
+
+	if was_closed {
+		sample_bezpath.close_path();
+	}
+
+	Some(sample_bezpath)
+}
+
+pub fn t_value_to_parametric(bezpath: &BezPath, t: f64, euclidian: bool, segments_length: Option<&[f64]>) -> (usize, f64) {
+	if euclidian {
+		let (segment_index, t) = bezpath_t_value_to_parametric(bezpath, BezPathTValue::GlobalEuclidean(t), segments_length);
+		let segment = bezpath.get_seg(segment_index + 1).unwrap();
+		return (segment_index, eval_pathseg_euclidean(segment, t, DEFAULT_ACCURACY));
+	}
+	bezpath_t_value_to_parametric(bezpath, BezPathTValue::GlobalParametric(t), segments_length)
+}
+
+/// Finds the t value of point on the given path segment i.e fractional distance along the segment's total length.
+/// It uses a binary search to find the value `t` such that the ratio `length_up_to_t / total_length` approximates the input `distance`.
+pub fn eval_pathseg_euclidean(path_segment: PathSeg, distance: f64, accuracy: f64) -> f64 {
+	let mut low_t = 0.;
+	let mut mid_t = 0.5;
+	let mut high_t = 1.;
+
+	let total_length = path_segment.perimeter(accuracy);
+
+	if !total_length.is_finite() || total_length <= f64::EPSILON {
+		return 0.;
+	}
+
+	let distance = distance.clamp(0., 1.);
+
+	while high_t - low_t > accuracy {
+		let current_length = path_segment.subsegment(0.0..mid_t).perimeter(accuracy);
+		let current_distance = current_length / total_length;
+
+		if current_distance > distance {
+			high_t = mid_t;
+		} else {
+			low_t = mid_t;
+		}
+		mid_t = (high_t + low_t) / 2.;
+	}
+
+	mid_t
+}
+
+/// Converts from a bezpath (composed of multiple segments) to a point along a certain segment represented.
+/// The returned tuple represents the segment index and the `t` value along that segment.
+/// Both the input global `t` value and the output `t` value are in euclidean space, meaning there is a constant rate of change along the arc length.
+fn global_euclidean_to_local_euclidean(bezpath: &BezPath, global_t: f64, lengths: &[f64], total_length: f64) -> (usize, f64) {
+	let mut accumulator = 0.;
+	for (index, length) in lengths.iter().enumerate() {
+		let length_ratio = length / total_length;
+		if (index == 0 || accumulator <= global_t) && global_t <= accumulator + length_ratio {
+			return (index, ((global_t - accumulator) / length_ratio).clamp(0., 1.));
+		}
+		accumulator += length_ratio;
+	}
+	(bezpath.segments().count() - 1, 1.)
+}
+
+enum BezPathTValue {
+	GlobalEuclidean(f64),
+	GlobalParametric(f64),
+}
+
+/// Convert a [BezPathTValue] to a parametric `(segment_index, t)` tuple.
+/// - Asserts that `t` values contained within the `SubpathTValue` argument lie in the range [0, 1].
+fn bezpath_t_value_to_parametric(bezpath: &BezPath, t: BezPathTValue, precomputed_segments_length: Option<&[f64]>) -> (usize, f64) {
+	let segment_count = bezpath.segments().count();
+	assert!(segment_count >= 1);
+
+	match t {
+		BezPathTValue::GlobalEuclidean(t) => {
+			let computed_segments_length;
+
+			let segments_length = if let Some(segments_length) = precomputed_segments_length {
+				segments_length
+			} else {
+				computed_segments_length = bezpath.segments().map(|segment| segment.perimeter(DEFAULT_ACCURACY)).collect::<Vec<f64>>();
+				computed_segments_length.as_slice()
+			};
+
+			let total_length = segments_length.iter().sum();
+
+			global_euclidean_to_local_euclidean(bezpath, t, segments_length, total_length)
+		}
+		BezPathTValue::GlobalParametric(global_t) => {
+			assert!((0.0..=1.).contains(&global_t));
+
+			if global_t == 1. {
+				return (segment_count - 1, 1.);
+			}
+
+			let scaled_t = global_t * segment_count as f64;
+			let segment_index = scaled_t.floor() as usize;
+			let t = scaled_t - segment_index as f64;
+
+			(segment_index, t)
+		}
+	}
+}
+
+/// Randomly places points across the filled surface of this subpath (which is assumed to be closed).
+/// The `separation_disk_diameter` determines the minimum distance between all points from one another.
+/// Conceptually, this works by "throwing a dart" at the subpath's bounding box and keeping the dart only if:
+/// - It's inside the shape
+/// - It's not closer than `separation_disk_diameter` to any other point from a previous accepted dart throw
+///
+/// This repeats until accepted darts fill all possible areas between one another.
+///
+/// While the conceptual process described above asymptotically slows down and is never guaranteed to produce a maximal set in finite time,
+/// this is implemented with an algorithm that produces a maximal set in O(n) time. The slowest part is actually checking if points are inside the subpath shape.
+pub fn poisson_disk_points(bezpath_index: usize, bezpaths: &[(BezPath, Rect)], separation_disk_diameter: f64, rng: impl FnMut() -> f64) -> Vec<DVec2> {
+	let (this_bezpath, this_bbox) = bezpaths[bezpath_index].clone();
+
+	if this_bezpath.elements().is_empty() {
+		return Vec::new();
+	}
+
+	let point_in_shape_checker = |point: DVec2| {
+		// Check against all paths the point is contained in to compute the correct winding number
+		let mut number = 0;
+
+		for (i, (shape, bbox)) in bezpaths.iter().enumerate() {
+			if bbox.x0 > point.x || bbox.y0 > point.y || bbox.x1 < point.x || bbox.y1 < point.y {
+				continue;
+			}
+
+			let winding = shape.winding(dvec2_to_point(point));
+			if winding == 0 && i == bezpath_index {
+				return false;
+			}
+			number += winding;
+		}
+
+		// Non-zero fill rule
+		number != 0
+	};
+
+	let line_intersect_shape_checker = |p0: (f64, f64), p1: (f64, f64)| {
+		for segment in this_bezpath.segments() {
+			if !segment.intersect_line(Line::new(p0, p1)).is_empty() {
+				return true;
+			}
+		}
+
+		false
+	};
+
+	let offset = DVec2::new(this_bbox.x0, this_bbox.y0);
+	let width = this_bbox.width();
+	let height = this_bbox.height();
+
+	poisson_disk_sample(offset, width, height, separation_disk_diameter, point_in_shape_checker, line_intersect_shape_checker, rng)
+}

node-graph/gvector-nodes/src/algorithms/merge_by_distance.rs

@@ -0,0 +1,214 @@
+use glam::{DAffine2, DVec2};
+use graphene_vector::{PointDomain, PointId, SegmentDomain, VectorData, VectorDataIndex};
+use petgraph::prelude::UnGraphMap;
+use rustc_hash::FxHashSet;
+
+pub trait MergeByDistanceExt {
+	/// Collapse all points with edges shorter than the specified distance
+	fn merge_by_distance_topological(&mut self, distance: f64);
+	fn merge_by_distance_spatial(&mut self, transform: DAffine2, distance: f64);
+}
+
+impl MergeByDistanceExt for VectorData {
+	fn merge_by_distance_topological(&mut self, distance: f64) {
+		// Treat self as an undirected graph
+		let indices = VectorDataIndex::build_from(self);
+
+		// TODO: We lose information on the winding order by using an undirected graph. Switch to a directed graph and fix the algorithm to handle that.
+		// Graph containing only short edges, referencing the data graph
+		let mut short_edges = UnGraphMap::new();
+
+		for segment_id in self.segment_ids().iter().copied() {
+			let length = indices.segment_chord_length(segment_id);
+			if length < distance {
+				let [start, end] = indices.segment_ends(segment_id);
+				let start = indices.point_graph.node_weight(start).unwrap().id;
+				let end = indices.point_graph.node_weight(end).unwrap().id;
+
+				short_edges.add_node(start);
+				short_edges.add_node(end);
+				short_edges.add_edge(start, end, segment_id);
+			}
+		}
+
+		// Group connected segments to collapse them into a single point
+		// TODO: there are a few possible algorithms for this - perhaps test empirically to find fastest
+		let collapse: Vec<FxHashSet<PointId>> = petgraph::algo::tarjan_scc(&short_edges).into_iter().map(|connected| connected.into_iter().collect()).collect();
+		let average_position = collapse
+			.iter()
+			.map(|collapse_set| {
+				let sum: DVec2 = collapse_set.iter().map(|&id| indices.point_position(id, self)).sum();
+				sum / collapse_set.len() as f64
+			})
+			.collect::<Vec<_>>();
+
+		// Collect points and segments to delete at the end to avoid invalidating indices
+		let mut points_to_delete = FxHashSet::default();
+		let mut segments_to_delete = FxHashSet::default();
+		for (mut collapse_set, average_pos) in collapse.into_iter().zip(average_position.into_iter()) {
+			// Remove any segments where both endpoints are in the collapse set
+			segments_to_delete.extend(self.segment_domain.iter().filter_map(|(id, start_offset, end_offset, _)| {
+				let start = self.point_domain.ids()[start_offset];
+				let end = self.point_domain.ids()[end_offset];
+				if collapse_set.contains(&start) && collapse_set.contains(&end) { Some(id) } else { None }
+			}));
+
+			// Delete all points but the first, set its position to the average, and update segments
+			let first_id = collapse_set.iter().copied().next().unwrap();
+			collapse_set.remove(&first_id);
+			let first_offset = indices.point_to_offset[&first_id];
+
+			// Look for segments with endpoints in `collapse_set` and replace them with the point we are collapsing to
+			for (_, start_offset, end_offset, handles) in self.segment_domain.iter_mut() {
+				let start_id = self.point_domain.ids()[*start_offset];
+				let end_id = self.point_domain.ids()[*end_offset];
+
+				// Update Bezier handles for moved points
+				if start_id == first_id {
+					let point_position = self.point_domain.position[*start_offset];
+					handles.move_start(average_pos - point_position);
+				}
+				if end_id == first_id {
+					let point_position = self.point_domain.position[*end_offset];
+					handles.move_end(average_pos - point_position);
+				}
+
+				// Replace removed points with the collapsed point
+				if collapse_set.contains(&start_id) {
+					let point_position = self.point_domain.position[*start_offset];
+					*start_offset = first_offset;
+					handles.move_start(average_pos - point_position);
+				}
+				if collapse_set.contains(&end_id) {
+					let point_position = self.point_domain.position[*end_offset];
+					*end_offset = first_offset;
+					handles.move_end(average_pos - point_position);
+				}
+			}
+
+			// Update the position of the collapsed point
+			self.point_domain.position[first_offset] = average_pos;
+
+			points_to_delete.extend(collapse_set)
+		}
+
+		// Remove faces whose start or end segments are removed
+		// TODO: Adjust faces and only delete if all (or all but one) segments are removed
+		self.region_domain
+			.retain_with_region(|_, segment_range| segments_to_delete.contains(segment_range.start()) || segments_to_delete.contains(segment_range.end()));
+		self.segment_domain.retain(|id| !segments_to_delete.contains(id), usize::MAX);
+		self.point_domain.retain(&mut self.segment_domain, |id| !points_to_delete.contains(id));
+	}
+
+	fn merge_by_distance_spatial(&mut self, transform: DAffine2, distance: f64) {
+		let point_count = self.point_domain.positions().len();
+
+		// Find min x and y for grid cell normalization
+		let mut min_x = f64::MAX;
+		let mut min_y = f64::MAX;
+
+		// Calculate mins without collecting all positions
+		for &pos in self.point_domain.positions() {
+			let transformed_pos = transform.transform_point2(pos);
+			min_x = min_x.min(transformed_pos.x);
+			min_y = min_y.min(transformed_pos.y);
+		}
+
+		// Create a spatial grid with cell size of 'distance'
+		use std::collections::HashMap;
+		let mut grid: HashMap<(i32, i32), Vec<usize>> = HashMap::new();
+
+		// Add points to grid cells without collecting all positions first
+		for i in 0..point_count {
+			let pos = transform.transform_point2(self.point_domain.positions()[i]);
+			let grid_x = ((pos.x - min_x) / distance).floor() as i32;
+			let grid_y = ((pos.y - min_y) / distance).floor() as i32;
+
+			grid.entry((grid_x, grid_y)).or_default().push(i);
+		}
+
+		// Create point index mapping for merged points
+		let mut point_index_map = vec![None; point_count];
+		let mut merged_positions = Vec::new();
+		let mut merged_indices = Vec::new();
+
+		// Process each point
+		for i in 0..point_count {
+			// Skip points that have already been processed
+			if point_index_map[i].is_some() {
+				continue;
+			}
+
+			let pos_i = transform.transform_point2(self.point_domain.positions()[i]);
+			let grid_x = ((pos_i.x - min_x) / distance).floor() as i32;
+			let grid_y = ((pos_i.y - min_y) / distance).floor() as i32;
+
+			let mut group = vec![i];
+
+			// Check only neighboring cells (3x3 grid around current cell)
+			for dx in -1..=1 {
+				for dy in -1..=1 {
+					let neighbor_cell = (grid_x + dx, grid_y + dy);
+
+					if let Some(indices) = grid.get(&neighbor_cell) {
+						for &j in indices {
+							if j > i && point_index_map[j].is_none() {
+								let pos_j = transform.transform_point2(self.point_domain.positions()[j]);
+								if pos_i.distance(pos_j) <= distance {
+									group.push(j);
+								}
+							}
+						}
+					}
+				}
+			}
+
+			// Create merged point - calculate positions as needed
+			let merged_position = group
+				.iter()
+				.map(|&idx| transform.transform_point2(self.point_domain.positions()[idx]))
+				.fold(DVec2::ZERO, |sum, pos| sum + pos)
+				/ group.len() as f64;
+
+			let merged_position = transform.inverse().transform_point2(merged_position);
+			let merged_index = merged_positions.len();
+
+			merged_positions.push(merged_position);
+			merged_indices.push(self.point_domain.ids()[group[0]]);
+
+			// Update mapping for all points in the group
+			for &idx in &group {
+				point_index_map[idx] = Some(merged_index);
+			}
+		}
+
+		// Create new point domain with merged points
+		let mut new_point_domain = PointDomain::new();
+		for (idx, pos) in merged_indices.into_iter().zip(merged_positions) {
+			new_point_domain.push(idx, pos);
+		}
+
+		// Update segment domain
+		let mut new_segment_domain = SegmentDomain::new();
+		for segment_idx in 0..self.segment_domain.ids().len() {
+			let id = self.segment_domain.ids()[segment_idx];
+			let start = self.segment_domain.start_point()[segment_idx];
+			let end = self.segment_domain.end_point()[segment_idx];
+			let handles = self.segment_domain.handles()[segment_idx];
+			let stroke = self.segment_domain.stroke()[segment_idx];
+
+			// Get new indices for start and end points
+			let new_start = point_index_map[start].unwrap();
+			let new_end = point_index_map[end].unwrap();
+
+			// Skip segments where start and end points were merged
+			if new_start != new_end {
+				new_segment_domain.push(id, new_start, new_end, handles, stroke);
+			}
+		}
+
+		// Create new vector data
+		self.point_domain = new_point_domain;
+		self.segment_domain = new_segment_domain;
+	}
+}

node-graph/gvector-nodes/src/algorithms/mod.rs

@@ -0,0 +1,5 @@
+pub mod bezpath_algorithms;
+pub mod merge_by_distance;
+pub mod offset_subpath;
+pub mod poisson_disk;
+pub mod spline;

node-graph/gvector-nodes/src/algorithms/offset_subpath.rs

@@ -0,0 +1,173 @@
+use bezier_rs::{Bezier, BezierHandles, Join, Subpath, TValue};
+use graphene_vector::PointId;
+
+/// Value to control smoothness and mathematical accuracy to offset a cubic Bezier.
+const CUBIC_REGULARIZATION_ACCURACY: f64 = 0.5;
+/// Accuracy of fitting offset curve to Bezier paths.
+const CUBIC_TO_BEZPATH_ACCURACY: f64 = 1e-3;
+/// Constant used to determine if `f64`s are equivalent.
+pub const MAX_ABSOLUTE_DIFFERENCE: f64 = 1e-3;
+
+fn segment_to_bezier(seg: kurbo::PathSeg) -> Bezier {
+	match seg {
+		kurbo::PathSeg::Line(line) => Bezier::from_linear_coordinates(line.p0.x, line.p0.y, line.p1.x, line.p1.y),
+		kurbo::PathSeg::Quad(quad_bez) => Bezier::from_quadratic_coordinates(quad_bez.p0.x, quad_bez.p0.y, quad_bez.p1.x, quad_bez.p1.y, quad_bez.p1.x, quad_bez.p1.y),
+		kurbo::PathSeg::Cubic(cubic_bez) => Bezier::from_cubic_coordinates(
+			cubic_bez.p0.x,
+			cubic_bez.p0.y,
+			cubic_bez.p1.x,
+			cubic_bez.p1.y,
+			cubic_bez.p2.x,
+			cubic_bez.p2.y,
+			cubic_bez.p3.x,
+			cubic_bez.p3.y,
+		),
+	}
+}
+
+// TODO: Replace the implementation to use only Kurbo API.
+/// Reduces the segments of the subpath into simple subcurves, then offset each subcurve a set `distance` away.
+/// The intersections of segments of the subpath are joined using the method specified by the `join` argument.
+pub fn offset_subpath(subpath: &Subpath<PointId>, distance: f64, join: Join) -> Subpath<PointId> {
+	// An offset at a distance 0 from the curve is simply the same curve.
+	// An offset of a single point is not defined.
+	if distance == 0. || subpath.len() <= 1 || subpath.len_segments() < 1 {
+		return subpath.clone();
+	}
+
+	let mut subpaths = subpath
+			.iter()
+			.filter(|bezier| !bezier.is_point())
+			.map(|bezier| bezier.to_cubic())
+			.map(|cubic| {
+				let Bezier { start, end, handles } = cubic;
+				let BezierHandles::Cubic { handle_start, handle_end } = handles else { unreachable!()};
+
+				let cubic_bez = kurbo::CubicBez::new((start.x, start.y), (handle_start.x, handle_start.y), (handle_end.x, handle_end.y), (end.x, end.y));
+				let cubic_offset = kurbo::offset::CubicOffset::new_regularized(cubic_bez, distance, CUBIC_REGULARIZATION_ACCURACY);
+				let offset_bezpath = kurbo::fit_to_bezpath(&cubic_offset, CUBIC_TO_BEZPATH_ACCURACY);
+
+				let beziers = offset_bezpath.segments().fold(Vec::new(), |mut acc, seg| {
+					acc.push(segment_to_bezier(seg));
+					acc
+				});
+
+				Subpath::from_beziers(&beziers, false)
+			})
+			.filter(|subpath| subpath.len() >= 2) // In some cases the reduced and scaled bézier is marked by is_point (so the subpath is empty).
+			.collect::<Vec<Subpath<PointId>>>();
+
+	let mut drop_common_point = vec![true; subpath.len()];
+
+	// Clip or join consecutive Subpaths
+	for i in 0..subpaths.len() - 1 {
+		let j = i + 1;
+		let subpath1 = &subpaths[i];
+		let subpath2 = &subpaths[j];
+
+		let last_segment = subpath1.get_segment(subpath1.len_segments() - 1).unwrap();
+		let first_segment = subpath2.get_segment(0).unwrap();
+
+		// If the anchors are approximately equal, there is no need to clip / join the segments
+		if last_segment.end().abs_diff_eq(first_segment.start(), MAX_ABSOLUTE_DIFFERENCE) {
+			continue;
+		}
+
+		// Calculate the angle formed between two consecutive Subpaths
+		let out_tangent = subpath.get_segment(i).unwrap().tangent(TValue::Parametric(1.));
+		let in_tangent = subpath.get_segment(j).unwrap().tangent(TValue::Parametric(0.));
+		let angle = out_tangent.angle_to(in_tangent);
+
+		// The angle is concave. The Subpath overlap and must be clipped
+		let mut apply_join = true;
+		if (angle > 0. && distance > 0.) || (angle < 0. && distance < 0.) {
+			// If the distance is large enough, there may still be no intersections. Also, if the angle is close enough to zero,
+			// subpath intersections may find no intersections. In this case, the points are likely close enough that we can approximate
+			// the points as being on top of one another.
+			if let Some((clipped_subpath1, clipped_subpath2)) = Subpath::clip_simple_subpaths(subpath1, subpath2) {
+				subpaths[i] = clipped_subpath1;
+				subpaths[j] = clipped_subpath2;
+				apply_join = false;
+			}
+		}
+		// The angle is convex. The Subpath must be joined using the specified join type
+		if apply_join {
+			drop_common_point[j] = false;
+			match join {
+				Join::Bevel => {}
+				Join::Miter(miter_limit) => {
+					let miter_manipulator_group = subpaths[i].miter_line_join(&subpaths[j], miter_limit);
+					if let Some(miter_manipulator_group) = miter_manipulator_group {
+						subpaths[i].manipulator_groups_mut().push(miter_manipulator_group);
+					}
+				}
+				Join::Round => {
+					let (out_handle, round_point, in_handle) = subpaths[i].round_line_join(&subpaths[j], subpath.manipulator_groups()[j].anchor);
+					let last_index = subpaths[i].manipulator_groups().len() - 1;
+					subpaths[i].manipulator_groups_mut()[last_index].out_handle = Some(out_handle);
+					subpaths[i].manipulator_groups_mut().push(round_point);
+					subpaths[j].manipulator_groups_mut()[0].in_handle = Some(in_handle);
+				}
+			}
+		}
+	}
+
+	// Clip any overlap in the last segment
+	if subpath.closed {
+		let out_tangent = subpath.get_segment(subpath.len_segments() - 1).unwrap().tangent(TValue::Parametric(1.));
+		let in_tangent = subpath.get_segment(0).unwrap().tangent(TValue::Parametric(0.));
+		let angle = out_tangent.angle_to(in_tangent);
+
+		let mut apply_join = true;
+		if (angle > 0. && distance > 0.) || (angle < 0. && distance < 0.) {
+			if let Some((clipped_subpath1, clipped_subpath2)) = Subpath::clip_simple_subpaths(&subpaths[subpaths.len() - 1], &subpaths[0]) {
+				// Merge the clipped subpaths
+				let last_index = subpaths.len() - 1;
+				subpaths[last_index] = clipped_subpath1;
+				subpaths[0] = clipped_subpath2;
+				apply_join = false;
+			}
+		}
+		if apply_join {
+			drop_common_point[0] = false;
+			match join {
+				Join::Bevel => {}
+				Join::Miter(miter_limit) => {
+					let last_subpath_index = subpaths.len() - 1;
+					let miter_manipulator_group = subpaths[last_subpath_index].miter_line_join(&subpaths[0], miter_limit);
+					if let Some(miter_manipulator_group) = miter_manipulator_group {
+						subpaths[last_subpath_index].manipulator_groups_mut().push(miter_manipulator_group);
+					}
+				}
+				Join::Round => {
+					let last_subpath_index = subpaths.len() - 1;
+					let (out_handle, round_point, in_handle) = subpaths[last_subpath_index].round_line_join(&subpaths[0], subpath.manipulator_groups()[0].anchor);
+					let last_index = subpaths[last_subpath_index].manipulator_groups().len() - 1;
+					subpaths[last_subpath_index].manipulator_groups_mut()[last_index].out_handle = Some(out_handle);
+					subpaths[last_subpath_index].manipulator_groups_mut().push(round_point);
+					subpaths[0].manipulator_groups_mut()[0].in_handle = Some(in_handle);
+				}
+			}
+		}
+	}
+
+	// Merge the subpaths. Drop points which overlap with one another.
+	let mut manipulator_groups = subpaths[0].manipulator_groups().to_vec();
+	for i in 1..subpaths.len() {
+		if drop_common_point[i] {
+			let last_group = manipulator_groups.pop().unwrap();
+			let mut manipulators_copy = subpaths[i].manipulator_groups().to_vec();
+			manipulators_copy[0].in_handle = last_group.in_handle;
+
+			manipulator_groups.append(&mut manipulators_copy);
+		} else {
+			manipulator_groups.append(&mut subpaths[i].manipulator_groups().to_vec());
+		}
+	}
+	if subpath.closed && drop_common_point[0] {
+		let last_group = manipulator_groups.pop().unwrap();
+		manipulator_groups[0].in_handle = last_group.in_handle;
+	}
+
+	Subpath::new(manipulator_groups, subpath.closed)
+}

node-graph/gvector-nodes/src/algorithms/poisson_disk.rs

@@ -0,0 +1,422 @@
+use glam::DVec2;
+use std::collections::HashMap;
+use std::f64;
+
+const DEEPEST_SUBDIVISION_LEVEL_BEFORE_DISCARDING: usize = 8;
+
+/// Fast (O(n) with respect to time and memory) algorithm for generating a maximal set of points using Poisson-disk sampling.
+/// Based on the paper:
+/// "Poisson Disk Point Sets by Hierarchical Dart Throwing"
+/// <https://scholarsarchive.byu.edu/facpub/237/>
+pub fn poisson_disk_sample(
+	offset: DVec2,
+	width: f64,
+	height: f64,
+	diameter: f64,
+	point_in_shape_checker: impl Fn(DVec2) -> bool,
+	line_intersect_shape_checker: impl Fn((f64, f64), (f64, f64)) -> bool,
+	rng: impl FnMut() -> f64,
+) -> Vec<DVec2> {
+	let mut rng = rng;
+	let diameter_squared = diameter.powi(2);
+
+	// Initialize a place to store the generated points within a spatial acceleration structure
+	let mut points_grid = AccelerationGrid::new(width, height, diameter);
+
+	// Pick a grid size for the base-level domain that's as large as possible, while also:
+	// - Dividing into an integer number of cells across the dartboard domain, to avoid wastefully throwing darts beyond the width and height of the dartboard domain
+	// - Being fully covered by the radius around a dart thrown anywhere in its area, where the worst-case is a corner which has a distance of sqrt(2) to the opposite corner
+	let greater_dimension = width.max(height);
+	let base_level_grid_size = greater_dimension / (greater_dimension * f64::consts::SQRT_2 / (diameter / 2.)).ceil();
+
+	// Initialize the problem by including all base-level squares in the active list since they're all part of the yet-to-be-targetted dartboard domain
+	let base_level = ActiveListLevel::new_filled(base_level_grid_size, offset, width, height, &point_in_shape_checker, &line_intersect_shape_checker);
+	// In the future, if necessary, this could be turned into a fixed-length array with worst-case length `f64::MANTISSA_DIGITS`
+	let mut active_list_levels = vec![base_level];
+
+	// Loop until all active squares have been processed, meaning all of the dartboard domain has been checked
+	while active_list_levels.iter().any(|active_list| active_list.not_empty()) {
+		// Randomly pick a square in the dartboard domain, with probability proportional to its area
+		let (active_square_level, active_square_index_in_level) = target_active_square(&active_list_levels, &mut rng);
+
+		// The level contains the list of all active squares at this target square's subdivision depth
+		let level = &mut active_list_levels[active_square_level];
+
+		// Take the targetted active square out of the list and get its size
+		let active_square = level.take_square(active_square_index_in_level);
+		let active_square_size = level.square_size();
+
+		// Skip this target square if it's within range of any current points, since more nearby points could have been added after this square was included in the active list
+		if !square_not_covered_by_poisson_points(active_square.top_left_corner(), active_square_size / 2., diameter_squared, &points_grid) {
+			continue;
+		}
+
+		// Throw a dart by picking a random point within this target square
+		let point = {
+			let active_top_left_corner = active_square.top_left_corner();
+			let x = active_top_left_corner.x + rng() * active_square_size;
+			let y = active_top_left_corner.y + rng() * active_square_size;
+			(x, y).into()
+		};
+
+		// If the dart hit a valid spot, save that point (we're now permanently done with this target square's region)
+		if point_not_covered_by_poisson_points(point, diameter_squared, &points_grid) {
+			// Silently reject the point if it lies outside the shape
+			if active_square.fully_in_shape() || point_in_shape_checker(point + offset) {
+				points_grid.insert(point);
+			}
+		}
+		// Otherwise, subdivide this target square and add valid sub-squares back to the active list for later targetting
+		else {
+			// Discard any targetable domain smaller than this limited number of subdivision levels since it's too small to matter
+			let next_level_deeper_level = active_square_level + 1;
+			if next_level_deeper_level > DEEPEST_SUBDIVISION_LEVEL_BEFORE_DISCARDING {
+				continue;
+			}
+
+			// If necessary for the following step, add another layer of depth to store squares at the next subdivision level
+			if active_list_levels.len() <= next_level_deeper_level {
+				active_list_levels.push(ActiveListLevel::new(active_square_size / 2.))
+			}
+
+			// Get the list of active squares at the level of depth beneath this target square's level
+			let next_level_deeper = &mut active_list_levels[next_level_deeper_level];
+
+			// Subdivide this target square into four sub-squares; running out of numerical precision will make this terminate at very small scales
+			let subdivided_size = active_square_size / 2.;
+			let active_top_left_corner = active_square.top_left_corner();
+			let subdivided = [
+				active_top_left_corner + DVec2::new(0., 0.),
+				active_top_left_corner + DVec2::new(subdivided_size, 0.),
+				active_top_left_corner + DVec2::new(0., subdivided_size),
+				active_top_left_corner + DVec2::new(subdivided_size, subdivided_size),
+			];
+
+			// Add the sub-squares which aren't within the radius of a nearby point to the sub-level's active list
+			let half_subdivided_size = subdivided_size / 2.;
+			let new_sub_squares = subdivided.into_iter().filter_map(|sub_square| {
+				// Any sub-squares within the radius of a nearby point are filtered out
+				if !square_not_covered_by_poisson_points(sub_square, half_subdivided_size, diameter_squared, &points_grid) {
+					return None;
+				}
+
+				// Fully inside the shape
+				if active_square.fully_in_shape() {
+					Some(ActiveSquare::new(sub_square, true))
+				}
+				// Intersecting the shape's border
+				else {
+					// The sub-square is fully inside the shape if its top-left corner is inside and its edges don't intersect the shape border
+					let point_with_offset = sub_square + offset;
+					let square_edges_intersect_shape = {
+						let min = point_with_offset;
+						let max = min + DVec2::splat(subdivided_size);
+
+						// Top edge line
+						line_intersect_shape_checker((min.x, min.y), (max.x, min.y)) ||
+						// Right edge line
+						line_intersect_shape_checker((max.x, min.y), (max.x, max.y)) ||
+						// Bottom edge line
+						line_intersect_shape_checker((max.x, max.y), (min.x, max.y)) ||
+						// Left edge line
+						line_intersect_shape_checker((min.x, max.y), (min.x, min.y))
+					};
+					let sub_square_fully_inside_shape = !square_edges_intersect_shape && point_in_shape_checker(point_with_offset) && point_in_shape_checker(point_with_offset + subdivided_size);
+
+					Some(ActiveSquare::new(sub_square, sub_square_fully_inside_shape))
+				}
+			});
+			next_level_deeper.add_squares(new_sub_squares);
+		}
+	}
+
+	points_grid.final_points(offset)
+}
+
+/// Randomly pick a square in the dartboard domain, with probability proportional to its area.
+/// Returns a tuple with the subdivision level depth and the square index at that depth.
+fn target_active_square(active_list_levels: &[ActiveListLevel], rng: &mut impl FnMut() -> f64) -> (usize, usize) {
+	let active_squares_total_area: f64 = active_list_levels.iter().map(|active_list| active_list.total_area()).sum();
+	let mut index_into_area = rng() * active_squares_total_area;
+
+	for (level, active_list_level) in active_list_levels.iter().enumerate() {
+		let subtracted = index_into_area - active_list_level.total_area();
+		if subtracted > 0. {
+			index_into_area = subtracted;
+			continue;
+		}
+
+		let active_square_index_in_level = (index_into_area / active_list_levels[level].square_area()).floor() as usize;
+		return (level, active_square_index_in_level);
+	}
+
+	panic!("index_into_area couldn't be be mapped to a square in any level of the active lists");
+}
+
+fn point_not_covered_by_poisson_points(point: DVec2, diameter_squared: f64, points_grid: &AccelerationGrid) -> bool {
+	points_grid.nearby_points(point).all(|nearby_point| {
+		let x_separation = nearby_point.x - point.x;
+		let y_separation = nearby_point.y - point.y;
+
+		x_separation.powi(2) + y_separation.powi(2) > diameter_squared
+	})
+}
+
+fn square_not_covered_by_poisson_points(point: DVec2, half_square_size: f64, diameter_squared: f64, points_grid: &AccelerationGrid) -> bool {
+	let square_center_x = point.x + half_square_size;
+	let square_center_y = point.y + half_square_size;
+
+	points_grid.nearby_points(point).all(|nearby_point| {
+		let x_distance = (square_center_x - nearby_point.x).abs() + half_square_size;
+		let y_distance = (square_center_y - nearby_point.y).abs() + half_square_size;
+
+		x_distance.powi(2) + y_distance.powi(2) > diameter_squared
+	})
+}
+
+#[inline(always)]
+fn cartesian_product<A, B>(a: A, b: B) -> impl Iterator<Item = (A::Item, B::Item)>
+where
+	A: Iterator + Clone,
+	B: Iterator + Clone,
+	A::Item: Clone,
+	B::Item: Clone,
+{
+	a.flat_map(move |i| (b.clone().map(move |j| (i.clone(), j))))
+}
+
+/// A square (represented by its top left corner position and width/height of `square_size`) that is currently a candidate for targetting by the dart throwing process.
+/// The positive sign bit encodes if the square is contained entirely within the masking shape, or negative if it's outside or intersects the shape path.
+pub struct ActiveSquare(DVec2);
+
+impl ActiveSquare {
+	pub fn new(top_left_corner: DVec2, fully_in_shape: bool) -> Self {
+		Self(if fully_in_shape { top_left_corner } else { -top_left_corner })
+	}
+
+	pub fn top_left_corner(&self) -> DVec2 {
+		self.0.abs()
+	}
+
+	pub fn fully_in_shape(&self) -> bool {
+		self.0.x.is_sign_positive()
+	}
+}
+
+pub struct ActiveListLevel {
+	/// List of all subdivided squares of the same size that are currently candidates for targetting by the dart throwing process
+	active_squares: Vec<ActiveSquare>,
+	/// Width and height of the squares in this level of subdivision
+	square_size: f64,
+	/// Current sum of the area in all active squares in this subdivision level
+	total_area: f64,
+}
+
+impl ActiveListLevel {
+	#[inline(always)]
+	pub fn new(square_size: f64) -> Self {
+		Self {
+			active_squares: Vec::new(),
+			square_size,
+			total_area: 0.,
+		}
+	}
+
+	pub fn new_filled(
+		square_size: f64,
+		offset: DVec2,
+		width: f64,
+		height: f64,
+		point_in_shape_checker: impl Fn(DVec2) -> bool,
+		line_intersect_shape_checker: impl Fn((f64, f64), (f64, f64)) -> bool,
+	) -> Self {
+		// These should divide evenly but rounding is to protect against small numerical imprecision errors
+		let x_squares = (width / square_size).round() as usize;
+		let y_squares = (height / square_size).round() as usize;
+
+		// Hashes based on the grid cell coordinates and direction of the line: (x, y, is_vertical)
+		let mut line_intersection_cache: HashMap<(usize, usize, bool), bool> = HashMap::new();
+
+		// Populate each square with its top-left corner coordinate
+		let active_squares: Vec<_> = cartesian_product(0..x_squares, 0..y_squares)
+			.filter_map(|(x, y)| {
+				let corner = DVec2::new(x as f64 * square_size, y as f64 * square_size);
+				let corner_with_offset = corner + offset;
+
+				// Lazily check (and cache) if the square's edges intersect the shape, which is an expensive operation
+				let mut square_edges_intersect_shape_value = None;
+				let mut square_edges_intersect_shape = || {
+					square_edges_intersect_shape_value.unwrap_or_else(|| {
+						let square_edges_intersect_shape = {
+							let min = corner_with_offset;
+							let max = min + DVec2::splat(square_size);
+
+							// Top edge line
+							*line_intersection_cache.entry((x, y, false)).or_insert_with(|| line_intersect_shape_checker((min.x, min.y), (max.x, min.y))) ||
+							// Right edge line
+							*line_intersection_cache.entry((x + 1, y, true)).or_insert_with(|| line_intersect_shape_checker((max.x, min.y), (max.x, max.y))) ||
+							// Bottom edge line
+							*line_intersection_cache.entry((x, y + 1, false)).or_insert_with(|| line_intersect_shape_checker((max.x, max.y), (min.x, max.y))) ||
+							// Left edge line
+							*line_intersection_cache.entry((x, y, true)).or_insert_with(|| line_intersect_shape_checker((min.x, max.y), (min.x, min.y)))
+						};
+						square_edges_intersect_shape_value = Some(square_edges_intersect_shape);
+						square_edges_intersect_shape
+					})
+				};
+
+				// Check if this cell's top-left corner is inside the shape
+				let point_in_shape = point_in_shape_checker(corner_with_offset);
+
+				// Determine if the square is inside the shape
+				let square_not_outside_shape = point_in_shape || square_edges_intersect_shape();
+				if square_not_outside_shape {
+					// Check if this cell's bottom-right corner is inside the shape
+					let opposite_corner_with_offset = DVec2::new((x + 1) as f64 * square_size, (y + 1) as f64 * square_size) + offset;
+					let opposite_corner_in_shape = point_in_shape_checker(opposite_corner_with_offset);
+
+					let square_in_shape = opposite_corner_in_shape && !square_edges_intersect_shape();
+					Some(ActiveSquare::new(corner, square_in_shape))
+				} else {
+					None
+				}
+			})
+			.collect();
+
+		// Sum every square's area to get the total
+		let total_area = square_size.powi(2) * active_squares.len() as f64;
+
+		Self {
+			active_squares,
+			square_size,
+			total_area,
+		}
+	}
+
+	#[must_use]
+	#[inline(always)]
+	pub fn take_square(&mut self, active_square_index: usize) -> ActiveSquare {
+		let targetted_square = self.active_squares.swap_remove(active_square_index);
+		self.total_area = self.square_size.powi(2) * self.active_squares.len() as f64;
+		targetted_square
+	}
+
+	#[inline(always)]
+	pub fn add_squares(&mut self, new_squares: impl Iterator<Item = ActiveSquare>) {
+		for new_square in new_squares {
+			self.active_squares.push(new_square);
+		}
+		self.total_area = self.square_size.powi(2) * self.active_squares.len() as f64;
+	}
+
+	#[inline(always)]
+	pub fn square_size(&self) -> f64 {
+		self.square_size
+	}
+
+	#[inline(always)]
+	pub fn square_area(&self) -> f64 {
+		self.square_size.powi(2)
+	}
+
+	#[inline(always)]
+	pub fn total_area(&self) -> f64 {
+		self.total_area
+	}
+
+	#[inline(always)]
+	pub fn not_empty(&self) -> bool {
+		!self.active_squares.is_empty()
+	}
+}
+
+#[derive(Clone, Default)]
+pub struct PointsList {
+	// The worst-case number of points in a 3x3 grid is 16 (one at each intersection of the four gridlines per axis)
+	storage_slots: [DVec2; 16],
+	length: usize,
+}
+
+impl PointsList {
+	#[inline(always)]
+	pub fn push(&mut self, point: DVec2) {
+		self.storage_slots[self.length] = point;
+		self.length += 1;
+	}
+
+	#[inline(always)]
+	pub fn list_cell_and_neighbors(&self) -> impl Iterator<Item = DVec2> {
+		// The negative bit is used to store whether a point belongs to a neighboring cell
+		self.storage_slots.into_iter().take(self.length).map(|point| (point.x.abs(), point.y.abs()).into())
+	}
+
+	#[inline(always)]
+	pub fn list_cell(&self) -> impl Iterator<Item = DVec2> {
+		// The negative bit is used to store whether a point belongs to a neighboring cell
+		self.storage_slots
+			.into_iter()
+			.take(self.length)
+			.filter(|point| point.x.is_sign_positive() && point.y.is_sign_positive())
+	}
+}
+
+pub struct AccelerationGrid {
+	size: f64,
+	dimension_x: usize,
+	dimension_y: usize,
+	cells: Vec<PointsList>,
+}
+
+impl AccelerationGrid {
+	#[inline(always)]
+	pub fn new(width: f64, height: f64, size: f64) -> Self {
+		let dimension_x = (width / size).ceil() as usize + 1;
+		let dimension_y = (height / size).ceil() as usize + 1;
+
+		Self {
+			size,
+			dimension_x,
+			dimension_y,
+			cells: vec![PointsList::default(); dimension_x * dimension_y],
+		}
+	}
+
+	#[inline(always)]
+	pub fn insert(&mut self, point: DVec2) {
+		let x = (point.x / self.size).floor() as usize;
+		let y = (point.y / self.size).floor() as usize;
+
+		// Insert this point at this cell and the surrounding cells in a 3x3 patch
+		for (x_offset, y_offset) in cartesian_product((-1)..=1, (-1)..=1) {
+			// Avoid going negative
+			let (x, y) = (x as isize + x_offset, y as isize + y_offset);
+			if x < 0 || y < 0 {
+				continue;
+			}
+			// Avoid going beyond the width or height
+			let (x, y) = (x as usize, y as usize);
+			if x > self.dimension_x - 1 || y > self.dimension_y - 1 {
+				continue;
+			}
+
+			// Get the cell corresponding to the (x, y) index
+			let cell = &mut self.cells[y * self.dimension_x + x];
+
+			// Store the given point in this grid cell, and use the negative bit to indicate if this belongs to a neighboring cell
+			cell.push(if x_offset == 0 && y_offset == 0 { point } else { -point });
+		}
+	}
+
+	#[inline(always)]
+	pub fn nearby_points(&self, point: DVec2) -> impl Iterator<Item = DVec2> {
+		let x = (point.x / self.size).floor() as usize;
+		let y = (point.y / self.size).floor() as usize;
+
+		self.cells[y * self.dimension_x + x].list_cell_and_neighbors()
+	}
+
+	#[inline(always)]
+	pub fn final_points(&self, offset: DVec2) -> Vec<DVec2> {
+		self.cells.iter().flat_map(|cell| cell.list_cell()).map(|point| point + offset).collect()
+	}
+}

node-graph/gvector-nodes/src/algorithms/spline.rs

@@ -0,0 +1,182 @@
+use glam::DVec2;
+
+/// Solve for the first handle of an open spline. (The opposite handle can be found by mirroring the result about the anchor.)
+pub fn solve_spline_first_handle_open(points: &[DVec2]) -> Vec<DVec2> {
+	let len_points = points.len();
+	if len_points == 0 {
+		return Vec::new();
+	}
+	if len_points == 1 {
+		return vec![points[0]];
+	}
+
+	// Matrix coefficients a, b and c (see https://mathworld.wolfram.com/CubicSpline.html).
+	// Because the `a` coefficients are all 1, they need not be stored.
+	// This algorithm does a variation of the above algorithm.
+	// Instead of using the traditional cubic (a + bt + ct^2 + dt^3), we use the bezier cubic.
+
+	let mut b = vec![DVec2::new(4., 4.); len_points];
+	b[0] = DVec2::new(2., 2.);
+	b[len_points - 1] = DVec2::new(2., 2.);
+
+	let mut c = vec![DVec2::new(1., 1.); len_points];
+
+	// 'd' is the the second point in a cubic bezier, which is what we solve for
+	let mut d = vec![DVec2::ZERO; len_points];
+
+	d[0] = DVec2::new(2. * points[1].x + points[0].x, 2. * points[1].y + points[0].y);
+	d[len_points - 1] = DVec2::new(3. * points[len_points - 1].x, 3. * points[len_points - 1].y);
+	for idx in 1..(len_points - 1) {
+		d[idx] = DVec2::new(4. * points[idx].x + 2. * points[idx + 1].x, 4. * points[idx].y + 2. * points[idx + 1].y);
+	}
+
+	// Solve with Thomas algorithm (see https://en.wikipedia.org/wiki/Tridiagonal_matrix_algorithm)
+	// Now we do row operations to eliminate `a` coefficients.
+	c[0] /= -b[0];
+	d[0] /= -b[0];
+	#[allow(clippy::assign_op_pattern)]
+	for i in 1..len_points {
+		b[i] += c[i - 1];
+		// For some reason this `+=` version makes the borrow checker mad:
+		// d[i] += d[i-1]
+		d[i] = d[i] + d[i - 1];
+		c[i] /= -b[i];
+		d[i] /= -b[i];
+	}
+
+	// At this point b[i] == -a[i + 1] and a[i] == 0.
+	// Now we do row operations to eliminate 'c' coefficients and solve.
+	d[len_points - 1] *= -1.;
+	#[allow(clippy::assign_op_pattern)]
+	for i in (0..len_points - 1).rev() {
+		d[i] = d[i] - (c[i] * d[i + 1]);
+		d[i] *= -1.; // d[i] /= b[i]
+	}
+
+	d
+}
+
+/// Solve for the first handle of a closed spline. (The opposite handle can be found by mirroring the result about the anchor.)
+/// If called with fewer than 3 points, this function will return an empty result.
+pub fn solve_spline_first_handle_closed(points: &[DVec2]) -> Vec<DVec2> {
+	let len_points = points.len();
+	if len_points < 3 {
+		return Vec::new();
+	}
+
+	// Matrix coefficients `a`, `b` and `c` (see https://mathworld.wolfram.com/CubicSpline.html).
+	// We don't really need to allocate them but it keeps the maths understandable.
+	let a = vec![DVec2::ONE; len_points];
+	let b = vec![DVec2::splat(4.); len_points];
+	let c = vec![DVec2::ONE; len_points];
+
+	let mut cmod = vec![DVec2::ZERO; len_points];
+	let mut u = vec![DVec2::ZERO; len_points];
+
+	// `x` is initially the output of the matrix multiplication, but is converted to the second value.
+	let mut x = vec![DVec2::ZERO; len_points];
+
+	for (i, point) in x.iter_mut().enumerate() {
+		let previous_i = i.checked_sub(1).unwrap_or(len_points - 1);
+		let next_i = (i + 1) % len_points;
+		*point = 3. * (points[next_i] - points[previous_i]);
+	}
+
+	// Solve using https://en.wikipedia.org/wiki/Tridiagonal_matrix_algorithm#Variants (the variant using periodic boundary conditions).
+	// This code below is based on the reference C language implementation provided in that section of the article.
+	let alpha = a[0];
+	let beta = c[len_points - 1];
+
+	// Arbitrary, but chosen such that division by zero is avoided.
+	let gamma = -b[0];
+
+	cmod[0] = alpha / (b[0] - gamma);
+	u[0] = gamma / (b[0] - gamma);
+	x[0] /= b[0] - gamma;
+
+	// Handle from from `1` to `len_points - 2` (inclusive).
+	for ix in 1..=(len_points - 2) {
+		let m = 1.0 / (b[ix] - a[ix] * cmod[ix - 1]);
+		cmod[ix] = c[ix] * m;
+		u[ix] = (0.0 - a[ix] * u[ix - 1]) * m;
+		x[ix] = (x[ix] - a[ix] * x[ix - 1]) * m;
+	}
+
+	// Handle `len_points - 1`.
+	let m = 1.0 / (b[len_points - 1] - alpha * beta / gamma - beta * cmod[len_points - 2]);
+	u[len_points - 1] = (alpha - a[len_points - 1] * u[len_points - 2]) * m;
+	x[len_points - 1] = (x[len_points - 1] - a[len_points - 1] * x[len_points - 2]) * m;
+
+	// Loop from `len_points - 2` to `0` (inclusive).
+	for ix in (0..=(len_points - 2)).rev() {
+		u[ix] = u[ix] - cmod[ix] * u[ix + 1];
+		x[ix] = x[ix] - cmod[ix] * x[ix + 1];
+	}
+
+	let fact = (x[0] + x[len_points - 1] * beta / gamma) / (1.0 + u[0] + u[len_points - 1] * beta / gamma);
+
+	for ix in 0..(len_points) {
+		x[ix] -= fact * u[ix];
+	}
+
+	let mut real = vec![DVec2::ZERO; len_points];
+	for i in 0..len_points {
+		let previous = i.checked_sub(1).unwrap_or(len_points - 1);
+		let next = (i + 1) % len_points;
+		real[i] = x[previous] * a[next] + x[i] * b[i] + x[next] * c[i];
+	}
+
+	// The matrix is now solved.
+
+	// Since we have computed the derivative, work back to find the start handle.
+	for i in 0..len_points {
+		x[i] = (x[i] / 3.) + points[i];
+	}
+
+	x
+}
+
+#[cfg(test)]
+mod tests {
+	use super::*;
+	#[test]
+	fn closed_spline() {
+		use graphene_vector::{dvec2_to_point, point_to_dvec2};
+		use kurbo::{BezPath, ParamCurve, ParamCurveDeriv};
+
+		// These points are just chosen arbitrary
+		let points = [DVec2::new(0., 0.), DVec2::new(0., 0.), DVec2::new(6., 5.), DVec2::new(7., 9.), DVec2::new(2., 3.)];
+
+		// List of first handle or second point in a cubic bezier curve.
+		let first_handles = solve_spline_first_handle_closed(&points);
+
+		// Construct the Subpath
+		let mut bezpath = BezPath::new();
+		bezpath.move_to(dvec2_to_point(points[0]));
+
+		for i in 0..first_handles.len() {
+			let next_i = i + 1;
+			let next_i = if next_i == first_handles.len() { 0 } else { next_i };
+
+			// First handle or second point of a cubic Bezier curve.
+			let p1 = dvec2_to_point(first_handles[i]);
+			// Second handle or third point of a cubic Bezier curve.
+			let p2 = dvec2_to_point(2. * points[next_i] - first_handles[next_i]);
+			// Endpoint or fourth point of a cubic Bezier curve.
+			let p3 = dvec2_to_point(points[next_i]);
+
+			bezpath.curve_to(p1, p2, p3);
+		}
+
+		// For each pair of bézier curves, ensure that the second derivative is continuous
+		for (bézier_a, bézier_b) in bezpath.segments().zip(bezpath.segments().skip(1).chain(bezpath.segments().take(1))) {
+			let derivative2_end_a = point_to_dvec2(bézier_a.to_cubic().deriv().eval(1.));
+			let derivative2_start_b = point_to_dvec2(bézier_b.to_cubic().deriv().eval(0.));
+
+			assert!(
+				derivative2_end_a.abs_diff_eq(derivative2_start_b, 1e-10),
+				"second derivative at the end of a {derivative2_end_a} is equal to the second derivative at the start of b {derivative2_start_b}"
+			);
+		}
+	}
+}

node-graph/gvector-nodes/src/generator_nodes.rs

@@ -0,0 +1,282 @@
+use super::misc::{ArcType, AsU64, GridType};
+use bezier_rs::Subpath;
+use glam::DVec2;
+use graphene_core::context::Ctx;
+use graphene_core::registry::types::{Angle, PixelSize};
+use graphene_vector::{HandleId, PointId, SegmentId, StrokeId, VectorData, VectorDataTable};
+
+trait CornerRadius {
+	fn generate(self, size: DVec2, clamped: bool) -> VectorDataTable;
+}
+impl CornerRadius for f64 {
+	fn generate(self, size: DVec2, clamped: bool) -> VectorDataTable {
+		let clamped_radius = if clamped { self.clamp(0., size.x.min(size.y).max(0.) / 2.) } else { self };
+		VectorDataTable::new(VectorData::from_subpath(Subpath::new_rounded_rect(size / -2., size / 2., [clamped_radius; 4])))
+	}
+}
+impl CornerRadius for [f64; 4] {
+	fn generate(self, size: DVec2, clamped: bool) -> VectorDataTable {
+		let clamped_radius = if clamped {
+			// Algorithm follows the CSS spec: <https://drafts.csswg.org/css-backgrounds/#corner-overlap>
+
+			let mut scale_factor: f64 = 1.;
+			for i in 0..4 {
+				let side_length = if i % 2 == 0 { size.x } else { size.y };
+				let adjacent_corner_radius_sum = self[i] + self[(i + 1) % 4];
+				if side_length < adjacent_corner_radius_sum {
+					scale_factor = scale_factor.min(side_length / adjacent_corner_radius_sum);
+				}
+			}
+			self.map(|x| x * scale_factor)
+		} else {
+			self
+		};
+		VectorDataTable::new(VectorData::from_subpath(Subpath::new_rounded_rect(size / -2., size / 2., clamped_radius)))
+	}
+}
+
+#[node_macro::node(category("Vector: Shape"))]
+fn circle(_: impl Ctx, _primary: (), #[default(50.)] radius: f64) -> VectorDataTable {
+	let radius = radius.abs();
+	VectorDataTable::new(VectorData::from_subpath(Subpath::new_ellipse(DVec2::splat(-radius), DVec2::splat(radius))))
+}
+
+#[node_macro::node(category("Vector: Shape"))]
+fn arc(
+	_: impl Ctx,
+	_primary: (),
+	#[default(50.)] radius: f64,
+	start_angle: Angle,
+	#[default(270.)]
+	#[range((0., 360.))]
+	sweep_angle: Angle,
+	arc_type: ArcType,
+) -> VectorDataTable {
+	VectorDataTable::new(VectorData::from_subpath(Subpath::new_arc(
+		radius,
+		start_angle / 360. * std::f64::consts::TAU,
+		sweep_angle / 360. * std::f64::consts::TAU,
+		match arc_type {
+			ArcType::Open => bezier_rs::ArcType::Open,
+			ArcType::Closed => bezier_rs::ArcType::Closed,
+			ArcType::PieSlice => bezier_rs::ArcType::PieSlice,
+		},
+	)))
+}
+
+#[node_macro::node(category("Vector: Shape"))]
+fn ellipse(_: impl Ctx, _primary: (), #[default(50)] radius_x: f64, #[default(25)] radius_y: f64) -> VectorDataTable {
+	let radius = DVec2::new(radius_x, radius_y);
+	let corner1 = -radius;
+	let corner2 = radius;
+
+	let mut ellipse = VectorData::from_subpath(Subpath::new_ellipse(corner1, corner2));
+
+	let len = ellipse.segment_domain.ids().len();
+	for i in 0..len {
+		ellipse
+			.colinear_manipulators
+			.push([HandleId::end(ellipse.segment_domain.ids()[i]), HandleId::primary(ellipse.segment_domain.ids()[(i + 1) % len])]);
+	}
+
+	VectorDataTable::new(ellipse)
+}
+
+#[node_macro::node(category("Vector: Shape"), properties("rectangle_properties"))]
+fn rectangle<T: CornerRadius>(
+	_: impl Ctx,
+	_primary: (),
+	#[default(100)] width: f64,
+	#[default(100)] height: f64,
+	_individual_corner_radii: bool, // TODO: Move this to the bottom once we have a migration capability
+	#[implementations(f64, [f64; 4])] corner_radius: T,
+	#[default(true)] clamped: bool,
+) -> VectorDataTable {
+	corner_radius.generate(DVec2::new(width, height), clamped)
+}
+
+#[node_macro::node(category("Vector: Shape"))]
+fn regular_polygon<T: AsU64>(
+	_: impl Ctx,
+	_primary: (),
+	#[default(6)]
+	#[hard_min(3.)]
+	#[implementations(u32, u64, f64)]
+	sides: T,
+	#[default(50)] radius: f64,
+) -> VectorDataTable {
+	let points = sides.as_u64();
+	let radius: f64 = radius * 2.;
+	VectorDataTable::new(VectorData::from_subpath(Subpath::new_regular_polygon(DVec2::splat(-radius), points, radius)))
+}
+
+#[node_macro::node(category("Vector: Shape"))]
+fn star<T: AsU64>(
+	_: impl Ctx,
+	_primary: (),
+	#[default(5)]
+	#[hard_min(2.)]
+	#[implementations(u32, u64, f64)]
+	sides: T,
+	#[default(50)] radius_1: f64,
+	#[default(25)] radius_2: f64,
+) -> VectorDataTable {
+	let points = sides.as_u64();
+	let diameter: f64 = radius_1 * 2.;
+	let inner_diameter = radius_2 * 2.;
+
+	VectorDataTable::new(VectorData::from_subpath(Subpath::new_star_polygon(DVec2::splat(-diameter), points, diameter, inner_diameter)))
+}
+
+#[node_macro::node(category("Vector: Shape"))]
+fn line(_: impl Ctx, _primary: (), #[default((0., -50.))] start: PixelSize, #[default((0., 50.))] end: PixelSize) -> VectorDataTable {
+	VectorDataTable::new(VectorData::from_subpath(Subpath::new_line(start, end)))
+}
+
+trait GridSpacing {
+	fn as_dvec2(&self) -> DVec2;
+}
+impl GridSpacing for f64 {
+	fn as_dvec2(&self) -> DVec2 {
+		DVec2::splat(*self)
+	}
+}
+impl GridSpacing for DVec2 {
+	fn as_dvec2(&self) -> DVec2 {
+		*self
+	}
+}
+
+#[node_macro::node(category("Vector: Shape"), properties("grid_properties"))]
+fn grid<T: GridSpacing>(
+	_: impl Ctx,
+	_primary: (),
+	grid_type: GridType,
+	#[hard_min(0.)]
+	#[default(10)]
+	#[implementations(f64, DVec2)]
+	spacing: T,
+	#[default(30., 30.)] angles: DVec2,
+	#[default(10)] columns: u32,
+	#[default(10)] rows: u32,
+) -> VectorDataTable {
+	let (x_spacing, y_spacing) = spacing.as_dvec2().into();
+	let (angle_a, angle_b) = angles.into();
+
+	let mut vector_data = VectorData::default();
+	let mut segment_id = SegmentId::ZERO;
+	let mut point_id = PointId::ZERO;
+
+	match grid_type {
+		GridType::Rectangular => {
+			// Create rectangular grid points and connect them with line segments
+			for y in 0..rows {
+				for x in 0..columns {
+					// Add current point to the grid
+					let current_index = vector_data.point_domain.ids().len();
+					vector_data.point_domain.push(point_id.next_id(), DVec2::new(x_spacing * x as f64, y_spacing * y as f64));
+
+					// Helper function to connect points with line segments
+					let mut push_segment = |to_index: Option<usize>| {
+						if let Some(other_index) = to_index {
+							vector_data
+								.segment_domain
+								.push(segment_id.next_id(), other_index, current_index, bezier_rs::BezierHandles::Linear, StrokeId::ZERO);
+						}
+					};
+
+					// Connect to the point to the left (horizontal connection)
+					push_segment((x > 0).then(|| current_index - 1));
+
+					// Connect to the point above (vertical connection)
+					push_segment(current_index.checked_sub(columns as usize));
+				}
+			}
+		}
+		GridType::Isometric => {
+			// Calculate isometric grid spacing based on angles
+			let tan_a = angle_a.to_radians().tan();
+			let tan_b = angle_b.to_radians().tan();
+			let spacing = DVec2::new(y_spacing / (tan_a + tan_b), y_spacing);
+
+			// Create isometric grid points and connect them with line segments
+			for y in 0..rows {
+				for x in 0..columns {
+					// Add current point to the grid with offset for odd columns
+					let current_index = vector_data.point_domain.ids().len();
+
+					let a_angles_eaten = x.div_ceil(2) as f64;
+					let b_angles_eaten = (x / 2) as f64;
+
+					let offset_y_fraction = b_angles_eaten * tan_b - a_angles_eaten * tan_a;
+
+					let position = DVec2::new(spacing.x * x as f64, spacing.y * y as f64 + offset_y_fraction * spacing.x);
+					vector_data.point_domain.push(point_id.next_id(), position);
+
+					// Helper function to connect points with line segments
+					let mut push_segment = |to_index: Option<usize>| {
+						if let Some(other_index) = to_index {
+							vector_data
+								.segment_domain
+								.push(segment_id.next_id(), other_index, current_index, bezier_rs::BezierHandles::Linear, StrokeId::ZERO);
+						}
+					};
+
+					// Connect to the point to the left
+					push_segment((x > 0).then(|| current_index - 1));
+
+					// Connect to the point directly above
+					push_segment(current_index.checked_sub(columns as usize));
+
+					// Additional diagonal connections for odd columns (creates hexagonal pattern)
+					if x % 2 == 1 {
+						// Connect to the point diagonally up-right (if not at right edge)
+						push_segment(current_index.checked_sub(columns as usize - 1).filter(|_| x + 1 < columns));
+
+						// Connect to the point diagonally up-left
+						push_segment(current_index.checked_sub(columns as usize + 1));
+					}
+				}
+			}
+		}
+	}
+
+	VectorDataTable::new(vector_data)
+}
+
+#[cfg(test)]
+mod tests {
+	use super::*;
+	#[test]
+	fn isometric_grid_test() {
+		// Doesn't crash with weird angles
+		grid((), (), GridType::Isometric, 0., (0., 0.).into(), 5, 5);
+		grid((), (), GridType::Isometric, 90., (90., 90.).into(), 5, 5);
+
+		// Works properly
+		let grid = grid((), (), GridType::Isometric, 10., (30., 30.).into(), 5, 5);
+		assert_eq!(grid.instance_ref_iter().next().unwrap().instance.point_domain.ids().len(), 5 * 5);
+		assert_eq!(grid.instance_ref_iter().next().unwrap().instance.segment_bezier_iter().count(), 4 * 5 + 4 * 9);
+		for (_, bezier, _, _) in grid.instance_ref_iter().next().unwrap().instance.segment_bezier_iter() {
+			assert_eq!(bezier.handles, bezier_rs::BezierHandles::Linear);
+			assert!(
+				((bezier.start - bezier.end).length() - 10.).abs() < 1e-5,
+				"Length of {} should be 10",
+				(bezier.start - bezier.end).length()
+			);
+		}
+	}
+
+	#[test]
+	fn skew_isometric_grid_test() {
+		let grid = grid((), (), GridType::Isometric, 10., (40., 30.).into(), 5, 5);
+		assert_eq!(grid.instance_ref_iter().next().unwrap().instance.point_domain.ids().len(), 5 * 5);
+		assert_eq!(grid.instance_ref_iter().next().unwrap().instance.segment_bezier_iter().count(), 4 * 5 + 4 * 9);
+		for (_, bezier, _, _) in grid.instance_ref_iter().next().unwrap().instance.segment_bezier_iter() {
+			assert_eq!(bezier.handles, bezier_rs::BezierHandles::Linear);
+			let vector = bezier.start - bezier.end;
+			let angle = (vector.angle_to(DVec2::X).to_degrees() + 180.) % 180.;
+			assert!([90., 150., 40.].into_iter().any(|target| (target - angle).abs() < 1e-10), "unexpected angle of {}", angle)
+		}
+	}
+}

node-graph/gvector-nodes/src/lib.rs

@@ -0,0 +1,5 @@
+pub mod algorithms;
+pub mod generator_nodes;
+pub mod misc;
+pub mod modification;
+pub mod vector_nodes;

node-graph/gvector-nodes/src/misc.rs

@@ -0,0 +1,98 @@
+use dyn_any::DynAny;
+use glam::DVec2;
+use kurbo::Point;
+
+/// Represents different ways of calculating the centroid.
+#[derive(Default, Debug, Clone, Copy, PartialEq, Eq, serde::Serialize, serde::Deserialize, Hash, DynAny, specta::Type, node_macro::ChoiceType)]
+#[widget(Radio)]
+pub enum CentroidType {
+	/// The center of mass for the area of a solid shape's interior, as if made out of an infinitely flat material.
+	#[default]
+	Area,
+	/// The center of mass for the arc length of a curved shape's perimeter, as if made out of an infinitely thin wire.
+	Length,
+}
+
+pub trait AsU64 {
+	fn as_u64(&self) -> u64;
+}
+impl AsU64 for u32 {
+	fn as_u64(&self) -> u64 {
+		*self as u64
+	}
+}
+impl AsU64 for u64 {
+	fn as_u64(&self) -> u64 {
+		*self
+	}
+}
+impl AsU64 for f64 {
+	fn as_u64(&self) -> u64 {
+		*self as u64
+	}
+}
+
+pub trait AsI64 {
+	fn as_i64(&self) -> i64;
+}
+impl AsI64 for u32 {
+	fn as_i64(&self) -> i64 {
+		*self as i64
+	}
+}
+impl AsI64 for u64 {
+	fn as_i64(&self) -> i64 {
+		*self as i64
+	}
+}
+impl AsI64 for f64 {
+	fn as_i64(&self) -> i64 {
+		*self as i64
+	}
+}
+
+#[derive(Default, Debug, Clone, Copy, PartialEq, Eq, serde::Serialize, serde::Deserialize, Hash, DynAny, specta::Type, node_macro::ChoiceType)]
+#[widget(Radio)]
+pub enum GridType {
+	#[default]
+	Rectangular,
+	Isometric,
+}
+
+#[repr(C)]
+#[derive(Debug, Clone, Copy, Default, PartialEq, Eq, serde::Serialize, serde::Deserialize, Hash, DynAny, specta::Type, node_macro::ChoiceType)]
+#[widget(Radio)]
+pub enum ArcType {
+	#[default]
+	Open,
+	Closed,
+	PieSlice,
+}
+
+#[repr(C)]
+#[derive(Debug, Clone, Copy, Default, PartialEq, Eq, serde::Serialize, serde::Deserialize, Hash, DynAny, specta::Type, node_macro::ChoiceType)]
+#[widget(Radio)]
+pub enum MergeByDistanceAlgorithm {
+	#[default]
+	Spatial,
+	Topological,
+}
+
+#[repr(C)]
+#[derive(Debug, Clone, Copy, Default, PartialEq, Eq, serde::Serialize, serde::Deserialize, Hash, DynAny, specta::Type, node_macro::ChoiceType)]
+#[widget(Radio)]
+pub enum PointSpacingType {
+	#[default]
+	/// The desired spacing distance between points.
+	Separation,
+	/// The exact number of points to span the path.
+	Quantity,
+}
+
+pub fn point_to_dvec2(point: Point) -> DVec2 {
+	DVec2 { x: point.x, y: point.y }
+}
+
+pub fn dvec2_to_point(value: DVec2) -> Point {
+	Point { x: value.x, y: value.y }
+}

node-graph/gvector-nodes/src/modification.rs

@@ -0,0 +1,728 @@
+use crate::misc::point_to_dvec2;
+use bezier_rs::BezierHandles;
+use dyn_any::DynAny;
+use glam::DVec2;
+use graphene_core::context::Ctx;
+use graphene_core::instances::Instance;
+use graphene_core::uuid::generate_uuid;
+use graphene_vector::{FillId, HandleId, HandleType, PointDomain, PointId, RegionDomain, RegionId, SegmentDomain, SegmentId, StrokeId, VectorData, VectorDataTable};
+use kurbo::{BezPath, PathEl, Point};
+use log::warn;
+use serde::de::{SeqAccess, Visitor};
+use serde::ser::SerializeSeq;
+use serde::{Deserialize, Deserializer, Serialize, Serializer};
+use std::collections::{HashMap, HashSet};
+use std::fmt;
+use std::hash::BuildHasher;
+use std::hash::Hash;
+
+/// Represents a procedural change to the [`PointDomain`] in [`VectorData`].
+#[derive(Clone, Debug, Default, PartialEq, serde::Serialize, serde::Deserialize)]
+pub struct PointModification {
+	add: Vec<PointId>,
+	remove: HashSet<PointId>,
+	#[serde(serialize_with = "serialize_hashmap", deserialize_with = "deserialize_hashmap")]
+	delta: HashMap<PointId, DVec2>,
+}
+
+impl Hash for PointModification {
+	fn hash<H: std::hash::Hasher>(&self, state: &mut H) {
+		generate_uuid().hash(state)
+	}
+}
+
+impl PointModification {
+	/// Apply this modification to the specified [`PointDomain`].
+	pub fn apply(&self, point_domain: &mut PointDomain, segment_domain: &mut SegmentDomain) {
+		point_domain.retain(segment_domain, |id| !self.remove.contains(id));
+
+		for (index, (id, position)) in point_domain.positions_mut().enumerate() {
+			let Some(&delta) = self.delta.get(&id) else { continue };
+			if !delta.is_finite() {
+				warn!("Invalid delta when applying a point modification");
+				continue;
+			}
+
+			*position += delta;
+
+			for (_, handles, start, end) in segment_domain.handles_mut() {
+				if start == index {
+					handles.move_start(delta);
+				}
+				if end == index {
+					handles.move_end(delta);
+				}
+			}
+		}
+
+		for &add_id in &self.add {
+			let Some(&position) = self.delta.get(&add_id) else { continue };
+			if !position.is_finite() {
+				warn!("Invalid position when applying a point modification");
+				continue;
+			}
+
+			point_domain.push(add_id, position);
+		}
+	}
+
+	/// Create a new modification that will convert an empty [`VectorData`] into the target [`VectorData`].
+	pub fn create_from_vector(vector_data: &VectorData) -> Self {
+		Self {
+			add: vector_data.point_domain.ids().to_vec(),
+			remove: HashSet::new(),
+			delta: vector_data.point_domain.ids().iter().copied().zip(vector_data.point_domain.positions().iter().cloned()).collect(),
+		}
+	}
+
+	fn push(&mut self, id: PointId, position: DVec2) {
+		self.add.push(id);
+		self.delta.insert(id, position);
+	}
+
+	fn remove(&mut self, id: PointId) {
+		self.remove.insert(id);
+		self.add.retain(|&add| add != id);
+		self.delta.remove(&id);
+	}
+}
+
+/// Represents a procedural change to the [`SegmentDomain`] in [`VectorData`].
+#[derive(Clone, Debug, Default, PartialEq, serde::Serialize, serde::Deserialize)]
+pub struct SegmentModification {
+	add: Vec<SegmentId>,
+	remove: HashSet<SegmentId>,
+	#[serde(serialize_with = "serialize_hashmap", deserialize_with = "deserialize_hashmap")]
+	start_point: HashMap<SegmentId, PointId>,
+	#[serde(serialize_with = "serialize_hashmap", deserialize_with = "deserialize_hashmap")]
+	end_point: HashMap<SegmentId, PointId>,
+	#[serde(serialize_with = "serialize_hashmap", deserialize_with = "deserialize_hashmap")]
+	handle_primary: HashMap<SegmentId, Option<DVec2>>,
+	#[serde(serialize_with = "serialize_hashmap", deserialize_with = "deserialize_hashmap")]
+	handle_end: HashMap<SegmentId, Option<DVec2>>,
+	#[serde(serialize_with = "serialize_hashmap", deserialize_with = "deserialize_hashmap")]
+	stroke: HashMap<SegmentId, StrokeId>,
+}
+
+impl SegmentModification {
+	/// Apply this modification to the specified [`SegmentDomain`].
+	pub fn apply(&self, segment_domain: &mut SegmentDomain, point_domain: &PointDomain) {
+		segment_domain.retain(|id| !self.remove.contains(id), point_domain.ids().len());
+
+		for (id, point) in segment_domain.start_point_mut() {
+			let Some(&new) = self.start_point.get(&id) else { continue };
+			let Some(index) = point_domain.resolve_id(new) else {
+				warn!("Invalid start ID when applying a segment modification");
+				continue;
+			};
+
+			*point = index;
+		}
+
+		for (id, point) in segment_domain.end_point_mut() {
+			let Some(&new) = self.end_point.get(&id) else { continue };
+			let Some(index) = point_domain.resolve_id(new) else {
+				warn!("Invalid end ID when applying a segment modification");
+				continue;
+			};
+
+			*point = index;
+		}
+
+		for (id, handles, start, end) in segment_domain.handles_mut() {
+			let Some(&start) = point_domain.positions().get(start) else { continue };
+			let Some(&end) = point_domain.positions().get(end) else { continue };
+
+			// Compute the actual start and end position based on the offset from the anchor
+			let start = self.handle_primary.get(&id).copied().map(|handle| handle.map(|handle| handle + start));
+			let end = self.handle_end.get(&id).copied().map(|handle| handle.map(|handle| handle + end));
+
+			if !start.unwrap_or_default().is_none_or(|start| start.is_finite()) || !end.unwrap_or_default().is_none_or(|end| end.is_finite()) {
+				warn!("Invalid handles when applying a segment modification");
+				continue;
+			}
+
+			match (start, end) {
+				// The new handles are fully specified by the modification
+				(Some(Some(handle_start)), Some(Some(handle_end))) => *handles = BezierHandles::Cubic { handle_start, handle_end },
+				(Some(Some(handle)), Some(None)) | (Some(None), Some(Some(handle))) => *handles = BezierHandles::Quadratic { handle },
+				(Some(None), Some(None)) => *handles = BezierHandles::Linear,
+				// Remove the end handle
+				(None, Some(None)) => {
+					if let BezierHandles::Cubic { handle_start, .. } = *handles {
+						*handles = BezierHandles::Quadratic { handle: handle_start }
+					}
+				}
+				// Change the end handle
+				(None, Some(Some(handle_end))) => match *handles {
+					BezierHandles::Linear => *handles = BezierHandles::Quadratic { handle: handle_end },
+					BezierHandles::Quadratic { handle: handle_start } => *handles = BezierHandles::Cubic { handle_start, handle_end },
+					BezierHandles::Cubic { handle_start, .. } => *handles = BezierHandles::Cubic { handle_start, handle_end },
+				},
+				// Remove the start handle
+				(Some(None), None) => *handles = BezierHandles::Linear,
+				// Change the start handle
+				(Some(Some(handle_start)), None) => match *handles {
+					BezierHandles::Linear => *handles = BezierHandles::Quadratic { handle: handle_start },
+					BezierHandles::Quadratic { .. } => *handles = BezierHandles::Quadratic { handle: handle_start },
+					BezierHandles::Cubic { handle_end, .. } => *handles = BezierHandles::Cubic { handle_start, handle_end },
+				},
+				// No change
+				(None, None) => {}
+			};
+		}
+
+		for (id, stroke) in segment_domain.stroke_mut() {
+			let Some(&new) = self.stroke.get(&id) else { continue };
+			*stroke = new;
+		}
+
+		for &add_id in &self.add {
+			let Some(&start) = self.start_point.get(&add_id) else { continue };
+			let Some(&end) = self.end_point.get(&add_id) else { continue };
+			let Some(&handle_start) = self.handle_primary.get(&add_id) else { continue };
+			let Some(&handle_end) = self.handle_end.get(&add_id) else { continue };
+			let Some(&stroke) = self.stroke.get(&add_id) else { continue };
+
+			let Some(start_index) = point_domain.resolve_id(start) else {
+				warn!("invalid start id: {:#?}", start);
+				continue;
+			};
+			let Some(end_index) = point_domain.resolve_id(end) else {
+				warn!("invalid end id: {:#?}", end);
+				continue;
+			};
+
+			let start_position = point_domain.positions()[start_index];
+			let end_position = point_domain.positions()[end_index];
+			let handles = match (handle_start, handle_end) {
+				(Some(handle_start), Some(handle_end)) => BezierHandles::Cubic {
+					handle_start: handle_start + start_position,
+					handle_end: handle_end + end_position,
+				},
+				(Some(handle), None) | (None, Some(handle)) => BezierHandles::Quadratic { handle: handle + start_position },
+				(None, None) => BezierHandles::Linear,
+			};
+
+			if !handles.is_finite() {
+				warn!("invalid handles");
+				continue;
+			}
+
+			segment_domain.push(add_id, start_index, end_index, handles, stroke);
+		}
+
+		assert!(
+			segment_domain.start_point().iter().all(|&index| index < point_domain.ids().len()),
+			"index should be in range {:#?}",
+			segment_domain
+		);
+		assert!(
+			segment_domain.end_point().iter().all(|&index| index < point_domain.ids().len()),
+			"index should be in range {:#?}",
+			segment_domain
+		);
+	}
+
+	/// Create a new modification that will convert an empty [`VectorData`] into the target [`VectorData`].
+	pub fn create_from_vector(vector_data: &VectorData) -> Self {
+		let point_id = |(&segment, &index)| (segment, vector_data.point_domain.ids()[index]);
+		Self {
+			add: vector_data.segment_domain.ids().to_vec(),
+			remove: HashSet::new(),
+			start_point: vector_data.segment_domain.ids().iter().zip(vector_data.segment_domain.start_point()).map(point_id).collect(),
+			end_point: vector_data.segment_domain.ids().iter().zip(vector_data.segment_domain.end_point()).map(point_id).collect(),
+			handle_primary: vector_data.segment_bezier_iter().map(|(id, b, _, _)| (id, b.handle_start().map(|handle| handle - b.start))).collect(),
+			handle_end: vector_data.segment_bezier_iter().map(|(id, b, _, _)| (id, b.handle_end().map(|handle| handle - b.end))).collect(),
+			stroke: vector_data.segment_domain.ids().iter().copied().zip(vector_data.segment_domain.stroke().iter().cloned()).collect(),
+		}
+	}
+
+	fn push(&mut self, id: SegmentId, points: [PointId; 2], handles: [Option<DVec2>; 2], stroke: StrokeId) {
+		self.remove.remove(&id);
+		self.add.push(id);
+		self.start_point.insert(id, points[0]);
+		self.end_point.insert(id, points[1]);
+		self.handle_primary.insert(id, handles[0]);
+		self.handle_end.insert(id, handles[1]);
+		self.stroke.insert(id, stroke);
+	}
+
+	fn remove(&mut self, id: SegmentId) {
+		self.remove.insert(id);
+		self.add.retain(|&add| add != id);
+		self.start_point.remove(&id);
+		self.end_point.remove(&id);
+		self.handle_primary.remove(&id);
+		self.handle_end.remove(&id);
+		self.stroke.remove(&id);
+	}
+}
+
+/// Represents a procedural change to the [`RegionDomain`] in [`VectorData`].
+#[derive(Clone, Debug, Default, PartialEq, serde::Serialize, serde::Deserialize)]
+pub struct RegionModification {
+	add: Vec<RegionId>,
+	remove: HashSet<RegionId>,
+	#[serde(serialize_with = "serialize_hashmap", deserialize_with = "deserialize_hashmap")]
+	segment_range: HashMap<RegionId, std::ops::RangeInclusive<SegmentId>>,
+	#[serde(serialize_with = "serialize_hashmap", deserialize_with = "deserialize_hashmap")]
+	fill: HashMap<RegionId, FillId>,
+}
+
+impl RegionModification {
+	/// Apply this modification to the specified [`RegionDomain`].
+	pub fn apply(&self, region_domain: &mut RegionDomain) {
+		region_domain.retain(|id| !self.remove.contains(id));
+
+		for (id, segment_range) in region_domain.segment_range_mut() {
+			let Some(new) = self.segment_range.get(&id) else { continue };
+			*segment_range = new.clone(); // Range inclusive is not copy
+		}
+
+		for (id, fill) in region_domain.fill_mut() {
+			let Some(&new) = self.fill.get(&id) else { continue };
+			*fill = new;
+		}
+
+		for &add_id in &self.add {
+			let Some(segment_range) = self.segment_range.get(&add_id) else { continue };
+			let Some(&fill) = self.fill.get(&add_id) else { continue };
+			region_domain.push(add_id, segment_range.clone(), fill);
+		}
+	}
+
+	/// Create a new modification that will convert an empty [`VectorData`] into the target [`VectorData`].
+	pub fn create_from_vector(vector_data: &VectorData) -> Self {
+		Self {
+			add: vector_data.region_domain.ids().to_vec(),
+			remove: HashSet::new(),
+			segment_range: vector_data.region_domain.ids().iter().copied().zip(vector_data.region_domain.segment_range().iter().cloned()).collect(),
+			fill: vector_data.region_domain.ids().iter().copied().zip(vector_data.region_domain.fill().iter().cloned()).collect(),
+		}
+	}
+}
+
+/// Represents a procedural change to the [`VectorData`].
+#[derive(Clone, Debug, Default, PartialEq, DynAny, serde::Serialize, serde::Deserialize)]
+pub struct VectorModification {
+	points: PointModification,
+	segments: SegmentModification,
+	regions: RegionModification,
+	add_g1_continuous: HashSet<[HandleId; 2]>,
+	remove_g1_continuous: HashSet<[HandleId; 2]>,
+}
+
+/// A modification type that can be added to a [`VectorModification`].
+#[derive(PartialEq, Clone, Debug, serde::Serialize, serde::Deserialize)]
+pub enum VectorModificationType {
+	InsertSegment { id: SegmentId, points: [PointId; 2], handles: [Option<DVec2>; 2] },
+	InsertPoint { id: PointId, position: DVec2 },
+
+	RemoveSegment { id: SegmentId },
+	RemovePoint { id: PointId },
+
+	SetG1Continuous { handles: [HandleId; 2], enabled: bool },
+	SetHandles { segment: SegmentId, handles: [Option<DVec2>; 2] },
+	SetPrimaryHandle { segment: SegmentId, relative_position: DVec2 },
+	SetEndHandle { segment: SegmentId, relative_position: DVec2 },
+	SetStartPoint { segment: SegmentId, id: PointId },
+	SetEndPoint { segment: SegmentId, id: PointId },
+
+	ApplyPointDelta { point: PointId, delta: DVec2 },
+	ApplyPrimaryDelta { segment: SegmentId, delta: DVec2 },
+	ApplyEndDelta { segment: SegmentId, delta: DVec2 },
+}
+
+impl VectorModification {
+	/// Apply this modification to the specified [`VectorData`].
+	pub fn apply(&self, vector_data: &mut VectorData) {
+		self.points.apply(&mut vector_data.point_domain, &mut vector_data.segment_domain);
+		self.segments.apply(&mut vector_data.segment_domain, &vector_data.point_domain);
+		self.regions.apply(&mut vector_data.region_domain);
+
+		let valid = |val: &[HandleId; 2]| vector_data.segment_domain.ids().contains(&val[0].segment) && vector_data.segment_domain.ids().contains(&val[1].segment);
+		vector_data
+			.colinear_manipulators
+			.retain(|val| !self.remove_g1_continuous.contains(val) && !self.remove_g1_continuous.contains(&[val[1], val[0]]) && valid(val));
+
+		for handles in &self.add_g1_continuous {
+			if !vector_data.colinear_manipulators.iter().any(|test| test == handles || test == &[handles[1], handles[0]]) && valid(handles) {
+				vector_data.colinear_manipulators.push(*handles);
+			}
+		}
+	}
+
+	/// Add a [`VectorModificationType`] to this modification.
+	pub fn modify(&mut self, vector_data_modification: &VectorModificationType) {
+		match vector_data_modification {
+			VectorModificationType::InsertSegment { id, points, handles } => self.segments.push(*id, *points, *handles, StrokeId::ZERO),
+			VectorModificationType::InsertPoint { id, position } => self.points.push(*id, *position),
+
+			VectorModificationType::RemoveSegment { id } => self.segments.remove(*id),
+			VectorModificationType::RemovePoint { id } => self.points.remove(*id),
+
+			VectorModificationType::SetG1Continuous { handles, enabled } => {
+				if *enabled {
+					if !self.add_g1_continuous.contains(&[handles[1], handles[0]]) {
+						self.add_g1_continuous.insert(*handles);
+					}
+					self.remove_g1_continuous.remove(handles);
+					self.remove_g1_continuous.remove(&[handles[1], handles[0]]);
+				} else {
+					if !self.remove_g1_continuous.contains(&[handles[1], handles[0]]) {
+						self.remove_g1_continuous.insert(*handles);
+					}
+					self.add_g1_continuous.remove(handles);
+					self.add_g1_continuous.remove(&[handles[1], handles[0]]);
+				}
+			}
+			VectorModificationType::SetHandles { segment, handles } => {
+				self.segments.handle_primary.insert(*segment, handles[0]);
+				self.segments.handle_end.insert(*segment, handles[1]);
+			}
+			VectorModificationType::SetPrimaryHandle { segment, relative_position } => {
+				self.segments.handle_primary.insert(*segment, Some(*relative_position));
+			}
+			VectorModificationType::SetEndHandle { segment, relative_position } => {
+				self.segments.handle_end.insert(*segment, Some(*relative_position));
+			}
+			VectorModificationType::SetStartPoint { segment, id } => {
+				self.segments.start_point.insert(*segment, *id);
+			}
+			VectorModificationType::SetEndPoint { segment, id } => {
+				self.segments.end_point.insert(*segment, *id);
+			}
+
+			VectorModificationType::ApplyPointDelta { point, delta } => {
+				*self.points.delta.entry(*point).or_default() += *delta;
+			}
+			VectorModificationType::ApplyPrimaryDelta { segment, delta } => {
+				let position = self.segments.handle_primary.entry(*segment).or_default();
+				*position = Some(position.unwrap_or_default() + *delta);
+			}
+			VectorModificationType::ApplyEndDelta { segment, delta } => {
+				let position = self.segments.handle_end.entry(*segment).or_default();
+				*position = Some(position.unwrap_or_default() + *delta);
+			}
+		}
+	}
+
+	/// Create a new modification that will convert an empty [`VectorData`] into the target [`VectorData`].
+	pub fn create_from_vector(vector_data: &VectorData) -> Self {
+		Self {
+			points: PointModification::create_from_vector(vector_data),
+			segments: SegmentModification::create_from_vector(vector_data),
+			regions: RegionModification::create_from_vector(vector_data),
+			add_g1_continuous: vector_data.colinear_manipulators.iter().copied().collect(),
+			remove_g1_continuous: HashSet::new(),
+		}
+	}
+}
+
+impl Hash for VectorModification {
+	fn hash<H: std::hash::Hasher>(&self, state: &mut H) {
+		generate_uuid().hash(state)
+	}
+}
+
+/// A node that applies a procedural modification to some [`VectorData`].
+#[node_macro::node(category(""))]
+async fn path_modify(_ctx: impl Ctx, mut vector_data: VectorDataTable, modification: Box<VectorModification>) -> VectorDataTable {
+	if vector_data.is_empty() {
+		vector_data.push(Instance::default());
+	}
+	let vector_data_instance = vector_data.get_mut(0).expect("push should give one item");
+	modification.apply(vector_data_instance.instance);
+	if vector_data.len() > 1 {
+		warn!("The path modify ran on {} instances of vector data. Only the first can be modified.", vector_data.len());
+	}
+	vector_data
+}
+
+// Do we want to enforce that all serialized/deserialized hashmaps are a vec of tuples?
+// TODO: Eventually remove this document upgrade code
+pub fn serialize_hashmap<K, V, S, H>(hashmap: &HashMap<K, V, H>, serializer: S) -> Result<S::Ok, S::Error>
+where
+	K: Serialize + Eq + Hash,
+	V: Serialize,
+	S: Serializer,
+	H: BuildHasher,
+{
+	let mut seq = serializer.serialize_seq(Some(hashmap.len()))?;
+	for (key, value) in hashmap {
+		seq.serialize_element(&(key, value))?;
+	}
+	seq.end()
+}
+
+pub fn deserialize_hashmap<'de, K, V, D, H>(deserializer: D) -> Result<HashMap<K, V, H>, D::Error>
+where
+	K: Deserialize<'de> + Eq + Hash,
+	V: Deserialize<'de>,
+	D: Deserializer<'de>,
+	H: BuildHasher + Default,
+{
+	struct HashMapVisitor<K, V, H> {
+		#[allow(clippy::type_complexity)]
+		marker: std::marker::PhantomData<fn() -> HashMap<K, V, H>>,
+	}
+
+	impl<'de, K, V, H> Visitor<'de> for HashMapVisitor<K, V, H>
+	where
+		K: Deserialize<'de> + Eq + Hash,
+		V: Deserialize<'de>,
+		H: BuildHasher + Default,
+	{
+		type Value = HashMap<K, V, H>;
+
+		fn expecting(&self, formatter: &mut fmt::Formatter) -> fmt::Result {
+			formatter.write_str("a sequence of tuples")
+		}
+
+		fn visit_seq<A>(self, mut seq: A) -> Result<Self::Value, A::Error>
+		where
+			A: SeqAccess<'de>,
+		{
+			let mut hashmap = HashMap::default();
+			while let Some((key, value)) = seq.next_element()? {
+				hashmap.insert(key, value);
+			}
+			Ok(hashmap)
+		}
+	}
+
+	let visitor = HashMapVisitor { marker: std::marker::PhantomData };
+	deserializer.deserialize_seq(visitor)
+}
+
+pub struct AppendBezpath<'a> {
+	first_point: Option<Point>,
+	last_point: Option<Point>,
+	first_point_index: Option<usize>,
+	last_point_index: Option<usize>,
+	first_segment_id: Option<SegmentId>,
+	last_segment_id: Option<SegmentId>,
+	point_id: PointId,
+	segment_id: SegmentId,
+	vector_data: &'a mut VectorData,
+}
+
+impl<'a> AppendBezpath<'a> {
+	fn new(vector_data: &'a mut VectorData) -> Self {
+		Self {
+			first_point: None,
+			last_point: None,
+			first_point_index: None,
+			last_point_index: None,
+			first_segment_id: None,
+			last_segment_id: None,
+			point_id: vector_data.point_domain.next_id(),
+			segment_id: vector_data.segment_domain.next_id(),
+			vector_data,
+		}
+	}
+
+	fn append_segment_and_close_path(&mut self, point: Point, handle: BezierHandles) {
+		let handle = if self.first_point.unwrap() != point {
+			// If the first point is not the same as the last point of the path then we append the segment
+			// with given handle and point and then close the path with linear handle.
+			self.append_segment(point, handle);
+			BezierHandles::Linear
+		} else {
+			// if the endpoints are the same then we close the path with given handle.
+			handle
+		};
+
+		// Create a new segment.
+		let next_segment_id = self.segment_id.next_id();
+		self.vector_data
+			.segment_domain
+			.push(next_segment_id, self.last_point_index.unwrap(), self.first_point_index.unwrap(), handle, StrokeId::ZERO);
+
+		// Create a new region.
+		let next_region_id = self.vector_data.region_domain.next_id();
+		let first_segment_id = self.first_segment_id.unwrap_or(next_segment_id);
+		let last_segment_id = next_segment_id;
+
+		self.vector_data.region_domain.push(next_region_id, first_segment_id..=last_segment_id, FillId::ZERO);
+	}
+
+	fn append_segment(&mut self, end_point: Point, handle: BezierHandles) {
+		// Append the point.
+		let next_point_index = self.vector_data.point_domain.ids().len();
+		let next_point_id = self.point_id.next_id();
+
+		self.vector_data.point_domain.push(next_point_id, point_to_dvec2(end_point));
+
+		// Append the segment.
+		let next_segment_id = self.segment_id.next_id();
+		self.vector_data
+			.segment_domain
+			.push(next_segment_id, self.last_point_index.unwrap(), next_point_index, handle, StrokeId::ZERO);
+
+		// Update the states.
+		self.last_point = Some(end_point);
+		self.last_point_index = Some(next_point_index);
+
+		self.first_segment_id = Some(self.first_segment_id.unwrap_or(next_segment_id));
+		self.last_segment_id = Some(next_segment_id);
+	}
+
+	fn append_first_point(&mut self, point: Point) {
+		self.first_point = Some(point);
+		self.last_point = Some(point);
+
+		// Append the first point.
+		let next_point_index = self.vector_data.point_domain.ids().len();
+		self.vector_data.point_domain.push(self.point_id.next_id(), point_to_dvec2(point));
+
+		// Update the state.
+		self.first_point_index = Some(next_point_index);
+		self.last_point_index = Some(next_point_index);
+	}
+
+	fn reset(&mut self) {
+		self.first_point = None;
+		self.last_point = None;
+		self.first_point_index = None;
+		self.last_point_index = None;
+		self.first_segment_id = None;
+		self.last_segment_id = None;
+	}
+
+	pub fn append_bezpath(vector_data: &'a mut VectorData, bezpath: BezPath) {
+		let mut this = Self::new(vector_data);
+		let mut elements = bezpath.elements().iter().peekable();
+
+		while let Some(element) = elements.next() {
+			let close_path = elements.peek().is_some_and(|elm| **elm == PathEl::ClosePath);
+
+			match *element {
+				PathEl::MoveTo(point) => this.append_first_point(point),
+				PathEl::LineTo(point) => {
+					let handle = BezierHandles::Linear;
+					if close_path {
+						this.append_segment_and_close_path(point, handle);
+					} else {
+						this.append_segment(point, handle);
+					}
+				}
+				PathEl::QuadTo(point, point1) => {
+					let handle = BezierHandles::Quadratic { handle: point_to_dvec2(point) };
+					if close_path {
+						this.append_segment_and_close_path(point1, handle);
+					} else {
+						this.append_segment(point1, handle);
+					}
+				}
+				PathEl::CurveTo(point, point1, point2) => {
+					let handle = BezierHandles::Cubic {
+						handle_start: point_to_dvec2(point),
+						handle_end: point_to_dvec2(point1),
+					};
+
+					if close_path {
+						this.append_segment_and_close_path(point2, handle);
+					} else {
+						this.append_segment(point2, handle);
+					}
+				}
+				PathEl::ClosePath => {
+					// Already handled using `append_segment_and_close_path()` hence we reset state and continue.
+					this.reset();
+				}
+			}
+		}
+	}
+}
+
+pub trait VectorDataExt {
+	/// Appends a Kurbo BezPath to the vector data.
+	fn append_bezpath(&mut self, bezpath: BezPath);
+}
+
+impl VectorDataExt for VectorData {
+	fn append_bezpath(&mut self, bezpath: BezPath) {
+		AppendBezpath::append_bezpath(self, bezpath);
+	}
+}
+
+pub trait HandleExt {
+	/// Set the handle's position relative to the anchor which is the start anchor for the primary handle and end anchor for the end handle.
+	#[must_use]
+	fn set_relative_position(self, relative_position: DVec2) -> VectorModificationType;
+}
+
+impl HandleExt for HandleId {
+	fn set_relative_position(self, relative_position: DVec2) -> VectorModificationType {
+		let Self { ty, segment } = self;
+		match ty {
+			HandleType::Primary => VectorModificationType::SetPrimaryHandle { segment, relative_position },
+			HandleType::End => VectorModificationType::SetEndHandle { segment, relative_position },
+		}
+	}
+}
+
+#[cfg(test)]
+mod tests {
+	use super::*;
+
+	#[test]
+	fn modify_new() {
+		let vector_data = VectorData::from_subpaths(
+			[bezier_rs::Subpath::new_ellipse(DVec2::ZERO, DVec2::ONE), bezier_rs::Subpath::new_rect(DVec2::NEG_ONE, DVec2::ZERO)],
+			false,
+		);
+
+		let modify = VectorModification::create_from_vector(&vector_data);
+
+		let mut new = VectorData::default();
+		modify.apply(&mut new);
+		assert_eq!(vector_data, new);
+	}
+
+	#[test]
+	fn modify_existing() {
+		use bezier_rs::{Bezier, Subpath};
+		let subpaths = [
+			Subpath::new_ellipse(DVec2::ZERO, DVec2::ONE),
+			Subpath::new_rect(DVec2::NEG_ONE, DVec2::ZERO),
+			Subpath::from_beziers(
+				&[
+					Bezier::from_quadratic_dvec2(DVec2::new(0., 0.), DVec2::new(5., 10.), DVec2::new(10., 0.)),
+					Bezier::from_quadratic_dvec2(DVec2::new(10., 0.), DVec2::new(15., 10.), DVec2::new(20., 0.)),
+				],
+				false,
+			),
+		];
+		let mut vector_data = VectorData::from_subpaths(subpaths, false);
+
+		let mut modify_new = VectorModification::create_from_vector(&vector_data);
+		let mut modify_original = VectorModification::default();
+
+		for modification in [&mut modify_new, &mut modify_original] {
+			let point = vector_data.point_domain.ids()[0];
+			modification.modify(&VectorModificationType::ApplyPointDelta { point, delta: DVec2::X * 0.5 });
+			let point = vector_data.point_domain.ids()[9];
+			modification.modify(&VectorModificationType::ApplyPointDelta { point, delta: DVec2::X });
+		}
+
+		let mut new = VectorData::default();
+		modify_new.apply(&mut new);
+
+		modify_original.apply(&mut vector_data);
+
+		assert_eq!(vector_data, new);
+		assert_eq!(vector_data.point_domain.positions()[0], DVec2::X);
+		assert_eq!(vector_data.point_domain.positions()[9], DVec2::new(11., 0.));
+		assert_eq!(
+			vector_data.segment_bezier_iter().nth(8).unwrap().1,
+			Bezier::from_quadratic_dvec2(DVec2::new(0., 0.), DVec2::new(5., 10.), DVec2::new(11., 0.))
+		);
+		assert_eq!(
+			vector_data.segment_bezier_iter().nth(9).unwrap().1,
+			Bezier::from_quadratic_dvec2(DVec2::new(11., 0.), DVec2::new(16., 10.), DVec2::new(20., 0.))
+		);
+	}
+}