1 /* 2 * Copyright (C) 2005, 2006 Apple Computer, Inc. All rights reserved. 3 * Copyright (C) 2009 Torch Mobile, Inc. 4 * Copyright (C) 2013 Google Inc. All rights reserved. 5 * 6 * Redistribution and use in source and binary forms, with or without 7 * modification, are permitted provided that the following conditions 8 * are met: 9 * 1. Redistributions of source code must retain the above copyright 10 * notice, this list of conditions and the following disclaimer. 11 * 2. Redistributions in binary form must reproduce the above copyright 12 * notice, this list of conditions and the following disclaimer in the 13 * documentation and/or other materials provided with the distribution. 14 * 15 * THIS SOFTWARE IS PROVIDED BY APPLE COMPUTER, INC. ``AS IS'' AND ANY 16 * EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 17 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR 18 * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL APPLE COMPUTER, INC. OR 19 * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, 20 * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, 21 * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR 22 * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY 23 * OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 24 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE 25 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 26 */ 27 28 #include "config.h" 29 #include "platform/transforms/TransformationMatrix.h" 30 31 #include "platform/geometry/FloatBox.h" 32 #include "platform/geometry/FloatQuad.h" 33 #include "platform/geometry/FloatRect.h" 34 #include "platform/geometry/IntRect.h" 35 #include "platform/geometry/LayoutRect.h" 36 #include "platform/transforms/AffineTransform.h" 37 38 #include "wtf/Assertions.h" 39 #include "wtf/MathExtras.h" 40 41 #if CPU(X86_64) 42 #include <emmintrin.h> 43 #endif 44 45 using namespace std; 46 47 namespace WebCore { 48 49 // 50 // Supporting Math Functions 51 // 52 // This is a set of function from various places (attributed inline) to do things like 53 // inversion and decomposition of a 4x4 matrix. They are used throughout the code 54 // 55 56 // 57 // Adapted from Matrix Inversion by Richard Carling, Graphics Gems <http://tog.acm.org/GraphicsGems/index.html>. 58 59 // EULA: The Graphics Gems code is copyright-protected. In other words, you cannot claim the text of the code 60 // as your own and resell it. Using the code is permitted in any program, product, or library, non-commercial 61 // or commercial. Giving credit is not required, though is a nice gesture. The code comes as-is, and if there 62 // are any flaws or problems with any Gems code, nobody involved with Gems - authors, editors, publishers, or 63 // webmasters - are to be held responsible. Basically, don't be a jerk, and remember that anything free comes 64 // with no guarantee. 65 66 // A clarification about the storage of matrix elements 67 // 68 // This class uses a 2 dimensional array internally to store the elements of the matrix. The first index into 69 // the array refers to the column that the element lies in; the second index refers to the row. 70 // 71 // In other words, this is the layout of the matrix: 72 // 73 // | m_matrix[0][0] m_matrix[1][0] m_matrix[2][0] m_matrix[3][0] | 74 // | m_matrix[0][1] m_matrix[1][1] m_matrix[2][1] m_matrix[3][1] | 75 // | m_matrix[0][2] m_matrix[1][2] m_matrix[2][2] m_matrix[3][2] | 76 // | m_matrix[0][3] m_matrix[1][3] m_matrix[2][3] m_matrix[3][3] | 77 78 typedef double Vector4[4]; 79 typedef double Vector3[3]; 80 81 const double SMALL_NUMBER = 1.e-8; 82 83 // inverse(original_matrix, inverse_matrix) 84 // 85 // calculate the inverse of a 4x4 matrix 86 // 87 // -1 88 // A = ___1__ adjoint A 89 // det A 90 91 // double = determinant2x2(double a, double b, double c, double d) 92 // 93 // calculate the determinant of a 2x2 matrix. 94 95 static double determinant2x2(double a, double b, double c, double d) 96 { 97 return a * d - b * c; 98 } 99 100 // double = determinant3x3(a1, a2, a3, b1, b2, b3, c1, c2, c3) 101 // 102 // Calculate the determinant of a 3x3 matrix 103 // in the form 104 // 105 // | a1, b1, c1 | 106 // | a2, b2, c2 | 107 // | a3, b3, c3 | 108 109 static double determinant3x3(double a1, double a2, double a3, double b1, double b2, double b3, double c1, double c2, double c3) 110 { 111 return a1 * determinant2x2(b2, b3, c2, c3) 112 - b1 * determinant2x2(a2, a3, c2, c3) 113 + c1 * determinant2x2(a2, a3, b2, b3); 114 } 115 116 // double = determinant4x4(matrix) 117 // 118 // calculate the determinant of a 4x4 matrix. 119 120 static double determinant4x4(const TransformationMatrix::Matrix4& m) 121 { 122 // Assign to individual variable names to aid selecting 123 // correct elements 124 125 double a1 = m[0][0]; 126 double b1 = m[0][1]; 127 double c1 = m[0][2]; 128 double d1 = m[0][3]; 129 130 double a2 = m[1][0]; 131 double b2 = m[1][1]; 132 double c2 = m[1][2]; 133 double d2 = m[1][3]; 134 135 double a3 = m[2][0]; 136 double b3 = m[2][1]; 137 double c3 = m[2][2]; 138 double d3 = m[2][3]; 139 140 double a4 = m[3][0]; 141 double b4 = m[3][1]; 142 double c4 = m[3][2]; 143 double d4 = m[3][3]; 144 145 return a1 * determinant3x3(b2, b3, b4, c2, c3, c4, d2, d3, d4) 146 - b1 * determinant3x3(a2, a3, a4, c2, c3, c4, d2, d3, d4) 147 + c1 * determinant3x3(a2, a3, a4, b2, b3, b4, d2, d3, d4) 148 - d1 * determinant3x3(a2, a3, a4, b2, b3, b4, c2, c3, c4); 149 } 150 151 // adjoint( original_matrix, inverse_matrix ) 152 // 153 // calculate the adjoint of a 4x4 matrix 154 // 155 // Let a denote the minor determinant of matrix A obtained by 156 // ij 157 // 158 // deleting the ith row and jth column from A. 159 // 160 // i+j 161 // Let b = (-1) a 162 // ij ji 163 // 164 // The matrix B = (b ) is the adjoint of A 165 // ij 166 167 static void adjoint(const TransformationMatrix::Matrix4& matrix, TransformationMatrix::Matrix4& result) 168 { 169 // Assign to individual variable names to aid 170 // selecting correct values 171 double a1 = matrix[0][0]; 172 double b1 = matrix[0][1]; 173 double c1 = matrix[0][2]; 174 double d1 = matrix[0][3]; 175 176 double a2 = matrix[1][0]; 177 double b2 = matrix[1][1]; 178 double c2 = matrix[1][2]; 179 double d2 = matrix[1][3]; 180 181 double a3 = matrix[2][0]; 182 double b3 = matrix[2][1]; 183 double c3 = matrix[2][2]; 184 double d3 = matrix[2][3]; 185 186 double a4 = matrix[3][0]; 187 double b4 = matrix[3][1]; 188 double c4 = matrix[3][2]; 189 double d4 = matrix[3][3]; 190 191 // Row column labeling reversed since we transpose rows & columns 192 result[0][0] = determinant3x3(b2, b3, b4, c2, c3, c4, d2, d3, d4); 193 result[1][0] = - determinant3x3(a2, a3, a4, c2, c3, c4, d2, d3, d4); 194 result[2][0] = determinant3x3(a2, a3, a4, b2, b3, b4, d2, d3, d4); 195 result[3][0] = - determinant3x3(a2, a3, a4, b2, b3, b4, c2, c3, c4); 196 197 result[0][1] = - determinant3x3(b1, b3, b4, c1, c3, c4, d1, d3, d4); 198 result[1][1] = determinant3x3(a1, a3, a4, c1, c3, c4, d1, d3, d4); 199 result[2][1] = - determinant3x3(a1, a3, a4, b1, b3, b4, d1, d3, d4); 200 result[3][1] = determinant3x3(a1, a3, a4, b1, b3, b4, c1, c3, c4); 201 202 result[0][2] = determinant3x3(b1, b2, b4, c1, c2, c4, d1, d2, d4); 203 result[1][2] = - determinant3x3(a1, a2, a4, c1, c2, c4, d1, d2, d4); 204 result[2][2] = determinant3x3(a1, a2, a4, b1, b2, b4, d1, d2, d4); 205 result[3][2] = - determinant3x3(a1, a2, a4, b1, b2, b4, c1, c2, c4); 206 207 result[0][3] = - determinant3x3(b1, b2, b3, c1, c2, c3, d1, d2, d3); 208 result[1][3] = determinant3x3(a1, a2, a3, c1, c2, c3, d1, d2, d3); 209 result[2][3] = - determinant3x3(a1, a2, a3, b1, b2, b3, d1, d2, d3); 210 result[3][3] = determinant3x3(a1, a2, a3, b1, b2, b3, c1, c2, c3); 211 } 212 213 // Returns false if the matrix is not invertible 214 static bool inverse(const TransformationMatrix::Matrix4& matrix, TransformationMatrix::Matrix4& result) 215 { 216 // Calculate the adjoint matrix 217 adjoint(matrix, result); 218 219 // Calculate the 4x4 determinant 220 // If the determinant is zero, 221 // then the inverse matrix is not unique. 222 double det = determinant4x4(matrix); 223 224 if (fabs(det) < SMALL_NUMBER) 225 return false; 226 227 // Scale the adjoint matrix to get the inverse 228 229 for (int i = 0; i < 4; i++) 230 for (int j = 0; j < 4; j++) 231 result[i][j] = result[i][j] / det; 232 233 return true; 234 } 235 236 // End of code adapted from Matrix Inversion by Richard Carling 237 238 // Perform a decomposition on the passed matrix, return false if unsuccessful 239 // From Graphics Gems: unmatrix.c 240 241 // Transpose rotation portion of matrix a, return b 242 static void transposeMatrix4(const TransformationMatrix::Matrix4& a, TransformationMatrix::Matrix4& b) 243 { 244 for (int i = 0; i < 4; i++) 245 for (int j = 0; j < 4; j++) 246 b[i][j] = a[j][i]; 247 } 248 249 // Multiply a homogeneous point by a matrix and return the transformed point 250 static void v4MulPointByMatrix(const Vector4 p, const TransformationMatrix::Matrix4& m, Vector4 result) 251 { 252 result[0] = (p[0] * m[0][0]) + (p[1] * m[1][0]) + 253 (p[2] * m[2][0]) + (p[3] * m[3][0]); 254 result[1] = (p[0] * m[0][1]) + (p[1] * m[1][1]) + 255 (p[2] * m[2][1]) + (p[3] * m[3][1]); 256 result[2] = (p[0] * m[0][2]) + (p[1] * m[1][2]) + 257 (p[2] * m[2][2]) + (p[3] * m[3][2]); 258 result[3] = (p[0] * m[0][3]) + (p[1] * m[1][3]) + 259 (p[2] * m[2][3]) + (p[3] * m[3][3]); 260 } 261 262 static double v3Length(Vector3 a) 263 { 264 return sqrt((a[0] * a[0]) + (a[1] * a[1]) + (a[2] * a[2])); 265 } 266 267 static void v3Scale(Vector3 v, double desiredLength) 268 { 269 double len = v3Length(v); 270 if (len != 0) { 271 double l = desiredLength / len; 272 v[0] *= l; 273 v[1] *= l; 274 v[2] *= l; 275 } 276 } 277 278 static double v3Dot(const Vector3 a, const Vector3 b) 279 { 280 return (a[0] * b[0]) + (a[1] * b[1]) + (a[2] * b[2]); 281 } 282 283 // Make a linear combination of two vectors and return the result. 284 // result = (a * ascl) + (b * bscl) 285 static void v3Combine(const Vector3 a, const Vector3 b, Vector3 result, double ascl, double bscl) 286 { 287 result[0] = (ascl * a[0]) + (bscl * b[0]); 288 result[1] = (ascl * a[1]) + (bscl * b[1]); 289 result[2] = (ascl * a[2]) + (bscl * b[2]); 290 } 291 292 // Return the cross product result = a cross b */ 293 static void v3Cross(const Vector3 a, const Vector3 b, Vector3 result) 294 { 295 result[0] = (a[1] * b[2]) - (a[2] * b[1]); 296 result[1] = (a[2] * b[0]) - (a[0] * b[2]); 297 result[2] = (a[0] * b[1]) - (a[1] * b[0]); 298 } 299 300 static bool decompose(const TransformationMatrix::Matrix4& mat, TransformationMatrix::DecomposedType& result) 301 { 302 TransformationMatrix::Matrix4 localMatrix; 303 memcpy(localMatrix, mat, sizeof(TransformationMatrix::Matrix4)); 304 305 // Normalize the matrix. 306 if (localMatrix[3][3] == 0) 307 return false; 308 309 int i, j; 310 for (i = 0; i < 4; i++) 311 for (j = 0; j < 4; j++) 312 localMatrix[i][j] /= localMatrix[3][3]; 313 314 // perspectiveMatrix is used to solve for perspective, but it also provides 315 // an easy way to test for singularity of the upper 3x3 component. 316 TransformationMatrix::Matrix4 perspectiveMatrix; 317 memcpy(perspectiveMatrix, localMatrix, sizeof(TransformationMatrix::Matrix4)); 318 for (i = 0; i < 3; i++) 319 perspectiveMatrix[i][3] = 0; 320 perspectiveMatrix[3][3] = 1; 321 322 if (determinant4x4(perspectiveMatrix) == 0) 323 return false; 324 325 // First, isolate perspective. This is the messiest. 326 if (localMatrix[0][3] != 0 || localMatrix[1][3] != 0 || localMatrix[2][3] != 0) { 327 // rightHandSide is the right hand side of the equation. 328 Vector4 rightHandSide; 329 rightHandSide[0] = localMatrix[0][3]; 330 rightHandSide[1] = localMatrix[1][3]; 331 rightHandSide[2] = localMatrix[2][3]; 332 rightHandSide[3] = localMatrix[3][3]; 333 334 // Solve the equation by inverting perspectiveMatrix and multiplying 335 // rightHandSide by the inverse. (This is the easiest way, not 336 // necessarily the best.) 337 TransformationMatrix::Matrix4 inversePerspectiveMatrix, transposedInversePerspectiveMatrix; 338 inverse(perspectiveMatrix, inversePerspectiveMatrix); 339 transposeMatrix4(inversePerspectiveMatrix, transposedInversePerspectiveMatrix); 340 341 Vector4 perspectivePoint; 342 v4MulPointByMatrix(rightHandSide, transposedInversePerspectiveMatrix, perspectivePoint); 343 344 result.perspectiveX = perspectivePoint[0]; 345 result.perspectiveY = perspectivePoint[1]; 346 result.perspectiveZ = perspectivePoint[2]; 347 result.perspectiveW = perspectivePoint[3]; 348 349 // Clear the perspective partition 350 localMatrix[0][3] = localMatrix[1][3] = localMatrix[2][3] = 0; 351 localMatrix[3][3] = 1; 352 } else { 353 // No perspective. 354 result.perspectiveX = result.perspectiveY = result.perspectiveZ = 0; 355 result.perspectiveW = 1; 356 } 357 358 // Next take care of translation (easy). 359 result.translateX = localMatrix[3][0]; 360 localMatrix[3][0] = 0; 361 result.translateY = localMatrix[3][1]; 362 localMatrix[3][1] = 0; 363 result.translateZ = localMatrix[3][2]; 364 localMatrix[3][2] = 0; 365 366 // Vector4 type and functions need to be added to the common set. 367 Vector3 row[3], pdum3; 368 369 // Now get scale and shear. 370 for (i = 0; i < 3; i++) { 371 row[i][0] = localMatrix[i][0]; 372 row[i][1] = localMatrix[i][1]; 373 row[i][2] = localMatrix[i][2]; 374 } 375 376 // Compute X scale factor and normalize first row. 377 result.scaleX = v3Length(row[0]); 378 v3Scale(row[0], 1.0); 379 380 // Compute XY shear factor and make 2nd row orthogonal to 1st. 381 result.skewXY = v3Dot(row[0], row[1]); 382 v3Combine(row[1], row[0], row[1], 1.0, -result.skewXY); 383 384 // Now, compute Y scale and normalize 2nd row. 385 result.scaleY = v3Length(row[1]); 386 v3Scale(row[1], 1.0); 387 result.skewXY /= result.scaleY; 388 389 // Compute XZ and YZ shears, orthogonalize 3rd row. 390 result.skewXZ = v3Dot(row[0], row[2]); 391 v3Combine(row[2], row[0], row[2], 1.0, -result.skewXZ); 392 result.skewYZ = v3Dot(row[1], row[2]); 393 v3Combine(row[2], row[1], row[2], 1.0, -result.skewYZ); 394 395 // Next, get Z scale and normalize 3rd row. 396 result.scaleZ = v3Length(row[2]); 397 v3Scale(row[2], 1.0); 398 result.skewXZ /= result.scaleZ; 399 result.skewYZ /= result.scaleZ; 400 401 // At this point, the matrix (in rows[]) is orthonormal. 402 // Check for a coordinate system flip. If the determinant 403 // is -1, then negate the matrix and the scaling factors. 404 v3Cross(row[1], row[2], pdum3); 405 if (v3Dot(row[0], pdum3) < 0) { 406 407 result.scaleX *= -1; 408 result.scaleY *= -1; 409 result.scaleZ *= -1; 410 411 for (i = 0; i < 3; i++) { 412 row[i][0] *= -1; 413 row[i][1] *= -1; 414 row[i][2] *= -1; 415 } 416 } 417 418 // Now, get the rotations out, as described in the gem. 419 420 // FIXME - Add the ability to return either quaternions (which are 421 // easier to recompose with) or Euler angles (rx, ry, rz), which 422 // are easier for authors to deal with. The latter will only be useful 423 // when we fix https://bugs.webkit.org/show_bug.cgi?id=23799, so I 424 // will leave the Euler angle code here for now. 425 426 // ret.rotateY = asin(-row[0][2]); 427 // if (cos(ret.rotateY) != 0) { 428 // ret.rotateX = atan2(row[1][2], row[2][2]); 429 // ret.rotateZ = atan2(row[0][1], row[0][0]); 430 // } else { 431 // ret.rotateX = atan2(-row[2][0], row[1][1]); 432 // ret.rotateZ = 0; 433 // } 434 435 double s, t, x, y, z, w; 436 437 t = row[0][0] + row[1][1] + row[2][2] + 1.0; 438 439 if (t > 1e-4) { 440 s = 0.5 / sqrt(t); 441 w = 0.25 / s; 442 x = (row[2][1] - row[1][2]) * s; 443 y = (row[0][2] - row[2][0]) * s; 444 z = (row[1][0] - row[0][1]) * s; 445 } else if (row[0][0] > row[1][1] && row[0][0] > row[2][2]) { 446 s = sqrt (1.0 + row[0][0] - row[1][1] - row[2][2]) * 2.0; // S=4*qx 447 x = 0.25 * s; 448 y = (row[0][1] + row[1][0]) / s; 449 z = (row[0][2] + row[2][0]) / s; 450 w = (row[2][1] - row[1][2]) / s; 451 } else if (row[1][1] > row[2][2]) { 452 s = sqrt (1.0 + row[1][1] - row[0][0] - row[2][2]) * 2.0; // S=4*qy 453 x = (row[0][1] + row[1][0]) / s; 454 y = 0.25 * s; 455 z = (row[1][2] + row[2][1]) / s; 456 w = (row[0][2] - row[2][0]) / s; 457 } else { 458 s = sqrt(1.0 + row[2][2] - row[0][0] - row[1][1]) * 2.0; // S=4*qz 459 x = (row[0][2] + row[2][0]) / s; 460 y = (row[1][2] + row[2][1]) / s; 461 z = 0.25 * s; 462 w = (row[1][0] - row[0][1]) / s; 463 } 464 465 result.quaternionX = x; 466 result.quaternionY = y; 467 result.quaternionZ = z; 468 result.quaternionW = w; 469 470 return true; 471 } 472 473 // Perform a spherical linear interpolation between the two 474 // passed quaternions with 0 <= t <= 1 475 static void slerp(double qa[4], const double qb[4], double t) 476 { 477 double ax, ay, az, aw; 478 double bx, by, bz, bw; 479 double cx, cy, cz, cw; 480 double angle; 481 double th, invth, scale, invscale; 482 483 ax = qa[0]; ay = qa[1]; az = qa[2]; aw = qa[3]; 484 bx = qb[0]; by = qb[1]; bz = qb[2]; bw = qb[3]; 485 486 angle = ax * bx + ay * by + az * bz + aw * bw; 487 488 if (angle < 0.0) { 489 ax = -ax; ay = -ay; 490 az = -az; aw = -aw; 491 angle = -angle; 492 } 493 494 if (angle + 1.0 > .05) { 495 if (1.0 - angle >= .05) { 496 th = acos (angle); 497 invth = 1.0 / sin (th); 498 scale = sin (th * (1.0 - t)) * invth; 499 invscale = sin (th * t) * invth; 500 } else { 501 scale = 1.0 - t; 502 invscale = t; 503 } 504 } else { 505 bx = -ay; 506 by = ax; 507 bz = -aw; 508 bw = az; 509 scale = sin(piDouble * (.5 - t)); 510 invscale = sin (piDouble * t); 511 } 512 513 cx = ax * scale + bx * invscale; 514 cy = ay * scale + by * invscale; 515 cz = az * scale + bz * invscale; 516 cw = aw * scale + bw * invscale; 517 518 qa[0] = cx; qa[1] = cy; qa[2] = cz; qa[3] = cw; 519 } 520 521 // End of Supporting Math Functions 522 523 TransformationMatrix::TransformationMatrix(const AffineTransform& t) 524 { 525 setMatrix(t.a(), t.b(), t.c(), t.d(), t.e(), t.f()); 526 } 527 528 TransformationMatrix& TransformationMatrix::scale(double s) 529 { 530 return scaleNonUniform(s, s); 531 } 532 533 TransformationMatrix& TransformationMatrix::rotateFromVector(double x, double y) 534 { 535 return rotate(rad2deg(atan2(y, x))); 536 } 537 538 TransformationMatrix& TransformationMatrix::flipX() 539 { 540 return scaleNonUniform(-1.0, 1.0); 541 } 542 543 TransformationMatrix& TransformationMatrix::flipY() 544 { 545 return scaleNonUniform(1.0, -1.0); 546 } 547 548 FloatPoint TransformationMatrix::projectPoint(const FloatPoint& p, bool* clamped) const 549 { 550 // This is basically raytracing. We have a point in the destination 551 // plane with z=0, and we cast a ray parallel to the z-axis from that 552 // point to find the z-position at which it intersects the z=0 plane 553 // with the transform applied. Once we have that point we apply the 554 // inverse transform to find the corresponding point in the source 555 // space. 556 // 557 // Given a plane with normal Pn, and a ray starting at point R0 and 558 // with direction defined by the vector Rd, we can find the 559 // intersection point as a distance d from R0 in units of Rd by: 560 // 561 // d = -dot (Pn', R0) / dot (Pn', Rd) 562 if (clamped) 563 *clamped = false; 564 565 if (m33() == 0) { 566 // In this case, the projection plane is parallel to the ray we are trying to 567 // trace, and there is no well-defined value for the projection. 568 return FloatPoint(); 569 } 570 571 double x = p.x(); 572 double y = p.y(); 573 double z = -(m13() * x + m23() * y + m43()) / m33(); 574 575 // FIXME: use multVecMatrix() 576 double outX = x * m11() + y * m21() + z * m31() + m41(); 577 double outY = x * m12() + y * m22() + z * m32() + m42(); 578 579 double w = x * m14() + y * m24() + z * m34() + m44(); 580 if (w <= 0) { 581 // Using int max causes overflow when other code uses the projected point. To 582 // represent infinity yet reduce the risk of overflow, we use a large but 583 // not-too-large number here when clamping. 584 const int largeNumber = 100000000 / kFixedPointDenominator; 585 outX = copysign(largeNumber, outX); 586 outY = copysign(largeNumber, outY); 587 if (clamped) 588 *clamped = true; 589 } else if (w != 1) { 590 outX /= w; 591 outY /= w; 592 } 593 594 return FloatPoint(static_cast<float>(outX), static_cast<float>(outY)); 595 } 596 597 FloatQuad TransformationMatrix::projectQuad(const FloatQuad& q, bool* clamped) const 598 { 599 FloatQuad projectedQuad; 600 601 bool clamped1 = false; 602 bool clamped2 = false; 603 bool clamped3 = false; 604 bool clamped4 = false; 605 606 projectedQuad.setP1(projectPoint(q.p1(), &clamped1)); 607 projectedQuad.setP2(projectPoint(q.p2(), &clamped2)); 608 projectedQuad.setP3(projectPoint(q.p3(), &clamped3)); 609 projectedQuad.setP4(projectPoint(q.p4(), &clamped4)); 610 611 if (clamped) 612 *clamped = clamped1 || clamped2 || clamped3 || clamped4; 613 614 // If all points on the quad had w < 0, then the entire quad would not be visible to the projected surface. 615 bool everythingWasClipped = clamped1 && clamped2 && clamped3 && clamped4; 616 if (everythingWasClipped) 617 return FloatQuad(); 618 619 return projectedQuad; 620 } 621 622 static float clampEdgeValue(float f) 623 { 624 ASSERT(!std::isnan(f)); 625 return min<float>(max<float>(f, (-LayoutUnit::max() / 2).toFloat()), (LayoutUnit::max() / 2).toFloat()); 626 } 627 628 LayoutRect TransformationMatrix::clampedBoundsOfProjectedQuad(const FloatQuad& q) const 629 { 630 FloatRect mappedQuadBounds = projectQuad(q).boundingBox(); 631 632 float left = clampEdgeValue(floorf(mappedQuadBounds.x())); 633 float top = clampEdgeValue(floorf(mappedQuadBounds.y())); 634 635 float right; 636 if (std::isinf(mappedQuadBounds.x()) && std::isinf(mappedQuadBounds.width())) 637 right = (LayoutUnit::max() / 2).toFloat(); 638 else 639 right = clampEdgeValue(ceilf(mappedQuadBounds.maxX())); 640 641 float bottom; 642 if (std::isinf(mappedQuadBounds.y()) && std::isinf(mappedQuadBounds.height())) 643 bottom = (LayoutUnit::max() / 2).toFloat(); 644 else 645 bottom = clampEdgeValue(ceilf(mappedQuadBounds.maxY())); 646 647 return LayoutRect(LayoutUnit::clamp(left), LayoutUnit::clamp(top), LayoutUnit::clamp(right - left), LayoutUnit::clamp(bottom - top)); 648 } 649 650 void TransformationMatrix::transformBox(FloatBox& box) const 651 { 652 FloatBox bounds; 653 bool firstPoint = true; 654 for (size_t i = 0; i < 2; ++i) { 655 for (size_t j = 0; j < 2; ++j) { 656 for (size_t k = 0; k < 2; ++k) { 657 FloatPoint3D point(box.x(), box.y(), box.z()); 658 point += FloatPoint3D(i * box.width(), j * box.height(), k * box.depth()); 659 point = mapPoint(point); 660 if (firstPoint) { 661 bounds.setOrigin(point); 662 firstPoint = false; 663 } else { 664 bounds.expandTo(point); 665 } 666 } 667 } 668 } 669 box = bounds; 670 } 671 672 FloatPoint TransformationMatrix::mapPoint(const FloatPoint& p) const 673 { 674 if (isIdentityOrTranslation()) 675 return FloatPoint(p.x() + static_cast<float>(m_matrix[3][0]), p.y() + static_cast<float>(m_matrix[3][1])); 676 677 return internalMapPoint(p); 678 } 679 680 FloatPoint3D TransformationMatrix::mapPoint(const FloatPoint3D& p) const 681 { 682 if (isIdentityOrTranslation()) 683 return FloatPoint3D(p.x() + static_cast<float>(m_matrix[3][0]), 684 p.y() + static_cast<float>(m_matrix[3][1]), 685 p.z() + static_cast<float>(m_matrix[3][2])); 686 687 return internalMapPoint(p); 688 } 689 690 IntRect TransformationMatrix::mapRect(const IntRect &rect) const 691 { 692 return enclosingIntRect(mapRect(FloatRect(rect))); 693 } 694 695 LayoutRect TransformationMatrix::mapRect(const LayoutRect& r) const 696 { 697 return enclosingLayoutRect(mapRect(FloatRect(r))); 698 } 699 700 FloatRect TransformationMatrix::mapRect(const FloatRect& r) const 701 { 702 if (isIdentityOrTranslation()) { 703 FloatRect mappedRect(r); 704 mappedRect.move(static_cast<float>(m_matrix[3][0]), static_cast<float>(m_matrix[3][1])); 705 return mappedRect; 706 } 707 708 FloatQuad result; 709 710 float maxX = r.maxX(); 711 float maxY = r.maxY(); 712 result.setP1(internalMapPoint(FloatPoint(r.x(), r.y()))); 713 result.setP2(internalMapPoint(FloatPoint(maxX, r.y()))); 714 result.setP3(internalMapPoint(FloatPoint(maxX, maxY))); 715 result.setP4(internalMapPoint(FloatPoint(r.x(), maxY))); 716 717 return result.boundingBox(); 718 } 719 720 FloatQuad TransformationMatrix::mapQuad(const FloatQuad& q) const 721 { 722 if (isIdentityOrTranslation()) { 723 FloatQuad mappedQuad(q); 724 mappedQuad.move(static_cast<float>(m_matrix[3][0]), static_cast<float>(m_matrix[3][1])); 725 return mappedQuad; 726 } 727 728 FloatQuad result; 729 result.setP1(internalMapPoint(q.p1())); 730 result.setP2(internalMapPoint(q.p2())); 731 result.setP3(internalMapPoint(q.p3())); 732 result.setP4(internalMapPoint(q.p4())); 733 return result; 734 } 735 736 TransformationMatrix& TransformationMatrix::scaleNonUniform(double sx, double sy) 737 { 738 m_matrix[0][0] *= sx; 739 m_matrix[0][1] *= sx; 740 m_matrix[0][2] *= sx; 741 m_matrix[0][3] *= sx; 742 743 m_matrix[1][0] *= sy; 744 m_matrix[1][1] *= sy; 745 m_matrix[1][2] *= sy; 746 m_matrix[1][3] *= sy; 747 return *this; 748 } 749 750 TransformationMatrix& TransformationMatrix::scale3d(double sx, double sy, double sz) 751 { 752 scaleNonUniform(sx, sy); 753 754 m_matrix[2][0] *= sz; 755 m_matrix[2][1] *= sz; 756 m_matrix[2][2] *= sz; 757 m_matrix[2][3] *= sz; 758 return *this; 759 } 760 761 TransformationMatrix& TransformationMatrix::rotate3d(double x, double y, double z, double angle) 762 { 763 // Normalize the axis of rotation 764 double length = sqrt(x * x + y * y + z * z); 765 if (length == 0) { 766 // A direction vector that cannot be normalized, such as [0, 0, 0], will cause the rotation to not be applied. 767 return *this; 768 } else if (length != 1) { 769 x /= length; 770 y /= length; 771 z /= length; 772 } 773 774 // Angles are in degrees. Switch to radians. 775 angle = deg2rad(angle); 776 777 double sinTheta = sin(angle); 778 double cosTheta = cos(angle); 779 780 TransformationMatrix mat; 781 782 // Optimize cases where the axis is along a major axis 783 if (x == 1.0 && y == 0.0 && z == 0.0) { 784 mat.m_matrix[0][0] = 1.0; 785 mat.m_matrix[0][1] = 0.0; 786 mat.m_matrix[0][2] = 0.0; 787 mat.m_matrix[1][0] = 0.0; 788 mat.m_matrix[1][1] = cosTheta; 789 mat.m_matrix[1][2] = sinTheta; 790 mat.m_matrix[2][0] = 0.0; 791 mat.m_matrix[2][1] = -sinTheta; 792 mat.m_matrix[2][2] = cosTheta; 793 mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; 794 mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; 795 mat.m_matrix[3][3] = 1.0; 796 } else if (x == 0.0 && y == 1.0 && z == 0.0) { 797 mat.m_matrix[0][0] = cosTheta; 798 mat.m_matrix[0][1] = 0.0; 799 mat.m_matrix[0][2] = -sinTheta; 800 mat.m_matrix[1][0] = 0.0; 801 mat.m_matrix[1][1] = 1.0; 802 mat.m_matrix[1][2] = 0.0; 803 mat.m_matrix[2][0] = sinTheta; 804 mat.m_matrix[2][1] = 0.0; 805 mat.m_matrix[2][2] = cosTheta; 806 mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; 807 mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; 808 mat.m_matrix[3][3] = 1.0; 809 } else if (x == 0.0 && y == 0.0 && z == 1.0) { 810 mat.m_matrix[0][0] = cosTheta; 811 mat.m_matrix[0][1] = sinTheta; 812 mat.m_matrix[0][2] = 0.0; 813 mat.m_matrix[1][0] = -sinTheta; 814 mat.m_matrix[1][1] = cosTheta; 815 mat.m_matrix[1][2] = 0.0; 816 mat.m_matrix[2][0] = 0.0; 817 mat.m_matrix[2][1] = 0.0; 818 mat.m_matrix[2][2] = 1.0; 819 mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; 820 mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; 821 mat.m_matrix[3][3] = 1.0; 822 } else { 823 // This case is the rotation about an arbitrary unit vector. 824 // 825 // Formula is adapted from Wikipedia article on Rotation matrix, 826 // http://en.wikipedia.org/wiki/Rotation_matrix#Rotation_matrix_from_axis_and_angle 827 // 828 // An alternate resource with the same matrix: http://www.fastgraph.com/makegames/3drotation/ 829 // 830 double oneMinusCosTheta = 1 - cosTheta; 831 mat.m_matrix[0][0] = cosTheta + x * x * oneMinusCosTheta; 832 mat.m_matrix[0][1] = y * x * oneMinusCosTheta + z * sinTheta; 833 mat.m_matrix[0][2] = z * x * oneMinusCosTheta - y * sinTheta; 834 mat.m_matrix[1][0] = x * y * oneMinusCosTheta - z * sinTheta; 835 mat.m_matrix[1][1] = cosTheta + y * y * oneMinusCosTheta; 836 mat.m_matrix[1][2] = z * y * oneMinusCosTheta + x * sinTheta; 837 mat.m_matrix[2][0] = x * z * oneMinusCosTheta + y * sinTheta; 838 mat.m_matrix[2][1] = y * z * oneMinusCosTheta - x * sinTheta; 839 mat.m_matrix[2][2] = cosTheta + z * z * oneMinusCosTheta; 840 mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; 841 mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; 842 mat.m_matrix[3][3] = 1.0; 843 } 844 multiply(mat); 845 return *this; 846 } 847 848 TransformationMatrix& TransformationMatrix::rotate3d(double rx, double ry, double rz) 849 { 850 // Angles are in degrees. Switch to radians. 851 rx = deg2rad(rx); 852 ry = deg2rad(ry); 853 rz = deg2rad(rz); 854 855 TransformationMatrix mat; 856 857 double sinTheta = sin(rz); 858 double cosTheta = cos(rz); 859 860 mat.m_matrix[0][0] = cosTheta; 861 mat.m_matrix[0][1] = sinTheta; 862 mat.m_matrix[0][2] = 0.0; 863 mat.m_matrix[1][0] = -sinTheta; 864 mat.m_matrix[1][1] = cosTheta; 865 mat.m_matrix[1][2] = 0.0; 866 mat.m_matrix[2][0] = 0.0; 867 mat.m_matrix[2][1] = 0.0; 868 mat.m_matrix[2][2] = 1.0; 869 mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; 870 mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; 871 mat.m_matrix[3][3] = 1.0; 872 873 TransformationMatrix rmat(mat); 874 875 sinTheta = sin(ry); 876 cosTheta = cos(ry); 877 878 mat.m_matrix[0][0] = cosTheta; 879 mat.m_matrix[0][1] = 0.0; 880 mat.m_matrix[0][2] = -sinTheta; 881 mat.m_matrix[1][0] = 0.0; 882 mat.m_matrix[1][1] = 1.0; 883 mat.m_matrix[1][2] = 0.0; 884 mat.m_matrix[2][0] = sinTheta; 885 mat.m_matrix[2][1] = 0.0; 886 mat.m_matrix[2][2] = cosTheta; 887 mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; 888 mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; 889 mat.m_matrix[3][3] = 1.0; 890 891 rmat.multiply(mat); 892 893 sinTheta = sin(rx); 894 cosTheta = cos(rx); 895 896 mat.m_matrix[0][0] = 1.0; 897 mat.m_matrix[0][1] = 0.0; 898 mat.m_matrix[0][2] = 0.0; 899 mat.m_matrix[1][0] = 0.0; 900 mat.m_matrix[1][1] = cosTheta; 901 mat.m_matrix[1][2] = sinTheta; 902 mat.m_matrix[2][0] = 0.0; 903 mat.m_matrix[2][1] = -sinTheta; 904 mat.m_matrix[2][2] = cosTheta; 905 mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0; 906 mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0; 907 mat.m_matrix[3][3] = 1.0; 908 909 rmat.multiply(mat); 910 911 multiply(rmat); 912 return *this; 913 } 914 915 TransformationMatrix& TransformationMatrix::translate(double tx, double ty) 916 { 917 m_matrix[3][0] += tx * m_matrix[0][0] + ty * m_matrix[1][0]; 918 m_matrix[3][1] += tx * m_matrix[0][1] + ty * m_matrix[1][1]; 919 m_matrix[3][2] += tx * m_matrix[0][2] + ty * m_matrix[1][2]; 920 m_matrix[3][3] += tx * m_matrix[0][3] + ty * m_matrix[1][3]; 921 return *this; 922 } 923 924 TransformationMatrix& TransformationMatrix::translate3d(double tx, double ty, double tz) 925 { 926 m_matrix[3][0] += tx * m_matrix[0][0] + ty * m_matrix[1][0] + tz * m_matrix[2][0]; 927 m_matrix[3][1] += tx * m_matrix[0][1] + ty * m_matrix[1][1] + tz * m_matrix[2][1]; 928 m_matrix[3][2] += tx * m_matrix[0][2] + ty * m_matrix[1][2] + tz * m_matrix[2][2]; 929 m_matrix[3][3] += tx * m_matrix[0][3] + ty * m_matrix[1][3] + tz * m_matrix[2][3]; 930 return *this; 931 } 932 933 TransformationMatrix& TransformationMatrix::translateRight(double tx, double ty) 934 { 935 if (tx != 0) { 936 m_matrix[0][0] += m_matrix[0][3] * tx; 937 m_matrix[1][0] += m_matrix[1][3] * tx; 938 m_matrix[2][0] += m_matrix[2][3] * tx; 939 m_matrix[3][0] += m_matrix[3][3] * tx; 940 } 941 942 if (ty != 0) { 943 m_matrix[0][1] += m_matrix[0][3] * ty; 944 m_matrix[1][1] += m_matrix[1][3] * ty; 945 m_matrix[2][1] += m_matrix[2][3] * ty; 946 m_matrix[3][1] += m_matrix[3][3] * ty; 947 } 948 949 return *this; 950 } 951 952 TransformationMatrix& TransformationMatrix::translateRight3d(double tx, double ty, double tz) 953 { 954 translateRight(tx, ty); 955 if (tz != 0) { 956 m_matrix[0][2] += m_matrix[0][3] * tz; 957 m_matrix[1][2] += m_matrix[1][3] * tz; 958 m_matrix[2][2] += m_matrix[2][3] * tz; 959 m_matrix[3][2] += m_matrix[3][3] * tz; 960 } 961 962 return *this; 963 } 964 965 TransformationMatrix& TransformationMatrix::skew(double sx, double sy) 966 { 967 // angles are in degrees. Switch to radians 968 sx = deg2rad(sx); 969 sy = deg2rad(sy); 970 971 TransformationMatrix mat; 972 mat.m_matrix[0][1] = tan(sy); // note that the y shear goes in the first row 973 mat.m_matrix[1][0] = tan(sx); // and the x shear in the second row 974 975 multiply(mat); 976 return *this; 977 } 978 979 TransformationMatrix& TransformationMatrix::applyPerspective(double p) 980 { 981 TransformationMatrix mat; 982 if (p != 0) 983 mat.m_matrix[2][3] = -1/p; 984 985 multiply(mat); 986 return *this; 987 } 988 989 TransformationMatrix TransformationMatrix::rectToRect(const FloatRect& from, const FloatRect& to) 990 { 991 ASSERT(!from.isEmpty()); 992 return TransformationMatrix(to.width() / from.width(), 993 0, 0, 994 to.height() / from.height(), 995 to.x() - from.x(), 996 to.y() - from.y()); 997 } 998 999 // this = mat * this. 1000 TransformationMatrix& TransformationMatrix::multiply(const TransformationMatrix& mat) 1001 { 1002 #if CPU(APPLE_ARMV7S) 1003 double* leftMatrix = &(m_matrix[0][0]); 1004 const double* rightMatrix = &(mat.m_matrix[0][0]); 1005 asm volatile (// First row of leftMatrix. 1006 "mov r3, %[leftMatrix]\n\t" 1007 "vld1.64 { d16-d19 }, [%[leftMatrix], :128]!\n\t" 1008 "vld1.64 { d0-d3}, [%[rightMatrix], :128]!\n\t" 1009 "vmul.f64 d4, d0, d16\n\t" 1010 "vld1.64 { d20-d23 }, [%[leftMatrix], :128]!\n\t" 1011 "vmla.f64 d4, d1, d20\n\t" 1012 "vld1.64 { d24-d27 }, [%[leftMatrix], :128]!\n\t" 1013 "vmla.f64 d4, d2, d24\n\t" 1014 "vld1.64 { d28-d31 }, [%[leftMatrix], :128]!\n\t" 1015 "vmla.f64 d4, d3, d28\n\t" 1016 1017 "vmul.f64 d5, d0, d17\n\t" 1018 "vmla.f64 d5, d1, d21\n\t" 1019 "vmla.f64 d5, d2, d25\n\t" 1020 "vmla.f64 d5, d3, d29\n\t" 1021 1022 "vmul.f64 d6, d0, d18\n\t" 1023 "vmla.f64 d6, d1, d22\n\t" 1024 "vmla.f64 d6, d2, d26\n\t" 1025 "vmla.f64 d6, d3, d30\n\t" 1026 1027 "vmul.f64 d7, d0, d19\n\t" 1028 "vmla.f64 d7, d1, d23\n\t" 1029 "vmla.f64 d7, d2, d27\n\t" 1030 "vmla.f64 d7, d3, d31\n\t" 1031 "vld1.64 { d0-d3}, [%[rightMatrix], :128]!\n\t" 1032 "vst1.64 { d4-d7 }, [r3, :128]!\n\t" 1033 1034 // Second row of leftMatrix. 1035 "vmul.f64 d4, d0, d16\n\t" 1036 "vmla.f64 d4, d1, d20\n\t" 1037 "vmla.f64 d4, d2, d24\n\t" 1038 "vmla.f64 d4, d3, d28\n\t" 1039 1040 "vmul.f64 d5, d0, d17\n\t" 1041 "vmla.f64 d5, d1, d21\n\t" 1042 "vmla.f64 d5, d2, d25\n\t" 1043 "vmla.f64 d5, d3, d29\n\t" 1044 1045 "vmul.f64 d6, d0, d18\n\t" 1046 "vmla.f64 d6, d1, d22\n\t" 1047 "vmla.f64 d6, d2, d26\n\t" 1048 "vmla.f64 d6, d3, d30\n\t" 1049 1050 "vmul.f64 d7, d0, d19\n\t" 1051 "vmla.f64 d7, d1, d23\n\t" 1052 "vmla.f64 d7, d2, d27\n\t" 1053 "vmla.f64 d7, d3, d31\n\t" 1054 "vld1.64 { d0-d3}, [%[rightMatrix], :128]!\n\t" 1055 "vst1.64 { d4-d7 }, [r3, :128]!\n\t" 1056 1057 // Third row of leftMatrix. 1058 "vmul.f64 d4, d0, d16\n\t" 1059 "vmla.f64 d4, d1, d20\n\t" 1060 "vmla.f64 d4, d2, d24\n\t" 1061 "vmla.f64 d4, d3, d28\n\t" 1062 1063 "vmul.f64 d5, d0, d17\n\t" 1064 "vmla.f64 d5, d1, d21\n\t" 1065 "vmla.f64 d5, d2, d25\n\t" 1066 "vmla.f64 d5, d3, d29\n\t" 1067 1068 "vmul.f64 d6, d0, d18\n\t" 1069 "vmla.f64 d6, d1, d22\n\t" 1070 "vmla.f64 d6, d2, d26\n\t" 1071 "vmla.f64 d6, d3, d30\n\t" 1072 1073 "vmul.f64 d7, d0, d19\n\t" 1074 "vmla.f64 d7, d1, d23\n\t" 1075 "vmla.f64 d7, d2, d27\n\t" 1076 "vmla.f64 d7, d3, d31\n\t" 1077 "vld1.64 { d0-d3}, [%[rightMatrix], :128]\n\t" 1078 "vst1.64 { d4-d7 }, [r3, :128]!\n\t" 1079 1080 // Fourth and last row of leftMatrix. 1081 "vmul.f64 d4, d0, d16\n\t" 1082 "vmla.f64 d4, d1, d20\n\t" 1083 "vmla.f64 d4, d2, d24\n\t" 1084 "vmla.f64 d4, d3, d28\n\t" 1085 1086 "vmul.f64 d5, d0, d17\n\t" 1087 "vmla.f64 d5, d1, d21\n\t" 1088 "vmla.f64 d5, d2, d25\n\t" 1089 "vmla.f64 d5, d3, d29\n\t" 1090 1091 "vmul.f64 d6, d0, d18\n\t" 1092 "vmla.f64 d6, d1, d22\n\t" 1093 "vmla.f64 d6, d2, d26\n\t" 1094 "vmla.f64 d6, d3, d30\n\t" 1095 1096 "vmul.f64 d7, d0, d19\n\t" 1097 "vmla.f64 d7, d1, d23\n\t" 1098 "vmla.f64 d7, d2, d27\n\t" 1099 "vmla.f64 d7, d3, d31\n\t" 1100 "vst1.64 { d4-d7 }, [r3, :128]\n\t" 1101 : [leftMatrix]"+r"(leftMatrix), [rightMatrix]"+r"(rightMatrix) 1102 : 1103 : "memory", "r3", "d0", "d1", "d2", "d3", "d4", "d5", "d6", "d7", "d16", "d17", "d18", "d19", "d20", "d21", "d22", "d23", "d24", "d25", "d26", "d27", "d28", "d29", "d30", "d31"); 1104 #elif defined(TRANSFORMATION_MATRIX_USE_X86_64_SSE2) 1105 // x86_64 has 16 XMM registers which is enough to do the multiplication fully in registers. 1106 __m128d matrixBlockA = _mm_load_pd(&(m_matrix[0][0])); 1107 __m128d matrixBlockC = _mm_load_pd(&(m_matrix[1][0])); 1108 __m128d matrixBlockE = _mm_load_pd(&(m_matrix[2][0])); 1109 __m128d matrixBlockG = _mm_load_pd(&(m_matrix[3][0])); 1110 1111 // First row. 1112 __m128d otherMatrixFirstParam = _mm_set1_pd(mat.m_matrix[0][0]); 1113 __m128d otherMatrixSecondParam = _mm_set1_pd(mat.m_matrix[0][1]); 1114 __m128d otherMatrixThirdParam = _mm_set1_pd(mat.m_matrix[0][2]); 1115 __m128d otherMatrixFourthParam = _mm_set1_pd(mat.m_matrix[0][3]); 1116 1117 // output00 and output01. 1118 __m128d accumulator = _mm_mul_pd(matrixBlockA, otherMatrixFirstParam); 1119 __m128d temp1 = _mm_mul_pd(matrixBlockC, otherMatrixSecondParam); 1120 __m128d temp2 = _mm_mul_pd(matrixBlockE, otherMatrixThirdParam); 1121 __m128d temp3 = _mm_mul_pd(matrixBlockG, otherMatrixFourthParam); 1122 1123 __m128d matrixBlockB = _mm_load_pd(&(m_matrix[0][2])); 1124 __m128d matrixBlockD = _mm_load_pd(&(m_matrix[1][2])); 1125 __m128d matrixBlockF = _mm_load_pd(&(m_matrix[2][2])); 1126 __m128d matrixBlockH = _mm_load_pd(&(m_matrix[3][2])); 1127 1128 accumulator = _mm_add_pd(accumulator, temp1); 1129 accumulator = _mm_add_pd(accumulator, temp2); 1130 accumulator = _mm_add_pd(accumulator, temp3); 1131 _mm_store_pd(&m_matrix[0][0], accumulator); 1132 1133 // output02 and output03. 1134 accumulator = _mm_mul_pd(matrixBlockB, otherMatrixFirstParam); 1135 temp1 = _mm_mul_pd(matrixBlockD, otherMatrixSecondParam); 1136 temp2 = _mm_mul_pd(matrixBlockF, otherMatrixThirdParam); 1137 temp3 = _mm_mul_pd(matrixBlockH, otherMatrixFourthParam); 1138 1139 accumulator = _mm_add_pd(accumulator, temp1); 1140 accumulator = _mm_add_pd(accumulator, temp2); 1141 accumulator = _mm_add_pd(accumulator, temp3); 1142 _mm_store_pd(&m_matrix[0][2], accumulator); 1143 1144 // Second row. 1145 otherMatrixFirstParam = _mm_set1_pd(mat.m_matrix[1][0]); 1146 otherMatrixSecondParam = _mm_set1_pd(mat.m_matrix[1][1]); 1147 otherMatrixThirdParam = _mm_set1_pd(mat.m_matrix[1][2]); 1148 otherMatrixFourthParam = _mm_set1_pd(mat.m_matrix[1][3]); 1149 1150 // output10 and output11. 1151 accumulator = _mm_mul_pd(matrixBlockA, otherMatrixFirstParam); 1152 temp1 = _mm_mul_pd(matrixBlockC, otherMatrixSecondParam); 1153 temp2 = _mm_mul_pd(matrixBlockE, otherMatrixThirdParam); 1154 temp3 = _mm_mul_pd(matrixBlockG, otherMatrixFourthParam); 1155 1156 accumulator = _mm_add_pd(accumulator, temp1); 1157 accumulator = _mm_add_pd(accumulator, temp2); 1158 accumulator = _mm_add_pd(accumulator, temp3); 1159 _mm_store_pd(&m_matrix[1][0], accumulator); 1160 1161 // output12 and output13. 1162 accumulator = _mm_mul_pd(matrixBlockB, otherMatrixFirstParam); 1163 temp1 = _mm_mul_pd(matrixBlockD, otherMatrixSecondParam); 1164 temp2 = _mm_mul_pd(matrixBlockF, otherMatrixThirdParam); 1165 temp3 = _mm_mul_pd(matrixBlockH, otherMatrixFourthParam); 1166 1167 accumulator = _mm_add_pd(accumulator, temp1); 1168 accumulator = _mm_add_pd(accumulator, temp2); 1169 accumulator = _mm_add_pd(accumulator, temp3); 1170 _mm_store_pd(&m_matrix[1][2], accumulator); 1171 1172 // Third row. 1173 otherMatrixFirstParam = _mm_set1_pd(mat.m_matrix[2][0]); 1174 otherMatrixSecondParam = _mm_set1_pd(mat.m_matrix[2][1]); 1175 otherMatrixThirdParam = _mm_set1_pd(mat.m_matrix[2][2]); 1176 otherMatrixFourthParam = _mm_set1_pd(mat.m_matrix[2][3]); 1177 1178 // output20 and output21. 1179 accumulator = _mm_mul_pd(matrixBlockA, otherMatrixFirstParam); 1180 temp1 = _mm_mul_pd(matrixBlockC, otherMatrixSecondParam); 1181 temp2 = _mm_mul_pd(matrixBlockE, otherMatrixThirdParam); 1182 temp3 = _mm_mul_pd(matrixBlockG, otherMatrixFourthParam); 1183 1184 accumulator = _mm_add_pd(accumulator, temp1); 1185 accumulator = _mm_add_pd(accumulator, temp2); 1186 accumulator = _mm_add_pd(accumulator, temp3); 1187 _mm_store_pd(&m_matrix[2][0], accumulator); 1188 1189 // output22 and output23. 1190 accumulator = _mm_mul_pd(matrixBlockB, otherMatrixFirstParam); 1191 temp1 = _mm_mul_pd(matrixBlockD, otherMatrixSecondParam); 1192 temp2 = _mm_mul_pd(matrixBlockF, otherMatrixThirdParam); 1193 temp3 = _mm_mul_pd(matrixBlockH, otherMatrixFourthParam); 1194 1195 accumulator = _mm_add_pd(accumulator, temp1); 1196 accumulator = _mm_add_pd(accumulator, temp2); 1197 accumulator = _mm_add_pd(accumulator, temp3); 1198 _mm_store_pd(&m_matrix[2][2], accumulator); 1199 1200 // Fourth row. 1201 otherMatrixFirstParam = _mm_set1_pd(mat.m_matrix[3][0]); 1202 otherMatrixSecondParam = _mm_set1_pd(mat.m_matrix[3][1]); 1203 otherMatrixThirdParam = _mm_set1_pd(mat.m_matrix[3][2]); 1204 otherMatrixFourthParam = _mm_set1_pd(mat.m_matrix[3][3]); 1205 1206 // output30 and output31. 1207 accumulator = _mm_mul_pd(matrixBlockA, otherMatrixFirstParam); 1208 temp1 = _mm_mul_pd(matrixBlockC, otherMatrixSecondParam); 1209 temp2 = _mm_mul_pd(matrixBlockE, otherMatrixThirdParam); 1210 temp3 = _mm_mul_pd(matrixBlockG, otherMatrixFourthParam); 1211 1212 accumulator = _mm_add_pd(accumulator, temp1); 1213 accumulator = _mm_add_pd(accumulator, temp2); 1214 accumulator = _mm_add_pd(accumulator, temp3); 1215 _mm_store_pd(&m_matrix[3][0], accumulator); 1216 1217 // output32 and output33. 1218 accumulator = _mm_mul_pd(matrixBlockB, otherMatrixFirstParam); 1219 temp1 = _mm_mul_pd(matrixBlockD, otherMatrixSecondParam); 1220 temp2 = _mm_mul_pd(matrixBlockF, otherMatrixThirdParam); 1221 temp3 = _mm_mul_pd(matrixBlockH, otherMatrixFourthParam); 1222 1223 accumulator = _mm_add_pd(accumulator, temp1); 1224 accumulator = _mm_add_pd(accumulator, temp2); 1225 accumulator = _mm_add_pd(accumulator, temp3); 1226 _mm_store_pd(&m_matrix[3][2], accumulator); 1227 #else 1228 Matrix4 tmp; 1229 1230 tmp[0][0] = (mat.m_matrix[0][0] * m_matrix[0][0] + mat.m_matrix[0][1] * m_matrix[1][0] 1231 + mat.m_matrix[0][2] * m_matrix[2][0] + mat.m_matrix[0][3] * m_matrix[3][0]); 1232 tmp[0][1] = (mat.m_matrix[0][0] * m_matrix[0][1] + mat.m_matrix[0][1] * m_matrix[1][1] 1233 + mat.m_matrix[0][2] * m_matrix[2][1] + mat.m_matrix[0][3] * m_matrix[3][1]); 1234 tmp[0][2] = (mat.m_matrix[0][0] * m_matrix[0][2] + mat.m_matrix[0][1] * m_matrix[1][2] 1235 + mat.m_matrix[0][2] * m_matrix[2][2] + mat.m_matrix[0][3] * m_matrix[3][2]); 1236 tmp[0][3] = (mat.m_matrix[0][0] * m_matrix[0][3] + mat.m_matrix[0][1] * m_matrix[1][3] 1237 + mat.m_matrix[0][2] * m_matrix[2][3] + mat.m_matrix[0][3] * m_matrix[3][3]); 1238 1239 tmp[1][0] = (mat.m_matrix[1][0] * m_matrix[0][0] + mat.m_matrix[1][1] * m_matrix[1][0] 1240 + mat.m_matrix[1][2] * m_matrix[2][0] + mat.m_matrix[1][3] * m_matrix[3][0]); 1241 tmp[1][1] = (mat.m_matrix[1][0] * m_matrix[0][1] + mat.m_matrix[1][1] * m_matrix[1][1] 1242 + mat.m_matrix[1][2] * m_matrix[2][1] + mat.m_matrix[1][3] * m_matrix[3][1]); 1243 tmp[1][2] = (mat.m_matrix[1][0] * m_matrix[0][2] + mat.m_matrix[1][1] * m_matrix[1][2] 1244 + mat.m_matrix[1][2] * m_matrix[2][2] + mat.m_matrix[1][3] * m_matrix[3][2]); 1245 tmp[1][3] = (mat.m_matrix[1][0] * m_matrix[0][3] + mat.m_matrix[1][1] * m_matrix[1][3] 1246 + mat.m_matrix[1][2] * m_matrix[2][3] + mat.m_matrix[1][3] * m_matrix[3][3]); 1247 1248 tmp[2][0] = (mat.m_matrix[2][0] * m_matrix[0][0] + mat.m_matrix[2][1] * m_matrix[1][0] 1249 + mat.m_matrix[2][2] * m_matrix[2][0] + mat.m_matrix[2][3] * m_matrix[3][0]); 1250 tmp[2][1] = (mat.m_matrix[2][0] * m_matrix[0][1] + mat.m_matrix[2][1] * m_matrix[1][1] 1251 + mat.m_matrix[2][2] * m_matrix[2][1] + mat.m_matrix[2][3] * m_matrix[3][1]); 1252 tmp[2][2] = (mat.m_matrix[2][0] * m_matrix[0][2] + mat.m_matrix[2][1] * m_matrix[1][2] 1253 + mat.m_matrix[2][2] * m_matrix[2][2] + mat.m_matrix[2][3] * m_matrix[3][2]); 1254 tmp[2][3] = (mat.m_matrix[2][0] * m_matrix[0][3] + mat.m_matrix[2][1] * m_matrix[1][3] 1255 + mat.m_matrix[2][2] * m_matrix[2][3] + mat.m_matrix[2][3] * m_matrix[3][3]); 1256 1257 tmp[3][0] = (mat.m_matrix[3][0] * m_matrix[0][0] + mat.m_matrix[3][1] * m_matrix[1][0] 1258 + mat.m_matrix[3][2] * m_matrix[2][0] + mat.m_matrix[3][3] * m_matrix[3][0]); 1259 tmp[3][1] = (mat.m_matrix[3][0] * m_matrix[0][1] + mat.m_matrix[3][1] * m_matrix[1][1] 1260 + mat.m_matrix[3][2] * m_matrix[2][1] + mat.m_matrix[3][3] * m_matrix[3][1]); 1261 tmp[3][2] = (mat.m_matrix[3][0] * m_matrix[0][2] + mat.m_matrix[3][1] * m_matrix[1][2] 1262 + mat.m_matrix[3][2] * m_matrix[2][2] + mat.m_matrix[3][3] * m_matrix[3][2]); 1263 tmp[3][3] = (mat.m_matrix[3][0] * m_matrix[0][3] + mat.m_matrix[3][1] * m_matrix[1][3] 1264 + mat.m_matrix[3][2] * m_matrix[2][3] + mat.m_matrix[3][3] * m_matrix[3][3]); 1265 1266 setMatrix(tmp); 1267 #endif 1268 return *this; 1269 } 1270 1271 void TransformationMatrix::multVecMatrix(double x, double y, double& resultX, double& resultY) const 1272 { 1273 resultX = m_matrix[3][0] + x * m_matrix[0][0] + y * m_matrix[1][0]; 1274 resultY = m_matrix[3][1] + x * m_matrix[0][1] + y * m_matrix[1][1]; 1275 double w = m_matrix[3][3] + x * m_matrix[0][3] + y * m_matrix[1][3]; 1276 if (w != 1 && w != 0) { 1277 resultX /= w; 1278 resultY /= w; 1279 } 1280 } 1281 1282 void TransformationMatrix::multVecMatrix(double x, double y, double z, double& resultX, double& resultY, double& resultZ) const 1283 { 1284 resultX = m_matrix[3][0] + x * m_matrix[0][0] + y * m_matrix[1][0] + z * m_matrix[2][0]; 1285 resultY = m_matrix[3][1] + x * m_matrix[0][1] + y * m_matrix[1][1] + z * m_matrix[2][1]; 1286 resultZ = m_matrix[3][2] + x * m_matrix[0][2] + y * m_matrix[1][2] + z * m_matrix[2][2]; 1287 double w = m_matrix[3][3] + x * m_matrix[0][3] + y * m_matrix[1][3] + z * m_matrix[2][3]; 1288 if (w != 1 && w != 0) { 1289 resultX /= w; 1290 resultY /= w; 1291 resultZ /= w; 1292 } 1293 } 1294 1295 bool TransformationMatrix::isInvertible() const 1296 { 1297 if (isIdentityOrTranslation()) 1298 return true; 1299 1300 double det = WebCore::determinant4x4(m_matrix); 1301 1302 if (fabs(det) < SMALL_NUMBER) 1303 return false; 1304 1305 return true; 1306 } 1307 1308 TransformationMatrix TransformationMatrix::inverse() const 1309 { 1310 if (isIdentityOrTranslation()) { 1311 // identity matrix 1312 if (m_matrix[3][0] == 0 && m_matrix[3][1] == 0 && m_matrix[3][2] == 0) 1313 return TransformationMatrix(); 1314 1315 // translation 1316 return TransformationMatrix(1, 0, 0, 0, 1317 0, 1, 0, 0, 1318 0, 0, 1, 0, 1319 -m_matrix[3][0], -m_matrix[3][1], -m_matrix[3][2], 1); 1320 } 1321 1322 TransformationMatrix invMat; 1323 bool inverted = WebCore::inverse(m_matrix, invMat.m_matrix); 1324 if (!inverted) 1325 return TransformationMatrix(); 1326 1327 return invMat; 1328 } 1329 1330 void TransformationMatrix::makeAffine() 1331 { 1332 m_matrix[0][2] = 0; 1333 m_matrix[0][3] = 0; 1334 1335 m_matrix[1][2] = 0; 1336 m_matrix[1][3] = 0; 1337 1338 m_matrix[2][0] = 0; 1339 m_matrix[2][1] = 0; 1340 m_matrix[2][2] = 1; 1341 m_matrix[2][3] = 0; 1342 1343 m_matrix[3][2] = 0; 1344 m_matrix[3][3] = 1; 1345 } 1346 1347 AffineTransform TransformationMatrix::toAffineTransform() const 1348 { 1349 return AffineTransform(m_matrix[0][0], m_matrix[0][1], m_matrix[1][0], 1350 m_matrix[1][1], m_matrix[3][0], m_matrix[3][1]); 1351 } 1352 1353 static inline void blendFloat(double& from, double to, double progress) 1354 { 1355 if (from != to) 1356 from = from + (to - from) * progress; 1357 } 1358 1359 void TransformationMatrix::blend(const TransformationMatrix& from, double progress) 1360 { 1361 if (from.isIdentity() && isIdentity()) 1362 return; 1363 1364 // decompose 1365 DecomposedType fromDecomp; 1366 DecomposedType toDecomp; 1367 if (!from.decompose(fromDecomp) || !decompose(toDecomp)) { 1368 if (progress < 0.5) 1369 *this = from; 1370 return; 1371 } 1372 1373 // interpolate 1374 blendFloat(fromDecomp.scaleX, toDecomp.scaleX, progress); 1375 blendFloat(fromDecomp.scaleY, toDecomp.scaleY, progress); 1376 blendFloat(fromDecomp.scaleZ, toDecomp.scaleZ, progress); 1377 blendFloat(fromDecomp.skewXY, toDecomp.skewXY, progress); 1378 blendFloat(fromDecomp.skewXZ, toDecomp.skewXZ, progress); 1379 blendFloat(fromDecomp.skewYZ, toDecomp.skewYZ, progress); 1380 blendFloat(fromDecomp.translateX, toDecomp.translateX, progress); 1381 blendFloat(fromDecomp.translateY, toDecomp.translateY, progress); 1382 blendFloat(fromDecomp.translateZ, toDecomp.translateZ, progress); 1383 blendFloat(fromDecomp.perspectiveX, toDecomp.perspectiveX, progress); 1384 blendFloat(fromDecomp.perspectiveY, toDecomp.perspectiveY, progress); 1385 blendFloat(fromDecomp.perspectiveZ, toDecomp.perspectiveZ, progress); 1386 blendFloat(fromDecomp.perspectiveW, toDecomp.perspectiveW, progress); 1387 1388 slerp(&fromDecomp.quaternionX, &toDecomp.quaternionX, progress); 1389 1390 // recompose 1391 recompose(fromDecomp); 1392 } 1393 1394 bool TransformationMatrix::decompose(DecomposedType& decomp) const 1395 { 1396 if (isIdentity()) { 1397 memset(&decomp, 0, sizeof(decomp)); 1398 decomp.perspectiveW = 1; 1399 decomp.scaleX = 1; 1400 decomp.scaleY = 1; 1401 decomp.scaleZ = 1; 1402 } 1403 1404 if (!WebCore::decompose(m_matrix, decomp)) 1405 return false; 1406 return true; 1407 } 1408 1409 void TransformationMatrix::recompose(const DecomposedType& decomp) 1410 { 1411 makeIdentity(); 1412 1413 // first apply perspective 1414 m_matrix[0][3] = decomp.perspectiveX; 1415 m_matrix[1][3] = decomp.perspectiveY; 1416 m_matrix[2][3] = decomp.perspectiveZ; 1417 m_matrix[3][3] = decomp.perspectiveW; 1418 1419 // now translate 1420 translate3d(decomp.translateX, decomp.translateY, decomp.translateZ); 1421 1422 // apply rotation 1423 double xx = decomp.quaternionX * decomp.quaternionX; 1424 double xy = decomp.quaternionX * decomp.quaternionY; 1425 double xz = decomp.quaternionX * decomp.quaternionZ; 1426 double xw = decomp.quaternionX * decomp.quaternionW; 1427 double yy = decomp.quaternionY * decomp.quaternionY; 1428 double yz = decomp.quaternionY * decomp.quaternionZ; 1429 double yw = decomp.quaternionY * decomp.quaternionW; 1430 double zz = decomp.quaternionZ * decomp.quaternionZ; 1431 double zw = decomp.quaternionZ * decomp.quaternionW; 1432 1433 // Construct a composite rotation matrix from the quaternion values 1434 TransformationMatrix rotationMatrix(1 - 2 * (yy + zz), 2 * (xy - zw), 2 * (xz + yw), 0, 1435 2 * (xy + zw), 1 - 2 * (xx + zz), 2 * (yz - xw), 0, 1436 2 * (xz - yw), 2 * (yz + xw), 1 - 2 * (xx + yy), 0, 1437 0, 0, 0, 1); 1438 1439 multiply(rotationMatrix); 1440 1441 // now apply skew 1442 if (decomp.skewYZ) { 1443 TransformationMatrix tmp; 1444 tmp.setM32(decomp.skewYZ); 1445 multiply(tmp); 1446 } 1447 1448 if (decomp.skewXZ) { 1449 TransformationMatrix tmp; 1450 tmp.setM31(decomp.skewXZ); 1451 multiply(tmp); 1452 } 1453 1454 if (decomp.skewXY) { 1455 TransformationMatrix tmp; 1456 tmp.setM21(decomp.skewXY); 1457 multiply(tmp); 1458 } 1459 1460 // finally, apply scale 1461 scale3d(decomp.scaleX, decomp.scaleY, decomp.scaleZ); 1462 } 1463 1464 bool TransformationMatrix::isIntegerTranslation() const 1465 { 1466 if (!isIdentityOrTranslation()) 1467 return false; 1468 1469 // Check for translate Z. 1470 if (m_matrix[3][2]) 1471 return false; 1472 1473 // Check for non-integer translate X/Y. 1474 if (static_cast<int>(m_matrix[3][0]) != m_matrix[3][0] || static_cast<int>(m_matrix[3][1]) != m_matrix[3][1]) 1475 return false; 1476 1477 return true; 1478 } 1479 1480 TransformationMatrix TransformationMatrix::to2dTransform() const 1481 { 1482 return TransformationMatrix(m_matrix[0][0], m_matrix[0][1], 0, m_matrix[0][3], 1483 m_matrix[1][0], m_matrix[1][1], 0, m_matrix[1][3], 1484 0, 0, 1, 0, 1485 m_matrix[3][0], m_matrix[3][1], 0, m_matrix[3][3]); 1486 } 1487 1488 void TransformationMatrix::toColumnMajorFloatArray(FloatMatrix4& result) const 1489 { 1490 result[0] = m11(); 1491 result[1] = m12(); 1492 result[2] = m13(); 1493 result[3] = m14(); 1494 result[4] = m21(); 1495 result[5] = m22(); 1496 result[6] = m23(); 1497 result[7] = m24(); 1498 result[8] = m31(); 1499 result[9] = m32(); 1500 result[10] = m33(); 1501 result[11] = m34(); 1502 result[12] = m41(); 1503 result[13] = m42(); 1504 result[14] = m43(); 1505 result[15] = m44(); 1506 } 1507 1508 bool TransformationMatrix::isBackFaceVisible() const 1509 { 1510 // Back-face visibility is determined by transforming the normal vector (0, 0, 1) and 1511 // checking the sign of the resulting z component. However, normals cannot be 1512 // transformed by the original matrix, they require being transformed by the 1513 // inverse-transpose. 1514 // 1515 // Since we know we will be using (0, 0, 1), and we only care about the z-component of 1516 // the transformed normal, then we only need the m33() element of the 1517 // inverse-transpose. Therefore we do not need the transpose. 1518 // 1519 // Additionally, if we only need the m33() element, we do not need to compute a full 1520 // inverse. Instead, knowing the inverse of a matrix is adjoint(matrix) / determinant, 1521 // we can simply compute the m33() of the adjoint (adjugate) matrix, without computing 1522 // the full adjoint. 1523 1524 double determinant = WebCore::determinant4x4(m_matrix); 1525 1526 // If the matrix is not invertible, then we assume its backface is not visible. 1527 if (fabs(determinant) < SMALL_NUMBER) 1528 return false; 1529 1530 double cofactor33 = determinant3x3(m11(), m12(), m14(), m21(), m22(), m24(), m41(), m42(), m44()); 1531 double zComponentOfTransformedNormal = cofactor33 / determinant; 1532 1533 return zComponentOfTransformedNormal < 0; 1534 } 1535 1536 SkMatrix44 TransformationMatrix::toSkMatrix44(const TransformationMatrix& matrix) 1537 { 1538 SkMatrix44 ret(SkMatrix44::kUninitialized_Constructor); 1539 ret.setDouble(0, 0, matrix.m11()); 1540 ret.setDouble(0, 1, matrix.m21()); 1541 ret.setDouble(0, 2, matrix.m31()); 1542 ret.setDouble(0, 3, matrix.m41()); 1543 ret.setDouble(1, 0, matrix.m12()); 1544 ret.setDouble(1, 1, matrix.m22()); 1545 ret.setDouble(1, 2, matrix.m32()); 1546 ret.setDouble(1, 3, matrix.m42()); 1547 ret.setDouble(2, 0, matrix.m13()); 1548 ret.setDouble(2, 1, matrix.m23()); 1549 ret.setDouble(2, 2, matrix.m33()); 1550 ret.setDouble(2, 3, matrix.m43()); 1551 ret.setDouble(3, 0, matrix.m14()); 1552 ret.setDouble(3, 1, matrix.m24()); 1553 ret.setDouble(3, 2, matrix.m34()); 1554 ret.setDouble(3, 3, matrix.m44()); 1555 return ret; 1556 } 1557 1558 } 1559