aboutsummaryrefslogtreecommitdiff
path: root/src/dirax.c
diff options
context:
space:
mode:
Diffstat (limited to 'src/dirax.c')
-rw-r--r--src/dirax.c99
1 files changed, 88 insertions, 11 deletions
diff --git a/src/dirax.c b/src/dirax.c
index c9c5f010..c1794764 100644
--- a/src/dirax.c
+++ b/src/dirax.c
@@ -32,6 +32,7 @@
#include "image.h"
#include "dirax.h"
+#include "utils.h"
typedef enum {
@@ -371,10 +372,13 @@ static int map_position(struct image *image, double x, double y,
}
+#define PEAK_WINDOW_SIZE (10)
+
static void search_peaks(struct image *image)
{
FILE *fh;
- int x, y;
+ int x, y, width;
+ int16_t *data;
fh = fopen("xfel.drx", "w");
if ( !fh ) {
@@ -383,23 +387,96 @@ static void search_peaks(struct image *image)
}
fprintf(fh, "%f\n", 0.5); /* Lie about the wavelength. */
- for ( x=0; x<image->width; x++ ) {
- for ( y=0; y<image->height; y++ ) {
+ data = image->data;
+ width = image->width;
+
+ for ( x=1; x<image->width-1; x++ ) {
+ for ( y=1; y<image->height-1; y++ ) {
+
+ double dx1, dx2, dy1, dy2;
+ double dxs, dys;
+ double grad;
+ int mask_x, mask_y;
+ int sx, sy;
+ double max;
+ unsigned int did_something = 1;
+
+ /* Get gradients */
+ dx1 = data[x+width*y] - data[(x+1)+width*y];
+ dx2 = data[(x-1)+width*y] - data[x+width*y];
+ dy1 = data[x+width*y] - data[(x+1)+width*(y+1)];
+ dy2 = data[x+width*(y-1)] - data[x+width*y];
+
+ /* Average gradient measurements from both sides */
+ dxs = ((dx1*dx1) + (dx2*dx2)) / 2;
+ dys = ((dy1*dy1) + (dy2*dy2)) / 2;
+
+ /* Calculate overall gradient */
+ grad = dxs + dys;
- int val;
+ if ( grad < 2000000 ) continue;
- val = image->data[y+image->height*(image->width-1-x)];
+ mask_x = x;
+ mask_y = y;
- if ( val > 1000 ) {
+ while ( (did_something)
+ && (distance(mask_x, mask_y, x, y)<50) ) {
- double rx, ry, rz;
+ max = data[mask_x+width*mask_y];
+ did_something = 0;
+
+ for ( sy=biggest(mask_y-PEAK_WINDOW_SIZE/2, 0);
+ sy<smallest(mask_y+PEAK_WINDOW_SIZE/2, image->height);
+ sy++ ) {
+ for ( sx=biggest(mask_x-PEAK_WINDOW_SIZE/2, 0);
+ sx<smallest(mask_x+PEAK_WINDOW_SIZE/2, width);
+ sx++ ) {
+
+ if ( data[sx+width*sy] > max ) {
+ max = data[sx+width*sy];
+ mask_x = sx;
+ mask_y = sy;
+ did_something = 1;
+ }
+
+ }
+ }
- /* Map and record reflection */
- map_position(image, x, y, &rx, &ry, &rz);
- fprintf(fh, "%10f %10f %10f %8f\n",
- rx/1e10, ry/1e10, rz/1e10, 1.0);
}
+ if ( !did_something ) {
+
+ //double d;
+ //int idx;
+
+ assert(mask_x<image->width);
+ assert(mask_y<image->height);
+ assert(mask_x>=0);
+ assert(mask_y>=0);
+
+ if ( distance(mask_x, mask_y, x, y) > 50.0 ) continue;
+
+ /* Check for a feature at exactly the
+ * same coordinates */
+ //image_feature_closest(flist, mask_x, mask_y,
+ // &d, &idx);
+
+ //if ( d > 1.0 ) {
+
+ double rx = 0.0;
+ double ry = 0.0;
+ double rz = 0.0;
+
+ /* Map and record reflection */
+ printf("%i %i\n", x, y);
+ map_position(image, x, y, &rx, &ry, &rz);
+ fprintf(fh, "%10f %10f %10f %8f\n",
+ rx/1e10, ry/1e10, rz/1e10, 1.0);
+ //}
+
+ }
+
+
}
}