| /* |
| * Copyright (C) 2005, 2006 Apple Computer, Inc. All rights reserved. |
| * Copyright (C) 2009 Torch Mobile, Inc. |
| * Copyright (C) 2013 Google Inc. All rights reserved. |
| * |
| * Redistribution and use in source and binary forms, with or without |
| * modification, are permitted provided that the following conditions |
| * are met: |
| * 1. Redistributions of source code must retain the above copyright |
| * notice, this list of conditions and the following disclaimer. |
| * 2. Redistributions in binary form must reproduce the above copyright |
| * notice, this list of conditions and the following disclaimer in the |
| * documentation and/or other materials provided with the distribution. |
| * |
| * THIS SOFTWARE IS PROVIDED BY APPLE COMPUTER, INC. ``AS IS'' AND ANY |
| * EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE |
| * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR |
| * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL APPLE COMPUTER, INC. OR |
| * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, |
| * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, |
| * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR |
| * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY |
| * OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT |
| * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE |
| * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. |
| */ |
| |
| #include "platform/transforms/TransformationMatrix.h" |
| #if HAVE(MIPS_MSA_INTRINSICS) |
| #include "platform/cpu/mips/CommonMacrosMSA.h" |
| #endif |
| |
| #include "platform/geometry/FloatBox.h" |
| #include "platform/geometry/FloatQuad.h" |
| #include "platform/geometry/FloatRect.h" |
| #include "platform/geometry/IntRect.h" |
| #include "platform/geometry/LayoutRect.h" |
| #include "platform/transforms/AffineTransform.h" |
| #include "platform/transforms/Rotation.h" |
| |
| #include "wtf/Assertions.h" |
| #include "wtf/MathExtras.h" |
| #include "wtf/text/WTFString.h" |
| |
| #include <cmath> |
| #include <cstdlib> |
| |
| #if CPU(X86_64) |
| #include <emmintrin.h> |
| #endif |
| |
| namespace blink { |
| |
| // |
| // Supporting Math Functions |
| // |
| // This is a set of function from various places (attributed inline) to do |
| // things like inversion and decomposition of a 4x4 matrix. They are used |
| // throughout the code |
| // |
| |
| // |
| // Adapted from Matrix Inversion by Richard Carling, Graphics Gems |
| // <http://tog.acm.org/GraphicsGems/index.html>. |
| |
| // EULA: The Graphics Gems code is copyright-protected. In other words, you |
| // cannot claim the text of the code as your own and resell it. Using the code |
| // is permitted in any program, product, or library, non-commercial or |
| // commercial. Giving credit is not required, though is a nice gesture. The code |
| // comes as-is, and if there are any flaws or problems with any Gems code, |
| // nobody involved with Gems - authors, editors, publishers, or webmasters - are |
| // to be held responsible. Basically, don't be a jerk, and remember that |
| // anything free comes with no guarantee. |
| |
| // A clarification about the storage of matrix elements |
| // |
| // This class uses a 2 dimensional array internally to store the elements of the |
| // matrix. The first index into the array refers to the column that the element |
| // lies in; the second index refers to the row. |
| // |
| // In other words, this is the layout of the matrix: |
| // |
| // | m_matrix[0][0] m_matrix[1][0] m_matrix[2][0] m_matrix[3][0] | |
| // | m_matrix[0][1] m_matrix[1][1] m_matrix[2][1] m_matrix[3][1] | |
| // | m_matrix[0][2] m_matrix[1][2] m_matrix[2][2] m_matrix[3][2] | |
| // | m_matrix[0][3] m_matrix[1][3] m_matrix[2][3] m_matrix[3][3] | |
| |
| typedef double Vector4[4]; |
| typedef double Vector3[3]; |
| |
| const double SmallNumber = 1.e-8; |
| |
| // inverse(original_matrix, inverse_matrix) |
| // |
| // calculate the inverse of a 4x4 matrix |
| // |
| // -1 |
| // A = ___1__ adjoint A |
| // det A |
| |
| // double = determinant2x2(double a, double b, double c, double d) |
| // |
| // calculate the determinant of a 2x2 matrix. |
| |
| static double determinant2x2(double a, double b, double c, double d) { |
| return a * d - b * c; |
| } |
| |
| // double = determinant3x3(a1, a2, a3, b1, b2, b3, c1, c2, c3) |
| // |
| // Calculate the determinant of a 3x3 matrix |
| // in the form |
| // |
| // | a1, b1, c1 | |
| // | a2, b2, c2 | |
| // | a3, b3, c3 | |
| |
| static double determinant3x3(double a1, |
| double a2, |
| double a3, |
| double b1, |
| double b2, |
| double b3, |
| double c1, |
| double c2, |
| double c3) { |
| return a1 * determinant2x2(b2, b3, c2, c3) - |
| b1 * determinant2x2(a2, a3, c2, c3) + |
| c1 * determinant2x2(a2, a3, b2, b3); |
| } |
| |
| // double = determinant4x4(matrix) |
| // |
| // calculate the determinant of a 4x4 matrix. |
| |
| static double determinant4x4(const TransformationMatrix::Matrix4& m) { |
| // Assign to individual variable names to aid selecting |
| // correct elements |
| |
| double a1 = m[0][0]; |
| double b1 = m[0][1]; |
| double c1 = m[0][2]; |
| double d1 = m[0][3]; |
| |
| double a2 = m[1][0]; |
| double b2 = m[1][1]; |
| double c2 = m[1][2]; |
| double d2 = m[1][3]; |
| |
| double a3 = m[2][0]; |
| double b3 = m[2][1]; |
| double c3 = m[2][2]; |
| double d3 = m[2][3]; |
| |
| double a4 = m[3][0]; |
| double b4 = m[3][1]; |
| double c4 = m[3][2]; |
| double d4 = m[3][3]; |
| |
| return a1 * determinant3x3(b2, b3, b4, c2, c3, c4, d2, d3, d4) - |
| b1 * determinant3x3(a2, a3, a4, c2, c3, c4, d2, d3, d4) + |
| c1 * determinant3x3(a2, a3, a4, b2, b3, b4, d2, d3, d4) - |
| d1 * determinant3x3(a2, a3, a4, b2, b3, b4, c2, c3, c4); |
| } |
| |
| #if !CPU(ARM64) |
| // adjoint( original_matrix, inverse_matrix ) |
| // |
| // calculate the adjoint of a 4x4 matrix |
| // |
| // Let a denote the minor determinant of matrix A obtained by |
| // ij |
| // |
| // deleting the ith row and jth column from A. |
| // |
| // i+j |
| // Let b = (-1) a |
| // ij ji |
| // |
| // The matrix B = (b ) is the adjoint of A |
| // ij |
| |
| static inline void adjoint(const TransformationMatrix::Matrix4& matrix, |
| TransformationMatrix::Matrix4& result) { |
| // Assign to individual variable names to aid |
| // selecting correct values |
| double a1 = matrix[0][0]; |
| double b1 = matrix[0][1]; |
| double c1 = matrix[0][2]; |
| double d1 = matrix[0][3]; |
| |
| double a2 = matrix[1][0]; |
| double b2 = matrix[1][1]; |
| double c2 = matrix[1][2]; |
| double d2 = matrix[1][3]; |
| |
| double a3 = matrix[2][0]; |
| double b3 = matrix[2][1]; |
| double c3 = matrix[2][2]; |
| double d3 = matrix[2][3]; |
| |
| double a4 = matrix[3][0]; |
| double b4 = matrix[3][1]; |
| double c4 = matrix[3][2]; |
| double d4 = matrix[3][3]; |
| |
| // Row column labeling reversed since we transpose rows & columns |
| result[0][0] = determinant3x3(b2, b3, b4, c2, c3, c4, d2, d3, d4); |
| result[1][0] = -determinant3x3(a2, a3, a4, c2, c3, c4, d2, d3, d4); |
| result[2][0] = determinant3x3(a2, a3, a4, b2, b3, b4, d2, d3, d4); |
| result[3][0] = -determinant3x3(a2, a3, a4, b2, b3, b4, c2, c3, c4); |
| |
| result[0][1] = -determinant3x3(b1, b3, b4, c1, c3, c4, d1, d3, d4); |
| result[1][1] = determinant3x3(a1, a3, a4, c1, c3, c4, d1, d3, d4); |
| result[2][1] = -determinant3x3(a1, a3, a4, b1, b3, b4, d1, d3, d4); |
| result[3][1] = determinant3x3(a1, a3, a4, b1, b3, b4, c1, c3, c4); |
| |
| result[0][2] = determinant3x3(b1, b2, b4, c1, c2, c4, d1, d2, d4); |
| result[1][2] = -determinant3x3(a1, a2, a4, c1, c2, c4, d1, d2, d4); |
| result[2][2] = determinant3x3(a1, a2, a4, b1, b2, b4, d1, d2, d4); |
| result[3][2] = -determinant3x3(a1, a2, a4, b1, b2, b4, c1, c2, c4); |
| |
| result[0][3] = -determinant3x3(b1, b2, b3, c1, c2, c3, d1, d2, d3); |
| result[1][3] = determinant3x3(a1, a2, a3, c1, c2, c3, d1, d2, d3); |
| result[2][3] = -determinant3x3(a1, a2, a3, b1, b2, b3, d1, d2, d3); |
| result[3][3] = determinant3x3(a1, a2, a3, b1, b2, b3, c1, c2, c3); |
| } |
| #endif |
| |
| // Returns false if the matrix is not invertible |
| static bool inverse(const TransformationMatrix::Matrix4& matrix, |
| TransformationMatrix::Matrix4& result) { |
| // Calculate the 4x4 determinant |
| // If the determinant is zero, |
| // then the inverse matrix is not unique. |
| double det = determinant4x4(matrix); |
| |
| if (fabs(det) < SmallNumber) |
| return false; |
| |
| #if CPU(ARM64) |
| double rdet = 1 / det; |
| const double* mat = &(matrix[0][0]); |
| double* pr = &(result[0][0]); |
| asm volatile( |
| // mat: v16 - v23 |
| // m11, m12, m13, m14 |
| // m21, m22, m23, m24 |
| // m31, m32, m33, m34 |
| // m41, m42, m43, m44 |
| "ld1 {v16.2d - v19.2d}, [%[mat]], 64 \n\t" |
| "ld1 {v20.2d - v23.2d}, [%[mat]] \n\t" |
| "ins v30.d[0], %[rdet] \n\t" |
| // Determinant: right mat2x2 |
| "trn1 v0.2d, v17.2d, v21.2d \n\t" |
| "trn2 v1.2d, v19.2d, v23.2d \n\t" |
| "trn2 v2.2d, v17.2d, v21.2d \n\t" |
| "trn1 v3.2d, v19.2d, v23.2d \n\t" |
| "trn2 v5.2d, v21.2d, v23.2d \n\t" |
| "trn1 v4.2d, v17.2d, v19.2d \n\t" |
| "trn2 v6.2d, v17.2d, v19.2d \n\t" |
| "trn1 v7.2d, v21.2d, v23.2d \n\t" |
| "trn2 v25.2d, v23.2d, v21.2d \n\t" |
| "trn1 v27.2d, v23.2d, v21.2d \n\t" |
| "fmul v0.2d, v0.2d, v1.2d \n\t" |
| "fmul v1.2d, v4.2d, v5.2d \n\t" |
| "fmls v0.2d, v2.2d, v3.2d \n\t" |
| "fmul v2.2d, v4.2d, v25.2d \n\t" |
| "fmls v1.2d, v6.2d, v7.2d \n\t" |
| "fmls v2.2d, v6.2d, v27.2d \n\t" |
| // Adjoint: |
| // v24: A11A12, v25: A13A14 |
| // v26: A21A22, v27: A23A24 |
| "fmul v3.2d, v18.2d, v0.d[1] \n\t" |
| "fmul v4.2d, v16.2d, v0.d[1] \n\t" |
| "fmul v5.2d, v16.2d, v1.d[1] \n\t" |
| "fmul v6.2d, v16.2d, v2.d[1] \n\t" |
| "fmls v3.2d, v20.2d, v1.d[1] \n\t" |
| "fmls v4.2d, v20.2d, v2.d[0] \n\t" |
| "fmls v5.2d, v18.2d, v2.d[0] \n\t" |
| "fmls v6.2d, v18.2d, v1.d[0] \n\t" |
| "fmla v3.2d, v22.2d, v2.d[1] \n\t" |
| "fmla v4.2d, v22.2d, v1.d[0] \n\t" |
| "fmla v5.2d, v22.2d, v0.d[0] \n\t" |
| "fmla v6.2d, v20.2d, v0.d[0] \n\t" |
| "fneg v3.2d, v3.2d \n\t" |
| "fneg v5.2d, v5.2d \n\t" |
| "trn1 v26.2d, v3.2d, v4.2d \n\t" |
| "trn1 v27.2d, v5.2d, v6.2d \n\t" |
| "trn2 v24.2d, v3.2d, v4.2d \n\t" |
| "trn2 v25.2d, v5.2d, v6.2d \n\t" |
| "fneg v24.2d, v24.2d \n\t" |
| "fneg v25.2d, v25.2d \n\t" |
| // Inverse |
| // v24: I11I12, v25: I13I14 |
| // v26: I21I22, v27: I23I24 |
| "fmul v24.2d, v24.2d, v30.d[0] \n\t" |
| "fmul v25.2d, v25.2d, v30.d[0] \n\t" |
| "fmul v26.2d, v26.2d, v30.d[0] \n\t" |
| "fmul v27.2d, v27.2d, v30.d[0] \n\t" |
| "st1 {v24.2d - v27.2d}, [%[pr]], 64 \n\t" |
| // Determinant: left mat2x2 |
| "trn1 v0.2d, v16.2d, v20.2d \n\t" |
| "trn2 v1.2d, v18.2d, v22.2d \n\t" |
| "trn2 v2.2d, v16.2d, v20.2d \n\t" |
| "trn1 v3.2d, v18.2d, v22.2d \n\t" |
| "trn2 v5.2d, v20.2d, v22.2d \n\t" |
| "trn1 v4.2d, v16.2d, v18.2d \n\t" |
| "trn2 v6.2d, v16.2d, v18.2d \n\t" |
| "trn1 v7.2d, v20.2d, v22.2d \n\t" |
| "trn2 v25.2d, v22.2d, v20.2d \n\t" |
| "trn1 v27.2d, v22.2d, v20.2d \n\t" |
| "fmul v0.2d, v0.2d, v1.2d \n\t" |
| "fmul v1.2d, v4.2d, v5.2d \n\t" |
| "fmls v0.2d, v2.2d, v3.2d \n\t" |
| "fmul v2.2d, v4.2d, v25.2d \n\t" |
| "fmls v1.2d, v6.2d, v7.2d \n\t" |
| "fmls v2.2d, v6.2d, v27.2d \n\t" |
| // Adjoint: |
| // v24: A31A32, v25: A33A34 |
| // v26: A41A42, v27: A43A44 |
| "fmul v3.2d, v19.2d, v0.d[1] \n\t" |
| "fmul v4.2d, v17.2d, v0.d[1] \n\t" |
| "fmul v5.2d, v17.2d, v1.d[1] \n\t" |
| "fmul v6.2d, v17.2d, v2.d[1] \n\t" |
| "fmls v3.2d, v21.2d, v1.d[1] \n\t" |
| "fmls v4.2d, v21.2d, v2.d[0] \n\t" |
| "fmls v5.2d, v19.2d, v2.d[0] \n\t" |
| "fmls v6.2d, v19.2d, v1.d[0] \n\t" |
| "fmla v3.2d, v23.2d, v2.d[1] \n\t" |
| "fmla v4.2d, v23.2d, v1.d[0] \n\t" |
| "fmla v5.2d, v23.2d, v0.d[0] \n\t" |
| "fmla v6.2d, v21.2d, v0.d[0] \n\t" |
| "fneg v3.2d, v3.2d \n\t" |
| "fneg v5.2d, v5.2d \n\t" |
| "trn1 v26.2d, v3.2d, v4.2d \n\t" |
| "trn1 v27.2d, v5.2d, v6.2d \n\t" |
| "trn2 v24.2d, v3.2d, v4.2d \n\t" |
| "trn2 v25.2d, v5.2d, v6.2d \n\t" |
| "fneg v24.2d, v24.2d \n\t" |
| "fneg v25.2d, v25.2d \n\t" |
| // Inverse |
| // v24: I31I32, v25: I33I34 |
| // v26: I41I42, v27: I43I44 |
| "fmul v24.2d, v24.2d, v30.d[0] \n\t" |
| "fmul v25.2d, v25.2d, v30.d[0] \n\t" |
| "fmul v26.2d, v26.2d, v30.d[0] \n\t" |
| "fmul v27.2d, v27.2d, v30.d[0] \n\t" |
| "st1 {v24.2d - v27.2d}, [%[pr]] \n\t" |
| : [mat] "+r"(mat), [pr] "+r"(pr) |
| : [rdet] "r"(rdet) |
| : "memory", "v0", "v1", "v2", "v3", "v4", "v5", "v6", "v7", "v16", "v17", |
| "v18", "v19", "v20", "v21", "v22", "v23", "24", "25", "v26", "v27", |
| "v28", "v29", "v30"); |
| #elif HAVE(MIPS_MSA_INTRINSICS) |
| const double rDet = 1 / det; |
| const double* mat = &(matrix[0][0]); |
| v2f64 mat0, mat1, mat2, mat3, mat4, mat5, mat6, mat7; |
| v2f64 rev2, rev3, rev4, rev5, rev6, rev7; |
| v2f64 tmp0, tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7; |
| v2f64 det0, det1, det2, tmp8, tmp9, tmp10, tmp11; |
| const v2f64 rdet = COPY_DOUBLE_TO_VECTOR(rDet); |
| // mat0 mat1 --> m00 m01 m02 m03 |
| // mat2 mat3 --> m10 m11 m12 m13 |
| // mat4 mat5 --> m20 m21 m22 m23 |
| // mat6 mat7 --> m30 m31 m32 m33 |
| LD_DP8(mat, 2, mat0, mat1, mat2, mat3, mat4, mat5, mat6, mat7); |
| |
| // Right half |
| rev3 = SLDI_D(mat3, mat3, 8); // m13 m12 |
| rev5 = SLDI_D(mat5, mat5, 8); // m23 m22 |
| rev7 = SLDI_D(mat7, mat7, 8); // m33 m32 |
| |
| // 2*2 Determinants |
| // for A00 & A01 |
| tmp0 = mat5 * rev7; |
| tmp1 = mat3 * rev7; |
| tmp2 = mat3 * rev5; |
| // for A10 & A11 |
| tmp3 = mat1 * rev7; |
| tmp4 = mat1 * rev5; |
| // for A20 & A21 |
| tmp5 = mat1 * rev3; |
| // for A30 & A31 |
| tmp6 = (v2f64)__msa_ilvr_d((v2i64)tmp1, (v2i64)tmp0); |
| tmp7 = (v2f64)__msa_ilvl_d((v2i64)tmp1, (v2i64)tmp0); |
| det0 = tmp6 - tmp7; |
| tmp6 = (v2f64)__msa_ilvr_d((v2i64)tmp3, (v2i64)tmp2); |
| tmp7 = (v2f64)__msa_ilvl_d((v2i64)tmp3, (v2i64)tmp2); |
| det1 = tmp6 - tmp7; |
| tmp6 = (v2f64)__msa_ilvr_d((v2i64)tmp5, (v2i64)tmp4); |
| tmp7 = (v2f64)__msa_ilvl_d((v2i64)tmp5, (v2i64)tmp4); |
| det2 = tmp6 - tmp7; |
| |
| // Co-factors |
| tmp0 = mat0 * (v2f64)__msa_splati_d((v2i64)det0, 0); |
| tmp1 = mat0 * (v2f64)__msa_splati_d((v2i64)det0, 1); |
| tmp2 = mat0 * (v2f64)__msa_splati_d((v2i64)det1, 0); |
| tmp3 = mat2 * (v2f64)__msa_splati_d((v2i64)det0, 0); |
| tmp4 = mat2 * (v2f64)__msa_splati_d((v2i64)det1, 1); |
| tmp5 = mat2 * (v2f64)__msa_splati_d((v2i64)det2, 0); |
| tmp6 = mat4 * (v2f64)__msa_splati_d((v2i64)det0, 1); |
| tmp7 = mat4 * (v2f64)__msa_splati_d((v2i64)det1, 1); |
| tmp8 = mat4 * (v2f64)__msa_splati_d((v2i64)det2, 1); |
| tmp9 = mat6 * (v2f64)__msa_splati_d((v2i64)det1, 0); |
| tmp10 = mat6 * (v2f64)__msa_splati_d((v2i64)det2, 0); |
| tmp11 = mat6 * (v2f64)__msa_splati_d((v2i64)det2, 1); |
| |
| tmp0 -= tmp7; |
| tmp1 -= tmp4; |
| tmp2 -= tmp5; |
| tmp3 -= tmp6; |
| tmp0 += tmp10; |
| tmp1 += tmp11; |
| tmp2 += tmp8; |
| tmp3 += tmp9; |
| |
| // Multiply with 1/det |
| tmp0 *= rdet; |
| tmp1 *= rdet; |
| tmp2 *= rdet; |
| tmp3 *= rdet; |
| |
| // Inverse: Upper half |
| result[0][0] = tmp3[1]; |
| result[0][1] = -tmp0[1]; |
| result[0][2] = tmp1[1]; |
| result[0][3] = -tmp2[1]; |
| result[1][0] = -tmp3[0]; |
| result[1][1] = tmp0[0]; |
| result[1][2] = -tmp1[0]; |
| result[1][3] = tmp2[0]; |
| // Left half |
| rev2 = SLDI_D(mat2, mat2, 8); // m13 m12 |
| rev4 = SLDI_D(mat4, mat4, 8); // m23 m22 |
| rev6 = SLDI_D(mat6, mat6, 8); // m33 m32 |
| |
| // 2*2 Determinants |
| // for A00 & A01 |
| tmp0 = mat4 * rev6; |
| tmp1 = mat2 * rev6; |
| tmp2 = mat2 * rev4; |
| // for A10 & A11 |
| tmp3 = mat0 * rev6; |
| tmp4 = mat0 * rev4; |
| // for A20 & A21 |
| tmp5 = mat0 * rev2; |
| // for A30 & A31 |
| tmp6 = (v2f64)__msa_ilvr_d((v2i64)tmp1, (v2i64)tmp0); |
| tmp7 = (v2f64)__msa_ilvl_d((v2i64)tmp1, (v2i64)tmp0); |
| det0 = tmp6 - tmp7; |
| tmp6 = (v2f64)__msa_ilvr_d((v2i64)tmp3, (v2i64)tmp2); |
| tmp7 = (v2f64)__msa_ilvl_d((v2i64)tmp3, (v2i64)tmp2); |
| det1 = tmp6 - tmp7; |
| tmp6 = (v2f64)__msa_ilvr_d((v2i64)tmp5, (v2i64)tmp4); |
| tmp7 = (v2f64)__msa_ilvl_d((v2i64)tmp5, (v2i64)tmp4); |
| det2 = tmp6 - tmp7; |
| |
| // Co-factors |
| tmp0 = mat3 * (v2f64)__msa_splati_d((v2i64)det0, 0); |
| tmp1 = mat1 * (v2f64)__msa_splati_d((v2i64)det0, 1); |
| tmp2 = mat1 * (v2f64)__msa_splati_d((v2i64)det0, 0); |
| tmp3 = mat1 * (v2f64)__msa_splati_d((v2i64)det1, 0); |
| tmp4 = mat3 * (v2f64)__msa_splati_d((v2i64)det1, 1); |
| tmp5 = mat3 * (v2f64)__msa_splati_d((v2i64)det2, 0); |
| tmp6 = mat5 * (v2f64)__msa_splati_d((v2i64)det0, 1); |
| tmp7 = mat5 * (v2f64)__msa_splati_d((v2i64)det1, 1); |
| tmp8 = mat5 * (v2f64)__msa_splati_d((v2i64)det2, 1); |
| tmp9 = mat7 * (v2f64)__msa_splati_d((v2i64)det1, 0); |
| tmp10 = mat7 * (v2f64)__msa_splati_d((v2i64)det2, 0); |
| tmp11 = mat7 * (v2f64)__msa_splati_d((v2i64)det2, 1); |
| tmp0 -= tmp6; |
| tmp1 -= tmp4; |
| tmp2 -= tmp7; |
| tmp3 -= tmp5; |
| tmp0 += tmp9; |
| tmp1 += tmp11; |
| tmp2 += tmp10; |
| tmp3 += tmp8; |
| |
| // Multiply with 1/det |
| tmp0 *= rdet; |
| tmp1 *= rdet; |
| tmp2 *= rdet; |
| tmp3 *= rdet; |
| |
| // Inverse: Lower half |
| result[2][0] = tmp0[1]; |
| result[2][1] = -tmp2[1]; |
| result[2][2] = tmp1[1]; |
| result[2][3] = -tmp3[1]; |
| result[3][0] = -tmp0[0]; |
| result[3][1] = tmp2[0]; |
| result[3][2] = -tmp1[0]; |
| result[3][3] = tmp3[0]; |
| #else |
| // Calculate the adjoint matrix |
| adjoint(matrix, result); |
| |
| // Scale the adjoint matrix to get the inverse |
| for (int i = 0; i < 4; i++) |
| for (int j = 0; j < 4; j++) |
| result[i][j] = result[i][j] / det; |
| #endif |
| return true; |
| } |
| |
| // End of code adapted from Matrix Inversion by Richard Carling |
| |
| // Perform a decomposition on the passed matrix, return false if unsuccessful |
| // From Graphics Gems: unmatrix.c |
| |
| // Transpose rotation portion of matrix a, return b |
| static void transposeMatrix4(const TransformationMatrix::Matrix4& a, |
| TransformationMatrix::Matrix4& b) { |
| for (int i = 0; i < 4; i++) |
| for (int j = 0; j < 4; j++) |
| b[i][j] = a[j][i]; |
| } |
| |
| // Multiply a homogeneous point by a matrix and return the transformed point |
| static void v4MulPointByMatrix(const Vector4 p, |
| const TransformationMatrix::Matrix4& m, |
| Vector4 result) { |
| result[0] = |
| (p[0] * m[0][0]) + (p[1] * m[1][0]) + (p[2] * m[2][0]) + (p[3] * m[3][0]); |
| result[1] = |
| (p[0] * m[0][1]) + (p[1] * m[1][1]) + (p[2] * m[2][1]) + (p[3] * m[3][1]); |
| result[2] = |
| (p[0] * m[0][2]) + (p[1] * m[1][2]) + (p[2] * m[2][2]) + (p[3] * m[3][2]); |
| result[3] = |
| (p[0] * m[0][3]) + (p[1] * m[1][3]) + (p[2] * m[2][3]) + (p[3] * m[3][3]); |
| } |
| |
| static double v3Length(Vector3 a) { |
| return std::sqrt((a[0] * a[0]) + (a[1] * a[1]) + (a[2] * a[2])); |
| } |
| |
| static void v3Scale(Vector3 v, double desiredLength) { |
| double len = v3Length(v); |
| if (len != 0) { |
| double l = desiredLength / len; |
| v[0] *= l; |
| v[1] *= l; |
| v[2] *= l; |
| } |
| } |
| |
| static double v3Dot(const Vector3 a, const Vector3 b) { |
| return (a[0] * b[0]) + (a[1] * b[1]) + (a[2] * b[2]); |
| } |
| |
| // Make a linear combination of two vectors and return the result. |
| // result = (a * ascl) + (b * bscl) |
| static void v3Combine(const Vector3 a, |
| const Vector3 b, |
| Vector3 result, |
| double ascl, |
| double bscl) { |
| result[0] = (ascl * a[0]) + (bscl * b[0]); |
| result[1] = (ascl * a[1]) + (bscl * b[1]); |
| result[2] = (ascl * a[2]) + (bscl * b[2]); |
| } |
| |
| // Return the cross product result = a cross b */ |
| static void v3Cross(const Vector3 a, const Vector3 b, Vector3 result) { |
| result[0] = (a[1] * b[2]) - (a[2] * b[1]); |
| result[1] = (a[2] * b[0]) - (a[0] * b[2]); |
| result[2] = (a[0] * b[1]) - (a[1] * b[0]); |
| } |
| |
| static bool decompose(const TransformationMatrix::Matrix4& mat, |
| TransformationMatrix::DecomposedType& result) { |
| TransformationMatrix::Matrix4 localMatrix; |
| memcpy(localMatrix, mat, sizeof(TransformationMatrix::Matrix4)); |
| |
| // Normalize the matrix. |
| if (localMatrix[3][3] == 0) |
| return false; |
| |
| int i, j; |
| for (i = 0; i < 4; i++) |
| for (j = 0; j < 4; j++) |
| localMatrix[i][j] /= localMatrix[3][3]; |
| |
| // perspectiveMatrix is used to solve for perspective, but it also provides |
| // an easy way to test for singularity of the upper 3x3 component. |
| TransformationMatrix::Matrix4 perspectiveMatrix; |
| memcpy(perspectiveMatrix, localMatrix, sizeof(TransformationMatrix::Matrix4)); |
| for (i = 0; i < 3; i++) |
| perspectiveMatrix[i][3] = 0; |
| perspectiveMatrix[3][3] = 1; |
| |
| if (determinant4x4(perspectiveMatrix) == 0) |
| return false; |
| |
| // First, isolate perspective. This is the messiest. |
| if (localMatrix[0][3] != 0 || localMatrix[1][3] != 0 || |
| localMatrix[2][3] != 0) { |
| // rightHandSide is the right hand side of the equation. |
| Vector4 rightHandSide; |
| rightHandSide[0] = localMatrix[0][3]; |
| rightHandSide[1] = localMatrix[1][3]; |
| rightHandSide[2] = localMatrix[2][3]; |
| rightHandSide[3] = localMatrix[3][3]; |
| |
| // Solve the equation by inverting perspectiveMatrix and multiplying |
| // rightHandSide by the inverse. (This is the easiest way, not |
| // necessarily the best.) |
| TransformationMatrix::Matrix4 inversePerspectiveMatrix, |
| transposedInversePerspectiveMatrix; |
| if (!inverse(perspectiveMatrix, inversePerspectiveMatrix)) |
| return false; |
| transposeMatrix4(inversePerspectiveMatrix, |
| transposedInversePerspectiveMatrix); |
| |
| Vector4 perspectivePoint; |
| v4MulPointByMatrix(rightHandSide, transposedInversePerspectiveMatrix, |
| perspectivePoint); |
| |
| result.perspectiveX = perspectivePoint[0]; |
| result.perspectiveY = perspectivePoint[1]; |
| result.perspectiveZ = perspectivePoint[2]; |
| result.perspectiveW = perspectivePoint[3]; |
| |
| // Clear the perspective partition |
| localMatrix[0][3] = localMatrix[1][3] = localMatrix[2][3] = 0; |
| localMatrix[3][3] = 1; |
| } else { |
| // No perspective. |
| result.perspectiveX = result.perspectiveY = result.perspectiveZ = 0; |
| result.perspectiveW = 1; |
| } |
| |
| // Next take care of translation (easy). |
| result.translateX = localMatrix[3][0]; |
| localMatrix[3][0] = 0; |
| result.translateY = localMatrix[3][1]; |
| localMatrix[3][1] = 0; |
| result.translateZ = localMatrix[3][2]; |
| localMatrix[3][2] = 0; |
| |
| // Vector4 type and functions need to be added to the common set. |
| Vector3 row[3], pdum3; |
| |
| // Now get scale and shear. |
| for (i = 0; i < 3; i++) { |
| row[i][0] = localMatrix[i][0]; |
| row[i][1] = localMatrix[i][1]; |
| row[i][2] = localMatrix[i][2]; |
| } |
| |
| // Compute X scale factor and normalize first row. |
| result.scaleX = v3Length(row[0]); |
| v3Scale(row[0], 1.0); |
| |
| // Compute XY shear factor and make 2nd row orthogonal to 1st. |
| result.skewXY = v3Dot(row[0], row[1]); |
| v3Combine(row[1], row[0], row[1], 1.0, -result.skewXY); |
| |
| // Now, compute Y scale and normalize 2nd row. |
| result.scaleY = v3Length(row[1]); |
| v3Scale(row[1], 1.0); |
| result.skewXY /= result.scaleY; |
| |
| // Compute XZ and YZ shears, orthogonalize 3rd row. |
| result.skewXZ = v3Dot(row[0], row[2]); |
| v3Combine(row[2], row[0], row[2], 1.0, -result.skewXZ); |
| result.skewYZ = v3Dot(row[1], row[2]); |
| v3Combine(row[2], row[1], row[2], 1.0, -result.skewYZ); |
| |
| // Next, get Z scale and normalize 3rd row. |
| result.scaleZ = v3Length(row[2]); |
| v3Scale(row[2], 1.0); |
| result.skewXZ /= result.scaleZ; |
| result.skewYZ /= result.scaleZ; |
| |
| // At this point, the matrix (in rows[]) is orthonormal. |
| // Check for a coordinate system flip. If the determinant |
| // is -1, then negate the matrix and the scaling factors. |
| v3Cross(row[1], row[2], pdum3); |
| if (v3Dot(row[0], pdum3) < 0) { |
| result.scaleX *= -1; |
| result.scaleY *= -1; |
| result.scaleZ *= -1; |
| |
| for (i = 0; i < 3; i++) { |
| row[i][0] *= -1; |
| row[i][1] *= -1; |
| row[i][2] *= -1; |
| } |
| } |
| |
| // Now, get the rotations out, as described in the gem. |
| |
| // FIXME - Add the ability to return either quaternions (which are |
| // easier to recompose with) or Euler angles (rx, ry, rz), which |
| // are easier for authors to deal with. The latter will only be useful |
| // when we fix https://bugs.webkit.org/show_bug.cgi?id=23799, so I |
| // will leave the Euler angle code here for now. |
| |
| // ret.rotateY = asin(-row[0][2]); |
| // if (cos(ret.rotateY) != 0) { |
| // ret.rotateX = atan2(row[1][2], row[2][2]); |
| // ret.rotateZ = atan2(row[0][1], row[0][0]); |
| // } else { |
| // ret.rotateX = atan2(-row[2][0], row[1][1]); |
| // ret.rotateZ = 0; |
| // } |
| |
| double s, t, x, y, z, w; |
| |
| t = row[0][0] + row[1][1] + row[2][2] + 1.0; |
| |
| if (t > 1e-4) { |
| s = 0.5 / std::sqrt(t); |
| w = 0.25 / s; |
| x = (row[2][1] - row[1][2]) * s; |
| y = (row[0][2] - row[2][0]) * s; |
| z = (row[1][0] - row[0][1]) * s; |
| } else if (row[0][0] > row[1][1] && row[0][0] > row[2][2]) { |
| s = std::sqrt(1.0 + row[0][0] - row[1][1] - row[2][2]) * 2.0; // S=4*qx |
| x = 0.25 * s; |
| y = (row[0][1] + row[1][0]) / s; |
| z = (row[0][2] + row[2][0]) / s; |
| w = (row[2][1] - row[1][2]) / s; |
| } else if (row[1][1] > row[2][2]) { |
| s = std::sqrt(1.0 + row[1][1] - row[0][0] - row[2][2]) * 2.0; // S=4*qy |
| x = (row[0][1] + row[1][0]) / s; |
| y = 0.25 * s; |
| z = (row[1][2] + row[2][1]) / s; |
| w = (row[0][2] - row[2][0]) / s; |
| } else { |
| s = std::sqrt(1.0 + row[2][2] - row[0][0] - row[1][1]) * 2.0; // S=4*qz |
| x = (row[0][2] + row[2][0]) / s; |
| y = (row[1][2] + row[2][1]) / s; |
| z = 0.25 * s; |
| w = (row[1][0] - row[0][1]) / s; |
| } |
| |
| result.quaternionX = x; |
| result.quaternionY = y; |
| result.quaternionZ = z; |
| result.quaternionW = w; |
| |
| return true; |
| } |
| |
| // Perform a spherical linear interpolation between the two |
| // passed quaternions with 0 <= t <= 1 |
| static void slerp(double qa[4], const double qb[4], double t) { |
| double ax, ay, az, aw; |
| double bx, by, bz, bw; |
| double cx, cy, cz, cw; |
| double product; |
| |
| ax = qa[0]; |
| ay = qa[1]; |
| az = qa[2]; |
| aw = qa[3]; |
| bx = qb[0]; |
| by = qb[1]; |
| bz = qb[2]; |
| bw = qb[3]; |
| |
| product = ax * bx + ay * by + az * bz + aw * bw; |
| |
| product = clampTo(product, -1.0, 1.0); |
| |
| const double epsilon = 1e-5; |
| if (std::abs(product - 1.0) < epsilon) { |
| // Result is qa, so just return |
| return; |
| } |
| |
| double denom = std::sqrt(1.0 - product * product); |
| double theta = std::acos(product); |
| double w = std::sin(t * theta) * (1.0 / denom); |
| |
| double scale1 = std::cos(t * theta) - product * w; |
| double scale2 = w; |
| |
| cx = ax * scale1 + bx * scale2; |
| cy = ay * scale1 + by * scale2; |
| cz = az * scale1 + bz * scale2; |
| cw = aw * scale1 + bw * scale2; |
| |
| qa[0] = cx; |
| qa[1] = cy; |
| qa[2] = cz; |
| qa[3] = cw; |
| } |
| |
| // End of Supporting Math Functions |
| |
| TransformationMatrix::TransformationMatrix(const AffineTransform& t) { |
| checkAlignment(); |
| setMatrix(t.a(), t.b(), t.c(), t.d(), t.e(), t.f()); |
| } |
| |
| TransformationMatrix& TransformationMatrix::scale(double s) { |
| return scaleNonUniform(s, s); |
| } |
| |
| FloatPoint TransformationMatrix::projectPoint(const FloatPoint& p, |
| bool* clamped) const { |
| // This is basically raytracing. We have a point in the destination |
| // plane with z=0, and we cast a ray parallel to the z-axis from that |
| // point to find the z-position at which it intersects the z=0 plane |
| // with the transform applied. Once we have that point we apply the |
| // inverse transform to find the corresponding point in the source |
| // space. |
| // |
| // Given a plane with normal Pn, and a ray starting at point R0 and |
| // with direction defined by the vector Rd, we can find the |
| // intersection point as a distance d from R0 in units of Rd by: |
| // |
| // d = -dot (Pn', R0) / dot (Pn', Rd) |
| if (clamped) |
| *clamped = false; |
| |
| if (m33() == 0) { |
| // In this case, the projection plane is parallel to the ray we are trying |
| // to trace, and there is no well-defined value for the projection. |
| return FloatPoint(); |
| } |
| |
| double x = p.x(); |
| double y = p.y(); |
| double z = -(m13() * x + m23() * y + m43()) / m33(); |
| |
| // FIXME: use multVecMatrix() |
| double outX = x * m11() + y * m21() + z * m31() + m41(); |
| double outY = x * m12() + y * m22() + z * m32() + m42(); |
| |
| double w = x * m14() + y * m24() + z * m34() + m44(); |
| if (w <= 0) { |
| // Using int max causes overflow when other code uses the projected point. |
| // To represent infinity yet reduce the risk of overflow, we use a large but |
| // not-too-large number here when clamping. |
| const int largeNumber = 100000000 / kFixedPointDenominator; |
| outX = copysign(largeNumber, outX); |
| outY = copysign(largeNumber, outY); |
| if (clamped) |
| *clamped = true; |
| } else if (w != 1) { |
| outX /= w; |
| outY /= w; |
| } |
| |
| return FloatPoint(static_cast<float>(outX), static_cast<float>(outY)); |
| } |
| |
| FloatQuad TransformationMatrix::projectQuad(const FloatQuad& q, |
| bool* clamped) const { |
| FloatQuad projectedQuad; |
| |
| bool clamped1 = false; |
| bool clamped2 = false; |
| bool clamped3 = false; |
| bool clamped4 = false; |
| |
| projectedQuad.setP1(projectPoint(q.p1(), &clamped1)); |
| projectedQuad.setP2(projectPoint(q.p2(), &clamped2)); |
| projectedQuad.setP3(projectPoint(q.p3(), &clamped3)); |
| projectedQuad.setP4(projectPoint(q.p4(), &clamped4)); |
| |
| if (clamped) |
| *clamped = clamped1 || clamped2 || clamped3 || clamped4; |
| |
| // If all points on the quad had w < 0, then the entire quad would not be |
| // visible to the projected surface. |
| bool everythingWasClipped = clamped1 && clamped2 && clamped3 && clamped4; |
| if (everythingWasClipped) |
| return FloatQuad(); |
| |
| return projectedQuad; |
| } |
| |
| static float clampEdgeValue(float f) { |
| ASSERT(!std::isnan(f)); |
| return clampTo(f, (-LayoutUnit::max() / 2).toFloat(), |
| (LayoutUnit::max() / 2).toFloat()); |
| } |
| |
| LayoutRect TransformationMatrix::clampedBoundsOfProjectedQuad( |
| const FloatQuad& q) const { |
| FloatRect mappedQuadBounds = projectQuad(q).boundingBox(); |
| |
| float left = clampEdgeValue(floorf(mappedQuadBounds.x())); |
| float top = clampEdgeValue(floorf(mappedQuadBounds.y())); |
| |
| float right; |
| if (std::isinf(mappedQuadBounds.x()) && std::isinf(mappedQuadBounds.width())) |
| right = (LayoutUnit::max() / 2).toFloat(); |
| else |
| right = clampEdgeValue(ceilf(mappedQuadBounds.maxX())); |
| |
| float bottom; |
| if (std::isinf(mappedQuadBounds.y()) && std::isinf(mappedQuadBounds.height())) |
| bottom = (LayoutUnit::max() / 2).toFloat(); |
| else |
| bottom = clampEdgeValue(ceilf(mappedQuadBounds.maxY())); |
| |
| return LayoutRect(LayoutUnit::clamp(left), LayoutUnit::clamp(top), |
| LayoutUnit::clamp(right - left), |
| LayoutUnit::clamp(bottom - top)); |
| } |
| |
| void TransformationMatrix::transformBox(FloatBox& box) const { |
| FloatBox bounds; |
| bool firstPoint = true; |
| for (size_t i = 0; i < 2; ++i) { |
| for (size_t j = 0; j < 2; ++j) { |
| for (size_t k = 0; k < 2; ++k) { |
| FloatPoint3D point(box.x(), box.y(), box.z()); |
| point += |
| FloatPoint3D(i * box.width(), j * box.height(), k * box.depth()); |
| point = mapPoint(point); |
| if (firstPoint) { |
| bounds.setOrigin(point); |
| firstPoint = false; |
| } else { |
| bounds.expandTo(point); |
| } |
| } |
| } |
| } |
| box = bounds; |
| } |
| |
| FloatPoint TransformationMatrix::mapPoint(const FloatPoint& p) const { |
| if (isIdentityOrTranslation()) |
| return FloatPoint(p.x() + static_cast<float>(m_matrix[3][0]), |
| p.y() + static_cast<float>(m_matrix[3][1])); |
| |
| return internalMapPoint(p); |
| } |
| |
| FloatPoint3D TransformationMatrix::mapPoint(const FloatPoint3D& p) const { |
| if (isIdentityOrTranslation()) |
| return FloatPoint3D(p.x() + static_cast<float>(m_matrix[3][0]), |
| p.y() + static_cast<float>(m_matrix[3][1]), |
| p.z() + static_cast<float>(m_matrix[3][2])); |
| |
| return internalMapPoint(p); |
| } |
| |
| IntRect TransformationMatrix::mapRect(const IntRect& rect) const { |
| return enclosingIntRect(mapRect(FloatRect(rect))); |
| } |
| |
| LayoutRect TransformationMatrix::mapRect(const LayoutRect& r) const { |
| return enclosingLayoutRect(mapRect(FloatRect(r))); |
| } |
| |
| FloatRect TransformationMatrix::mapRect(const FloatRect& r) const { |
| if (isIdentityOrTranslation()) { |
| FloatRect mappedRect(r); |
| mappedRect.move(static_cast<float>(m_matrix[3][0]), |
| static_cast<float>(m_matrix[3][1])); |
| return mappedRect; |
| } |
| |
| FloatQuad result; |
| |
| float maxX = r.maxX(); |
| float maxY = r.maxY(); |
| result.setP1(internalMapPoint(FloatPoint(r.x(), r.y()))); |
| result.setP2(internalMapPoint(FloatPoint(maxX, r.y()))); |
| result.setP3(internalMapPoint(FloatPoint(maxX, maxY))); |
| result.setP4(internalMapPoint(FloatPoint(r.x(), maxY))); |
| |
| return result.boundingBox(); |
| } |
| |
| FloatQuad TransformationMatrix::mapQuad(const FloatQuad& q) const { |
| if (isIdentityOrTranslation()) { |
| FloatQuad mappedQuad(q); |
| mappedQuad.move(static_cast<float>(m_matrix[3][0]), |
| static_cast<float>(m_matrix[3][1])); |
| return mappedQuad; |
| } |
| |
| FloatQuad result; |
| result.setP1(internalMapPoint(q.p1())); |
| result.setP2(internalMapPoint(q.p2())); |
| result.setP3(internalMapPoint(q.p3())); |
| result.setP4(internalMapPoint(q.p4())); |
| return result; |
| } |
| |
| TransformationMatrix& TransformationMatrix::scaleNonUniform(double sx, |
| double sy) { |
| m_matrix[0][0] *= sx; |
| m_matrix[0][1] *= sx; |
| m_matrix[0][2] *= sx; |
| m_matrix[0][3] *= sx; |
| |
| m_matrix[1][0] *= sy; |
| m_matrix[1][1] *= sy; |
| m_matrix[1][2] *= sy; |
| m_matrix[1][3] *= sy; |
| return *this; |
| } |
| |
| TransformationMatrix& TransformationMatrix::scale3d(double sx, |
| double sy, |
| double sz) { |
| scaleNonUniform(sx, sy); |
| |
| m_matrix[2][0] *= sz; |
| m_matrix[2][1] *= sz; |
| m_matrix[2][2] *= sz; |
| m_matrix[2][3] *= sz; |
| return *this; |
| } |
| |
| TransformationMatrix& TransformationMatrix::rotate3d(const Rotation& rotation) { |
| return rotate3d(rotation.axis.x(), rotation.axis.y(), rotation.axis.z(), |
| rotation.angle); |
| } |
| |
| TransformationMatrix& TransformationMatrix::rotate3d(double x, |
| double y, |
| double z, |
| double angle) { |
| // Normalize the axis of rotation |
| double length = std::sqrt(x * x + y * y + z * z); |
| if (length == 0) { |
| // A direction vector that cannot be normalized, such as [0, 0, 0], will |
| // cause the rotation to not be applied. |
| return *this; |
| } else if (length != 1) { |
| x /= length; |
| y /= length; |
| z /= length; |
| } |
| |
| // Angles are in degrees. Switch to radians. |
| angle = deg2rad(angle); |
| |
| double sinTheta = std::sin(angle); |
| double cosTheta = std::cos(angle); |
| |
| TransformationMatrix mat; |
| |
| // Optimize cases where the axis is along a major axis |
| if (x == 1.0 && y == 0.0 && z == 0.0) { |
| mat.m_matrix[0][0] = 1.0; |
| mat.m_matrix[0][1] = 0.0; |
| mat.m_matrix[0][2] = 0.0; |
| mat.m_matrix[1][0] = 0.0; |
| mat.m_matrix[1][1] = cosTheta; |
| mat.m_matrix[1][2] = sinTheta; |
| mat.m_matrix[2][0] = 0.0; |
| mat.m_matrix[2][1] = -sinTheta; |
| mat.m_matrix[2][2] = cosTheta; |
| mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; |
| mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; |
| mat.m_matrix[3][3] = 1.0; |
| } else if (x == 0.0 && y == 1.0 && z == 0.0) { |
| mat.m_matrix[0][0] = cosTheta; |
| mat.m_matrix[0][1] = 0.0; |
| mat.m_matrix[0][2] = -sinTheta; |
| mat.m_matrix[1][0] = 0.0; |
| mat.m_matrix[1][1] = 1.0; |
| mat.m_matrix[1][2] = 0.0; |
| mat.m_matrix[2][0] = sinTheta; |
| mat.m_matrix[2][1] = 0.0; |
| mat.m_matrix[2][2] = cosTheta; |
| mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; |
| mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; |
| mat.m_matrix[3][3] = 1.0; |
| } else if (x == 0.0 && y == 0.0 && z == 1.0) { |
| mat.m_matrix[0][0] = cosTheta; |
| mat.m_matrix[0][1] = sinTheta; |
| mat.m_matrix[0][2] = 0.0; |
| mat.m_matrix[1][0] = -sinTheta; |
| mat.m_matrix[1][1] = cosTheta; |
| mat.m_matrix[1][2] = 0.0; |
| mat.m_matrix[2][0] = 0.0; |
| mat.m_matrix[2][1] = 0.0; |
| mat.m_matrix[2][2] = 1.0; |
| mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; |
| mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; |
| mat.m_matrix[3][3] = 1.0; |
| } else { |
| // This case is the rotation about an arbitrary unit vector. |
| // |
| // Formula is adapted from Wikipedia article on Rotation matrix, |
| // http://en.wikipedia.org/wiki/Rotation_matrix#Rotation_matrix_from_axis_and_angle |
| // |
| // An alternate resource with the same matrix: |
| // http://www.fastgraph.com/makegames/3drotation/ |
| // |
| double oneMinusCosTheta = 1 - cosTheta; |
| mat.m_matrix[0][0] = cosTheta + x * x * oneMinusCosTheta; |
| mat.m_matrix[0][1] = y * x * oneMinusCosTheta + z * sinTheta; |
| mat.m_matrix[0][2] = z * x * oneMinusCosTheta - y * sinTheta; |
| mat.m_matrix[1][0] = x * y * oneMinusCosTheta - z * sinTheta; |
| mat.m_matrix[1][1] = cosTheta + y * y * oneMinusCosTheta; |
| mat.m_matrix[1][2] = z * y * oneMinusCosTheta + x * sinTheta; |
| mat.m_matrix[2][0] = x * z * oneMinusCosTheta + y * sinTheta; |
| mat.m_matrix[2][1] = y * z * oneMinusCosTheta - x * sinTheta; |
| mat.m_matrix[2][2] = cosTheta + z * z * oneMinusCosTheta; |
| mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; |
| mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; |
| mat.m_matrix[3][3] = 1.0; |
| } |
| multiply(mat); |
| return *this; |
| } |
| |
| TransformationMatrix& TransformationMatrix::rotate3d(double rx, |
| double ry, |
| double rz) { |
| // Angles are in degrees. Switch to radians. |
| rx = deg2rad(rx); |
| ry = deg2rad(ry); |
| rz = deg2rad(rz); |
| |
| TransformationMatrix mat; |
| |
| double sinTheta = std::sin(rz); |
| double cosTheta = std::cos(rz); |
| |
| mat.m_matrix[0][0] = cosTheta; |
| mat.m_matrix[0][1] = sinTheta; |
| mat.m_matrix[0][2] = 0.0; |
| mat.m_matrix[1][0] = -sinTheta; |
| mat.m_matrix[1][1] = cosTheta; |
| mat.m_matrix[1][2] = 0.0; |
| mat.m_matrix[2][0] = 0.0; |
| mat.m_matrix[2][1] = 0.0; |
| mat.m_matrix[2][2] = 1.0; |
| mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; |
| mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; |
| mat.m_matrix[3][3] = 1.0; |
| |
| TransformationMatrix rmat(mat); |
| |
| sinTheta = std::sin(ry); |
| cosTheta = std::cos(ry); |
| |
| mat.m_matrix[0][0] = cosTheta; |
| mat.m_matrix[0][1] = 0.0; |
| mat.m_matrix[0][2] = -sinTheta; |
| mat.m_matrix[1][0] = 0.0; |
| mat.m_matrix[1][1] = 1.0; |
| mat.m_matrix[1][2] = 0.0; |
| mat.m_matrix[2][0] = sinTheta; |
| mat.m_matrix[2][1] = 0.0; |
| mat.m_matrix[2][2] = cosTheta; |
| mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; |
| mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; |
| mat.m_matrix[3][3] = 1.0; |
| |
| rmat.multiply(mat); |
| |
| sinTheta = std::sin(rx); |
| cosTheta = std::cos(rx); |
| |
| mat.m_matrix[0][0] = 1.0; |
| mat.m_matrix[0][1] = 0.0; |
| mat.m_matrix[0][2] = 0.0; |
| mat.m_matrix[1][0] = 0.0; |
| mat.m_matrix[1][1] = cosTheta; |
| mat.m_matrix[1][2] = sinTheta; |
| mat.m_matrix[2][0] = 0.0; |
| mat.m_matrix[2][1] = -sinTheta; |
| mat.m_matrix[2][2] = cosTheta; |
| mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; |
| mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; |
| mat.m_matrix[3][3] = 1.0; |
| |
| rmat.multiply(mat); |
| |
| multiply(rmat); |
| return *this; |
| } |
| |
| TransformationMatrix& TransformationMatrix::translate(double tx, double ty) { |
| m_matrix[3][0] += tx * m_matrix[0][0] + ty * m_matrix[1][0]; |
| m_matrix[3][1] += tx * m_matrix[0][1] + ty * m_matrix[1][1]; |
| m_matrix[3][2] += tx * m_matrix[0][2] + ty * m_matrix[1][2]; |
| m_matrix[3][3] += tx * m_matrix[0][3] + ty * m_matrix[1][3]; |
| return *this; |
| } |
| |
| TransformationMatrix& TransformationMatrix::translate3d(double tx, |
| double ty, |
| double tz) { |
| m_matrix[3][0] += |
| tx * m_matrix[0][0] + ty * m_matrix[1][0] + tz * m_matrix[2][0]; |
| m_matrix[3][1] += |
| tx * m_matrix[0][1] + ty * m_matrix[1][1] + tz * m_matrix[2][1]; |
| m_matrix[3][2] += |
| tx * m_matrix[0][2] + ty * m_matrix[1][2] + tz * m_matrix[2][2]; |
| m_matrix[3][3] += |
| tx * m_matrix[0][3] + ty * m_matrix[1][3] + tz * m_matrix[2][3]; |
| return *this; |
| } |
| |
| TransformationMatrix& TransformationMatrix::translateRight(double tx, |
| double ty) { |
| if (tx != 0) { |
| m_matrix[0][0] += m_matrix[0][3] * tx; |
| m_matrix[1][0] += m_matrix[1][3] * tx; |
| m_matrix[2][0] += m_matrix[2][3] * tx; |
| m_matrix[3][0] += m_matrix[3][3] * tx; |
| } |
| |
| if (ty != 0) { |
| m_matrix[0][1] += m_matrix[0][3] * ty; |
| m_matrix[1][1] += m_matrix[1][3] * ty; |
| m_matrix[2][1] += m_matrix[2][3] * ty; |
| m_matrix[3][1] += m_matrix[3][3] * ty; |
| } |
| |
| return *this; |
| } |
| |
| TransformationMatrix& TransformationMatrix::translateRight3d(double tx, |
| double ty, |
| double tz) { |
| translateRight(tx, ty); |
| if (tz != 0) { |
| m_matrix[0][2] += m_matrix[0][3] * tz; |
| m_matrix[1][2] += m_matrix[1][3] * tz; |
| m_matrix[2][2] += m_matrix[2][3] * tz; |
| m_matrix[3][2] += m_matrix[3][3] * tz; |
| } |
| |
| return *this; |
| } |
| |
| TransformationMatrix& TransformationMatrix::skew(double sx, double sy) { |
| // angles are in degrees. Switch to radians |
| sx = deg2rad(sx); |
| sy = deg2rad(sy); |
| |
| TransformationMatrix mat; |
| mat.m_matrix[0][1] = |
| std::tan(sy); // note that the y shear goes in the first row |
| mat.m_matrix[1][0] = std::tan(sx); // and the x shear in the second row |
| |
| multiply(mat); |
| return *this; |
| } |
| |
| TransformationMatrix& TransformationMatrix::applyPerspective(double p) { |
| TransformationMatrix mat; |
| if (p != 0) |
| mat.m_matrix[2][3] = -1 / p; |
| |
| multiply(mat); |
| return *this; |
| } |
| |
| TransformationMatrix& TransformationMatrix::applyTransformOrigin(double x, |
| double y, |
| double z) { |
| translateRight3d(x, y, z); |
| translate3d(-x, -y, -z); |
| return *this; |
| } |
| |
| TransformationMatrix& TransformationMatrix::zoom(double zoomFactor) { |
| m_matrix[0][3] /= zoomFactor; |
| m_matrix[1][3] /= zoomFactor; |
| m_matrix[2][3] /= zoomFactor; |
| m_matrix[3][0] *= zoomFactor; |
| m_matrix[3][1] *= zoomFactor; |
| m_matrix[3][2] *= zoomFactor; |
| return *this; |
| } |
| |
| // Calculates *this = *this * mat. |
| // Note: A * B means that the transforms represented by A happen first, and |
| // then the transforms represented by B. That is, the matrix A * B corresponds |
| // to a CSS transform list <transform-function-A> <transform-function-B>. |
| // Some branches of this function may make use of the fact that |
| // transpose(A * B) == transpose(B) * transpose(A); remember that |
| // m_matrix[a][b] is matrix element row b, col a. |
| // FIXME: As of 2016-05-04, the ARM64 branch is NOT triggered by tests on the CQ |
| // bots, see crbug.com/477892 and crbug.com/584508. |
| TransformationMatrix& TransformationMatrix::multiply( |
| const TransformationMatrix& mat) { |
| #if CPU(ARM64) |
| double* rightMatrix = &(m_matrix[0][0]); |
| const double* leftMatrix = &(mat.m_matrix[0][0]); |
| asm volatile( |
| // Load mat.m_matrix to v16 - v23. |
| // Load this.m_matrix to v24 - v31. |
| // Result: this = mat * this |
| // | v0, v1 | | v16, v17 | | v24, v25 | |
| // | v2, v3 | = | v18, v19 | * | v26, v27 | |
| // | v4, v5 | | v20, v21 | | v28, v29 | |
| // | v6, v7 | | v22, v23 | | v30, v31 | |
| "mov x9, %[rightMatrix] \t\n" |
| "ld1 {v16.2d - v19.2d}, [%[leftMatrix]], 64 \t\n" |
| "ld1 {v20.2d - v23.2d}, [%[leftMatrix]] \t\n" |
| "ld1 {v24.2d - v27.2d}, [%[rightMatrix]], 64 \t\n" |
| "ld1 {v28.2d - v31.2d}, [%[rightMatrix]] \t\n" |
| |
| "fmul v0.2d, v24.2d, v16.d[0] \t\n" |
| "fmul v1.2d, v25.2d, v16.d[0] \t\n" |
| "fmul v2.2d, v24.2d, v18.d[0] \t\n" |
| "fmul v3.2d, v25.2d, v18.d[0] \t\n" |
| "fmul v4.2d, v24.2d, v20.d[0] \t\n" |
| "fmul v5.2d, v25.2d, v20.d[0] \t\n" |
| "fmul v6.2d, v24.2d, v22.d[0] \t\n" |
| "fmul v7.2d, v25.2d, v22.d[0] \t\n" |
| |
| "fmla v0.2d, v26.2d, v16.d[1] \t\n" |
| "fmla v1.2d, v27.2d, v16.d[1] \t\n" |
| "fmla v2.2d, v26.2d, v18.d[1] \t\n" |
| "fmla v3.2d, v27.2d, v18.d[1] \t\n" |
| "fmla v4.2d, v26.2d, v20.d[1] \t\n" |
| "fmla v5.2d, v27.2d, v20.d[1] \t\n" |
| "fmla v6.2d, v26.2d, v22.d[1] \t\n" |
| "fmla v7.2d, v27.2d, v22.d[1] \t\n" |
| |
| "fmla v0.2d, v28.2d, v17.d[0] \t\n" |
| "fmla v1.2d, v29.2d, v17.d[0] \t\n" |
| "fmla v2.2d, v28.2d, v19.d[0] \t\n" |
| "fmla v3.2d, v29.2d, v19.d[0] \t\n" |
| "fmla v4.2d, v28.2d, v21.d[0] \t\n" |
| "fmla v5.2d, v29.2d, v21.d[0] \t\n" |
| "fmla v6.2d, v28.2d, v23.d[0] \t\n" |
| "fmla v7.2d, v29.2d, v23.d[0] \t\n" |
| |
| "fmla v0.2d, v30.2d, v17.d[1] \t\n" |
| "fmla v1.2d, v31.2d, v17.d[1] \t\n" |
| "fmla v2.2d, v30.2d, v19.d[1] \t\n" |
| "fmla v3.2d, v31.2d, v19.d[1] \t\n" |
| "fmla v4.2d, v30.2d, v21.d[1] \t\n" |
| "fmla v5.2d, v31.2d, v21.d[1] \t\n" |
| "fmla v6.2d, v30.2d, v23.d[1] \t\n" |
| "fmla v7.2d, v31.2d, v23.d[1] \t\n" |
| |
| "st1 {v0.2d - v3.2d}, [x9], 64 \t\n" |
| "st1 {v4.2d - v7.2d}, [x9] \t\n" |
| : [leftMatrix] "+r"(leftMatrix), [rightMatrix] "+r"(rightMatrix) |
| : |
| : "memory", "x9", "v16", "v17", "v18", "v19", "v20", "v21", "v22", "v23", |
| "v24", "v25", "v26", "v27", "v28", "v29", "v30", "v31", "v0", "v1", |
| "v2", "v3", "v4", "v5", "v6", "v7"); |
| #elif HAVE(MIPS_MSA_INTRINSICS) |
| v2f64 vleftM0, vleftM1, vleftM2, vleftM3, vleftM4, vleftM5, vleftM6, vleftM7; |
| v2f64 vRightM0, vRightM1, vRightM2, vRightM3, vRightM4, vRightM5, vRightM6, |
| vRightM7; |
| v2f64 vTmpM0, vTmpM1, vTmpM2, vTmpM3; |
| |
| vRightM0 = LD_DP(&(m_matrix[0][0])); |
| vRightM1 = LD_DP(&(m_matrix[0][2])); |
| vRightM2 = LD_DP(&(m_matrix[1][0])); |
| vRightM3 = LD_DP(&(m_matrix[1][2])); |
| vRightM4 = LD_DP(&(m_matrix[2][0])); |
| vRightM5 = LD_DP(&(m_matrix[2][2])); |
| vRightM6 = LD_DP(&(m_matrix[3][0])); |
| vRightM7 = LD_DP(&(m_matrix[3][2])); |
| |
| vleftM0 = LD_DP(&(mat.m_matrix[0][0])); |
| vleftM2 = LD_DP(&(mat.m_matrix[0][2])); |
| vleftM4 = LD_DP(&(mat.m_matrix[1][0])); |
| vleftM6 = LD_DP(&(mat.m_matrix[1][2])); |
| |
| vleftM1 = (v2f64)__msa_splati_d((v2i64)vleftM0, 1); |
| vleftM0 = (v2f64)__msa_splati_d((v2i64)vleftM0, 0); |
| vleftM3 = (v2f64)__msa_splati_d((v2i64)vleftM2, 1); |
| vleftM2 = (v2f64)__msa_splati_d((v2i64)vleftM2, 0); |
| vleftM5 = (v2f64)__msa_splati_d((v2i64)vleftM4, 1); |
| vleftM4 = (v2f64)__msa_splati_d((v2i64)vleftM4, 0); |
| vleftM7 = (v2f64)__msa_splati_d((v2i64)vleftM6, 1); |
| vleftM6 = (v2f64)__msa_splati_d((v2i64)vleftM6, 0); |
| |
| vTmpM0 = vleftM0 * vRightM0; |
| vTmpM1 = vleftM0 * vRightM1; |
| vTmpM0 += vleftM1 * vRightM2; |
| vTmpM1 += vleftM1 * vRightM3; |
| vTmpM0 += vleftM2 * vRightM4; |
| vTmpM1 += vleftM2 * vRightM5; |
| vTmpM0 += vleftM3 * vRightM6; |
| vTmpM1 += vleftM3 * vRightM7; |
| |
| vTmpM2 = vleftM4 * vRightM0; |
| vTmpM3 = vleftM4 * vRightM1; |
| vTmpM2 += vleftM5 * vRightM2; |
| vTmpM3 += vleftM5 * vRightM3; |
| vTmpM2 += vleftM6 * vRightM4; |
| vTmpM3 += vleftM6 * vRightM5; |
| vTmpM2 += vleftM7 * vRightM6; |
| vTmpM3 += vleftM7 * vRightM7; |
| |
| vleftM0 = LD_DP(&(mat.m_matrix[2][0])); |
| vleftM2 = LD_DP(&(mat.m_matrix[2][2])); |
| vleftM4 = LD_DP(&(mat.m_matrix[3][0])); |
| vleftM6 = LD_DP(&(mat.m_matrix[3][2])); |
| |
| ST_DP(vTmpM0, &(m_matrix[0][0])); |
| ST_DP(vTmpM1, &(m_matrix[0][2])); |
| ST_DP(vTmpM2, &(m_matrix[1][0])); |
| ST_DP(vTmpM3, &(m_matrix[1][2])); |
| |
| vleftM1 = (v2f64)__msa_splati_d((v2i64)vleftM0, 1); |
| vleftM0 = (v2f64)__msa_splati_d((v2i64)vleftM0, 0); |
| vleftM3 = (v2f64)__msa_splati_d((v2i64)vleftM2, 1); |
| vleftM2 = (v2f64)__msa_splati_d((v2i64)vleftM2, 0); |
| vleftM5 = (v2f64)__msa_splati_d((v2i64)vleftM4, 1); |
| vleftM4 = (v2f64)__msa_splati_d((v2i64)vleftM4, 0); |
| vleftM7 = (v2f64)__msa_splati_d((v2i64)vleftM6, 1); |
| vleftM6 = (v2f64)__msa_splati_d((v2i64)vleftM6, 0); |
| |
| vTmpM0 = vleftM0 * vRightM0; |
| vTmpM1 = vleftM0 * vRightM1; |
| vTmpM0 += vleftM1 * vRightM2; |
| vTmpM1 += vleftM1 * vRightM3; |
| vTmpM0 += vleftM2 * vRightM4; |
| vTmpM1 += vleftM2 * vRightM5; |
| vTmpM0 += vleftM3 * vRightM6; |
| vTmpM1 += vleftM3 * vRightM7; |
| |
| vTmpM2 = vleftM4 * vRightM0; |
| vTmpM3 = vleftM4 * vRightM1; |
| vTmpM2 += vleftM5 * vRightM2; |
| vTmpM3 += vleftM5 * vRightM3; |
| vTmpM2 += vleftM6 * vRightM4; |
| vTmpM3 += vleftM6 * vRightM5; |
| vTmpM2 += vleftM7 * vRightM6; |
| vTmpM3 += vleftM7 * vRightM7; |
| |
| ST_DP(vTmpM0, &(m_matrix[2][0])); |
| ST_DP(vTmpM1, &(m_matrix[2][2])); |
| ST_DP(vTmpM2, &(m_matrix[3][0])); |
| ST_DP(vTmpM3, &(m_matrix[3][2])); |
| #elif defined(TRANSFORMATION_MATRIX_USE_X86_64_SSE2) |
| // x86_64 has 16 XMM registers which is enough to do the multiplication fully |
| // in registers. |
| __m128d matrixBlockA = _mm_load_pd(&(m_matrix[0][0])); |
| __m128d matrixBlockC = _mm_load_pd(&(m_matrix[1][0])); |
| __m128d matrixBlockE = _mm_load_pd(&(m_matrix[2][0])); |
| __m128d matrixBlockG = _mm_load_pd(&(m_matrix[3][0])); |
| |
| // First row. |
| __m128d otherMatrixFirstParam = _mm_set1_pd(mat.m_matrix[0][0]); |
| __m128d otherMatrixSecondParam = _mm_set1_pd(mat.m_matrix[0][1]); |
| __m128d otherMatrixThirdParam = _mm_set1_pd(mat.m_matrix[0][2]); |
| __m128d otherMatrixFourthParam = _mm_set1_pd(mat.m_matrix[0][3]); |
| |
| // output00 and output01. |
| __m128d accumulator = _mm_mul_pd(matrixBlockA, otherMatrixFirstParam); |
| __m128d temp1 = _mm_mul_pd(matrixBlockC, otherMatrixSecondParam); |
| __m128d temp2 = _mm_mul_pd(matrixBlockE, otherMatrixThirdParam); |
| __m128d temp3 = _mm_mul_pd(matrixBlockG, otherMatrixFourthParam); |
| |
| __m128d matrixBlockB = _mm_load_pd(&(m_matrix[0][2])); |
| __m128d matrixBlockD = _mm_load_pd(&(m_matrix[1][2])); |
| __m128d matrixBlockF = _mm_load_pd(&(m_matrix[2][2])); |
| __m128d matrixBlockH = _mm_load_pd(&(m_matrix[3][2])); |
| |
| accumulator = _mm_add_pd(accumulator, temp1); |
| accumulator = _mm_add_pd(accumulator, temp2); |
| accumulator = _mm_add_pd(accumulator, temp3); |
| _mm_store_pd(&m_matrix[0][0], accumulator); |
| |
| // output02 and output03. |
| accumulator = _mm_mul_pd(matrixBlockB, otherMatrixFirstParam); |
| temp1 = _mm_mul_pd(matrixBlockD, otherMatrixSecondParam); |
| temp2 = _mm_mul_pd(matrixBlockF, otherMatrixThirdParam); |
| temp3 = _mm_mul_pd(matrixBlockH, otherMatrixFourthParam); |
| |
| accumulator = _mm_add_pd(accumulator, temp1); |
| accumulator = _mm_add_pd(accumulator, temp2); |
| accumulator = _mm_add_pd(accumulator, temp3); |
| _mm_store_pd(&m_matrix[0][2], accumulator); |
| |
| // Second row. |
| otherMatrixFirstParam = _mm_set1_pd(mat.m_matrix[1][0]); |
| otherMatrixSecondParam = _mm_set1_pd(mat.m_matrix[1][1]); |
| otherMatrixThirdParam = _mm_set1_pd(mat.m_matrix[1][2]); |
| otherMatrixFourthParam = _mm_set1_pd(mat.m_matrix[1][3]); |
| |
| // output10 and output11. |
| accumulator = _mm_mul_pd(matrixBlockA, otherMatrixFirstParam); |
| temp1 = _mm_mul_pd(matrixBlockC, otherMatrixSecondParam); |
| temp2 = _mm_mul_pd(matrixBlockE, otherMatrixThirdParam); |
| temp3 = _mm_mul_pd(matrixBlockG, otherMatrixFourthParam); |
| |
| accumulator = _mm_add_pd(accumulator, temp1); |
| accumulator = _mm_add_pd(accumulator, temp2); |
| accumulator = _mm_add_pd(accumulator, temp3); |
| _mm_store_pd(&m_matrix[1][0], accumulator); |
| |
| // output12 and output13. |
| accumulator = _mm_mul_pd(matrixBlockB, otherMatrixFirstParam); |
| temp1 = _mm_mul_pd(matrixBlockD, otherMatrixSecondParam); |
| temp2 = _mm_mul_pd(matrixBlockF, otherMatrixThirdParam); |
| temp3 = _mm_mul_pd(matrixBlockH, otherMatrixFourthParam); |
| |
| accumulator = _mm_add_pd(accumulator, temp1); |
| accumulator = _mm_add_pd(accumulator, temp2); |
| accumulator = _mm_add_pd(accumulator, temp3); |
| _mm_store_pd(&m_matrix[1][2], accumulator); |
| |
| // Third row. |
| otherMatrixFirstParam = _mm_set1_pd(mat.m_matrix[2][0]); |
| otherMatrixSecondParam = _mm_set1_pd(mat.m_matrix[2][1]); |
| otherMatrixThirdParam = _mm_set1_pd(mat.m_matrix[2][2]); |
| otherMatrixFourthParam = _mm_set1_pd(mat.m_matrix[2][3]); |
| |
| // output20 and output21. |
| accumulator = _mm_mul_pd(matrixBlockA, otherMatrixFirstParam); |
| temp1 = _mm_mul_pd(matrixBlockC, otherMatrixSecondParam); |
| temp2 = _mm_mul_pd(matrixBlockE, otherMatrixThirdParam); |
| temp3 = _mm_mul_pd(matrixBlockG, otherMatrixFourthParam); |
| |
| accumulator = _mm_add_pd(accumulator, temp1); |
| accumulator = _mm_add_pd(accumulator, temp2); |
| accumulator = _mm_add_pd(accumulator, temp3); |
| _mm_store_pd(&m_matrix[2][0], accumulator); |
| |
| // output22 and output23. |
| accumulator = _mm_mul_pd(matrixBlockB, otherMatrixFirstParam); |
| temp1 = _mm_mul_pd(matrixBlockD, otherMatrixSecondParam); |
| temp2 = _mm_mul_pd(matrixBlockF, otherMatrixThirdParam); |
| temp3 = _mm_mul_pd(matrixBlockH, otherMatrixFourthParam); |
| |
| accumulator = _mm_add_pd(accumulator, temp1); |
| accumulator = _mm_add_pd(accumulator, temp2); |
| accumulator = _mm_add_pd(accumulator, temp3); |
| _mm_store_pd(&m_matrix[2][2], accumulator); |
| |
| // Fourth row. |
| otherMatrixFirstParam = _mm_set1_pd(mat.m_matrix[3][0]); |
| otherMatrixSecondParam = _mm_set1_pd(mat.m_matrix[3][1]); |
| otherMatrixThirdParam = _mm_set1_pd(mat.m_matrix[3][2]); |
| otherMatrixFourthParam = _mm_set1_pd(mat.m_matrix[3][3]); |
| |
| // output30 and output31. |
| accumulator = _mm_mul_pd(matrixBlockA, otherMatrixFirstParam); |
| temp1 = _mm_mul_pd(matrixBlockC, otherMatrixSecondParam); |
| temp2 = _mm_mul_pd(matrixBlockE, otherMatrixThirdParam); |
| temp3 = _mm_mul_pd(matrixBlockG, otherMatrixFourthParam); |
| |
| accumulator = _mm_add_pd(accumulator, temp1); |
| accumulator = _mm_add_pd(accumulator, temp2); |
| accumulator = _mm_add_pd(accumulator, temp3); |
| _mm_store_pd(&m_matrix[3][0], accumulator); |
| |
| // output32 and output33. |
| accumulator = _mm_mul_pd(matrixBlockB, otherMatrixFirstParam); |
| temp1 = _mm_mul_pd(matrixBlockD, otherMatrixSecondParam); |
| temp2 = _mm_mul_pd(matrixBlockF, otherMatrixThirdParam); |
| temp3 = _mm_mul_pd(matrixBlockH, otherMatrixFourthParam); |
| |
| accumulator = _mm_add_pd(accumulator, temp1); |
| accumulator = _mm_add_pd(accumulator, temp2); |
| accumulator = _mm_add_pd(accumulator, temp3); |
| _mm_store_pd(&m_matrix[3][2], accumulator); |
| #else |
| Matrix4 tmp; |
| |
| tmp[0][0] = (mat.m_matrix[0][0] * m_matrix[0][0] + |
| mat.m_matrix[0][1] * m_matrix[1][0] + |
| mat.m_matrix[0][2] * m_matrix[2][0] + |
| mat.m_matrix[0][3] * m_matrix[3][0]); |
| tmp[0][1] = (mat.m_matrix[0][0] * m_matrix[0][1] + |
| mat.m_matrix[0][1] * m_matrix[1][1] + |
| mat.m_matrix[0][2] * m_matrix[2][1] + |
| mat.m_matrix[0][3] * m_matrix[3][1]); |
| tmp[0][2] = (mat.m_matrix[0][0] * m_matrix[0][2] + |
| mat.m_matrix[0][1] * m_matrix[1][2] + |
| mat.m_matrix[0][2] * m_matrix[2][2] + |
| mat.m_matrix[0][3] * m_matrix[3][2]); |
| tmp[0][3] = (mat.m_matrix[0][0] * m_matrix[0][3] + |
| mat.m_matrix[0][1] * m_matrix[1][3] + |
| mat.m_matrix[0][2] * m_matrix[2][3] + |
| mat.m_matrix[0][3] * m_matrix[3][3]); |
| |
| tmp[1][0] = (mat.m_matrix[1][0] * m_matrix[0][0] + |
| mat.m_matrix[1][1] * m_matrix[1][0] + |
| mat.m_matrix[1][2] * m_matrix[2][0] + |
| mat.m_matrix[1][3] * m_matrix[3][0]); |
| tmp[1][1] = (mat.m_matrix[1][0] * m_matrix[0][1] + |
| mat.m_matrix[1][1] * m_matrix[1][1] + |
| mat.m_matrix[1][2] * m_matrix[2][1] + |
| mat.m_matrix[1][3] * m_matrix[3][1]); |
| tmp[1][2] = (mat.m_matrix[1][0] * m_matrix[0][2] + |
| mat.m_matrix[1][1] * m_matrix[1][2] + |
| mat.m_matrix[1][2] * m_matrix[2][2] + |
| mat.m_matrix[1][3] * m_matrix[3][2]); |
| tmp[1][3] = (mat.m_matrix[1][0] * m_matrix[0][3] + |
| mat.m_matrix[1][1] * m_matrix[1][3] + |
| mat.m_matrix[1][2] * m_matrix[2][3] + |
| mat.m_matrix[1][3] * m_matrix[3][3]); |
| |
| tmp[2][0] = (mat.m_matrix[2][0] * m_matrix[0][0] + |
| mat.m_matrix[2][1] * m_matrix[1][0] + |
| mat.m_matrix[2][2] * m_matrix[2][0] + |
| mat.m_matrix[2][3] * m_matrix[3][0]); |
| tmp[2][1] = (mat.m_matrix[2][0] * m_matrix[0][1] + |
| mat.m_matrix[2][1] * m_matrix[1][1] + |
| mat.m_matrix[2][2] * m_matrix[2][1] + |
| mat.m_matrix[2][3] * m_matrix[3][1]); |
| tmp[2][2] = (mat.m_matrix[2][0] * m_matrix[0][2] + |
| mat.m_matrix[2][1] * m_matrix[1][2] + |
| mat.m_matrix[2][2] * m_matrix[2][2] + |
| mat.m_matrix[2][3] * m_matrix[3][2]); |
| tmp[2][3] = (mat.m_matrix[2][0] * m_matrix[0][3] + |
| mat.m_matrix[2][1] * m_matrix[1][3] + |
| mat.m_matrix[2][2] * m_matrix[2][3] + |
| mat.m_matrix[2][3] * m_matrix[3][3]); |
| |
| tmp[3][0] = (mat.m_matrix[3][0] * m_matrix[0][0] + |
| mat.m_matrix[3][1] * m_matrix[1][0] + |
| mat.m_matrix[3][2] * m_matrix[2][0] + |
| mat.m_matrix[3][3] * m_matrix[3][0]); |
| tmp[3][1] = (mat.m_matrix[3][0] * m_matrix[0][1] + |
| mat.m_matrix[3][1] * m_matrix[1][1] + |
| mat.m_matrix[3][2] * m_matrix[2][1] + |
| mat.m_matrix[3][3] * m_matrix[3][1]); |
| tmp[3][2] = (mat.m_matrix[3][0] * m_matrix[0][2] + |
| mat.m_matrix[3][1] * m_matrix[1][2] + |
| mat.m_matrix[3][2] * m_matrix[2][2] + |
| mat.m_matrix[3][3] * m_matrix[3][2]); |
| tmp[3][3] = (mat.m_matrix[3][0] * m_matrix[0][3] + |
| mat.m_matrix[3][1] * m_matrix[1][3] + |
| mat.m_matrix[3][2] * m_matrix[2][3] + |
| mat.m_matrix[3][3] * m_matrix[3][3]); |
| |
| setMatrix(tmp); |
| #endif |
| return *this; |
| } |
| |
| void TransformationMatrix::multVecMatrix(double x, |
| double y, |
| double& resultX, |
| double& resultY) const { |
| resultX = m_matrix[3][0] + x * m_matrix[0][0] + y * m_matrix[1][0]; |
| resultY = m_matrix[3][1] + x * m_matrix[0][1] + y * m_matrix[1][1]; |
| double w = m_matrix[3][3] + x * m_matrix[0][3] + y * m_matrix[1][3]; |
| if (w != 1 && w != 0) { |
| resultX /= w; |
| resultY /= w; |
| } |
| } |
| |
| void TransformationMatrix::multVecMatrix(double x, |
| double y, |
| double z, |
| double& resultX, |
| double& resultY, |
| double& resultZ) const { |
| resultX = m_matrix[3][0] + x * m_matrix[0][0] + y * m_matrix[1][0] + |
| z * m_matrix[2][0]; |
| resultY = m_matrix[3][1] + x * m_matrix[0][1] + y * m_matrix[1][1] + |
| z * m_matrix[2][1]; |
| resultZ = m_matrix[3][2] + x * m_matrix[0][2] + y * m_matrix[1][2] + |
| z * m_matrix[2][2]; |
| double w = m_matrix[3][3] + x * m_matrix[0][3] + y * m_matrix[1][3] + |
| z * m_matrix[2][3]; |
| if (w != 1 && w != 0) { |
| resultX /= w; |
| resultY /= w; |
| resultZ /= w; |
| } |
| } |
| |
| bool TransformationMatrix::isInvertible() const { |
| if (isIdentityOrTranslation()) |
| return true; |
| |
| double det = blink::determinant4x4(m_matrix); |
| |
| if (fabs(det) < SmallNumber) |
| return false; |
| |
| return true; |
| } |
| |
| TransformationMatrix TransformationMatrix::inverse() const { |
| if (isIdentityOrTranslation()) { |
| // identity matrix |
| if (m_matrix[3][0] == 0 && m_matrix[3][1] == 0 && m_matrix[3][2] == 0) |
| return TransformationMatrix(); |
| |
| // translation |
| return TransformationMatrix(1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, |
| -m_matrix[3][0], -m_matrix[3][1], |
| -m_matrix[3][2], 1); |
| } |
| |
| TransformationMatrix invMat; |
| bool inverted = blink::inverse(m_matrix, invMat.m_matrix); |
| if (!inverted) |
| return TransformationMatrix(); |
| |
| return invMat; |
| } |
| |
| void TransformationMatrix::makeAffine() { |
| m_matrix[0][2] = 0; |
| m_matrix[0][3] = 0; |
| |
| m_matrix[1][2] = 0; |
| m_matrix[1][3] = 0; |
| |
| m_matrix[2][0] = 0; |
| m_matrix[2][1] = 0; |
| m_matrix[2][2] = 1; |
| m_matrix[2][3] = 0; |
| |
| m_matrix[3][2] = 0; |
| m_matrix[3][3] = 1; |
| } |
| |
| AffineTransform TransformationMatrix::toAffineTransform() const { |
| return AffineTransform(m_matrix[0][0], m_matrix[0][1], m_matrix[1][0], |
| m_matrix[1][1], m_matrix[3][0], m_matrix[3][1]); |
| } |
| |
| void TransformationMatrix::flattenTo2d() { |
| m_matrix[2][0] = 0; |
| m_matrix[2][1] = 0; |
| m_matrix[0][2] = 0; |
| m_matrix[1][2] = 0; |
| m_matrix[2][2] = 1; |
| m_matrix[3][2] = 0; |
| m_matrix[2][3] = 0; |
| } |
| |
| static inline void blendFloat(double& from, double to, double progress) { |
| if (from != to) |
| from = from + (to - from) * progress; |
| } |
| |
| void TransformationMatrix::blend(const TransformationMatrix& from, |
| double progress) { |
| if (from.isIdentity() && isIdentity()) |
| return; |
| |
| // decompose |
| DecomposedType fromDecomp; |
| DecomposedType toDecomp; |
| if (!from.decompose(fromDecomp) || !decompose(toDecomp)) { |
| if (progress < 0.5) |
| *this = from; |
| return; |
| } |
| |
| // interpolate |
| blendFloat(fromDecomp.scaleX, toDecomp.scaleX, progress); |
| blendFloat(fromDecomp.scaleY, toDecomp.scaleY, progress); |
| blendFloat(fromDecomp.scaleZ, toDecomp.scaleZ, progress); |
| blendFloat(fromDecomp.skewXY, toDecomp.skewXY, progress); |
| blendFloat(fromDecomp.skewXZ, toDecomp.skewXZ, progress); |
| blendFloat(fromDecomp.skewYZ, toDecomp.skewYZ, progress); |
| blendFloat(fromDecomp.translateX, toDecomp.translateX, progress); |
| blendFloat(fromDecomp.translateY, toDecomp.translateY, progress); |
| blendFloat(fromDecomp.translateZ, toDecomp.translateZ, progress); |
| blendFloat(fromDecomp.perspectiveX, toDecomp.perspectiveX, progress); |
| blendFloat(fromDecomp.perspectiveY, toDecomp.perspectiveY, progress); |
| blendFloat(fromDecomp.perspectiveZ, toDecomp.perspectiveZ, progress); |
| blendFloat(fromDecomp.perspectiveW, toDecomp.perspectiveW, progress); |
| |
| slerp(&fromDecomp.quaternionX, &toDecomp.quaternionX, progress); |
| |
| // recompose |
| recompose(fromDecomp); |
| } |
| |
| bool TransformationMatrix::decompose(DecomposedType& decomp) const { |
| if (isIdentity()) { |
| memset(&decomp, 0, sizeof(decomp)); |
| decomp.perspectiveW = 1; |
| decomp.scaleX = 1; |
| decomp.scaleY = 1; |
| decomp.scaleZ = 1; |
| } |
| |
| if (!blink::decompose(m_matrix, decomp)) |
| return false; |
| return true; |
| } |
| |
| void TransformationMatrix::recompose(const DecomposedType& decomp) { |
| makeIdentity(); |
| |
| // first apply perspective |
| m_matrix[0][3] = decomp.perspectiveX; |
| m_matrix[1][3] = decomp.perspectiveY; |
| m_matrix[2][3] = decomp.perspectiveZ; |
| m_matrix[3][3] = decomp.perspectiveW; |
| |
| // now translate |
| translate3d(decomp.translateX, decomp.translateY, decomp.translateZ); |
| |
| // apply rotation |
| double xx = decomp.quaternionX * decomp.quaternionX; |
| double xy = decomp.quaternionX * decomp.quaternionY; |
| double xz = decomp.quaternionX * decomp.quaternionZ; |
| double xw = decomp.quaternionX * decomp.quaternionW; |
| double yy = decomp.quaternionY * decomp.quaternionY; |
| double yz = decomp.quaternionY * decomp.quaternionZ; |
| double yw = decomp.quaternionY * decomp.quaternionW; |
| double zz = decomp.quaternionZ * decomp.quaternionZ; |
| double zw = decomp.quaternionZ * decomp.quaternionW; |
| |
| // Construct a composite rotation matrix from the quaternion values |
| TransformationMatrix rotationMatrix( |
| 1 - 2 * (yy + zz), 2 * (xy - zw), 2 * (xz + yw), 0, 2 * (xy + zw), |
| 1 - 2 * (xx + zz), 2 * (yz - xw), 0, 2 * (xz - yw), 2 * (yz + xw), |
| 1 - 2 * (xx + yy), 0, 0, 0, 0, 1); |
| |
| multiply(rotationMatrix); |
| |
| // now apply skew |
| if (decomp.skewYZ) { |
| TransformationMatrix tmp; |
| tmp.setM32(decomp.skewYZ); |
| multiply(tmp); |
| } |
| |
| if (decomp.skewXZ) { |
| TransformationMatrix tmp; |
| tmp.setM31(decomp.skewXZ); |
| multiply(tmp); |
| } |
| |
| if (decomp.skewXY) { |
| TransformationMatrix tmp; |
| tmp.setM21(decomp.skewXY); |
| multiply(tmp); |
| } |
| |
| // finally, apply scale |
| scale3d(decomp.scaleX, decomp.scaleY, decomp.scaleZ); |
| } |
| |
| bool TransformationMatrix::isIntegerTranslation() const { |
| if (!isIdentityOrTranslation()) |
| return false; |
| |
| // Check for translate Z. |
| if (m_matrix[3][2]) |
| return false; |
| |
| // Check for non-integer translate X/Y. |
| if (static_cast<int>(m_matrix[3][0]) != m_matrix[3][0] || |
| static_cast<int>(m_matrix[3][1]) != m_matrix[3][1]) |
| return false; |
| |
| return true; |
| } |
| |
| FloatSize TransformationMatrix::to2DTranslation() const { |
| ASSERT(isIdentityOr2DTranslation()); |
| return FloatSize(m_matrix[3][0], m_matrix[3][1]); |
| } |
| |
| void TransformationMatrix::toColumnMajorFloatArray(FloatMatrix4& result) const { |
| result[0] = m11(); |
| result[1] = m12(); |
| result[2] = m13(); |
| result[3] = m14(); |
| result[4] = m21(); |
| result[5] = m22(); |
| result[6] = m23(); |
| result[7] = m24(); |
| result[8] = m31(); |
| result[9] = m32(); |
| result[10] = m33(); |
| result[11] = m34(); |
| result[12] = m41(); |
| result[13] = m42(); |
| result[14] = m43(); |
| result[15] = m44(); |
| } |
| |
| SkMatrix44 TransformationMatrix::toSkMatrix44( |
| const TransformationMatrix& matrix) { |
| SkMatrix44 ret(SkMatrix44::kUninitialized_Constructor); |
| ret.setDouble(0, 0, matrix.m11()); |
| ret.setDouble(0, 1, matrix.m21()); |
| ret.setDouble(0, 2, matrix.m31()); |
| ret.setDouble(0, 3, matrix.m41()); |
| ret.setDouble(1, 0, matrix.m12()); |
| ret.setDouble(1, 1, matrix.m22()); |
| ret.setDouble(1, 2, matrix.m32()); |
| ret.setDouble(1, 3, matrix.m42()); |
| ret.setDouble(2, 0, matrix.m13()); |
| ret.setDouble(2, 1, matrix.m23()); |
| ret.setDouble(2, 2, matrix.m33()); |
| ret.setDouble(2, 3, matrix.m43()); |
| ret.setDouble(3, 0, matrix.m14()); |
| ret.setDouble(3, 1, matrix.m24()); |
| ret.setDouble(3, 2, matrix.m34()); |
| ret.setDouble(3, 3, matrix.m44()); |
| return ret; |
| } |
| |
| String TransformationMatrix::toString(bool asMatrix) const { |
| if (asMatrix) { |
| // Return as a matrix in row-major order. |
| return String::format( |
| "[%lg,%lg,%lg,%lg,\n%lg,%lg,%lg,%lg,\n%lg,%lg,%lg,%lg,\n%lg,%lg,%lg,%" |
| "lg]", |
| m11(), m21(), m31(), m41(), m12(), m22(), m32(), m42(), m13(), m23(), |
| m33(), m43(), m14(), m24(), m34(), m44()); |
| } |
| |
| TransformationMatrix::DecomposedType decomposition; |
| if (!decompose(decomposition)) |
| return toString(true) + " (degenerate)"; |
| |
| if (isIdentityOrTranslation()) { |
| if (decomposition.translateX == 0 && decomposition.translateY == 0 && |
| decomposition.translateZ == 0) |
| return "identity"; |
| return String::format("translation(%lg,%lg,%lg)", decomposition.translateX, |
| decomposition.translateY, decomposition.translateZ); |
| } |
| |
| return String::format( |
| "translation(%lg,%lg,%lg), scale(%lg,%lg,%lg), skew(%lg,%lg,%lg), " |
| "quaternion(%lg,%lg,%lg,%lg), perspective(%lg,%lg,%lg,%lg)", |
| decomposition.translateX, decomposition.translateY, |
| decomposition.translateZ, decomposition.scaleX, decomposition.scaleY, |
| decomposition.scaleZ, decomposition.skewXY, decomposition.skewXZ, |
| decomposition.skewYZ, decomposition.quaternionX, |
| decomposition.quaternionY, decomposition.quaternionZ, |
| decomposition.quaternionW, decomposition.perspectiveX, |
| decomposition.perspectiveY, decomposition.perspectiveZ, |
| decomposition.perspectiveW); |
| } |
| |
| } // namespace blink |