Perform intensity extraction on the image as it was before noise filtering

5 files changed, 51 insertions, 29 deletions

diff --git a/src/filters.c b/src/filters.c
index 172be222..17792b85 100644
--- a/src/filters.c
+++ b/src/filters.c
@@ -61,7 +61,25 @@ static void clean_panel(struct image *image, int sx, int sy)
 }
 
 
-static void noise_filter(struct image *image)
+/* Pre-processing to make life easier */
+void filter_cm(struct image *image)
+{
+	int px, py;
+
+	if ( (image->width != 1024) || (image->height != 1024) ) return;
+
+	for ( px=0; px<2; px++ ) {
+	for ( py=0; py<8; py++ ) {
+
+		clean_panel(image, 512*px, 128*py);
+
+	}
+	}
+
+}
+
+
+void filter_noise(struct image *image, float *old)
 {
 	int x, y;
 
@@ -71,6 +89,10 @@ static void noise_filter(struct image *image)
 		int dx, dy;
 		int val = image->data[x+image->width*y];
 
+		if ( old != NULL ) old[x+image->width*y] = val;
+
+		/* FIXME: This isn't really the right thing to do
+		 * at the edges. */
 		if ( (x==0) || (x==image->width-1)
 		  || (y==0) || (y==image->height-1) ) {
 			if ( val < 0 ) val = 0;
@@ -94,27 +116,3 @@ static void noise_filter(struct image *image)
 	}
 	}
 }
-
-
-/* Pre-processing to make life easier */
-void filter_cm(struct image *image)
-{
-	int px, py;
-
-	if ( (image->width != 1024) || (image->height != 1024) ) return;
-
-	for ( px=0; px<2; px++ ) {
-	for ( py=0; py<8; py++ ) {
-
-		clean_panel(image, 512*px, 128*py);
-
-	}
-	}
-
-}
-
-
-void filter_noise(struct image *image)
-{
-	noise_filter(image);
-}
diff --git a/src/filters.h b/src/filters.h
index 97bcf5a2..6b35d7e8 100644
--- a/src/filters.h
+++ b/src/filters.h
@@ -19,7 +19,7 @@
 
 
 extern void filter_cm(struct image *image);
-extern void filter_noise(struct image *image);
+extern void filter_noise(struct image *image, float *old);
 
 
 #endif	/* FILTERS_H */
diff --git a/src/image.h b/src/image.h
index dfd02d05..5f2cf3d0 100644
--- a/src/image.h
+++ b/src/image.h
@@ -70,7 +70,7 @@ struct rvec
 /* Structure describing an image */
 struct image {
 
-	float			*data;     /* Integer counts after bloom */
+	float			*data;
 	double			*twotheta;
 	struct molecule		*molecule;
 	UnitCell		*indexed_cell;
diff --git a/src/indexamajig.c b/src/indexamajig.c
index 9958369d..9bdfa5e0 100644
--- a/src/indexamajig.c
+++ b/src/indexamajig.c
@@ -280,6 +280,8 @@ int main(int argc, char *argv[])
 		char line[1024];
 		struct hdfile *hdfile;
 		struct image *simage;
+		float *data_for_measurement;
+		size_t data_size;
 
 		rval = fgets(line, 1023, fh);
 		if ( rval == NULL ) continue;
@@ -309,6 +311,27 @@ int main(int argc, char *argv[])
 			filter_cm(&image);
 		}
 
+		/* Take snapshot of image after CM subtraction but before
+		 * the aggressive noise filter. */
+		data_size = image.width*image.height*sizeof(float);
+		data_for_measurement = malloc(data_size);
+
+		if ( config_noisefilter ) {
+			filter_noise(&image, data_for_measurement);
+		} else {
+
+			int x, y;
+
+			for ( x=0; x<image.width; x++ ) {
+			for ( y=0; y<image.height; y++ ) {
+				float val;
+				val = image.data[x+image.width*y];
+				data_for_measurement[x+image.width*y] = val;
+			}
+			}
+
+		}
+
 		/* Perform 'fine' peak search */
 		search_peaks(&image);
 
@@ -338,7 +361,7 @@ int main(int argc, char *argv[])
 		/* Measure intensities if requested */
 		if ( config_nearbragg ) {
 			/* Use original data (temporarily) */
-			simage->data = image.data;
+			simage->data = data_for_measurement;
 			output_intensities(simage, image.indexed_cell);
 			simage->data = NULL;
 		}
@@ -365,6 +388,7 @@ done:
 		free(image.data);
 		free(image.det.panels);
 		image_feature_list_free(image.features);
+		free(data_for_measurement);
 		hdfile_close(hdfile);
 		H5close();
 
diff --git a/src/render.c b/src/render.c
index a4be5a14..206ef9b3 100644
--- a/src/render.c
+++ b/src/render.c
@@ -78,7 +78,7 @@ float *render_get_image_binned(DisplayWindow *dw, int binning, float *max)
 		hdf5_read(dw->hdfile, image);
 		dw->image_dirty = 0;
 		if ( dw->cmfilter ) filter_cm(image);
-		if ( dw->noisefilter ) filter_noise(image);
+		if ( dw->noisefilter ) filter_noise(image, NULL);
 
 		/* Deal with the old image, if existing */
 		if ( dw->image != NULL ) {