[Piglit] [PATCH 02/23] util: Add matrix inverse and transpose functions.

[Fabian Bieler]

---
 tests/util/piglit-matrix.c | 90 ++++++++++++++++++++++++++++++++++++++++++++++
 tests/util/piglit-matrix.h |  5 +++
 2 files changed, 95 insertions(+)

diff --git a/tests/util/piglit-matrix.c b/tests/util/piglit-matrix.c
index ce90c0a..3c295e6 100644
--- a/tests/util/piglit-matrix.c
+++ b/tests/util/piglit-matrix.c
@@ -22,6 +22,7 @@
  */
 
 
+#include <assert.h>
 #include <math.h>
 #include <stdio.h>
 #include "piglit-matrix.h"
@@ -451,6 +452,95 @@ piglit_print_matrix(const float mat[16])
 }
 
 
+/**
+ * Invert matrix using cramer's rule.
+ * This assumes that the matrix is non-singular (or non-near-singular).
+ */
+void
+piglit_matrix_inverse(float inv[16], const float m[16])
+{
+	inv[0] = m[5] * m[10] * m[15] - m[5] * m[11] * m[14] -
+		 m[9] * m[6] * m[15] + m[9] * m[7] * m[14] +
+		 m[13] * m[6] * m[11] - m[13] * m[7] * m[10];
+
+	inv[4] = -m[4] * m[10] * m[15] + m[4] * m[11] * m[14] +
+		 m[8] * m[6] * m[15] - m[8] * m[7] * m[14] -
+		 m[12] * m[6] * m[11] + m[12] * m[7] * m[10];
+
+	inv[8] = m[4] * m[9] * m[15] - m[4] * m[11] * m[13] -
+		 m[8] * m[5] * m[15] + m[8] * m[7] * m[13] +
+		 m[12] * m[5] * m[11] - m[12] * m[7] * m[9];
+
+	inv[12] = -m[4] * m[9] * m[14] + m[4] * m[10] * m[13] +
+		  m[8] * m[5] * m[14] - m[8] * m[6] * m[13] -
+		  m[12] * m[5] * m[10] + m[12] * m[6] * m[9];
+
+	inv[1] = -m[1] * m[10] * m[15] + m[1] * m[11] * m[14] +
+		 m[9] * m[2] * m[15] - m[9] * m[3] * m[14] -
+		 m[13] * m[2] * m[11] + m[13] * m[3] * m[10];
+
+	inv[5] = m[0] * m[10] * m[15] - m[0] * m[11] * m[14] -
+		 m[8] * m[2] * m[15] + m[8] * m[3] * m[14] +
+		 m[12] * m[2] * m[11] - m[12] * m[3] * m[10];
+
+	inv[9] = -m[0] * m[9] * m[15] + m[0] * m[11] * m[13] +
+		 m[8] * m[1] * m[15] - m[8] * m[3] * m[13] -
+		 m[12] * m[1] * m[11] + m[12] * m[3] * m[9];
+
+	inv[13] = m[0] * m[9] * m[14] - m[0] * m[10] * m[13] -
+		  m[8] * m[1] * m[14] + m[8] * m[2] * m[13] +
+		  m[12] * m[1] * m[10] - m[12] * m[2] * m[9];
+
+	inv[2] = m[1] * m[6] * m[15] - m[1] * m[7] * m[14] -
+		 m[5] * m[2] * m[15] + m[5] * m[3] * m[14] +
+		 m[13] * m[2] * m[7] - m[13] * m[3] * m[6];
+
+	inv[6] = -m[0] * m[6] * m[15] + m[0] * m[7] * m[14] +
+		 m[4] * m[2] * m[15] - m[4] * m[3] * m[14] -
+		 m[12] * m[2] * m[7] + m[12] * m[3] * m[6];
+
+	inv[10] = m[0] * m[5] * m[15] - m[0] * m[7] * m[13] -
+		  m[4] * m[1] * m[15] + m[4] * m[3] * m[13] +
+		  m[12] * m[1] * m[7] - m[12] * m[3] * m[5];
+
+	inv[14] = -m[0] * m[5] * m[14] + m[0] * m[6] * m[13] +
+		  m[4] * m[1] * m[14] - m[4] * m[2] * m[13] -
+		  m[12] * m[1] * m[6] + m[12] * m[2] * m[5];
+
+	inv[3] = -m[1] * m[6] * m[11] + m[1] * m[7] * m[10] +
+		 m[5] * m[2] * m[11] - m[5] * m[3] * m[10] -
+		 m[9] * m[2] * m[7] + m[9] * m[3] * m[6];
+
+	inv[7] = m[0] * m[6] * m[11] - m[0] * m[7] * m[10] -
+		 m[4] * m[2] * m[11] + m[4] * m[3] * m[10] +
+		 m[8] * m[2] * m[7] - m[8] * m[3] * m[6];
+
+	inv[11] = -m[0] * m[5] * m[11] + m[0] * m[7] * m[9] +
+		  m[4] * m[1] * m[11] - m[4] * m[3] * m[9] -
+		  m[8] * m[1] * m[7] + m[8] * m[3] * m[5];
+
+	inv[15] = m[0] * m[5] * m[10] - m[0] * m[6] * m[9] -
+		  m[4] * m[1] * m[10] + m[4] * m[2] * m[9] +
+		  m[8] * m[1] * m[6] - m[8] * m[2] * m[5];
+
+	float det = m[0] * inv[0] + m[1] * inv[4] + m[2] * inv[8] + m[3] * inv[12];
+
+	assert(fabsf(det) > 1e-10f);
+
+	for (int i = 0; i < 16; i++)
+		inv[i] = inv[i] / det;
+}
+
+
+void
+piglit_matrix_transpose(float out[16], const float m[16])
+{
+	for (int i = 0; i < 4; ++i)
+		for (int j = 0; j < 4; ++j)
+			out[i+4*j] = m[4*i+j];
+}
+
+
 /*
  * XXX to do items:
  * We could add a simple set of matrix stack functions.
diff --git a/tests/util/piglit-matrix.h b/tests/util/piglit-matrix.h
index 5a51d91..31c80d3 100644
--- a/tests/util/piglit-matrix.h
+++ b/tests/util/piglit-matrix.h
@@ -83,6 +83,11 @@ piglit_project_to_window(float win[3],
 void
 piglit_print_matrix(const float mat[16]);
 
+void
+piglit_matrix_inverse(float inv[16], const float m[16]);
+
+void
+piglit_matrix_transpose(float out[16], const float m[16]);
 
 #ifdef __cplusplus
 } /* end extern "C" */
-- 
2.7.4