From ee2a9703c723d807c4e590906ab6a45868651a81 Mon Sep 17 00:00:00 2001
From: Thomas White <taw@physics.org>
Date: Mon, 3 Mar 2014 15:44:26 +0100
Subject: Include background and peak height in stream

---
 libcrystfel/src/integration.c | 18 ++++++++++++++++++
 1 file changed, 18 insertions(+)

(limited to 'libcrystfel/src/integration.c')

diff --git a/libcrystfel/src/integration.c b/libcrystfel/src/integration.c
index b86399e9..bcd0a88b 100644
--- a/libcrystfel/src/integration.c
+++ b/libcrystfel/src/integration.c
@@ -216,6 +216,9 @@ struct peak_box
 	double sigma;
 	double J;  /* Profile scaling factor */
 
+	/* Highest (non-ignored) value in peak */
+	double peak;
+
 	/* Offsets to final observed position */
 	double offs_fs;
 	double offs_ss;
@@ -841,6 +844,7 @@ static int check_box(struct intcontext *ic, struct peak_box *bx, int *sat)
 	                   &cdx, &cdy, &cdz);
 	get_indices(bx->refl, &hr, &kr, &lr);
 
+	bx->peak = -INFINITY;
 	for ( p=0; p<ic->w; p++ ) {
 	for ( q=0; q<ic->w; q++ ) {
 
@@ -874,6 +878,12 @@ static int check_box(struct intcontext *ic, struct peak_box *bx, int *sat)
 			if ( sat != NULL ) *sat = 1;
 		}
 
+		if ( (bx->bm[p+ic->w*q] != BM_IG)
+		  && (bx->bm[p+ic->w*q] != BM_BH)
+		  && (boxi(ic, bx, p, q) > bx->peak) ) {
+			bx->peak = boxi(ic, bx, p, q);
+		}
+
 		/* Ignore if this pixel is closer to the next reciprocal lattice
 		 * point */
 		dv = get_q_for_panel(bx->p, fs, ss, NULL, ic->k);
@@ -1342,9 +1352,15 @@ static void integrate_prof2d_once(struct intcontext *ic, struct peak_box *bx)
 	if ( bg_ok(bx) ) {
 
 		double pfs, pss;
+		double bgmean;
+		double sig2_bg;  /* unused */
+
+		mean_var_area(ic, bx, BM_BG, &bgmean, &sig2_bg);
 
 		set_intensity(bx->refl, bx->intensity);
 		set_esd_intensity(bx->refl, bx->sigma);
+		set_peak(bx->refl, bx->peak);
+		set_mean_bg(bx->refl, bgmean);
 		set_redundancy(bx->refl, 1);
 
 		/* Update position */
@@ -1600,6 +1616,8 @@ static void integrate_rings_once(Reflection *refl, struct image *image,
 	set_intensity(refl, intensity);
 	set_esd_intensity(refl, sigma);
 	set_redundancy(refl, 1);
+	set_mean_bg(refl, bgmean);
+	set_peak(refl, bx->peak);
 
 	get_indices(refl, &h, &k, &l);
 	one_over_d = resolution(cell, h, k, l);
-- 
cgit v1.2.3