Attempt to match multiple cells

3 files changed, 98 insertions, 22 deletions

diff --git a/src/dirax.c b/src/dirax.c
index a9694f46..d18d1fa6 100644
--- a/src/dirax.c
+++ b/src/dirax.c
@@ -41,6 +41,8 @@
 
 #define DIRAX_VERBOSE 0
 
+#define MAX_DIRAX_CELL_CANDIDATES (5)
+
 
 typedef enum {
 	DIRAX_INPUT_NONE,
@@ -51,7 +53,8 @@ typedef enum {
 
 static void dirax_parseline(const char *line, struct image *image)
 {
-	int rf, i;
+	int rf, i, di, acl, acl_nh;
+	float d;
 
 	#if DIRAX_VERBOSE
 	char *copy;
@@ -77,7 +80,7 @@ static void dirax_parseline(const char *line, struct image *image)
 		if ( line[i] == 'R' ) rf = 1;
 		if ( (line[i] == 'D') && rf ) {
 			image->dirax_read_cell = 1;
-			image->indexed_cell = cell_new();
+			image->candidate_cells[image->ncells] = cell_new();
 			return;
 		}
 		i++;
@@ -86,31 +89,49 @@ static void dirax_parseline(const char *line, struct image *image)
 	/* Parse unit cell vectors as appropriate */
 	if ( image->dirax_read_cell == 1 ) {
 		/* First row of unit cell values */
-		float ax, ay, az, d;
+		float ax, ay, az;
 		sscanf(line, "%f %f %f %f %f %f", &d, &d, &d, &ax, &ay, &az);
-		cell_set_cartesian_a(image->indexed_cell,
+		cell_set_cartesian_a(image->candidate_cells[image->ncells],
 		                     ax*1e-10, ay*1e-10, az*1e-10);
 		image->dirax_read_cell++;
 		return;
 	} else if ( image->dirax_read_cell == 2 ) {
 		/* First row of unit cell values */
-		float bx, by, bz, d;
+		float bx, by, bz;
 		sscanf(line, "%f %f %f %f %f %f", &d, &d, &d, &bx, &by, &bz);
-		cell_set_cartesian_b(image->indexed_cell,
+		cell_set_cartesian_b(image->candidate_cells[image->ncells],
 		                     bx*1e-10, by*1e-10, bz*1e-10);
 		image->dirax_read_cell++;
 		return;
 	} else if ( image->dirax_read_cell == 3 ) {
 		/* First row of unit cell values */
-		float cx, cy, cz, d;
+		float cx, cy, cz;
 		sscanf(line, "%f %f %f %f %f %f", &d, &d, &d, &cx, &cy, &cz);
-		cell_set_cartesian_c(image->indexed_cell,
+		cell_set_cartesian_c(image->candidate_cells[image->ncells++],
 		                     cx*1e-10, cy*1e-10, cz*1e-10);
 		image->dirax_read_cell = 0;
 		return;
 	}
 
 	image->dirax_read_cell = 0;
+
+	if ( sscanf(line, "%i %i %f %f %f %f %f %f %i", &acl, &acl_nh,
+	                               &d, &d, &d, &d, &d, &d, &di) == 9 ) {
+		if ( acl_nh > image->best_acl_nh ) {
+
+			int i, found = 0;
+
+			for ( i=0; i<image->n_acls_tried; i++ ) {
+				if ( image->acls_tried[i] == acl ) found = 1;
+			}
+
+			if ( !found ) {
+				image->best_acl_nh = acl_nh;
+				image->best_acl = acl;
+			}
+
+		}
+	}
 }
 
 
@@ -135,6 +156,8 @@ static void dirax_sendline(const char *line, struct image *image)
 
 static void dirax_send_next(struct image *image)
 {
+	char tmp[32];
+
 	switch ( image->dirax_step ) {
 
 	case 1 :
@@ -166,14 +189,32 @@ static void dirax_send_next(struct image *image)
 		break;
 
 	case 8 :
-		/* Skip DirAx's 'acl' prompt */
-		dirax_sendline("\n", image);
+		if ( image->best_acl_nh == 0 ) {
+			STATUS("No more cells to try.\n");
+			dirax_sendline("exit", image);
+			break;
+		}
+		snprintf(tmp, 31, "acl %i\n", image->best_acl);
+		image->acls_tried[image->n_acls_tried++] = image->best_acl;
+		dirax_sendline(tmp, image);
 		break;
 
 	case 9 :
 		dirax_sendline("cell\n", image);
 		break;
 
+	case 10 :
+		if ( image->n_acls_tried == MAX_DIRAX_CELL_CANDIDATES ) {
+			STATUS("That's enough cells.\n");
+			dirax_sendline("exit\n", image);
+		} else {
+			/* Go back round for another cell */
+			image->best_acl_nh = 0;
+			image->dirax_step = 6;
+			dirax_sendline("acl\n", image);
+		}
+		break;
+
 	default:
 		dirax_sendline("exit\n", image);
 		return;
@@ -349,6 +390,8 @@ void run_dirax(struct image *image)
 
 	image->dirax_step = 1;	/* This starts the "initialisation" procedure */
 	image->dirax_read_cell = 0;
+	image->n_acls_tried = 0;
+	image->best_acl_nh = 0;
 
 	do {
 
diff --git a/src/image.h b/src/image.h
index 5715d570..dfd02d05 100644
--- a/src/image.h
+++ b/src/image.h
@@ -26,6 +26,9 @@
 #include "detector.h"
 
 
+#define MAX_CELL_CANDIDATES (32)
+
+
 /* Structure describing a feature in an image */
 struct imagefeature {
 
@@ -71,6 +74,8 @@ struct image {
 	double			*twotheta;
 	struct molecule		*molecule;
 	UnitCell		*indexed_cell;
+	UnitCell		*candidate_cells[MAX_CELL_CANDIDATES];
+	int			ncells;
 	struct detector		det;
 
 	struct quaternion	orientation;
@@ -94,6 +99,10 @@ struct image {
 	/* DirAx auto-indexing high-level stuff */
 	int			dirax_step;
 	int			dirax_read_cell;
+	int			best_acl;
+	int			best_acl_nh;
+	int			acls_tried[MAX_CELL_CANDIDATES];
+	int			n_acls_tried;
 
 };
 
diff --git a/src/index.c b/src/index.c
index 713380e7..155e11e4 100644
--- a/src/index.c
+++ b/src/index.c
@@ -61,7 +61,6 @@ void index_pattern(struct image *image, IndexingMethod indm, int no_match,
                    int verbose)
 {
 	int i;
-	UnitCell *new_cell = NULL;
 	int nc = 0;
 
 	/* Map positions to 3D */
@@ -83,27 +82,52 @@ void index_pattern(struct image *image, IndexingMethod indm, int no_match,
 
 	write_drx(image);
 
+	image->ncells = 0;
+
 	/* Index (or not) as appropriate */
 	if ( indm == INDEXING_NONE ) return;
 	if ( indm == INDEXING_DIRAX ) run_dirax(image);
 
-	if ( image->indexed_cell == NULL ) {
-		STATUS("No cell found.\n");
+	if ( image->ncells == 0 ) {
+		STATUS("No candidate cells found.\n");
+		return;
+	}
+
+	if ( no_match ) {
+		image->indexed_cell = image->candidate_cells[0];
+		if ( verbose ) {
+			STATUS("--------------------\n");
+			STATUS("The indexed cell (matching not performed):\n");
+			cell_print(image->indexed_cell);
+			STATUS("--------------------\n");
+		}
 		return;
-	} else if ( verbose ) {
-		STATUS("--------------------\n");
-		STATUS("The indexed cell (before matching):\n");
-		cell_print(image->indexed_cell);
-		STATUS("--------------------\n");
 	}
 
-	if ( !no_match ) {
-		new_cell = match_cell(image->indexed_cell,
+	for ( i=0; i<image->ncells; i++ ) {
+
+		UnitCell *new_cell = NULL;
+
+		if ( verbose ) {
+			STATUS("--------------------\n");
+			STATUS("Candidate cell %i (before matching):\n", i);
+			cell_print(image->candidate_cells[i]);
+			STATUS("--------------------\n");
+		}
+
+		new_cell = match_cell(image->candidate_cells[i],
 			              image->molecule->cell, verbose);
-		free(image->indexed_cell);
 		image->indexed_cell = new_cell;
 		if ( new_cell == NULL ) {
-			STATUS("Cell found, but not matched.\n");
+			STATUS("Cell %i not matched.\n", i);
+		} else {
+			goto done;
 		}
+
+	}
+
+done:
+	for ( i=0; i<image->ncells; i++ ) {
+		free(image->candidate_cells[i]);
 	}
 }