From 9c54c3d3f8824f5817c716eba5990c46489f9df1 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 3 Dec 2010 17:42:29 +0100 Subject: Fix integrate_peak() problems which caused incorrect peak positions --- src/peaks.c | 44 ++++++++++++++++++++++---------------------- 1 file changed, 22 insertions(+), 22 deletions(-) (limited to 'src/peaks.c') diff --git a/src/peaks.c b/src/peaks.c index 1d548aa0..263040dd 100644 --- a/src/peaks.c +++ b/src/peaks.c @@ -49,12 +49,6 @@ /* Window size for Zaefferer peak detection */ #define PEAK_WINDOW_SIZE (10) -/* Integration radius for peaks */ -#define INTEGRATION_RADIUS (10) - -/* Integration radius for background estimation */ -#define OUT_INTEGRATION_RADIUS (15) - static int in_streak(int x, int y) { @@ -188,26 +182,33 @@ static int cull_peaks(struct image *image) int integrate_peak(struct image *image, int xp, int yp, float *xc, float *yc, float *intensity, double *pbg, double *pmax, - int do_polar, int do_sa) + int do_polar, int do_sa, int centroid) { signed int x, y; - const int lim = INTEGRATION_RADIUS * INTEGRATION_RADIUS; - const int out_lim = OUT_INTEGRATION_RADIUS * OUT_INTEGRATION_RADIUS; + int lim, out_lim; double total = 0.0; int xct = 0; int yct = 0; double noise = 0.0; int noise_counts = 0; double max = 0.0; + struct panel *p = NULL; + + p = find_panel(image->det, xp, yp); + if ( p == NULL ) return 1; + if ( p->no_index ) return 1; + + lim = p->peak_sep/2; + out_lim = lim + 1; - for ( x=-OUT_INTEGRATION_RADIUS; x<+OUT_INTEGRATION_RADIUS; x++ ) { - for ( y=-OUT_INTEGRATION_RADIUS; y<+OUT_INTEGRATION_RADIUS; y++ ) { + for ( x=-out_lim; x<+out_lim; x++ ) { + for ( y=-out_lim; y<+out_lim; y++ ) { - struct panel *p = NULL; double val, sa, pix_area, Lsq, dsq, proj_area; float tt = 0.0; double phi, pa, pb, pol; uint16_t flags; + struct panel *p2; /* Outer mask radius */ if ( x*x + y*y > out_lim ) continue; @@ -215,6 +216,10 @@ int integrate_peak(struct image *image, int xp, int yp, if ( ((x+xp)>=image->width) || ((x+xp)<0) ) continue; if ( ((y+yp)>=image->height) || ((y+yp)<0) ) continue; + /* Strayed off one panel? */ + p2 = find_panel(image->det, x+xp, y+yp); + if ( p2 != p ) return 1; + /* Veto this peak if we tried to integrate in a bad region */ if ( image->flags != NULL ) { flags = image->flags[(x+xp)+image->width*(y+yp)]; @@ -233,11 +238,6 @@ int integrate_peak(struct image *image, int xp, int yp, if ( val > max ) max = val; - p = find_panel(image->det, x+xp, y+yp); - if ( p == NULL ) return 1; - - if ( p->no_index ) return 1; - if ( do_sa ) { /* Area of one pixel */ @@ -281,7 +281,7 @@ int integrate_peak(struct image *image, int xp, int yp, } /* The centroid is excitingly undefined if there is no intensity */ - if ( total != 0 ) { + if ( centroid && (total != 0) ) { *xc = (float)xct / total; *yc = (float)yct / total; *intensity = total; @@ -420,7 +420,7 @@ void search_peaks(struct image *image, float threshold, float min_gradient) * Don't bother doing polarisation/SA correction, because the * intensity of this peak is only an estimate at this stage. */ r = integrate_peak(image, mask_x, mask_y, - &fx, &fy, &intensity, NULL, NULL, 0, 0); + &fx, &fy, &intensity, NULL, NULL, 0, 0, 1); if ( r ) { /* Bad region - don't detect peak */ nrej_bad++; @@ -766,7 +766,7 @@ void output_intensities(struct image *image, UnitCell *cell, * revised coordinates. */ r = integrate_peak(image, f->x, f->y, &x, &y, &intensity, &bg, &max, - polar, sa); + polar, sa, 1); if ( r ) { /* The original peak (which also went * through integrate_peak(), but with @@ -788,7 +788,7 @@ void output_intensities(struct image *image, UnitCell *cell, image->cpeaks[i].x, image->cpeaks[i].y, &x, &y, &intensity, &bg, &max, - polar, sa); + polar, sa, 1); if ( r ) { /* Plain old ordinary peak veto */ n_veto++; @@ -809,7 +809,7 @@ void output_intensities(struct image *image, UnitCell *cell, image->cpeaks[i].x, image->cpeaks[i].y, &x, &y, &intensity, &bg, &max, polar, - sa); + sa, 0); if ( r ) { /* Plain old ordinary peak veto */ n_veto++; -- cgit v1.2.3