aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--src/hrs-scaling.c4
-rw-r--r--src/partialator.c3
-rw-r--r--src/post-refinement.c16
3 files changed, 17 insertions, 6 deletions
diff --git a/src/hrs-scaling.c b/src/hrs-scaling.c
index fcf7efae..5c5ec8f0 100644
--- a/src/hrs-scaling.c
+++ b/src/hrs-scaling.c
@@ -263,6 +263,10 @@ static double iterate_scale(struct image *images, int n,
e_val = gsl_vector_alloc(n);
e_vec = gsl_matrix_alloc(n, n);
val = gsl_eigen_symmv(M, e_val, e_vec, work);
+ if ( val ) {
+ ERROR("Couldn't diagonalise matrix.\n");
+ return 0.0;
+ }
gsl_eigen_symmv_free(work);
val = gsl_eigen_symmv_sort(e_val, e_vec, GSL_EIGEN_SORT_ABS_DESC);
diff --git a/src/partialator.c b/src/partialator.c
index dc29d0d6..a410c507 100644
--- a/src/partialator.c
+++ b/src/partialator.c
@@ -22,6 +22,7 @@
#include <getopt.h>
#include <assert.h>
#include <pthread.h>
+#include <gsl/gsl_errno.h>
#include "utils.h"
#include "hdf5-file.h"
@@ -274,6 +275,8 @@ int main(int argc, char *argv[])
}
STATUS("There are %i patterns to process\n", n_total_patterns);
+ gsl_set_error_handler_off();
+
images = malloc(n_total_patterns * sizeof(struct image));
if ( images == NULL ) {
ERROR("Couldn't allocate memory for images.\n");
diff --git a/src/post-refinement.c b/src/post-refinement.c
index db8d4990..0945affd 100644
--- a/src/post-refinement.c
+++ b/src/post-refinement.c
@@ -320,13 +320,17 @@ static double pr_iterate(struct image *image, const RefList *full,
//show_matrix_eqn(M, v, NUM_PARAMS);
shifts = gsl_vector_alloc(NUM_PARAMS);
- gsl_linalg_HH_solve(M, v, shifts);
-
max_shift = 0.0;
- for ( param=0; param<NUM_PARAMS; param++ ) {
- double shift = gsl_vector_get(shifts, param);
- apply_shift(image, param, shift);
- if ( fabs(shift) > max_shift ) max_shift = fabs(shift);
+ if ( gsl_linalg_HH_solve(M, v, shifts) ) {
+ ERROR("Couldn't solve normal equations!\n");
+ } else {
+
+ for ( param=0; param<NUM_PARAMS; param++ ) {
+ double shift = gsl_vector_get(shifts, param);
+ apply_shift(image, param, shift);
+ if ( fabs(shift) > max_shift ) max_shift = fabs(shift);
+ }
+
}
gsl_matrix_free(M);