Commit f69a0db1 authored by Thomas White's avatar Thomas White
Browse files

Add intensity output

parent 1a834921
......@@ -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 \
......
/*
* 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");
}
/*
* 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 */
......@@ -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++;
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment