Nooooo! Use RADIANS internally...

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

10 files changed, 37 insertions, 48 deletions

diff --git a/src/displaywindow.c b/src/displaywindow.c
index 434508c..af17934 100644
--- a/src/displaywindow.c
+++ b/src/displaywindow.c
@@ -38,6 +38,7 @@
 #include "imagedisplay.h"
 #include "intensities.h"
 #include "glbits.h"
+#include "utils.h"
 
 static gint displaywindow_changeview(GtkWidget *widget, GtkRadioAction *action, DisplayWindow *dw) {
 
@@ -242,7 +243,7 @@ static gint displaywindow_setaxis_response(GtkWidget *widget, gint response, Dis
 		offset = gtk_entry_get_text(GTK_ENTRY(dw->tiltaxis_entry));
 		sscanf(offset, "%f", &off);
 		
-		mapping_adjust_axis(dw->ctx, off);
+		mapping_adjust_axis(dw->ctx, deg2rad(off));
 		
 	}
 	
diff --git a/src/glbits.c b/src/glbits.c
index bf523fa..5aed16b 100644
--- a/src/glbits.c
+++ b/src/glbits.c
@@ -27,6 +27,7 @@
 #include "trackball.h"
 #include "reflections.h"
 #include "image.h"
+#include "utils.h"
 
 /* Utility function to load and compile a shader, checking the info log */
 static GLhandleARB glbits_load_shader(const char *filename, GLenum type) {
@@ -588,7 +589,7 @@ void glbits_prepare(DisplayWindow *dw) {
 		 *	so rotate tilt axis anticlockwise by omega.
 		 *	Since the rotation is about +z, this is already anticlockwise
 		 *	when looking down z. */
-		glRotatef(ctx->images->images[0].omega, 0.0, 0.0, 1.0);
+		glRotatef(rad2deg(ctx->images->images[0].omega), 0.0, 0.0, 1.0);
 		glScalef(50.0, 1.0, 1.0);
 		glTranslatef(-0.5, 0.0, 0.0);
 		glMaterialfv(GL_FRONT_AND_BACK, GL_EMISSION, yellow);
diff --git a/src/image.h b/src/image.h
index c617330..6b2f4e7 100644
--- a/src/image.h
+++ b/src/image.h
@@ -42,9 +42,9 @@ typedef struct {
 typedef struct imagerecord_struct {
 
 	uint16_t		*image;
-	double			tilt;		/* Degrees.  Defines where the pattern lies in reciprocal space */
-	double			omega;		/* Degrees.  Defines where the pattern lies in reciprocal space */
-	double			slop;		/* Degrees.  Defines where the pattern lies "on the negative" */
+	double			tilt;		/* Radians.  Defines where the pattern lies in reciprocal space */
+	double			omega;		/* Radians.  Defines where the pattern lies in reciprocal space */
+	double			slop;		/* Radians.  Defines where the pattern lies "on the negative" */
 	
 	FormulationMode		fmode;
 	double			pixel_size;
diff --git a/src/imagedisplay.c b/src/imagedisplay.c
index bc2fa53..00b7285 100644
--- a/src/imagedisplay.c
+++ b/src/imagedisplay.c
@@ -192,13 +192,13 @@ static gboolean imagedisplay_redraw(GtkWidget *drawingarea, GdkEventExpose *even
 			imagedisplay->imagerecord.y_centre * scale,
 			(imagedisplay->imagerecord.x_centre + imagedisplay->imagerecord.width) * scale,
 			(imagedisplay->imagerecord.y_centre + imagedisplay->imagerecord.width
-					* tan(deg2rad(imagedisplay->imagerecord.omega))) * scale);
+					* tan(imagedisplay->imagerecord.omega)) * scale);
 		imagedisplay_draw_line(imagedisplay->gc_tiltaxis,
 			imagedisplay->imagerecord.x_centre * scale,
 			imagedisplay->imagerecord.y_centre * scale,
 			(imagedisplay->imagerecord.x_centre - imagedisplay->imagerecord.width) * scale,
 			(imagedisplay->imagerecord.y_centre - imagedisplay->imagerecord.width
-					* tan(deg2rad(imagedisplay->imagerecord.omega))) * scale);
+					* tan(imagedisplay->imagerecord.omega)) * scale);
 	}
 	
 	/* NB This calls the function above, which sorts out stuff */
