Fix handling of overlapping reflections

5 files changed, 33 insertions, 43 deletions

diff --git a/libcrystfel/src/integration.c b/libcrystfel/src/integration.c
index 4d9a1219..07a3e393 100644
--- a/libcrystfel/src/integration.c
+++ b/libcrystfel/src/integration.c
@@ -159,6 +159,7 @@ struct intcontext
 	int w;
 	enum boxmask_val *bm;  /* Box mask */
 	struct image *image;
+	int **masks;  /* Peak location mask from make_BgMask() */
 
 	struct peak_box *boxes;
 	int n_boxes;
@@ -872,6 +873,20 @@ static int check_box(struct intcontext *ic, struct peak_box *bx, int *sat)
 			bx->bm[p+ic->w*q] = BM_BH;
 		}
 
+		/* If this is a background pixel, it shouldn't contain any
+		 * pixels which are in the peak region of ANY reflection */
+		if ( (ic->masks != NULL) && (bx->bm[p+ic->w*q] == BM_BG)
+		  && (ic->masks[bx->pn][fs + bx->p->w*ss] > 0) ) {
+			bx->bm[p+ic->w*q] = BM_BH;
+		}
+
+		/* If this is a peak pixel, it shouldn't contain any pixels
+		 * which are in the peak region of more than one reflection */
+		if ( (ic->masks != NULL) && (bx->bm[p+ic->w*q] == BM_PK)
+		  && (ic->masks[bx->pn][fs + bx->p->w*ss] > 1) ) {
+			bx->bm[p+ic->w*q] = BM_BH;
+		}
+
 		if ( (bx->bm[p+ic->w*q] != BM_IG)
 		  && (bx->bm[p+ic->w*q] != BM_BH)
 		  && (boxi(ic, bx, p, q) > bx->p->max_adu) ) {
@@ -1476,7 +1491,7 @@ static void integrate_prof2d(IntegrationMethod meth, PartialityModel pmodel,
                              Crystal *cr, struct image *image, IntDiag int_diag,
                              signed int idh, signed int idk, signed int idl,
                              double ir_inn, double ir_mid, double ir_out,
-                             int results_pipe)
+                             int results_pipe, int **masks)
 {
 	RefList *list;
 	UnitCell *cell;
@@ -1500,6 +1515,7 @@ static void integrate_prof2d(IntegrationMethod meth, PartialityModel pmodel,
 	ic.int_diag_h = idh;
 	ic.int_diag_k = idk;
 	ic.int_diag_l = idl;
+	ic.masks = masks;
 	if ( init_intcontext(&ic) ) {
 		ERROR("Failed to initialise integration.\n");
 		return;
@@ -1733,7 +1749,7 @@ static void integrate_rings(IntegrationMethod meth, PartialityModel pmodel,
                             Crystal *cr, struct image *image, IntDiag int_diag,
                             signed int idh, signed int idk, signed int idl,
                             double ir_inn, double ir_mid, double ir_out,
-                            int results_pipe)
+                            int results_pipe, int **masks)
 {
 	RefList *list;
 	Reflection *refl;
@@ -1762,6 +1778,7 @@ static void integrate_rings(IntegrationMethod meth, PartialityModel pmodel,
 	ic.int_diag_k = idk;
 	ic.int_diag_l = idl;
 	ic.meth = meth;
+	ic.masks = masks;
 	if ( init_intcontext(&ic) ) {
 		ERROR("Failed to initialise integration.\n");
 		return;
@@ -1818,6 +1835,11 @@ void integrate_all_4(struct image *image, IntegrationMethod meth,
                      int results_pipe)
 {
 	int i;
+	int *masks[image->det->n_panels];
+
+	for ( i=0; i<image->det->n_panels; i++ ) {
+		masks[i] = make_BgMask(image, &image->det->panels[i], ir_inn);
+	}
 
 	for ( i=0; i<image->n_crystals; i++ ) {
 
@@ -1833,7 +1855,7 @@ void integrate_all_4(struct image *image, IntegrationMethod meth,
 			integrate_rings(meth, pmodel, cr, image,
 			                int_diag, idh, idk, idl,
 			                ir_inn, ir_mid, ir_out,
-			                results_pipe);
+			                results_pipe, masks);
 			res = estimate_resolution(crystal_get_cell(cr),
 			                          image->features);
 			break;
@@ -1842,7 +1864,7 @@ void integrate_all_4(struct image *image, IntegrationMethod meth,
 			integrate_prof2d(meth, pmodel, cr, image,
 			                 int_diag, idh, idk, idl,
 			                 ir_inn, ir_mid, ir_out,
-			                 results_pipe);
+			                 results_pipe, masks);
 			res = estimate_resolution(crystal_get_cell(cr),
 			                          image->features);
 			break;
@@ -1859,6 +1881,10 @@ void integrate_all_4(struct image *image, IntegrationMethod meth,
 		}
 
 	}
+
+	for ( i=0; i<image->det->n_panels; i++ ) {
+		free(masks[i]);
+	}
 }
 
 
@@ -1947,37 +1973,3 @@ IntegrationMethod integration_method(const char *str, int *err)
 
 }
 
-
-/**
- * flag_overlaps:
- * @image: An image structure
- * @ir_inn: The radius of the peak region to use
- *
- * Flags, in the bad pixel mask for @image, every pixel for which more than one
- * reflection is predicted to have its peak region.
- *
- */
-void flag_overlaps(struct image *image, double ir_inn)
-{
-	int i;
-
-	for ( i=0; i<image->det->n_panels; i++ ) {
-
-		int *mask;
-		int fs, ss;
-		struct panel *p = &image->det->panels[i];
-
-		mask = make_BgMask(image, p, ir_inn);
-
-		for ( ss=0; ss<p->h; ss++ ) {
-		for ( fs=0; fs<p->w; fs++ ) {
-			if ( mask[fs+p->w*ss] > 1 ) {
-				image->bad[i][fs+p->w*ss] = 1;
-			}
-		}
-		}
-
-		free(mask);
-
-	}
-}
diff --git a/libcrystfel/src/integration.h b/libcrystfel/src/integration.h
index 6e80405c..9e1c8ee1 100644
--- a/libcrystfel/src/integration.h
+++ b/libcrystfel/src/integration.h
@@ -112,8 +112,6 @@ extern void integrate_all_4(struct image *image, IntegrationMethod meth,
                             signed int idh, signed int idk, signed int idl,
                             int results_pipe);
 
-extern void flag_overlaps(struct image *image, double ir_inn);
-
 
 #ifdef __cplusplus
 }
diff --git a/src/process_image.c b/src/process_image.c
index 77c5ece7..96752832 100644
--- a/src/process_image.c
+++ b/src/process_image.c
@@ -181,8 +181,6 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs,
 		crystal_set_mosaicity(image.crystals[i], 0.0);  /* radians */
 	}
 
-	flag_overlaps(&image, iargs->ir_inn);
-
 	/* Integrate all the crystals at once - need all the crystals so that
 	 * overlaps can be detected. */
 	integrate_all_4(&image, iargs->int_meth, PMODEL_SPHERE, iargs->push_res,
diff --git a/tests/integration_check.c b/tests/integration_check.c
index b5dbea67..06bdd7d8 100644
--- a/tests/integration_check.c
+++ b/tests/integration_check.c
@@ -140,6 +140,7 @@ int main(int argc, char *argv[])
 		ic.ir_out = ir_out;
 		ic.meth = INTEGRATION_RINGS;
 		ic.int_diag = INTDIAG_NONE;
+		ic.masks = NULL;
 		if ( init_intcontext(&ic) ) {
 			ERROR("Failed to initialise integration.\n");
 			return 1;
diff --git a/tests/prof2d_check.c b/tests/prof2d_check.c
index a1c89d36..3352db49 100644
--- a/tests/prof2d_check.c
+++ b/tests/prof2d_check.c
@@ -177,7 +177,8 @@ int main(int argc, char *argv[])
 
 	reflist_free(list);  /* integrate_prof2d() will predict again */
 	integrate_prof2d(INTEGRATION_PROF2D, PMODEL_SPHERE, cr, &image,
-	                 INTDIAG_NONE, 0, 0, 0, ir_inn, ir_mid, ir_out, 0);
+	                 INTDIAG_NONE, 0, 0, 0, ir_inn, ir_mid, ir_out, 0,
+	                 NULL);
 	list = crystal_get_reflections(cr);
 
 	printf("Weak reflections:\n");