Take nominal photon energy from file

12 files changed, 57 insertions, 18 deletions

diff --git a/src/calibrate_detector.c b/src/calibrate_detector.c
index 5ca468c6..b20b1619 100644
--- a/src/calibrate_detector.c
+++ b/src/calibrate_detector.c
@@ -189,7 +189,8 @@ static void add_image(void *args, int cookie)
 		return;
 	}
 
-	hdf5_read(hdfile, &image, 1);
+	/* FIXME: Nominal photon energy */
+	hdf5_read(hdfile, &image, 1, 2000.0);
 
 	if ( pargs->config_cmfilter ) {
 		filter_cm(&image);
diff --git a/src/cubeit.c b/src/cubeit.c
index ad78704e..b613a3fa 100644
--- a/src/cubeit.c
+++ b/src/cubeit.c
@@ -237,7 +237,8 @@ static void sum_image(void *pg, int cookie)
 		return;
 	}
 
-	hdf5_read(hdfile, &image, 1);
+	/* FIXME: Nominal photon energy */
+	hdf5_read(hdfile, &image, 1, 2000.0);
 
 	cell_get_cartesian(apargs->cell, &ax, &ay, &az, &bx, &by,
 	                                &bz, &cx, &cy, &cz);
diff --git a/src/hdf5-file.c b/src/hdf5-file.c
index bcd51f99..a24a663a 100644
--- a/src/hdf5-file.c
+++ b/src/hdf5-file.c
@@ -390,7 +390,8 @@ static void debodge_saturation(struct hdfile *f, struct image *image)
 }
 
 
-int hdf5_read(struct hdfile *f, struct image *image, int satcorr)
+int hdf5_read(struct hdfile *f, struct image *image, int satcorr,
+              double nominal_photon_energy)
 {
 	herr_t r;
 	float *buf;
@@ -432,8 +433,8 @@ int hdf5_read(struct hdfile *f, struct image *image, int satcorr)
 	/* Read wavelength from file */
 	image->lambda = get_wavelength(f);
 	if ( image->lambda < 0.0 ) {
-		ERROR("Couldn't read wavelength - using 2 keV.\n");
-		image->lambda = ph_en_to_lambda(eV_to_J(2000.0));
+		ERROR("Couldn't read wavelength - using nominal value.\n");
+		image->lambda = ph_en_to_lambda(eV_to_J(nominal_photon_energy));
 	}
 
 	image->f0 = get_f0(f);
diff --git a/src/hdf5-file.h b/src/hdf5-file.h
index 145c9cac..4f849638 100644
--- a/src/hdf5-file.h
+++ b/src/hdf5-file.h
@@ -27,7 +27,8 @@ struct hdfile;
 extern int hdf5_write(const char *filename, const void *data,
                       int width, int height, int type);
 
-extern int hdf5_read(struct hdfile *f, struct image *image, int satcorr);
+extern int hdf5_read(struct hdfile *f, struct image *image, int satcorr,
+                     double nominal_photon_energy);
 
 extern struct hdfile *hdfile_open(const char *filename);
 extern int hdfile_set_image(struct hdfile *f, const char *path);
diff --git a/src/index.c b/src/index.c
index 97cf3887..0c656042 100644
--- a/src/index.c
+++ b/src/index.c
@@ -42,7 +42,8 @@ static IndexingPrivate *indexing_private(IndexingMethod indm)
 
 
 IndexingPrivate *prepare_indexing(IndexingMethod indm, UnitCell *cell,
-                                  const char *filename, struct detector *det)
+                                  const char *filename, struct detector *det,
+                                  double nominal_photon_energy)
 {
 	switch ( indm ) {
 	case INDEXING_NONE :
@@ -50,7 +51,8 @@ IndexingPrivate *prepare_indexing(IndexingMethod indm, UnitCell *cell,
 	case INDEXING_DIRAX :
 		return indexing_private(indm);
 	case INDEXING_TEMPLATE :
-		return generate_templates(cell, filename, det);
+		return generate_templates(cell, filename, det,
+		                          nominal_photon_energy);
 	}
 	return 0;
 }
diff --git a/src/index.h b/src/index.h
index 259d3eb6..0b740995 100644
--- a/src/index.h
+++ b/src/index.h
@@ -33,7 +33,8 @@ typedef struct _indexingprivate IndexingPrivate;
 
 extern IndexingPrivate *prepare_indexing(IndexingMethod indm, UnitCell *cell,
                                          const char *filename,
-                                         struct detector *det);
+                                         struct detector *det,
+                                         double nominal_photon_energy);
 
 extern void map_all_peaks(struct image *image);
 
