linebender / kurbo

A Rust library for manipulating curves
Apache License 2.0
727 stars 71 forks source link

`<QuadBez as ParamCurveNearest>::nearest` can panic #315

Closed richard-uk1 closed 1 year ago

richard-uk1 commented 1 year ago

I had a panic in our production (which uses kurbo for the analytic nearst point calculation). Looking at the code, the only way it can panic is if the cubic solver returned no solutions. This can happen if the cubic is in fact a quadratic or only a constant, which theoretically cannot happen but in practice can if a number is rounded to 0.

So I propose that we add a check so that if there are no roots then we check the ends of the curve. Because this branch will be very rarely hit in practice, iterating over a lot of curves shouldn't be slowed (due to branch prediction) although of course it depends on what the hardware does and so benchmarking is probably in order.

Example modification that ensures r_best is Some:

impl ParamCurveNearest for QuadBez {
    /// Find the nearest point, using analytical algorithm based on cubic root finding.
    fn nearest(&self, p: Point, _accuracy: f64) -> Nearest {
        fn eval_t(p: Point, t_best: &mut f64, r_best: &mut Option<f64>, t: f64, p0: Point) {
            let r = (p0 - p).hypot2();
            if r_best.map(|r_best| r < r_best).unwrap_or(true) {
                *r_best = Some(r);
                *t_best = t;
            }
        }
        fn try_t(
            q: &QuadBez,
            p: Point,
            t_best: &mut f64,
            r_best: &mut Option<f64>,
            t: f64,
        ) -> bool {
            if !(0.0..=1.0).contains(&t) {
                return true;
            }
            eval_t(p, t_best, r_best, t, q.eval(t));
            false
        }
        let d0 = self.p1 - self.p0;
        let d1 = self.p0.to_vec2() + self.p2.to_vec2() - 2.0 * self.p1.to_vec2();
        let d = self.p0 - p;
        let c0 = d.dot(d0);
        let c1 = 2.0 * d0.hypot2() + d.dot(d1);
        let c2 = 3.0 * d1.dot(d0);
        let c3 = d1.hypot2();
        let roots = solve_cubic(c0, c1, c2, c3);
        let mut r_best = None;
        let mut t_best = 0.0;
        let mut need_ends = false;
        // added
        if roots.is_empty() {
            need_ends = true;
        }
        for &t in &roots {
            need_ends |= try_t(self, p, &mut t_best, &mut r_best, t);
        }
        if need_ends {
            eval_t(p, &mut t_best, &mut r_best, 0.0, self.p0);
            eval_t(p, &mut t_best, &mut r_best, 1.0, self.p2);
        }

        Nearest {
            t: t_best,
            distance_sq: r_best.unwrap(),
        }
    }
}
richard-uk1 commented 1 year ago

I haven't found a failing case yet, and would like to before I make a PR.