Switch to a more robust intersection algorithm based on implicitization

This commit is contained in:
Patrick Walton 2017-09-20 16:18:55 -07:00
parent dafb18fd1d
commit 4ac11e9010
6 changed files with 114 additions and 88 deletions

View File

@ -35,6 +35,10 @@ const POINT_LABEL_FONT_SIZE: number = 12.0;
const POINT_LABEL_OFFSET: glmatrix.vec2 = glmatrix.vec2.fromValues(12.0, 12.0);
const POINT_RADIUS: number = 2.0;
const LIGHT_STROKE_STYLE: string = "rgb(192, 192, 192)";
const LINE_STROKE_STYLE: string = "rgb(0, 128, 0)";
const CURVE_STROKE_STYLE: string = "rgb(128, 0, 0)";
const BUILTIN_URIS = {
font: BUILTIN_FONT_URI,
svg: BUILTIN_SVG_URI,
@ -214,7 +218,7 @@ class MeshDebuggerView extends PathfinderView {
context.scale(this.camera.scale, this.camera.scale);
const invScaleFactor = window.devicePixelRatio / this.camera.scale;
context.font = `${12.0 * invScaleFactor}px ${POINT_LABEL_FONT}`;
context.font = `12px ${POINT_LABEL_FONT}`;
context.lineWidth = invScaleFactor;
const bQuads = new Uint32Array(meshes.bQuads);
@ -268,28 +272,42 @@ class MeshDebuggerView extends PathfinderView {
context.beginPath();
context.moveTo(upperLeftPosition[0], upperLeftPosition[1]);
if (upperControlPointPosition != null) {
context.strokeStyle = CURVE_STROKE_STYLE;
context.quadraticCurveTo(upperControlPointPosition[0],
upperControlPointPosition[1],
upperRightPosition[0],
upperRightPosition[1]);
} else {
context.strokeStyle = LINE_STROKE_STYLE;
context.lineTo(upperRightPosition[0], upperRightPosition[1]);
}
context.stroke();
context.strokeStyle = LIGHT_STROKE_STYLE;
context.beginPath();
context.moveTo(upperRightPosition[0], upperRightPosition[1]);
context.lineTo(lowerRightPosition[0], lowerRightPosition[1]);
context.stroke();
context.beginPath();
context.moveTo(lowerRightPosition[0], lowerRightPosition[1]);
if (lowerControlPointPosition != null) {
context.strokeStyle = CURVE_STROKE_STYLE;
context.quadraticCurveTo(lowerControlPointPosition[0],
lowerControlPointPosition[1],
lowerLeftPosition[0],
lowerLeftPosition[1]);
} else {
context.strokeStyle = LINE_STROKE_STYLE;
context.lineTo(lowerLeftPosition[0], lowerLeftPosition[1]);
}
context.stroke();
context.closePath();
context.strokeStyle = LIGHT_STROKE_STYLE;
context.beginPath();
context.moveTo(lowerLeftPosition[0], lowerLeftPosition[1]);
context.lineTo(upperLeftPosition[0], upperLeftPosition[1]);
context.stroke();
}
@ -323,9 +341,12 @@ function drawVertexIfNecessary(context: CanvasRenderingContext2D,
context.arc(position[0], position[1], POINT_RADIUS * invScaleFactor, 0, 2.0 * Math.PI);
context.fill();
context.save();
context.scale(invScaleFactor, invScaleFactor);
context.fillText("" + vertexIndex,
position[0] + POINT_LABEL_OFFSET[0] * invScaleFactor,
position[1] + POINT_LABEL_OFFSET[1] * invScaleFactor);
position[0] / invScaleFactor + POINT_LABEL_OFFSET[0],
position[1] / invScaleFactor + POINT_LABEL_OFFSET[1]);
context.restore();
}
function main() {

View File

@ -14,7 +14,6 @@ use geometry::{self, SubdividedQuadraticBezier};
use log::LogLevel;
use pathfinder_path_utils::PathBuffer;
use pathfinder_path_utils::curve::Curve;
use pathfinder_path_utils::intersection::Intersection;
use pathfinder_path_utils::line::Line;
use std::collections::BinaryHeap;
use std::cmp::Ordering;
@ -846,10 +845,7 @@ impl<'a> Partitioner<'a> {
&upper_right_endpoint_position);
let lower_line = Line::new(lower_left_vertex_position,
lower_right_endpoint_position);
Intersection::calculate(&upper_curve, &lower_line).map(|intersection| {
lower_line.sample(intersection.t_b)
})
upper_curve.intersect(&lower_line)
}
(u32::MAX, lower_control_point_vertex_index) => {
@ -860,10 +856,7 @@ impl<'a> Partitioner<'a> {
&lower_right_endpoint_position);
let upper_line = Line::new(upper_left_vertex_position,
upper_right_endpoint_position);
Intersection::calculate(&upper_line, &lower_curve).map(|intersection| {
upper_line.sample(intersection.t_a)
})
lower_curve.intersect(&upper_line)
}
(upper_control_point_vertex_index, lower_control_point_vertex_index) => {
@ -877,10 +870,7 @@ impl<'a> Partitioner<'a> {
let lower_curve = Curve::new(&lower_left_vertex_position,
&lower_control_point,
&lower_right_endpoint_position);
Intersection::calculate(&upper_curve, &lower_curve).map(|intersection| {
upper_curve.sample(intersection.t_a)
})
upper_curve.intersect(&lower_curve)
}
}
}

View File

@ -15,6 +15,7 @@ use euclid::Point2D;
use std::f32;
use PathSegment;
use intersection::{Intersect, Side};
pub struct Curve {
pub endpoints: [Point2D<f32>; 2],
@ -72,4 +73,9 @@ impl Curve {
None
}
}
#[inline]
pub fn intersect<T>(&self, other: &T) -> Option<Point2D<f32>> where T: Side {
<Curve as Intersect>::intersect(self, other)
}
}

View File

@ -10,95 +10,81 @@
//! Intersections of two segments.
use euclid::{Point2D, Rect};
use euclid::approxeq::ApproxEq;
use euclid::Point2D;
use curve::Curve;
use line::Line;
use lerp;
use {det2x2, det3x3, lerp};
const SUBDIVISION_TOLERANCE: f32 = 0.0001;
const MAX_SUBDIVISIONS: u32 = 1000;
pub struct Intersection {
pub t_a: f32,
pub t_b: f32,
pub trait Side {
fn side(&self, point: &Point2D<f32>) -> f32;
}
impl Intersection {
pub(crate) trait Intersect {
fn sample(&self, t: f32) -> Point2D<f32>;
/// Requires that any curves be monotonic. (See the `monotonic` module for that.)
///
/// This should work for line segments, but it is inefficient.
///
/// See T.W. Sederberg, "Computer Aided Geometric Design Course Notes" § 7.6.
pub fn calculate<A, B>(a: &A, b: &B) -> Option<Intersection> where A: Intersect, B: Intersect {
let (mut a_lower_t, mut a_upper_t) = (0.0, 1.0);
let (mut b_lower_t, mut b_upper_t) = (0.0, 1.0);
/// See T.W. Sederberg, "Computer Aided Geometric Design Course Notes" § 17.8.
fn intersect<T>(&self, other: &T) -> Option<Point2D<f32>> where T: Side {
let (mut t_min, mut t_max) = (0.0, 1.0);
let mut iteration = 0;
while t_max - t_min > f32::approx_epsilon() {
let (p_min, p_max) = (self.sample(t_min), self.sample(t_max));
for _ in 0..MAX_SUBDIVISIONS {
let a_lower_point = a.sample(a_lower_t);
let a_upper_point = a.sample(a_upper_t);
let b_lower_point = b.sample(b_lower_t);
let b_upper_point = b.sample(b_upper_t);
let side_min = other.side(&p_min).signum();
let side_max = other.side(&p_max).signum();
if iteration == 0 && side_min == side_max {
return None
}
let a_distance = (a_upper_point - a_lower_point).length();
let b_distance = (b_upper_point - b_lower_point).length();
let t_mid = lerp(t_min, t_max, 0.5);
let p_mid = self.sample(t_mid);
let side_mid = other.side(&p_mid).signum();
let need_to_subdivide_a = a_distance >= SUBDIVISION_TOLERANCE;
let need_to_subdivide_b = b_distance >= SUBDIVISION_TOLERANCE;
if !need_to_subdivide_b && !need_to_subdivide_a {
if side_mid == side_min {
t_min = t_mid
} else if side_mid == side_max {
t_max = t_mid
} else {
break
}
let a_rect;
if need_to_subdivide_a {
let a_middle_t = lerp(a_lower_t, a_upper_t, 0.5);
let a_middle_point = a.sample(a_middle_t);
let a_lower_rect =
Rect::from_points(&[a_lower_point, a_middle_point]);
let a_upper_rect =
Rect::from_points(&[a_middle_point, a_upper_point]);
let b_rect = Rect::from_points(&[b_lower_point, b_upper_point]);
if a_lower_rect.intersects(&b_rect) {
a_upper_t = a_middle_t;
a_rect = a_lower_rect;
} else if a_upper_rect.intersects(&b_rect) {
a_lower_t = a_middle_t;
a_rect = a_upper_rect;
} else {
return None
}
} else {
a_rect = Rect::from_points(&[a_lower_point, a_upper_point])
}
if need_to_subdivide_b {
let b_middle_t = lerp(b_lower_t, b_upper_t, 0.5);
let b_middle_point = b.sample(b_middle_t);
let b_lower_rect = Rect::from_points(&[b_lower_point, b_middle_point]);
let b_upper_rect = Rect::from_points(&[b_middle_point, b_upper_point]);
if b_lower_rect.intersects(&a_rect) {
b_upper_t = b_middle_t
} else if b_upper_rect.intersects(&a_rect) {
b_lower_t = b_middle_t
} else {
return None
}
}
iteration += 1
}
Some(Intersection {
t_a: lerp(a_lower_t, a_upper_t, 0.5),
t_b: lerp(b_lower_t, b_upper_t, 0.5),
})
Some(self.sample(lerp(t_min, t_max, 0.5)))
}
}
pub trait Intersect {
fn sample(&self, t: f32) -> Point2D<f32>;
impl Side for Line {
#[inline]
fn side(&self, point: &Point2D<f32>) -> f32 {
self.to_vector().cross(*point - self.endpoints[0])
}
}
impl Side for Curve {
/// See T.W. Sederberg, "Computer Aided Geometric Design Course Notes" § 17.6.1.
fn side(&self, point: &Point2D<f32>) -> f32 {
fn l(factor: f32, point: &Point2D<f32>, point_i: &Point2D<f32>, point_j: &Point2D<f32>)
-> f32 {
factor * det3x3(&[
point.x, point.y, 1.0,
point_i.x, point_i.y, 1.0,
point_j.x, point_j.y, 1.0,
])
}
let l20 = l(1.0 * 1.0, point, &self.endpoints[1], &self.endpoints[0]);
det2x2(&[
l(2.0 * 1.0, point, &self.endpoints[1], &self.control_point), l20,
l20, l(2.0 * 1.0, point, &self.control_point, &self.endpoints[0]),
])
}
}
impl Intersect for Line {

View File

@ -199,6 +199,17 @@ impl<I> Iterator for Transform2DPathStream<I> where I: Iterator<Item = PathSegme
}
#[inline]
pub fn lerp(a: f32, b: f32, t: f32) -> f32 {
pub(crate) fn lerp(a: f32, b: f32, t: f32) -> f32 {
a + (b - a) * t
}
#[inline]
pub(crate) fn det2x2(m: &[f32; 4]) -> f32 {
m[0] * m[3] - m[1] * m[2]
}
pub(crate) fn det3x3(m: &[f32; 9]) -> f32 {
m[0] * det2x2(&[m[4], m[5], m[7], m[8]]) -
m[1] * det2x2(&[m[3], m[5], m[6], m[8]]) +
m[2] * det2x2(&[m[3], m[4], m[6], m[7]])
}

View File

@ -10,7 +10,9 @@
//! Geometry utilities for straight line segments.
use euclid::Point2D;
use euclid::{Point2D, Vector2D};
use intersection::{Intersect, Side};
pub struct Line {
pub endpoints: [Point2D<f32>; 2],
@ -28,4 +30,14 @@ impl Line {
pub fn sample(&self, t: f32) -> Point2D<f32> {
self.endpoints[0].lerp(self.endpoints[1], t)
}
#[inline]
pub(crate) fn to_vector(&self) -> Vector2D<f32> {
self.endpoints[1] - self.endpoints[0]
}
#[inline]
pub fn intersect<T>(&self, other: &T) -> Option<Point2D<f32>> where T: Side {
<Line as Intersect>::intersect(self, other)
}
}