commit 6959d45972e2034b8c22d86a22dd232f549cd578
author Mike Klein <mtklein@chromium.org> Thu Apr 05 15:24:00 2018 -0400
committer Skia Commit-Bot <skia-commit-bot@chromium.org> Thu Apr 05 19:59:38 2018 +0000
tree ae95914db21f45fc1d8948ef7fede9e7090c10e6
parent d75edb2e3bfebe530ffb4388e0d51d57260377b0

better matrix 3x3 invert

I didn't think hard about this... just folded through constants.

Change-Id: I994773009964e183691cf5f10861756acfe4b2ae
Reviewed-on: https://skia-review.googlesource.com/118986
Commit-Queue: Mike Klein <mtklein@chromium.org>
Reviewed-by: Brian Osman <brianosman@google.com>

src/LinearAlgebra.c

diff --git a/src/LinearAlgebra.c b/src/LinearAlgebra.c
index a872613..e7e335e 100644
--- a/src/LinearAlgebra.c
+++ b/src/LinearAlgebra.c
@@ -9,16 +9,6 @@
 #include "LinearAlgebra.h"
 #include "PortableMath.h"
 
-static bool Matrix4x4_isfinite(const skcms_Matrix4x4* m) {
-    for (int r = 0; r < 4; ++r)
-    for (int c = 0; c < 4; ++c) {
-        if (!isfinitef_(m->vals[r][c])) {
-            return false;
-        }
-    }
-    return true;
-}
-
 bool skcms_Matrix4x4_invert(const skcms_Matrix4x4* src, skcms_Matrix4x4* dst) {
     double a00 = (double)src->vals[0][0],
            a01 = (double)src->vals[1][0],
@@ -96,28 +86,69 @@
     dst->vals[2][3] = (float)( a31*b01 - a30*b03 - a32*b00 );
     dst->vals[3][3] = (float)( a20*b03 - a21*b01 + a22*b00 );
 
-    return Matrix4x4_isfinite(dst);
+    for (int r = 0; r < 4; ++r)
+    for (int c = 0; c < 4; ++c) {
+        if (!isfinitef_(dst->vals[r][c])) {
+            return false;
+        }
+    }
+    return true;
 }
 
 bool skcms_Matrix3x3_invert(const skcms_Matrix3x3* src, skcms_Matrix3x3* dst) {
-    // TODO: I am being _very_ lazy.
-    skcms_Matrix4x4 m = {{
-        { src->vals[0][0], src->vals[0][1], src->vals[0][2], 0.0f },
-        { src->vals[1][0], src->vals[1][1], src->vals[1][2], 0.0f },
-        { src->vals[2][0], src->vals[2][1], src->vals[2][2], 0.0f },
-        {            0.0f,            0.0f,            0.0f, 1.0f },
-    }};
+    double a00 = (double)src->vals[0][0],
+           a01 = (double)src->vals[1][0],
+           a02 = (double)src->vals[2][0],
+           a10 = (double)src->vals[0][1],
+           a11 = (double)src->vals[1][1],
+           a12 = (double)src->vals[2][1],
+           a20 = (double)src->vals[0][2],
+           a21 = (double)src->vals[1][2],
+           a22 = (double)src->vals[2][2];
 
-    skcms_Matrix4x4 inv;
-    if (!skcms_Matrix4x4_invert(&m, &inv)) {
+    double b0 = a00*a11 - a01*a10,
+           b1 = a00*a12 - a02*a10,
+           b2 = a01*a12 - a02*a11,
+           b3 = a20,
+           b4 = a21,
+           b5 = a22;
+
+    double determinant = b0*b5
+                       - b1*b4
+                       + b2*b3;
+
+    if (determinant == 0) {
         return false;
     }
 
-    *dst = (skcms_Matrix3x3){{
-        { inv.vals[0][0], inv.vals[0][1], inv.vals[0][2] },
-        { inv.vals[1][0], inv.vals[1][1], inv.vals[1][2] },
-        { inv.vals[2][0], inv.vals[2][1], inv.vals[2][2] },
-    }};
+    double invdet = 1.0 / determinant;
+    if (!isfinitef_((float)invdet)) {
+        return false;
+    }
+
+    b0 *= invdet;
+    b1 *= invdet;
+    b2 *= invdet;
+    b3 *= invdet;
+    b4 *= invdet;
+    b5 *= invdet;
+
+    dst->vals[0][0] = (float)( a11*b5 - a12*b4 );
+    dst->vals[1][0] = (float)( a02*b4 - a01*b5 );
+    dst->vals[2][0] = (float)(        +     b2 );
+    dst->vals[0][1] = (float)( a12*b3 - a10*b5 );
+    dst->vals[1][1] = (float)( a00*b5 - a02*b3 );
+    dst->vals[2][1] = (float)(        -     b1 );
+    dst->vals[0][2] = (float)( a10*b4 - a11*b3 );
+    dst->vals[1][2] = (float)( a01*b3 - a00*b4 );
+    dst->vals[2][2] = (float)(        +     b0 );
+
+    for (int r = 0; r < 3; ++r)
+    for (int c = 0; c < 3; ++c) {
+        if (!isfinitef_(dst->vals[r][c])) {
+            return false;
+        }
+    }
     return true;
 }