4

I have an ascending cubic Bézier curve. ($x_0 \leq x_1 \leq x_2 \leq x_3$) Considering this property, there is always one and only one $y$ value per $x$ value.

The point ($x, y$) along the curve is determined by the following equation, from Wikipedia:

$B(t) = (1-t)^3 P_0 + 3t (1-t)^2 P_1 + 3t^2 (1-t) P_2 + t^3 P_3$

However, this equation takes $t$ as a parameter, which is the progression along the curve, or, as Wikipedia explains,

The $t$ in the function for a linear Bézier curve can be thought of as describing how far $B(t)$ is from $P_0$ to $P_n$. For example when $t=0.25$, $B(t)$ is one quarter of the way from point $P_0$ to $P_n$.

enter image description here

I would like to use the $B(t)$ function with a $x$ parameter instead. Something like $B(x)$, or anything that can convert a given $x$ to a $t$.

Edit: I realized that I could solve $B(t)$ for $t$:

$B(t) = (1-t)^3 P_0 + 3t (1-t)^2 P_1 + 3t^2 (1-t) P_2 + t^3 P_3$

$x(t) = (1-t)^3 x_0 + 3t (1-t)^2 x_1 + 3t^2 (1-t) x_2 + t^3 x_3$

$x(t) = (-x_0+3x_1-3x_2+x_3)t^3 + (3x_0-6x_1+3x_2)t^2 + (-3x_0+3x_1)t + x_0$

$t(x) = ...?$

Edit 2: With the Newton-Raphson method, I managed to reach my goal. I'm attaching a screenshot.

Newton-Raphson curve

The red dots represent my interpolation (and surprisingly, extrapolation as well. This was an unexpected but great surprise!).

The gray dots represent a small optimization I made to the method: the number of iterations based on the difficulty to approximate. I increase the precision in an invertly proportional relation to the absolute value of the derivative: $$p\ \alpha\ \frac{1}{|B'(t_n)|}$$

I have to tweak the parameters, but so far, this method is efficient and pretty safe, even at a nearly straight slope. Thanks for the help!

P.S.: The more explicit solution, without the precision adjustment (which is unfortunately only well represented in procedural code):

$$t\approx x-\sum_{n=0}^p (\frac{B(t)_x-x}{B'(t)_x})$$

(This takes advantage of the fact that $-x$ is nullified in the derivative Also note that $t$ is reassigned at each step of the sum).

P.P.S.: The accompanying code:

p = 5;
p_max = 100;

t = x;

for (n = 0; n < min(p, p_max); n++)
{
    d = dB(t).x;

    if (p < p_max) p += min(0.02 / d^2), 1);

    t -= (B(t).x - x) / d;
}

return t;
Lazlo
  • 189
  • 7
  • This is actually a duplicate of http://math.stackexchange.com/questions/26846/is-there-an-explicit-form-for-cubic-bezier-curves, but Henning Makholm's solution to use Newton's method is better (imho) than the binary tree search accepted in the other question. – Lazlo Dec 27 '11 at 03:36
  • Actually, it is a false statement that my curve is ascending. Fortunately, the solution (Newton-Raphson) works in any case. – Lazlo Dec 27 '11 at 06:26

2 Answers2

2

There are closed formulas for solving general cubic equations, but as Wolfram Alpha showed you, they are so horrible (involving complex cube roots even if your coefficients are real) that you might as well not bother, and instead find your $t$ numerically -- such as with bisection or Newton-Raphson.

  • Considering I'm writing code, I don't really know when $t$ has "enough precision". I managed to rearrange the equation by coefficients, and it seems there could be some common factors, but I'm bad at trial and error. – Lazlo Dec 27 '11 at 03:10
  • @Lazlo: You need to decide for yourself how much precision you want -- or take it as a parameter to your function and let the caller worry about it :-) – hmakholm left over Monica Dec 27 '11 at 03:13
  • Sadly, I'm the caller. I guess I'm just a bit stubborn about it. There should be some way of finding $t(x)$ in a Bézier, no? – Lazlo Dec 27 '11 at 03:16
  • 1
    Rearranging the equation is not going to help you in the general case, since every third degree polynomial can be gotten out of the Bézier construction by choosing appropriate control points. Because the inequalities for the control points you have do not reduce the degrees of freedom of the problem, I wouldn't expect them to be enough to allow a tractable closed-form solution. – hmakholm left over Monica Dec 27 '11 at 03:17
  • 1
    @Lazlo: "There should be some way" -- there is, and that way is numerical approximation. If you're the caller, you must know what the solution is going to be used for, and thence how much of an error is acceptable for that purpose. (Note that even something simple such as the quadratic formula comes with rounding errors when implemented on a machine; you just don't get to choose your precision in that case). – hmakholm left over Monica Dec 27 '11 at 03:19
  • In the so-called casus irreducibilis case (cubic with all roots real), there are trigonometric formulae for the roots of a cubic formula. I still consider them to be rather unwieldy. I'd instead use one of those specialized numerical methods for Bézier curves myself. – J. M. ain't a mathematician Dec 27 '11 at 03:29
  • @J.M.: There are specialized numerical methods for Bézier curves? Sounds like it might be worth an answer of its own. – hmakholm left over Monica Dec 27 '11 at 03:50
  • I linked to one such method in the comments to the question ("Bézier clipping", they call it, and @Lazlo can now use that phrase for his future searches); I'd treat it as an embellishment of Newton-Raphson that exploits the beautiful properties of the Bernstein polynomials. I'm not really sure what else to say that should be an answer instead of a comment, so I'm content with this. :) – J. M. ain't a mathematician Dec 27 '11 at 03:54
0

Lazlo Bonin. Well, it's an equation, solved only with parameter t. In this equation P0, P1, P2, P3 are the points with (X,Y) (for 2D) parametrs. So you need functions that multiply Point and double (f.e. (t^3)*P3 ) Also you need function to sum points. So, after you make all calculations you get point. t is a parameter, in space [0,1]. So you need to make an accuracy or step of increasing parametr t. For example, step 0.01, you get 100 points between your points. Here is how a function of calaculating cube splain can be made:

private List<Point> cube(Point p0, Point p1, Point p2, Point p3)
    {
        List<Point> result = new List<Point>(accuracy + 1);
        double t = 0;
        while (t <= 1)
        {
            Point p = p3.Multiply(t * t * t);
            p = p.Add(p2.Multiply((1 - t) * 3 * t * t));
            p = p.Add(p1.Multiply(3 * t * (1 - t) * (1 - t)));
            p = p.Add(p0.Multiply(Math.Pow(1 - t, 3)));
            result.Add(p);
            t += 1.0 / accuracy;
        }
        return result;
    }
  • I know how to use $P=B(t)$. I want to solve $P_x=B(t)$ for $t$, thus finding $P_x$, or $x$. – Lazlo Dec 27 '11 at 03:00