diff --git a/src/iprtest.c b/src/iprtest.c
index 4a06a4a..976ad51 100644
--- a/src/iprtest.c
+++ b/src/iprtest.c
@@ -111,27 +111,6 @@ int main(int argc, char *argv[]) {
 	GtkWidget *window;
 	GtkWidget *drawing;
 	IPRTestContext *ctx;
-	double dlx, dly, dlz;
-	
-	mapping_rotate(1.0, 0.0, 0.0, &dlx, &dly, &dlz, 0.0, 0.0);
-	printf("omega=0.0, tilt=0.0   => image 'x' = { %f %f %f }\n", dlx, dly, dlz);
-	mapping_rotate(0.0, 1.0, 0.0, &dlx, &dly, &dlz, 0.0, 0.0);
-	printf("omega=0.0, tilt=0.0   => image 'y' = { %f %f %f }\n", dlx, dly, dlz);
-
-	mapping_rotate(1.0, 0.0, 0.0, &dlx, &dly, &dlz, deg2rad(90.0), 0.0);
-	printf("omega=90.0, tilt=0.0   => image 'x' = { %f %f %f }\n", dlx, dly, dlz);
-	mapping_rotate(0.0, 1.0, 0.0, &dlx, &dly, &dlz, deg2rad(90.0), 0.0);
-	printf("omega=90.0, tilt=0.0   => image 'y' = { %f %f %f }\n", dlx, dly, dlz);
-
-	mapping_rotate(1.0, 0.0, 0.0, &dlx, &dly, &dlz, 0.0, deg2rad(90.0));
-	printf("omega=0.0, tilt=90.0   => image 'x' = { %f %f %f }\n", dlx, dly, dlz);
-	mapping_rotate(0.0, 1.0, 0.0, &dlx, &dly, &dlz, 0.0, deg2rad(90.0));
-	printf("omega=0.0, tilt=90.0   => image 'y' = { %f %f %f }\n", dlx, dly, dlz);
-
-	mapping_rotate(1.0, 0.0, 0.0, &dlx, &dly, &dlz, deg2rad(90.0), deg2rad(90.0));
-	printf("omega=90.0, tilt=90.0   => image 'x' = { %f %f %f }\n", dlx, dly, dlz);
-	mapping_rotate(0.0, 1.0, 0.0, &dlx, &dly, &dlz, deg2rad(90.0), deg2rad(90.0));
-	printf("omega=90.0, tilt=90.0   => image 'y' = { %f %f %f }\n", dlx, dly, dlz);
 	
 	gtk_init(&argc, &argv);
 	
@@ -139,8 +118,8 @@ int main(int argc, char *argv[]) {
 	
 	ctx->image = malloc(sizeof(ImageRecord));
 	ctx->image->image = NULL;
-	ctx->image->tilt = deg2rad(90.0);
-	ctx->image->omega = deg2rad(90.0);
+	ctx->image->tilt = deg2rad(45.0);
+	ctx->image->omega = deg2rad(0.0);
 	ctx->image->slop = deg2rad(0.0);
 	ctx->image->fmode = FORMULATION_PIXELSIZE;
 	ctx->image->pixel_size = 0.6e8;
@@ -156,15 +135,15 @@ int main(int argc, char *argv[]) {
 	
 	ctx->cell = malloc(sizeof(Basis));
 	/* The "true" cell */
-	ctx->cell->a.x = 5.0e9;  ctx->cell->a.y = 0.0;  ctx->cell->a.z = 0.0;
-	ctx->cell->b.x = 0.0;  ctx->cell->b.y = 5.0e9;  ctx->cell->b.z = 0.0;
-	ctx->cell->c.x = 0.0;  ctx->cell->c.y = 0.0;  ctx->cell->c.z = 5.0e9;
+	ctx->cell->a.x = 5.0e9;  ctx->cell->a.y = 0.0e9;  ctx->cell->a.z = 0.0e9;
+	ctx->cell->b.x = 0.0e9;  ctx->cell->b.y = 5.0e9;  ctx->cell->b.z = 0.0e9;
+	ctx->cell->c.x = 0.0e9;  ctx->cell->c.y = 0.0e9;  ctx->cell->c.z = 5.0e9;
 	ctx->reflections = reflection_list_from_cell(ctx->cell);
 	ctx->image->features = reproject_get_reflections(ctx->image, ctx->reflections);
 	/* The "model" cell */
-	ctx->cell->a.x = 5.2e9;  ctx->cell->a.y = 0.1e9;  ctx->cell->a.z = 0.0;
-	ctx->cell->b.x = 0.1e9;  ctx->cell->b.y = 4.5e9;  ctx->cell->b.z = 0.0;
-	ctx->cell->c.x = 0.0;  ctx->cell->c.y = 0.0;  ctx->cell->c.z = 5.0e9;
+	ctx->cell->a.x = 5.2e9;  ctx->cell->a.y = 0.1e9;  ctx->cell->a.z = 0.0e9;
+	ctx->cell->b.x = 0.2e9;  ctx->cell->b.y = 4.8e9;  ctx->cell->b.z = 0.0e9;
+	ctx->cell->c.x = 0.0e9;  ctx->cell->c.y = 0.0e9;  ctx->cell->c.z = 5.0e9;
 	ctx->gen = reflection_list_from_cell(ctx->cell);
 	ctx->image->rflist = reproject_get_reflections(ctx->image, ctx->gen);
 	
diff --git a/src/mapping.c b/src/mapping.c
index bc6e062..c9872f9 100644
--- a/src/mapping.c
+++ b/src/mapping.c
@@ -19,6 +19,7 @@
 #include "image.h"
 #include "displaywindow.h"
 #include "cache.h"
+#include "utils.h"
 
 void mapping_rotate(double x, double y, double z, double *ddx, double *ddy, double *ddz, double omega, double tilt) {
 
@@ -66,8 +67,8 @@ int mapping_map_to_space(ImageFeature *refl, double *ddx, double *ddy, double *d
 	
 	imagerecord = refl->parent;
 	x = refl->x - imagerecord->x_centre;	y = refl->y - imagerecord->y_centre;
-	tilt = 2*M_PI*(imagerecord->tilt/360);	/* Convert to Radians */
-	omega = 2*M_PI*(imagerecord->omega/360);	/* Likewise */
+	tilt = imagerecord->tilt;
+	omega = imagerecord->omega;
 	k = 1/imagerecord->lambda;
 	
 	/* Calculate an angular description of the reflection */
@@ -162,7 +163,8 @@ void mapping_adjust_axis(ControlContext *ctx, double offset) {
 	int i;
 
 	for ( i=0; i<ctx->images->n_images; i++ ) {
-		printf("Image #%3i: old omega=%f deg, new omega=%f deg\n", i, ctx->images->images[i].omega, ctx->images->images[i].omega+offset);
+		printf("Image #%3i: old omega=%f deg, new omega=%f deg\n", i, rad2deg(ctx->images->images[i].omega),
+							rad2deg(ctx->images->images[i].omega+offset));
 		ctx->images->images[i].omega += offset;
 	}
 	
diff --git a/src/mrc.c b/src/mrc.c
index b68166f..57df51c 100644
--- a/src/mrc.c
+++ b/src/mrc.c
@@ -79,7 +79,7 @@ int mrc_read(ControlContext *ctx) {
 		} else {
 			ctx->lambda = lambda(1000*ext[i].voltage);
 		}
-		ctx->omega = ext[i].tilt_axis;
+		ctx->omega = deg2rad(ext[i].tilt_axis);
 		ctx->pixel_size = ext[i].pixel_size;
 		
 		fseek(fh, mrc.next + sizeof(MRCHeader) + mrc.nx*mrc.ny*2*i, SEEK_SET);
@@ -93,7 +93,7 @@ int mrc_read(ControlContext *ctx) {
 		}
 		free(image);
 		
-		image_add(ctx->images, uimage, mrc.nx, mrc.ny, ext[i].a_tilt, ctx);
+		image_add(ctx->images, uimage, mrc.nx, mrc.ny, deg2rad(ext[i].a_tilt), ctx);
 		
 	}
 
diff --git a/src/qdrp.c b/src/qdrp.c
index 5c7df98..aae0c2e 100644
--- a/src/qdrp.c
+++ b/src/qdrp.c
@@ -20,6 +20,7 @@
 #include "readpng.h"
 #include "control.h"
 #include "reflections.h"
+#include "utils.h"
 
 static void qdrp_chomp(char *str) {
 
@@ -65,9 +66,9 @@ static int qdrp_parsefileline(ControlContext *ctx, const char *line) {
 	skip_chars;
 	file = strndup(line+mark, pos-mark);
 					
-	tilt = strtod(tilt_s, NULL);
+	tilt = deg2rad(strtod(tilt_s, NULL));
 	free(tilt_s);
-	printf("QD: Reading file: Tilt = %f deg, File='%s'\n", tilt, file);
+	printf("QD: Reading file: Tilt = %f deg, File='%s'\n", rad2deg(tilt), file);
 	
 	if ( readpng_read(file, tilt, ctx) ) {
 		printf("Reconstruction failed.\n");
@@ -102,7 +103,7 @@ static int qdrp_parseline(ControlContext *ctx, const char *line) {
 		}
 					
 		if ( !ctx->omega_set ) {
-			printf("QD: Parameter 'omega' not specified, and not using tilt-axis-search.\n");
+			printf("QD: Parameter 'omega' not specified.\n");
 			return 1;
 		}
 				
@@ -175,10 +176,10 @@ static int qdrp_parseline(ControlContext *ctx, const char *line) {
 			mark = pos;
 			skip_chars;
 			omega_s = strndup(line+mark, pos-mark);
-			ctx->omega = strtod(omega_s, NULL);
+			ctx->omega = deg2rad(strtod(omega_s, NULL));
 			free(omega_s);
 			ctx->omega_set = 1;
-			printf("QD: omega = %f deg\n", ctx->omega);
+			printf("QD: omega = %f deg\n", rad2deg(ctx->omega));
 				
 		}
 		
diff --git a/src/refine.c b/src/refine.c
index 3a01a9d..fe932d0 100644
--- a/src/refine.c
+++ b/src/refine.c
@@ -192,6 +192,11 @@ ImageFeature *refine_fit_image(Basis *cell, ImageRecord *image) {
 		
 		/* Prepare for solving */
 		refine_fix_unconstrained(M);
+		gsl_matrix_set(M, 2, 0, 0.0);
+		gsl_matrix_set(M, 2, 1, 0.0);
+		gsl_matrix_set(M, 2, 2, 1.0);
+		gsl_matrix_set(M, 0, 2, 0.0);
+		gsl_matrix_set(M, 1, 2, 0.0);
 		matrix_vector_show(M, p, "RF: ");
 		
 		/* Calculate and display the determinant */
diff --git a/src/reproject.c b/src/reproject.c
index 41ec8c5..4036863 100644
--- a/src/reproject.c
+++ b/src/reproject.c
@@ -58,8 +58,8 @@ ImageFeatureList *reproject_get_reflections(ImageRecord *image, ReflectionList *
 			
 	flist = image_feature_list_new();
 	
-	tilt = deg2rad(image->tilt);
-	omega = deg2rad(image->omega);
+	tilt = image->tilt;
+	omega = image->omega;
 	
 	/* Calculate the (normalised) incident electron wavevector
 	 *	n is rotated "into" the reconstruction, so only one omega step. */