diff --git a/src/indexamajig.c b/src/indexamajig.c
index 54d5e630..348e3e8b 100644
--- a/src/indexamajig.c
+++ b/src/indexamajig.c
@@ -36,6 +36,7 @@
 #include "filters.h"
 #include "reflections.h"
 #include "thread-pool.h"
+#include "beam-parameters.h"
 
 
 enum {
@@ -71,6 +72,7 @@ struct static_index_args
 	const double *intensities;
 	struct gpu_context *gctx;
 	int peaks;
+	double nominal_photon_energy;
 
 	/* Output stream */
 	pthread_mutex_t *output_mutex;  /* Protects the output stream */
@@ -120,6 +122,9 @@ static void show_help(const char *s)
 "                           dirax    : invoke DirAx\n"
 "                           template : index by template matching\n"
 " -g. --geometry=<file>    Get detector geometry from file.\n"
+" -b, --beam=<file>        Get beam parameters from file (provides nominal\n"
+"                           wavelength value if no per-shot value is found in\n"
+"                           the HDF5 files.\n"
 " -p, --pdb=<file>         PDB file from which to get the unit cell to match.\n"
 "                           Default: 'molecule.pdb'.\n"
 " -x, --prefix=<p>         Prefix filenames from input file with <p>.\n"
@@ -340,7 +345,8 @@ static void process_image(void *pp, int cookie)
 		return;
 	}
 
