diff options
| -rw-r--r-- | src/Makefile.am | 2 | ||||
| -rw-r--r-- | src/intensities.c | 139 | ||||
| -rw-r--r-- | src/intensities.h | 23 | ||||
| -rw-r--r-- | src/pattern_sim.c | 5 |
4 files changed, 168 insertions, 1 deletions
diff --git a/src/Makefile.am b/src/Makefile.am index 7a8e9e43..47f6f002 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -3,7 +3,7 @@ bin_PROGRAMS = pattern_sim integr_sim AM_CFLAGS = -Wall pattern_sim_SOURCES = pattern_sim.c diffraction.c utils.c image.c cell.c \ - hdf5-file.c ewald.c detector.c sfac.c + hdf5-file.c ewald.c detector.c sfac.c intensities.c pattern_sim_LDADD = @LIBS@ integr_sim_SOURCES = integr_sim.c diffraction.c utils.c ewald.c sfac.c \ diff --git a/src/intensities.c b/src/intensities.c new file mode 100644 index 00000000..7630bcd8 --- /dev/null +++ b/src/intensities.c @@ -0,0 +1,139 @@ +/* + * intensities.c + * + * Extract intensities from patterns + * + * (c) 2007-2009 Thomas White <thomas.white@desy.de> + * + * pattern_sim - Simulate diffraction patterns from small crystals + * + */ + + +#include <stdlib.h> +#include <math.h> +#include <stdio.h> +#include <assert.h> + +#include "image.h" +#include "intensities.h" +#include "cell.h" +#include "sfac.h" + + +#define MAX_HITS (1024) + + +struct reflhit { + signed int h; + signed int k; + signed int l; + double min_distance; + int x; + int y; +}; + + +static int sum_nearby_points(uint16_t *data, int width, int x, int y) +{ + int dx, dy; + int intensity = 0; + + for ( dx=-3; dx<=3; dx++ ) { + for ( dy=-3; dy<=3; dy++ ) { + intensity += data[(x+dx) + width*(y+dy)]; + data[(x+dx) + width*(y+dy)] = 0; + } + } + + return intensity; +} + + +void output_intensities(struct image *image) +{ + int x, y; + double ax, ay, az; + double bx, by, bz; + double cx, cy, cz; + struct reflhit hits[MAX_HITS]; + int n_hits = 0; + int i; + + cell_get_cartesian(image->molecule->cell, &ax, &ay, &az, + &bx, &by, &bz, + &cx, &cy, &cz); + + for ( x=0; x<image->width; x++ ) { + for ( y=0; y<image->height; y++ ) { + + double hd, kd, ld; + signed int h, k, l; + struct threevec q; + double dist; + int found = 0; + int j; + + q = image->qvecs[x + image->width*y]; + + hd = q.u * ax + q.v * ay + q.w * az; + kd = q.u * bx + q.v * by + q.w * bz; + ld = q.u * cx + q.v * cy + q.w * cz; + + h = (signed int)rint(hd); + k = (signed int)rint(kd); + l = (signed int)rint(ld); + + dist = sqrt(pow(fmod(hd, 1.0), 2.0) + pow(fmod(kd, 1.0), 2.0) + + pow(fmod(ld, 1.0), 2.0)); + + if ( dist > 0.2 ) continue; + + for ( j=0; j<n_hits; j++ ) { + if ( (hits[j].h == h) && (hits[j].k == k) + && (hits[j].l == l) ) { + + if ( dist < hits[j].min_distance ) { + hits[j].min_distance = dist; + hits[j].x = x; + hits[j].y = y; + } + + found = 1; + + } + } + + if ( !found ) { + hits[n_hits].min_distance = dist; + hits[n_hits].x = x; + hits[n_hits].y = y; + hits[n_hits].h = h; + hits[n_hits].k = k; + hits[n_hits].l = l; + n_hits++; + assert(n_hits < MAX_HITS); + } + + } + } + + STATUS("Found %i reflections\n", n_hits); + + /* Explicit printf() used here (not normally allowed) because + * we really want to output to stdout */ + for ( i=0; i<n_hits; i++ ) { + + int intensity; + intensity = sum_nearby_points(image->data, image->width, + hits[i].x, hits[i].y); + + printf("%3i %3i %3i %6i (at %i,%i)\n", + hits[i].h, hits[i].k, hits[i].l, intensity, + hits[i].x, hits[i].y); + + } + + /* Blank line at end */ + printf("\n"); +} diff --git a/src/intensities.h b/src/intensities.h new file mode 100644 index 00000000..1a788f53 --- /dev/null +++ b/src/intensities.h @@ -0,0 +1,23 @@ +/* + * intensities.h + * + * Extract intensities from patterns + * + * (c) 2007-2009 Thomas White <thomas.white@desy.de> + * + * pattern_sim - Simulate diffraction patterns from small crystals + * + */ + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +#ifndef INTENSITIES_H +#define INTENSITIES_H + +#include "image.h" + +extern void output_intensities(struct image *image); + +#endif /* INTENSITIES_H */ diff --git a/src/pattern_sim.c b/src/pattern_sim.c index 3f3d93ea..e3c5ec8a 100644 --- a/src/pattern_sim.c +++ b/src/pattern_sim.c @@ -25,6 +25,7 @@ #include "utils.h" #include "hdf5-file.h" #include "detector.h" +#include "intensities.h" static void show_help(const char *s) @@ -224,6 +225,10 @@ int main(int argc, char *argv[]) get_diffraction(&image); record_image(&image); + if ( config_nearbragg ) { + output_intensities(&image); + } + snprintf(filename, 1023, "results/sim-%i.h5", number); number++; |