Make background subtraction optional

4 files changed, 25 insertions, 12 deletions

diff --git a/src/indexamajig.c b/src/indexamajig.c
index b7ec6345..eef9a381 100644
--- a/src/indexamajig.c
+++ b/src/indexamajig.c
@@ -65,6 +65,7 @@ struct static_index_args
 	int config_satcorr;
 	int config_closer;
 	int config_insane;
+	int config_bgsub;
 	float threshold;
 	float min_gradient;
 	struct detector *det;
@@ -179,6 +180,8 @@ static void show_help(const char *s)
 "     --unpolarized        Don't correct for the polarisation of the X-rays.\n"
 "     --no-sat-corr        Don't correct values of saturated peaks using a\n"
 "                           table included in the HDF5 file.\n"
+"     --bg-sub             Subtract local background estimates from\n"
+"                           integrated intensities.\n"
 "     --threshold=<n>      Only accept peaks above <n> ADU.  Default: 800.\n"
 "     --min-gradient=<n>   Minimum gradient for Zaefferer peak search.\n"
 "                           Default: 100,000.\n"
@@ -333,7 +336,8 @@ static void process_image(void *pp, int cookie)
 			                                 0, 0.1);
 
 		integrate_reflections(&image, config_polar,
-			              pargs->static_args.config_closer);
+			              pargs->static_args.config_closer,
+			              pargs->static_args.config_bgsub);
 
 		/* OR */
 
@@ -505,6 +509,7 @@ int main(int argc, char *argv[])
 	int config_checkprefix = 1;
 	int config_closer = 1;
 	int config_insane = 0;
+	int config_bgsub = 0;
 	int config_basename = 0;
 	float threshold = 800.0;
 	float min_gradient = 100000.0;
@@ -565,7 +570,8 @@ int main(int argc, char *argv[])
 		{"record",             1, NULL,                5},
 		{"cpus",               1, NULL,                6},
 		{"cpugroup",           1, NULL,                7},
-		{"cpuoffset",         1, NULL,                 8},
+		{"cpuoffset",          1, NULL,                8},
+		{"bg-sub",             1, &config_bgsub,       0},
 		{0, 0, NULL, 0}
 	};
 
@@ -872,6 +878,7 @@ int main(int argc, char *argv[])
 	qargs.static_args.config_satcorr = config_satcorr;
 	qargs.static_args.config_closer = config_closer;
 	qargs.static_args.config_insane = config_insane;
+	qargs.static_args.config_bgsub = config_bgsub;
 	qargs.static_args.cellr = cellr;
 	qargs.static_args.threshold = threshold;
 	qargs.static_args.min_gradient = min_gradient;
diff --git a/src/pattern_sim.c b/src/pattern_sim.c
index ba760e32..f1fbce7e 100644
--- a/src/pattern_sim.c
+++ b/src/pattern_sim.c
@@ -558,7 +558,7 @@ int main(int argc, char *argv[])
 			image.reflections = find_projected_peaks(&image,
 			                            image.indexed_cell, 0, 0.1);
 
-			integrate_reflections(&image, 0, 0);
+			integrate_reflections(&image, 0, 0, 0);
 
 			/* OR */
 
diff --git a/src/peaks.c b/src/peaks.c
index b838b810..4b87ba08 100644
--- a/src/peaks.c
+++ b/src/peaks.c
@@ -144,7 +144,7 @@ static int cull_peaks(struct image *image)
 int integrate_peak(struct image *image, int cfs, int css,
                    double *pfs, double *pss, double *intensity,
                    double *pbg, double *pmax, double *sigma,
-                   int do_polar, int centroid)
+                   int do_polar, int centroid, int bgsub)
 {
 	signed int fs, ss;
 	double lim, out_lim;
@@ -246,11 +246,14 @@ int integrate_peak(struct image *image, int cfs, int css,
 	if ( centroid && (total != 0) ) {
 		*pfs = (double)fsct / total;
 		*pss = (double)ssct / total;
-		*intensity = total - pixel_counts*noise_mean;
 	} else {
 		*pfs = (double)cfs;
 		*pss = (double)css;
+	}
+	if ( bgsub ) {
 		*intensity = total - pixel_counts*noise_mean;
+	} else {
+		*intensity = total;
 	}
 
 	if ( in_bad_region(image->det, *pfs, *pss) ) return 1;
@@ -380,7 +383,7 @@ static void search_peaks_in_panel(struct image *image, float threshold,
 		 * intensity of this peak is only an estimate at this stage. */
 		r = integrate_peak(image, mask_fs, mask_ss,
 		                   &f_fs, &f_ss, &intensity,
-		                   NULL, NULL, NULL, 0, 1);
+		                   NULL, NULL, NULL, 0, 1, 0);
 		if ( r ) {
 			/* Bad region - don't detect peak */
 			nrej_bad++;
@@ -592,7 +595,8 @@ int peak_sanity_check(struct image *image, UnitCell *cell,
 
 
 /* Integrate the list of predicted reflections in "image" */
-void integrate_reflections(struct image *image, int polar, int use_closer)
+void integrate_reflections(struct image *image, int polar, int use_closer,
+                           int bgsub)
 {
 	Reflection *refl;
 	RefListIterator *iter;
@@ -635,7 +639,8 @@ void integrate_reflections(struct image *image, int polar, int use_closer)
 		}
 
 		r = integrate_peak(image, pfs, pss, &fs, &ss,
-		                   &intensity, &bg, &max, &sigma, polar, 0);
+		                   &intensity, &bg, &max, &sigma, polar, 0,
+		                   bgsub);
 
 		/* Record intensity and set redundancy to 1 on success */
 		if ( r == 0 ) {
diff --git a/src/peaks.h b/src/peaks.h
index 2faa2b5e..ec473fd1 100644
--- a/src/peaks.h
+++ b/src/peaks.h
@@ -25,7 +25,7 @@ extern void search_peaks(struct image *image, float threshold,
                          float min_gradient);
 
 extern void integrate_reflections(struct image *image,
-                                  int polar, int use_closer);
+                                  int polar, int use_closer, int bgsub);
 
 extern int peak_sanity_check(struct image *image, UnitCell *cell,
                              int circular_domain, double domain_r);
@@ -33,9 +33,10 @@ extern int peak_sanity_check(struct image *image, UnitCell *cell,
 extern RefList *find_projected_peaks(struct image *image, UnitCell *cell,
                                      int circular_domain, double domain_r);
 
-extern int integrate_peak(struct image *image, int xp, int yp,
-                          double *xc, double *yc, double *intensity,
+/* Exported so it can be poked by integration_check */
+extern int integrate_peak(struct image *image, int cfs, int css,
+                          double *pfs, double *pss, double *intensity,
                           double *pbg, double *pmax, double *sigma,
-                          int do_polar, int centroid);
+                          int do_polar, int centroid, int bgsub);
 
 #endif	/* PEAKS_H */