From 652872f38e4f1738366b16c350f27efd8cf7c99d Mon Sep 17 00:00:00 2001
From: Thomas White <taw@physics.org>
Date: Mon, 15 Nov 2021 16:14:46 +0100
Subject: pair_peaks: Use reciprocal space distance instead of pixel distance

The reciprocal space distance limit has been set as one third of the
smallest inter-Bragg spacing.

Fixes: https://gitlab.desy.de/thomas.white/crystfel/-/issues/38
---
 libcrystfel/src/predict-refine.c | 26 ++++++++++++++++++++++----
 1 file changed, 22 insertions(+), 4 deletions(-)

(limited to 'libcrystfel')

diff --git a/libcrystfel/src/predict-refine.c b/libcrystfel/src/predict-refine.c
index 38aff911..3c1fabdd 100644
--- a/libcrystfel/src/predict-refine.c
+++ b/libcrystfel/src/predict-refine.c
@@ -179,11 +179,14 @@ static int pair_peaks(struct image *image, Crystal *cr,
 	double cx, cy, cz;
 	double dx, dy;
 	RefList *all_reflist;
+	double lowest_one_over_d;
 
 	all_reflist = reflist_new();
 	cell_get_cartesian(crystal_get_cell(cr),
 	                   &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz);
 
+	lowest_one_over_d = lowest_reflection(crystal_get_cell(cr));
+
 	crystal_get_det_shift(cr, &dx, &dy);
 
 	/* First, create a RefList containing the most likely indices for each
@@ -255,17 +258,32 @@ static int pair_peaks(struct image *image, Crystal *cr,
 	 * good pairings */
 	for ( i=0; i<n; i++ ) {
 
-		double fs, ss, pd;
+		double fs, ss;
 		signed int h, k, l;
+		int pnl;
+		double refl_r[3];
+		double pk_r[3];
 		Reflection *refl = rps[i].refl;
 
 		get_indices(refl, &h, &k, &l);
 
 		/* Is the supposed reflection anywhere near the peak? */
 		get_detector_pos(refl, &fs, &ss);
-		pd = pow(fs - rps[i].peak->fs, 2.0)
-		   + pow(ss - rps[i].peak->ss, 2.0);
-		if ( pd > 10.0 * 10.0 ) continue; /* FIXME Hardcoded distance (GitLab #38) */
+
+		pnl = get_panel_number(refl);
+		detgeom_transform_coords(&image->detgeom->panels[pnl],
+		                         fs, ss,
+		                         image->lambda, dx, dy, refl_r);
+		detgeom_transform_coords(&image->detgeom->panels[pnl],
+		                         rps[i].peak->fs, rps[i].peak->ss,
+		                         image->lambda, dx, dy, pk_r);
+
+		if ( modulus(refl_r[0] - pk_r[0],
+		             refl_r[1] - pk_r[1],
+		             refl_r[2] - pk_r[2]) > lowest_one_over_d / 3.0 )
+		{
+			continue;
+		}
 
 		rps[n_acc] = rps[i];
 		rps[n_acc].refl = reflection_new(h, k, l);
-- 
cgit v1.2.3