diff options
Diffstat (limited to 'src')
| -rw-r--r-- | src/Makefile.am | 8 | ||||
| -rw-r--r-- | src/cell.c | 148 | ||||
| -rw-r--r-- | src/cell.h | 55 | ||||
| -rw-r--r-- | src/image.c | 174 | ||||
| -rw-r--r-- | src/image.h | 121 | ||||
| -rw-r--r-- | src/main.c | 2 | ||||
| -rw-r--r-- | src/relrod.c | 180 | ||||
| -rw-r--r-- | src/relrod.h | 24 | ||||
| -rw-r--r-- | src/sim-main.c | 82 | ||||
| -rw-r--r-- | src/utils.c | 124 | ||||
| -rw-r--r-- | src/utils.h | 42 |
11 files changed, 958 insertions, 2 deletions
diff --git a/src/Makefile.am b/src/Makefile.am index ac173bb1..b480cedd 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -1,4 +1,8 @@ -bin_PROGRAMS = template_index -template_index_SOURCES = main.c +bin_PROGRAMS = template_index simulate_patterns + +template_index_SOURCES = main.c relrod.c utils.c image.c cell.c template_index_LDADD = @LIBS@ -lm AM_CFLAGS = -Wall -g + +simulate_patterns_SOURCES = sim-main.c relrod.c utils.c image.c cell.c +simulate_patterns_LDADD = @LIBS@ -lm diff --git a/src/cell.c b/src/cell.c new file mode 100644 index 00000000..d2495b0f --- /dev/null +++ b/src/cell.c @@ -0,0 +1,148 @@ +/* + * cell.c + * + * Unit Cell Calculations + * + * (c) 2007-2009 Thomas White <thomas.white@desy.de> + * + * template_index - Indexing diffraction patterns by template matching + * + */ + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +#include <math.h> +#include <stdlib.h> +#include <stdio.h> + +#include "cell.h" +#include "utils.h" + + +/* Update the cartesian representation from the crystallographic one */ +static void cell_update_cartesian(UnitCell *cell) +{ + double tmp, V, cosalphastar, cstar; + + if ( !cell ) return; + + /* a in terms of x, y and z + * +a (cryst) is defined to lie along +x (cart) */ + cell->ax = cell->a; + cell->ay = 0.0; + cell->az = 0.0; + + /* b in terms of x, y and z + * b (cryst) is defined to lie in the xy (cart) plane */ + cell->bx = cell->b*cos(cell->gamma); + cell->by = cell->b*sin(cell->gamma); + cell->bz = 0.0; + + tmp = cos(cell->alpha)*cos(cell->alpha) + - cos(cell->beta)*cos(cell->beta) + - cos(cell->gamma)*cos(cell->gamma) + + 2.0*cos(cell->alpha)*cos(cell->beta)*cos(cell->gamma); + V = cell->a * cell->b * cell->c * sqrt(1.0 - tmp); + + cosalphastar = cos(cell->beta)*cos(cell->gamma) - cos(cell->alpha); + cosalphastar /= sin(cell->beta)*sin(cell->gamma); + + cstar = (cell->a * cell->b * sin(cell->gamma))/V; + + /* c in terms of x, y and z */ + cell->cx = cell->c*cos(cell->beta); + cell->cy = -cell->c*sin(cell->beta)*cosalphastar; + cell->cz = 1.0/cstar; +} + + +/* Update the crystallographic representation from the cartesian one */ +static void cell_update_crystallographic(UnitCell *cell) +{ + if ( !cell ) return; + + cell->a = modulus(cell->ax, cell->ay, cell->az); + cell->b = modulus(cell->bx, cell->by, cell->bz); + cell->c = modulus(cell->cx, cell->cy, cell->cz); + + cell->alpha = angle_between(cell->bx, cell->by, cell->bz, + cell->cx, cell->cy, cell->cz); + cell->beta = angle_between(cell->ax, cell->ay, cell->az, + cell->cx, cell->cy, cell->cz); + cell->gamma = angle_between(cell->ax, cell->ay, cell->az, + cell->bx, cell->by, cell->bz); + + printf("a=%f nm\n", cell->a); + printf("b=%f nm\n", cell->b); + printf("c=%f nm\n", cell->c); + printf("alpha = %f deg\n", rad2deg(cell->alpha)); + printf(" beta = %f deg\n", rad2deg(cell->beta)); + printf("gamma = %f deg\n", rad2deg(cell->gamma)); +} + + +UnitCell *cell_new() +{ + UnitCell *cell; + + cell = malloc(sizeof(UnitCell)); + if ( !cell ) return NULL; + + cell->a = 1.0; + cell->b = 1.0; + cell->c = 1.0; + cell->alpha = M_PI_2; + cell->beta = M_PI_2; + cell->gamma = M_PI_2; + + cell_update_cartesian(cell); + + return cell; +} + + +void cell_set_parameters(UnitCell *cell, double a, double b, double c, + double alpha, double beta, double gamma) +{ + if ( !cell ) return; + + cell->a = a; + cell->b = b; + cell->c = c; + cell->alpha = alpha; + cell->beta = beta; + cell->gamma = gamma; + + cell_update_cartesian(cell); +} + + +void cell_set_cartesian(UnitCell *cell, + double ax, double ay, double az, + double bx, double by, double bz, + double cx, double cy, double cz) +{ + if ( !cell ) return; + + cell->ax = ax; cell->ay = ay; cell->az = az; + cell->bx = bx; cell->by = by; cell->bz = bz; + cell->cx = cx; cell->cy = cy; cell->cz = cz; + + cell_update_crystallographic(cell); +} + + +UnitCell *cell_new_from_parameters(double a, double b, double c, + double alpha, double beta, double gamma) +{ + UnitCell *cell; + + cell = cell_new(); + if ( !cell ) return NULL; + + cell_set_parameters(cell, a, b, c, alpha, beta, gamma); + + return cell; +} diff --git a/src/cell.h b/src/cell.h new file mode 100644 index 00000000..e99c88df --- /dev/null +++ b/src/cell.h @@ -0,0 +1,55 @@ +/* + * cell.h + * + * Unit Cell Calculations + * + * (c) 2007-2009 Thomas White <thomas.white@desy.de> + * + * template_index - Indexing diffraction patterns by template matching + * + */ + +#ifndef CELL_H +#define CELL_H + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +typedef struct { + + /* Crystallographic representation */ + double a; /* nm */ + double b; /* nm */ + double c; /* nm */ + double alpha; /* Radians */ + double beta; /* Radians */ + double gamma; /* Radians */ + + /* Cartesian representation */ + double ax; double bx; double cx; + double ay; double by; double cy; + double az; double bz; double cz; + + /* Cartesian representation of reciprocal axes */ + double axs; double bxs; double cxs; + double ays; double bys; double cys; + double azs; double bzs; double czs; + +} UnitCell; + +extern UnitCell *cell_new(void); + +/* Lengths in nm, angles in radians */ +extern UnitCell *cell_new_from_parameters(double a, double b, double c, + double alpha, double beta, double gamma); + +extern void cell_set_cartesian(UnitCell *cell, + double ax, double ay, double az, + double bx, double by, double bz, + double cx, double cy, double cz); + +extern void cell_set_parameters(UnitCell *cell, double a, double b, double c, + double alpha, double beta, double gamma); + +#endif /* CELL_H */ diff --git a/src/image.c b/src/image.c new file mode 100644 index 00000000..b3825b38 --- /dev/null +++ b/src/image.c @@ -0,0 +1,174 @@ +/* + * image.c + * + * Handle images and image features + * + * (c) 2006-2009 Thomas White <thomas.white@desy.de> + * + * template_index - Indexing diffraction patterns by template matching + * + */ + + +#define _GNU_SOURCE 1 +#include <stdlib.h> +#include <assert.h> +#include <math.h> + +#include "image.h" +#include "utils.h" + + +int image_add(ImageList *list, struct image *image) +{ + if ( list->images ) { + list->images = realloc(list->images, + (list->n_images+1)*sizeof(struct image)); + } else { + assert(list->n_images == 0); + list->images = malloc(sizeof(struct image)); + } + + /* Copy the metadata */ + list->images[list->n_images] = *image; + + /* Fill out some extra fields */ + list->images[list->n_images].features = NULL; + list->images[list->n_images].rflist = NULL; + + list->n_images++; + + return list->n_images - 1; +} + + +ImageList *image_list_new() +{ + ImageList *list; + + list = malloc(sizeof(ImageList)); + + list->n_images = 0; + list->images = NULL; + + return list; +} + + +void image_add_feature_reflection(ImageFeatureList *flist, double x, double y, + struct image *parent, double intensity) +{ + if ( flist->features ) { + flist->features = realloc(flist->features, + (flist->n_features+1)*sizeof(ImageFeature)); + } else { + assert(flist->n_features == 0); + flist->features = malloc(sizeof(ImageFeature)); + } + + flist->features[flist->n_features].x = x; + flist->features[flist->n_features].y = y; + flist->features[flist->n_features].intensity = intensity; + flist->features[flist->n_features].parent = parent; + flist->features[flist->n_features].partner = NULL; + flist->features[flist->n_features].partner_d = 0.0; + + flist->n_features++; + +} + + +void image_add_feature(ImageFeatureList *flist, double x, double y, + struct image *parent, double intensity) +{ + image_add_feature_reflection(flist, x, y, parent, intensity); +} + + +ImageFeatureList *image_feature_list_new() +{ + ImageFeatureList *flist; + + flist = malloc(sizeof(ImageFeatureList)); + + flist->n_features = 0; + flist->features = NULL; + + return flist; +} + + +void image_feature_list_free(ImageFeatureList *flist) +{ + if ( !flist ) return; + + if ( flist->features ) free(flist->features); + free(flist); +} + + +ImageFeature *image_feature_closest(ImageFeatureList *flist, double x, double y, + double *d, int *idx) +{ + int i; + double dmin = +HUGE_VAL; + int closest = 0; + + for ( i=0; i<flist->n_features; i++ ) { + + double d; + + d = distance(flist->features[i].x, flist->features[i].y, x, y); + + if ( d < dmin ) { + dmin = d; + closest = i; + } + + } + + if ( dmin < +HUGE_VAL ) { + *d = dmin; + *idx = closest; + return &flist->features[closest]; + } + + *d = +INFINITY; + return NULL; +} + + +ImageFeature *image_feature_second_closest(ImageFeatureList *flist, + double x, double y, double *d, + int *idx) +{ + int i; + double dmin = +HUGE_VAL; + int closest = 0; + double dfirst; + int idxfirst; + + image_feature_closest(flist, x, y, &dfirst, &idxfirst); + + for ( i=0; i<flist->n_features; i++ ) { + + double d; + + d = distance(flist->features[i].x, flist->features[i].y, x, y); + + if ( (d < dmin) && (i != idxfirst) ) { + dmin = d; + closest = i; + } + + } + + if ( dmin < +HUGE_VAL ) { + *d = dmin; + *idx = closest; + return &flist->features[closest]; + } + + return NULL; + +} diff --git a/src/image.h b/src/image.h new file mode 100644 index 00000000..f53f7d92 --- /dev/null +++ b/src/image.h @@ -0,0 +1,121 @@ +/* + * image.h + * + * Handle images and image features + * + * (c) 2006-2009 Thomas White <thomas.white@desy.de> + * + * template_index - Indexing diffraction patterns by template matching + * + */ + + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +#ifndef IMAGE_H +#define IMAGE_H + +#include <stdint.h> + + +typedef enum { + FORMULATION_CLEN, + FORMULATION_PIXELSIZE +} FormulationMode; + + +typedef struct imagefeature_struct { + + struct image *parent; + double x; + double y; + double intensity; + + /* Partner for this feature (in another feature list) or NULL */ + struct imagefeature_struct *partner; + + /* Distance between this feature and its partner, if any. */ + double partner_d; + + /* The reflection this was projected from, if any */ + struct reflection_struct *reflection; + +} ImageFeature; + + +typedef struct { + + ImageFeature *features; + int n_features; + +} ImageFeatureList; + + +struct image { + + uint16_t *data; + + /* Radians. Defines where the pattern lies in reciprocal space */ + double tilt; + + /* Radians. Defines where the pattern lies in reciprocal space */ + double omega; + + /* Image parameters can be given as camera length or pixel size. + * If FORMULATION_CLEN, then camera_len and resolution must be given. + * If FORMULATION_PIXELSIZE, then pixel_size only is needed.*/ + FormulationMode fmode; + double pixel_size; + double camera_len; + double resolution; /* pixels per metre */ + + /* Wavelength must always be given */ + double lambda; + + int width; + int height; + double x_centre; + double y_centre; + + ImageFeatureList *features; /* "Experimental" features */ + ImageFeatureList *rflist; /* "Predicted" features */ + +}; + + +typedef struct imagelist_struct { + + int n_images; + struct image *images; + +} ImageList; + + +extern ImageList *image_list_new(void); + +extern int image_add(ImageList *list, struct image *image); + +extern ImageFeatureList *image_feature_list_new(void); + +extern void image_feature_list_free(ImageFeatureList *flist); + +extern void image_add_feature(ImageFeatureList *flist, double x, double y, + struct image *parent, double intensity); + +extern void image_add_feature_reflection(ImageFeatureList *flist, + double x, double y, + struct image *parent, + double intensity); + +extern ImageFeature *image_feature_closest(ImageFeatureList *flist, + double x, double y, double *d, + int *idx); + +extern ImageFeature *image_feature_second_closest(ImageFeatureList *flist, + double x, double y, double *d, + int *idx); + + +#endif /* IMAGE_H */ @@ -60,6 +60,8 @@ int main(int argc, char *argv[]) return 1; } + printf("Generating templates...\n"); + printf("Input files (%i):\n", nin); for ( i=0; i<nin; i++ ) { printf("%6i: %s\n", i+1, in_files[i]); diff --git a/src/relrod.c b/src/relrod.c new file mode 100644 index 00000000..674f7971 --- /dev/null +++ b/src/relrod.c @@ -0,0 +1,180 @@ +/* + * relrod.c + * + * Calculate reflection positions via line-sphere intersection test + * + * (c) 2007-2009 Thomas White <thomas.white@desy.de> + * + * template_index - Indexing diffraction patterns by template matching + * + */ + + +#include <stdlib.h> +#include <math.h> +#include <stdio.h> + +#include "image.h" +#include "utils.h" +#include "cell.h" + + +static void mapping_rotate(double x, double y, double z, + double *ddx, double *ddy, double *ddz, + double omega, double tilt) +{ + double nx, ny, nz; + double x_temp, y_temp, z_temp; + + /* First: rotate image clockwise until tilt axis is aligned + * horizontally. */ + nx = x*cos(omega) + y*sin(omega); + ny = -x*sin(omega) + y*cos(omega); + nz = z; + + /* Now, tilt about the x-axis ANTICLOCKWISE around +x, i.e. the + * "wrong" way. This is because the crystal is rotated in the + * experiment, not the Ewald sphere. */ + x_temp = nx; y_temp = ny; z_temp = nz; + nx = x_temp; + ny = cos(tilt)*y_temp + sin(tilt)*z_temp; + nz = -sin(tilt)*y_temp + cos(tilt)*z_temp; + + /* Finally, reverse the omega rotation to restore the location of the + * image in 3D space */ + x_temp = nx; y_temp = ny; z_temp = nz; + nx = x_temp*cos(-omega) + y_temp*sin(-omega); + ny = -x_temp*sin(-omega) + y_temp*cos(-omega); + nz = z_temp; + + *ddx = nx; + *ddy = ny; + *ddz = nz; +} + + +void get_reflections(struct image *image, UnitCell *cell) +{ + ImageFeatureList *flist; + double smax = 0.2e9; + double tilt, omega, k; + double nx, ny, nz; /* "normal" vector */ + double kx, ky, kz; /* Electron wavevector ("normal" times 1/lambda */ + double ux, uy, uz; /* "up" vector */ + double rx, ry, rz; /* "right" vector */ + + flist = image_feature_list_new(); + + tilt = image->tilt; + omega = image->omega; + k = 1.0/image->lambda; + + /* Calculate the (normalised) incident electron wavevector */ + mapping_rotate(0.0, 0.0, 1.0, &nx, &ny, &nz, omega, tilt); + kx = nx / image->lambda; + ky = ny / image->lambda; + kz = nz / image->lambda; /* This is the centre of the Ewald sphere */ + + /* Determine where "up" is */ + mapping_rotate(0.0, 1.0, 0.0, &ux, &uy, &uz, omega, tilt); + + /* Determine where "right" is */ + mapping_rotate(1.0, 0.0, 0.0, &rx, &ry, &rz, omega, tilt); + + do { + + double xl, yl, zl; + double a, b, c; + double s1, s2, s, t; + double g_sq, gn; + + /* Get the coordinates of the reciprocal lattice point */ + // xl = reflection->x; + // yl = reflection->y; + // zl = reflection->z; + g_sq = modulus_squared(xl, yl, zl); + gn = xl*nx + yl*ny + zl*nz; + + /* Next, solve the relrod equation to calculate + * the excitation error */ + a = 1.0; + b = 2.0*(gn - k); + c = -2.0*gn*k + g_sq; + t = -0.5*(b + sign(b)*sqrt(b*b - 4.0*a*c)); + s1 = t/a; + s2 = c/t; + if ( fabs(s1) < fabs(s2) ) s = s1; else s = s2; + + /* Skip this reflection if s is large */ + if ( fabs(s) <= smax ) { + + double xi, yi, zi; + double gx, gy, gz; + double theta; + double x, y; + + /* Determine the intersection point */ + xi = xl + s*nx; yi = yl + s*ny; zi = zl + s*nz; + + /* Calculate Bragg angle */ + gx = xi - kx; + gy = yi - ky; + gz = zi - kz; /* This is the vector from the centre of + * the sphere to the intersection */ + theta = angle_between(-kx, -ky, -kz, gx, gy, gz); + + if ( theta > 0.0 ) { + + double dx, dy, psi; + + /* Calculate azimuth of point in image + * (anticlockwise from +x) */ + dx = xi*rx + yi*ry + zi*rz; + dy = xi*ux + yi*uy + zi*uz; + psi = atan2(dy, dx); + + /* Get image coordinates from polar + * representation */ + if ( image->fmode == FORMULATION_CLEN ) { + x = image->camera_len*tan(theta)*cos(psi); + y = image->camera_len*tan(theta)*sin(psi); + x *= image->resolution; + y *= image->resolution; + } else if ( image->fmode==FORMULATION_PIXELSIZE ) { + x = tan(theta)*cos(psi) / image->lambda; + y = tan(theta)*sin(psi) / image->lambda; + x /= image->pixel_size; + y /= image->pixel_size; + } else { + fprintf(stderr, + "Unrecognised formulation mode " + "in reproject_get_reflections\n"); + return NULL; + } + + x += image->x_centre; + y += image->y_centre; + + /* Sanity check */ + if ( (x>=0) && (x<image->width) + && (y>=0) && (y<image->height) ) { + + /* Record the reflection + * Intensity should be multiplied by + * relrod spike function except + * reprojected reflections aren't used + * quantitatively for anything + */ + image_add_feature_reflection(flist, x, + y, image, 1.0); + + } /* else it's outside the picture somewhere */ + + } /* else this is the central beam */ + + } + + } while ( 1 ); + + image->rflist = flist; +} diff --git a/src/relrod.h b/src/relrod.h new file mode 100644 index 00000000..1dac6057 --- /dev/null +++ b/src/relrod.h @@ -0,0 +1,24 @@ +/* + * relrod.h + * + * Calculate reflection positions via line-sphere intersection test + * + * (c) 2007-2009 Thomas White <thomas.white@desy.de> + * + * template_index - Indexing diffraction patterns by template matching + * + */ + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +#ifndef RELROD_H +#define RELROD_H + +#include "image.h" +#include "cell.h" + +extern void get_reflections(struct image *image, UnitCell *cell); + +#endif /* RELROD_H */ diff --git a/src/sim-main.c b/src/sim-main.c new file mode 100644 index 00000000..38ad45bc --- /dev/null +++ b/src/sim-main.c @@ -0,0 +1,82 @@ +/* + * sim-main.c + * + * Simulate test data + * + * (c) 2006-2009 Thomas White <thomas.white@desy.de> + * + * template_index - Indexing diffraction patterns by template matching + * + */ + + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +#include <stdarg.h> +#include <stdlib.h> +#include <stdio.h> +#include <string.h> +#include <unistd.h> + +#include "image.h" +#include "relrod.h" +#include "cell.h" +#include "utils.h" + + +static void main_show_help(const char *s) +{ + printf("Syntax: %s <file1.h5> <file2.h5> [...]\n\n", s); + printf("Index diffraction patterns\n\n"); + printf(" -h Display this help message\n"); +} + + +int main(int argc, char *argv[]) +{ + int c; + ImageList *list; + UnitCell *cell; + struct image image; + + while ((c = getopt(argc, argv, "h")) != -1) { + + switch ( c ) { + + case 'h' : { + main_show_help(argv[0]); + return 0; + } + + } + + } + + printf("Generating test data...\n"); + list = image_list_new(); + image.width = 512; + image.height = 512; + image.tilt = 0.0; + image.omega = 0.0; + image.fmode = FORMULATION_CLEN; + image.x_centre = 128; + image.y_centre = 128; + image.camera_len = 1.0; /* one metre */ + image.resolution = 5120; + image.lambda = 0.6e-9; + image.data = malloc(512*512*2); + image_add(list, &image); + + cell = cell_new_from_parameters(1.0, + 1.0, + 1.0, + deg2rad(90.0), + deg2rad(90.0), + deg2rad(90.0)); + + get_reflections(&image, cell); + + return 0; +} diff --git a/src/utils.c b/src/utils.c new file mode 100644 index 00000000..0c7d7ad4 --- /dev/null +++ b/src/utils.c @@ -0,0 +1,124 @@ +/* + * utils.c + * + * Utility stuff + * + * (c) 2006-2009 Thomas White <thomas.white@desy.de> + * + * template_index - Indexing diffraction patterns by template matching + * + */ + +#include <math.h> +#include <string.h> + +#include "utils.h" + + +/* Return the MOST POSITIVE of two numbers */ +unsigned int biggest(signed int a, signed int b) +{ + if ( a>b ) { + return a; + } + return b; +} + + +/* Return the LEAST POSITIVE of two numbers */ +unsigned int smallest(signed int a, signed int b) +{ + if ( a<b ) { + return a; + } + return b; +} + + +double distance(double x1, double y1, double x2, double y2) +{ + return sqrt((x2-x1)*(x2-x1) + (y2-y1)*(y2-y1)); +} + + +double modulus(double x, double y, double z) +{ + return sqrt(x*x + y*y + z*z); +} + + +double modulus_squared(double x, double y, double z) { + return x*x + y*y + z*z; +} + + +double distance3d(double x1, double y1, double z1, + double x2, double y2, double z2) +{ + return modulus(x1-x2, y1-y2, z1-z2); +} + + +/* Angle between two vectors. Answer in radians */ +double angle_between(double x1, double y1, double z1, + double x2, double y2, double z2) +{ + double mod1 = modulus(x1, y1, z1); + double mod2 = modulus(x2, y2, z2); + return acos( (x1*x2 + y1*y2 + z1*z2) / (mod1*mod2) ); +} + + +/* As above, answer in degrees */ +double angle_between_d(double x1, double y1, double z1, + double x2, double y2, double z2) +{ + return rad2deg(angle_between(x1, y1, z1, x2, y2, z2)); +} + + +/* Wavelength of an electron (in m) given accelerating potential (in V) */ +double lambda(double V) +{ + double m = 9.110E-31; + double h = 6.625E-34; + double e = 1.60E-19; + double c = 2.998E8; + + return h / sqrt(2*m*e*V*(1+((e*V) / (2*m*c*c)))); +} + + +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; + } + } +} + + +int sign(double a) +{ + if ( a < 0 ) return -1; + if ( a > 0 ) return +1; + return 0; +} diff --git a/src/utils.h b/src/utils.h new file mode 100644 index 00000000..ef37bc73 --- /dev/null +++ b/src/utils.h @@ -0,0 +1,42 @@ +/* + * utils.h + * + * Utility stuff + * + * (c) 2006-2009 Thomas White <thomas.white@desy.de> + * + * template_index - Indexing diffraction patterns by template matching + * + */ + +#ifndef UTILS_H +#define UTILS_H + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +#include <math.h> + +extern unsigned int biggest(signed int a, signed int b); +extern unsigned int smallest(signed int a, signed int b); +extern double distance(double x1, double y1, double x2, double y2); +extern double modulus(double x, double y, double z); +extern double modulus_squared(double x, double y, double z); +extern double angle_between(double x1, double y1, double z1, + double x2, double y2, double z2); +extern double angle_between_d(double x1, double y1, double z1, + double x2, double y2, double z2); +extern double lambda(double voltage); +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); +extern int sign(double a); + +#define rad2deg(a) ((a)*180/M_PI) +#define deg2rad(a) ((a)*M_PI/180) + +#define is_odd(a) ((a)%2==1) + +#endif /* UTILS_H */ |