aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--src/Makefile.am2
-rw-r--r--src/intensities.c139
-rw-r--r--src/intensities.h23
-rw-r--r--src/pattern_sim.c5
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++;