views:

1558

answers:

5
+6  Q: 

Bezier clipping

I'm trying to find/make an algorithm to compute the intersection (a new filled object) of two arbitrary filled 2D objects. The objects are defined using either lines or cubic beziers and may have holes or self-intersect. I'm aware of several existing algorithms doing the same with polygons, listed here. However, I'd like to support beziers without subdividing them into polygons, and the output should have roughly the same control points as the input in areas where there are no intersections.

This is for an interactive program to do some CSG but the clipping doesn't need to be real-time. I've searched for a while but haven't found good starting points.

+1  A: 

There are a number of academic research papers on doing bezier clipping:

http://www.andrew.cmu.edu/user/sowen/abstracts/Se306.html

http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.61.6669

http://www.dm.unibo.it/~casciola/html/research_rr.html

I recommend the interval methods because as you describe, you don't have to divide down to polygons, and you can get guaranteed results as well as define your own arbitrary precision for the resultset. For more information on interval rendering, you may also refer to http://www.sunfishstudio.com

devinmoore
+3  A: 

It's just a root-finding problem in the end.

Consider the case that you just want to check if two bezier curves intersect in the x-axis only. The equation simply becomes:

0 = bezier1(t1).x - bezier2(t2).x;

That's one equation and two unkowns (t1 and t2). Now add the second dimension:

0 = bezier1(t1).x - bezier2(t2).x;
0 = bezier1(t1).y - bezier2(t2).y;

That's two equations and two unkowns, e.g. a system that can be solved.

The bezier equations itself can be expressed several ways. I favor the polynomial form for a reason I'll give you later. The polynomial form is just:

  bezier(t).x = a.x * t^3 + b.x * t^2 + c.x * t + d.x;

a,b,c and d are the polynomial coefficients. You can get them with a bit of math from the control points (if you need help with this part write me a comment, I have that part in my source-code somewhere).

A polynomial of degree 3 may have at most three roots (or zeros). Subtracting two polynomials of the same degree yields another polynomial with that degree. In other words: The degree does not grow, and that's great because finding roots of a polynomial of that degree can be done with a closed form solution, so you don't have to use iterative algorithms. Code to do so can be found en masse on the net.

Since your beziers are only defined in the interval 0 <= t <= 1 you should ignore all roots outside this range.

Next step: You have an intersection between the curves only If you find the same roots for your X and Y polynom. All others can be ignored. At this point you split your curves and insert an intersection point. (If you need help on this part again write a comment).

Here's the reason why I like to do these things in the polynomial form so much: Finding out if two polynomials have a root at the same position at all is cheaper to calculate the roots first and search for a duplicate. To do so you build the sylvester determinant of the two polynomials and calculate the determinant from it.

If this determinant is zero, the two polynomials don't only have a root and you can examine the next pair of possible intersecting beziers.

The sylvester matrix for two polynomials of degree three looks like this (assume a,b,c and d are the polynomial coefficients again):

  |  a.x  b.x  c.x  d.x    0    0  |
  |    0  a.x  b.x  c.x  d.x    0  |
  |    0    0  a.x  b.x  c.x  d.x  |
  |  a.y  b.y  c.y  d.y    0    0  |
  |    0  a.y  b.y  c.y  d.y    0  |
  |    0    0  a.y  b.y  c.y  d.y  |

That seems like a hell of a matrix (and it is), but calculating the determinant of it is cheaper than root-finding because of it's structure.

Since the first column contains only two non-zero elements you can get the determinant via cofactor expansion of exactly these elements.

Btw - if yoou work with only quadratic beziers (e.g. those that only have a single control-point) things become a bit easier of course.


Useful links:

Sylvester Matrix on Wikipedia

Cofactor Expansion on Wikipedia