-	hdf5_read(hdfile, &image, pargs->static_args.config_satcorr);
+	hdf5_read(hdfile, &image, pargs->static_args.config_satcorr,
+	          pargs->static_args.nominal_photon_energy);
 
 	if ( config_cmfilter ) {
 		filter_cm(&image);
@@ -531,6 +537,8 @@ int main(int argc, char *argv[])
 	char prepare_filename[1024];
 	IndexingPrivate *ipriv;
 	struct queue_args qargs;
+	struct beam_params *beam = NULL;
+	double nominal_photon_energy;
 
 	/* Long options */
 	const struct option longopts[] = {
@@ -545,6 +553,7 @@ int main(int argc, char *argv[])
 		{"write-drx",          0, &config_writedrx,    1},
 		{"indexing",           1, NULL,               'z'},
 		{"geometry",           1, NULL,               'g'},
+		{"beam",               1, NULL,               'b'},
 		{"simulate",           0, &config_simulate,    1},
 		{"filter-cm",          0, &config_cmfilter,    1},
 		{"filter-noise",       0, &config_noisefilter, 1},
@@ -566,7 +575,7 @@ int main(int argc, char *argv[])
 	};
 
 	/* Short options */
-	while ((c = getopt_long(argc, argv, "hi:wp:j:x:g:t:o:",
+	while ((c = getopt_long(argc, argv, "hi:wp:j:x:g:t:o:b:",
 	                        longopts, NULL)) != -1) {
 
 		switch (c) {
@@ -610,6 +619,15 @@ int main(int argc, char *argv[])
 			threshold = strtof(optarg, NULL);
 			break;
 
+		case 'b' :
+			beam = get_beam_parameters(optarg);
+			if ( beam == NULL ) {
+				ERROR("Failed to load beam parameters"
+				      " from '%s'\n", optarg);
+				return 1;
+			}
+			break;
+
 		case 2 :
 			speaks = strdup(optarg);
 			break;
@@ -742,6 +760,14 @@ int main(int argc, char *argv[])
 	fprintf(ofh, "\n");
 	fflush(ofh);
 
+	if ( beam != NULL ) {
+		nominal_photon_energy = beam->photon_energy;
+	} else {
+		STATUS("No beam parameters file was given, so I'm taking the"
+		       " nominal photon energy to be 2 keV.\n");
+		nominal_photon_energy = 2000.0;
+	}
+
 	/* Get first filename and use it to set up the indexing */
 	if ( fh != stdin ) {
 		rval = fgets(prepare_line, 1023, fh);
@@ -758,7 +784,8 @@ int main(int argc, char *argv[])
 		STATUS("Stuff might break.\n");
 		prepare_filename[0] = '\0';
 	}
-	ipriv = prepare_indexing(indm, cell, prepare_filename, det);
+	ipriv = prepare_indexing(indm, cell, prepare_filename, det,
+	                         nominal_photon_energy);
 	if ( ipriv == NULL ) {
 		ERROR("Failed to prepare indexing.\n");
 		return 1;
@@ -792,6 +819,7 @@ int main(int argc, char *argv[])
 	qargs.static_args.peaks = peaks;
 	qargs.static_args.output_mutex = &output_mutex;
 	qargs.static_args.ofh = ofh;
+	qargs.static_args.nominal_photon_energy = nominal_photon_energy;
 
 	qargs.fh = fh;
 	qargs.prefix = prefix;
diff --git a/src/powder_plot.c b/src/powder_plot.c
index d9049cc3..f8038274 100644
--- a/src/powder_plot.c
+++ b/src/powder_plot.c
@@ -102,7 +102,7 @@ int main(int argc, char *argv[])
 
 	hdfile = hdfile_open(filename);
 	hdfile_set_image(hdfile, "/data/data");
-	hdf5_read(hdfile, &image, 1);
+	hdf5_read(hdfile, &image, 1, 2000.0);
 
 	for ( x=0; x<image.width; x++ ) {
 	for ( y=0; y<image.height; y++ ) {
diff --git a/src/reintegrate.c b/src/reintegrate.c
index c9de3338..220cad26 100644
--- a/src/reintegrate.c
+++ b/src/reintegrate.c
@@ -129,7 +129,8 @@ static void process_image(void *pg, int cookie)
 		return;
 	}
 
-	hdf5_read(hdfile, &image, pargs->config_satcorr);
+	/* FIXME: Nominal photon energy */
+	hdf5_read(hdfile, &image, pargs->config_satcorr, 2000.0);
 
 	map_all_peaks(&image);
 
diff --git a/src/render.c b/src/render.c
index 7b694a25..128ede3a 100644
--- a/src/render.c
+++ b/src/render.c
@@ -78,7 +78,8 @@ float *render_get_image_binned(DisplayWindow *dw, int binning, float *max)
 		image->features = NULL;
 		image->data = NULL;
 
-		hdf5_read(dw->hdfile, image, 1);
+		/* We don't care about the photon energy here */
+		hdf5_read(dw->hdfile, image, 1, 0.0);
 		dw->image_dirty = 0;
 		if ( dw->cmfilter ) filter_cm(image);
 		if ( dw->noisefilter ) filter_noise(image, NULL);
diff --git a/src/templates.c b/src/templates.c
index 105fa755..e9bd7b72 100644
--- a/src/templates.c
+++ b/src/templates.c
@@ -109,7 +109,8 @@ UnitCell *rotate_cell(UnitCell *in, double omega, double phi, double rot)
 
 /* Generate templates for the given cell using a representative image */
 IndexingPrivate *generate_templates(UnitCell *cell, const char *filename,
-                                    struct detector *det)
+                                    struct detector *det,
+                                    double nominal_photon_energy)
 {
 	struct _indexingprivate_template *priv;
 	const char *holo;
@@ -129,7 +130,7 @@ IndexingPrivate *generate_templates(UnitCell *cell, const char *filename,
 		ERROR("Couldn't select path\n");
 		return NULL;
 	}
-	hdf5_read(hdfile, &image, 0);
+	hdf5_read(hdfile, &image, 0, nominal_photon_energy);
 	hdfile_close(hdfile);
 	image.det = det;
 
diff --git a/src/templates.h b/src/templates.h
index 50aa7443..33e6c75e 100644
--- a/src/templates.h
+++ b/src/templates.h
@@ -22,7 +22,8 @@
 #include "cell.h"
 
 extern IndexingPrivate *generate_templates(UnitCell *cell, const char *filename,
-                                           struct detector *det);
+                                           struct detector *det,
+                                           double nominal_photon_energy);
 
 
 extern void match_templates(struct image *image, IndexingPrivate *ipriv);