1
1

Address review comments

This commit is contained in:
2022-12-17 01:26:59 +01:00
parent 7c56b50c3b
commit f3e93e9ace
3 changed files with 47 additions and 49 deletions

View File

@@ -774,17 +774,17 @@ template<typename T, int NumCol, int NumRow>
MatBase<T, NumCol, NumRow> from_rotation(const EulerXYZ<T> &rotation)
{
using MatT = MatBase<T, NumCol, NumRow>;
using IntermediateType = typename TypeTraits<T>::IntermediateType;
IntermediateType ci = math::cos(rotation.x);
IntermediateType cj = math::cos(rotation.y);
IntermediateType ch = math::cos(rotation.z);
IntermediateType si = math::sin(rotation.x);
IntermediateType sj = math::sin(rotation.y);
IntermediateType sh = math::sin(rotation.z);
IntermediateType cc = ci * ch;
IntermediateType cs = ci * sh;
IntermediateType sc = si * ch;
IntermediateType ss = si * sh;
using DoublePrecision = typename TypeTraits<T>::DoublePrecision;
DoublePrecision ci = math::cos(rotation.x);
DoublePrecision cj = math::cos(rotation.y);
DoublePrecision ch = math::cos(rotation.z);
DoublePrecision si = math::sin(rotation.x);
DoublePrecision sj = math::sin(rotation.y);
DoublePrecision sh = math::sin(rotation.z);
DoublePrecision cc = ci * ch;
DoublePrecision cs = ci * sh;
DoublePrecision sc = si * ch;
DoublePrecision ss = si * sh;
MatT mat;
mat[0][0] = T(cj * ch);
@@ -805,21 +805,21 @@ template<typename T, int NumCol, int NumRow>
MatBase<T, NumCol, NumRow> from_rotation(const Quaternion<T> &rotation)
{
using MatT = MatBase<T, NumCol, NumRow>;
using IntermediateType = typename TypeTraits<T>::IntermediateType;
IntermediateType q0 = M_SQRT2 * IntermediateType(rotation.x);
IntermediateType q1 = M_SQRT2 * IntermediateType(rotation.y);
IntermediateType q2 = M_SQRT2 * IntermediateType(rotation.z);
IntermediateType q3 = M_SQRT2 * IntermediateType(rotation.w);
using DoublePrecision = typename TypeTraits<T>::DoublePrecision;
DoublePrecision q0 = M_SQRT2 * DoublePrecision(rotation.x);
DoublePrecision q1 = M_SQRT2 * DoublePrecision(rotation.y);
DoublePrecision q2 = M_SQRT2 * DoublePrecision(rotation.z);
DoublePrecision q3 = M_SQRT2 * DoublePrecision(rotation.w);
IntermediateType qda = q0 * q1;
IntermediateType qdb = q0 * q2;
IntermediateType qdc = q0 * q3;
IntermediateType qaa = q1 * q1;
IntermediateType qab = q1 * q2;
IntermediateType qac = q1 * q3;
IntermediateType qbb = q2 * q2;
IntermediateType qbc = q2 * q3;
IntermediateType qcc = q3 * q3;
DoublePrecision qda = q0 * q1;
DoublePrecision qdb = q0 * q2;
DoublePrecision qdc = q0 * q3;
DoublePrecision qaa = q1 * q1;
DoublePrecision qab = q1 * q2;
DoublePrecision qac = q1 * q3;
DoublePrecision qbb = q2 * q2;
DoublePrecision qbc = q2 * q3;
DoublePrecision qcc = q3 * q3;
MatT mat;
mat[0][0] = T(1.0 - qbb - qcc);

View File

@@ -135,6 +135,7 @@ template<typename T> struct AxisAngle {
/**
* Create a rotation from an axis and an angle.
* \note `axis` does not have to be normalized.
* Use `AxisAngleNormalized` instead to skip normalization cost.
*/
AxisAngle(const vec3_type &axis, T angle);
@@ -217,13 +218,13 @@ template<typename T> struct AxisAngleNormalized : public AxisAngle<T> {
/**
* Intermediate Types.
*
* Some functions need to have higher precision than standard floats.
* Some functions need to have higher precision than standard floats for some operations.
*/
template<typename T> struct TypeTraits {
using IntermediateType = T;
using DoublePrecision = T;
};
template<> struct TypeTraits<float> {
using IntermediateType = double;
using DoublePrecision = double;
};
}; // namespace detail

View File

@@ -12,10 +12,6 @@
#include <Eigen/Dense>
#include <Eigen/Eigenvalues>
using Eigen::Map;
using Eigen::Matrix;
using Eigen::Stride;
/* -------------------------------------------------------------------- */
/** \name Matrix multiplication
* \{ */
@@ -74,7 +70,7 @@ namespace blender::math {
template<typename T, int Size> T determinant(const MatBase<T, Size, Size> &mat)
{
return Map<const Matrix<T, Size, Size>>(mat.base_ptr()).determinant();
return Eigen::Map<const Eigen::Matrix<T, Size, Size>>(mat.base_ptr()).determinant();
}
template float determinant(const float2x2 &mat);
@@ -86,7 +82,8 @@ template double determinant(const double4x4 &mat);
template<typename T> bool is_negative(const MatBase<T, 4, 4> &mat)
{
return Map<const Matrix<T, 3, 3>, 0, Stride<4, 1>>(mat.base_ptr()).determinant() < T(0);
return Eigen::Map<const Eigen::Matrix<T, 3, 3>, 0, Eigen::Stride<4, 1>>(mat.base_ptr())
.determinant() < T(0);
}
template bool is_negative(const float4x4 &mat);
@@ -139,8 +136,8 @@ template<typename T, int Size>
MatBase<T, Size, Size> invert(const MatBase<T, Size, Size> &mat, bool &r_success)
{
MatBase<T, Size, Size> result;
Map<const Matrix<T, Size, Size>> M(mat.base_ptr());
Map<Matrix<T, Size, Size>> R(result.base_ptr());
Eigen::Map<const Eigen::Matrix<T, Size, Size>> M(mat.base_ptr());
Eigen::Map<Eigen::Matrix<T, Size, Size>> R(result.base_ptr());
M.computeInverseWithCheck(R, r_success, 0.0f);
if (!r_success) {
R = R.Zero();
@@ -177,21 +174,21 @@ MatBase<T, Size, Size> pseudo_invert(const MatBase<T, Size, Size> &mat, T epsilo
{
using namespace Eigen;
using MatrixT = Matrix<T, Size, Size>;
using VectorT = Matrix<T, Size, 1>;
using MatrixT = Eigen::Matrix<T, Size, Size>;
using VectorT = Eigen::Matrix<T, Size, 1>;
/* Blender and Eigen matrices are both column-major.
* Since our matrix is squared, we can use thinU/V. */
/** WORKAROUND:
* (ComputeThinU | ComputeThinV) must be set as runtime parameters in Eigen < 3.4.0.
* But this requires the matrix type to be dynamic to avoid an assert.
*/
using MatrixDynamicT = Matrix<T, Dynamic, Dynamic>;
using MatrixDynamicT = Eigen::Matrix<T, Dynamic, Dynamic>;
JacobiSVD<MatrixDynamicT, NoQRPreconditioner> svd(
Map<const MatrixDynamicT>(mat.base_ptr(), Size, Size), ComputeThinU | ComputeThinV);
Eigen::Map<const MatrixDynamicT>(mat.base_ptr(), Size, Size), ComputeThinU | ComputeThinV);
(Map<MatrixT>(U.base_ptr())) = svd.matrixU();
(Map<VectorT>(S_val)) = svd.singularValues();
(Map<MatrixT>(V.base_ptr())) = svd.matrixV();
(Eigen::Map<MatrixT>(U.base_ptr())) = svd.matrixU();
(Eigen::Map<VectorT>(S_val)) = svd.singularValues();
(Eigen::Map<MatrixT>(V.base_ptr())) = svd.matrixV();
}
/* Invert or nullify component based on epsilon comparison. */
@@ -237,20 +234,20 @@ static void polar_decompose(const MatBase<T, 3, 3> &mat3,
{
using namespace Eigen;
using MatrixT = Matrix<T, 3, 3>;
using VectorT = Matrix<T, 3, 1>;
using MatrixT = Eigen::Matrix<T, 3, 3>;
using VectorT = Eigen::Matrix<T, 3, 1>;
/* Blender and Eigen matrices are both column-major.
* Since our matrix is squared, we can use thinU/V. */
/** WORKAROUND:
* (ComputeThinU | ComputeThinV) must be set as runtime parameters in Eigen < 3.4.0.
* But this requires the matrix type to be dynamic to avoid an assert.
*/
using MatrixDynamicT = Matrix<T, Dynamic, Dynamic>;
using MatrixDynamicT = Eigen::Matrix<T, Dynamic, Dynamic>;
JacobiSVD<MatrixDynamicT, NoQRPreconditioner> svd(
Map<const MatrixDynamicT>(mat3.base_ptr(), 3, 3), ComputeThinU | ComputeThinV);
Eigen::Map<const MatrixDynamicT>(mat3.base_ptr(), 3, 3), ComputeThinU | ComputeThinV);
(Map<MatrixT>(W.base_ptr())) = svd.matrixU();
(Map<VectorT>(S_val)) = svd.singularValues();
(Eigen::Map<MatrixT>(W.base_ptr())) = svd.matrixU();
(Eigen::Map<VectorT>(S_val)) = svd.singularValues();
(Map<MatrixT>(V.base_ptr())) = svd.matrixV();
}