Let's say I have a Bezier curve B(u), if I increment u parameter at a constant rate I don't obtain a costant speed movement along the curve, because the relation between u parameter and the point obtained evaluating the curve is not linear.
I've read and implemented a David Eberly's article. It explains how to move at constant speed along a parametric curve.
Suppose I have a function F(t) that takes as input a time value t and a speed function sigma that returns the speed value at time t, I can obtain a costant speed movement along the curve, varying t parameter at a constant rate: B(F(t))
The core of the article I'm using is the following function:
float umin, umax; // The curve parameter interval [umin,umax].
Point Y (float u); // The position Y(u), umin <= u <= umax.
Point DY (float u); // The derivative dY(u)/du, umin <= u <= umax.
float LengthDY (float u) { return Length(DY(u)); }
float ArcLength (float t) { return Integral(umin,u,LengthDY()); }
float L = ArcLength(umax); // The total length of the curve.
float tmin, tmax; // The user-specified time interval [tmin,tmax]
float Sigma (float t); // The user-specified speed at time t.
float GetU (float t) // tmin <= t <= tmax
{
float h = (t - tmin)/n; // step size, `n' is application-specified
float u = umin; // initial condition
t = tmin; // initial condition
for (int i = 1; i <= n; i++)
{
// The divisions here might be a problem if the divisors are
// nearly zero.
float k1 = h*Sigma(t)/LengthDY(u);
float k2 = h*Sigma(t + h/2)/LengthDY(u + k1/2);
float k3 = h*Sigma(t + h/2)/LengthDY(u + k2/2);
float k4 = h*Sigma(t + h)/LengthDY(u + k3);
t += h;
u += (k1 + 2*(k2 + k3) + k4)/6;
}
return u;
}
It allows me to get the curve parameter u calculated using the supplied time t and sigma function.
Now the function works fine when the speed sigma is costant. If sigma is represent a uniform accelartion I'm getting wrong values from it.
Here's an example, of a straight Bezier curve, where P0 and P1 are the control points, T0 T1 the tangent. The curve's defined:
[x,y,z]= B(u) =(1–u)3P0 + 3(1–u)2uT0 + 3(1–u)u2T1 + u3P2
Let's say I want to know the position along the curve at time t = 3.
If I a constant velocity:
float sigma(float t)
{
return 1f;
}
and the following data:
V0 = 1;
V1 = 1;
t0 = 0;
L = 10;
I can analitically calculate the position:
px = v0 * t = 1 * 3 = 3
If I solve the same equation using my Bezier spline and the algorithm above with n =5 I obtain:
px = 3.002595;
Considering numerically approximation the value is quite precise (I did a lot of test on that. I omit details but Bezier my curves implementation is fine and the length of the curve itself is calculated quite precisely using Gaussian Quadrature).
Now If I try to define sigma as a uniform acceleration function, I get bad results. Consider the following data:
V0 = 1;
V1 = 2;
t0 = 0;
L = 10;
I can calculate the time a particle will reach the P1 using linear motion equations:
L = 0.5 * (V0 + V1) * t1 =>
t1 = 2 * L / (V1 + V0) = 2 * 10 / 3 = 6.6666666
Having t I can calculate acceleration:
a = (V1 - V0) / (t1 - t0) = (2 - 1) / 6.6666666 = 0.15
I have all data to define my sigma function:
float sigma (float t)
{
float speed = V0 + a * t;
}
If I analitically solve this I'd expect the following velocity of a particle after time t =3 :
Vx = V0 + a * t = 1 + 0.15 * 3 = 1.45
and the position will be:
px = 0.5 * (V0 + Vx) * t = 0.5 * (1 + 1.45) * 3 = 3.675
But if I calculate it with the alorithm above, the position results:
px = 4.358587
that is quite different from what I'm expecting.
Sorry for the long post, if anyone has enough patience to read it, I'd be glad.
Do you have any suggestion?What am I missing?Anyone can tell me what I'm doing wrong?
EDIT: I'm trying with 3D Bezier curve. Defined this way:
public Vector3 Bezier(float t)
{
float a = 1f - t;
float a_2 = a * a;
float a_3 = a_2 *a;
float t_2 = t * t;
Vector3 point = (P0 * a_3) + (3f * a_2 * t * T0) + (3f * a * t_2 * T1) + t_2 * t * P1 ;
return point;
}
and the derivative:
public Vector3 Derivative(float t)
{
float a = 1f - t;
float a_2 = a * a;
float t_2 = t * t;
float t6 = 6f*t;
Vector3 der = -3f * a_2 * P0 + 3f * a_2 * T0 - t6 * a * T0 - 3f* t_2 * T1 + t6 * a * T1 + 3f * t_2 * P1;
return der;
}