1// Copyright (C) 2017 Klaralvdalens Datakonsult AB (KDAB).
2// SPDX-License-Identifier: LicenseRef-Qt-Commercial OR LGPL-3.0-only OR GPL-2.0-only OR GPL-3.0-only
3
4#include "bezierevaluator_p.h"
5#include <private/keyframe_p.h>
6#include <QtCore/qglobal.h>
7#include <QtCore/qdebug.h>
8
9#include <cmath>
10
11QT_BEGIN_NAMESPACE
12
13namespace {
14
15inline double qCbrt(double x)
16{
17 // Android is just broken and doesn't define cbrt in std namespace
18#if defined(Q_OS_ANDROID)
19 if (x > 0.0)
20 return std::pow(x, 1.0 / 3.0);
21 else if (x < 0.0)
22 return -std::pow(-x, 1.0 / 3.0);
23 else
24 return 0.0;
25#else
26 return std::cbrt(x: x);
27#endif
28}
29
30} // anonymous
31
32namespace Qt3DAnimation {
33namespace Animation {
34
35/*!
36 \internal
37
38 Evaluates the value of the cubic bezier at time \a time.
39 This requires first finding the value of the bezier parameter, u,
40 corresponding to the requested time which should itself be
41 sandwiched by the provided times and keyframes.
42
43 Once u is found, substitute this back into the cubic Bezier
44 equation using the y components of the keyframe control points.
45 */
46float BezierEvaluator::valueForTime(float time) const
47{
48 const float u = parameterForTime(time);
49
50 // Calculate powers of u and (1-u) that we need
51 const float u2 = u * u;
52 const float u3 = u2 * u;
53 const float mu = 1.0f - u;
54 const float mu2 = mu * mu;
55 const float mu3 = mu2 * mu;
56
57 // The cubic Bezier control points
58 const float p0 = m_keyframe0.value;
59 const float p1 = m_keyframe0.rightControlPoint.y();
60 const float p2 = m_keyframe1.leftControlPoint.y();
61 const float p3 = m_keyframe1.value;
62
63 // Evaluate the cubic Bezier function
64 return p0 * mu3 + 3.0f * p1 * mu2 * u + 3.0f * p2 * mu * u2 + p3 * u3;
65}
66
67/*!
68 \internal
69
70 Calculates the value of the Bezier parameter, u, for the
71 requested time which is the x coordinate of the Keyframes.
72
73 Given 4 ordered control points p0, p1, p2, and p3, the cubic
74 Bezier equation is:
75
76 x(u) = (1-u)^3 p0 + 3 (1-u)^2 u p1 + 3 (1-u) u^2 p2 + u^3 p3
77
78 To find the value of u that corresponds with a given x
79 value (time in the case of keyframes), we can expand the
80 above equation, and then collect terms to arrive at:
81
82 0 = a u^3 + b u^2 + c u + d
83
84 where
85
86 a = p3 - p0 + 3 (p1 - p2)
87 b = 3 (p0 - 2 p1 + p2)
88 c = 3 (p1 - p0)
89 d = p0 - x(u)
90
91 We can then use findCubicRoots to locate the single root of
92 this cubic equation found in the range [0,1] used for this
93 section of the FCurve. This works because the FCurve ensures
94 that the function it represents via the Bezier control points
95 in the Keyframes is single valued. (as a function of time).
96 Time, therefore must be single valued on the interval and
97 therefore have a single root for any given time in the interval
98 covered by the Keyframes.
99 */
100float BezierEvaluator::parameterForTime(float time) const
101{
102 Q_ASSERT(time >= m_time0);
103 Q_ASSERT(time <= m_time1);
104
105 const float p0 = m_time0;
106 const float p1 = m_keyframe0.rightControlPoint.x();
107 const float p2 = m_keyframe1.leftControlPoint.x();
108 const float p3 = m_time1;
109
110 const float coeffs[4] = {
111 p0 - time, // d
112 3.0f * (p1 - p0), // c
113 3.0f * (p0 - 2.0f * p1 + p2), // b
114 p3 - p0 + 3.0f * (p1 - p2) // a
115 };
116
117 float roots[3];
118 const int numberOfRoots = findCubicRoots(coefficients: coeffs, roots);
119 for (int i = 0; i < numberOfRoots; ++i) {
120 if (roots[i] >= -0.01f && roots[i] <= 1.01f)
121 return qMin(a: qMax(a: roots[i], b: 0.0f), b: 1.0f);
122 }
123
124 qWarning() << "Failed to find root of cubic bezier at time" << time
125 << "with coeffs: a =" << coeffs[3] << "b =" << coeffs[2]
126 << "c =" << coeffs[1] << "d =" << coeffs[0];
127 return 0.0f;
128}
129
130bool almostZero(float value, float threshold=1e-3f)
131{
132 // 1e-3 might seem excessively fuzzy, but any smaller value will make the
133 // factors a, b, and c large enough to knock out the cubic solver.
134 return value > -threshold && value < threshold;
135}
136
137/*!
138 \internal
139
140 Finds the roots of the cubic equation ax^3 + bx^2 + cx + d = 0 for
141 real coefficients and returns the number of roots. The roots are
142 put into the \a roots array. The coefficients should be passed in
143 as coeffs[0] = d, coeffs[1] = c, coeffs[2] = b, coeffs[3] = a.
144 */
145int BezierEvaluator::findCubicRoots(const float coeffs[4], float roots[3])
146{
147 const float a = coeffs[3];
148 const float b = coeffs[2];
149 const float c = coeffs[1];
150 const float d = coeffs[0];
151
152 // Simple cases with linear, quadratic or invalid equations
153 if (almostZero(value: a)) {
154 if (almostZero(value: b)) {
155 if (almostZero(value: c))
156 return 0;
157
158 roots[0] = -d / c;
159 return 1;
160 }
161 const float discriminant = c * c - 4.f * b * d;
162 if (discriminant < 0.f)
163 return 0;
164
165 if (discriminant == 0.f) {
166 roots[0] = -c / (2.f * b);
167 return 1;
168 }
169
170 roots[0] = (-c + std::sqrt(x: discriminant)) / (2.f * b);
171 roots[1] = (-c - std::sqrt(x: discriminant)) / (2.f * b);
172 return 2;
173 }
174
175 // See https://en.wikipedia.org/wiki/Cubic_function#General_solution_to_the_cubic_equation_with_real_coefficients
176 // for a description. We depress the general cubic to a form that can more easily be solved. Solve it and then
177 // substitue the results back to get the roots of the original cubic.
178 int numberOfRoots = 0;
179 const double oneThird = 1.0 / 3.0;
180 const double piByThree = M_PI / 3.0;
181
182 // Put cubic into normal format: x^3 + Ax^2 + Bx + C = 0
183 const double A = double(b / a);
184 const double B = double(c / a);
185 const double C = double(d / a);
186
187 // Substitute x = y - A/3 to eliminate quadratic term (depressed form):
188 // x^3 + px + q = 0
189 const double Asq = A * A;
190 const double p = oneThird * (-oneThird * Asq + B);
191 const double q = 1.0 / 2.0 * (2.0 / 27.0 * A * Asq - oneThird * A * B + C);
192
193 // Use Cardano's formula
194 const double pCubed = p * p * p;
195 const double discriminant = q * q + pCubed;
196
197 if (almostZero(value: discriminant, threshold: 1e-6f)) {
198 if (qIsNull(d: q)) {
199 // One repeated triple root
200 roots[0] = 0.0;
201 numberOfRoots = 1;
202 } else {
203 // One single and one double root
204 double u = qCbrt(x: -q);
205 roots[0] = 2.0 * u;
206 roots[1] = -u;
207 numberOfRoots = 2;
208 }
209 } else if (discriminant < 0) {
210 // Three real solutions
211 double phi = oneThird * std::acos(x: -q / std::sqrt(x: -pCubed));
212 double t = 2.0 * std::sqrt(x: -p);
213
214 roots[0] = t * std::cos(x: phi);
215 roots[1] = -t * std::cos(x: phi + piByThree);
216 roots[2] = -t * std::cos(x: phi - piByThree);
217 numberOfRoots = 3;
218 } else {
219 // One real solution
220 double sqrtDisc = std::sqrt(x: discriminant);
221 double u = qCbrt(x: sqrtDisc - q);
222 double v = -qCbrt(x: sqrtDisc + q);
223
224 roots[0] = u + v;
225 numberOfRoots = 1;
226 }
227
228 // Substitute back in
229 const double sub = oneThird * A;
230 for (int i = 0; i < numberOfRoots; ++i) {
231 roots[i] -= sub;
232 // Take care of cases where we are close to zero or one
233 if (almostZero(value: roots[i], threshold: 1e-6f))
234 roots[i] = 0.f;
235 if (almostZero(value: roots[i] - 1.f, threshold: 1e-6f))
236 roots[i] = 1.f;
237 }
238
239 return numberOfRoots;
240}
241
242} // namespace Animation
243} // namespace Qt3DAnimation
244
245QT_END_NAMESPACE
246

source code of qt3d/src/animation/backend/bezierevaluator.cpp