1 files changed, 114 insertions, 192 deletions

diff --git a/src/process_image.c b/src/process_image.c
index b785705a..d3355f38 100644
--- a/src/process_image.c
+++ b/src/process_image.c
@@ -3,11 +3,11 @@
  *
  * The processing pipeline for one image
  *
- * Copyright © 2012-2020 Deutsches Elektronen-Synchrotron DESY,
+ * Copyright © 2012-2021 Deutsches Elektronen-Synchrotron DESY,
  *                       a research centre of the Helmholtz Association.
  *
  * Authors:
- *   2010-2019 Thomas White <taw@physics.org>
+ *   2010-2020 Thomas White <taw@physics.org>
  *   2014-2017 Valerio Mariani <valerio.mariani@desy.de>
  *   2017      Stijn de Graaf
  *
@@ -33,32 +33,31 @@
 #endif
 
 #include <stdlib.h>
-#include <hdf5.h>
 #include <gsl/gsl_errno.h>
 #include <gsl/gsl_statistics_double.h>
 #include <gsl/gsl_sort.h>
 #include <unistd.h>
 #include <sys/stat.h>
 
-#include "utils.h"
-#include "hdf5-file.h"
-#include "index.h"
-#include "peaks.h"
-#include "detector.h"
-#include "filters.h"
-#include "thread-pool.h"
-#include "geometry.h"
-#include "stream.h"
-#include "reflist-utils.h"
+#include <utils.h>
+#include <index.h>
+#include <peaks.h>
+#include <filters.h>
+#include <thread-pool.h>
+#include <geometry.h>
+#include <stream.h>
+#include <reflist-utils.h>
+#include <integration.h>
+#include <detgeom.h>
+#include <image-msgpack.h>
+
 #include "process_image.h"
-#include "integration.h"
 #include "predict-refine.h"
 #include "im-sandbox.h"
-#include "time-accounts.h"
 #include "im-zmq.h"
+#include "time-accounts.h"
 
-
-static float **backup_image_data(float **dp, struct detector *det)
+static float **backup_image_data(float **dp, struct detgeom *det)
 {
 	float **bu;
 	int i;
@@ -86,7 +85,7 @@ static float **backup_image_data(float **dp, struct detector *det)
 }
 
 
-static void restore_image_data(float **dp, struct detector *det, float **bu)
+static void restore_image_data(float **dp, struct detgeom *det, float **bu)
 {
 	int i;
 
@@ -100,16 +99,22 @@ static void restore_image_data(float **dp, struct detector *det, float **bu)
 }
 
 
