add mat3 stuff

src/matrix-math.c

@@ -213,6 +213,156 @@ void mat4BuildPerspective(mat4 out, GLfloat fovy, GLfloat aspect, GLfloat n, GLf
 	mat4BuildProjection(out, r, l, t, b, n, f);
 }
 
+/**
+ * linearly interpolates between a and b, storing the result in out
+ */
+void mat4Interpolate(mat4 out, mat4 a, mat4 b, GLfloat f) {
+	for (int i = 0; i < 16; i++) {
+		out[i] = (1 - f) * a[i] + f * b[i];
+	}
+}
+
+void mat4From3(mat4 out, mat3 in) {
+	mat4Identity(out);
+	memcpy(&out[0], &in[0], sizeof(GLfloat) * 3);
+	memcpy(&out[4], &in[3], sizeof(GLfloat) * 3);
+	memcpy(&out[8], &in[6], sizeof(GLfloat) * 3);
+}
+
+/**
+ * https://stackoverflow.com/questions/1148309/inverting-a-4x4-matrix
+ */
+void mat4Inverse(mat4 out, mat4 m) {
+	GLfloat inv[16], det;
+	int i;
+
+	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];
+
+	det = m[0] * inv[0] + m[1] * inv[4] + m[2] * inv[8] + m[3] * inv[12];
+
+	if (det == 0) {
+		mat4Copy(m, out); // singular matrix, can't invert
+		return;
+	}
+
+	det = 1.0 / det;
+
+	for (i = 0; i < 16; i++)
+		out[i] = inv[i] * det;
+
+	return;
+}
+
 
 /**
  * subtracts b from a, storing the result in out
@@ -269,3 +419,87 @@ void vec3Set(vec3 out, GLfloat x, GLfloat y, GLfloat z) {
 	out[1] = y;
 	out[2] = z;
 }
+
+
+void mat3Copy(mat3 src, mat3 dst) {
+	for (int i = 0; i < 9; i++) {
+		dst[i] = src[i];
+	}
+}
+
+void mat3From4(mat3 out, mat4 in) {
+	memcpy(&out[0], &in[0], sizeof(GLfloat) * 3);
+	memcpy(&out[3], &in[4], sizeof(GLfloat) * 3);
+	memcpy(&out[6], &in[8], sizeof(GLfloat) * 3);
+}
+
+void mat3Transpose(mat3 out, mat3 in) {
+	mat3 result;
+
+	for (int i = 0; i < 3; i++) {
+		for (int j = 0; j < 3; j++) {
+			result[i * 3 + j] = in[j * 3 + i];
+		}
+	}
+
+	mat3Copy(result, out);
+}
+
+void mat3Minor(mat3 out, mat3 in) {
+	mat3 result;
+
+	// TODO: check if this is correct
+	result[0] = in[4] * in[8] - in[5] * in[7];
+	result[1] = in[3] * in[8] - in[5] * in[6];
+	result[2] = in[3] * in[7] - in[4] * in[6];
+
+	result[3] = in[1] * in[8] - in[2] * in[7];
+	result[4] = in[0] * in[8] - in[2] * in[6];
+	result[5] = in[0] * in[7] - in[1] * in[6];
+
+
+	result[6] = in[1] * in[5] - in[2] * in[4];
+	result[7] = in[0] * in[5] - in[2] * in[3];
+	result[8] = in[0] * in[4] - in[1] * in[3];
+
+	mat3Copy(result, out);
+}
+
+void mat3Cofactor(mat3 out, mat3 in) {
+	mat3Minor(out, in);
+
+	out[1] *= -1;
+	out[3] *= -1;
+	out[5] *= -1;
+	out[7] *= -1;
+}
+
+void mat3Adjoint(mat3 out, mat3 in) {
+	mat3Cofactor(out, in);
+	mat3Transpose(out, in);
+}
+
+void mat3MultiplyScalar(mat3 out, mat3 in, GLfloat x) {
+	for (int i = 0; i < 9; i++) {
+		out[i] = in[i] * x;
+	}
+}
+
+GLfloat mat3Determinant(mat3 M) {
+	return
+		+ M[0] * M[4] * M[8]
+		+ M[3] * M[7] * M[2]
+		+ M[6] * M[1] * M[5]
+
+		- M[2] * M[4] * M[6]
+		- M[5] * M[7] * M[0]
+		- M[8] * M[1] * M[3]
+	;
+}
+
+void mat3Inverse(mat3 out, mat3 in) {
+	mat3 result;
+	mat3Adjoint(result, in);
+	mat3MultiplyScalar(result, result, 1 / mat3Determinant(in));
+	mat3Copy(result, out);
+}

src/matrix-math.h

@@ -11,6 +11,7 @@ typedef GLfloat vec4[4];
 typedef GLfloat vec3[3];
 
 typedef GLfloat mat4[16];
+typedef GLfloat mat3[9];
 
 extern void mat4Identity(mat4 mat);
 extern void mat4Copy(mat4 src, mat4 dst);
@@ -25,6 +26,9 @@ extern void mat4RotateZ(mat4 out, mat4 in, GLfloat a);
 extern void mat4BuildLookAt(mat4 out, vec3 eye, vec3 center, vec3 up);
 extern void mat4BuildProjection(mat4 out, GLfloat r, GLfloat l, GLfloat t, GLfloat b, GLfloat n, GLfloat f);
 extern void mat4BuildPerspective(mat4 out, GLfloat fovy, GLfloat aspect, GLfloat n, GLfloat f);
+extern void mat4Interpolate(mat4 out, mat4 a, mat4 b, GLfloat f);
+extern void mat4From3(mat4 out, mat3 in);
+extern void mat4Inverse(mat4 m, mat4 out);
 
 extern void vec3Subtract(vec3 out, vec3 a, vec3 b);
 extern void vec3CrossProduct(vec3 out, vec3 a, vec3 b);
@@ -33,4 +37,14 @@ extern GLfloat vec3Length(vec3 in);
 extern GLfloat vec3DotProduct(vec3 a, vec3 b);
 extern void vec3Set(vec3 out, GLfloat x, GLfloat y, GLfloat z);
 
+extern void mat3Copy(mat3 src, mat3 dst);
+extern void mat3Inverse(mat3 out, mat3 in);
+extern GLfloat mat3Determinant(mat3 m);
+extern void mat3MultiplyScalar(mat3 out, mat3 in, GLfloat x);
+extern void mat3Adjoint(mat3 out, mat3 in);
+extern void mat3Cofactor(mat3 out, mat3 in);
+extern void mat3Minor(mat3 out, mat3 in);
+extern void mat3Transpose(mat3 out, mat3 in);
+extern void mat3From4(mat3 out, mat4 in);
+
 #endif