DirAx stuff

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

14 files changed, 641 insertions, 104 deletions

diff --git a/src/Makefile.am b/src/Makefile.am
index bdfe197..0be30fc 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -1,7 +1,8 @@
 bin_PROGRAMS = dtr
 dtr_SOURCES = main.c displaywindow.c trackball.c reflections.c readpng.c mrc.c imagedisplay.c utils.c itrans.c qdrp.c ipr.c cache.c \
-		itrans-threshold.c itrans-zaefferer.c itrans-lsq.c itrans-stat.c control.c mapping.c reproject.c prealign.c basis.c
-dtr_LDADD = @LIBS@ @GTK_LIBS@ -lm @GTKGLEXT_LIBS@ -lgsl -lgslcblas
+		itrans-threshold.c itrans-zaefferer.c itrans-lsq.c itrans-stat.c control.c mapping.c reproject.c prealign.c basis.c \
+		dirax.c
+dtr_LDADD = @LIBS@ @GTK_LIBS@ -lm @GTKGLEXT_LIBS@ -lgsl -lgslcblas -lutil
 AM_CFLAGS = -Wall -g @CFLAGS@ @GTK_CFLAGS@ @GTKGLEXT_CFLAGS@
 AM_CPPFLAGS = -DDATADIR=\""$(datadir)"\"
 
