Improve snapping performance (#3067)

* Use lyon_geom for intersection calculation of bezier segments

* Implement approximate nearest point calculation and loosen bounding boxes

* Add algorithm explanation

* Update editor/src/messages/tool/common_functionality/snapping/layer_snapper.rs

Co-authored-by: Keavon Chambers <keavon@keavon.com>

---------

Co-authored-by: Keavon Chambers <keavon@keavon.com>

Author: TrueDoctor

Cargo.lock

@@ -159,6 +159,7 @@ syn = { version = "2.0", default-features = false, features = [
 	"proc-macro",
 ] }
 kurbo = { version = "0.11.0", features = ["serde"] }
+lyon_geom = "1.0"
 petgraph = { version = "0.7.1", default-features = false, features = [
 	"graphmap",
 ] }

Cargo.toml

@@ -161,6 +161,17 @@ impl DocumentMetadata {
 			.reduce(Quad::combine_bounds)
 	}
 
+	/// Get the loose bounding box of the click target of the specified layer in the specified transform space
+	pub fn loose_bounding_box_with_transform(&self, layer: LayerNodeIdentifier, transform: DAffine2) -> Option<[DVec2; 2]> {
+		self.click_targets(layer)?
+			.iter()
+			.filter_map(|click_target| match click_target.target_type() {
+				ClickTargetType::Subpath(subpath) => subpath.loose_bounding_box_with_transform(transform),
+				ClickTargetType::FreePoint(_) => click_target.bounding_box_with_transform(transform),
+			})
+			.reduce(Quad::combine_bounds)
+	}
+
 	/// Calculate the corners of the bounding box but with a nonzero size.
 	///
 	/// If the layer bounds are `0` in either axis then they are changed to be `1`.

editor/src/messages/portfolio/document/utility_types/document_metadata.rs

@@ -332,7 +332,8 @@ impl SnapManager {
 			}
 			return;
 		}
