Keep track of the "un-centering" matrix, as well as the "centering"

This makes it easy to reverse the transformation, if required, which it
is when comparing centered cells.

1 files changed, 25 insertions, 0 deletions

diff --git a/libcrystfel/src/rational.c b/libcrystfel/src/rational.c
index 58178559..4f02c8b1 100644
--- a/libcrystfel/src/rational.c
+++ b/libcrystfel/src/rational.c
@@ -546,6 +546,31 @@ void rtnl_mtx_print(const RationalMatrix *m)
 }
 
 
+void rtnl_mtx_mtxmult(const RationalMatrix *A, const RationalMatrix *B,
+                      RationalMatrix *ans)
+{
+	int i, j;
+
+	assert(ans->cols == A->cols);
+	assert(ans->rows == B->rows);
+	assert(A->cols == B->rows);
+
+	for ( i=0; i<ans->rows; i++ ) {
+		for ( j=0; j<ans->cols; j++ ) {
+			int k;
+			Rational sum = rtnl_zero();
+			for ( k=0; k<A->rows; k++ ) {
+				Rational add;
+				add = rtnl_mul(rtnl_mtx_get(A, i, k),
+				               rtnl_mtx_get(B, k, j));
+				sum = rtnl_add(sum, add);
+			}
+			rtnl_mtx_set(ans, i, j, sum);
+		}
+	}
+}
+
+
 void rtnl_mtx_mult(const RationalMatrix *m, const Rational *vec, Rational *ans)
 {
 	int i, j;