diff --git a/src/basis.c b/src/basis.c
index 19004ee..8a282e3 100644
--- a/src/basis.c
+++ b/src/basis.c
@@ -87,7 +87,7 @@ static int basis_lfom(ControlContext *ctx, double vx, double vy, double vz) {
 		
 			Reflection *check;
 			
-			check = reflection_find_nearest(ctx->reflectionlist, tcentre->x+vx*j, tcentre->y+vy*j, tcentre->z+vz*j);
+			check = reflectionlist_find_nearest(ctx->reflectionlist, tcentre->x+vx*j, tcentre->y+vy*j, tcentre->z+vz*j);
 			if ( check && (distance3d(check->x, check->y, check->z, tcentre->x+vx*j, tcentre->y+vy*j, tcentre->z+vz*j) < tol) ) {
 				lfom++;
 			}
@@ -135,14 +135,14 @@ static ReflectionList *basis_find_seeds(ControlContext *ctx) {
 	z = z_temp;
 	
 	/* Find the point in the middle of the "wedge" */
-	scale = reflection_largest_g(ctx->reflectionlist)/6;
+	scale = reflectionlist_largest_g(ctx->reflectionlist)/6;
 	x *= scale;
 	y *= scale;
 	z *= scale;
 	reflection_add(ctx->reflectionlist, x, y, z, 1.0, REFLECTION_VECTOR_MARKER_2);
 	
 	/* Find an "origin" reflection */
-	centre = reflection_find_nearest(ctx->reflectionlist, x, y, z);
+	centre = reflectionlist_find_nearest(ctx->reflectionlist, x, y, z);
 	if ( !centre ) return NULL;
 	centre->found = 1;
 	reflection_add(ctx->reflectionlist, centre->x, centre->y, centre->z, 1.0, REFLECTION_GENERATED);
@@ -161,7 +161,7 @@ static ReflectionList *basis_find_seeds(ControlContext *ctx) {
 			accept = 1;
 		
 			/* Find a "candidate vector" reflection */
-			vector = reflection_find_nearest_longer_unknown(ctx->reflectionlist, centre->x, centre->y, centre->z, 1e9); /* 0.5 A^-1 */
+			vector = reflectionlist_find_nearest_longer_unknown(ctx->reflectionlist, centre->x, centre->y, centre->z, 1e9);
 			if ( !vector ) {
 				printf("BS: Only found %i seeds\n", i);
 				return seeds;
@@ -174,7 +174,7 @@ static ReflectionList *basis_find_seeds(ControlContext *ctx) {
 			vz = vector->z - centre->z;
 			
 			/* Proximity test: don't duplicate seeds */
-			check = reflection_find_nearest(seeds, vx, vy, vz);
+			check = reflectionlist_find_nearest(seeds, vx, vy, vz);
 			if ( check ) {
 				if ( distance3d(vx, vy, vz, check->x, check->y, check->z) < 1e9 ) {
 					/* Too close to another seed */
@@ -209,7 +209,7 @@ static Basis *basis_assemble_seeds(ReflectionList *seeds) {
 	int i, j, b;
 	ReflectionList *seeds_sorted;
 
-	seeds_sorted = reflection_sort_lfom(seeds);
+	seeds_sorted = reflectionlist_sort_lfom(seeds);
 	
 	basis = malloc(10*sizeof(Basis));
 	
diff --git a/src/cache.c b/src/cache.c
index 078fa72..ad47ca2 100644
--- a/src/cache.c
+++ b/src/cache.c
@@ -41,10 +41,10 @@ ReflectionList *cache_load(const char *filename) {
 	reflectionlist = reflectionlist_new();
 	f = fopen(filename, "rb");
 	if ( !f ) {
-		fprintf(stderr, "Couldn't read reflections from cache\n");
+		fprintf(stderr, "Couldn't open cache file\n");
 	}
 	if ( fread(&ch, sizeof(CacheHeader), 1, f) == 0 ) {
-		fprintf(stderr, "Couldn't read reflections from cache\n");
+		fprintf(stderr, "Couldn't read cache header\n");
 		fclose(f);
 		return NULL;
 	}
@@ -58,7 +58,7 @@ ReflectionList *cache_load(const char *filename) {
 			fprintf(stderr, "Couldn't read reflections from cache\n");
 			fclose(f);
 			free(cr);
-			reflection_clear(reflectionlist);
+			reflectionlist_clear(reflectionlist);
 			return NULL;
 		}
 			
@@ -129,9 +129,13 @@ unsigned int cache_is_cachefile(const char *filename) {
 	size_t nread;
 
 	fh = fopen(filename, "rb");
+	if ( !fh ) {
+		printf("Couldn't open file '%s'\n", filename);
+		return 0;
+	}
 	nread = fread(&ch, sizeof(CacheHeader), 1, fh);
 	fclose(fh);
-
+	
 	if ( nread != 1 ) {
 		return 0;
 	}
@@ -140,11 +144,6 @@ unsigned int cache_is_cachefile(const char *filename) {
 		return 1;
 	}
 	
-	/* Backwards compatability */
-	if ( strncmp(ch.top, "DTR-CACHE", 9) == 0 ) {
-		return 1;
-	}
-	
 	return 0;
 
 }
diff --git a/src/control.c b/src/control.c
index 270be6f..032700d 100644
--- a/src/control.c
+++ b/src/control.c
@@ -53,6 +53,7 @@ ControlContext *control_ctx_new() {
 	ctx->x_centre = 0;
 	ctx->y_centre = 0;
 	ctx->have_centres = 0;
+	ctx->cell = NULL;
 	
 	return ctx;
 
diff --git a/src/control.h b/src/control.h
index 5b0d187..bb42395 100644
--- a/src/control.h
+++ b/src/control.h
@@ -26,7 +26,8 @@ typedef enum ift_enum {
 	INPUT_NONE,
 	INPUT_QDRP,
 	INPUT_MRC,
-	INPUT_CACHE
+	INPUT_CACHE,
+	INPUT_DRX
 } InputFileType;
 
 typedef enum {
@@ -44,8 +45,9 @@ typedef enum {
 } PeakSearchMode;
 
 typedef enum {
-	RECONSTRUCTION_MAPPING,
-	RECONSTRUCTION_PREDICTION
+	RECONSTRUCTION_MAPPING		= 0,
+	RECONSTRUCTION_PREDICTION	= 1<<0,
+	RECONSTRUCTION_DIRAX		= 1<<1
 } ReconstructionMode;
 
 typedef struct imagerecord_struct {
@@ -85,6 +87,7 @@ typedef struct cctx_struct {
 	unsigned int			prealign;
 	unsigned int			savecache;
 	unsigned int			have_centres;
+	unsigned int			use_dirax;
 	
 	/* Input filename */
 	char 				*filename;
@@ -113,18 +116,31 @@ typedef struct cctx_struct {
 	/* Output */
 	struct reflectionlist_struct	*reflectionlist;
 	struct dw_struct		*dw;
+	struct basis_struct		*cell;		
 	
 	/* GTK bits */
 	GtkWidget			*combo_algorithm;
 	GtkWidget			*combo_peaksearch;
 	GtkWidget			*checkbox_prealign;
 	GtkWidget			*checkbox_savecache;
+	GtkWidget			*checkbox_dirax;
 	
 	/* IPR stuff */
 	int				ipr_cur_image;
 	struct imagedisplay_struct	*ipr_id;
 	struct reflectionlist_struct	*ipr_lat;
 	struct basis_struct		*ipr_basis;
+	
+	/* DirAx low-level stuff */
+	GIOChannel			*dirax;
+	int				dirax_pty;
+	char				*dirax_rbuffer;
+	int				dirax_rbufpos;
+	int				dirax_rbuflen;
+	
+	/* DirAx high-level stuff */
+	int				dirax_step;
+	int				dirax_read_cell;
 
 } ControlContext;
 
diff --git a/src/dirax.c b/src/dirax.c
new file mode 100644
index 0000000..4f4d135
--- /dev/null
+++ b/src/dirax.c
@@ -0,0 +1,437 @@
+/*
+ * dirax.c
+ *
+ * Invoke the DirAx auto-indexing program
+ * also: handle DirAx input files
+ *
+ * (c) 2007 Thomas White <taw27@cam.ac.uk>
+ *
+ *  dtr - Diffraction Tomography Reconstruction
+ *
+ */
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include <glib.h>
+#include <stdlib.h>
+#include <stdio.h>
+#include <math.h>
+#include <string.h>
+#include <pty.h>
+#include <unistd.h>
+#include <sys/wait.h>
+#include <fcntl.h>
+#include <assert.h>
+#include <sys/ioctl.h>
+#include <termio.h>
+	
+#include "control.h"
+#include "reflections.h"
+#include "utils.h"
+#include "basis.h"
+#include "displaywindow.h"
+
+typedef enum {
+	DIRAX_INPUT_NONE,
+	DIRAX_INPUT_LINE,
+	DIRAX_INPUT_PROMPT
+} DirAxInputType;
+
+static void dirax_parseline(const char *line, ControlContext *ctx) {
+	
+	int i, rf;
+	char *copy;
+	
+	copy = strdup(line);
+	for ( i=0; i<strlen(copy); i++ ) {
+		if ( copy[i] == '\r' ) copy[i]='r';
+		if ( copy[i] == '\n' ) copy[i]='n';
+	}
+	printf("DX: DirAx: %s\n", copy);
+	free(copy);
+	
+	/* Is this the first line of a unit cell specification? */
+	rf = 0; i = 0;
+	while ( (i<strlen(line)) && ((line[i] == 'R') || (line[i] == 'D') || (line[i] == ' ')) ) {
+		if ( line[i] == 'R' ) rf = 1;
+		if ( (line[i] == 'D') && rf ) {
+			ctx->dirax_read_cell = 1;
+			if ( ctx->cell ) {
+				free(ctx->cell);
+			}
+			ctx->cell = malloc(sizeof(Basis));
+			ctx->cell->a.x = 0.0;  ctx->cell->a.y = 0.0;  ctx->cell->a.z = 0.0;
+			ctx->cell->b.x = 0.0;  ctx->cell->b.y = 0.0;  ctx->cell->b.z = 0.0;
+			ctx->cell->c.x = 0.0;  ctx->cell->c.y = 0.0;  ctx->cell->c.z = 0.0;
+			return;
+		}
+		i++;
+	}
+	
+	/* Parse unit cell vectors as appropriate */
+	if ( ctx->dirax_read_cell == 1 ) {
+		/* First row of unit cell values */
+		float x, y, z;
+		sscanf(line, "%f %f %f", &x, &y, &z);
+		ctx->cell->a.x = x*1e10;  ctx->cell->a.y = y*1e10;  ctx->cell->a.z = z*1e10;
+		ctx->dirax_read_cell++;
+		return;
+	} else if ( ctx->dirax_read_cell == 2 ) {
+		/* First row of unit cell values */
+		float x, y, z;
+		sscanf(line, "%f %f %f", &x, &y, &z);
+		ctx->cell->b.x = x*1e10;  ctx->cell->b.y = y*1e10;  ctx->cell->b.z = z*1e10;
+		ctx->dirax_read_cell++;
+		return;
+	} else if ( ctx->dirax_read_cell == 3 ) {
+		/* First row of unit cell values */
+		float x, y, z;
+		sscanf(line, "%f %f %f", &x, &y, &z);
+		ctx->cell->c.x = x*1e10;  ctx->cell->c.y = y*1e10;  ctx->cell->c.z = z*1e10;
+		printf("DX: Read a reciprocal unit cell\n");
+		displaywindow_update(ctx->dw);
+		ctx->dirax_read_cell = 0;
+		return;
+	}
+	
+	ctx->dirax_read_cell = 0;
+
+}
+
+static void dirax_sendline(const char *line, ControlContext *ctx) {
+	
+	char *copy;
+	int i;
+	
+	write(ctx->dirax_pty, line, strlen(line));
+
+	copy = strdup(line);
+	for ( i=0; i<strlen(copy); i++ ) {
+		if ( copy[i] == '\r' ) copy[i]='\0';
+		if ( copy[i] == '\n' ) copy[i]='\0';
+	}
+	printf("sent '%s'", copy);	/* No newline here */
+	free(copy);
+	
+}
+
+static void dirax_send_next(ControlContext *ctx) {
+
+	switch ( ctx->dirax_step ) {
+	
+		case 0 : {
+			dirax_sendline("read dtr.drx\n", ctx);
+			ctx->dirax_step++;
+			break;
+		}
+		
+		case 1 : {
+			dirax_sendline("go\n", ctx);
+			ctx->dirax_step++;
+			break;
+		}
+		
+		case 2 : {
+			dirax_sendline("cell\n", ctx);
+			ctx->dirax_step++;
+			break;
+		}
+		
+		default: {
+			printf("idle");
+		}
+		
+	}
+	
+	fsync(ctx->dirax_pty);
+	
+}
+
+static gboolean dirax_readable(GIOChannel *dirax, GIOCondition condition, ControlContext *ctx) {
+	
+	int rval;
+	
+	rval = read(ctx->dirax_pty, ctx->dirax_rbuffer+ctx->dirax_rbufpos, ctx->dirax_rbuflen-ctx->dirax_rbufpos);
+	
+	if ( (rval == -1) || (rval == 0) ) {
+	
+		printf("DX: Lost connection to DirAx\n");
+		g_io_channel_shutdown(ctx->dirax, FALSE, NULL);
+		ctx->dirax = NULL;
+		return FALSE;
+		
+	} else {
+	
+		int no_string = 0;
+		
+		ctx->dirax_rbufpos += rval;
+		assert(ctx->dirax_rbufpos <= ctx->dirax_rbuflen);
+		
+		while ( (!no_string) && (ctx->dirax_rbufpos > 0) ) {
+	
+			int i;
+			int block_ready = 0;
+			DirAxInputType type = DIRAX_INPUT_NONE;
+			
+			/* See if there's a full line in the buffer yet */
+			for ( i=0; i<ctx->dirax_rbufpos; i++ ) {	/* Means the last value looked at is roffset-1 */
+			
+				if ( (ctx->dirax_rbuffer[i] == '\n') ) { 
+					block_ready = 1;
+					type = DIRAX_INPUT_LINE;
+					break;
+				}
+				
+				if ( i <= ctx->dirax_rbufpos-7 ) {
+					if ( strncmp(ctx->dirax_rbuffer+i, "Dirax> ", 7) == 0 ) {
+						block_ready = 1;
+						type = DIRAX_INPUT_PROMPT;
+						break;
+					}
+				}
+			
+			}
+			
+			if ( block_ready == 1 ) {
+			
+				switch ( type ) {
+				
+					case DIRAX_INPUT_LINE : {
+					
+						char *block_buffer = NULL;
+						unsigned int new_rbuflen;
+						unsigned int endbit_length;
+						
+						block_buffer = malloc(i+1);
+						memcpy(block_buffer, ctx->dirax_rbuffer, i);
+						block_buffer[i] = '\0';
+						dirax_parseline(block_buffer, ctx);
+						free(block_buffer);
+						
+						/* Now the block's been parsed, it should be forgotten about */
+						endbit_length = i+1;
+						memmove(ctx->dirax_rbuffer, ctx->dirax_rbuffer + endbit_length,
+										ctx->dirax_rbuflen - endbit_length);
+						/* Subtract the number of bytes removed */
+						ctx->dirax_rbufpos = ctx->dirax_rbufpos - endbit_length;
+						new_rbuflen = ctx->dirax_rbuflen - endbit_length;
+						if ( new_rbuflen == 0 ) {
+							new_rbuflen = 256;
+						}
+						ctx->dirax_rbuffer = realloc(ctx->dirax_rbuffer, new_rbuflen);
+						ctx->dirax_rbuflen = new_rbuflen;
+						
+						break;
+						
+					}
+					
+					case DIRAX_INPUT_PROMPT : {
+						
+						memmove(ctx->dirax_rbuffer+i, ctx->dirax_rbuffer+i+7, 7);
+						ctx->dirax_rbufpos -= 7;
+						printf("DX: Got prompt: ");
+						dirax_send_next(ctx);
+						printf("\n");
+						
+						break;
+					
+					}
+					
+					default : {
+						printf("DX: Unrecognised input mode (this never happens!)\n");
+						abort();
+					}
+				
+				}
+			
+			} else {
+			
+				if ( ctx->dirax_rbufpos == ctx->dirax_rbuflen ) {
+			
+					/* More buffer space is needed */
+					ctx->dirax_rbuffer = realloc(ctx->dirax_rbuffer, ctx->dirax_rbuflen + 256);
+					ctx->dirax_rbuflen = ctx->dirax_rbuflen + 256;
+					/* The new space gets used at the next read, shortly... */
+			
+				}
+				no_string = 1;
+			
+			}
+			
+		}
+		
+	}
+	
+	return TRUE;
+	
+}
+
+void dirax_invoke(ControlContext *ctx) {
+
+	FILE *fh;
+	Reflection *ref;
+	int pty;
+	unsigned int opts;
+	pid_t dirax_pid;
+	
+	printf("DX: Starting DirAx...\n");
+	
+	fh = fopen("dtr.drx", "w");
+	if ( !fh ) {
+		printf("DX: Couldn't open temporary file dtr.drx\n");
+		return;
+	}
+	fprintf(fh, "%f\n", 0.5);  /* Lie about the wavelength.  DirAx can't handle the truth. */
+	ref = ctx->reflectionlist->reflections;
+	while ( ref ) {
+		fprintf(fh, "%f %f %f %f\n", ref->x/1e10, ref->y/1e10, ref->z/1e10, ref->intensity);
+		ref = ref->next;
+	}
+	fclose(fh);
+	
+	struct termios tty;
+	if ( ioctl(0, TCGETS, &tty) ) {
+		printf("DX: Failed to read pty attributes\n");
+		return;
+	}
+	tty.c_lflag &= ~ECHO;
+	dirax_pid = forkpty(&pty, NULL, &tty, NULL);
+	if ( dirax_pid == -1 ) {
+		printf("DX: Failed to fork.\n");
+		return;
+	}
+	if ( dirax_pid == 0 ) {
+		
+		/* Child process: invoke DirAx */
+		execlp("dirax", "", (char *)NULL);
+		printf("(from the Other Side) Failed to invoke DirAx.\n");
+		_exit(0);
+		
+	}
+	
+	ctx->dirax_pty = pty;
+	ctx->dirax_rbuffer = malloc(256);
+	ctx->dirax_rbuflen = 256;
+	ctx->dirax_rbufpos = 0;
+	
+	/* Set non-blocking */
+	opts = fcntl(pty, F_GETFL);
+	fcntl(pty, F_SETFL, opts | O_NONBLOCK);
+	
+	ctx->dirax_step = 0;
+	ctx->dirax_read_cell = 0;
+	ctx->dirax = g_io_channel_unix_new(pty);
+	g_io_add_watch(ctx->dirax, G_IO_IN | G_IO_HUP, (GIOFunc)dirax_readable, ctx);
+	
+	return;
+	
+}
+
+/* Despite being part of the same module, this has very little to do with invoking DirAx */
+ReflectionList *dirax_load(const char *filename) {
+
+	FILE *fh;
+	char line[256];
+	ReflectionList *list;
+	int lambda_set = 0;
+	
+	fh = fopen(filename, "r");
+	if ( !fh ) {
+		printf("Couldn't open file '%s'\n", filename);
+		return 0;
+	}
+	
+	list = reflectionlist_new();
+	
+	while ( !feof(fh) ) {
+	
+		size_t ptr;
+		float lambda, theta, chib, phib;
+		
+		fgets(line, 255, fh);
+		ptr = skipspace(line);
+		if ( line[ptr] == '!' ) continue;
+		if ( line[ptr] == '\n' ) continue;
+		if ( line[ptr] == '\r' ) continue;
+		if ( sscanf(line+ptr, "%f %f %f\n", &theta, &phib, &chib) == 3 ) {
+			
+			double s, x, y, z;
+			float blah, intensity;
+			
+			/* Try to find an intensity value. Use dummy value if it fails */
+			if ( sscanf(line+ptr, "%f %f %f %f\n", &blah, &blah, &blah, &intensity) != 4 ) {
+				intensity = 1.0;
+			}
+				
+			if ( !lambda_set ) {
+				printf("DX: Wavelength not specified\n");
+				continue;
+			}
+			
+			chib = deg2rad(chib);
+			phib = deg2rad(phib);
+			theta = deg2rad(theta);
+			s = 2*(sin(theta)/lambda);
+			x = -s*cos(chib)*sin(phib);
+			y = +s*cos(chib)*cos(phib);
+			z = +s*sin(chib);
+			reflection_add(list, x, y, z, 1.0, REFLECTION_NORMAL);
+			
+			continue;
+			
+		}
+		if ( sscanf(line+ptr, "%f\n", &lambda) == 1 ) {
+			if ( lambda_set ) {
+				printf("DX: Warning: Found something which looks like a second wavelength\n");
+			}
+			lambda /= 1e10; /* Convert from A to m */
+			lambda_set = 1;
+		}
+		
+	}
+	
+	fclose(fh);
+	
+	return list;
+	
+}
+
+int dirax_is_drxfile(const char *filename) {
+
+	FILE *fh;
+	float lambda;
+	char line[256];
+	
+	fh = fopen(filename, "r");
+	if ( !fh ) {
+		printf("Couldn't open file '%s'\n", filename);
+		return 0;
+	}
+	
+	while ( !feof(fh) ) {
+	
+		size_t ptr;
+		
+		fgets(line, 255, fh);
+		ptr = skipspace(line);
+		if ( line[ptr] == '!' ) continue;
+		if ( line[ptr] == '\n' ) continue;
+		if ( line[ptr] == '\r' ) continue;
+		fscanf(fh, "%f\n", &lambda);
+		fclose(fh);
+		if ( lambda > 0.5 ) {
+			return 1;
+		} else {
+			return 0;
+		}
+		
+	}
+	
+	fclose(fh);
+	
+	return 0;
+	
+}
+    
diff --git a/src/dirax.h b/src/dirax.h
new file mode 100644
index 0000000..b06c6c9
--- /dev/null
+++ b/src/dirax.h
@@ -0,0 +1,28 @@
+/*
+ * dirax.h
+ *
+ * Invoke the DirAx auto-indexing program
+ * also: handle DirAx input files
+ *
+ * (c) 2007 Thomas White <taw27@cam.ac.uk>
+ *
+ *  dtr - Diffraction Tomography Reconstruction
+ *
+ */
+ 
+#ifndef DIRAX_H
+#define DIRAX_H
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include "control.h"
+#include "reflections.h"
+
+extern void dirax_invoke(ControlContext *ctx);
+extern ReflectionList *dirax_load(const char *filename);
+extern int dirax_is_drxfile(const char *filename);
+
+#endif	/* DIRAX_H */
+
diff --git a/src/displaywindow.c b/src/displaywindow.c
index b38007d..866f63f 100644
--- a/src/displaywindow.c
+++ b/src/displaywindow.c
@@ -31,13 +31,14 @@
 #include "ipr.h"
 #include "displaywindow.h"
 #include "basis.h"
+#include "dirax.h"
 
 enum {
 	DW_ORTHO,
 	DW_PERSPECTIVE
 };
 
-static void displaywindow_about() {
+static void displaywindow_about(GtkWidget *widget, DisplayWindow *dw) {
 
 	GtkWidget *window;
 
@@ -70,10 +71,14 @@ static void displaywindow_about() {
 
 }
 
-static void displaywindow_close() {
+static void displaywindow_close(GtkWidget *widget, DisplayWindow *dw) {
 	gtk_exit(0);
 }
 
+static void displaywindow_dirax(GtkWidget *widget, DisplayWindow *dw) {
+	dirax_invoke(dw->ctx);
+}
+
 static void displaywindow_gl_set_ortho(DisplayWindow *dw, GLfloat w, GLfloat h) {
 
 	glMatrixMode(GL_PROJECTION);
@@ -169,6 +174,7 @@ static void displaywindow_addmenubar(DisplayWindow *dw) {
 		{ "ViewAction", NULL, "_View", NULL, NULL, NULL },
 
 		{ "ToolsAction", NULL, "_Tools", NULL, NULL, NULL },
+		{ "DirAxAction", GTK_STOCK_EXECUTE, "Invoke _DirAx", "<Ctrl>D", NULL, G_CALLBACK(displaywindow_dirax) },
 		
 		{ "HelpAction", NULL, "_Help", NULL, NULL, NULL },
 		{ "AboutAction", GTK_STOCK_ABOUT, "_About DTR...", NULL, NULL, G_CALLBACK(displaywindow_about) },
@@ -343,7 +349,6 @@ static void displaywindow_gl_create_list(DisplayWindow *dw) {
 	GLfloat green[] = { 0.0, 1.0, 0.0, 1.0 };
 	GLfloat yellow[] = { 1.0, 1.0, 0.0, 1.0 };
 	GLfloat glass[] = { 0.2, 0.0, 0.8, 000.1 };
-	GLfloat glass_spec[] = { 0.8, 0.8, 0.8, 1.0 };
 	GLfloat black[] = { 0.0, 0.0, 0.0, 1.0 };
 	Reflection *reflection;
 	GLUquadricObj *quadric;
@@ -360,15 +365,15 @@ static void displaywindow_gl_create_list(DisplayWindow *dw) {
 	}
 	reflection = ctx->reflectionlist->reflections;
 	i = 0;
-	do {
+	while ( reflection != NULL ) {
 		if ( reflection->type == REFLECTION_NORMAL ) i++;
 		reflection = reflection->next;
-	} while ( reflection != NULL );
+	};
 	dw->gl_ref_num_vertices = i;
 	i = 0;
 	reflection = ctx->reflectionlist->reflections;
 	vertices = malloc(3*dw->gl_ref_num_vertices*sizeof(GLfloat));
-	do {
+	while ( reflection != NULL ) {
 		if ( reflection->type == REFLECTION_NORMAL ) {
 			vertices[3*i + 0] = reflection->x/1e9;
 			vertices[3*i + 1] = reflection->y/1e9;
@@ -376,7 +381,7 @@ static void displaywindow_gl_create_list(DisplayWindow *dw) {
 			i++;
 		}
 		reflection = reflection->next;
-	} while ( reflection != NULL );	
+	};	
 	if ( dw->gl_use_buffers ) {
 		glBindBufferARB(GL_ARRAY_BUFFER, dw->gl_ref_vertex_buffer);
 		glBufferDataARB(GL_ARRAY_BUFFER, 3*dw->gl_ref_num_vertices*sizeof(GLfloat), vertices, GL_STATIC_DRAW);
@@ -392,21 +397,21 @@ static void displaywindow_gl_create_list(DisplayWindow *dw) {
 	}
 	reflection = ctx->reflectionlist->reflections;
 	i = 0;
-	do {
+	while ( reflection != NULL ) {
 		if ( reflection->type == REFLECTION_MARKER ) i++;
 		reflection = reflection->next;
-	} while ( reflection != NULL );
+	};
 	dw->gl_marker_num_vertices = i*VERTICES_IN_A_BLOB;
 	i = 0;
 	reflection = ctx->reflectionlist->reflections;
 	vertices = malloc(3*dw->gl_marker_num_vertices*sizeof(GLfloat));
 	normals = malloc(3*dw->gl_marker_num_vertices*sizeof(GLfloat));
-	do {
+	while ( reflection != NULL ) {
 		if ( reflection->type == REFLECTION_MARKER ) {
 			DRAW_BLOB
 		}
 		reflection = reflection->next;
-	} while ( reflection != NULL );	
+	};	
 	if ( dw->gl_use_buffers ) {
 		glBindBufferARB(GL_ARRAY_BUFFER, dw->gl_marker_vertex_buffer);
 		glBufferDataARB(GL_ARRAY_BUFFER, 3*dw->gl_marker_num_vertices*sizeof(GLfloat), vertices, GL_STATIC_DRAW);
@@ -426,21 +431,21 @@ static void displaywindow_gl_create_list(DisplayWindow *dw) {
 	}
 	reflection = ctx->reflectionlist->reflections;
 	i = 0;
-	do {
+	while ( reflection != NULL ) {
 		if ( reflection->type == REFLECTION_GENERATED ) i++;
 		reflection = reflection->next;
-	} while ( reflection != NULL );
+	};
 	dw->gl_gen_num_vertices = i*VERTICES_IN_A_BLOB;
 	i = 0;
 	reflection = ctx->reflectionlist->reflections;
 	vertices = malloc(3*dw->gl_gen_num_vertices*sizeof(GLfloat));
 	normals = malloc(3*dw->gl_gen_num_vertices*sizeof(GLfloat));
-	do {
+	while ( reflection != NULL ) {
 		if ( reflection->type == REFLECTION_GENERATED ) {
 			DRAW_BLOB
 		}
 		reflection = reflection->next;
-	} while ( reflection != NULL );	
+	};	
 	if ( dw->gl_use_buffers ) {
 		glBindBufferARB(GL_ARRAY_BUFFER, dw->gl_gen_vertex_buffer);
 		glBufferDataARB(GL_ARRAY_BUFFER, 3*dw->gl_gen_num_vertices*sizeof(GLfloat), vertices, GL_STATIC_DRAW);
@@ -519,7 +524,7 @@ static void displaywindow_gl_create_list(DisplayWindow *dw) {
 	/* Plot the other reflections */
 	reflection = ctx->reflectionlist->reflections;
 	quadric = gluNewQuadric();
-	do {
+	while ( reflection != NULL ) {
 	
 		if ( reflection->type == REFLECTION_CENTRAL ) {
 		
@@ -583,61 +588,61 @@ static void displaywindow_gl_create_list(DisplayWindow *dw) {
 		
 		reflection = reflection->next;
 		
-	} while ( reflection != NULL );
+	};
 	
-	/* If this is an iterative prediction-refinement reconstruction, draw the unit cell */
-	if ( ctx->rmode == RECONSTRUCTION_PREDICTION ) {
+	/* Draw the reciprocal unit cell if one is available */
+	if ( ctx->cell ) {
 		glBegin(GL_LINE_STRIP);
 			glMaterialfv(GL_FRONT_AND_BACK, GL_AMBIENT_AND_DIFFUSE, red);
 			glMaterialfv(GL_FRONT_AND_BACK, GL_SPECULAR, black);
 			glMaterialf(GL_FRONT_AND_BACK, GL_SHININESS, 0.0);
 			glVertex3f(0.0, 0.0, 0.0);
-			glVertex3f(ctx->ipr_basis->a.x/1e9, ctx->ipr_basis->a.y/1e9, ctx->ipr_basis->a.z/1e9);
-			glVertex3f(ctx->ipr_basis->a.x/1e9 + ctx->ipr_basis->b.x/1e9,
-					ctx->ipr_basis->a.y/1e9 + ctx->ipr_basis->b.y/1e9,
-					ctx->ipr_basis->a.z/1e9 + ctx->ipr_basis->b.z/1e9);
-			glVertex3f(ctx->ipr_basis->b.x/1e9, ctx->ipr_basis->b.y/1e9, ctx->ipr_basis->b.z/1e9);
+			glVertex3f(ctx->cell->a.x/1e9, ctx->cell->a.y/1e9, ctx->cell->a.z/1e9);
+			glVertex3f(ctx->cell->a.x/1e9 + ctx->cell->b.x/1e9,
+					ctx->cell->a.y/1e9 + ctx->cell->b.y/1e9,
+					ctx->cell->a.z/1e9 + ctx->cell->b.z/1e9);
+			glVertex3f(ctx->cell->b.x/1e9, ctx->cell->b.y/1e9, ctx->cell->b.z/1e9);
 			glVertex3f(0.0, 0.0, 0.0);
-			glVertex3f(ctx->ipr_basis->c.x/1e9, ctx->ipr_basis->c.y/1e9, ctx->ipr_basis->c.z/1e9);
-			glVertex3f(ctx->ipr_basis->c.x/1e9 + ctx->ipr_basis->a.x/1e9,
-					ctx->ipr_basis->c.y/1e9 + ctx->ipr_basis->a.y/1e9,
-					ctx->ipr_basis->c.z/1e9 + ctx->ipr_basis->a.z/1e9);
-			glVertex3f(ctx->ipr_basis->c.x/1e9 + ctx->ipr_basis->a.x/1e9 + ctx->ipr_basis->b.x/1e9,
-					ctx->ipr_basis->c.y/1e9 + ctx->ipr_basis->a.y/1e9 + ctx->ipr_basis->b.y/1e9,
-					ctx->ipr_basis->c.z/1e9 + ctx->ipr_basis->a.z/1e9 + ctx->ipr_basis->b.z/1e9);
-			glVertex3f(ctx->ipr_basis->c.x/1e9 + ctx->ipr_basis->b.x/1e9,
-					ctx->ipr_basis->c.y/1e9 + ctx->ipr_basis->b.y/1e9,
-					ctx->ipr_basis->c.z/1e9 + ctx->ipr_basis->b.z/1e9);
-			glVertex3f(ctx->ipr_basis->c.x/1e9, ctx->ipr_basis->c.y/1e9, ctx->ipr_basis->c.z/1e9);
+			glVertex3f(ctx->cell->c.x/1e9, ctx->cell->c.y/1e9, ctx->cell->c.z/1e9);
+			glVertex3f(ctx->cell->c.x/1e9 + ctx->cell->a.x/1e9,
+					ctx->cell->c.y/1e9 + ctx->cell->a.y/1e9,
+					ctx->cell->c.z/1e9 + ctx->cell->a.z/1e9);
+			glVertex3f(ctx->cell->c.x/1e9 + ctx->cell->a.x/1e9 + ctx->cell->b.x/1e9,
+					ctx->cell->c.y/1e9 + ctx->cell->a.y/1e9 + ctx->cell->b.y/1e9,
+					ctx->cell->c.z/1e9 + ctx->cell->a.z/1e9 + ctx->cell->b.z/1e9);
+			glVertex3f(ctx->cell->c.x/1e9 + ctx->cell->b.x/1e9,
+					ctx->cell->c.y/1e9 + ctx->cell->b.y/1e9,
+					ctx->cell->c.z/1e9 + ctx->cell->b.z/1e9);
+			glVertex3f(ctx->cell->c.x/1e9, ctx->cell->c.y/1e9, ctx->cell->c.z/1e9);
 			
-			glVertex3f(ctx->ipr_basis->c.x/1e9 + ctx->ipr_basis->b.x/1e9,
-					ctx->ipr_basis->c.y/1e9 + ctx->ipr_basis->b.y/1e9,
-					ctx->ipr_basis->c.z/1e9 + ctx->ipr_basis->b.z/1e9);
-			glVertex3f(ctx->ipr_basis->b.x/1e9, ctx->ipr_basis->b.y/1e9, ctx->ipr_basis->b.z/1e9);
-			glVertex3f(ctx->ipr_basis->c.x/1e9 + ctx->ipr_basis->b.x/1e9,
-					ctx->ipr_basis->c.y/1e9 + ctx->ipr_basis->b.y/1e9,
-					ctx->ipr_basis->c.z/1e9 + ctx->ipr_basis->b.z/1e9);
-			glVertex3f(ctx->ipr_basis->c.x/1e9 + ctx->ipr_basis->a.x/1e9 + ctx->ipr_basis->b.x/1e9,
-					ctx->ipr_basis->c.y/1e9 + ctx->ipr_basis->a.y/1e9 + ctx->ipr_basis->b.y/1e9,
-					ctx->ipr_basis->c.z/1e9 + ctx->ipr_basis->a.z/1e9 + ctx->ipr_basis->b.z/1e9);
-			glVertex3f(ctx->ipr_basis->a.x/1e9 + ctx->ipr_basis->b.x/1e9,
-					ctx->ipr_basis->a.y/1e9 + ctx->ipr_basis->b.y/1e9,
-					ctx->ipr_basis->a.z/1e9 + ctx->ipr_basis->b.z/1e9);
-			glVertex3f(ctx->ipr_basis->c.x/1e9 + ctx->ipr_basis->a.x/1e9 + ctx->ipr_basis->b.x/1e9,
-					ctx->ipr_basis->c.y/1e9 + ctx->ipr_basis->a.y/1e9 + ctx->ipr_basis->b.y/1e9,
-					ctx->ipr_basis->c.z/1e9 + ctx->ipr_basis->a.z/1e9 + ctx->ipr_basis->b.z/1e9);
-			glVertex3f(ctx->ipr_basis->c.x/1e9 + ctx->ipr_basis->a.x/1e9,
-					ctx->ipr_basis->c.y/1e9 + ctx->ipr_basis->a.y/1e9,
-					ctx->ipr_basis->c.z/1e9 + ctx->ipr_basis->a.z/1e9);
-			glVertex3f(ctx->ipr_basis->a.x/1e9, ctx->ipr_basis->a.y/1e9, ctx->ipr_basis->a.z/1e9);
+			glVertex3f(ctx->cell->c.x/1e9 + ctx->cell->b.x/1e9,
+					ctx->cell->c.y/1e9 + ctx->cell->b.y/1e9,
+					ctx->cell->c.z/1e9 + ctx->cell->b.z/1e9);
+			glVertex3f(ctx->cell->b.x/1e9, ctx->cell->b.y/1e9, ctx->cell->b.z/1e9);
+			glVertex3f(ctx->cell->c.x/1e9 + ctx->cell->b.x/1e9,
+					ctx->cell->c.y/1e9 + ctx->cell->b.y/1e9,
+					ctx->cell->c.z/1e9 + ctx->cell->b.z/1e9);
+			glVertex3f(ctx->cell->c.x/1e9 + ctx->cell->a.x/1e9 + ctx->cell->b.x/1e9,
+					ctx->cell->c.y/1e9 + ctx->cell->a.y/1e9 + ctx->cell->b.y/1e9,
+					ctx->cell->c.z/1e9 + ctx->cell->a.z/1e9 + ctx->cell->b.z/1e9);
+			glVertex3f(ctx->cell->a.x/1e9 + ctx->cell->b.x/1e9,
+					ctx->cell->a.y/1e9 + ctx->cell->b.y/1e9,
+					ctx->cell->a.z/1e9 + ctx->cell->b.z/1e9);
+			glVertex3f(ctx->cell->c.x/1e9 + ctx->cell->a.x/1e9 + ctx->cell->b.x/1e9,
+					ctx->cell->c.y/1e9 + ctx->cell->a.y/1e9 + ctx->cell->b.y/1e9,
+					ctx->cell->c.z/1e9 + ctx->cell->a.z/1e9 + ctx->cell->b.z/1e9);
+			glVertex3f(ctx->cell->c.x/1e9 + ctx->cell->a.x/1e9,
+					ctx->cell->c.y/1e9 + ctx->cell->a.y/1e9,
+					ctx->cell->c.z/1e9 + ctx->cell->a.z/1e9);
+			glVertex3f(ctx->cell->a.x/1e9, ctx->cell->a.y/1e9, ctx->cell->a.z/1e9);
 		glEnd();
 	}
 
 	/* Zero plane */
 	glBegin(GL_QUADS);
 		glMaterialfv(GL_FRONT_AND_BACK, GL_AMBIENT_AND_DIFFUSE, glass);
-		glMaterialfv(GL_FRONT_AND_BACK, GL_SPECULAR, glass_spec);
-		glMaterialf(GL_FRONT_AND_BACK, GL_SHININESS, 50.0);
+		glMaterialfv(GL_FRONT_AND_BACK, GL_SPECULAR, black);
+		glMaterialf(GL_FRONT_AND_BACK, GL_SHININESS, 0.0);
 		glVertex3f(50, 50, 0.0);
 		glVertex3f(50, -50, 0.0);
 		glVertex3f(-50, -50, 0.0);
diff --git a/src/ipr.c b/src/ipr.c
index b74ecca..e16d5d8 100644
--- a/src/ipr.c
+++ b/src/ipr.c
@@ -95,7 +95,7 @@ static gint ipr_clicked(GtkWidget *widget, GdkEventButton *event, ControlContext
 	if ( ctx->ipr_cur_image == ctx->n_images ) ctx->ipr_cur_image = 0;
 	
 	imagedisplay_clear_circles(ctx->ipr_id);
-	reflection_clear_markers(ctx->reflectionlist);
+	reflectionlist_clear_markers(ctx->reflectionlist);
 	
 	refl = reproject_get_reflections(ctx->images[ctx->ipr_cur_image], &n, ctx->ipr_lat, ctx);
 	for ( j=0; j<n; j++ ) {
diff --git a/src/main.c b/src/main.c
index 578c2fe..f80da7e 100644
--- a/src/main.c
+++ b/src/main.c
@@ -33,12 +33,13 @@
 #include "mapping.h"
 #include "prealign.h"
 #include "control.h"
+#include "dirax.h"
 
 void main_do_reconstruction(ControlContext *ctx) {
 
 	int val = 0;
 	
-	if ( ctx->inputfiletype != INPUT_CACHE ) {
+	if ( (ctx->inputfiletype != INPUT_CACHE) && (ctx->inputfiletype != INPUT_DRX) ) {
 		prealign_sum_stack(ctx);
 		mapping_create(ctx);
 	}
@@ -64,7 +65,7 @@ static gint main_method_window_response(GtkWidget *method_window, gint response,
 
 	if ( response == GTK_RESPONSE_OK ) {
 	
-		int val = -1;
+		int val = 0;
 		
 		switch ( gtk_combo_box_get_active(GTK_COMBO_BOX(ctx->combo_algorithm)) ) {
 			case 0 : ctx->rmode = RECONSTRUCTION_MAPPING; break;
@@ -103,7 +104,10 @@ static gint main_method_window_response(GtkWidget *method_window, gint response,
 			val = mrc_read(ctx);
 		} else if ( ctx->inputfiletype == INPUT_CACHE ) {
 			ctx->reflectionlist = cache_load(ctx->filename);
-			val=0;
+			if ( !ctx->reflectionlist ) val = 1;
+		} else if ( ctx->inputfiletype == INPUT_DRX ) {
+			ctx->reflectionlist = dirax_load(ctx->filename);
+			if ( !ctx->reflectionlist ) val = 1;
 		}
 		
 		if ( val ) {
@@ -174,7 +178,12 @@ void main_method_dialog_open(ControlContext *ctx) {
 	gtk_table_attach_defaults(GTK_TABLE(table), ctx->checkbox_savecache, 1, 3, 4, 5);
 	
 	if ( ctx->inputfiletype == INPUT_CACHE ) {
-		gtk_combo_box_append_text(GTK_COMBO_BOX(ctx->combo_peaksearch), "Get from cache file");
+		gtk_combo_box_append_text(GTK_COMBO_BOX(ctx->combo_peaksearch), "3D coordinates from cache file");
+		gtk_widget_set_sensitive(GTK_WIDGET(ctx->combo_peaksearch), FALSE);
+		gtk_combo_box_set_active(GTK_COMBO_BOX(ctx->combo_peaksearch), 5);
+	}
+	if ( ctx->inputfiletype == INPUT_DRX ) {
+		gtk_combo_box_append_text(GTK_COMBO_BOX(ctx->combo_peaksearch), "3D coordinates from DirAx file");
 		gtk_widget_set_sensitive(GTK_WIDGET(ctx->combo_peaksearch), FALSE);
 		gtk_combo_box_set_active(GTK_COMBO_BOX(ctx->combo_peaksearch), 5);
 	}
@@ -194,6 +203,7 @@ int main(int argc, char *argv[]) {
 	ControlContext *ctx;
 	InputFileType type;
 	struct stat stat_buffer;
+	FILE *fh;
 
 	gtk_init(&argc, &argv);
 	
@@ -212,10 +222,17 @@ int main(int argc, char *argv[]) {
 		return 1;
 	}
 	
-	filename = basename(argv[1]);
+	filename = argv[1];
 	ctx = control_ctx_new();
 	type = INPUT_NONE;
 	
+	fh = fopen(filename, "r");
+	if ( !fh ) {
+		printf("Couldn't open file '%s'\n", filename);
+		return 1;
+	}
+	fclose(fh);
+	
 	if ( qdrp_is_qdrprc(filename) ) {
 		printf("QDRP input file detected.\n");
 		ctx->inputfiletype = INPUT_QDRP;
@@ -225,6 +242,9 @@ int main(int argc, char *argv[]) {
 	} else if ( cache_is_cachefile(filename) ) {
 		printf("Cached reflection file detected.\n");
 		ctx->inputfiletype = INPUT_CACHE;
+	} else if ( dirax_is_drxfile(filename) ) {
+		printf("Dirax input file detected.\n");
+		ctx->inputfiletype = INPUT_DRX;
 	} else {
 		fprintf(stderr, "Unrecognised input file type\n");
 		return 1;
diff --git a/src/reflections.c b/src/reflections.c
index 56871f0..dc3ce2f 100644
--- a/src/reflections.c
+++ b/src/reflections.c
@@ -39,7 +39,7 @@ ReflectionList *reflectionlist_new() {
 
 }
 
-void reflection_clear_markers(ReflectionList *reflectionlist) {
+void reflectionlist_clear_markers(ReflectionList *reflectionlist) {
 
 	Reflection *reflection = reflectionlist->reflections;
 	Reflection *prev = NULL;
@@ -70,7 +70,7 @@ void reflection_clear_markers(ReflectionList *reflectionlist) {
 	reflectionlist->last_reflection = prev;
 	
 }
-void reflection_clear(ReflectionList *reflectionlist) {
+void reflectionlist_clear(ReflectionList *reflectionlist) {
 
 	Reflection *reflection = reflectionlist->reflections;
 	do {
@@ -83,8 +83,8 @@ void reflection_clear(ReflectionList *reflectionlist) {
 
 }
 
-void reflection_free(ReflectionList *reflectionlist) {
-	reflection_clear(reflectionlist);
+void reflectionlist_free(ReflectionList *reflectionlist) {
+	reflectionlist_clear(reflectionlist);
 	free(reflectionlist);
 }
 
@@ -101,7 +101,7 @@ Reflection *reflection_add(ReflectionList *reflectionlist, double x, double y, d
 		reflectionlist->list_capped = 1;
 	}
 	
-	nearest = reflection_find_nearest_type(reflectionlist, x, y, z, type);
+	nearest = reflectionlist_find_nearest_type(reflectionlist, x, y, z, type);
 	if ( nearest && distance3d(x, y, z, nearest->x, nearest->y, nearest->z) < 0.1e9 ) return NULL;
 	
 	new_reflection = malloc(sizeof(Reflection));
@@ -226,7 +226,7 @@ void reflection_add_from_reflection(ReflectionList *reflectionlist, Reflection *
 	reflectionlist->n_reflections++;
 }
 
-double reflection_largest_g(ReflectionList *reflectionlist) {
+double reflectionlist_largest_g(ReflectionList *reflectionlist) {
 
 	double max = 0.0;
 	Reflection *reflection;
@@ -245,7 +245,7 @@ double reflection_largest_g(ReflectionList *reflectionlist) {
 
 }
 
-Reflection *reflection_find_nearest(ReflectionList *reflectionlist, double x, double y, double z) {
+Reflection *reflectionlist_find_nearest(ReflectionList *reflectionlist, double x, double y, double z) {
 
 	double max = +INFINITY;
 	Reflection *reflection;
@@ -268,7 +268,7 @@ Reflection *reflection_find_nearest(ReflectionList *reflectionlist, double x, do
 
 }
 
-Reflection *reflection_find_nearest_longer_unknown(ReflectionList *reflectionlist, double x, double y, double z, double min_distance) {
+Reflection *reflectionlist_find_nearest_longer_unknown(ReflectionList *reflectionlist, double x, double y, double z, double min_distance) {
 
 	double max = +INFINITY;
 	Reflection *reflection;
@@ -291,7 +291,7 @@ Reflection *reflection_find_nearest_longer_unknown(ReflectionList *reflectionlis
 
 }
 
-Reflection *reflection_find_nearest_type(ReflectionList *reflectionlist, double x, double y, double z, ReflectionType type) {
+Reflection *reflectionlist_find_nearest_type(ReflectionList *reflectionlist, double x, double y, double z, ReflectionType type) {
 
 	double max = +INFINITY;
 	Reflection *reflection;
@@ -316,7 +316,7 @@ Reflection *reflection_find_nearest_type(ReflectionList *reflectionlist, double
 }
 
 /* This destroys the lfom values in the input list */
-ReflectionList *reflection_sort_lfom(ReflectionList *in) {
+ReflectionList *reflectionlist_sort_lfom(ReflectionList *in) {
 
 	ReflectionList *out;
 	Reflection *best;
diff --git a/src/reflections.h b/src/reflections.h
index c753430..11fbdcf 100644
--- a/src/reflections.h
+++ b/src/reflections.h
@@ -59,22 +59,21 @@ typedef struct reflectionlist_struct {
 } ReflectionList;
 
 extern ReflectionList *reflectionlist_new(void);
-extern void reflection_clear(ReflectionList *reflectionlist);
-extern void reflection_clear_markers(ReflectionList *reflectionlist);
-extern void reflection_free(ReflectionList *reflectionlist);
-
-extern Reflection *reflection_add(ReflectionList *reflectionlist, double x, double y, double z, double intensity, ReflectionType type);
+extern void reflectionlist_clear(ReflectionList *reflectionlist);
+extern void reflectionlist_clear_markers(ReflectionList *reflectionlist);
+extern void reflectionlist_free(ReflectionList *reflectionlist);
 
 #include "control.h"
-
+extern Reflection *reflection_add(ReflectionList *reflectionlist, double x, double y, double z, double intensity, ReflectionType type);
 extern void reflection_add_from_dp(ControlContext *ctx, double x, double y, ImageRecord *imagerecord, double intensity);
 extern void reflection_add_from_reflection(ReflectionList *reflectionlist, Reflection *r);
-extern double reflection_largest_g(ReflectionList *reflectionlist);
 extern int reflection_map_to_space(ImageReflection *refl, double *ddx, double *ddy, double *ddz, double *twotheta);
-extern Reflection *reflection_find_nearest(ReflectionList *reflectionlist, double x, double y, double z);
-extern Reflection *reflection_find_nearest_longer_unknown(ReflectionList *reflectionlist, double x, double y, double z, double min_distance);
-extern Reflection *reflection_find_nearest_type(ReflectionList *reflectionlist, double x, double y, double z, ReflectionType type);
-extern ReflectionList *reflection_sort_lfom(ReflectionList *in);
+
+extern double reflectionlist_largest_g(ReflectionList *reflectionlist);
+extern Reflection *reflectionlist_find_nearest(ReflectionList *reflectionlist, double x, double y, double z);
+extern Reflection *reflectionlist_find_nearest_longer_unknown(ReflectionList *reflectionlist, double x, double y, double z, double min_distance);
+extern Reflection *reflectionlist_find_nearest_type(ReflectionList *reflectionlist, double x, double y, double z, ReflectionType type);
+extern ReflectionList *reflectionlist_sort_lfom(ReflectionList *in);
 
 #endif	/* REFLECTION_H */
 
diff --git a/src/utils.c b/src/utils.c
index 4618b44..1e750a8 100644
--- a/src/utils.c
+++ b/src/utils.c
@@ -14,6 +14,7 @@
 #include <gsl/gsl_matrix.h>
 
 #include "utils.h"
+#include "string.h"
 
 /* Return the MOST POSITIVE of two numbers */
 unsigned int biggest(signed int a, signed int b) {
@@ -81,3 +82,30 @@ void matrix_renormalise(gsl_matrix *m) {
 
 }
 
+size_t skipspace(const char *s) {
+
+	size_t i;
+	
+	for ( i=0; i<strlen(s); i++ ) {
+		if ( (s[i] != ' ') && (s[i] != '\t') ) return i;
+	}
+	
+	return strlen(s);
+
+}
+
+void chomp(char *s) {
+
+	size_t i;
+	
+	if ( !s ) return;
+	
+	for ( i=0; i<strlen(s); i++ ) {
+		if ( (s[i] == '\n') || (s[i] == '\r') ) {
+			s[i] = '\0';
+			return;
+		}
+	}
+	
+}
+
diff --git a/src/utils.h b/src/utils.h
index 1d268e6..9f8d636 100644
--- a/src/utils.h
+++ b/src/utils.h
@@ -18,6 +18,7 @@
 #endif
 
 #include <gsl/gsl_matrix.h>
+#include <math.h>
 
 extern unsigned int biggest(signed int a, signed int b);
 extern unsigned int smallest(signed int a, signed int b);
@@ -28,6 +29,8 @@ extern double angle_between_d(double x1, double y1, double z1, double x2, double
 extern double lambda(double voltage);
 extern void matrix_renormalise(gsl_matrix *m);
 extern double distance3d(double x1, double y1, double z1, double x2, double y2, double z2);
+extern size_t skipspace(const char *s);
+extern void chomp(char *s);
 
 #define rad2deg(a) ((a)*180/M_PI)
 #define deg2rad(a) ((a)*M_PI/180)