Nils Pipenbrinck
Thanks for the explanation! Have you implemented this already? I'm puzzled by this root finding approach as elsewhere I've seen mentioned that it's a 9th degree polynomial to intersect 2 cubic beziers. Drawing 2 loops on paper that intersect each other...
jjrv
... (see 2nd example in http://truetex.com/bezint.htm) clearly produces 6 intersections between the curves and 2 self-intersections. I'd guess to find the 6 roots you need at least a 6th degree polynomial? Can you squeeze 6 roots out of your solution?
jjrv
Excuse the digging, but when you say "If this determinant is zero, the two polynomials don't only have a root and you can examine the next pair of possible intersecting beziers.", do you mean that when the determinant is zero there *is* an intersection or there *isn't*?
zneak
Hi zneak: If the determinant is zero there is is a root in X and Y at *exactly the same value of t, so there is an intersection of the two curves. (the t may not be in the interval 0..1 though). If the determinant is <> zero you can be sure that the curves don't intersect anywhere.
Nils Pipenbrinck
A: 

I wrote code to do this a long, long time ago. The project I was working on defined 2D objects using piecewise Bezier boundaries that were generated as PostScipt paths.

The approach I used was:

Let curves p, q, be defined by Bezier control points. Do they intersect?

Compute the bounding boxes of the control points. If they don't overlap, then the two curves don't intersect. Otherwise:

p.x(t), p.y(t), q.x(u), q.y(u) are cubic polynomials on 0 <= t,u <= 1.0. The distance squared (p.x - q.x) ** 2 + (p.y - q.y) ** 2 is a polynomial on (t,u). Use Newton-Raphson to try and solve that for zero. Discard any solutions outside 0 <= t,u <= 1.0

N-R may or may not converge. The curves might not intersect, or N-R can just blow up when the two curves are nearly parallel. So cut off N-R if it's not converging after after some arbitrary number of iterations. This can be a fairly small number.

If N-R doesn't converge on a solution, split one curve (say, p) into two curves pa, pb at t = 0.5. This is easy, it's just computing midpoints, as the linked article shows. Then recursively test (q, pa) and (q, pb) for intersections. (Note that in the next layer of recursion that q has become p, so that p and q are alternately split on each ply of the recursion.)

Most of the recursive calls return quickly because the bounding boxes are non-overlapping.

You will have to cut off the recursion at some arbitrary depth, to handle weird cases where the two curves are parallel and don't quite touch, but the distance between them is arbitrarily small -- perhaps only 1 ULP of difference.

When you do find an intersection, you're not done, because cubic curves can have multiple crossings. So you have to split each curve at the intersecting point, and recursively check for more interections between (pa, qa), (pa, qb), (pb, qa), (pb, qb).

Die in Sente
A: 

I know I'm at risk of being redundant, but I was investigating the same issue and found a solution that I'd read in academic papers but hadn't found a working solution for.

You can rewrite the bezier curves as a set of two bi-variate cubic equations like this:

  • ∆x = ax₁*t₁^3 + bx₁*t₁^2 + cx₁*t₁ + dx₁ - ax₂*t₂^3 - bx₂*t₂^2 - cx₂*t₂ - dx₂
  • ∆y = ay₁*t₁^3 + by₁*t₁^2 + cy₁*t₁ + dy₁ - ay₂*t₂^3 - by₂*t₂^2 - cy₂*t₂ - dy₂

Obviously, the curves intersect for values of (t₁,t₂) where ∆x = ∆y = 0. Unfortunately, it's complicated by the fact that it is difficult to find roots in two dimensions, and approximate approaches tend to (as another writer put it) blow up.

But if you're using integers or rational numbers for your control points, then you can use Groebner bases to rewrite your bi-variate order-3 polynomials into a (possibly-up-to-order-9-thus-your-nine-possible-intersections) monovariate polynomial. After that you just need to find your roots (for, say t₂) in one dimension, and plug your results back into one of your original equations to find the other dimension.

Burchburger has a layman-friendly introduction to Groebner Bases called "Gröbner Bases: A Short Introduction for Systems Theorists" that was very helpful for me. Google it. The other paper that was helpful was one called "Fast, precise flattening of cubic Bézier path and offset curves" by TF Hain, which has lots of utility equations for bezier curves, including how to find the polynomial coefficients for the x and y equations.

As for whether the Bezier clipping will help with this particular method, I doubt it, but bezier clipping is a method for narrowing down where intersections might be, not for finding a final (though possibly approximate) answer of where it is. A lot of time with this method will be spent in finding the mono-variate equation, and that task doesn't get any easier with clipping. Finding the roots is by comparison trivial.

However, one of the advantages of this method is that it doesn't depend on recursively subdividing the curve, and the problem becomes a simple one-dimensional root-finding problem, which is not easy, but well documented. The major disadvantage is that computing Groebner bases is costly and becomes very unwieldy if you're dealing with floating point polynomials or using higher order Bezier curves.

If you want some finished code in Haskell to find the intersections, let me know.

Adam
A: 

I found the following publication to be the best of information regarding Bezier Clipping:

T. W. Sederberg, BYU, Computer Aided Geometric Design Course Notes

Chapter 7 that talks about Curve Intersection is available online. It outlines 4 different approaches to find intersections and describes Bezier Clipping in detail:

http://www.tsplines.com/technology/edu/CurveIntersection.pdf

lehni