diff --git a/doc/changelog.dox b/doc/changelog.dox index 512623a10..0b8f5e2b4 100644 --- a/doc/changelog.dox +++ b/doc/changelog.dox @@ -718,6 +718,14 @@ See also: on ES2 builds by mistake --- the @gl_extension{EXT,texture_sRGB_R8} and @gl_extension{EXT,texture_sRGB_RG8} extensions require OpenGL ES 3.0 at least +- @ref Math::Complex::fromMatrix() and @ref Math::Quaternion::fromMatrix() + now additionaly assert that the input matrix is a pure rotation without any + reflections. Before it only asserted for orthogonality, but that led to + arbitrary or even invalid quaternions when a reflection matrix was passed, + in case of complex numbers the reflection information was just lost in the + process. Existing code that calls these with unsanitized inputs now + additionally needs to account for reflection as suggested in the + documentation. - @ref SceneGraph::Object::addChild() no longer requires the type constructor to have the last parameter a parent object object pointer, as that was quite limiting. Instead it's calling @ref SceneGraph::Object::setParent() diff --git a/doc/snippets/MagnumMath.cpp b/doc/snippets/MagnumMath.cpp index aa2dbd499..0377acd1d 100644 --- a/doc/snippets/MagnumMath.cpp +++ b/doc/snippets/MagnumMath.cpp @@ -41,6 +41,8 @@ #include "Magnum/Math/StrictWeakOrdering.h" #include "Magnum/Math/Swizzle.h" +#define DOXYGEN_IGNORE(...) __VA_ARGS__ + using namespace Magnum; using namespace Magnum::Math::Literals; @@ -1033,6 +1035,17 @@ Matrix3 transformation = static_cast(transformation); } +{ +/* [Matrix3-rotation-extract-reflection] */ +Matrix3 transformation = DOXYGEN_IGNORE({}); +Matrix2x2 rotation = transformation.rotation(); +Vector2 scaling = transformation.scaling(); +if(rotation.determinant() < 0.0f) { + rotation[0] *= -1.0f; + scaling[0] *= -1.0f; +} +/* [Matrix3-rotation-extract-reflection] */ +} { /* [Matrix4-usage] */ using namespace Math::Literals; @@ -1045,6 +1058,18 @@ Matrix4 transformation = static_cast(transformation); } +{ +/* [Matrix4-rotation-extract-reflection] */ +Matrix4 transformation = DOXYGEN_IGNORE({}); +Matrix3x3 rotation = transformation.rotation(); +Vector3 scaling = transformation.scaling(); +if(rotation.determinant() < 0.0f) { + rotation[0] *= -1.0f; + scaling[0] *= -1.0f; +} +/* [Matrix4-rotation-extract-reflection] */ +} + { /* [Quaternion-fromEuler] */ Rad x, y, z; diff --git a/src/Magnum/Math/Complex.h b/src/Magnum/Math/Complex.h index fed3febeb..d5b0b7039 100644 --- a/src/Magnum/Math/Complex.h +++ b/src/Magnum/Math/Complex.h @@ -115,15 +115,14 @@ template class Complex { /** * @brief Create complex number from rotation matrix * - * Expects that the matrix is orthogonal (i.e. pure rotation). + * Expects that the matrix is a pure rotation, i.e. orthogonal and + * without any reflection. See @ref Matrix3::rotation() const for an + * example of how to extract rotation, reflection and scaling + * components from a 2D transformation matrix. * @see @ref toMatrix(), @ref DualComplex::fromMatrix(), - * @ref Matrix::isOrthogonal() + * @ref Matrix::isOrthogonal(), @ref Matrix::determinant() */ - static Complex fromMatrix(const Matrix2x2& matrix) { - CORRADE_ASSERT(matrix.isOrthogonal(), - "Math::Complex::fromMatrix(): the matrix is not orthogonal:" << Corrade::Utility::Debug::newline << matrix, {}); - return Implementation::complexFromMatrix(matrix); - } + static Complex fromMatrix(const Matrix2x2& matrix); /** * @brief Default constructor @@ -620,6 +619,21 @@ template inline Complex slerp(const Complex& normalizedA, const C return (std::sin((T(1) - t)*a)*normalizedA + std::sin(t*a)*normalizedB)/std::sin(a); } +template inline Complex Complex::fromMatrix(const Matrix2x2& matrix) { + /* Checking for determinant equal to 1 ensures we have a pure rotation + without shear or reflections. + + Assuming a column of an identity matrix is allowed to have a length of + 1 ± ε, the determinant would then be (1 ± ε)^2. Which is + 1 ± 2ε + e^2, and given that higher powers of ε are unrepresentable, the + fuzzy comparison should be 1 ± 2ε. This is similar to + Vector::isNormalized(), which compares the dot product (length squared) + to 1 ± 2ε. */ + CORRADE_ASSERT(std::abs(matrix.determinant() - T(1)) < T(2)*TypeTraits::epsilon(), + "Math::Complex::fromMatrix(): the matrix is not a rotation:" << Corrade::Utility::Debug::newline << matrix, {}); + return Implementation::complexFromMatrix(matrix); +} + #ifndef CORRADE_NO_DEBUG /** @debugoperator{Complex} */ template Corrade::Utility::Debug& operator<<(Corrade::Utility::Debug& debug, const Complex& value) { diff --git a/src/Magnum/Math/Matrix.h b/src/Magnum/Math/Matrix.h index 3b328e32d..ffdc564f2 100644 --- a/src/Magnum/Math/Matrix.h +++ b/src/Magnum/Math/Matrix.h @@ -232,8 +232,9 @@ template class Matrix: public RectangularMatrixLaplace's formula: @f[ * \det \boldsymbol{A} = \sum_{j=1}^n (-1)^{i+j} a_{i,j} \det \boldsymbol{A}_{i,j} * @f] @f$ \boldsymbol{A}_{i,j} @f$ is @f$ \boldsymbol{A} @f$ without diff --git a/src/Magnum/Math/Matrix3.h b/src/Magnum/Math/Matrix3.h index 3afa7599b..e89c87a3a 100644 --- a/src/Magnum/Math/Matrix3.h +++ b/src/Magnum/Math/Matrix3.h @@ -324,9 +324,10 @@ template class Matrix3: public Matrix3x3 { * represent shear and reflection. Especially when non-uniform scaling * is involved, decomposition of the result into primary linear * transformations may have multiple equivalent solutions. See - * @ref Algorithms::svd() and @ref Algorithms::qr() for further info. - * See also @ref rotationShear(), @ref rotation() const and - * @ref scaling() const for extracting further properties. + * @ref rotation() const, @ref Algorithms::svd() and + * @ref Algorithms::qr() for further info. See also + * @ref rotationShear() and @ref scaling() const for extracting further + * properties. * * @see @ref from(const Matrix2x2&, const Vector2&), * @ref rotation(Rad), @ref Matrix4::rotationScaling() @@ -337,7 +338,7 @@ template class Matrix3: public Matrix3x3 { } /** - * @brief 2D rotation and shear part of the matrix + * @brief 2D rotation, reflection and shear part of the matrix * * Normalized upper-left 2x2 part of the matrix. Assuming the following * matrix, with the upper-left 2x2 part represented by column vectors @@ -362,7 +363,9 @@ template class Matrix3: public Matrix3x3 { * not require orthogonal input. See also @ref rotationScaling() and * @ref scaling() const for extracting other properties. The * @ref Algorithms::svd() and @ref Algorithms::qr() can be used to - * separate the rotation / shear properties. + * separate the rotation / shear components; see @ref rotation() const + * for an example of decomposing a rotation + reflection matrix into a + * pure rotation and signed scaling. * * @see @ref from(const Matrix2x2&, const Vector2&), * @ref rotation(Rad), @ref Matrix4::rotationShear() const @@ -373,7 +376,7 @@ template class Matrix3: public Matrix3x3 { } /** - * @brief 2D rotation part of the matrix + * @brief 2D rotation and reflection part of the matrix * * Normalized upper-left 2x2 part of the matrix. Expects that the * normalized part is orthogonal. Assuming the following matrix, with @@ -399,12 +402,23 @@ template class Matrix3: public Matrix3x3 { * added orthogonality requirement. See also @ref rotationScaling() and * @ref scaling() const for extracting other properties. * - * @note Extracting rotation part of a matrix this way may cause - * assertions in case you have unsanitized input (for example a - * model transformation loaded from an external source) or when - * you accumulate many transformations together (for example when - * controlling a FPS camera). To mitigate this, either first - * reorthogonalize the matrix using + * There's usually several solutions for decomposing the matrix into a + * rotation @f$ \boldsymbol{R} @f$ and a scaling @f$ \boldsymbol{S} @f$ + * that satisfy @f$ \boldsymbol{R} \boldsymbol{S} = \boldsymbol{M} @f$. + * One possibility that gives you always a pure rotation matrix without + * reflections (which can then be fed to @ref Complex::fromMatrix(), + * for example) is to flip an arbitrary column of the 2x2 part if its + * @ref determinant() is negative, and apply the sign flip to the + * corresponding scaling component instead: + * + * @snippet MagnumMath.cpp Matrix3-rotation-extract-reflection + * + * @note Extracting rotation part of a matrix with this function may + * cause assertions in case you have unsanitized input (for + * example a model transformation loaded from an external source) + * or when you accumulate many transformations together (for + * example when controlling a FPS camera). To mitigate this, + * either first reorthogonalize the matrix using * @ref Algorithms::gramSchmidtOrthogonalize(), decompose it to * basic linear transformations using @ref Algorithms::svd() or * @ref Algorithms::qr() or use a different transformation @@ -420,10 +434,10 @@ template class Matrix3: public Matrix3x3 { Matrix2x2 rotation() const; /** - * @brief 2D rotation part of the matrix assuming there is no scaling + * @brief 2D rotation and reflection part of the matrix assuming there is no scaling * - * Similar to @ref rotation(), but expects that the rotation part is - * orthogonal, saving the extra renormalization. Assuming the + * Similar to @ref rotation() const, but expects that the rotation part + * is orthogonal, saving the extra renormalization. Assuming the * following matrix, with the upper-left 2x2 part represented by column * vectors @f$ \boldsymbol{a} @f$ and @f$ \boldsymbol{b} @f$: @f[ * \begin{pmatrix} @@ -511,13 +525,10 @@ template class Matrix3: public Matrix3x3 { * @f] * * Note that the returned vector is sign-less and the signs are instead - * contained in @ref rotation() const / @ref rotationShear() const in - * order to ensure @f$ \boldsymbol{R} \boldsymbol{S} = \boldsymbol{M} @f$ - * for @f$ \boldsymbol{R} @f$ and @f$ \boldsymbol{S} @f$ extracted out - * of @f$ \boldsymbol{M} @f$. The signs can be extracted for example by - * applying @ref Math::sign() on a @ref diagonal(), but keep in mind - * that the signs can be negative even for pure rotation matrices. - * + * contained in @ref rotation() const / @ref rotationShear() const, + * meaning these contain rotation together with a potential reflection. + * See @ref rotation() const for an example of decomposing a rotation + + * reflection matrix into a pure rotation and signed scaling. * @see @ref scalingSquared(), @ref uniformScaling(), * @ref rotation() const, @ref Matrix4::scaling() const */ diff --git a/src/Magnum/Math/Matrix4.h b/src/Magnum/Math/Matrix4.h index 959f39ff1..348ad4151 100644 --- a/src/Magnum/Math/Matrix4.h +++ b/src/Magnum/Math/Matrix4.h @@ -575,9 +575,10 @@ template class Matrix4: public Matrix4x4 { * represent shear and reflection. Especially when non-uniform scaling * is involved, decomposition of the result into primary linear * transformations may have multiple equivalent solutions. See - * @ref Algorithms::svd() and @ref Algorithms::qr() for further info. - * See also @ref rotationShear(), @ref rotation() const and - * @ref scaling() const for extracting further properties. + * @ref rotation() const, @ref Algorithms::svd() and + * @ref Algorithms::qr() for further info. See also + * @ref rotationShear() and @ref scaling() const for extracting further + * properties. * * @see @ref from(const Matrix3x3&, const Vector3&), * @ref rotation(Rad, const Vector3&), @@ -590,7 +591,7 @@ template class Matrix4: public Matrix4x4 { } /** - * @brief 3D rotation and shear part of the matrix + * @brief 3D rotation, reflection and shear part of the matrix * * Normalized upper-left 3x3 part of the matrix. Assuming the following * matrix, with the upper-left 3x3 part represented by column vectors @@ -618,7 +619,9 @@ template class Matrix4: public Matrix4x4 { * not require orthogonal input. See also @ref rotationScaling() and * @ref scaling() const for extracting other properties. The * @ref Algorithms::svd() and @ref Algorithms::qr() can be used to - * separate the rotation / shear properties. + * separate the rotation / shear components; see @ref rotation() const + * for an example of decomposing a rotation + reflection matrix into a + * pure rotation and signed scaling. * * @see @ref from(const Matrix3x3&, const Vector3&), * @ref rotation(Rad, const Vector3&), @@ -631,7 +634,7 @@ template class Matrix4: public Matrix4x4 { } /** - * @brief 3D rotation part of the matrix + * @brief 3D rotation and reflection part of the matrix * * Normalized upper-left 3x3 part of the matrix. Expects that the * normalized part is orthogonal. Assuming the following matrix, with @@ -660,12 +663,23 @@ template class Matrix4: public Matrix4x4 { * added orthogonality requirement. See also @ref rotationScaling() and * @ref scaling() const for extracting other properties. * - * @note Extracting rotation part of a matrix this way may cause - * assertions in case you have unsanitized input (for example a - * model transformation loaded from an external source) or when - * you accumulate many transformations together (for example when - * controlling a FPS camera). To mitigate this, either first - * reorthogonalize the matrix using + * There's usually several solutions for decomposing the matrix into a + * rotation @f$ \boldsymbol{R} @f$ and a scaling @f$ \boldsymbol{S} @f$ + * that satisfy @f$ \boldsymbol{R} \boldsymbol{S} = \boldsymbol{M} @f$. + * One possibility that gives you always a pure rotation matrix without + * reflections (which can then be fed to @ref Quaternion::fromMatrix(), + * for example) is to flip an arbitrary column of the 3x3 part if its + * @ref determinant() is negative, and apply the sign flip to the + * corresponding scaling component instead: + * + * @snippet MagnumMath.cpp Matrix4-rotation-extract-reflection + * + * @note Extracting rotation part of a matrix with this function may + * cause assertions in case you have unsanitized input (for + * example a model transformation loaded from an external source) + * or when you accumulate many transformations together (for + * example when controlling a FPS camera). To mitigate this, + * either first reorthogonalize the matrix using * @ref Algorithms::gramSchmidtOrthogonalize(), decompose it to * basic linear transformations using @ref Algorithms::svd() or * @ref Algorithms::qr() or use a different transformation @@ -682,10 +696,10 @@ template class Matrix4: public Matrix4x4 { Matrix3x3 rotation() const; /** - * @brief 3D rotation part of the matrix assuming there is no scaling + * @brief 3D rotation and reflection part of the matrix assuming there is no scaling * - * Similar to @ref rotation(), but expects that the rotation part is - * orthogonal, saving the extra renormalization. Assuming the + * Similar to @ref rotation() const, but expects that the rotation part + * is orthogonal, saving the extra renormalization. Assuming the * following matrix, with the upper-left 3x3 part represented by column * vectors @f$ \boldsymbol{a} @f$, @f$ \boldsymbol{b} @f$ and * @f$ \boldsymbol{c} @f$: @f[ @@ -783,13 +797,10 @@ template class Matrix4: public Matrix4x4 { * @f] * * Note that the returned vector is sign-less and the signs are instead - * contained in @ref rotation() const / @ref rotationShear() const in - * order to ensure @f$ \boldsymbol{R} \boldsymbol{S} = \boldsymbol{M} @f$ - * for @f$ \boldsymbol{R} @f$ and @f$ \boldsymbol{S} @f$ extracted out - * of @f$ \boldsymbol{M} @f$. The signs can be extracted for example by - * applying @ref Math::sign() on a @ref diagonal(), but keep in mind - * that the signs can be negative even for pure rotation matrices. - * + * contained in @ref rotation() const / @ref rotationShear() const, + * meaning these contain rotation together with a potential reflection. + * See @ref rotation() const for an example of decomposing a rotation + + * reflection matrix into a pure rotation and signed scaling. * @see @ref scalingSquared(), @ref uniformScaling(), * @ref rotation() const, @ref Matrix3::scaling() const */ diff --git a/src/Magnum/Math/Quaternion.h b/src/Magnum/Math/Quaternion.h index f4b0708ef..de92fa39a 100644 --- a/src/Magnum/Math/Quaternion.h +++ b/src/Magnum/Math/Quaternion.h @@ -293,9 +293,12 @@ template class Quaternion { /** * @brief Create a quaternion from a rotation matrix * - * Expects that the matrix is orthogonal (i.e. pure rotation). + * Expects that the matrix is a pure rotation, i.e. orthogonal and + * without any reflection. See @ref Matrix4::rotation() const for an + * example of how to extract rotation, reflection and scaling + * components from a 3D transformation matrix. * @see @ref toMatrix(), @ref DualComplex::fromMatrix(), - * @ref Matrix::isOrthogonal() + * @ref Matrix::isOrthogonal(), @ref Matrix::determinant() */ static Quaternion fromMatrix(const Matrix3x3& matrix); @@ -750,8 +753,17 @@ template inline Quaternion Quaternion::rotation(const Rad angl } template inline Quaternion Quaternion::fromMatrix(const Matrix3x3& matrix) { - CORRADE_ASSERT(matrix.isOrthogonal(), - "Math::Quaternion::fromMatrix(): the matrix is not orthogonal:" << Corrade::Utility::Debug::newline << matrix, {}); + /* Checking for determinant equal to 1 ensures we have a pure rotation + without shear or reflections. + + Assuming a column of an identity matrix is allowed to have a length of + 1 ± ε, the determinant would then be (1 ± ε)^3. Which is + 1 ± 3ε + 3e^2 + ε^3, and given that higher powers of ε are + unrepresentable, the fuzzy comparison should be 1 ± 3ε. This is similar + to Vector::isNormalized(), which compares the dot product (length + squared) to 1 ± 2ε. */ + CORRADE_ASSERT(std::abs(matrix.determinant() - T(1)) < T(3)*TypeTraits::epsilon(), + "Math::Quaternion::fromMatrix(): the matrix is not a rotation:" << Corrade::Utility::Debug::newline << matrix, {}); return Implementation::quaternionFromMatrix(matrix); } diff --git a/src/Magnum/Math/Test/ComplexTest.cpp b/src/Magnum/Math/Test/ComplexTest.cpp index 1b1a91f25..18a472068 100644 --- a/src/Magnum/Math/Test/ComplexTest.cpp +++ b/src/Magnum/Math/Test/ComplexTest.cpp @@ -94,7 +94,7 @@ struct ComplexTest: Corrade::TestSuite::Tester { void angleNotNormalized(); void rotation(); void matrix(); - void matrixNotOrthogonal(); + void matrixNotRotation(); void lerp(); void lerpNotNormalized(); void slerp(); @@ -148,7 +148,7 @@ ComplexTest::ComplexTest() { &ComplexTest::angleNotNormalized, &ComplexTest::rotation, &ComplexTest::matrix, - &ComplexTest::matrixNotOrthogonal, + &ComplexTest::matrixNotRotation, &ComplexTest::lerp, &ComplexTest::lerpNotNormalized, &ComplexTest::slerp, @@ -499,21 +499,35 @@ void ComplexTest::matrix() { CORRADE_COMPARE(a.toMatrix(), m); CORRADE_COMPARE(Complex::fromMatrix(m), a); + + /* One reflection is bad (asserts in the test below), but two are fine */ + CORRADE_COMPARE(Complex::fromMatrix(( + Matrix3::scaling({-1.0f, -1.0f})*Matrix3::rotation(37.0_degf) + ).rotationScaling()), Complex::rotation(180.0_degf + 37.0_degf)); } -void ComplexTest::matrixNotOrthogonal() { +void ComplexTest::matrixNotRotation() { #ifdef CORRADE_NO_ASSERT CORRADE_SKIP("CORRADE_NO_ASSERT defined, can't test assertions"); #endif std::ostringstream out; Error redirectError{&out}; - - Complex::fromMatrix(Matrix3::rotation(37.0_degf).rotationScaling()*2); + /* Shear, using rotation() instead of rotationScaling() as that isn't + supposed to "fix" the shear */ + Complex::fromMatrix((Matrix3::scaling({2.0f, 1.0f})* + Matrix3::rotation(37.0_degf)).rotation()); + /* Reflection, using rotation() instead of rotationScaling() as that isn't + supposed to "fix" the reflection either */ + Complex::fromMatrix((Matrix3::scaling({-1.0f, 1.0f})* + Matrix3::rotation(37.0_degf)).rotation()); CORRADE_COMPARE(out.str(), - "Math::Complex::fromMatrix(): the matrix is not orthogonal:\n" - "Matrix(1.59727, -1.20363,\n" - " 1.20363, 1.59727)\n"); + "Math::Complex::fromMatrix(): the matrix is not a rotation:\n" + "Matrix(0.935781, -0.833258,\n" + " 0.352581, 0.552885)\n" + "Math::Complex::fromMatrix(): the matrix is not a rotation:\n" + "Matrix(-0.798635, 0.601815,\n" + " 0.601815, 0.798635)\n"); } void ComplexTest::lerp() { diff --git a/src/Magnum/Math/Test/QuaternionTest.cpp b/src/Magnum/Math/Test/QuaternionTest.cpp index 21ccd4ae8..274cb9b3c 100644 --- a/src/Magnum/Math/Test/QuaternionTest.cpp +++ b/src/Magnum/Math/Test/QuaternionTest.cpp @@ -98,7 +98,7 @@ struct QuaternionTest: Corrade::TestSuite::Tester { void angleNormalizedButOver1(); void angleNotNormalized(); void matrix(); - void matrixNotOrthogonal(); + void matrixNotRotation(); void euler(); void eulerNotNormalized(); @@ -181,7 +181,7 @@ QuaternionTest::QuaternionTest() { &QuaternionTest::angleNormalizedButOver1, &QuaternionTest::angleNotNormalized, &QuaternionTest::matrix, - &QuaternionTest::matrixNotOrthogonal, + &QuaternionTest::matrixNotRotation, &QuaternionTest::euler, &QuaternionTest::eulerNotNormalized, @@ -601,24 +601,37 @@ void QuaternionTest::matrix() { Math::max(m4.diagonal()[1], m4.diagonal()[2]), Corrade::TestSuite::Compare::Greater); CORRADE_COMPARE(Quaternion::fromMatrix(m4), q4); + + /* One reflection is bad (asserts in the test below), but two are fine */ + CORRADE_COMPARE(Quaternion::fromMatrix(( + Matrix4::scaling({-1.0f, -1.0f, 1.0f})*Matrix4::rotationZ(37.0_degf) + ).rotation()), Quaternion::rotation(180.0_degf + 37.0_degf, Vector3::zAxis())); } -void QuaternionTest::matrixNotOrthogonal() { +void QuaternionTest::matrixNotRotation() { #ifdef CORRADE_NO_ASSERT CORRADE_SKIP("CORRADE_NO_ASSERT defined, can't test assertions"); #endif std::ostringstream out; Error redirectError{&out}; - - Vector3 axis = Vector3(-3.0f, 1.0f, 5.0f).normalized(); - Matrix3x3 m = Matrix4::rotation(37.0_degf, axis).rotationScaling(); - Quaternion::fromMatrix(m*2); + /* Shear, using rotation() instead of rotationScaling() as that isn't + supposed to "fix" the shear */ + Quaternion::fromMatrix((Matrix4::scaling({2.0f, 1.0f, 1.0f})* + Matrix4::rotationZ(37.0_degf)).rotation()); + /* Reflection, using rotation() instead of rotationScaling() as that isn't + supposed to "fix" the reflection either */ + Quaternion::fromMatrix((Matrix4::scaling({-1.0f, 1.0f, 1.0f})* + Matrix4::rotationZ(37.0_degf)).rotation()); CORRADE_COMPARE(out.str(), - "Math::Quaternion::fromMatrix(): the matrix is not orthogonal:\n" - "Matrix(1.70083, -1.05177, 0.0308525,\n" - " 0.982733, 1.60878, 0.667885,\n" - " -0.376049, -0.552819, 1.88493)\n"); + "Math::Quaternion::fromMatrix(): the matrix is not a rotation:\n" + "Matrix(0.935781, -0.833258, 0,\n" + " 0.352581, 0.552885, 0,\n" + " 0, 0, 1)\n" + "Math::Quaternion::fromMatrix(): the matrix is not a rotation:\n" + "Matrix(-0.798635, 0.601815, 0,\n" + " 0.601815, 0.798635, 0,\n" + " 0, 0, 1)\n"); } void QuaternionTest::euler() {