Скажем, у меня есть Bezier curveB(u)
, если приращение u
параметра с постоянной скоростью не получить скорость перемещения Костанта вдоль кривой, так как соотношение между u
параметром и точкой полученная оценка кривой не является линейной.Безье кубических кривых: перемещение с равномерным ускорением
Я прочитал и внедрил article Дэвида Эберли. В нем объясняется, как двигаться с постоянной скоростью по параметрической кривой.
Предположим, что у меня есть функция F(t)
, который принимает в качестве входных данных значение времени t
и функция скорости sigma
, которая возвращает значение скорости в момент времени t
, можно получить скорость перемещения Костанта вдоль кривой, изменяя параметр т с постоянной скоростью : B(F(t))
ядро статьи я использую следующие функции:
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;
}
это позволяет мне получить параметр u
кривой вычисленной с помощью прилагаемого времени t
и сигма. Теперь функция работает нормально, когда сигма скорости является дорогостоящей. Если сигма представляет собой единое ускорение, я получаю от него неправильные значения.
Вот пример прямой кривой Безье, где P0 и P1 - контрольные точки, T0 T1 - касательная. это определяется кривая:
[x,y,z]= B(u) =(1–u)3P0 + 3(1–u)2uT0 + 3(1–u)u2T1 + u3P2
Допустим, я хочу знать позицию вдоль кривой во время t = 3
. Если я постоянная скорость:
float sigma(float t)
{
return 1f;
}
и следующие данные:
V0 = 1;
V1 = 1;
t0 = 0;
L = 10;
Я могу аналитически рассчитать положение:
px = v0 * t = 1 * 3 = 3
Если я решаю то же самое уравнение, используя мой Безье сплайн и выше приведенный выше алгоритм с n =5
. Я получаю:
px = 3.002595;
Учитывая численное приближение, значение является довольно точным (я много сделал для этого. Я опускаю детали, но Безье моя реализация кривых прекрасна, и сама длина кривой вычисляется достаточно точно, используя Gaussian Quadrature).
Теперь, если я попытаюсь определить сигму как функцию равномерного ускорения, я получаю плохие результаты. Рассмотрим следующие данные:
V0 = 1;
V1 = 2;
t0 = 0;
L = 10;
можно вычислить время частица достигнет Р1 с использованием линейных уравнений движения:
L = 0.5 * (V0 + V1) * t1 =>
t1 = 2 * L/(V1 + V0) = 2 * 10/3 = 6.6666666
Имея t
можно рассчитать ускорение:
a = (V1 - V0)/(t1 - t0) = (2 - 1)/6.6666666 = 0.15
У меня есть все данные для определения моей сигма-функции:
float sigma (float t)
{
float speed = V0 + a * t;
}
Если я решить эту проблему аналитически я ожидаю следующую скорость частицы после времени t =3
:
Vx = V0 + a * t = 1 + 0.15 * 3 = 1.45
и позиция будет:
px = 0.5 * (V0 + Vx) * t = 0.5 * (1 + 1.45) * 3 = 3.675
Но если рассчитать его с alorithm выше, результаты позиции:
px = 4.358587
, что совсем другое fr о, чего я ожидаю.
Извините за длинный пост, если у кого есть достаточно терпения, чтобы прочитать его, я был бы рад.
У вас есть предложения? Что мне не хватает? Кто-нибудь может сказать мне, что я делаю неправильно?
EDIT: Я пытаюсь с кривой 3D Безье. Определяется следующим образом:
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;
}
и производное:
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;
}
и что алгоритм дает вам для t = 6.6666 ...? Это значение 10, то есть L, или другое? – lmsteffan