diff options
| -rw-r--r-- | libcrystfel/src/integration.c | 2 | ||||
| -rw-r--r-- | libcrystfel/src/utils.c | 21 | ||||
| -rw-r--r-- | libcrystfel/src/utils.h | 2 | ||||
| -rw-r--r-- | 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++ ) { |