From f7d118ec0bc9e7dd964bbb9a8bc4fad46b8412d7 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Wed, 13 Jan 2010 10:46:44 +0100 Subject: Make DirAx use optional --- Makefile.am | 2 +- src/Makefile.am | 2 +- src/dirax.c | 59 ++------------------------------- src/dirax.h | 2 +- src/image.h | 5 +++ src/index.c | 98 +++++++++++++++++++++++++++++++++++++++++++++++++++++++ src/index.h | 25 ++++++++++++++ src/indexamajig.c | 15 ++++++--- 8 files changed, 145 insertions(+), 63 deletions(-) create mode 100644 src/index.c create mode 100644 src/index.h diff --git a/Makefile.am b/Makefile.am index 0406c2f4..c4954775 100644 --- a/Makefile.am +++ b/Makefile.am @@ -2,5 +2,5 @@ EXTRA_DIST = configure src/cell.h src/hdf5-file.h src/image.h src/relrod.h \ src/utils.h src/diffraction.h src/detector.h src/ewald.h \ src/sfac.h src/intensities.h src/reflections.h src/list_tmp.h \ src/statistics.h src/displaywindow.h src/render.h src/hdfsee.h \ - data/displaywindow.ui src/dirax.h src/peaks.h + data/displaywindow.ui src/dirax.h src/peaks.h src/index.h SUBDIRS = src data diff --git a/src/Makefile.am b/src/Makefile.am index 948d2fbf..8dd6af85 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -17,7 +17,7 @@ process_hkl_SOURCES = process_hkl.c sfac.c statistics.c cell.c utils.c \ process_hkl_LDADD = @LIBS@ indexamajig_SOURCES = indexamajig.c hdf5-file.c utils.c dirax.c cell.c image.c \ - intensities.c ewald.c peaks.c + intensities.c ewald.c peaks.c index.c indexamajig_LDADD = @LIBS@ if HAVE_GTK diff --git a/src/dirax.c b/src/dirax.c index 2db591f4..991a3813 100644 --- a/src/dirax.c +++ b/src/dirax.c @@ -323,63 +323,13 @@ static gboolean dirax_readable(GIOChannel *dirax, GIOCondition condition, } -static int map_position(struct image *image, double x, double y, - double *rx, double *ry, double *rz) -{ - /* "Input" space */ - double d; - - /* Angular description of reflection */ - double theta, psi, k; - - x -= image->x_centre; - y -= image->y_centre; - k = 1.0 / image->lambda; - - /* FIXME: Don't process lower CCD for now */ - if ( y < 0 ) return 0; - - if ( image->fmode == FORMULATION_CLEN ) { - - /* Convert pixels to metres */ - x /= image->resolution; - y /= image->resolution; /* Convert pixels to metres */ - d = sqrt((x*x) + (y*y)); - theta = atan2(d, image->camera_len); - - } else if (image->fmode == FORMULATION_PIXELSIZE ) { - - /* Convert pixels to metres^-1 */ - x *= image->pixel_size; - y *= image->pixel_size; /* Convert pixels to metres^-1 */ - d = sqrt((x*x) + (y*y)); - theta = atan2(d, k); - - } else { - ERROR("Unrecognised formulation mode in mapping_scale.\n"); - return -1; - } - - psi = atan2(y, x); - - *rx = k*sin(theta)*cos(psi); - *ry = k*sin(theta)*sin(psi); - *rz = k - k*cos(theta); - - return 0; -} - - -void index_pattern(struct image *image, int no_index, int dump_peaks) +void run_dirax(struct image *image, int no_index) { unsigned int opts; int saved_stderr; FILE *fh; int i; - /* Do peak search and splurge out 'xfel.drx' */ - search_peaks(image, dump_peaks); - fh = fopen("xfel.drx", "w"); if ( !fh ) { ERROR("Couldn't open temporary file xfel.drx\n"); @@ -390,14 +340,11 @@ void index_pattern(struct image *image, int no_index, int dump_peaks) for ( i=0; ifeatures); i++ ) { struct imagefeature *f; - double rx = 0.0; - double ry = 0.0; - double rz = 0.0; f = image_get_feature(image->features, i); - map_position(image, f->x, f->y, &rx, &ry, &rz); + fprintf(fh, "%10f %10f %10f %8f\n", - rx/1e10, ry/1e10, rz/1e10, 1.0); + f->rx/1e10, f->ry/1e10, f->rz/1e10, 1.0); } fclose(fh); diff --git a/src/dirax.h b/src/dirax.h index 4ebad677..58a0416e 100644 --- a/src/dirax.h +++ b/src/dirax.h @@ -18,7 +18,7 @@ #endif -extern void index_pattern(struct image *image, int no_index, int dump_peaks); +extern void run_dirax(struct image *image, int no_index); #endif /* DIRAX_H */ diff --git a/src/image.h b/src/image.h index 9649fe8e..7cf47761 100644 --- a/src/image.h +++ b/src/image.h @@ -46,6 +46,11 @@ struct imagefeature { /* Distance between this feature and its partner, if any. */ double partner_d; + /* Reciprocal space coordinates (m^-1 of course) of this feature */ + double rx; + double ry; + double rz; + }; /* An opaque type representing a list of image features */ diff --git a/src/index.c b/src/index.c new file mode 100644 index 00000000..c8d5dbe6 --- /dev/null +++ b/src/index.c @@ -0,0 +1,98 @@ +/* + * index.c + * + * Perform indexing (somehow) + * + * (c) 2006-2009 Thomas White + * + * Part of CrystFEL - crystallography with a FEL + * + */ + + +#ifdef HAVE_CONFIG_H +#include +#endif + +#include +#include +#include +#include +#include + +#include "image.h" +#include "utils.h" +#include "peaks.h" +#include "dirax.h" + + +static int map_position(struct image *image, double x, double y, + double *rx, double *ry, double *rz) +{ + /* "Input" space */ + double d; + + /* Angular description of reflection */ + double theta, psi, k; + + x -= image->x_centre; + y -= image->y_centre; + k = 1.0 / image->lambda; + + /* FIXME: Don't process lower CCD for now */ + if ( y < 0 ) return 0; + + if ( image->fmode == FORMULATION_CLEN ) { + + /* Convert pixels to metres */ + x /= image->resolution; + y /= image->resolution; /* Convert pixels to metres */ + d = sqrt((x*x) + (y*y)); + theta = atan2(d, image->camera_len); + + } else if (image->fmode == FORMULATION_PIXELSIZE ) { + + /* Convert pixels to metres^-1 */ + x *= image->pixel_size; + y *= image->pixel_size; /* Convert pixels to metres^-1 */ + d = sqrt((x*x) + (y*y)); + theta = atan2(d, k); + + } else { + ERROR("Unrecognised formulation mode in mapping_scale.\n"); + return -1; + } + + psi = atan2(y, x); + + *rx = k*sin(theta)*cos(psi); + *ry = k*sin(theta)*sin(psi); + *rz = k - k*cos(theta); + + return 0; +} + + +void index_pattern(struct image *image, int no_index, int dump_peaks, + int use_dirax) +{ + int i; + + /* Perform 'fine' peak search */ + search_peaks(image, dump_peaks); + + /* Map positions to 3D. FIXME: Handle lower detector */ + for ( i=0; ifeatures); i++ ) { + + struct imagefeature *f; + + f = image_get_feature(image->features, i); + map_position(image, f->x, f->y, &f->rx, &f->ry, &f->rz); + + } + + if ( use_dirax ) { + run_dirax(image, no_index); + return; + } +} diff --git a/src/index.h b/src/index.h new file mode 100644 index 00000000..ea0b4485 --- /dev/null +++ b/src/index.h @@ -0,0 +1,25 @@ +/* + * index.h + * + * Perform indexing (somehow) + * + * (c) 2006-2009 Thomas White + * + * Part of CrystFEL - crystallography with a FEL + * + */ + + +#ifndef INDEX_H +#define INDEX_H + +#ifdef HAVE_CONFIG_H +#include +#endif + + +extern void index_pattern(struct image *image, int no_index, int dump_peaks, + int use_dirax); + + +#endif /* INDEX_H */ diff --git a/src/indexamajig.c b/src/indexamajig.c index 7c2e1b5d..7f863a65 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -23,7 +23,7 @@ #include "utils.h" #include "hdf5-file.h" -#include "dirax.h" +#include "index.h" #include "intensities.h" #include "ewald.h" #include "peaks.h" @@ -39,8 +39,10 @@ static void show_help(const char *s) "\n" " -i, --input= Specify file containing list of images to process.\n" " '-' means stdin, which is the default.\n" -" --no-index Do everything else (including fine peak search),\n" -" but don't invoke the indexing program.\n" +" --no-index Do everything else (including fine peak search and\n" +" writing 'xfel.drx' if DirAx is being used), but\n" +" don't actually index.\n" +" --dirax Use DirAx for indexing.\n" "\n"); } @@ -55,6 +57,7 @@ int main(int argc, char *argv[]) int n_hits; int config_noindex = 0; int config_dumpfound = 0; + int config_dirax = 0; /* Long options */ const struct option longopts[] = { @@ -62,6 +65,7 @@ int main(int argc, char *argv[]) {"input", 1, NULL, 'i'}, {"no-index", 0, &config_noindex, 1}, {"dump-found-peaks", 0, &config_dumpfound, 1}, + {"dirax", 0, &config_dirax, 1}, {0, 0, NULL, 0} }; @@ -135,7 +139,10 @@ int main(int argc, char *argv[]) if ( fom > 0 ) { /* Calculate orientation matrix (by magic) */ - index_pattern(&image, config_noindex, config_dumpfound); + index_pattern(&image, config_noindex, config_dumpfound, + config_dirax); + + if ( image.molecule == NULL ) continue; /* View head-on (unit cell is tilted) */ image.orientation.x = 0.0; -- cgit v1.2.3