Refinement test framework

git-svn-id: svn://cook.msm.cam.ac.uk:745/diff-tomo/dtr@271 bf6ca9ba-c028-0410-8290-897cf20841d1

1 files changed, 21 insertions, 17 deletions

diff --git a/src/refine.c b/src/refine.c
index b56cc05..7bb2585 100644
--- a/src/refine.c
+++ b/src/refine.c
@@ -93,11 +93,11 @@ void refine_do_sequence(ControlContext *ctx) {
 	gtk_widget_show_all(window_fit);
 	
 	/* Perform final refinement */
-	refine_do_cell(ctx);
 	printf("Best omega offset = %f deg (%f nm^-1)\n", rad2deg(omega_offs_best), fit_best/1e9);
 	for ( j=0; j<ctx->images->n_images; j++ ) {
 		ctx->images->images[j].omega += omega_offs_best;
 	}
+	refine_do_cell(ctx);
 	reproject_lattice_changed(ctx);
 	mapping_adjust_axis(ctx, omega_offs_best);
 
@@ -165,13 +165,11 @@ static void refine_simplex_transform(SimplexVertex *s, int v_worst, double fac)
 	
 }
 
-static double refine_iteration(SimplexVertex *s, Deviation *d, int nf) {
+static double refine_iteration(SimplexVertex *s, Deviation *d, int nf, int debug) {
 
 	int v_worst, v_best, v_second_worst, i;
 	double fom_worst, fom_new, fom_best, fom_second_worst;
-
-//	printf("---\n");
-
+	
 	/* Find the least favourable vertex of the simplex */
 	v_worst = 0;
 	fom_worst = 0.0;
@@ -185,7 +183,7 @@ static double refine_iteration(SimplexVertex *s, Deviation *d, int nf) {
 		
 		fom = refine_mean_dev(d, nf, s, i);
 		
-		//printf("Mean deviation at simplex vertex %i = %f nm^-1\n", i, fom/1e9);
+		if ( debug ) printf("Mean deviation at simplex vertex %i = %f nm^-1\n", i, fom/1e9);
 		if ( fom > fom_worst ) {
 			v_second_worst = v_worst;
 			fom_second_worst = fom_worst;
@@ -198,27 +196,28 @@ static double refine_iteration(SimplexVertex *s, Deviation *d, int nf) {
 		}
 		
 	}
+	if ( debug ) printf("The worst vertex is number %i\n", v_worst);
 	
 	/* Reflect this vertex across the opposite face */
 	refine_simplex_transform(s, v_worst, -1.0);
 	
 	/* Is the worst vertex any better? */
 	fom_new = refine_mean_dev(d, nf, s, v_worst);
-//	printf("New mean deviation for the worst vertex after reflection is %f nm^-1\n", fom_new/1e9);
+	if ( debug ) printf("New mean deviation for the worst vertex after reflection is %f nm^-1\n", fom_new/1e9);
 	if ( fom_new > fom_worst ) {
 	
 		double fom_new_new;
 		
 		/* It's worse than before.  Contract in 1D and see if that helps. */
-//		printf("Worse.  Trying a 1D contraction...\n");
+		if ( debug ) printf("Worse.  Trying a 1D contraction...\n");
 		refine_simplex_transform(s, v_worst, 0.5);
 		fom_new_new = refine_mean_dev(d, nf, s, v_worst);
-//		printf("Mean deviation after 1D contraction is %f nm^-1\n", fom_new_new/1e9);
+		if ( debug ) printf("Mean deviation after 1D contraction is %f nm^-1\n", fom_new_new/1e9);
 		if ( fom_new_new > fom_second_worst ) {
 			
 			int i;
 			
-//			printf("Not as good as the second worst vertex: contracting around the best vertex (%i)\n", v_best);
+			if ( debug ) printf("Not as good as the second worst vertex: contracting around the best vertex (%i)\n", v_best);
 			for ( i=0; i<10; i++ ) {
 				if ( i != v_best ) refine_simplex_transform(s, i, 0.5);
 			}
@@ -227,22 +226,24 @@ static double refine_iteration(SimplexVertex *s, Deviation *d, int nf) {
 		
 	} else {
 	
+		#if 0
 		/* It's better.  Try to expand in this direction */
 		double fom_new_new;
 		SimplexVertex save;
 		
-//		printf("This is better.  Trying to expand...\n");
+		if ( debug ) printf("This is better.  Trying to expand...\n");
 		
 		save = s[v_worst];
 		refine_simplex_transform(s, v_worst, 2);
 		/* Better? */
 		fom_new_new = refine_mean_dev(d, nf, s, v_worst);
-//		printf("Mean deviation after expansion is %f nm^-1\n", fom_new_new/1e9);
+		if ( debug ) printf("Mean deviation after expansion is %f nm^-1\n", fom_new_new/1e9);
 		if ( fom_new_new > fom_new ) {
 			/* "Got too confident" */
 			s[v_worst] = save;
-//			printf("Got too confident - reverting\n");
+			if ( debug ) printf("Got too confident - reverting\n");
 		} /* else yay. */
+		#endif
 	
 	}
 	
@@ -257,6 +258,7 @@ double refine_do_cell(ControlContext *ctx) {
 	double delta;
 	int i, nf, f, it, maxiter;
 	const double tol = 0.001e9;	/* Stopping condition */
+	int debug = 1;
 	
 	if ( !ctx->cell_lattice ) {
 		displaywindow_error("No reciprocal unit cell has been found.", ctx->dw);
@@ -282,7 +284,7 @@ double refine_do_cell(ControlContext *ctx) {
 		}
 		
 	}
-	printf("RF: There are %i partnered features in total\n", nf);
+	if ( debug ) printf("RF: There are %i partnered features in total\n", nf);
 	
 	/* Initialise the 'deviation table' */
 	d = malloc(nf*sizeof(Deviation));
@@ -355,13 +357,15 @@ double refine_do_cell(ControlContext *ctx) {
 	
 		double conv;
 		
-//		printf("Simplex method iteration %i\n", it);
-		conv = refine_iteration(s, d, nf);
+		if ( debug ) printf("------------------- Simplex method iteration %i -------------------\n", it);
+		conv = refine_iteration(s, d, nf, debug);
 		if ( conv < tol ) {
-			printf("Converged after %i iterations (%f nm^-1)\n", it, conv/1e9);
+			if ( debug ) printf("RF: Converged after %i iterations (%f nm^-1)\n", it, conv/1e9);
 			break;
 		}
+		
 	}
+	if ( it == maxiter ) printf("RF: Did not converge.\n");
 	
 	/* Apply the final values to the cell */
 	ctx->cell->a.x += s[0].dax;	ctx->cell->b.x += s[0].dbx;	ctx->cell->c.x += s[0].dcx;