aboutsummaryrefslogtreecommitdiff
path: root/src
diff options
context:
space:
mode:
authorThomas White <taw@physics.org>2010-05-26 18:51:23 +0200
committerThomas White <taw@physics.org>2010-05-26 18:51:23 +0200
commit8941c1ba4ba4bc4ed00b41db371fb9d75ac137ca (patch)
tree6e33386eee5bd51ec52dd533e8af26234b48d277 /src
parent2ea2488a3752ccb12bae10ae5bf8446f2f3f99f3 (diff)
Add peak sanity check
Diffstat (limited to 'src')
-rw-r--r--src/image.h12
-rw-r--r--src/indexamajig.c16
-rw-r--r--src/peaks.c77
-rw-r--r--src/peaks.h1
4 files changed, 86 insertions, 20 deletions
diff --git a/src/image.h b/src/image.h
index 0479debc..85265028 100644
--- a/src/image.h
+++ b/src/image.h
@@ -67,6 +67,16 @@ struct rvec
};
+struct reflhit {
+ signed int h;
+ signed int k;
+ signed int l;
+ double min_distance;
+ int x;
+ int y;
+};
+
+
/* Structure describing an image */
struct image {
@@ -77,6 +87,8 @@ struct image {
int ncells;
struct detector det;
const char *filename;
+ struct reflhit *hits;
+ int n_hits;
int id;
diff --git a/src/indexamajig.c b/src/indexamajig.c
index c21dc22f..db433369 100644
--- a/src/indexamajig.c
+++ b/src/indexamajig.c
@@ -67,6 +67,7 @@ struct process_args
struct process_result
{
int hit;
+ int peaks_sane;
};
@@ -194,6 +195,10 @@ static struct image *get_simage(struct image *template, int alternate)
image->indexed_cell = template->indexed_cell;
image->f0 = template->f0;
+ /* Prevent muppetry */
+ image->hits = NULL;
+ image->n_hits = 0;
+
return image;
}
@@ -252,6 +257,8 @@ static void *process_image(void *pargsv)
image.indexed_cell = NULL;
image.id = pargs->id;
image.filename = filename;
+ image.hits = NULL;
+ image.n_hits = 0;
STATUS("Processing '%s'\n", image.filename);
@@ -318,6 +325,15 @@ static void *process_image(void *pargsv)
/* No cell at this point? Then we're done. */
if ( image.indexed_cell == NULL ) goto done;
+ /* Sanity check */
+ if ( !peak_sanity_check(&image, image.indexed_cell) ) {
+ STATUS("Failed peak sanity check.\n");
+ result->peaks_sane = 0;
+ goto done;
+ } else {
+ result->peaks_sane = 1;
+ }
+
/* Measure intensities if requested */
if ( config_nearbragg ) {
/* Use original data (temporarily) */
diff --git a/src/peaks.c b/src/peaks.c
index 20779bb3..fe5dcbea 100644
--- a/src/peaks.c
+++ b/src/peaks.c
@@ -44,15 +44,6 @@
/* Degree of polarisation of X-ray beam */
#define POL (1.0)
-struct reflhit {
- signed int h;
- signed int k;
- signed int l;
- double min_distance;
- int x;
- int y;
-};
-
#define PEAK_WINDOW_SIZE (10)
#define MAX_PEAKS (2048)
@@ -385,23 +376,17 @@ void dump_peaks(struct image *image, pthread_mutex_t *mutex)
}
-void output_intensities(struct image *image, UnitCell *cell,
- pthread_mutex_t *mutex, int unpolar)
+static int find_projected_peaks(struct image *image, UnitCell *cell)
{
int x, y;
double ax, ay, az;
double bx, by, bz;
double cx, cy, cz;
- double a, b, c, al, be, ga;
- double asx, asy, asz;
- double bsx, bsy, bsz;
- double csx, csy, csz;
- struct reflhit hits[MAX_HITS];
+ struct reflhit *hits;
int n_hits = 0;
- int i;
- int n_found;
- int n_indclose = 0;
- int n_foundclose = 0;
+
+ hits = malloc(sizeof(struct reflhit)*MAX_HITS);
+ if ( hits == NULL ) return 0;
cell_get_cartesian(cell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz);
@@ -462,6 +447,58 @@ void output_intensities(struct image *image, UnitCell *cell,
}
STATUS("Found %i reflections\n", n_hits);
+ image->hits = hits;
+ image->n_hits = n_hits;
+
+ return n_hits;
+}
+
+
+int peak_sanity_check(struct image *image, UnitCell *cell)
+{
+ int i;
+ const int n_hits = image->n_hits;
+ const struct reflhit *hits = image->hits;
+ int n_sane = 0;
+
+ find_projected_peaks(image, cell);
+ if ( image->n_hits == 0 ) return 0; /* Failed sanity check: no peaks */
+
+ for ( i=0; i<n_hits; i++ ) {
+
+ double d;
+ int idx;
+ struct imagefeature *f;
+
+ f = image_feature_closest(image->features, hits[i].x, hits[i].y,
+ &d, &idx);
+ if ( (f != NULL) && (d < PEAK_CLOSE) ) {
+ n_sane++;
+ }
+
+ }
+
+ if ( (float)n_sane / (float)n_hits < 0.8 ) return 0;
+
+ return 1;
+}
+
+
+void output_intensities(struct image *image, UnitCell *cell,
+ pthread_mutex_t *mutex, int unpolar)
+{
+ int i;
+ int n_found;
+ int n_indclose = 0;
+ int n_foundclose = 0;
+ double asx, asy, asz;
+ double bsx, bsy, bsz;
+ double csx, csy, csz;
+ double a, b, c, al, be, ga;
+ int n_hits = image->n_hits;
+ struct reflhit *hits = image->hits;
+
+ if ( image->n_hits == 0 ) return;
/* Get exclusive access to the output stream if necessary */
if ( mutex != NULL ) pthread_mutex_lock(mutex);
diff --git a/src/peaks.h b/src/peaks.h
index fa912455..10cf5b0a 100644
--- a/src/peaks.h
+++ b/src/peaks.h
@@ -23,5 +23,6 @@ extern void search_peaks(struct image *image);
extern void dump_peaks(struct image *image, pthread_mutex_t *mutex);
extern void output_intensities(struct image *image, UnitCell *cell,
pthread_mutex_t *mutex, int unpolar);
+extern int peak_sanity_check(struct image *image, UnitCell *cell);
#endif /* PEAKS_H */