Fall back to Newton's method instead of dividing by zero when finding critical
points during curve monotonic conversion. Closes #146.
This commit is contained in:
parent
5568e6f64f
commit
3310b15826
|
@ -503,6 +503,8 @@ impl Contour {
|
||||||
}
|
}
|
||||||
|
|
||||||
fn make_monotonic(&mut self) {
|
fn make_monotonic(&mut self) {
|
||||||
|
debug!("--- make_monotonic() ---");
|
||||||
|
|
||||||
let contour = self.take();
|
let contour = self.take();
|
||||||
self.bounds = contour.bounds;
|
self.bounds = contour.bounds;
|
||||||
|
|
||||||
|
@ -550,6 +552,8 @@ impl Contour {
|
||||||
}
|
}
|
||||||
|
|
||||||
fn handle_cubic(contour: &mut Contour, segment: Segment) {
|
fn handle_cubic(contour: &mut Contour, segment: Segment) {
|
||||||
|
debug!("handle_cubic({:?})", segment);
|
||||||
|
|
||||||
match segment.as_cubic_segment().y_extrema() {
|
match segment.as_cubic_segment().y_extrema() {
|
||||||
(Some(t0), Some(t1)) => {
|
(Some(t0), Some(t1)) => {
|
||||||
let (segments_01, segment_2) = segment.as_cubic_segment().split(t1);
|
let (segments_01, segment_2) = segment.as_cubic_segment().split(t1);
|
||||||
|
|
|
@ -12,8 +12,11 @@
|
||||||
|
|
||||||
use crate::basic::line_segment::LineSegmentF32;
|
use crate::basic::line_segment::LineSegmentF32;
|
||||||
use crate::basic::point::Point2DF32;
|
use crate::basic::point::Point2DF32;
|
||||||
|
use crate::util::{self, EPSILON};
|
||||||
use pathfinder_simd::default::F32x4;
|
use pathfinder_simd::default::F32x4;
|
||||||
|
|
||||||
|
const MAX_NEWTON_ITERATIONS: u32 = 32;
|
||||||
|
|
||||||
#[derive(Clone, Copy, Debug, PartialEq)]
|
#[derive(Clone, Copy, Debug, PartialEq)]
|
||||||
pub struct Segment {
|
pub struct Segment {
|
||||||
pub baseline: LineSegmentF32,
|
pub baseline: LineSegmentF32,
|
||||||
|
@ -302,12 +305,35 @@ impl<'s> CubicSegment<'s> {
|
||||||
let pxv0v1v2 = p0p1p2p3 - pxp0p1p2;
|
let pxv0v1v2 = p0p1p2p3 - pxp0p1p2;
|
||||||
let (v0, v1, v2) = (pxv0v1v2[1], pxv0v1v2[2], pxv0v1v2[3]);
|
let (v0, v1, v2) = (pxv0v1v2[1], pxv0v1v2[2], pxv0v1v2[3]);
|
||||||
|
|
||||||
|
let (t0, t1);
|
||||||
let (v0_to_v1, v2_to_v1) = (v0 - v1, v2 - v1);
|
let (v0_to_v1, v2_to_v1) = (v0 - v1, v2 - v1);
|
||||||
let discrim = f32::sqrt(v1 * v1 - v0 * v2);
|
let denom = v0_to_v1 + v2_to_v1;
|
||||||
let denom = 1.0 / (v0_to_v1 + v2_to_v1);
|
|
||||||
|
|
||||||
let t0 = (v0_to_v1 + discrim) * denom;
|
if util::approx_eq(denom, 0.0) {
|
||||||
let t1 = (v0_to_v1 - discrim) * denom;
|
// Let's not divide by zero (issue #146). Fall back to Newton's method.
|
||||||
|
// FIXME(pcwalton): Can we have two roots here?
|
||||||
|
let mut t = 0.5;
|
||||||
|
for _ in 0..MAX_NEWTON_ITERATIONS {
|
||||||
|
let dydt = 3.0 * ((denom * t - v0_to_v1 - v0_to_v1) * t + v0);
|
||||||
|
if f32::abs(dydt) <= EPSILON {
|
||||||
|
break
|
||||||
|
}
|
||||||
|
let d2ydt2 = 6.0 * (denom * t - v0_to_v1);
|
||||||
|
t -= dydt / d2ydt2;
|
||||||
|
}
|
||||||
|
t0 = t;
|
||||||
|
t1 = 0.0;
|
||||||
|
debug!("... t=(newton) {}", t);
|
||||||
|
} else {
|
||||||
|
// Algebraically compute the values for t.
|
||||||
|
let discrim = f32::sqrt(v1 * v1 - v0 * v2);
|
||||||
|
let denom_recip = 1.0 / denom;
|
||||||
|
|
||||||
|
t0 = (v0_to_v1 + discrim) * denom_recip;
|
||||||
|
t1 = (v0_to_v1 - discrim) * denom_recip;
|
||||||
|
|
||||||
|
debug!("... t=({} +/- {})/{} t0={} t1={}", v0_to_v1, discrim, denom, t0, t1);
|
||||||
|
}
|
||||||
|
|
||||||
return match (
|
return match (
|
||||||
t0 > EPSILON && t0 < 1.0 - EPSILON,
|
t0 > EPSILON && t0 < 1.0 - EPSILON,
|
||||||
|
@ -318,8 +344,6 @@ impl<'s> CubicSegment<'s> {
|
||||||
(false, true) => (Some(t1), None),
|
(false, true) => (Some(t1), None),
|
||||||
(true, true) => (Some(f32::min(t0, t1)), Some(f32::max(t0, t1))),
|
(true, true) => (Some(f32::min(t0, t1)), Some(f32::max(t0, t1))),
|
||||||
};
|
};
|
||||||
|
|
||||||
const EPSILON: f32 = 0.001;
|
|
||||||
}
|
}
|
||||||
|
|
||||||
#[inline]
|
#[inline]
|
||||||
|
|
|
@ -12,6 +12,14 @@
|
||||||
|
|
||||||
use std::f32;
|
use std::f32;
|
||||||
|
|
||||||
|
pub const EPSILON: f32 = 0.001;
|
||||||
|
|
||||||
|
/// Approximate equality.
|
||||||
|
#[inline]
|
||||||
|
pub fn approx_eq(a: f32, b: f32) -> bool {
|
||||||
|
f32::abs(a - b) <= EPSILON
|
||||||
|
}
|
||||||
|
|
||||||
/// Linear interpolation.
|
/// Linear interpolation.
|
||||||
#[inline]
|
#[inline]
|
||||||
pub fn lerp(a: f32, b: f32, t: f32) -> f32 {
|
pub fn lerp(a: f32, b: f32, t: f32) -> f32 {
|
||||||
|
|
Loading…
Reference in New Issue