-static int file_wait_open_read(struct sb_shm *sb_shared, struct image *image,
-                               TimeAccounts *taccs, char *last_task,
-                               signed int wait_for_file, int cookie,
-                               struct imagefile **pimfile)
+static struct image *file_wait_open_read(const char *filename,
+                                         const char *event,
+                                         DataTemplate *dtempl,
+                                         struct sb_shm *sb_shared,
+                                         TimeAccounts *taccs,
+                                         char *last_task,
+                                         signed int wait_for_file,
+                                         int cookie,
+                                         int no_image_data,
+                                         int no_mask_data)
 {
 	signed int file_wait_time = wait_for_file;
 	int wait_message_done = 0;
 	int read_retry_done = 0;
 	int r;
-	struct imagefile *imfile;
+	struct image *image;
 
 	time_accounts_set(taccs, TACC_WAITFILE);
 	set_last_task(last_task, "wait for file");
@@ -119,14 +124,14 @@ static int file_wait_open_read(struct sb_shm *sb_shared, struct image *image,
 		struct stat statbuf;
 
 		sb_shared->pings[cookie]++;
-		r = stat(image->filename, &statbuf);
+		r = stat(filename, &statbuf);
 		if ( r ) {
 
 			if ( (wait_for_file != 0) && (file_wait_time != 0) ) {
 
 				if ( !wait_message_done ) {
 					STATUS("Waiting for '%s'\n",
-					       image->filename);
+					       filename);
 					wait_message_done = 1;
 				}
 
@@ -138,8 +143,8 @@ static int file_wait_open_read(struct sb_shm *sb_shared, struct image *image,
 
 			}
 
-			ERROR("File %s not found\n", image->filename);
-			return 1;
+			ERROR("File %s not found\n", filename);
+			return NULL;
 		}
 
 	} while ( r );
@@ -149,44 +154,29 @@ static int file_wait_open_read(struct sb_shm *sb_shared, struct image *image,
 	sb_shared->pings[cookie]++;
 
 	do {
-		imfile = imagefile_open(image->filename);
-		if ( imfile == NULL ) {
-			if ( wait_for_file && !read_retry_done ) {
-				read_retry_done = 1;
-				r = 1;
-				STATUS("File '%s' exists but could not be opened."
-				       "  Trying again after 10 seconds.\n",
-				       image->filename);
-				sleep(10);
-				continue;
-			}
-			ERROR("Couldn't open file: %s\n", image->filename);
-			return 1;
-		}
 
 		time_accounts_set(taccs, TACC_HDF5READ);
 		set_last_task(last_task, "read file");
 		sb_shared->pings[cookie]++;
 
-		r = imagefile_read(imfile, image, image->event);
-		if ( r ) {
+		image = image_read(dtempl, filename, event,
+		                   no_image_data, no_mask_data);
+		if ( image == NULL ) {
 			if ( wait_for_file && !read_retry_done ) {
 				read_retry_done = 1;
-				imagefile_close(imfile);
 				STATUS("File '%s' exists but could not be read."
 				       "  Trying again after 10 seconds.\n",
-				       image->filename);
+				       filename);
 				sleep(10);
 				continue;
 			}
-			ERROR("Couldn't open file: %s\n", image->filename);
-			return 1;
+			ERROR("Couldn't open file: %s\n", filename);
+			return NULL;
 		}
 
 	} while ( r );
 
-	*pimfile = imfile;
-	return 0;
+	return image;
 }
 
 
@@ -195,8 +185,7 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs,
                    int serial, struct sb_shm *sb_shared, TimeAccounts *taccs,
                    char *last_task)
 {
-	struct imagefile *imfile = NULL;
-	struct image image;
+	struct image *image;
 	int i;
 	int r;
 	int ret;
@@ -204,27 +193,21 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs,
 	float **prefilter;
 	int any_crystals;
 
-	image.features = NULL;
-	image.copyme = iargs->copyme;
-	image.id = cookie;
-	image.beam = iargs->beam;
-	image.det = copy_geom(iargs->det);
-	image.crystals = NULL;
-	image.n_crystals = 0;
-	image.serial = serial;
-	image.indexed_by = INDEXING_NONE;
-
-	image.filename = pargs->filename_p_e->filename;
-	image.event = pargs->filename_p_e->ev;
 	if ( pargs->msgpack_obj != NULL ) {
 		set_last_task(last_task, "unpacking messagepack object");
-		if ( unpack_msgpack_data(pargs->msgpack_obj, &image,
-		                         iargs->no_image_data) ) return;
+		image = image_msgpack_read(iargs->dtempl,
+		                           pargs->msgpack_obj,
+		                           iargs->no_image_data);
+		if ( image == NULL ) return;
 	} else {
-		if ( file_wait_open_read(sb_shared, &image, taccs, last_task,
-		                         iargs->wait_for_file, cookie,
-		                         &imfile) )
-		{
+		image = file_wait_open_read(pargs->filename, pargs->event,
+		                            iargs->dtempl,
+		                            sb_shared, taccs, last_task,
+		                            iargs->wait_for_file,
+		                            cookie,
+		                            iargs->no_image_data,
+		                            iargs->no_mask_data);
+		if ( image == NULL ) {
 			if ( iargs->wait_for_file != 0 ) {
 				pthread_mutex_lock(&sb_shared->totals_lock);
 				sb_shared->should_shutdown = 1;
@@ -234,62 +217,43 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs,
 		}
 	}
 
-	image.bw = iargs->beam->bandwidth;
+	image->serial = serial;
 
 	/* Take snapshot of image before applying horrible noise filters */
 	time_accounts_set(taccs, TACC_FILTER);
 	set_last_task(last_task, "image filter");
 	sb_shared->pings[cookie]++;
-	prefilter = backup_image_data(image.dp, image.det);
+	prefilter = backup_image_data(image->dp, image->detgeom);
 
 	if ( iargs->median_filter > 0 ) {
-		filter_median(&image, iargs->median_filter);
+		filter_median(image, iargs->median_filter);
 	}
 
 	if ( iargs->noisefilter ) {
-		filter_noise(&image);
+		filter_noise(image);
 	}
 
 	time_accounts_set(taccs, TACC_RESRANGE);
 	set_last_task(last_task, "resolution range");
 	sb_shared->pings[cookie]++;
-	mark_resolution_range_as_bad(&image, iargs->highres, +INFINITY);
+	mark_resolution_range_as_bad(image, iargs->highres, +INFINITY);
 
 	time_accounts_set(taccs, TACC_PEAKSEARCH);
 	sb_shared->pings[cookie]++;
 	switch ( iargs->peaks ) {
 
-		struct hdfile *hdfile;
-
 		case PEAK_HDF5:
-		set_last_task(last_task, "peaksearch:hdf5");
-		hdfile = imagefile_get_hdfile(imfile);
-		if ( (hdfile == NULL)
-		  || (get_peaks_2(&image, hdfile, iargs->hdf5_peak_path,
-		                  iargs->half_pixel_shift)) )
-		{
-			ERROR("Failed to get peaks from HDF5 file.\n");
-		}
-		if ( !iargs->no_revalidate ) {
-			validate_peaks(&image, iargs->min_snr,
-				       iargs->pk_inn, iargs->pk_mid,
-		                       iargs->pk_out, iargs->use_saturated,
-				       iargs->check_hdf5_snr);
-		}
-		break;
-
 		case PEAK_CXI:
-		set_last_task(last_task, "peaksearch:cxi");
-		hdfile = imagefile_get_hdfile(imfile);
-		if ( (hdfile == NULL)
-		 ||  (get_peaks_cxi_2(&image, hdfile, iargs->hdf5_peak_path,
-		                      pargs->filename_p_e,
-		                      iargs->half_pixel_shift)) )
-		{
-			ERROR("Failed to get peaks from CXI file.\n");
+		set_last_task(last_task, "peaksearch:hdf5orcxi");
+		image->features = image_read_peaks(iargs->dtempl,
+		                                   pargs->filename,
+		                                   pargs->event,
+		                                   iargs->half_pixel_shift);
+		if ( image->features == NULL ) {
+			ERROR("Failed to get peaks from HDF5 file.\n");
 		}
 		if ( !iargs->no_revalidate ) {
-			validate_peaks(&image, iargs->min_snr,
+			validate_peaks(image, iargs->min_snr,
 				       iargs->pk_inn, iargs->pk_mid,
 		                       iargs->pk_out, iargs->use_saturated,
 				       iargs->check_hdf5_snr);
@@ -298,7 +262,7 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs,
 
 		case PEAK_ZAEF:
 		set_last_task(last_task, "peaksearch:zaef");
-		search_peaks(&image, iargs->threshold,
+		search_peaks(image, iargs->threshold,
 		             iargs->min_sq_gradient, iargs->min_snr,
 		             iargs->pk_inn, iargs->pk_mid, iargs->pk_out,
 		             iargs->use_saturated);
@@ -306,30 +270,24 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs,
 
 		case PEAK_PEAKFINDER8:
 		set_last_task(last_task, "peaksearch:pf8");
-		if ( search_peaks_peakfinder8(&image, 2048,
-		                               iargs->threshold,
-		                               iargs->min_snr,
-		                               iargs->min_pix_count,
-		                               iargs->max_pix_count,
-		                               iargs->local_bg_radius,
-		                               iargs->min_res,
-		                               iargs->max_res,
-		                               iargs->use_saturated) ) {
-			if ( image.event != NULL ) {
-				ERROR("Failed to find peaks in image %s"
-				      "(event %s).\n", image.filename,
-				      get_event_string(image.event));
-			} else {
-				ERROR("Failed to find peaks in image %s.",
-				      image.filename);
-			}
-
+		if ( search_peaks_peakfinder8(image, 2048,
+		                              iargs->threshold,
+		                              iargs->min_snr,
+		                              iargs->min_pix_count,
+		                              iargs->max_pix_count,
+		                              iargs->local_bg_radius,
+		                              iargs->min_res,
+		                              iargs->max_res,
+		                              iargs->use_saturated) ) {
+			ERROR("Failed to find peaks in image %s"
+			      "(event %s).\n",
+			      image->filename, image->ev);
 		}
 		break;
 
 		case PEAK_PEAKFINDER9:
 		set_last_task(last_task, "peaksearch:pf9");
-		if ( search_peaks_peakfinder9(&image,
+		if ( search_peaks_peakfinder9(image,
 		                              iargs->min_snr_biggest_pix,
 		                              iargs->min_snr_peak_pix,
 		                              iargs->min_snr,
@@ -337,31 +295,29 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs,
 		                              iargs->min_peak_over_neighbour,
 		                              iargs->local_bg_radius) )
 		{
-			if ( image.event != NULL ) {
-				ERROR("Failed to find peaks in image %s"
-				      "(event %s).\n", image.filename,
-				      get_event_string(image.event));
-			} else {
-				ERROR("Failed to find peaks in image %s.",
-				      image.filename);
-			}
+			ERROR("Failed to find peaks in image %s"
+			      "(event %s).\n",
+			      image->filename, image->ev);
 		}
 		break;
 
 		case PEAK_MSGPACK:
-		get_peaks_msgpack(pargs->msgpack_obj, &image,
-		                  iargs->half_pixel_shift);
+		image->features = image_msgpack_read_peaks(iargs->dtempl,
+		                                           pargs->msgpack_obj,
+		                                           iargs->half_pixel_shift);
 		break;
 
 		case PEAK_NONE:
+		case PEAK_ERROR:
 		break;
 
 	}
 
-	image.peak_resolution = estimate_peak_resolution(image.features,
-	                                                 image.lambda);
+	image->peak_resolution = estimate_peak_resolution(image->features,
+	                                                  image->lambda,
+	                                                  image->detgeom);
 
-	restore_image_data(image.dp, image.det, prefilter);
+	restore_image_data(image->dp, image->detgeom, prefilter);
 
 	rn = getcwd(NULL, 0);
 
@@ -369,33 +325,24 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs,
 	if ( r ) {
 		ERROR("Failed to chdir to temporary folder: %s\n",
 		      strerror(errno));
-		imagefile_close(imfile);
 		return;
 	}
 
 	/* Set beam parameters */
 	if ( iargs->fix_divergence >= 0.0 ) {
-		image.div = iargs->fix_divergence;
+		image->div = iargs->fix_divergence;
 	} else {
-		image.div = 0.0;
-	}
-	/* Set beam spectrum */
-	if ( iargs->spectrum != NULL ) {
-		image.spectrum = iargs->spectrum;
-	} else {
-		image.spectrum = spectrum_generate_gaussian(image.lambda,
-		                                            image.bw);
+		image->div = 0.0;
 	}
 
-	if ( image_feature_count(image.features) < iargs->min_peaks ) {
+	if ( image_feature_count(image->features) < iargs->min_peaks ) {
 		r = chdir(rn);
 		if ( r ) {
 			ERROR("Failed to chdir: %s\n", strerror(errno));
-			imagefile_close(imfile);
 			return;
 		}
 		free(rn);
-		image.hit = 0;
+		image->hit = 0;
 
 		if ( iargs->stream_nonhits ) {
 			goto streamwrite;
@@ -403,18 +350,17 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs,
 			goto out;
 		}
 	}
-	image.hit = 1;
+	image->hit = 1;
 
 	/* Index the pattern */
 	time_accounts_set(taccs, TACC_INDEXING);
 	set_last_task(last_task, "indexing");
-	index_pattern_3(&image, iargs->ipriv, &sb_shared->pings[cookie],
+	index_pattern_3(image, iargs->ipriv, &sb_shared->pings[cookie],
 	                last_task);
 
 	r = chdir(rn);
 	if ( r ) {
 		ERROR("Failed to chdir: %s\n", strerror(errno));
-		imagefile_close(imfile);
 		return;
 	}
 	free(rn);
@@ -423,22 +369,22 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs,
 	time_accounts_set(taccs, TACC_PREDPARAMS);
 	set_last_task(last_task, "prediction params");
 	if ( iargs->fix_profile_r >= 0.0 ) {
-		for ( i=0; i<image.n_crystals; i++ ) {
-			crystal_set_profile_radius(image.crystals[i],
+		for ( i=0; i<image->n_crystals; i++ ) {
+			crystal_set_profile_radius(image->crystals[i],
 			                           iargs->fix_profile_r);
-			crystal_set_mosaicity(image.crystals[i], 0.0);
+			crystal_set_mosaicity(image->crystals[i], 0.0);
 		}
 	} else {
-		for ( i=0; i<image.n_crystals; i++ ) {
-			crystal_set_profile_radius(image.crystals[i], 0.02e9);
-			crystal_set_mosaicity(image.crystals[i], 0.0);
+		for ( i=0; i<image->n_crystals; i++ ) {
+			crystal_set_profile_radius(image->crystals[i], 0.02e9);
+			crystal_set_mosaicity(image->crystals[i], 0.0);
 		}
 	}
 
 
 	if ( iargs->fix_profile_r < 0.0 ) {
-		for ( i=0; i<image.n_crystals; i++ ) {
-			if ( refine_radius(image.crystals[i], &image) ) {
+		for ( i=0; i<image->n_crystals; i++ ) {
+			if ( refine_radius(image->crystals[i], image) ) {
 				ERROR("WARNING: Radius determination failed\n");
 			}
 		}
@@ -448,7 +394,7 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs,
 	time_accounts_set(taccs, TACC_INTEGRATION);
 	set_last_task(last_task, "integration");
 	sb_shared->pings[cookie]++;
-	integrate_all_5(&image, iargs->int_meth, PMODEL_XSPHERE,
+	integrate_all_5(image, iargs->int_meth, PMODEL_XSPHERE,
 	                iargs->push_res,
 	                iargs->ir_inn, iargs->ir_mid, iargs->ir_out,
 	                iargs->int_diag, iargs->int_diag_h,
@@ -459,21 +405,18 @@ streamwrite:
 	time_accounts_set(taccs, TACC_WRITESTREAM);
 	set_last_task(last_task, "stream write");
 	sb_shared->pings[cookie]++;
-	ret = write_chunk(st, &image, imfile,
-	                  iargs->stream_peaks, iargs->stream_refls,
-	                  pargs->filename_p_e->ev);
+	ret = stream_write_chunk(st, image, iargs->stream_flags);
 	if ( ret != 0 ) {
 		ERROR("Error writing stream file.\n");
 	}
 
 	int n = 0;
-	for ( i=0; i<image.n_crystals; i++ ) {
-		n += crystal_get_num_implausible_reflections(image.crystals[i]);
+	for ( i=0; i<image->n_crystals; i++ ) {
+		n += crystal_get_num_implausible_reflections(image->crystals[i]);
 	}
 	if ( n > 0 ) {
 		STATUS("WARNING: %i implausibly negative reflection%s in %s "
-		       "%s\n", n, n>1?"s":"", image.filename,
-		       get_event_string(image.event));
+		       "%s\n", n, n>1?"s":"", image->filename, image->ev);
 	}
 
 out:
@@ -483,40 +426,19 @@ out:
 	sb_shared->pings[cookie]++;
 	pthread_mutex_lock(&sb_shared->totals_lock);
 	any_crystals = 0;
-	for ( i=0; i<image.n_crystals; i++ ) {
-		if ( crystal_get_user_flag(image.crystals[i]) == 0 ) {
+	for ( i=0; i<image->n_crystals; i++ ) {
+		if ( crystal_get_user_flag(image->crystals[i]) == 0 ) {
 			sb_shared->n_crystals++;
 			any_crystals = 1;
 		}
 	}
 	sb_shared->n_processed++;
-	sb_shared->n_hits += image.hit;
+	sb_shared->n_hits += image->hit;
 	sb_shared->n_hadcrystals += any_crystals;
 	pthread_mutex_unlock(&sb_shared->totals_lock);
 
-	for ( i=0; i<image.n_crystals; i++ ) {
-		cell_free(crystal_get_cell(image.crystals[i]));
-		reflist_free(crystal_get_reflections(image.crystals[i]));
-		crystal_free(image.crystals[i]);
-	}
-	free(image.crystals);
-
-	/* Free spectrum only if we generated it for this image */
-	if ( iargs->spectrum == NULL ) {
-		spectrum_free(image.spectrum);
-	}
-
-	for ( i=0; i<image.det->n_panels; i++ ) {
-		free(image.dp[i]);
-		free(image.bad[i]);
-		free(image.sat[i]);
-	}
-	free(image.dp);
-	free(image.bad);
-	free(image.sat);
+	/* Free image (including detgeom) */
+	image_free(image);
 
-	image_feature_list_free(image.features);
-	free_detector_geometry(image.det);
-	if ( imfile != NULL ) imagefile_close(imfile);
 	set_last_task(last_task, "sandbox");
 }