-		let Some(bounds) = document.metadata().bounding_box_with_transform(layer, DAffine2::IDENTITY) else {
+		// We use a loose bounding box here since these are potential candidates which will be filtered later anyway
+		let Some(bounds) = document.metadata().loose_bounding_box_with_transform(layer, DAffine2::IDENTITY) else {
 			return;
 		};
 		let layer_bounds = document.metadata().transform_to_document(layer) * Quad::from_box(bounds);

editor/src/messages/tool/common_functionality/snapping.rs

@@ -3,7 +3,7 @@ use crate::consts::HIDE_HANDLE_DISTANCE;
 use crate::messages::portfolio::document::utility_types::document_metadata::LayerNodeIdentifier;
 use crate::messages::portfolio::document::utility_types::misc::*;
 use crate::messages::prelude::*;
-use glam::{DAffine2, DVec2};
+use glam::{DAffine2, DVec2, FloatExt};
 use graphene_std::math::math_ext::QuadExt;
 use graphene_std::renderer::Quad;
 use graphene_std::subpath::pathseg_points;
@@ -13,7 +13,7 @@ use graphene_std::vector::algorithms::bezpath_algorithms::{pathseg_normals_to_po
 use graphene_std::vector::algorithms::intersection::filtered_segment_intersections;
 use graphene_std::vector::misc::dvec2_to_point;
 use graphene_std::vector::misc::point_to_dvec2;
-use kurbo::{Affine, DEFAULT_ACCURACY, Nearest, ParamCurve, ParamCurveNearest, PathSeg};
+use kurbo::{Affine, ParamCurve, PathSeg};
 
 #[derive(Clone, Debug, Default)]
 pub struct LayerSnapper {
@@ -107,9 +107,11 @@ impl LayerSnapper {
 			if path.document_curve.start().distance_squared(path.document_curve.end()) < tolerance * tolerance * 2. {
 				continue;
 			}
-			let Nearest { distance_sq, t } = path.document_curve.nearest(dvec2_to_point(point.document_point), DEFAULT_ACCURACY);
-			let snapped_point_document = point_to_dvec2(path.document_curve.eval(t));
-			let distance = distance_sq.sqrt();
+			let Some((distance_squared, closest)) = path.approx_nearest_point(point.document_point, 10) else {
+				continue;
+			};
+			let snapped_point_document = point_to_dvec2(closest);
+			let distance = distance_squared.sqrt();
 
 			if distance < tolerance {
 				snap_results.curves.push(SnappedCurve {
@@ -322,6 +324,99 @@ struct SnapCandidatePath {
 	bounds: Option<Quad>,
 }
 
+impl SnapCandidatePath {
+	/// Calculates the point on the curve which lies closest to `point`.
+	///
+	/// ## Algorithm:
+	/// 1. We first perform a coarse scan of the path segment to find the most promising starting point.
+	/// 2. Afterwards we refine this point by performing a binary search to either side assuming that the segment contains at most one extremal point.
+	/// 3. The smaller of the two resulting distances is returned.
+	///
+	/// ## Visualization:
+	/// ```text
+	///        Query Point (×)
+	///             ×
+	///            /|\
+	///           / | \  distance checks
+	///          /  |  \
+	///         v   v   v
+	///     ●---●---●---●---●  <- Curve with coarse scan points
+	///     0  0.25 0.5 0.75 1  (parameter t values)
+	///         ^       ^
+	///         |   |   |
+	///        min mid max
+	///    Find closest scan point
+	///
+	///    Refine left region using binary search:
+	///
+	///            ●------●------●
+	///           0.25  0.375   0.5
+	///
+	///    Result:           |     (=0.4)
+	///    And the right region:
+	///
+	///            ●------●------●
+	///           0.5   0.625   0.75
+	///    Result: |               (=0.5)
+	///
+	///    The t value with minimal dist is thus 0.4
+	///    Return: (dist_closest, point_on_curve)
+	/// ```
+	pub fn approx_nearest_point(&self, point: DVec2, lut_steps: usize) -> Option<(f64, kurbo::Point)> {
+		let point = dvec2_to_point(point);
+
+		let time_values = (0..lut_steps).map(|x| x as f64 / lut_steps as f64);
+		let points = time_values.map(|t| (t, self.document_curve.eval(t)));
+		let points_with_distances = points.map(|(t, p)| (t, p.distance_squared(point), p));
+		let (t, _, _) = points_with_distances.min_by(|(_, a, _), (_, b, _)| a.partial_cmp(b).unwrap_or(std::cmp::Ordering::Equal))?;
+
+		let min_t = (t - (lut_steps as f64).recip()).max(0.);
+		let max_t = (t + (lut_steps as f64).recip()).min(1.);
+		let left = self.refine_nearest_point(point, min_t, t);
+		let right = self.refine_nearest_point(point, t, max_t);
+
+		if left.0 < right.0 { Some(left) } else { Some(right) }
+	}
+
+	/// Refines the nearest point search within a given parameter range using binary search.
+	///
+	/// This method performs iterative refinement by:
+	/// 1. Evaluating the midpoint of the current parameter range
+	/// 2. Comparing distances at the endpoints and midpoint
+	/// 3. Narrowing the search range to the side with the shorter distance
+	/// 4. Continuing until convergence (when the range becomes very small)
+	///
+	/// Returns a tuple of (parameter_t, closest_point) where parameter_t is in the range [min_t, max_t].
+	fn refine_nearest_point(&self, point: kurbo::Point, mut min_t: f64, mut max_t: f64) -> (f64, kurbo::Point) {
+		let mut min_dist = self.document_curve.eval(min_t).distance_squared(point);
+		let mut max_dist = self.document_curve.eval(max_t).distance_squared(point);
+		let mut mid_t = max_t.lerp(min_t, 0.5);
+		let mut mid_point = self.document_curve.eval(mid_t);
+		let mut mid_dist = mid_point.distance_squared(point);
+
+		for _ in 0..10 {
+			if (min_dist - max_dist).abs() < 1e-3 {
+				return (mid_dist, mid_point);
+			}
+			if mid_dist > min_dist && mid_dist > max_dist {
+				return (mid_dist, mid_point);
+			}
+			if max_dist > min_dist {
+				max_t = mid_t;
+				max_dist = mid_dist;
+			} else {
+				min_t = mid_t;
+				min_dist = mid_dist;
+			}
+			mid_t = max_t.lerp(min_t, 0.5);
+			mid_point = self.document_curve.eval(mid_t);
+			mid_dist = mid_point.distance_squared(point);
+		}
+
+		(mid_dist, mid_point)
+	}
+}
+
 #[derive(Clone, Debug, Default)]
 pub struct SnapCandidatePoint {
 	pub document_point: DVec2,

editor/src/messages/tool/common_functionality/snapping/layer_snapper.rs

@@ -35,6 +35,7 @@ tinyvec = { workspace = true }
 parley = { workspace = true }
 skrifa = { workspace = true }
 kurbo = { workspace = true }
+lyon_geom = { workspace = true }
 log = { workspace = true }
 base64 = { workspace = true }
 poly-cool = { workspace = true }

node-graph/gcore/Cargo.toml

@@ -1,5 +1,24 @@
 use super::contants::MIN_SEPARATION_VALUE;
 use kurbo::{BezPath, DEFAULT_ACCURACY, ParamCurve, PathSeg, Shape};
+use lyon_geom::{CubicBezierSegment, Point};
+
+/// Converts a kurbo cubic bezier to a lyon_geom CubicBezierSegment
+fn kurbo_cubic_to_lyon(cubic: kurbo::CubicBez) -> CubicBezierSegment<f64> {
+	CubicBezierSegment {
+		from: Point::new(cubic.p0.x, cubic.p0.y),
+		ctrl1: Point::new(cubic.p1.x, cubic.p1.y),
+		ctrl2: Point::new(cubic.p2.x, cubic.p2.y),
+		to: Point::new(cubic.p3.x, cubic.p3.y),
+	}
+}
+
+/// Fast cubic-cubic intersection using lyon_geom's analytical approach
+fn cubic_cubic_intersections_lyon(cubic1: kurbo::CubicBez, cubic2: kurbo::CubicBez) -> Vec<(f64, f64)> {
+	let lyon_cubic1 = kurbo_cubic_to_lyon(cubic1);
+	let lyon_cubic2 = kurbo_cubic_to_lyon(cubic2);
+
+	lyon_cubic1.cubic_intersections_t(&lyon_cubic2).to_vec()
+}
 
 /// Calculates the intersection points the bezpath has with a given segment and returns a list of `(usize, f64)` tuples,
 /// where the `usize` represents the index of the segment in the bezpath, and the `f64` represents the `t`-value local to
@@ -37,6 +56,8 @@ pub fn segment_intersections(segment1: PathSeg, segment2: PathSeg, accuracy: Opt
 	match (segment1, segment2) {
 		(PathSeg::Line(line), segment2) => segment2.intersect_line(line).iter().map(|i| (i.line_t, i.segment_t)).collect(),
 		(segment1, PathSeg::Line(line)) => segment1.intersect_line(line).iter().map(|i| (i.segment_t, i.line_t)).collect(),
+		// Fast path for cubic-cubic intersections using lyon_geom
+		(PathSeg::Cubic(cubic1), PathSeg::Cubic(cubic2)) => cubic_cubic_intersections_lyon(cubic1, cubic2),
 		(segment1, segment2) => {
 			let mut intersections = Vec::new();
 			segment_intersections_inner(segment1, 0., 1., segment2, 0., 1., accuracy, &mut intersections);
@@ -51,6 +72,21 @@ pub fn subsegment_intersections(segment1: PathSeg, min_t1: f64, max_t1: f64, seg
 	match (segment1, segment2) {
 		(PathSeg::Line(line), segment2) => segment2.intersect_line(line).iter().map(|i| (i.line_t, i.segment_t)).collect(),
 		(segment1, PathSeg::Line(line)) => segment1.intersect_line(line).iter().map(|i| (i.segment_t, i.line_t)).collect(),
+		// Fast path for cubic-cubic intersections using lyon_geom with subsegment parameters
+		(PathSeg::Cubic(cubic1), PathSeg::Cubic(cubic2)) => {
+			let sub_cubic1 = cubic1.subsegment(min_t1..max_t1);
+			let sub_cubic2 = cubic2.subsegment(min_t2..max_t2);
+
+			cubic_cubic_intersections_lyon(sub_cubic1, sub_cubic2)
+				.into_iter()
+				// Convert subsegment t-values back to original segment t-values
+				.map(|(t1, t2)| {
+					let original_t1 = min_t1 + t1 * (max_t1 - min_t1);
+					let original_t2 = min_t2 + t2 * (max_t2 - min_t2);
+					(original_t1, original_t2)
+				})
+				.collect()
+		}
 		(segment1, segment2) => {
 			let mut intersections = Vec::new();
 			segment_intersections_inner(segment1, min_t1, max_t1, segment2, min_t2, max_t2, accuracy, &mut intersections);
@@ -59,12 +95,33 @@ pub fn subsegment_intersections(segment1: PathSeg, min_t1: f64, max_t1: f64, seg
 	}
 }
 
+fn approx_bounding_box(path_seg: PathSeg) -> kurbo::Rect {
+	use kurbo::Rect;
+	match path_seg {
+		PathSeg::Line(line) => kurbo::Rect::from_points(line.p0, line.p1),
+		PathSeg::Quad(quad_bez) => {
+			let r1 = Rect::from_points(quad_bez.p0, quad_bez.p1);
+			let r2 = Rect::from_points(quad_bez.p1, quad_bez.p2);
+			r1.union(r2)
+		}
+		PathSeg::Cubic(cubic_bez) => {
+			let r1 = Rect::from_points(cubic_bez.p0, cubic_bez.p1);
+			let r2 = Rect::from_points(cubic_bez.p2, cubic_bez.p3);
+			r1.union(r2)
+		}
+	}
+}
+
 /// Implements [https://pomax.github.io/bezierinfo/#curveintersection] to find intersection between two Bezier segments
 /// by splitting the segment recursively until the size of the subsegment's bounding box is smaller than the accuracy.
 #[allow(clippy::too_many_arguments)]
 fn segment_intersections_inner(segment1: PathSeg, min_t1: f64, max_t1: f64, segment2: PathSeg, min_t2: f64, max_t2: f64, accuracy: f64, intersections: &mut Vec<(f64, f64)>) {
-	let bbox1 = segment1.subsegment(min_t1..max_t1).bounding_box();
-	let bbox2 = segment2.subsegment(min_t2..max_t2).bounding_box();
+	let bbox1 = approx_bounding_box(segment1.subsegment(min_t1..max_t1));
+	let bbox2 = approx_bounding_box(segment2.subsegment(min_t2..max_t2));
+
+	if intersections.len() > 50 {
+		return;
+	}
 
 	let mid_t1 = (min_t1 + max_t1) / 2.;
 	let mid_t2 = (min_t2 + max_t2) / 2.;