From ecfe23ca70540e2f8ee6efb600c7b577d95f236f Mon Sep 17 00:00:00 2001
From: Thomas White <taw@physics.org>
Date: Fri, 9 Aug 2013 17:52:36 +0200
Subject: Improve show_matrix_eqn()

---
 libcrystfel/src/integration.c |  2 +-
 libcrystfel/src/utils.c       | 21 ++++++++++++++-------
 libcrystfel/src/utils.h       |  2 +-
 src/post-refinement.c         |  8 +++++---
 4 files changed, 21 insertions(+), 12 deletions(-)

diff --git a/libcrystfel/src/integration.c b/libcrystfel/src/integration.c
index 0b82d9e5..2816784f 100644
--- a/libcrystfel/src/integration.c
+++ b/libcrystfel/src/integration.c
@@ -411,7 +411,7 @@ static void fit_bg(struct intcontext *ic, struct peak_box *bx)
 	}
 
 	if ( bx->verbose ) {
-		show_matrix_eqn(bx->bgm, v, 3);
+		show_matrix_eqn(bx->bgm, v);
 	}
 
 	/* SVD is massive overkill here */
diff --git a/libcrystfel/src/utils.c b/libcrystfel/src/utils.c
index 697c773f..6c731a83 100644
--- a/libcrystfel/src/utils.c
+++ b/libcrystfel/src/utils.c
@@ -57,22 +57,29 @@
  * show_matrix_eqn:
  * @M: A matrix
  * @v: A vector
- * @r: The number of elements in @v and the side length of @M.
  *
  * Displays a matrix equation of the form @M.a = @v.
- *
- * @M must be square.
  **/
-void show_matrix_eqn(gsl_matrix *M, gsl_vector *v, int r)
+void show_matrix_eqn(gsl_matrix *M, gsl_vector *v)
 {
 	int i, j;
 
-	for ( i=0; i<r; i++ ) {
+	if ( M->size1 != v->size ) {
+		ERROR("Matrix and vector sizes don't agree.\n");
+		return;
+	}
+
+	for ( i=0; i<M->size1; i++ ) {
 		STATUS("[ ");
-		for ( j=0; j<r; j++ ) {
+		for ( j=0; j<M->size2; j++ ) {
 			STATUS("%+9.3e ", gsl_matrix_get(M, i, j));
 		}
-		STATUS("][ a%2i ] = [ %+9.3e ]\n", i, gsl_vector_get(v, i));
+		if ( i < M->size2 ) {
+			STATUS("][ a%2i ] = [ %+9.3e ]\n", i,
+			       gsl_vector_get(v, i));
+		} else {
+			STATUS("]        = [ +%9.3e ]\n", gsl_vector_get(v, i));
+		}
 	}
 }
 
diff --git a/libcrystfel/src/utils.h b/libcrystfel/src/utils.h
index a206ccfd..f4477717 100644
--- a/libcrystfel/src/utils.h
+++ b/libcrystfel/src/utils.h
@@ -96,7 +96,7 @@ extern struct rvec quat_rot(struct rvec q, struct quaternion z);
 
 /* --------------------------- Useful functions ----------------------------- */
 
-extern void show_matrix_eqn(gsl_matrix *M, gsl_vector *v, int r);
+extern void show_matrix_eqn(gsl_matrix *M, gsl_vector *v);
 extern void show_matrix(gsl_matrix *M);
 extern size_t notrail(char *s);
 extern void chomp(char *s);
diff --git a/src/post-refinement.c b/src/post-refinement.c
index 628e33d2..94c8a816 100644
--- a/src/post-refinement.c
+++ b/src/post-refinement.c
@@ -361,7 +361,8 @@ static int check_eigen(gsl_vector *e_val)
 }
 
 
-static gsl_vector *solve_svd(gsl_vector *v, gsl_matrix *M, int *n_filt)
+static gsl_vector *solve_svd(gsl_vector *v, gsl_matrix *M, int *n_filt,
+                             int verbose)
 {
 	gsl_matrix *s_vec;
 	gsl_vector *s_val;
@@ -452,6 +453,7 @@ static double pr_iterate(Crystal *cr, const RefList *full,
 	RefList *reflections;
 	double max_shift;
 	int nref = 0;
+	const int verbose = 1;
 
 	reflections = crystal_get_reflections(cr);
 
@@ -533,7 +535,7 @@ static double pr_iterate(Crystal *cr, const RefList *full,
 
 		nref++;
 	}
-	//show_matrix_eqn(M, v, NUM_PARAMS);
+	if ( verbose ) show_matrix_eqn(M, v);
 
 	//STATUS("%i reflections went into the equations.\n", nref);
 	if ( nref == 0 ) {
@@ -544,7 +546,7 @@ static double pr_iterate(Crystal *cr, const RefList *full,
 	}
 
 	max_shift = 0.0;
-	shifts = solve_svd(v, M, &prdata->n_filtered);
+	shifts = solve_svd(v, M, &prdata->n_filtered, verbose);
 	if ( shifts != NULL ) {
 
 		for ( param=0; param<NUM_PARAMS; param++ ) {
-- 
cgit v1.2.3