diff options
| -rw-r--r-- | src/Makefile.am | 10 | ||||
| -rw-r--r-- | src/index.c | 9 | ||||
| -rw-r--r-- | src/index.h | 1 | ||||
| -rw-r--r-- | src/indexamajig.c | 3 | ||||
| -rw-r--r-- | src/mosflm.c | 470 | ||||
| -rw-r--r-- | src/mosflm.h | 28 |
6 files changed, 516 insertions, 5 deletions
diff --git a/src/Makefile.am b/src/Makefile.am index 9d0100d5..dfd82d1f 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -25,8 +25,8 @@ process_hkl_LDADD = @LIBS@ indexamajig_SOURCES = indexamajig.c hdf5-file.c utils.c cell.c image.c \ peaks.c index.c filters.c diffraction.c detector.c \ - sfac.c dirax.c reflections.c templates.c symmetry.c \ - geometry.c thread-pool.c beam-parameters.c + sfac.c dirax.c mosflm.c reflections.c templates.c \ + symmetry.c geometry.c thread-pool.c beam-parameters.c indexamajig_LDADD = @LIBS@ if HAVE_OPENCL indexamajig_SOURCES += diffraction-gpu.c cl-utils.c @@ -76,8 +76,8 @@ endif reintegrate_SOURCES = reintegrate.c cell.c hdf5-file.c utils.c detector.c \ peaks.c image.c stream.c \ - index.c dirax.c templates.c geometry.c symmetry.c \ - thread-pool.c + index.c dirax.c mosflm.c templates.c geometry.c \ + symmetry.c thread-pool.c reintegrate_LDADD = @LIBS@ estimate_background_SOURCES = estimate_background.c stream.c utils.c cell.c @@ -87,7 +87,7 @@ INCLUDES = "-I$(top_srcdir)/data" EXTRA_DIST = cell.h hdf5-file.h image.h utils.h diffraction.h detector.h \ sfac.h reflections.h list_tmp.h statistics.h displaywindow.h \ - render.h hdfsee.h dirax.h peaks.h index.h filters.h \ + render.h hdfsee.h dirax.h mosflm.h peaks.h index.h filters.h \ diffraction-gpu.h cl-utils.h symmetry.h \ povray.h index-priv.h geometry.h templates.h render_hkl.h \ stream.h thread-pool.h beam-parameters.h post-refinement.h \ diff --git a/src/index.c b/src/index.c index bb725fcb..3e0a09b1 100644 --- a/src/index.c +++ b/src/index.c @@ -24,6 +24,7 @@ #include "utils.h" #include "peaks.h" #include "dirax.h" +#include "mosflm.h" #include "sfac.h" #include "detector.h" #include "index.h" @@ -50,6 +51,8 @@ IndexingPrivate *prepare_indexing(IndexingMethod indm, UnitCell *cell, return indexing_private(indm); case INDEXING_DIRAX : return indexing_private(indm); + case INDEXING_MOSFLM : + return indexing_private(indm); case INDEXING_TEMPLATE : return generate_templates(cell, filename, det, nominal_photon_energy); @@ -67,6 +70,9 @@ void cleanup_indexing(IndexingPrivate *priv) case INDEXING_DIRAX : free(priv); break; + case INDEXING_MOSFLM : + free(priv); + break; case INDEXING_TEMPLATE : free_templates(priv); } @@ -140,6 +146,9 @@ void index_pattern(struct image *image, UnitCell *cell, IndexingMethod indm, case INDEXING_DIRAX : run_dirax(image); break; + case INDEXING_MOSFLM : + run_mosflm(image); + break; case INDEXING_TEMPLATE : match_templates(image, ipriv); break; diff --git a/src/index.h b/src/index.h index 81549e00..29be9dcc 100644 --- a/src/index.h +++ b/src/index.h @@ -27,6 +27,7 @@ typedef enum { INDEXING_NONE, INDEXING_DIRAX, + INDEXING_MOSFLM, INDEXING_TEMPLATE } IndexingMethod; diff --git a/src/indexamajig.c b/src/indexamajig.c index 5a33e74d..8d762fa4 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -127,6 +127,7 @@ static void show_help(const char *s) " --indexing=<method> Use 'method' for indexing. Choose from:\n" " none : no indexing (default)\n" " dirax : invoke DirAx\n" +" mosflm : invoke MOSFLM (DPS)\n" " template : index by template matching\n" " -g. --geometry=<file> Get detector geometry from file.\n" " -b, --beam=<file> Get beam parameters from file (provides nominal\n" @@ -788,6 +789,8 @@ int main(int argc, char *argv[]) indm = INDEXING_NONE; } else if ( strcmp(indm_str, "dirax") == 0) { indm = INDEXING_DIRAX; + } else if ( strcmp(indm_str, "mosflm") == 0) { + indm = INDEXING_MOSFLM; } else if ( strcmp(indm_str, "template") == 0) { indm = INDEXING_TEMPLATE; } else { diff --git a/src/mosflm.c b/src/mosflm.c new file mode 100644 index 00000000..aed22720 --- /dev/null +++ b/src/mosflm.c @@ -0,0 +1,470 @@ +/* + * mosflm.c + * + * Invoke the DPS auto-indexing algorithm through MOSFLM + * + * This will actuaully run DirAx for now... will be fixed soon... + * + * (c) 2006-2010 Thomas White <taw@physics.org> + * + * Part of CrystFEL - crystallography with a FEL + * + */ + + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + + +#include <stdlib.h> +#include <stdio.h> +#include <math.h> +#include <string.h> +#include <unistd.h> +#include <sys/wait.h> +#include <fcntl.h> +#include <assert.h> +#include <sys/ioctl.h> + +#if HAVE_FORKPTY_LINUX +#include <pty.h> +#elif HAVE_FORKPTY_BSD +#include <util.h> +#endif + + +#include "image.h" +#include "mosflm.h" +#include "utils.h" +#include "sfac.h" +#include "peaks.h" + + + + + + +#define DIRAX_VERBOSE 0 + +#define MAX_DIRAX_CELL_CANDIDATES (5) + + +typedef enum { + DIRAX_INPUT_NONE, + DIRAX_INPUT_LINE, + DIRAX_INPUT_PROMPT, + DIRAX_INPUT_ACL +} DirAxInputType; + + +static void dirax_parseline(const char *line, struct image *image) +{ + int rf, i, di, acl, acl_nh; + float d; + + #if DIRAX_VERBOSE + 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]='\0'; + } + STATUS("DirAx: %s\n", copy); + free(copy); + #endif + + if ( strstr(line, "reflections from file") ) { + ERROR("DirAx can't understand this data.\n"); + return; + } + + /* 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 ) { + image->dirax_read_cell = 1; + image->candidate_cells[image->ncells] = cell_new(); + return; + } + i++; + } + + /* Parse unit cell vectors as appropriate */ + if ( image->dirax_read_cell == 1 ) { + /* First row of unit cell values */ + float ax, ay, az; + int r; + r = sscanf(line, "%f %f %f %f %f %f", + &d, &d, &d, &ax, &ay, &az); + if ( r != 6 ) { + ERROR("Couldn't understand cell line\n"); + image->dirax_read_cell = 0; + free(image->candidate_cells[image->ncells]); + return; + } + 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 ) { + /* Second row of unit cell values */ + float bx, by, bz; + int r; + r = sscanf(line, "%f %f %f %f %f %f", + &d, &d, &d, &bx, &by, &bz); + if ( r != 6 ) { + ERROR("Couldn't understand cell line\n"); + image->dirax_read_cell = 0; + free(image->candidate_cells[image->ncells]); + return; + } + 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 ) { + /* Third row of unit cell values */ + float cx, cy, cz; + int r; + r = sscanf(line, "%f %f %f %f %f %f", + &d, &d, &d, &cx, &cy, &cz); + if ( r != 6 ) { + ERROR("Couldn't understand cell line\n"); + image->dirax_read_cell = 0; + free(image->candidate_cells[image->ncells]); + return; + } + 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; + } + + } + } +} + + +static void dirax_sendline(const char *line, struct image *image) +{ + #if DIRAX_VERBOSE + char *copy; + int i; + + copy = strdup(line); + for ( i=0; i<strlen(copy); i++ ) { + if ( copy[i] == '\r' ) copy[i]='\0'; + if ( copy[i] == '\n' ) copy[i]='\0'; + } + STATUS("To DirAx: '%s'\n", copy); + free(copy); + #endif + + write(image->dirax_pty, line, strlen(line)); +} + + +static void dirax_send_next(struct image *image) +{ + char tmp[32]; + + switch ( image->dirax_step ) { + + case 1 : + dirax_sendline("\\echo off\n", image); + break; + + case 2 : + snprintf(tmp, 31, "read xfel-%i.drx\n", image->id); + dirax_sendline(tmp, image); + break; + + case 3 : + dirax_sendline("dmax 1000\n", image); + break; + + case 4 : + dirax_sendline("indexfit 2\n", image); + break; + + case 5 : + dirax_sendline("levelfit 1000\n", image); + break; + + case 6 : + dirax_sendline("go\n", image); + break; + + case 7 : + dirax_sendline("acl\n", image); + break; + + case 8 : + if ( image->best_acl_nh == 0 ) { + STATUS("No more cells to try.\n"); + /* At this point, DirAx is presenting its ACL prompt + * and waiting for a single number. Use an extra + * newline to choose automatic ACL selection before + * exiting. */ + dirax_sendline("\nexit\n", image); + break; + } + snprintf(tmp, 31, "%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 ) { + dirax_sendline("exit\n", image); + } else { + /* Go back round for another cell */ + image->best_acl_nh = 0; + image->dirax_step = 7; + dirax_sendline("acl\n", image); + } + break; + + default: + dirax_sendline("exit\n", image); + return; + + } + + image->dirax_step++; +} + + +static int dirax_readable(struct image *image) +{ + int rval; + int no_string = 0; + + rval = read(image->dirax_pty, image->dirax_rbuffer+image->dirax_rbufpos, + image->dirax_rbuflen-image->dirax_rbufpos); + + if ( (rval == -1) || (rval == 0) ) return 1; + + image->dirax_rbufpos += rval; + assert(image->dirax_rbufpos <= image->dirax_rbuflen); + + while ( (!no_string) && (image->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<image->dirax_rbufpos-1; i++ ) { + /* Means the last value looked at is rbufpos-2 */ + + /* Is there a prompt in the buffer? */ + if ( (i+7 <= image->dirax_rbufpos) + && (!strncmp(image->dirax_rbuffer+i, "Dirax> ", 7)) ) { + block_ready = 1; + type = DIRAX_INPUT_PROMPT; + break; + } + + /* How about an ACL prompt? */ + if ( (i+10 <= image->dirax_rbufpos) + && (!strncmp(image->dirax_rbuffer+i, "acl/auto [", 10)) ) { + block_ready = 1; + type = DIRAX_INPUT_ACL; + break; + } + + if ( (image->dirax_rbuffer[i] == '\r') + && (image->dirax_rbuffer[i+1] == '\n') ) { + block_ready = 1; + type = DIRAX_INPUT_LINE; + break; + } + + } + + if ( block_ready ) { + + unsigned int new_rbuflen; + unsigned int endbit_length; + char *block_buffer = NULL; + + switch ( type ) { + + case DIRAX_INPUT_LINE : + + block_buffer = malloc(i+1); + memcpy(block_buffer, image->dirax_rbuffer, i); + block_buffer[i] = '\0'; + + if ( block_buffer[0] == '\r' ) { + memmove(block_buffer, block_buffer+1, i); + } + + dirax_parseline(block_buffer, image); + free(block_buffer); + endbit_length = i+2; + break; + + case DIRAX_INPUT_PROMPT : + + dirax_send_next(image); + endbit_length = i+7; + break; + + case DIRAX_INPUT_ACL : + + dirax_send_next(image); + endbit_length = i+10; + break; + + default : + + /* Obviously, this never happens :) */ + ERROR("Unrecognised DirAx input mode! " + "I don't know how to understand DirAx\n"); + return 1; + + } + + /* Now the block's been parsed, it should be + * forgotten about */ + memmove(image->dirax_rbuffer, + image->dirax_rbuffer + endbit_length, + image->dirax_rbuflen - endbit_length); + + /* Subtract the number of bytes removed */ + image->dirax_rbufpos = image->dirax_rbufpos + - endbit_length; + new_rbuflen = image->dirax_rbuflen - endbit_length; + if ( new_rbuflen == 0 ) new_rbuflen = 256; + image->dirax_rbuffer = realloc(image->dirax_rbuffer, + new_rbuflen); + image->dirax_rbuflen = new_rbuflen; + + } else { + + if ( image->dirax_rbufpos==image->dirax_rbuflen ) { + + /* More buffer space is needed */ + image->dirax_rbuffer = realloc( + image->dirax_rbuffer, + image->dirax_rbuflen + 256); + image->dirax_rbuflen = image->dirax_rbuflen + + 256; + /* The new space gets used at the next + * read, shortly... */ + + } + no_string = 1; + + } + + } + + return 0; +} + + + + + +void run_mosflm(struct image *image) +{ + unsigned int opts; + int status; + int rval; + + printf("Mosflm is not yet implemented. Using DirAx for now.\n"); + + image->dirax_pid = forkpty(&image->dirax_pty, NULL, NULL, NULL); + if ( image->dirax_pid == -1 ) { + ERROR("Failed to fork for DirAx\n"); + return; + } + if ( image->dirax_pid == 0 ) { + + /* Child process: invoke DirAx */ + struct termios t; + + /* Turn echo off */ + tcgetattr(STDIN_FILENO, &t); + t.c_lflag &= ~(ECHO | ECHOE | ECHOK | ECHONL); + tcsetattr(STDIN_FILENO, TCSANOW, &t); + + execlp("dirax", "", (char *)NULL); + ERROR("Failed to invoke DirAx.\n"); + _exit(0); + + } + + image->dirax_rbuffer = malloc(256); + image->dirax_rbuflen = 256; + image->dirax_rbufpos = 0; + + /* Set non-blocking */ + opts = fcntl(image->dirax_pty, F_GETFL); + fcntl(image->dirax_pty, F_SETFL, opts | O_NONBLOCK); + + 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 { + + fd_set fds; + struct timeval tv; + int sval; + + FD_ZERO(&fds); + FD_SET(image->dirax_pty, &fds); + + tv.tv_sec = 10; + tv.tv_usec = 0; + + sval = select(image->dirax_pty+1, &fds, NULL, NULL, &tv); + + if ( sval == -1 ) { + ERROR("select() failed.\n"); + rval = 1; + } else if ( sval != 0 ) { + rval = dirax_readable(image); + } else { + ERROR("No response from DirAx..\n"); + rval = 1; + } + + } while ( !rval ); + + close(image->dirax_pty); + free(image->dirax_rbuffer); + wait(&status); + + return; +} diff --git a/src/mosflm.h b/src/mosflm.h new file mode 100644 index 00000000..8d09d8d5 --- /dev/null +++ b/src/mosflm.h @@ -0,0 +1,28 @@ +/* + * mosflm.h + * + * Invoke the DPS auto-indexing algorithm through MOSFLM + * + * This will actuaully run DirAx for now... will be fixed soon... + * + * (c) 2006-2010 Thomas White <taw@physics.org> + * + * Part of CrystFEL - crystallography with a FEL + * + */ + + +#ifndef MOSFLM_H +#define MOSFLM_H + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +#include "utils.h" + + +extern void run_mosflm(struct image *image); + + +#endif /* MOSFLM_H */ |