Welcome to mirror list, hosted at ThFree Co, Russian Federation.

git.blender.org/blender.git - Unnamed repository; edit this file 'description' to name the repository.
summaryrefslogtreecommitdiff
path: root/source/blender/blenlib/intern/math_rotation.c
diff options
context:
space:
mode:
Diffstat (limited to 'source/blender/blenlib/intern/math_rotation.c')
-rw-r--r--source/blender/blenlib/intern/math_rotation.c199
1 files changed, 117 insertions, 82 deletions
diff --git a/source/blender/blenlib/intern/math_rotation.c b/source/blender/blenlib/intern/math_rotation.c
index f0bfc7c21e1..ff45bbee5c9 100644
--- a/source/blender/blenlib/intern/math_rotation.c
+++ b/source/blender/blenlib/intern/math_rotation.c
@@ -176,7 +176,7 @@ void quat_to_compatible_quat(float q[4], const float a[4], const float old[4])
}
}
-/* skip error check, currently only needed by mat3_to_quat_is_ok */
+/* Skip error check, currently only needed by #mat3_to_quat_legacy. */
static void quat_to_mat3_no_error(float m[3][3], const float q[4])
{
double q0, q1, q2, q3, qda, qdb, qdc, qaa, qab, qac, qbb, qbc, qcc;
@@ -269,9 +269,11 @@ void quat_to_mat4(float m[4][4], const float q[4])
m[3][3] = 1.0f;
}
-void mat3_normalized_to_quat(float q[4], const float mat[3][3])
+void mat3_normalized_to_quat_fast(float q[4], const float mat[3][3])
{
BLI_ASSERT_UNIT_M3(mat);
+ /* Caller must ensure matrices aren't negative for valid results, see: T24291, T94231. */
+ BLI_assert(!is_negative_m3(mat));
/* Check the trace of the matrix - bad precision if close to -1. */
const float trace = mat[0][0] + mat[1][1] + mat[2][2];
@@ -332,34 +334,54 @@ void mat3_normalized_to_quat(float q[4], const float mat[3][3])
normalize_qt(q);
}
-void mat3_to_quat(float q[4], const float m[3][3])
-{
- float unit_mat[3][3];
- /* work on a copy */
- /* this is needed AND a 'normalize_qt' in the end */
- normalize_m3_m3(unit_mat, m);
- mat3_normalized_to_quat(q, unit_mat);
+static void mat3_normalized_to_quat_with_checks(float q[4], float mat[3][3])
+{
+ const float det = determinant_m3_array(mat);
+ if (UNLIKELY(!isfinite(det))) {
+ unit_m3(mat);
+ }
+ else if (UNLIKELY(det < 0.0f)) {
+ negate_m3(mat);
+ }
+ mat3_normalized_to_quat_fast(q, mat);
}
-void mat4_normalized_to_quat(float q[4], const float m[4][4])
+void mat3_normalized_to_quat(float q[4], const float mat[3][3])
{
- float mat3[3][3];
+ float unit_mat_abs[3][3];
+ copy_m3_m3(unit_mat_abs, mat);
+ mat3_normalized_to_quat_with_checks(q, unit_mat_abs);
+}
- copy_m3_m4(mat3, m);
- mat3_normalized_to_quat(q, mat3);
+void mat3_to_quat(float q[4], const float mat[3][3])
+{
+ float unit_mat_abs[3][3];
+ normalize_m3_m3(unit_mat_abs, mat);
+ mat3_normalized_to_quat_with_checks(q, unit_mat_abs);
}
-void mat4_to_quat(float q[4], const float m[4][4])
+void mat4_normalized_to_quat(float q[4], const float mat[4][4])
{
- float mat3[3][3];
+ float unit_mat_abs[3][3];
+ copy_m3_m4(unit_mat_abs, mat);
+ mat3_normalized_to_quat_with_checks(q, unit_mat_abs);
+}
- copy_m3_m4(mat3, m);
- mat3_to_quat(q, mat3);
+void mat4_to_quat(float q[4], const float mat[4][4])
+{
+ float unit_mat_abs[3][3];
+ copy_m3_m4(unit_mat_abs, mat);
+ normalize_m3(unit_mat_abs);
+ mat3_normalized_to_quat_with_checks(q, unit_mat_abs);
}
-void mat3_to_quat_is_ok(float q[4], const float wmat[3][3])
+void mat3_to_quat_legacy(float q[4], const float wmat[3][3])
{
+ /* Legacy version of #mat3_to_quat which has slightly different behavior.
+ * Keep for particle-system & boids since replacing this will make subtle changes
+ * that impact hair in existing files. See: D15772. */
+
float mat[3][3], matr[3][3], matn[3][3], q1[4], q2[4], angle, si, co, nor[3];
/* work on a copy */
@@ -498,7 +520,10 @@ void rotation_between_quats_to_quat(float q[4], const float q1[4], const float q
mul_qt_qtqt(q, tquat, q2);
}
-float quat_split_swing_and_twist(const float q_in[4], int axis, float r_swing[4], float r_twist[4])
+float quat_split_swing_and_twist(const float q_in[4],
+ const int axis,
+ float r_swing[4],
+ float r_twist[4])
{
BLI_assert(axis >= 0 && axis <= 2);
@@ -915,53 +940,63 @@ float tri_to_quat(float q[4], const float a[3], const float b[3], const float c[
return len;
}
-void sin_cos_from_fraction(const int numerator, const int denominator, float *r_sin, float *r_cos)
+void sin_cos_from_fraction(int numerator, int denominator, float *r_sin, float *r_cos)
{
- BLI_assert((numerator <= denominator) && (denominator > 0));
- if ((denominator & 3) == 0) {
- const int denominator_4 = denominator / 4;
- if (numerator <= denominator_4) {
- /* Fall through. */
- }
- else {
- if (numerator <= denominator_4 * 2) {
- const float phi = (float)(2.0 * M_PI) *
- ((float)(numerator - denominator_4) / (float)denominator);
- *r_sin = cosf(phi);
- *r_cos = -sinf(phi);
- }
- else if (numerator <= denominator_4 * 3) {
- const float phi = (float)(2.0 * M_PI) *
- ((float)(numerator - (denominator_4 * 2)) / (float)denominator);
- *r_sin = -sinf(phi);
- *r_cos = -cosf(phi);
- }
- else {
- const float phi = (float)(2.0 * M_PI) *
- ((float)(numerator - (denominator_4 * 3)) / (float)denominator);
- *r_cos = sinf(phi);
- *r_sin = -cosf(phi);
- }
- return;
- }
- }
- else if ((denominator & 1) == 0) {
- const int denominator_2 = denominator / 2;
- if (numerator <= denominator_2) {
- /* Fall through. */
- }
- else {
- const float phi = (float)(2.0 * M_PI) *
- ((float)(numerator - denominator_2) / (float)denominator);
- *r_sin = -sinf(phi);
- *r_cos = -cosf(phi);
- return;
- }
+ /* By default, creating a circle from an integer: calling #sinf & #cosf on the fraction doesn't
+ * create symmetrical values (because floats can't represent Pi exactly).
+ * Resolve this when the rotation is calculated from a fraction by mapping the `numerator`
+ * to lower values so X/Y values for points around a circle are exactly symmetrical, see T87779.
+ *
+ * Multiply both the `numerator` and `denominator` by eight to ensure we can divide the circle
+ * into 8 octants. For each octant, we then use symmetry and negation to bring the `numerator`
+ * closer to the origin where precision is highest.
+ *
+ * Cases 2, 4, 5 and 7, use the trigonometric identity sin(-x) == -sin(x).
+ * Cases 1, 2, 5 and 6, swap the pointers `r_sin` and `r_cos`.
+ */
+ BLI_assert(0 <= numerator);
+ BLI_assert(numerator <= denominator);
+ BLI_assert(denominator > 0);
+
+ numerator *= 8; /* Multiply numerator the same as denominator. */
+ const int octant = numerator / denominator; /* Determine the octant. */
+ denominator *= 8; /* Ensure denominator is a multiple of eight. */
+ float cos_sign = 1.0f; /* Either 1.0f or -1.0f. */
+
+ switch (octant) {
+ case 0:
+ /* Primary octant, nothing to do. */
+ break;
+ case 1:
+ case 2:
+ numerator = (denominator / 4) - numerator;
+ SWAP(float *, r_sin, r_cos);
+ break;
+ case 3:
+ case 4:
+ numerator = (denominator / 2) - numerator;
+ cos_sign = -1.0f;
+ break;
+ case 5:
+ case 6:
+ numerator = numerator - (denominator * 3 / 4);
+ SWAP(float *, r_sin, r_cos);
+ cos_sign = -1.0f;
+ break;
+ case 7:
+ numerator = numerator - denominator;
+ break;
+ default:
+ BLI_assert_unreachable();
}
- const float phi = (float)(2.0 * M_PI) * ((float)numerator / (float)denominator);
- *r_sin = sinf(phi);
- *r_cos = cosf(phi);
+ BLI_assert(-denominator / 4 <= numerator); /* Numerator may be negative. */
+ BLI_assert(numerator <= denominator / 4);
+ BLI_assert(cos_sign == -1.0f || cos_sign == 1.0f);
+
+ const float angle = (float)(2.0 * M_PI) * ((float)numerator / (float)denominator);
+ *r_sin = sinf(angle);
+ *r_cos = cosf(angle) * cos_sign;
}
void print_qt(const char *str, const float q[4])
@@ -1371,10 +1406,10 @@ void mat4_normalized_to_eul(float eul[3], const float m[4][4])
copy_m3_m4(mat3, m);
mat3_normalized_to_eul(eul, mat3);
}
-void mat4_to_eul(float eul[3], const float m[4][4])
+void mat4_to_eul(float eul[3], const float mat[4][4])
{
float mat3[3][3];
- copy_m3_m4(mat3, m);
+ copy_m3_m4(mat3, mat);
mat3_to_eul(eul, mat3);
}
@@ -1409,7 +1444,7 @@ void eul_to_quat(float quat[4], const float eul[3])
quat[3] = cj * cs - sj * sc;
}
-void rotate_eul(float beul[3], const char axis, const float ang)
+void rotate_eul(float beul[3], const char axis, const float angle)
{
float eul[3], mat1[3][3], mat2[3][3], totmat[3][3];
@@ -1417,13 +1452,13 @@ void rotate_eul(float beul[3], const char axis, const float ang)
eul[0] = eul[1] = eul[2] = 0.0f;
if (axis == 'X') {
- eul[0] = ang;
+ eul[0] = angle;
}
else if (axis == 'Y') {
- eul[1] = ang;
+ eul[1] = angle;
}
else {
- eul[2] = ang;
+ eul[2] = angle;
}
eul_to_mat3(mat1, eul);
@@ -1442,7 +1477,7 @@ void compatible_eul(float eul[3], const float oldrot[3])
const float pi_x2 = (2.0f * (float)M_PI);
float deul[3];
- unsigned int i;
+ uint i;
/* correct differences of about 360 degrees first */
for (i = 0; i < 3; i++) {
@@ -1779,23 +1814,23 @@ void mat3_to_compatible_eulO(float eul[3],
void mat4_normalized_to_compatible_eulO(float eul[3],
const float oldrot[3],
const short order,
- const float m[4][4])
+ const float mat[4][4])
{
float mat3[3][3];
/* for now, we'll just do this the slow way (i.e. copying matrices) */
- copy_m3_m4(mat3, m);
+ copy_m3_m4(mat3, mat);
mat3_normalized_to_compatible_eulO(eul, oldrot, order, mat3);
}
void mat4_to_compatible_eulO(float eul[3],
const float oldrot[3],
const short order,
- const float m[4][4])
+ const float mat[4][4])
{
float mat3[3][3];
/* for now, we'll just do this the slow way (i.e. copying matrices) */
- copy_m3_m4(mat3, m);
+ copy_m3_m4(mat3, mat);
normalize_m3(mat3);
mat3_normalized_to_compatible_eulO(eul, oldrot, order, mat3);
}
@@ -1814,7 +1849,7 @@ void quat_to_compatible_eulO(float eul[3],
/* rotate the given euler by the given angle on the specified axis */
/* NOTE: is this safe to do with different axis orders? */
-void rotate_eulO(float beul[3], const short order, char axis, float ang)
+void rotate_eulO(float beul[3], const short order, const char axis, const float angle)
{
float eul[3], mat1[3][3], mat2[3][3], totmat[3][3];
@@ -1823,13 +1858,13 @@ void rotate_eulO(float beul[3], const short order, char axis, float ang)
zero_v3(eul);
if (axis == 'X') {
- eul[0] = ang;
+ eul[0] = angle;
}
else if (axis == 'Y') {
- eul[1] = ang;
+ eul[1] = angle;
}
else {
- eul[2] = ang;
+ eul[2] = angle;
}
eulO_to_mat3(mat1, eul, order);
@@ -2143,8 +2178,8 @@ void quat_apply_track(float quat[4], short axis, short upflag)
* up axis is used X->Y, Y->X, Z->X, if this first up axis isn't used then rotate 90d
* the strange bit shift below just find the low axis {X:Y, Y:X, Z:X} */
if (upflag != (2 - axis) >> 1) {
- float q[4] = {sqrt_1_2, 0.0, 0.0, 0.0}; /* assign 90d rotation axis */
- q[axis + 1] = ((axis == 1)) ? sqrt_1_2 : -sqrt_1_2; /* flip non Y axis */
+ float q[4] = {sqrt_1_2, 0.0, 0.0, 0.0}; /* assign 90d rotation axis */
+ q[axis + 1] = (axis == 1) ? sqrt_1_2 : -sqrt_1_2; /* flip non Y axis */
mul_qt_qtqt(quat, quat, q);
}
}
@@ -2325,8 +2360,8 @@ bool mat3_from_axis_conversion(
value = ((src_forward << (0 * 3)) | (src_up << (1 * 3)) | (dst_forward << (2 * 3)) |
(dst_up << (3 * 3)));
- for (uint i = 0; i < (ARRAY_SIZE(_axis_convert_matrix)); i++) {
- for (uint j = 0; j < (ARRAY_SIZE(*_axis_convert_lut)); j++) {
+ for (uint i = 0; i < ARRAY_SIZE(_axis_convert_matrix); i++) {
+ for (uint j = 0; j < ARRAY_SIZE(*_axis_convert_lut); j++) {
if (_axis_convert_lut[i][j] == value) {
copy_m3_m3(r_mat, _axis_convert_matrix[i]);
return true;
generated by cgit v1.2.3 (git 2.25.1) at 2024-06-21 21:56:22 +0300