• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1 /*
2  * Copyright (C) 2005, 2006 Apple Computer, Inc.  All rights reserved.
3  * Copyright (C) 2009 Torch Mobile, Inc.
4  *
5  * Redistribution and use in source and binary forms, with or without
6  * modification, are permitted provided that the following conditions
7  * are met:
8  * 1. Redistributions of source code must retain the above copyright
9  *    notice, this list of conditions and the following disclaimer.
10  * 2. Redistributions in binary form must reproduce the above copyright
11  *    notice, this list of conditions and the following disclaimer in the
12  *    documentation and/or other materials provided with the distribution.
13  *
14  * THIS SOFTWARE IS PROVIDED BY APPLE COMPUTER, INC. ``AS IS'' AND ANY
15  * EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
16  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
17  * PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL APPLE COMPUTER, INC. OR
18  * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
19  * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
20  * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
21  * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY
22  * OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
23  * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
24  * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
25  */
26 
27 #include "config.h"
28 #include "TransformationMatrix.h"
29 
30 #include "FloatPoint3D.h"
31 #include "FloatRect.h"
32 #include "FloatQuad.h"
33 #include "IntRect.h"
34 
35 #include <wtf/Assertions.h>
36 #include <wtf/MathExtras.h>
37 
38 namespace WebCore {
39 
40 //
41 // Supporting Math Functions
42 //
43 // This is a set of function from various places (attributed inline) to do things like
44 // inversion and decomposition of a 4x4 matrix. They are used throughout the code
45 //
46 
47 //
48 // Adapted from Matrix Inversion by Richard Carling, Graphics Gems <http://tog.acm.org/GraphicsGems/index.html>.
49 
50 // EULA: The Graphics Gems code is copyright-protected. In other words, you cannot claim the text of the code
51 // as your own and resell it. Using the code is permitted in any program, product, or library, non-commercial
52 // or commercial. Giving credit is not required, though is a nice gesture. The code comes as-is, and if there
53 // are any flaws or problems with any Gems code, nobody involved with Gems - authors, editors, publishers, or
54 // webmasters - are to be held responsible. Basically, don't be a jerk, and remember that anything free comes
55 // with no guarantee.
56 
57 typedef double Vector4[4];
58 typedef double Vector3[3];
59 
60 const double SMALL_NUMBER = 1.e-8;
61 
62 // inverse(original_matrix, inverse_matrix)
63 //
64 // calculate the inverse of a 4x4 matrix
65 //
66 // -1
67 // A  = ___1__ adjoint A
68 //       det A
69 
70 //  double = determinant2x2(double a, double b, double c, double d)
71 //
72 //  calculate the determinant of a 2x2 matrix.
73 
determinant2x2(double a,double b,double c,double d)74 static double determinant2x2(double a, double b, double c, double d)
75 {
76     return a * d - b * c;
77 }
78 
79 //  double = determinant3x3(a1, a2, a3, b1, b2, b3, c1, c2, c3)
80 //
81 //  Calculate the determinant of a 3x3 matrix
82 //  in the form
83 //
84 //      | a1,  b1,  c1 |
85 //      | a2,  b2,  c2 |
86 //      | a3,  b3,  c3 |
87 
determinant3x3(double a1,double a2,double a3,double b1,double b2,double b3,double c1,double c2,double c3)88 static double determinant3x3(double a1, double a2, double a3, double b1, double b2, double b3, double c1, double c2, double c3)
89 {
90     return a1 * determinant2x2(b2, b3, c2, c3)
91          - b1 * determinant2x2(a2, a3, c2, c3)
92          + c1 * determinant2x2(a2, a3, b2, b3);
93 }
94 
95 //  double = determinant4x4(matrix)
96 //
97 //  calculate the determinant of a 4x4 matrix.
98 
determinant4x4(const TransformationMatrix::Matrix4 & m)99 static double determinant4x4(const TransformationMatrix::Matrix4& m)
100 {
101     // Assign to individual variable names to aid selecting
102     // correct elements
103 
104     double a1 = m[0][0];
105     double b1 = m[0][1];
106     double c1 = m[0][2];
107     double d1 = m[0][3];
108 
109     double a2 = m[1][0];
110     double b2 = m[1][1];
111     double c2 = m[1][2];
112     double d2 = m[1][3];
113 
114     double a3 = m[2][0];
115     double b3 = m[2][1];
116     double c3 = m[2][2];
117     double d3 = m[2][3];
118 
119     double a4 = m[3][0];
120     double b4 = m[3][1];
121     double c4 = m[3][2];
122     double d4 = m[3][3];
123 
124     return a1 * determinant3x3(b2, b3, b4, c2, c3, c4, d2, d3, d4)
125          - b1 * determinant3x3(a2, a3, a4, c2, c3, c4, d2, d3, d4)
126          + c1 * determinant3x3(a2, a3, a4, b2, b3, b4, d2, d3, d4)
127          - d1 * determinant3x3(a2, a3, a4, b2, b3, b4, c2, c3, c4);
128 }
129 
130 // adjoint( original_matrix, inverse_matrix )
131 //
132 //   calculate the adjoint of a 4x4 matrix
133 //
134 //    Let  a   denote the minor determinant of matrix A obtained by
135 //         ij
136 //
137 //    deleting the ith row and jth column from A.
138 //
139 //                  i+j
140 //   Let  b   = (-1)    a
141 //        ij            ji
142 //
143 //  The matrix B = (b  ) is the adjoint of A
144 //                   ij
145 
adjoint(const TransformationMatrix::Matrix4 & matrix,TransformationMatrix::Matrix4 & result)146 static void adjoint(const TransformationMatrix::Matrix4& matrix, TransformationMatrix::Matrix4& result)
147 {
148     // Assign to individual variable names to aid
149     // selecting correct values
150     double a1 = matrix[0][0];
151     double b1 = matrix[0][1];
152     double c1 = matrix[0][2];
153     double d1 = matrix[0][3];
154 
155     double a2 = matrix[1][0];
156     double b2 = matrix[1][1];
157     double c2 = matrix[1][2];
158     double d2 = matrix[1][3];
159 
160     double a3 = matrix[2][0];
161     double b3 = matrix[2][1];
162     double c3 = matrix[2][2];
163     double d3 = matrix[2][3];
164 
165     double a4 = matrix[3][0];
166     double b4 = matrix[3][1];
167     double c4 = matrix[3][2];
168     double d4 = matrix[3][3];
169 
170     // Row column labeling reversed since we transpose rows & columns
171     result[0][0]  =   determinant3x3(b2, b3, b4, c2, c3, c4, d2, d3, d4);
172     result[1][0]  = - determinant3x3(a2, a3, a4, c2, c3, c4, d2, d3, d4);
173     result[2][0]  =   determinant3x3(a2, a3, a4, b2, b3, b4, d2, d3, d4);
174     result[3][0]  = - determinant3x3(a2, a3, a4, b2, b3, b4, c2, c3, c4);
175 
176     result[0][1]  = - determinant3x3(b1, b3, b4, c1, c3, c4, d1, d3, d4);
177     result[1][1]  =   determinant3x3(a1, a3, a4, c1, c3, c4, d1, d3, d4);
178     result[2][1]  = - determinant3x3(a1, a3, a4, b1, b3, b4, d1, d3, d4);
179     result[3][1]  =   determinant3x3(a1, a3, a4, b1, b3, b4, c1, c3, c4);
180 
181     result[0][2]  =   determinant3x3(b1, b2, b4, c1, c2, c4, d1, d2, d4);
182     result[1][2]  = - determinant3x3(a1, a2, a4, c1, c2, c4, d1, d2, d4);
183     result[2][2]  =   determinant3x3(a1, a2, a4, b1, b2, b4, d1, d2, d4);
184     result[3][2]  = - determinant3x3(a1, a2, a4, b1, b2, b4, c1, c2, c4);
185 
186     result[0][3]  = - determinant3x3(b1, b2, b3, c1, c2, c3, d1, d2, d3);
187     result[1][3]  =   determinant3x3(a1, a2, a3, c1, c2, c3, d1, d2, d3);
188     result[2][3]  = - determinant3x3(a1, a2, a3, b1, b2, b3, d1, d2, d3);
189     result[3][3]  =   determinant3x3(a1, a2, a3, b1, b2, b3, c1, c2, c3);
190 }
191 
192 // Returns false if the matrix is not invertible
inverse(const TransformationMatrix::Matrix4 & matrix,TransformationMatrix::Matrix4 & result)193 static bool inverse(const TransformationMatrix::Matrix4& matrix, TransformationMatrix::Matrix4& result)
194 {
195     // Calculate the adjoint matrix
196     adjoint(matrix, result);
197 
198     // Calculate the 4x4 determinant
199     // If the determinant is zero,
200     // then the inverse matrix is not unique.
201     double det = determinant4x4(matrix);
202 
203     if (fabs(det) < SMALL_NUMBER)
204         return false;
205 
206     // Scale the adjoint matrix to get the inverse
207 
208     for (int i = 0; i < 4; i++)
209         for (int j = 0; j < 4; j++)
210             result[i][j] = result[i][j] / det;
211 
212     return true;
213 }
214 
215 // End of code adapted from Matrix Inversion by Richard Carling
216 
217 // Perform a decomposition on the passed matrix, return false if unsuccessful
218 // From Graphics Gems: unmatrix.c
219 
220 // Transpose rotation portion of matrix a, return b
transposeMatrix4(const TransformationMatrix::Matrix4 & a,TransformationMatrix::Matrix4 & b)221 static void transposeMatrix4(const TransformationMatrix::Matrix4& a, TransformationMatrix::Matrix4& b)
222 {
223     for (int i = 0; i < 4; i++)
224         for (int j = 0; j < 4; j++)
225             b[i][j] = a[j][i];
226 }
227 
228 // Multiply a homogeneous point by a matrix and return the transformed point
v4MulPointByMatrix(const Vector4 p,const TransformationMatrix::Matrix4 & m,Vector4 result)229 static void v4MulPointByMatrix(const Vector4 p, const TransformationMatrix::Matrix4& m, Vector4 result)
230 {
231     result[0] = (p[0] * m[0][0]) + (p[1] * m[1][0]) +
232                 (p[2] * m[2][0]) + (p[3] * m[3][0]);
233     result[1] = (p[0] * m[0][1]) + (p[1] * m[1][1]) +
234                 (p[2] * m[2][1]) + (p[3] * m[3][1]);
235     result[2] = (p[0] * m[0][2]) + (p[1] * m[1][2]) +
236                 (p[2] * m[2][2]) + (p[3] * m[3][2]);
237     result[3] = (p[0] * m[0][3]) + (p[1] * m[1][3]) +
238                 (p[2] * m[2][3]) + (p[3] * m[3][3]);
239 }
240 
v3Length(Vector3 a)241 static double v3Length(Vector3 a)
242 {
243     return sqrt((a[0] * a[0]) + (a[1] * a[1]) + (a[2] * a[2]));
244 }
245 
v3Scale(Vector3 v,double desiredLength)246 static void v3Scale(Vector3 v, double desiredLength)
247 {
248     double len = v3Length(v);
249     if (len != 0) {
250         double l = desiredLength / len;
251         v[0] *= l;
252         v[1] *= l;
253         v[2] *= l;
254     }
255 }
256 
v3Dot(const Vector3 a,const Vector3 b)257 static double v3Dot(const Vector3 a, const Vector3 b)
258 {
259     return (a[0] * b[0]) + (a[1] * b[1]) + (a[2] * b[2]);
260 }
261 
262 // Make a linear combination of two vectors and return the result.
263 // result = (a * ascl) + (b * bscl)
v3Combine(const Vector3 a,const Vector3 b,Vector3 result,double ascl,double bscl)264 static void v3Combine(const Vector3 a, const Vector3 b, Vector3 result, double ascl, double bscl)
265 {
266     result[0] = (ascl * a[0]) + (bscl * b[0]);
267     result[1] = (ascl * a[1]) + (bscl * b[1]);
268     result[2] = (ascl * a[2]) + (bscl * b[2]);
269 }
270 
271 // Return the cross product result = a cross b */
v3Cross(const Vector3 a,const Vector3 b,Vector3 result)272 static void v3Cross(const Vector3 a, const Vector3 b, Vector3 result)
273 {
274     result[0] = (a[1] * b[2]) - (a[2] * b[1]);
275     result[1] = (a[2] * b[0]) - (a[0] * b[2]);
276     result[2] = (a[0] * b[1]) - (a[1] * b[0]);
277 }
278 
decompose(const TransformationMatrix::Matrix4 & mat,TransformationMatrix::DecomposedType & result)279 static bool decompose(const TransformationMatrix::Matrix4& mat, TransformationMatrix::DecomposedType& result)
280 {
281     TransformationMatrix::Matrix4 localMatrix;
282     memcpy(localMatrix, mat, sizeof(TransformationMatrix::Matrix4));
283 
284     // Normalize the matrix.
285     if (localMatrix[3][3] == 0)
286         return false;
287 
288     int i, j;
289     for (i = 0; i < 4; i++)
290         for (j = 0; j < 4; j++)
291             localMatrix[i][j] /= localMatrix[3][3];
292 
293     // perspectiveMatrix is used to solve for perspective, but it also provides
294     // an easy way to test for singularity of the upper 3x3 component.
295     TransformationMatrix::Matrix4 perspectiveMatrix;
296     memcpy(perspectiveMatrix, localMatrix, sizeof(TransformationMatrix::Matrix4));
297     for (i = 0; i < 3; i++)
298         perspectiveMatrix[i][3] = 0;
299     perspectiveMatrix[3][3] = 1;
300 
301     if (determinant4x4(perspectiveMatrix) == 0)
302         return false;
303 
304     // First, isolate perspective.  This is the messiest.
305     if (localMatrix[0][3] != 0 || localMatrix[1][3] != 0 || localMatrix[2][3] != 0) {
306         // rightHandSide is the right hand side of the equation.
307         Vector4 rightHandSide;
308         rightHandSide[0] = localMatrix[0][3];
309         rightHandSide[1] = localMatrix[1][3];
310         rightHandSide[2] = localMatrix[2][3];
311         rightHandSide[3] = localMatrix[3][3];
312 
313         // Solve the equation by inverting perspectiveMatrix and multiplying
314         // rightHandSide by the inverse.  (This is the easiest way, not
315         // necessarily the best.)
316         TransformationMatrix::Matrix4 inversePerspectiveMatrix, transposedInversePerspectiveMatrix;
317         inverse(perspectiveMatrix, inversePerspectiveMatrix);
318         transposeMatrix4(inversePerspectiveMatrix, transposedInversePerspectiveMatrix);
319 
320         Vector4 perspectivePoint;
321         v4MulPointByMatrix(rightHandSide, transposedInversePerspectiveMatrix, perspectivePoint);
322 
323         result.perspectiveX = perspectivePoint[0];
324         result.perspectiveY = perspectivePoint[1];
325         result.perspectiveZ = perspectivePoint[2];
326         result.perspectiveW = perspectivePoint[3];
327 
328         // Clear the perspective partition
329         localMatrix[0][3] = localMatrix[1][3] = localMatrix[2][3] = 0;
330         localMatrix[3][3] = 1;
331     } else {
332         // No perspective.
333         result.perspectiveX = result.perspectiveY = result.perspectiveZ = 0;
334         result.perspectiveW = 1;
335     }
336 
337     // Next take care of translation (easy).
338     result.translateX = localMatrix[3][0];
339     localMatrix[3][0] = 0;
340     result.translateY = localMatrix[3][1];
341     localMatrix[3][1] = 0;
342     result.translateZ = localMatrix[3][2];
343     localMatrix[3][2] = 0;
344 
345     // Vector4 type and functions need to be added to the common set.
346     Vector3 row[3], pdum3;
347 
348     // Now get scale and shear.
349     for (i = 0; i < 3; i++) {
350         row[i][0] = localMatrix[i][0];
351         row[i][1] = localMatrix[i][1];
352         row[i][2] = localMatrix[i][2];
353     }
354 
355     // Compute X scale factor and normalize first row.
356     result.scaleX = v3Length(row[0]);
357     v3Scale(row[0], 1.0);
358 
359     // Compute XY shear factor and make 2nd row orthogonal to 1st.
360     result.skewXY = v3Dot(row[0], row[1]);
361     v3Combine(row[1], row[0], row[1], 1.0, -result.skewXY);
362 
363     // Now, compute Y scale and normalize 2nd row.
364     result.scaleY = v3Length(row[1]);
365     v3Scale(row[1], 1.0);
366     result.skewXY /= result.scaleY;
367 
368     // Compute XZ and YZ shears, orthogonalize 3rd row.
369     result.skewXZ = v3Dot(row[0], row[2]);
370     v3Combine(row[2], row[0], row[2], 1.0, -result.skewXZ);
371     result.skewYZ = v3Dot(row[1], row[2]);
372     v3Combine(row[2], row[1], row[2], 1.0, -result.skewYZ);
373 
374     // Next, get Z scale and normalize 3rd row.
375     result.scaleZ = v3Length(row[2]);
376     v3Scale(row[2], 1.0);
377     result.skewXZ /= result.scaleZ;
378     result.skewYZ /= result.scaleZ;
379 
380     // At this point, the matrix (in rows[]) is orthonormal.
381     // Check for a coordinate system flip.  If the determinant
382     // is -1, then negate the matrix and the scaling factors.
383     v3Cross(row[1], row[2], pdum3);
384     if (v3Dot(row[0], pdum3) < 0) {
385         for (i = 0; i < 3; i++) {
386             result.scaleX *= -1;
387             row[i][0] *= -1;
388             row[i][1] *= -1;
389             row[i][2] *= -1;
390         }
391     }
392 
393     // Now, get the rotations out, as described in the gem.
394 
395     // FIXME - Add the ability to return either quaternions (which are
396     // easier to recompose with) or Euler angles (rx, ry, rz), which
397     // are easier for authors to deal with. The latter will only be useful
398     // when we fix https://bugs.webkit.org/show_bug.cgi?id=23799, so I
399     // will leave the Euler angle code here for now.
400 
401     // ret.rotateY = asin(-row[0][2]);
402     // if (cos(ret.rotateY) != 0) {
403     //     ret.rotateX = atan2(row[1][2], row[2][2]);
404     //     ret.rotateZ = atan2(row[0][1], row[0][0]);
405     // } else {
406     //     ret.rotateX = atan2(-row[2][0], row[1][1]);
407     //     ret.rotateZ = 0;
408     // }
409 
410     double s, t, x, y, z, w;
411 
412     t = row[0][0] + row[1][1] + row[2][2] + 1.0;
413 
414     if (t > 1e-4) {
415         s = 0.5 / sqrt(t);
416         w = 0.25 / s;
417         x = (row[2][1] - row[1][2]) * s;
418         y = (row[0][2] - row[2][0]) * s;
419         z = (row[1][0] - row[0][1]) * s;
420     } else if (row[0][0] > row[1][1] && row[0][0] > row[2][2]) {
421         s = sqrt (1.0 + row[0][0] - row[1][1] - row[2][2]) * 2.0; // S=4*qx
422         x = 0.25 * s;
423         y = (row[0][1] + row[1][0]) / s;
424         z = (row[0][2] + row[2][0]) / s;
425         w = (row[2][1] - row[1][2]) / s;
426     } else if (row[1][1] > row[2][2]) {
427         s = sqrt (1.0 + row[1][1] - row[0][0] - row[2][2]) * 2.0; // S=4*qy
428         x = (row[0][1] + row[1][0]) / s;
429         y = 0.25 * s;
430         z = (row[1][2] + row[2][1]) / s;
431         w = (row[0][2] - row[2][0]) / s;
432     } else {
433         s = sqrt(1.0 + row[2][2] - row[0][0] - row[1][1]) * 2.0; // S=4*qz
434         x = (row[0][2] + row[2][0]) / s;
435         y = (row[1][2] + row[2][1]) / s;
436         z = 0.25 * s;
437         w = (row[1][0] - row[0][1]) / s;
438     }
439 
440     result.quaternionX = x;
441     result.quaternionY = y;
442     result.quaternionZ = z;
443     result.quaternionW = w;
444 
445     return true;
446 }
447 
448 // Perform a spherical linear interpolation between the two
449 // passed quaternions with 0 <= t <= 1
slerp(double qa[4],const double qb[4],double t)450 static void slerp(double qa[4], const double qb[4], double t)
451 {
452     double ax, ay, az, aw;
453     double bx, by, bz, bw;
454     double cx, cy, cz, cw;
455     double angle;
456     double th, invth, scale, invscale;
457 
458     ax = qa[0]; ay = qa[1]; az = qa[2]; aw = qa[3];
459     bx = qb[0]; by = qb[1]; bz = qb[2]; bw = qb[3];
460 
461     angle = ax * bx + ay * by + az * bz + aw * bw;
462 
463     if (angle < 0.0) {
464         ax = -ax; ay = -ay;
465         az = -az; aw = -aw;
466         angle = -angle;
467     }
468 
469     if (angle + 1.0 > .05) {
470         if (1.0 - angle >= .05) {
471             th = acos (angle);
472             invth = 1.0 / sin (th);
473             scale = sin (th * (1.0 - t)) * invth;
474             invscale = sin (th * t) * invth;
475         } else {
476             scale = 1.0 - t;
477             invscale = t;
478         }
479     } else {
480         bx = -ay;
481         by = ax;
482         bz = -aw;
483         bw = az;
484         scale = sin(piDouble * (.5 - t));
485         invscale = sin (piDouble * t);
486     }
487 
488     cx = ax * scale + bx * invscale;
489     cy = ay * scale + by * invscale;
490     cz = az * scale + bz * invscale;
491     cw = aw * scale + bw * invscale;
492 
493     qa[0] = cx; qa[1] = cy; qa[2] = cz; qa[3] = cw;
494 }
495 
496 // End of Supporting Math Functions
497 
scale(double s)498 TransformationMatrix& TransformationMatrix::scale(double s)
499 {
500     return scaleNonUniform(s, s);
501 }
502 
rotateFromVector(double x,double y)503 TransformationMatrix& TransformationMatrix::rotateFromVector(double x, double y)
504 {
505     return rotate(rad2deg(atan2(y, x)));
506 }
507 
flipX()508 TransformationMatrix& TransformationMatrix::flipX()
509 {
510     return scaleNonUniform(-1.0f, 1.0f);
511 }
512 
flipY()513 TransformationMatrix& TransformationMatrix::flipY()
514 {
515     return scaleNonUniform(1.0f, -1.0f);
516 }
517 
projectPoint(const FloatPoint & p) const518 FloatPoint TransformationMatrix::projectPoint(const FloatPoint& p) const
519 {
520     // This is basically raytracing. We have a point in the destination
521     // plane with z=0, and we cast a ray parallel to the z-axis from that
522     // point to find the z-position at which it intersects the z=0 plane
523     // with the transform applied. Once we have that point we apply the
524     // inverse transform to find the corresponding point in the source
525     // space.
526     //
527     // Given a plane with normal Pn, and a ray starting at point R0 and
528     // with direction defined by the vector Rd, we can find the
529     // intersection point as a distance d from R0 in units of Rd by:
530     //
531     // d = -dot (Pn', R0) / dot (Pn', Rd)
532 
533     double x = p.x();
534     double y = p.y();
535     double z = -(m13() * x + m23() * y + m43()) / m33();
536 
537     double outX = x * m11() + y * m21() + z * m31() + m41();
538     double outY = x * m12() + y * m22() + z * m32() + m42();
539 
540     double w = x * m14() + y * m24() + z * m34() + m44();
541     if (w != 1 && w != 0) {
542         outX /= w;
543         outY /= w;
544     }
545 
546     return FloatPoint(static_cast<float>(outX), static_cast<float>(outY));
547 }
548 
projectQuad(const FloatQuad & q) const549 FloatQuad TransformationMatrix::projectQuad(const FloatQuad& q) const
550 {
551     FloatQuad projectedQuad;
552     projectedQuad.setP1(projectPoint(q.p1()));
553     projectedQuad.setP2(projectPoint(q.p2()));
554     projectedQuad.setP3(projectPoint(q.p3()));
555     projectedQuad.setP4(projectPoint(q.p4()));
556     return projectedQuad;
557 }
558 
mapPoint(const FloatPoint & p) const559 FloatPoint TransformationMatrix::mapPoint(const FloatPoint& p) const
560 {
561     if (isIdentityOrTranslation())
562         return FloatPoint(p.x() + static_cast<float>(m_matrix[3][0]), p.y() + static_cast<float>(m_matrix[3][1]));
563 
564     double x, y;
565     multVecMatrix(p.x(), p.y(), x, y);
566     return FloatPoint(static_cast<float>(x), static_cast<float>(y));
567 }
568 
mapPoint(const FloatPoint3D & p) const569 FloatPoint3D TransformationMatrix::mapPoint(const FloatPoint3D& p) const
570 {
571     if (isIdentityOrTranslation())
572         return FloatPoint3D(p.x() + static_cast<float>(m_matrix[3][0]),
573                             p.y() + static_cast<float>(m_matrix[3][1]),
574                             p.z() + static_cast<float>(m_matrix[3][2]));
575 
576     double x, y, z;
577     multVecMatrix(p.x(), p.y(), p.z(), x, y, z);
578     return FloatPoint3D(static_cast<float>(x), static_cast<float>(y), static_cast<float>(z));
579 }
580 
mapRect(const IntRect & rect) const581 IntRect TransformationMatrix::mapRect(const IntRect &rect) const
582 {
583     return enclosingIntRect(mapRect(FloatRect(rect)));
584 }
585 
mapRect(const FloatRect & r) const586 FloatRect TransformationMatrix::mapRect(const FloatRect& r) const
587 {
588     if (isIdentityOrTranslation()) {
589         FloatRect mappedRect(r);
590         mappedRect.move(static_cast<float>(m_matrix[3][0]), static_cast<float>(m_matrix[3][1]));
591         return mappedRect;
592     }
593 
594     FloatQuad resultQuad = mapQuad(FloatQuad(r));
595     return resultQuad.boundingBox();
596 }
597 
mapQuad(const FloatQuad & q) const598 FloatQuad TransformationMatrix::mapQuad(const FloatQuad& q) const
599 {
600     if (isIdentityOrTranslation()) {
601         FloatQuad mappedQuad(q);
602         mappedQuad.move(static_cast<float>(m_matrix[3][0]), static_cast<float>(m_matrix[3][1]));
603         return mappedQuad;
604     }
605 
606     FloatQuad result;
607     result.setP1(mapPoint(q.p1()));
608     result.setP2(mapPoint(q.p2()));
609     result.setP3(mapPoint(q.p3()));
610     result.setP4(mapPoint(q.p4()));
611     return result;
612 }
613 
scaleNonUniform(double sx,double sy)614 TransformationMatrix& TransformationMatrix::scaleNonUniform(double sx, double sy)
615 {
616     TransformationMatrix mat;
617     mat.m_matrix[0][0] = sx;
618     mat.m_matrix[1][1] = sy;
619 
620     multLeft(mat);
621     return *this;
622 }
623 
scale3d(double sx,double sy,double sz)624 TransformationMatrix& TransformationMatrix::scale3d(double sx, double sy, double sz)
625 {
626     TransformationMatrix mat;
627     mat.m_matrix[0][0] = sx;
628     mat.m_matrix[1][1] = sy;
629     mat.m_matrix[2][2] = sz;
630 
631     multLeft(mat);
632     return *this;
633 }
634 
rotate3d(double x,double y,double z,double angle)635 TransformationMatrix& TransformationMatrix::rotate3d(double x, double y, double z, double angle)
636 {
637     // angles are in degrees. Switch to radians
638     angle = deg2rad(angle);
639 
640     angle /= 2.0f;
641     double sinA = sin(angle);
642     double cosA = cos(angle);
643     double sinA2 = sinA * sinA;
644 
645     // normalize
646     double length = sqrt(x * x + y * y + z * z);
647     if (length == 0) {
648         // bad vector, just use something reasonable
649         x = 0;
650         y = 0;
651         z = 1;
652     } else if (length != 1) {
653         x /= length;
654         y /= length;
655         z /= length;
656     }
657 
658     TransformationMatrix mat;
659 
660     // optimize case where axis is along major axis
661     if (x == 1.0f && y == 0.0f && z == 0.0f) {
662         mat.m_matrix[0][0] = 1.0f;
663         mat.m_matrix[0][1] = 0.0f;
664         mat.m_matrix[0][2] = 0.0f;
665         mat.m_matrix[1][0] = 0.0f;
666         mat.m_matrix[1][1] = 1.0f - 2.0f * sinA2;
667         mat.m_matrix[1][2] = 2.0f * sinA * cosA;
668         mat.m_matrix[2][0] = 0.0f;
669         mat.m_matrix[2][1] = -2.0f * sinA * cosA;
670         mat.m_matrix[2][2] = 1.0f - 2.0f * sinA2;
671         mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0f;
672         mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0f;
673         mat.m_matrix[3][3] = 1.0f;
674     } else if (x == 0.0f && y == 1.0f && z == 0.0f) {
675         mat.m_matrix[0][0] = 1.0f - 2.0f * sinA2;
676         mat.m_matrix[0][1] = 0.0f;
677         mat.m_matrix[0][2] = -2.0f * sinA * cosA;
678         mat.m_matrix[1][0] = 0.0f;
679         mat.m_matrix[1][1] = 1.0f;
680         mat.m_matrix[1][2] = 0.0f;
681         mat.m_matrix[2][0] = 2.0f * sinA * cosA;
682         mat.m_matrix[2][1] = 0.0f;
683         mat.m_matrix[2][2] = 1.0f - 2.0f * sinA2;
684         mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0f;
685         mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0f;
686         mat.m_matrix[3][3] = 1.0f;
687     } else if (x == 0.0f && y == 0.0f && z == 1.0f) {
688         mat.m_matrix[0][0] = 1.0f - 2.0f * sinA2;
689         mat.m_matrix[0][1] = 2.0f * sinA * cosA;
690         mat.m_matrix[0][2] = 0.0f;
691         mat.m_matrix[1][0] = -2.0f * sinA * cosA;
692         mat.m_matrix[1][1] = 1.0f - 2.0f * sinA2;
693         mat.m_matrix[1][2] = 0.0f;
694         mat.m_matrix[2][0] = 0.0f;
695         mat.m_matrix[2][1] = 0.0f;
696         mat.m_matrix[2][2] = 1.0f;
697         mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0f;
698         mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0f;
699         mat.m_matrix[3][3] = 1.0f;
700     } else {
701         double x2 = x*x;
702         double y2 = y*y;
703         double z2 = z*z;
704 
705         mat.m_matrix[0][0] = 1.0f - 2.0f * (y2 + z2) * sinA2;
706         mat.m_matrix[0][1] = 2.0f * (x * y * sinA2 + z * sinA * cosA);
707         mat.m_matrix[0][2] = 2.0f * (x * z * sinA2 - y * sinA * cosA);
708         mat.m_matrix[1][0] = 2.0f * (y * x * sinA2 - z * sinA * cosA);
709         mat.m_matrix[1][1] = 1.0f - 2.0f * (z2 + x2) * sinA2;
710         mat.m_matrix[1][2] = 2.0f * (y * z * sinA2 + x * sinA * cosA);
711         mat.m_matrix[2][0] = 2.0f * (z * x * sinA2 + y * sinA * cosA);
712         mat.m_matrix[2][1] = 2.0f * (z * y * sinA2 - x * sinA * cosA);
713         mat.m_matrix[2][2] = 1.0f - 2.0f * (x2 + y2) * sinA2;
714         mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0f;
715         mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0f;
716         mat.m_matrix[3][3] = 1.0f;
717     }
718     multLeft(mat);
719     return *this;
720 }
721 
rotate3d(double rx,double ry,double rz)722 TransformationMatrix& TransformationMatrix::rotate3d(double rx, double ry, double rz)
723 {
724     // angles are in degrees. Switch to radians
725     rx = deg2rad(rx);
726     ry = deg2rad(ry);
727     rz = deg2rad(rz);
728 
729     TransformationMatrix mat;
730 
731     rz /= 2.0f;
732     double sinA = sin(rz);
733     double cosA = cos(rz);
734     double sinA2 = sinA * sinA;
735 
736     mat.m_matrix[0][0] = 1.0f - 2.0f * sinA2;
737     mat.m_matrix[0][1] = 2.0f * sinA * cosA;
738     mat.m_matrix[0][2] = 0.0f;
739     mat.m_matrix[1][0] = -2.0f * sinA * cosA;
740     mat.m_matrix[1][1] = 1.0f - 2.0f * sinA2;
741     mat.m_matrix[1][2] = 0.0f;
742     mat.m_matrix[2][0] = 0.0f;
743     mat.m_matrix[2][1] = 0.0f;
744     mat.m_matrix[2][2] = 1.0f;
745     mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0f;
746     mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0f;
747     mat.m_matrix[3][3] = 1.0f;
748 
749     TransformationMatrix rmat(mat);
750 
751     ry /= 2.0f;
752     sinA = sin(ry);
753     cosA = cos(ry);
754     sinA2 = sinA * sinA;
755 
756     mat.m_matrix[0][0] = 1.0f - 2.0f * sinA2;
757     mat.m_matrix[0][1] = 0.0f;
758     mat.m_matrix[0][2] = -2.0f * sinA * cosA;
759     mat.m_matrix[1][0] = 0.0f;
760     mat.m_matrix[1][1] = 1.0f;
761     mat.m_matrix[1][2] = 0.0f;
762     mat.m_matrix[2][0] = 2.0f * sinA * cosA;
763     mat.m_matrix[2][1] = 0.0f;
764     mat.m_matrix[2][2] = 1.0f - 2.0f * sinA2;
765     mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0f;
766     mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0f;
767     mat.m_matrix[3][3] = 1.0f;
768 
769     rmat.multLeft(mat);
770 
771     rx /= 2.0f;
772     sinA = sin(rx);
773     cosA = cos(rx);
774     sinA2 = sinA * sinA;
775 
776     mat.m_matrix[0][0] = 1.0f;
777     mat.m_matrix[0][1] = 0.0f;
778     mat.m_matrix[0][2] = 0.0f;
779     mat.m_matrix[1][0] = 0.0f;
780     mat.m_matrix[1][1] = 1.0f - 2.0f * sinA2;
781     mat.m_matrix[1][2] = 2.0f * sinA * cosA;
782     mat.m_matrix[2][0] = 0.0f;
783     mat.m_matrix[2][1] = -2.0f * sinA * cosA;
784     mat.m_matrix[2][2] = 1.0f - 2.0f * sinA2;
785     mat.m_matrix[0][3] = mat.m_matrix[1][3] = mat.m_matrix[2][3] = 0.0f;
786     mat.m_matrix[3][0] = mat.m_matrix[3][1] = mat.m_matrix[3][2] = 0.0f;
787     mat.m_matrix[3][3] = 1.0f;
788 
789     rmat.multLeft(mat);
790 
791     multLeft(rmat);
792     return *this;
793 }
794 
translate(double tx,double ty)795 TransformationMatrix& TransformationMatrix::translate(double tx, double ty)
796 {
797     m_matrix[3][0] += tx * m_matrix[0][0] + ty * m_matrix[1][0];
798     m_matrix[3][1] += tx * m_matrix[0][1] + ty * m_matrix[1][1];
799     m_matrix[3][2] += tx * m_matrix[0][2] + ty * m_matrix[1][2];
800     m_matrix[3][3] += tx * m_matrix[0][3] + ty * m_matrix[1][3];
801     return *this;
802 }
803 
translate3d(double tx,double ty,double tz)804 TransformationMatrix& TransformationMatrix::translate3d(double tx, double ty, double tz)
805 {
806     m_matrix[3][0] += tx * m_matrix[0][0] + ty * m_matrix[1][0] + tz * m_matrix[2][0];
807     m_matrix[3][1] += tx * m_matrix[0][1] + ty * m_matrix[1][1] + tz * m_matrix[2][1];
808     m_matrix[3][2] += tx * m_matrix[0][2] + ty * m_matrix[1][2] + tz * m_matrix[2][2];
809     m_matrix[3][3] += tx * m_matrix[0][3] + ty * m_matrix[1][3] + tz * m_matrix[2][3];
810     return *this;
811 }
812 
translateRight(double tx,double ty)813 TransformationMatrix& TransformationMatrix::translateRight(double tx, double ty)
814 {
815     if (tx != 0) {
816         m_matrix[0][0] +=  m_matrix[0][3] * tx;
817         m_matrix[1][0] +=  m_matrix[1][3] * tx;
818         m_matrix[2][0] +=  m_matrix[2][3] * tx;
819         m_matrix[3][0] +=  m_matrix[3][3] * tx;
820     }
821 
822     if (ty != 0) {
823         m_matrix[0][1] +=  m_matrix[0][3] * ty;
824         m_matrix[1][1] +=  m_matrix[1][3] * ty;
825         m_matrix[2][1] +=  m_matrix[2][3] * ty;
826         m_matrix[3][1] +=  m_matrix[3][3] * ty;
827     }
828 
829     return *this;
830 }
831 
translateRight3d(double tx,double ty,double tz)832 TransformationMatrix& TransformationMatrix::translateRight3d(double tx, double ty, double tz)
833 {
834     translateRight(tx, ty);
835     if (tz != 0) {
836         m_matrix[0][2] +=  m_matrix[0][3] * tz;
837         m_matrix[1][2] +=  m_matrix[1][3] * tz;
838         m_matrix[2][2] +=  m_matrix[2][3] * tz;
839         m_matrix[3][2] +=  m_matrix[3][3] * tz;
840     }
841 
842     return *this;
843 }
844 
skew(double sx,double sy)845 TransformationMatrix& TransformationMatrix::skew(double sx, double sy)
846 {
847     // angles are in degrees. Switch to radians
848     sx = deg2rad(sx);
849     sy = deg2rad(sy);
850 
851     TransformationMatrix mat;
852     mat.m_matrix[0][1] = tan(sy); // note that the y shear goes in the first row
853     mat.m_matrix[1][0] = tan(sx); // and the x shear in the second row
854 
855     multLeft(mat);
856     return *this;
857 }
858 
applyPerspective(double p)859 TransformationMatrix& TransformationMatrix::applyPerspective(double p)
860 {
861     TransformationMatrix mat;
862     if (p != 0)
863         mat.m_matrix[2][3] = -1/p;
864 
865     multLeft(mat);
866     return *this;
867 }
868 
rectToRect(const FloatRect & from,const FloatRect & to)869 TransformationMatrix TransformationMatrix::rectToRect(const FloatRect& from, const FloatRect& to)
870 {
871     ASSERT(!from.isEmpty());
872     return TransformationMatrix(to.width() / from.width(),
873                                 0, 0,
874                                 to.height() / from.height(),
875                                 to.x() - from.x(),
876                                 to.y() - from.y());
877 }
878 
879 //
880 // *this = mat * *this
881 //
multLeft(const TransformationMatrix & mat)882 TransformationMatrix& TransformationMatrix::multLeft(const TransformationMatrix& mat)
883 {
884     Matrix4 tmp;
885 
886     tmp[0][0] = (mat.m_matrix[0][0] * m_matrix[0][0] + mat.m_matrix[0][1] * m_matrix[1][0]
887                + mat.m_matrix[0][2] * m_matrix[2][0] + mat.m_matrix[0][3] * m_matrix[3][0]);
888     tmp[0][1] = (mat.m_matrix[0][0] * m_matrix[0][1] + mat.m_matrix[0][1] * m_matrix[1][1]
889                + mat.m_matrix[0][2] * m_matrix[2][1] + mat.m_matrix[0][3] * m_matrix[3][1]);
890     tmp[0][2] = (mat.m_matrix[0][0] * m_matrix[0][2] + mat.m_matrix[0][1] * m_matrix[1][2]
891                + mat.m_matrix[0][2] * m_matrix[2][2] + mat.m_matrix[0][3] * m_matrix[3][2]);
892     tmp[0][3] = (mat.m_matrix[0][0] * m_matrix[0][3] + mat.m_matrix[0][1] * m_matrix[1][3]
893                + mat.m_matrix[0][2] * m_matrix[2][3] + mat.m_matrix[0][3] * m_matrix[3][3]);
894 
895     tmp[1][0] = (mat.m_matrix[1][0] * m_matrix[0][0] + mat.m_matrix[1][1] * m_matrix[1][0]
896                + mat.m_matrix[1][2] * m_matrix[2][0] + mat.m_matrix[1][3] * m_matrix[3][0]);
897     tmp[1][1] = (mat.m_matrix[1][0] * m_matrix[0][1] + mat.m_matrix[1][1] * m_matrix[1][1]
898                + mat.m_matrix[1][2] * m_matrix[2][1] + mat.m_matrix[1][3] * m_matrix[3][1]);
899     tmp[1][2] = (mat.m_matrix[1][0] * m_matrix[0][2] + mat.m_matrix[1][1] * m_matrix[1][2]
900                + mat.m_matrix[1][2] * m_matrix[2][2] + mat.m_matrix[1][3] * m_matrix[3][2]);
901     tmp[1][3] = (mat.m_matrix[1][0] * m_matrix[0][3] + mat.m_matrix[1][1] * m_matrix[1][3]
902                + mat.m_matrix[1][2] * m_matrix[2][3] + mat.m_matrix[1][3] * m_matrix[3][3]);
903 
904     tmp[2][0] = (mat.m_matrix[2][0] * m_matrix[0][0] + mat.m_matrix[2][1] * m_matrix[1][0]
905                + mat.m_matrix[2][2] * m_matrix[2][0] + mat.m_matrix[2][3] * m_matrix[3][0]);
906     tmp[2][1] = (mat.m_matrix[2][0] * m_matrix[0][1] + mat.m_matrix[2][1] * m_matrix[1][1]
907                + mat.m_matrix[2][2] * m_matrix[2][1] + mat.m_matrix[2][3] * m_matrix[3][1]);
908     tmp[2][2] = (mat.m_matrix[2][0] * m_matrix[0][2] + mat.m_matrix[2][1] * m_matrix[1][2]
909                + mat.m_matrix[2][2] * m_matrix[2][2] + mat.m_matrix[2][3] * m_matrix[3][2]);
910     tmp[2][3] = (mat.m_matrix[2][0] * m_matrix[0][3] + mat.m_matrix[2][1] * m_matrix[1][3]
911                + mat.m_matrix[2][2] * m_matrix[2][3] + mat.m_matrix[2][3] * m_matrix[3][3]);
912 
913     tmp[3][0] = (mat.m_matrix[3][0] * m_matrix[0][0] + mat.m_matrix[3][1] * m_matrix[1][0]
914                + mat.m_matrix[3][2] * m_matrix[2][0] + mat.m_matrix[3][3] * m_matrix[3][0]);
915     tmp[3][1] = (mat.m_matrix[3][0] * m_matrix[0][1] + mat.m_matrix[3][1] * m_matrix[1][1]
916                + mat.m_matrix[3][2] * m_matrix[2][1] + mat.m_matrix[3][3] * m_matrix[3][1]);
917     tmp[3][2] = (mat.m_matrix[3][0] * m_matrix[0][2] + mat.m_matrix[3][1] * m_matrix[1][2]
918                + mat.m_matrix[3][2] * m_matrix[2][2] + mat.m_matrix[3][3] * m_matrix[3][2]);
919     tmp[3][3] = (mat.m_matrix[3][0] * m_matrix[0][3] + mat.m_matrix[3][1] * m_matrix[1][3]
920                + mat.m_matrix[3][2] * m_matrix[2][3] + mat.m_matrix[3][3] * m_matrix[3][3]);
921 
922     setMatrix(tmp);
923     return *this;
924 }
925 
multVecMatrix(double x,double y,double & resultX,double & resultY) const926 void TransformationMatrix::multVecMatrix(double x, double y, double& resultX, double& resultY) const
927 {
928     resultX = m_matrix[3][0] + x * m_matrix[0][0] + y * m_matrix[1][0];
929     resultY = m_matrix[3][1] + x * m_matrix[0][1] + y * m_matrix[1][1];
930     double w = m_matrix[3][3] + x * m_matrix[0][3] + y * m_matrix[1][3];
931     if (w != 1 && w != 0) {
932         resultX /= w;
933         resultY /= w;
934     }
935 }
936 
multVecMatrix(double x,double y,double z,double & resultX,double & resultY,double & resultZ) const937 void TransformationMatrix::multVecMatrix(double x, double y, double z, double& resultX, double& resultY, double& resultZ) const
938 {
939     resultX = m_matrix[3][0] + x * m_matrix[0][0] + y * m_matrix[1][0] + z * m_matrix[2][0];
940     resultY = m_matrix[3][1] + x * m_matrix[0][1] + y * m_matrix[1][1] + z * m_matrix[2][1];
941     resultZ = m_matrix[3][2] + x * m_matrix[0][2] + y * m_matrix[1][2] + z * m_matrix[2][2];
942     double w = m_matrix[3][3] + x * m_matrix[0][3] + y * m_matrix[1][3] + z * m_matrix[2][3];
943     if (w != 1 && w != 0) {
944         resultX /= w;
945         resultY /= w;
946         resultZ /= w;
947     }
948 }
949 
isInvertible() const950 bool TransformationMatrix::isInvertible() const
951 {
952     if (isIdentityOrTranslation())
953         return true;
954 
955     double det = WebCore::determinant4x4(m_matrix);
956 
957     if (fabs(det) < SMALL_NUMBER)
958         return false;
959 
960     return true;
961 }
962 
inverse() const963 TransformationMatrix TransformationMatrix::inverse() const
964 {
965     if (isIdentityOrTranslation()) {
966         // identity matrix
967         if (m_matrix[3][0] == 0 && m_matrix[3][1] == 0 && m_matrix[3][2] == 0)
968             return TransformationMatrix();
969 
970         // translation
971         return TransformationMatrix(1, 0, 0, 0,
972                                     0, 1, 0, 0,
973                                     0, 0, 1, 0,
974                                     -m_matrix[3][0], -m_matrix[3][1], -m_matrix[3][2], 1);
975     }
976 
977     TransformationMatrix invMat;
978     bool inverted = WebCore::inverse(m_matrix, invMat.m_matrix);
979     if (!inverted)
980         return TransformationMatrix();
981 
982     return invMat;
983 }
984 
makeAffine()985 void TransformationMatrix::makeAffine()
986 {
987     m_matrix[0][2] = 0;
988     m_matrix[0][3] = 0;
989 
990     m_matrix[1][2] = 0;
991     m_matrix[1][3] = 0;
992 
993     m_matrix[2][0] = 0;
994     m_matrix[2][1] = 0;
995     m_matrix[2][2] = 1;
996     m_matrix[2][3] = 0;
997 
998     m_matrix[3][2] = 0;
999     m_matrix[3][3] = 1;
1000 }
1001 
blendFloat(double & from,double to,double progress)1002 static inline void blendFloat(double& from, double to, double progress)
1003 {
1004     if (from != to)
1005         from = from + (to - from) * progress;
1006 }
1007 
blend(const TransformationMatrix & from,double progress)1008 void TransformationMatrix::blend(const TransformationMatrix& from, double progress)
1009 {
1010     if (from.isIdentity() && isIdentity())
1011         return;
1012 
1013     // decompose
1014     DecomposedType fromDecomp;
1015     DecomposedType toDecomp;
1016     from.decompose(fromDecomp);
1017     decompose(toDecomp);
1018 
1019     // interpolate
1020     blendFloat(fromDecomp.scaleX, toDecomp.scaleX, progress);
1021     blendFloat(fromDecomp.scaleY, toDecomp.scaleY, progress);
1022     blendFloat(fromDecomp.scaleZ, toDecomp.scaleZ, progress);
1023     blendFloat(fromDecomp.skewXY, toDecomp.skewXY, progress);
1024     blendFloat(fromDecomp.skewXZ, toDecomp.skewXZ, progress);
1025     blendFloat(fromDecomp.skewYZ, toDecomp.skewYZ, progress);
1026     blendFloat(fromDecomp.translateX, toDecomp.translateX, progress);
1027     blendFloat(fromDecomp.translateY, toDecomp.translateY, progress);
1028     blendFloat(fromDecomp.translateZ, toDecomp.translateZ, progress);
1029     blendFloat(fromDecomp.perspectiveX, toDecomp.perspectiveX, progress);
1030     blendFloat(fromDecomp.perspectiveY, toDecomp.perspectiveY, progress);
1031     blendFloat(fromDecomp.perspectiveZ, toDecomp.perspectiveZ, progress);
1032     blendFloat(fromDecomp.perspectiveW, toDecomp.perspectiveW, progress);
1033 
1034     slerp(&fromDecomp.quaternionX, &toDecomp.quaternionX, progress);
1035 
1036     // recompose
1037     recompose(fromDecomp);
1038 }
1039 
decompose(DecomposedType & decomp) const1040 bool TransformationMatrix::decompose(DecomposedType& decomp) const
1041 {
1042     if (isIdentity()) {
1043         memset(&decomp, 0, sizeof(decomp));
1044         decomp.perspectiveW = 1;
1045         decomp.scaleX = 1;
1046         decomp.scaleY = 1;
1047         decomp.scaleZ = 1;
1048     }
1049 
1050     if (!WebCore::decompose(m_matrix, decomp))
1051         return false;
1052     return true;
1053 }
1054 
recompose(const DecomposedType & decomp)1055 void TransformationMatrix::recompose(const DecomposedType& decomp)
1056 {
1057     makeIdentity();
1058 
1059     // first apply perspective
1060     m_matrix[0][3] = (float) decomp.perspectiveX;
1061     m_matrix[1][3] = (float) decomp.perspectiveY;
1062     m_matrix[2][3] = (float) decomp.perspectiveZ;
1063     m_matrix[3][3] = (float) decomp.perspectiveW;
1064 
1065     // now translate
1066     translate3d((float) decomp.translateX, (float) decomp.translateY, (float) decomp.translateZ);
1067 
1068     // apply rotation
1069     double xx = decomp.quaternionX * decomp.quaternionX;
1070     double xy = decomp.quaternionX * decomp.quaternionY;
1071     double xz = decomp.quaternionX * decomp.quaternionZ;
1072     double xw = decomp.quaternionX * decomp.quaternionW;
1073     double yy = decomp.quaternionY * decomp.quaternionY;
1074     double yz = decomp.quaternionY * decomp.quaternionZ;
1075     double yw = decomp.quaternionY * decomp.quaternionW;
1076     double zz = decomp.quaternionZ * decomp.quaternionZ;
1077     double zw = decomp.quaternionZ * decomp.quaternionW;
1078 
1079     // Construct a composite rotation matrix from the quaternion values
1080     TransformationMatrix rotationMatrix(1 - 2 * (yy + zz), 2 * (xy - zw), 2 * (xz + yw), 0,
1081                            2 * (xy + zw), 1 - 2 * (xx + zz), 2 * (yz - xw), 0,
1082                            2 * (xz - yw), 2 * (yz + xw), 1 - 2 * (xx + yy), 0,
1083                            0, 0, 0, 1);
1084 
1085     multLeft(rotationMatrix);
1086 
1087     // now apply skew
1088     if (decomp.skewYZ) {
1089         TransformationMatrix tmp;
1090         tmp.setM32((float) decomp.skewYZ);
1091         multLeft(tmp);
1092     }
1093 
1094     if (decomp.skewXZ) {
1095         TransformationMatrix tmp;
1096         tmp.setM31((float) decomp.skewXZ);
1097         multLeft(tmp);
1098     }
1099 
1100     if (decomp.skewXY) {
1101         TransformationMatrix tmp;
1102         tmp.setM21((float) decomp.skewXY);
1103         multLeft(tmp);
1104     }
1105 
1106     // finally, apply scale
1107     scale3d((float) decomp.scaleX, (float) decomp.scaleY, (float) decomp.scaleZ);
1108 }
1109 
1110 }
1111